Merge remote-tracking branch 'git/feature/detailed-heart-simulation'
This commit is contained in:
+25
-18
@@ -1,35 +1,42 @@
|
||||
#include <iostream>
|
||||
#include <memory>
|
||||
#include <thread>
|
||||
#include <chrono>
|
||||
#include <vector>
|
||||
#include <string>
|
||||
#include "MedicalLib/MedicalLib.h"
|
||||
#include "MedicalLib/Patient.h"
|
||||
#include "MedicalLib/Organ.h"
|
||||
#include "MedicalLib/Heart.h"
|
||||
#include "MedicalLib/Lungs.h"
|
||||
#include "MedicalLib/Brain.h"
|
||||
#include "MedicalLib/Liver.h"
|
||||
#include "MedicalLib/Kidneys.h"
|
||||
#include "MedicalLib/Stomach.h"
|
||||
|
||||
int main() {
|
||||
// Initialize a new patient
|
||||
Patient patient = initializePatient(1);
|
||||
// Initialize a new patient with a 12-lead heart
|
||||
Patient patient = initializePatient(1, 12);
|
||||
std::cout << "Patient created with ID: " << patient.patientId << std::endl;
|
||||
|
||||
// Simulate some time passing
|
||||
updatePatient(patient, 60.0);
|
||||
std::cout << "\nPatient state updated after 60 seconds." << std::endl;
|
||||
// Simulate time passing and print a live summary
|
||||
const double simulationTime_s = 10.0;
|
||||
const double deltaTime_s = 0.05; // 20 Hz update rate for display
|
||||
const int numSteps = static_cast<int>(simulationTime_s / deltaTime_s);
|
||||
|
||||
// Get a summary for a specific organ
|
||||
std::cout << "\nHeart Summary:\n" << getOrganSummary(patient, "Heart") << std::endl;
|
||||
std::cout << "\n--- Simulating " << simulationTime_s << " seconds of heart activity... ---" << std::endl;
|
||||
|
||||
// Get a summary for all organs
|
||||
std::cout << "\nFull Patient Summary:\n" << getPatientSummary(patient) << std::endl;
|
||||
for (int i = 0; i < numSteps; ++i) {
|
||||
// Clear console on systems that support ANSI escape codes
|
||||
#if defined(__linux__) || defined(__APPLE__)
|
||||
std::cout << "\033[2J\033[1;1H";
|
||||
#endif
|
||||
|
||||
// Get a specific organ and call a method on it
|
||||
if (const Heart* heart = getOrgan<Heart>(patient)) {
|
||||
std::cout << "\nSuccessfully retrieved Heart organ." << std::endl;
|
||||
std::cout << "Direct access to heart rate: " << heart->getHeartRate() << " bpm" << std::endl;
|
||||
updatePatient(patient, deltaTime_s);
|
||||
|
||||
std::cout << "Time: " << i * deltaTime_s << "s / " << simulationTime_s << "s\n" << std::endl;
|
||||
std::cout << getOrganSummary(patient, "Heart") << std::endl;
|
||||
|
||||
std::this_thread::sleep_for(std::chrono::milliseconds(50));
|
||||
}
|
||||
|
||||
std::cout << "\n--- Simulation Complete. Final State: ---\n" << std::endl;
|
||||
std::cout << getPatientSummary(patient) << std::endl;
|
||||
|
||||
return 0;
|
||||
}
|
||||
|
||||
+65
-10
@@ -1,37 +1,92 @@
|
||||
#pragma once
|
||||
|
||||
#include "Organ.h"
|
||||
#include <vector>
|
||||
#include <string>
|
||||
#include <map>
|
||||
#include <deque>
|
||||
|
||||
/**
|
||||
* @brief Represents the Heart organ.
|
||||
* @brief Represents the state of a heart valve.
|
||||
*/
|
||||
enum class ValveStatus { OPEN, CLOSED };
|
||||
|
||||
/**
|
||||
* @brief Represents a single heart valve and its potential pathologies.
|
||||
*/
|
||||
struct Valve {
|
||||
std::string name;
|
||||
ValveStatus status = ValveStatus::CLOSED;
|
||||
double stenosis = 0.0; // Degree of narrowing [0, 1]
|
||||
double regurgitation = 0.0; // Degree of leakage [0, 1]
|
||||
};
|
||||
|
||||
/**
|
||||
* @brief Represents the state of a heart chamber.
|
||||
*/
|
||||
enum class ChamberState { SYSTOLE, DIASTOLE };
|
||||
|
||||
/**
|
||||
* @brief Represents a single chamber of the heart.
|
||||
*/
|
||||
struct Chamber {
|
||||
std::string name;
|
||||
ChamberState state = ChamberState::DIASTOLE;
|
||||
double volume_mL = 0.0;
|
||||
double pressure_mmHg = 0.0;
|
||||
double endDiastolicVolume_mL = 120.0;
|
||||
double endSystolicVolume_mL = 50.0;
|
||||
};
|
||||
|
||||
/**
|
||||
* @brief Represents the Heart organ, with detailed mechanical and electrical simulation.
|
||||
*/
|
||||
class MEDICAL_LIB_API Heart : public Organ {
|
||||
public:
|
||||
/**
|
||||
* @brief Constructor for the Heart class.
|
||||
* @param id The ID of the organ.
|
||||
* @param numLeads The number of EKG leads to simulate (e.g., 3, 5, or 12).
|
||||
*/
|
||||
Heart(int id);
|
||||
Heart(int id, int numLeads = 12);
|
||||
|
||||
/**
|
||||
* @brief Updates the heart's state over time.
|
||||
* @brief Updates the heart's state over time, simulating the cardiac cycle.
|
||||
* @param deltaTime_s The time elapsed in seconds.
|
||||
*/
|
||||
void update(double deltaTime_s) override;
|
||||
|
||||
/**
|
||||
* @brief Gets a string summary of the heart's vitals.
|
||||
* @brief Gets a string summary of the heart's vitals, including EKG and mechanical data.
|
||||
* @return A string containing the heart's vital signs.
|
||||
*/
|
||||
std::string getSummary() const override;
|
||||
|
||||
// Specific getters for Heart properties
|
||||
// --- Electrical Properties ---
|
||||
double getHeartRate() const;
|
||||
double getBloodPressureSystolic() const;
|
||||
double getBloodPressureDiastolic() const;
|
||||
const std::map<std::string, std::deque<double>>& getEkgData() const;
|
||||
|
||||
// --- Mechanical Properties ---
|
||||
double getEjectionFraction() const;
|
||||
double getAorticPressure() const; // Represents systemic blood pressure
|
||||
|
||||
private:
|
||||
double heartRate;
|
||||
double bloodPressureSystolic;
|
||||
double bloodPressureDiastolic;
|
||||
// --- Electrical Simulation ---
|
||||
double simulateEkgWaveform(double timeInCycle);
|
||||
double heartRate; // Underlying target heart rate (bpm)
|
||||
double measuredHeartRate;
|
||||
int numLeads;
|
||||
std::vector<std::string> leadNames;
|
||||
std::map<std::string, std::deque<double>> ekgData;
|
||||
size_t ekgHistorySize;
|
||||
// Cycle timing
|
||||
double totalTime_s;
|
||||
double cardiacCyclePosition_s;
|
||||
double lastRPeakTime_s;
|
||||
bool rPeakDetectedInCycle;
|
||||
|
||||
// --- Mechanical Simulation ---
|
||||
Chamber leftAtrium, rightAtrium, leftVentricle, rightVentricle;
|
||||
Valve mitralValve, tricuspidValve, aorticValve, pulmonaryValve;
|
||||
double ejectionFraction; // Percentage
|
||||
};
|
||||
|
||||
@@ -27,12 +27,20 @@ struct Patient {
|
||||
#endif
|
||||
|
||||
/**
|
||||
* @brief Initializes a new patient with baseline vital signs.
|
||||
* @brief Initializes a new patient with baseline vital signs and a 12-lead heart.
|
||||
* @param patientId The ID for the new patient.
|
||||
* @return A Patient struct with default healthy values.
|
||||
*/
|
||||
MEDICAL_LIB_API Patient initializePatient(int patientId);
|
||||
|
||||
/**
|
||||
* @brief Initializes a new patient with a specific number of heart leads.
|
||||
* @param patientId The ID for the new patient.
|
||||
* @param numHeartLeads The number of EKG leads for the heart.
|
||||
* @return A Patient struct with default healthy values.
|
||||
*/
|
||||
MEDICAL_LIB_API Patient initializePatient(int patientId, int numHeartLeads);
|
||||
|
||||
/**
|
||||
* @brief Updates the patient's vital signs based on the time elapsed.
|
||||
* @param patient The patient to update.
|
||||
|
||||
+170
-19
@@ -2,6 +2,16 @@
|
||||
#include <random>
|
||||
#include <algorithm>
|
||||
#include <sstream>
|
||||
#include <vector>
|
||||
#include <string>
|
||||
#include <cmath>
|
||||
#include <numeric>
|
||||
#include <deque>
|
||||
|
||||
// Helper for Gaussian function to model EKG waves
|
||||
static double gaussian(double x, double mu, double sigma) {
|
||||
return exp(-0.5 * pow((x - mu) / sigma, 2));
|
||||
}
|
||||
|
||||
// Helper function for random fluctuations
|
||||
static double getFluctuation(double stddev) {
|
||||
@@ -11,33 +21,174 @@ static double getFluctuation(double stddev) {
|
||||
return d(gen);
|
||||
}
|
||||
|
||||
Heart::Heart(int id) : Organ(id, "Heart"), heartRate(75.0), bloodPressureSystolic(120.0), bloodPressureDiastolic(80.0) {}
|
||||
Heart::Heart(int id, int numLeads)
|
||||
: Organ(id, "Heart"),
|
||||
heartRate(75.0),
|
||||
measuredHeartRate(75.0),
|
||||
numLeads(numLeads),
|
||||
ekgHistorySize(200),
|
||||
totalTime_s(0.0),
|
||||
cardiacCyclePosition_s(0.0),
|
||||
lastRPeakTime_s(-1.0),
|
||||
rPeakDetectedInCycle(false),
|
||||
ejectionFraction(0.55) {
|
||||
|
||||
// Initialize Chambers
|
||||
leftAtrium.name = "Left Atrium";
|
||||
rightAtrium.name = "Right Atrium";
|
||||
leftVentricle.name = "Left Ventricle";
|
||||
rightVentricle.name = "Right Ventricle";
|
||||
|
||||
// Initialize Valves
|
||||
mitralValve.name = "Mitral Valve";
|
||||
tricuspidValve.name = "Tricuspid Valve";
|
||||
aorticValve.name = "Aortic Valve";
|
||||
pulmonaryValve.name = "Pulmonary Valve";
|
||||
|
||||
const std::vector<std::string> allLeadNames = {
|
||||
"I", "II", "III", "aVR", "aVL", "aVF",
|
||||
"V1", "V2", "V3", "V4", "V5", "V6"
|
||||
};
|
||||
|
||||
int leadsToCreate = std::min((int)allLeadNames.size(), numLeads);
|
||||
for (int i = 0; i < leadsToCreate; ++i) {
|
||||
leadNames.push_back(allLeadNames[i]);
|
||||
ekgData[allLeadNames[i]] = std::deque<double>();
|
||||
}
|
||||
}
|
||||
|
||||
void Heart::update(double deltaTime_s) {
|
||||
const double baseline_hr = 75.0;
|
||||
const double baseline_bp_systolic = 120.0;
|
||||
const double baseline_bp_diastolic = 80.0;
|
||||
const double theta = 0.1; // Mean reversion speed
|
||||
const double hr_stddev = 0.1;
|
||||
const double bp_stddev = 0.1;
|
||||
// --- Electrical Simulation Update ---
|
||||
totalTime_s += deltaTime_s;
|
||||
heartRate += getFluctuation(0.01); // Slow variation in underlying rate
|
||||
heartRate = std::max(60.0, std::min(heartRate, 100.0));
|
||||
double cycleDuration_s = 60.0 / heartRate;
|
||||
|
||||
heartRate += theta * (baseline_hr - heartRate) * deltaTime_s + getFluctuation(hr_stddev * deltaTime_s);
|
||||
bloodPressureSystolic += theta * (baseline_bp_systolic - bloodPressureSystolic) * deltaTime_s + getFluctuation(bp_stddev * deltaTime_s);
|
||||
bloodPressureDiastolic += theta * (baseline_bp_diastolic - bloodPressureDiastolic) * deltaTime_s + getFluctuation(bp_stddev * deltaTime_s);
|
||||
double oldCyclePosition = cardiacCyclePosition_s;
|
||||
cardiacCyclePosition_s += deltaTime_s;
|
||||
|
||||
heartRate = std::clamp(heartRate, 60.0, 100.0);
|
||||
bloodPressureSystolic = std::clamp(bloodPressureSystolic, 90.0, 120.0);
|
||||
bloodPressureDiastolic = std::clamp(bloodPressureDiastolic, 60.0, 80.0);
|
||||
// R-peak detection for measured heart rate
|
||||
double rPeakCycleTime = 0.22 * cycleDuration_s;
|
||||
if (oldCyclePosition < rPeakCycleTime && cardiacCyclePosition_s >= rPeakCycleTime && !rPeakDetectedInCycle) {
|
||||
if (lastRPeakTime_s > 0) {
|
||||
measuredHeartRate = 60.0 / (totalTime_s - lastRPeakTime_s);
|
||||
}
|
||||
lastRPeakTime_s = totalTime_s;
|
||||
rPeakDetectedInCycle = true;
|
||||
}
|
||||
|
||||
if (cardiacCyclePosition_s > cycleDuration_s) {
|
||||
cardiacCyclePosition_s -= cycleDuration_s;
|
||||
rPeakDetectedInCycle = false;
|
||||
}
|
||||
|
||||
double timeInCycle = cardiacCyclePosition_s / cycleDuration_s;
|
||||
double baseVoltage = simulateEkgWaveform(timeInCycle);
|
||||
for(const auto& leadName : leadNames) {
|
||||
double leadModifier = 1.0 - (std::distance(leadNames.begin(), std::find(leadNames.begin(), leadNames.end(), leadName)) * 0.1);
|
||||
ekgData[leadName].push_front(baseVoltage * leadModifier);
|
||||
if (ekgData[leadName].size() > ekgHistorySize) ekgData[leadName].pop_back();
|
||||
}
|
||||
|
||||
// --- Mechanical Simulation Update ---
|
||||
// Define pressures outside the heart
|
||||
double venousPressure = 5.0; // CVP
|
||||
double pulmonaryArteryPressure = 20.0;
|
||||
double aorticPressure = getAorticPressure();
|
||||
|
||||
// Determine chamber states based on EKG-timed cycle
|
||||
if (timeInCycle >= 0.0 && timeInCycle < 0.15) { // Atrial Systole
|
||||
leftAtrium.state = ChamberState::SYSTOLE;
|
||||
rightAtrium.state = ChamberState::SYSTOLE;
|
||||
// End of ventricular filling, capture EDV
|
||||
if(oldCyclePosition < 0.15 && timeInCycle >= 0.15) {
|
||||
leftVentricle.endDiastolicVolume_mL = leftVentricle.volume_mL;
|
||||
}
|
||||
} else {
|
||||
leftAtrium.state = ChamberState::DIASTOLE;
|
||||
rightAtrium.state = ChamberState::DIASTOLE;
|
||||
}
|
||||
|
||||
if (timeInCycle >= 0.20 && timeInCycle < 0.5) { // Ventricular Systole
|
||||
leftVentricle.state = ChamberState::SYSTOLE;
|
||||
rightVentricle.state = ChamberState::SYSTOLE;
|
||||
// End of ventricular ejection, capture ESV and calculate EF
|
||||
if(oldCyclePosition < 0.5 && timeInCycle >= 0.5) {
|
||||
leftVentricle.endSystolicVolume_mL = leftVentricle.volume_mL;
|
||||
if (leftVentricle.endDiastolicVolume_mL > 0) {
|
||||
ejectionFraction = (leftVentricle.endDiastolicVolume_mL - leftVentricle.endSystolicVolume_mL) / leftVentricle.endDiastolicVolume_mL;
|
||||
}
|
||||
}
|
||||
} else {
|
||||
leftVentricle.state = ChamberState::DIASTOLE;
|
||||
rightVentricle.state = ChamberState::DIASTOLE;
|
||||
}
|
||||
|
||||
// Update chamber pressures based on state (very simplified model)
|
||||
leftAtrium.pressure_mmHg = (leftAtrium.state == ChamberState::SYSTOLE) ? 10.0 : 5.0;
|
||||
rightAtrium.pressure_mmHg = (rightAtrium.state == ChamberState::SYSTOLE) ? 7.0 : 2.0;
|
||||
leftVentricle.pressure_mmHg = (leftVentricle.state == ChamberState::SYSTOLE) ? 125.0 * sin((timeInCycle - 0.2) / 0.3 * M_PI) : 5.0;
|
||||
rightVentricle.pressure_mmHg = (rightVentricle.state == ChamberState::SYSTOLE) ? 25.0 * sin((timeInCycle - 0.2) / 0.3 * M_PI) : 2.0;
|
||||
|
||||
// Update Valve Status
|
||||
tricuspidValve.status = (rightAtrium.pressure_mmHg > rightVentricle.pressure_mmHg) ? ValveStatus::OPEN : ValveStatus::CLOSED;
|
||||
mitralValve.status = (leftAtrium.pressure_mmHg > leftVentricle.pressure_mmHg) ? ValveStatus::OPEN : ValveStatus::CLOSED;
|
||||
pulmonaryValve.status = (rightVentricle.pressure_mmHg > pulmonaryArteryPressure) ? ValveStatus::OPEN : ValveStatus::CLOSED;
|
||||
aorticValve.status = (leftVentricle.pressure_mmHg > aorticPressure) ? ValveStatus::OPEN : ValveStatus::CLOSED;
|
||||
|
||||
// Update Chamber Volumes (simplified blood flow)
|
||||
double flowRate = 500.0 * deltaTime_s; // mL/s
|
||||
if (mitralValve.status == ValveStatus::OPEN) leftVentricle.volume_mL += flowRate;
|
||||
if (tricuspidValve.status == ValveStatus::OPEN) rightVentricle.volume_mL += flowRate;
|
||||
if (aorticValve.status == ValveStatus::OPEN) leftVentricle.volume_mL -= flowRate * 1.5;
|
||||
if (pulmonaryValve.status == ValveStatus::OPEN) rightVentricle.volume_mL -= flowRate * 1.5;
|
||||
|
||||
// Clamp volumes to realistic values
|
||||
leftVentricle.volume_mL = std::max(40.0, std::min(leftVentricle.volume_mL, 130.0));
|
||||
rightVentricle.volume_mL = std::max(40.0, std::min(rightVentricle.volume_mL, 130.0));
|
||||
}
|
||||
|
||||
double Heart::simulateEkgWaveform(double timeInCycle) {
|
||||
double p_time = 0.1, q_time = 0.2, r_time = 0.22, s_time = 0.24, t_time = 0.4;
|
||||
double p_amp = 0.15, q_amp = -0.1, r_amp = 1.0, s_amp = -0.25, t_amp = 0.3;
|
||||
double p_sigma = 0.04, qrs_sigma = 0.02, t_sigma = 0.06;
|
||||
double voltage = 0.0;
|
||||
voltage += p_amp * gaussian(timeInCycle, p_time, p_sigma);
|
||||
voltage += q_amp * gaussian(timeInCycle, q_time, qrs_sigma);
|
||||
voltage += r_amp * gaussian(timeInCycle, r_time, qrs_sigma);
|
||||
voltage += s_amp * gaussian(timeInCycle, s_time, qrs_sigma);
|
||||
voltage += t_amp * gaussian(timeInCycle, t_time, t_sigma);
|
||||
return voltage;
|
||||
}
|
||||
|
||||
std::string Heart::getSummary() const {
|
||||
std::stringstream ss;
|
||||
ss << "Type: " << organType << " (ID: " << organId << ")\n"
|
||||
<< " Heart Rate: " << heartRate << " bpm\n"
|
||||
<< " Blood Pressure: " << bloodPressureSystolic << "/" << bloodPressureDiastolic << " mmHg";
|
||||
ss.precision(2);
|
||||
ss << std::fixed;
|
||||
ss << "--- Heart Summary ---\n"
|
||||
<< "Heart Rate (Measured): " << getHeartRate() << " bpm\n"
|
||||
<< "Ejection Fraction: " << getEjectionFraction() * 100.0 << "%\n"
|
||||
<< "Aortic Pressure: " << getAorticPressure() << " mmHg\n\n"
|
||||
<< "--- Chambers ---\n"
|
||||
<< " LV Volume: " << leftVentricle.volume_mL << " mL\n"
|
||||
<< " LV Pressure: " << leftVentricle.pressure_mmHg << " mmHg\n"
|
||||
<< " RV Volume: " << rightVentricle.volume_mL << " mL\n"
|
||||
<< " RV Pressure: " << rightVentricle.pressure_mmHg << " mmHg\n\n"
|
||||
<< "--- Valves ---\n"
|
||||
<< " Aortic Valve: " << (aorticValve.status == ValveStatus::OPEN ? "OPEN" : "CLOSED") << "\n"
|
||||
<< " Mitral Valve: " << (mitralValve.status == ValveStatus::OPEN ? "OPEN" : "CLOSED") << "\n";
|
||||
return ss.str();
|
||||
}
|
||||
|
||||
double Heart::getHeartRate() const { return heartRate; }
|
||||
double Heart::getBloodPressureSystolic() const { return bloodPressureSystolic; }
|
||||
double Heart::getBloodPressureDiastolic() const { return bloodPressureDiastolic; }
|
||||
double Heart::getHeartRate() const { return measuredHeartRate; }
|
||||
const std::map<std::string, std::deque<double>>& Heart::getEkgData() const { return ekgData; }
|
||||
double Heart::getEjectionFraction() const { return ejectionFraction; }
|
||||
double Heart::getAorticPressure() const {
|
||||
// Aortic pressure decays over time, but gets "pumped up" by ventricular ejection
|
||||
// This is a placeholder for a more complex arterial model.
|
||||
if(aorticValve.status == ValveStatus::OPEN) {
|
||||
return leftVentricle.pressure_mmHg;
|
||||
}
|
||||
// Simplified diastolic pressure decay
|
||||
return 80.0 + 40.0 * exp(-cardiacCyclePosition_s);
|
||||
}
|
||||
|
||||
+14
-3
@@ -15,16 +15,17 @@
|
||||
#include <memory>
|
||||
|
||||
/**
|
||||
* @brief Initializes a new patient with baseline vital signs and a standard set of organs.
|
||||
* @brief Initializes a new patient with a specific number of heart leads.
|
||||
* @param patientId The ID for the new patient.
|
||||
* @param numHeartLeads The number of EKG leads for the heart.
|
||||
* @return A Patient struct with default healthy values.
|
||||
*/
|
||||
Patient initializePatient(int patientId) {
|
||||
Patient initializePatient(int patientId, int numHeartLeads) {
|
||||
Patient patient;
|
||||
patient.patientId = patientId;
|
||||
|
||||
// Initialize default organs
|
||||
patient.organs.push_back(std::make_unique<Heart>(1));
|
||||
patient.organs.push_back(std::make_unique<Heart>(1, numHeartLeads));
|
||||
patient.organs.push_back(std::make_unique<Lungs>(2));
|
||||
patient.organs.push_back(std::make_unique<Brain>(3));
|
||||
patient.organs.push_back(std::make_unique<Liver>(4));
|
||||
@@ -41,6 +42,16 @@ Patient initializePatient(int patientId) {
|
||||
return patient;
|
||||
}
|
||||
|
||||
/**
|
||||
* @brief Initializes a new patient with baseline vital signs and a 12-lead heart.
|
||||
* @param patientId The ID for the new patient.
|
||||
* @return A Patient struct with default healthy values.
|
||||
*/
|
||||
Patient initializePatient(int patientId) {
|
||||
return initializePatient(patientId, 12);
|
||||
}
|
||||
|
||||
|
||||
/**
|
||||
* @brief Updates the patient's state by updating the state of all their organs.
|
||||
* @param patient The patient to update.
|
||||
|
||||
Reference in New Issue
Block a user