Changeset 19765


Ignore:
Timestamp:
Oct 12, 2019, 8:06:55 PM (6 days ago)
Author:
tbretz
Message:
fHeight is not used right now... (note that the current default coeffieints also don't match). The loops for ozone and aerosols stopped one entry to early which produced a step in aerosol contents as for all hogher values the last entry was used, aerosol and ozone was only properly initialized if numericall by chance the counting hit the bin border, InitOzone errornesously checked for valid Aerosols, prepared a more reaosnable default consztrucotr: the current one did more or less nothing as the fObLevel was not yet initialized.
Location:
trunk/Mars/msim
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • trunk/Mars/msim/MAtmosphere.cc

    r19763 r19765  
    137137
    138138MAtmRayleigh::MAtmRayleigh() : fR(MCorsikaRunHeader::fgEarthRadius),
    139     fHeight{0, 300000, 1e+06, 5e+06, 1e+07},
     139    //fHeight{0, 300000, 1e+06, 5e+06, 1e+07},
    140140    //fAtmA{-144.838, -124.071, 0.360027, -0.000824761, 0.00221589},
    141141    fAtmB{1192.34, 1173.98, 1412.08, 810.682, 1},
     
    171171        fObsLevel,                     // fObsLevel,                   //   0km
    172172        TMath::Max(fObsLevel, 7.75e5), // TMath::Max(fObsLevel, 4e5),  //   4km
    173         16.5e5,                       //  10e5,                       //  10km
    174         50.0e5,                       //  40e5,                       //  40km
     173        16.5e5,                        //  10e5,                       //  10km
     174        50.0e5,                        //  40e5,                       //  40km
    175175        105.0e5,                       // 100e5                        // 100km
    176176    };
     
    414414
    415415        Double_t path_slant = 0;
    416         for (Double_t h = fObsLevel; h <= 50e5; h += dh)
     416        for (Double_t h=fObsLevel; h<50e5+dh/2; h+=dh)
    417417        {
    418418            // h is the true height vertical above ground
    419             if (fmod(h,1e4) == 0)
    420                 ozone_path[(int)(h/1e4)][j] = path_slant;
     419            //if (fmod(h,1e4) == 0)
     420            {
     421                ozone_path[TMath::FloorNint(h/1e4)][j] = path_slant;
     422            }
    421423
    422424            const Double_t km  = h/1e5;
     
    450452
    451453        Double_t path_slant = 0;
    452         for (Double_t h = fObsLevel; h <= 30e5; h += dh)
     454        for (Double_t h=fObsLevel; h<=30e5+dh/2; h+=dh)
    453455        {
    454456            // h is the true height vertical above ground
    455             if (fmod(h,1e4) == 0)
    456                 aerosol_path[(int)(h/1e4)][j] = path_slant;
     457            //if (fmod(h,1e4) == 0)
     458            {
     459                aerosol_path[TMath::FloorNint(h/1e4)][j] = path_slant;
     460            }
    457461
    458462            const Double_t km  = h/1e5;
     
    482486    }
    483487
    484     if (!HasValidAerosol())
     488    if (!HasValidOzone())
    485489        return kFALSE;
    486490
     
    539543
    540544    //******* Ozone absorption *******
    541     if (h > 50.e5)
    542         h = 50.e5;
    543545
    544546    // Vigroux Ozone absorption coefficient a.s.l. through interpolation:
     
    549551    // Now use the pre-calculated values of the path integral
    550552    // for h and theta
    551     const UInt_t H = TMath::Nint(h/1e4);
    552     const UInt_t T = TMath::Min(89, TMath::Nint(theta/STEPTHETA));
     553    const UInt_t H = TMath::Min(500, TMath::Nint(h/1e4));
     554    const UInt_t T = TMath::Min( 89, TMath::Nint(theta/STEPTHETA));
    553555
    554556    const Double_t path = ozone_path[H][T];
     
    563565
    564566    //******* Mie (aerosol) *******
    565     if (h > 30.e5)
    566         h = 30.e5;
    567567
    568568    //const float aero_betap[15] = {0.27, 0.26, 0.25, 0.24, 0.24, 0.23, 0.20, 0.180, 0.167, 0.158, 0.150, 0.142, 0.135, 0.127, 0.120};
     
    572572    // Now use the pre-calculated values of the path integral
    573573    // for h and theta
    574     const UInt_t H = TMath::Nint(h/1e4);
    575     const UInt_t T = TMath::Min(89, TMath::Nint(theta/STEPTHETA));
    576 
     574    const UInt_t H = TMath::Min(300, TMath::Nint(h/1e4));
     575    const UInt_t T = TMath::Min( 89, TMath::Nint(theta/STEPTHETA));
    577576
    578577    const Double_t path = aerosol_path[H][T];
  • trunk/Mars/msim/MAtmosphere.h

    r19763 r19765  
    9898        Init(h);//, "ozone.txt", "aerosols.txt");
    9999    }
     100    MAtmosphere(Double_t obs=-1, const char *name1=0, const char *name2=0) : fAbsCoeffOzone(0), fAbsCoeffAerosols(0)
     101    {
     102        MAtmRayleigh::Init(obs);
    100103
    101     MAtmosphere(const char *name1=0, const char *name2=0) : fAbsCoeffOzone(0), fAbsCoeffAerosols(0)
    102     {
    103104        if (name1)
    104105            InitOzone(name1);
Note: See TracChangeset for help on using the changeset viewer.