// ESAF : Euso Simulation and Analysis Framework
// $Id: LowtranAtmosphere.cc,v 1.12 2005/02/10 17:59:23 moreggia Exp $
// Corinne Berat created Apr, 2 2004
/*****************************************************************************
* ESAF: Euso Simulation and Analysis Framework *
* *
* Id: LowtranAtmosphere *
* Package: atmosphere *
* Coordinator: S. Moreggia, C. Berat *
* *
*****************************************************************************/
//_____________________________________________________________________________
//
// LowtranAtmosphere
//
// <extensive class description>
//
// Config file parameters
// ======================
//
// <parameter name>: <parameter description>
// -Valid options: <available options>
//
#include "LowtranAtmosphere.hh"
#include "Config.hh"
#include "LowtranAtmosphereData.hh"
#include <TH1F.h>
#include <TROOT.h>
#include <math.h>
ClassImp(LowtranAtmosphere)
//_____________________________________________________________________________
LowtranAtmosphere::LowtranAtmosphere(string type) : Atmosphere(), fData(0) {
//
// ctor
//
Msg(EsafMsg::Info) << "Building with the configuration:"<< MsgDispatch;
fType = type;
Build();
}
//_____________________________________________________________________________
LowtranAtmosphere::~LowtranAtmosphere() {
//
// ctor
//
SafeDelete(fData);
}
//___________________________________________________________________________________
void LowtranAtmosphere::CreateInstance() {
//
// Create the single atmosphere instance as being a LowtranAtmosphere object
//
if(!Atmosphere::fChild) Atmosphere::fChild = new LowtranAtmosphere("lowtran");
else Atmosphere::fChild->Msg(EsafMsg::Warning) << "Trying to create two atmospheres !!"<< MsgDispatch;
}
//_____________________________________________________________________________
void LowtranAtmosphere::Build() {
//
// Build atmosphere according to config parameters
//
fData = new LowtranAtmosphereData();
if(!fData) Msg(EsafMsg::Panic) << "No new LowtranData, memory pb" << MsgDispatch ;
}
//_____________________________________________________________________________
Double_t LowtranAtmosphere::Interpolate(string& val,string& mod,Double_t h) const {
//
// Calculate interpolation
//
const Double_t* data;
const Double_t* dalt;
dalt = fData->GetAltitudeTable();
if(val == "press")
data = fData->GetPressureTable();
else if(val == "temp")
data = fData->GetTemperatureTable();
else if(val == "Airdens")
data = fData->GetAir_DensityTable();
else if(val == "Odens")
data = fData->GetO_DensityTable();
else if(val == "O2dens")
data = fData->GetO2_DensityTable();
else if(val == "N2dens")
data = fData->GetN2_DensityTable();
else
Msg(EsafMsg::Panic)<< "Wrong val argument in LOWTRAN_00Atmosphere::Interpolate"
<< MsgDispatch ;
if(h <0) return data[0];
Int_t nb = fData->NumberOfElements();
if(h > 120*km) return data[nb-1];
Double_t val1, val2;
Double_t h1;
Double_t step;
Int_t i = 0;
if(h < 25*km) {
step = 1*km;
i = Int_t(floor(h/step));
}
else if(h < 50*km) {
step = 2.5*km;
i = 25 + Int_t( floor((h-25*km)/step) );
}
else {
step = 5*km;
i = 35 + Int_t( floor((h-50*km)/step) );
}
h1 = dalt[i];
val1 = data[i];
val2 = data[i+1];
if(mod == "linear" || val1 == 0 || val2 == 0)
return (val2 - val1)*(h - h1)/step + val1;
else if(mod == "expon")
return val1 * exp(-(h - h1) * log(val1/val2)/step);
else
Msg(EsafMsg::Panic)<<"Wrong mod argument in LOWTRAN_00Atmosphere::Interpolate"
<< MsgDispatch ;
return 0;
}
//_____________________________________________________________________________
Double_t LowtranAtmosphere::Pressure(Double_t h) const {
//
// Pressure at a given altitude
//
if(h > 100*km)
return 0.;
// linear interpolation between table values
string val = "press";
string mod = "expon";
return Interpolate(val,mod,h);
}
//_____________________________________________________________________________
Double_t LowtranAtmosphere::Temperature(Double_t h) const {
//
// Temperature at a given altitude
//
if(h > 100*km)
return 0.;
// linear interpolation between table values
string val = "temp";
string mod = "linear";
return Interpolate(val,mod,h);
}
//_____________________________________________________________________________
Double_t LowtranAtmosphere::Air_Density(Double_t h) const {
//
// Atmospheric air density at a given altitude
//
if(h > 100*km)
return 0.;
// exponential interpolation between table values
string val = "Airdens";
string mod = "expon";
return Interpolate(val,mod,h);
}
//_____________________________________________________________________________
Double_t LowtranAtmosphere::O_Density(Double_t h) const {
//
// Oxygen density at a given altitude
//
if(h > 100*km)
return 0.;
// exponential interpolation between table values
string val = "Odens";
string mod = "linear";
return Interpolate(val,mod,h);
}
//_____________________________________________________________________________
Double_t LowtranAtmosphere::O2_Density(Double_t h) const {
//
// Di-oxygen density at a given altitude
//
if(h > 100*km)
return 0.;
// exponential interpolation between table values
string val = "O2dens";
string mod = "linear";
return Interpolate(val,mod,h);
}
//_____________________________________________________________________________
Double_t LowtranAtmosphere::O3_Density(Double_t h) const {
//
// Ozone density at a given altitude
//
if(h > 100*km)
return 0.;
// exponential interpolation between table values
string val = "O3dens";
string mod = "linear";
return Interpolate(val,mod,h);
}
//_____________________________________________________________________________
Double_t LowtranAtmosphere::N2_Density(Double_t h) const {
//
// Di-nitrogen density at a given altitude
//
if(h > 100*km)
return 0.;
// exponential interpolation between table values
string val = "N2dens";
string mod = "linear";
return Interpolate(val,mod,h);
}