/* Distributed under the Apache License, Version 2.0.
   See accompanying NOTICE file for details.*/

#include "engine/human_adult/whole_body/controller/Controller.h"

#include "cdm/patient/SEPatient.h"
#include "cdm/properties/SEScalar0To1.h"
#include "cdm/properties/SEScalarArea.h"
#include "cdm/properties/SEScalarFrequency.h"
#include "cdm/properties/SEScalarHeatCapacitance.h"
#include "cdm/properties/SEScalarHeatCapacitancePerMass.h"
#include "cdm/properties/SEScalarLength.h"
#include "cdm/properties/SEScalarMass.h"
#include "cdm/properties/SEScalarMassPerMass.h"
#include "cdm/properties/SEScalarMassPerVolume.h"
#include "cdm/properties/SEScalarPower.h"
#include "cdm/properties/SEScalarPressure.h"
#include "cdm/properties/SEScalarPressurePerVolume.h"
#include "cdm/properties/SEScalarTime.h"
#include "cdm/properties/SEScalarVolume.h"

namespace pulse { namespace human_adult_whole_body
{
  bool Controller::SetupPatient(SEPatient const& patient)
  {
    bool err = false;
    std::stringstream ss;
    m_InitialPatient->Copy(patient);
    m_CurrentPatient->Clear();

    //Sex is the only thing we absolutely need to be defined, the CDM assumes male if not provided

    //AGE ---------------------------------------------------------------
    double age_yr;
    double ageMin_yr = 18.0;
    double ageMax_yr = 65.0;
    double ageStandard_yr = 44.0;
    if (!m_InitialPatient->HasAge())
    {
      m_InitialPatient->GetAge().SetValue(ageStandard_yr, TimeUnit::yr);
      ss << "No patient age set. Using the standard value of " << ageStandard_yr << " years.";
      Info(ss);
    }
    age_yr = m_InitialPatient->GetAge().GetValue(TimeUnit::yr);
    if (age_yr < ageMin_yr)
    {
      ss << "Patient age of " << age_yr << " years is too young. We do not model pediatrics. Minimum age allowed is " << ageMin_yr << " years.";
      Error(ss);
      err = true;
    }
    else if (age_yr > ageMax_yr)
    {
      ss << "Patient age of " << age_yr << " years is too old. We do not model geriatrics. Maximum age allowed is " << ageMax_yr << " years.";
      Error(ss);
      err = true;
    }

    //HEIGHT ---------------------------------------------------------------
    //From CDC values for 20 year olds
    //Mins are 3rd percentile, Maxs are 97th percentile, and standard is 50th percentile
    /// \cite Centers2016clinical
    double heightMinMale_cm = 163.0;
    double heightMaxMale_cm = 190.0;
    double heightStandardMale_cm = 177.0;
    double heightMinFemale_cm = 151.0;
    double heightMaxFemale_cm = 175.5;
    double heightStandardFemale_cm = 163.0;
    //Male
    double heightMin_cm = heightMinMale_cm;
    double heightMax_cm = heightMaxMale_cm;
    double heightStandard_cm = heightStandardMale_cm;
    if (m_InitialPatient->GetSex() == ePatient_Sex::Female)
    {
      //Female
      heightMin_cm = heightMinFemale_cm;
      heightMax_cm = heightMaxFemale_cm;
      heightStandard_cm = heightStandardFemale_cm;
    }
    if (!m_InitialPatient->HasHeight())
    {
      m_InitialPatient->GetHeight().SetValue(heightStandard_cm, LengthUnit::cm);
      ss << "No patient height set. Using the standard value of " << heightStandard_cm << " cm.";
      Info(ss);
    }
    double height_cm = m_InitialPatient->GetHeight().GetValue(LengthUnit::cm);
    double height_ft = Convert(height_cm, LengthUnit::cm, LengthUnit::ft);
    //Check for outrageous values
    if (height_ft < 4.5 || height_ft > 7.0)
    {
      Error("Patient height setting is outrageous. It must be between 4.5 and 7.0 ft");
      err = true;
    }
    if (height_cm < heightMin_cm)
    {
      ss << "Patient height of " << height_cm << " cm is outside of typical ranges - below 3rd percentile (" << heightMax_cm << " cm). No guarantees of model validity.";
      Warning(ss);
    }
    else if (height_cm > heightMax_cm)
    {
      ss << "Patient height of " << height_cm << " cm is outside of typical ranges - above 97th percentile(" << heightMin_cm << " cm). No guarantees of model validity.";
      Warning(ss);
    }

    //IDEAL BODY WEIGHT ----------------------------------------------------
    /// \cite green2017green
    //page 295
    //Devine Formula
    if (m_InitialPatient->HasIdealBodyWeight())
      Warning("Ignorning provided patient ideal body weight. It is determined by weight and body fat fraction.");
    double height_in = m_InitialPatient->GetHeight().GetValue(LengthUnit::in);
    double idealWeight_kg = 0.0;
    if (m_InitialPatient->GetSex() == ePatient_Sex::Female)
    {
      //Female
      idealWeight_kg = 45.5 + 2.3 * (height_in - 60.0);
    }
    else
    {
      //Male
      idealWeight_kg = 50.0 + 2.3 * (height_in - 60.0);
    }
    m_InitialPatient->GetIdealBodyWeight().SetValue(idealWeight_kg, MassUnit::kg);
    ss << "Patient ideal body weight computed and set to " << idealWeight_kg << " kg.";
    Info(ss);

    //WEIGHT ---------------------------------------------------------------
    /// \cite World2006bmi
    double weight_kg;
    double BMI_kg_per_m2;
    double BMIStandard_kg_per_m2 = 21.75;
    double BMIObese_kg_per_m2 = 30.0;
    double BMIOverweight_kg_per_m2 = 25.0;
    double BMIUnderweight_kg_per_m2 = 18.5;
    double BMISeverelyUnderweight_kg_per_m2 = 16.0;
    if (!m_InitialPatient->HasWeight())
    {
      weight_kg = BMIStandard_kg_per_m2 * pow(m_InitialPatient->GetHeight().GetValue(LengthUnit::m), 2);
      m_InitialPatient->GetWeight().SetValue(weight_kg, MassUnit::kg);
      ss << "No patient weight set. Using the standard BMI value of 21.75 kg/m^2, resulting in a weight of " << weight_kg << " kg.";
      Info(ss);
    }
    weight_kg = m_InitialPatient->GetWeight(MassUnit::kg);
    BMI_kg_per_m2 = weight_kg / pow(m_InitialPatient->GetHeight().GetValue(LengthUnit::m), 2);
    if (BMI_kg_per_m2 > BMIObese_kg_per_m2)
    {
      ss << "Patient Body Mass Index (BMI) of " << BMI_kg_per_m2 << "  kg/m^2 is too high. Obese patients must be modeled by adding/using a condition. Maximum BMI allowed is " << BMIObese_kg_per_m2 << " kg/m^2.";
      Error(ss);
      err = true;
    }
    if (BMI_kg_per_m2 > BMIOverweight_kg_per_m2)
    {
      ss << "Patient Body Mass Index (BMI) of " << BMI_kg_per_m2 << " kg/m^2 is overweight. No guarantees of model validity.";
      Warning(ss);
    }
    if (BMI_kg_per_m2 < BMIUnderweight_kg_per_m2)
    {
      ss << "Patient Body Mass Index (BMI) of " << BMI_kg_per_m2 << " kg/m^2 is underweight. No guarantees of model validity.";
      Warning(ss);
    }
    if (BMI_kg_per_m2 < BMISeverelyUnderweight_kg_per_m2)
    {
      ss << "Patient Body Mass Index (BMI) of " << BMI_kg_per_m2 << " kg/m^2 is too low. Severly underweight patients must be modeled by adding/using a condition. Maximum BMI allowed is " << BMISeverelyUnderweight_kg_per_m2 << " kg/m^2.";
      Error(ss);
      err = true;
    }

    //BODY FAT FRACTION ---------------------------------------------------------------
    //From American Council on Exercise
    /// \cite muth2009what
    double fatFraction;
    double fatFractionStandardMale = 0.21;
    double fatFractionStandardFemale = 0.28;
    double fatFractionMaxMale = 0.25; //Obese
    double fatFractionMaxFemale = 0.32; //Obese
    double fatFractionMinMale = 0.02; //Essential fat
    double fatFractionMinFemale = 0.10; //Essential fat  
    //Male
    double fatFractionMin = fatFractionMinMale;
    double fatFractionMax = fatFractionMaxMale;
    double fatFractionStandard = fatFractionStandardMale;
    if (m_InitialPatient->GetSex() == ePatient_Sex::Female)
    {
      //Female
      fatFractionMin = fatFractionMinFemale;
      fatFractionMax = fatFractionMaxFemale;
      fatFractionStandard = fatFractionStandardFemale;
    }

    if (!m_InitialPatient->HasBodyFatFraction())
    {
      fatFraction = fatFractionStandard;
      m_InitialPatient->GetBodyFatFraction().SetValue(fatFraction);
      ss << "No patient body fat fraction set. Using the standard value of " << fatFraction << ".";
      Info(ss);
    }
    fatFraction = m_InitialPatient->GetBodyFatFraction().GetValue();
    if (fatFraction > fatFractionMax)
    {
      ss << "Patient body fat fraction of " << fatFraction << " is too high. Obese patients must be modeled by adding/using a condition. Maximum body fat fraction allowed is " << fatFractionMax << ".";
      Error(ss);
      err = true;
    }
    else if (fatFraction < fatFractionMin)
    {
      ss << "Patient body fat fraction  " << fatFraction << " is too low. Patients must have essential fat. Minimum body fat fraction allowed is " << fatFractionMin << ".";
      Error(ss);
      err = true;
    }

    //Lean Body Mass ---------------------------------------------------------------
    if (m_InitialPatient->HasLeanBodyMass())
    {
      ss << "Patient lean body mass cannot be set. It is determined by weight and body fat fraction.";
      Error(ss);
      err = true;
    }
    double leanBodyMass_kg = weight_kg * (1.0 - fatFraction);
    m_InitialPatient->GetLeanBodyMass().SetValue(leanBodyMass_kg, MassUnit::kg);
    ss << "Patient lean body mass computed and set to " << leanBodyMass_kg << " kg.";
    Info(ss);

    //Body Density ---------------------------------------------------------------
    if (m_InitialPatient->HasBodyDensity())
    {
      ss << "Patient body density cannot be set. It is determined using body fat fraction.";
      Error(ss);
      err = true;
    }
    //Using the average of Siri and Brozek formulas
    /// \cite siri1961body 
    /// \cite brovzek1963densitometric
    double SiriBodyDensity_g_Per_cm3 = 4.95 / (fatFraction + 4.50);
    double BrozekBodyDensity_g_Per_cm3 = 4.57 / (fatFraction + 4.142);
    double bodyDensity_g_Per_cm3 = (SiriBodyDensity_g_Per_cm3 + BrozekBodyDensity_g_Per_cm3) / 2.0;
    m_InitialPatient->GetBodyDensity().SetValue(bodyDensity_g_Per_cm3, MassPerVolumeUnit::g_Per_cm3);
    ss << "Patient body density computed and set to " << bodyDensity_g_Per_cm3 << " g/cm^3.";
    Info(ss);

    //Heart Rate ---------------------------------------------------------------
    double heartRate_bpm;
    double heartStandard_bpm = 72.0;
    double heartRateMax_bpm = 100.0;
    double heartRateTachycardia_bpm = 110;
    double heartRateMin_bpm = 60.0;
    double heartRateBradycardia_bpm = 50;
    if (!m_InitialPatient->HasHeartRateBaseline())
    {
      heartRate_bpm = heartStandard_bpm;
      m_InitialPatient->GetHeartRateBaseline().SetValue(heartRate_bpm, FrequencyUnit::Per_min);
      ss << "No patient heart rate baseline set. Using the standard value of " << heartRate_bpm << " bpm.";
      Info(ss);
    }
    heartRate_bpm = m_InitialPatient->GetHeartRateBaseline(FrequencyUnit::Per_min);
    if (heartRate_bpm > heartRateMax_bpm)
    {
      if (heartRate_bpm <= heartRateTachycardia_bpm)
      {
        ss << "Patient heart rate baseline of " << heartRate_bpm << " bpm is tachycardic. Tachycardia heart rate  is [" << heartRateMax_bpm << "," << heartRateTachycardia_bpm << "] bpm.";
        Info(ss);
      }
      else
      {
        ss << "Patient heart rate baseline of " << heartRate_bpm << " bpm is too high. Maximum heart rate baseline allowed is " << heartRateTachycardia_bpm << " bpm.";
        Error(ss);
        err = true;
      }
    }
    else if (heartRate_bpm < heartRateMin_bpm)
    {
      if (heartRate_bpm <= heartRateTachycardia_bpm)
      {
        ss << "Patient heart rate baseline of " << heartRate_bpm << " bpm is bradycardic. Bradycardia heart rate  is [" << heartRateBradycardia_bpm << "," << heartRateMin_bpm << "] bpm.";
        Info(ss);
      }
      else
      {
        ss << "Patient heart rate baseline of " << heartRate_bpm << " is too low. Minimum heart rate baseline allowed is " << heartRateBradycardia_bpm << " bpm.";
        Error(ss);
        err = true;
      }
    }

    //Tanaka H, Monahan KD, Seals DR (January 2001). "Age-predicted maximal heart rate revisited". J. Am. Coll. Cardiol. 37(1): 153�6. doi:10.1016/S0735-1097(00)01054-8.PMID 11153730.
    double computedHeartRateMaximum_bpm = 208.0 - (0.7 * m_InitialPatient->GetAge(TimeUnit::yr));
    if (!m_InitialPatient->HasHeartRateMaximum())
    {
      m_InitialPatient->GetHeartRateMaximum().SetValue(computedHeartRateMaximum_bpm, FrequencyUnit::Per_min);
      ss << "No patient heart rate maximum set. Using a computed value of " << computedHeartRateMaximum_bpm << " bpm.";
      Info(ss);
    }
    else
    {
      if (m_InitialPatient->GetHeartRateMaximum(FrequencyUnit::Per_min) < heartRate_bpm)
      {
        ss << "Patient heart rate maximum must be greater than the baseline heart rate.";
        Error(ss);
        err = true;
      }
      ss << "Specified patient heart rate maximum of " << m_InitialPatient->GetHeartRateMaximum(FrequencyUnit::Per_min) << " bpm differs from computed value of " << computedHeartRateMaximum_bpm << " bpm. No guarantees of model validity.";
      Warning(ss);
    }
    if (!m_InitialPatient->HasHeartRateMinimum())
    {
      m_InitialPatient->GetHeartRateMinimum().SetValue(0.001, FrequencyUnit::Per_min);
      ss << "No patient heart rate minimum set. Using a default value of " << 0.001 << " bpm.";
      Info(ss);
    }
    if (m_InitialPatient->GetHeartRateMinimum(FrequencyUnit::Per_min) > heartRate_bpm)
    {
      ss << "Patient heart rate minimum must be less than the baseline heart rate.";
      Error(ss);
      err = true;
    }

    //Arterial Pressures ---------------------------------------------------------------
    double systolic_mmHg;
    double diastolic_mmHg;
    double systolicStandard_mmHg = 114.0;
    double diastolicStandard_mmHg = 73.5;
    double systolicMax_mmHg = 120.0; //Hypertension
    double diastolicMax_mmHg = 80.0; //Hypertension
    double systolicMin_mmHg = 90.0; //Hypotension
    double diastolicMin_mmHg = 60.0; //Hypotension
    double narrowestPulseFactor = 0.75; //From Wikipedia: Pulse Pressure
    if (!m_InitialPatient->HasSystolicArterialPressureBaseline())
    {
      systolic_mmHg = systolicStandard_mmHg;
      m_InitialPatient->GetSystolicArterialPressureBaseline().SetValue(systolic_mmHg, PressureUnit::mmHg);
      ss << "No patient systolic pressure baseline set. Using the standard value of " << systolic_mmHg << " mmHg.";
      Info(ss);
    }
    systolic_mmHg = m_InitialPatient->GetSystolicArterialPressureBaseline(PressureUnit::mmHg);
    if (systolic_mmHg < systolicMin_mmHg)
    {
      ss << "Patient systolic pressure baseline of " << systolic_mmHg << " mmHg is too low. Hypotension must be modeled by adding/using a condition. Minimum systolic pressure baseline allowed is " << systolicMin_mmHg << " mmHg.";
      Error(ss);
      err = true;
    }
    else if (systolic_mmHg > systolicMax_mmHg)
    {
      ss << "Patient systolic pressure baseline of " << systolic_mmHg << " mmHg is too high. Hypertension must be modeled by adding/using a condition. Maximum systolic pressure baseline allowed is " << systolicMax_mmHg << " mmHg.";
      Error(ss);
      err = true;
    }

    if (!m_InitialPatient->HasDiastolicArterialPressureBaseline())
    {
      diastolic_mmHg = diastolicStandard_mmHg;
      m_InitialPatient->GetDiastolicArterialPressureBaseline().SetValue(diastolic_mmHg, PressureUnit::mmHg);
      ss << "No patient diastolic pressure baseline set. Using the standard value of " << diastolic_mmHg << " mmHg.";
      Info(ss);
    }
    diastolic_mmHg = m_InitialPatient->GetDiastolicArterialPressureBaseline(PressureUnit::mmHg);
    if (diastolic_mmHg < diastolicMin_mmHg)
    {
      ss << "Patient diastolic pressure baseline of " << diastolic_mmHg << " mmHg is too low. Hypotension must be modeled by adding/using a condition. Minimum diastolic pressure baseline allowed is " << diastolicMin_mmHg << " mmHg.";
      Error(ss);
      err = true;
    }
    else if (diastolic_mmHg > diastolicMax_mmHg)
    {
      ss << "Patient diastolic pressure baseline of " << diastolic_mmHg << " mmHg is too high. Hypertension must be modeled by adding/using a condition. Maximum diastolic pressure baseline allowed is " << diastolicMax_mmHg << " mmHg.";
      Error(ss);
      err = true;
    }

    if (diastolic_mmHg > 0.75 * systolic_mmHg)
    {
      ss << "Patient baseline pulse pressure (systolic vs. diastolic pressure fraction) of " << diastolic_mmHg / systolic_mmHg << " is abnormally narrow. Minimum fraction allowed is " << narrowestPulseFactor << " .";
      Error(ss);
      err = true;
    }

    if (m_InitialPatient->HasMeanArterialPressureBaseline())
    {
      ss << "Patient mean arterial pressure baseline cannot be set. It is determined through homeostatic simulation.";
      Error(ss);
      err = true;
    }
    double MAP_mmHg = 1.0 / 3.0 * systolic_mmHg + 2.0 / 3.0 * diastolic_mmHg;
    m_InitialPatient->GetMeanArterialPressureBaseline().SetValue(MAP_mmHg, PressureUnit::mmHg);

    //Blood Volume ---------------------------------------------------------------
    /// \cite Morgan2006Clinical
    double bloodVolume_mL;
    double computedBloodVolume_mL = 65.6 * pow(weight_kg, 1.02);
    double bloodVolumeMin_mL = computedBloodVolume_mL * 0.85; //Stage 1 Hypovolemia
    double bloodVolumeMax_mL = computedBloodVolume_mL * 1.15; //Just go the same distance on the other side
    if (!m_InitialPatient->HasBloodVolumeBaseline())
    {
      bloodVolume_mL = computedBloodVolume_mL;
      m_InitialPatient->GetBloodVolumeBaseline().SetValue(bloodVolume_mL, VolumeUnit::mL);
      ss << "No patient blood volume baseline set. Using a computed value of " << computedBloodVolume_mL << " mL.";
      Info(ss);
    }
    bloodVolume_mL = m_InitialPatient->GetBloodVolumeBaseline(VolumeUnit::mL);
    if (bloodVolume_mL != computedBloodVolume_mL)
    {
      ss << "Specified patient blood volume baseline of " << bloodVolume_mL << " mL differs from computed value of " << computedBloodVolume_mL << " mL. No guarantees of model validity and there is a good chance the patient will not reach a starting homeostatic point.";
      Warning(ss);
    }
    if (bloodVolume_mL < bloodVolumeMin_mL)
    {
      ss << "Patient blood volume baseline of " << bloodVolume_mL << " mL is too low. Hypovolemia must be modeled by adding/using a condition. Minimum blood volume baseline allowed is " << bloodVolumeMin_mL << " mL.";
      Error(ss);
      err = true;
    }
    else if (bloodVolume_mL > bloodVolumeMax_mL)
    {
      ss << "Patient blood volume baseline of " << bloodVolume_mL << " mL is too high. Excessive volume must be modeled by adding/using a condition. Maximum blood volume baseline allowed is " << bloodVolumeMax_mL << " mL.";
      Error(ss);
      err = true;
    }

    //Respiration Rate ---------------------------------------------------------------
    //Note: This is overwritten after stabilization
    /// \cite green2017green
    double respirationRate_bpm;
    double respirationRateStandard_bpm = 12.0;
    double respirationRateMax_bpm = 20.0;
    double respirationRateMin_bpm = 8.0;
    if (!m_InitialPatient->HasRespirationRateBaseline())
    {
      respirationRate_bpm = respirationRateStandard_bpm;
      m_InitialPatient->GetRespirationRateBaseline().SetValue(respirationRate_bpm, FrequencyUnit::Per_min);
      ss << "No patient respiration rate baseline set. Using the standard value of " << respirationRate_bpm << " bpm.";
      Info(ss);
    }
    respirationRate_bpm = m_InitialPatient->GetRespirationRateBaseline(FrequencyUnit::Per_min);
    if (respirationRate_bpm > respirationRateMax_bpm)
    {
      ss << "Patient respiration rate baseline of " << respirationRate_bpm << " bpm is too high. Non-healthy values must be modeled by adding/using a condition. Maximum respiration rate baseline allowed is " << respirationRateMax_bpm << " bpm.";
      Error(ss);
      err = true;
    }
    else if (respirationRate_bpm < respirationRateMin_bpm)
    {
      ss << "Patient respiration rate baseline of " << respirationRate_bpm << " bpm is too low. Non-healthy values must be modeled by adding/using a condition. Minimum respiration rate baseline allowed is " << respirationRateMin_bpm << " bpm.";
      Error(ss);
      err = true;
    }

    //Right Lung Ratio ---------------------------------------------------------------
    double rightLungRatio;
    double rightLungRatioStandard = 0.525;
    double rightLungRatioMax = 0.60;
    double rightLungRatioMin = 0.50;
    if (!m_InitialPatient->HasRightLungRatio())
    {
      rightLungRatio = rightLungRatioStandard;
      m_InitialPatient->GetRightLungRatio().SetValue(rightLungRatio);
      ss << "No patient right lung ratio set. Using the standard value of " << rightLungRatio << ".";
      Info(ss);
    }
    rightLungRatio = m_InitialPatient->GetRightLungRatio().GetValue();
    if (rightLungRatio > rightLungRatioMax)
    {
      ss << "Patient right lung ratio of " << rightLungRatio << " is too high. Non-healthy values must be modeled by adding/using a condition. Maximum right lung ratio allowed is " << rightLungRatioMax << ".";
      Error(ss);
      err = true;
    }
    else if (rightLungRatio < rightLungRatioMin)
    {
      ss << "Patient right lung ratio of " << rightLungRatio << " is too low. Non-healthy values must be modeled by adding/using a condition. Minimum right lung ratio allowed is " << rightLungRatioMin << ".";
      Error(ss);
      err = true;
    }

    //Respiratory Volumes ---------------------------------------------------------------
    //These are based on weight
    /// \cite ganong1995review
    double totalLungCapacity_L;
    double computedTotalLungCapacity_L = 80.0 * idealWeight_kg / 1000.0;
    if (!m_InitialPatient->HasTotalLungCapacity())
    {
      totalLungCapacity_L = computedTotalLungCapacity_L;
      m_InitialPatient->GetTotalLungCapacity().SetValue(totalLungCapacity_L, VolumeUnit::L);
      ss << "No patient total lung capacity set. Using a computed value of " << computedTotalLungCapacity_L << " L.";
      Info(ss);
    }
    totalLungCapacity_L = m_InitialPatient->GetTotalLungCapacity(VolumeUnit::L);
    if (totalLungCapacity_L != computedTotalLungCapacity_L)
    {
      ss << "Specified total lung capacity of " << totalLungCapacity_L << " L differs from computed value of " << computedTotalLungCapacity_L << " L. No guarantees of model validity.";
      Warning(ss);
    }

    double functionalResidualCapacity_L;
    double computedFunctionalResidualCapacity_L = 30.0 * idealWeight_kg / 1000.0;
    if (!m_InitialPatient->HasFunctionalResidualCapacity())
    {
      functionalResidualCapacity_L = computedFunctionalResidualCapacity_L;
      m_InitialPatient->GetFunctionalResidualCapacity().SetValue(functionalResidualCapacity_L, VolumeUnit::L);
      ss << "No patient functional residual capacity set. Using a computed value of " << computedFunctionalResidualCapacity_L << " L.";
      Info(ss);
    }
    functionalResidualCapacity_L = m_InitialPatient->GetFunctionalResidualCapacity(VolumeUnit::L);
    if (functionalResidualCapacity_L != computedFunctionalResidualCapacity_L)
    {
      ss << "Specified functional residual capacity of " << functionalResidualCapacity_L << " L differs from computed value of " << computedFunctionalResidualCapacity_L << " L. No guarantees of model validity.";
      Warning(ss);
    }

    double residualVolume_L;
    double computRedesidualVolume_L = 16.0 * idealWeight_kg / 1000.0;
    if (!m_InitialPatient->HasResidualVolume())
    {
      residualVolume_L = computRedesidualVolume_L;
      m_InitialPatient->GetResidualVolume().SetValue(residualVolume_L, VolumeUnit::L);
      ss << "No patient residual volume set. Using a computed value of " << computRedesidualVolume_L << " L.";
      Info(ss);
    }
    residualVolume_L = m_InitialPatient->GetResidualVolume(VolumeUnit::L);
    if (residualVolume_L != computRedesidualVolume_L)
    {
      ss << "Specified residual volume of " << residualVolume_L << " L differs from computed value of " << computRedesidualVolume_L << " L. No guarantees of model validity.";
      Warning(ss);
    }

    if (m_InitialPatient->HasTidalVolumeBaseline())
    {
      ss << "Patient tidal volume baseline cannot be set. It is determined through homeostatic simulation.";
      Error(ss);
      err = true;
    }
    if (m_InitialPatient->HasVitalCapacity())
    {
      ss << "Patient vital capacity cannot be set. It is directly computed via other lung volume patient parameters.";
      Error(ss);
      err = true;
    }
    if (m_InitialPatient->HasExpiratoryReserveVolume())
    {
      ss << "Patient expiratory reserve volume cannot be set. It is directly computed via other lung volume patient parameters.";
      Error(ss);
      err = true;
    }
    if (m_InitialPatient->HasInspiratoryReserveVolume())
    {
      ss << "Patient inspiratory reserve volume cannot be set. It is directly computed via other lung volume patient parameters.";
      Error(ss);
      err = true;
    }
    if (m_InitialPatient->HasInspiratoryCapacity())
    {
      ss << "Patient inspiratory capacity cannot be set. It is directly computed via other lung volume patient parameters.";
      Error(ss);
      err = true;
    }

    double tidalVolume_L = 37.0 * idealWeight_kg / 1000.0 - functionalResidualCapacity_L;
    double vitalCapacity = totalLungCapacity_L - residualVolume_L;
    double expiratoryReserveVolume = functionalResidualCapacity_L - residualVolume_L;
    double inspiratoryReserveVolume = totalLungCapacity_L - functionalResidualCapacity_L - tidalVolume_L;
    double inspiratoryCapacity = totalLungCapacity_L - functionalResidualCapacity_L;
    //No negative volumes
    if (totalLungCapacity_L < 0.0 || functionalResidualCapacity_L < 0.0 || residualVolume_L < 0.0 || tidalVolume_L < 0.0 ||
      vitalCapacity < 0.0 || expiratoryReserveVolume < 0.0 || inspiratoryReserveVolume < 0.0 || inspiratoryCapacity < 0.0)
    {
      ss << "All patient lung volumes must be positive.";
      Error(ss);
      err = true;
    }
    m_InitialPatient->GetTidalVolumeBaseline().SetValue(tidalVolume_L, VolumeUnit::L); //This is overwritten after stabilization  
    ss << "Patient tidal volume computed and set to " << tidalVolume_L << " L.";
    Info(ss);

    m_InitialPatient->GetVitalCapacity().SetValue(vitalCapacity, VolumeUnit::L);
    ss << "Patient vital capacity computed and set to " << vitalCapacity << " L.";
    Info(ss);

    m_InitialPatient->GetExpiratoryReserveVolume().SetValue(expiratoryReserveVolume, VolumeUnit::L);
    ss << "Patient expiratory reserve volume computed and set to " << expiratoryReserveVolume << " L.";
    Info(ss);

    m_InitialPatient->GetInspiratoryReserveVolume().SetValue(inspiratoryReserveVolume, VolumeUnit::L);
    ss << "Patient inspiratory reserve volume computed and set to " << inspiratoryReserveVolume << " L.";
    Info(ss);

    m_InitialPatient->GetInspiratoryCapacity().SetValue(inspiratoryCapacity, VolumeUnit::L);
    ss << "Patient inspiratory capacity computed and set to " << inspiratoryCapacity << " L.";
    Info(ss);

    //Alveoli Surface Area ---------------------------------------------------------------
    /// \cite roberts2000gaseous
    double standardAlveoliSurfaceArea_m2 = 70.0;
    double alveoliSurfaceArea_m2;
    //Scale the alveoli surface area based on the size of the patient�s lungs  
    /// \cite ganong1995review
    double standardTotalLungCapacity_L = 6.17; //This is the Total Lung Capacity of our standard patient  
    double computedAlveoliSurfaceArea_m2 = totalLungCapacity_L / standardTotalLungCapacity_L * standardAlveoliSurfaceArea_m2;
    if (!m_InitialPatient->HasAlveoliSurfaceArea())
    {
      alveoliSurfaceArea_m2 = computedAlveoliSurfaceArea_m2;
      m_InitialPatient->GetAlveoliSurfaceArea().SetValue(alveoliSurfaceArea_m2, AreaUnit::m2);
      ss << "No patient alveoli surface area set. Using a computed value of " << computedAlveoliSurfaceArea_m2 << " m^2.";
      Info(ss);
    }
    alveoliSurfaceArea_m2 = m_InitialPatient->GetAlveoliSurfaceArea(AreaUnit::m2);
    if (alveoliSurfaceArea_m2 != computedAlveoliSurfaceArea_m2)
    {
      ss << "Specified alveoli surface area of " << alveoliSurfaceArea_m2 << " m^2 differs from computed value of " << computedAlveoliSurfaceArea_m2 << " m^2. No guarantees of model validity.";
      Warning(ss);
    }

    //Skin Surface Area ---------------------------------------------------------------
    /// \cite du1989formula
    double skinSurfaceArea_m2;
    double computSkinSurfaceArea_m2 = 0.20247 * pow(weight_kg, 0.425) * pow(Convert(height_cm, LengthUnit::cm, LengthUnit::m), 0.725);
    if (!m_InitialPatient->HasSkinSurfaceArea())
    {
      skinSurfaceArea_m2 = computSkinSurfaceArea_m2;
      m_InitialPatient->GetSkinSurfaceArea().SetValue(skinSurfaceArea_m2, AreaUnit::m2);
      ss << "No patient skin surface area set. Using a computed value of " << computSkinSurfaceArea_m2 << " m^2.";
      Info(ss);
    }
    skinSurfaceArea_m2 = m_InitialPatient->GetSkinSurfaceArea(AreaUnit::m2);
    if (skinSurfaceArea_m2 != computSkinSurfaceArea_m2)
    {
      ss << "Specified skin surface area of " << skinSurfaceArea_m2 << " cm differs from computed value of " << computSkinSurfaceArea_m2 << " cm. No guarantees of model validity.";
      Warning(ss);
    }

    //Basal Metabolic Rate ---------------------------------------------------------------
    //The basal metabolic rate is determined from the Harris-Benedict formula, with differences dependent on sex, age, height and mass
    /// \cite roza1984metabolic
    double BMR_kcal_Per_day;
    double computBMR_kcal_Per_day = 88.632 + 13.397 * weight_kg + 4.799 * height_cm - 5.677 * age_yr; //Male
    if (m_InitialPatient->GetSex() == ePatient_Sex::Female)
    {
      computBMR_kcal_Per_day = 447.593 + 9.247 * weight_kg + 3.098 * height_cm - 4.330 * age_yr; //Female
    }
    if (!m_InitialPatient->HasBasalMetabolicRate())
    {
      BMR_kcal_Per_day = computBMR_kcal_Per_day;
      m_InitialPatient->GetBasalMetabolicRate().SetValue(BMR_kcal_Per_day, PowerUnit::kcal_Per_day);

      ss << "No patient basal metabolic rate set. Using a computed value of " << computBMR_kcal_Per_day << " kcal/day.";
      Info(ss);
    }
    BMR_kcal_Per_day = m_InitialPatient->GetBasalMetabolicRate(PowerUnit::kcal_Per_day);
    if (BMR_kcal_Per_day != computBMR_kcal_Per_day)
    {
      ss << "Specified basal metabolic rate of " << BMR_kcal_Per_day << " kcal/day differs from computed value of " << computBMR_kcal_Per_day << " kcal/day. No guarantees of model validity.";
      Warning(ss);
    }

    // Copy the initial patient to the current patient
    m_CurrentPatient->Copy(*m_InitialPatient);

    if (err)
      return false;
    return true;
  }
END_NAMESPACE_EX
