Changeset 19332

Timestamp:

10/30/18 15:39:36 (6 years ago)

Author:

tbretz

Message:

Implemented variable block length to support thinning option, read more data values from the headers and store and print them.

Location:

trunk/Mars/mcorsika

Files:

• MCorsikaRunHeader.cc (modified) (9 diffs)
• MCorsikaRunHeader.h (modified) (4 diffs)

trunk/Mars/mcorsika/MCorsikaRunHeader.cc

@@ -47 +47 @@
 //  + UInt_t  fNumReuse
 //
+// Class Version 4:
+// ----------------
+//  + UInt_t  fCerenkovFileOption
+//  + Float_t fEnergyCutoffHadrons
+//  + Float_t fEnergyCutoffMuons
+//  + Float_t fEnergyCutoffElectrons
+//  + Float_t fEnergyCutoffPhotons
+//  + Float_t fThinningEnergyFractionH
+//  + Float_t fThinningEnergyFractionEM
+//  + Float_t fThinningWeightLimitH
+//  + Float_t fThinningWeightLimitEM
+//  + Float_t fThinningMaxRadius
+//
 ////////////////////////////////////////////////////////////////////////////
 
@@ -83 +96 @@
 // Read in one run header from the binary file
 //
-Bool_t MCorsikaRunHeader::ReadEvt(MCorsikaFormat * fInFormat)
-{
-    Float_t f[272];
-    if (!fInFormat->Read(f, 272 * sizeof(Float_t)))
+Bool_t MCorsikaRunHeader::ReadEvt(MCorsikaFormat * fInFormat, const uint32_t &blockLength)
+{
+    vector<Float_t> f(blockLength);
+    if (!fInFormat->Read(f.data(), blockLength))
         return kFALSE;
 
@@ -104 +117 @@
 
     memset(fObsLevel, 0, 10*4);
-    memcpy(fObsLevel, f+4, fNumObsLevel*4);
+    memcpy(fObsLevel, f.data()+4, fNumObsLevel*4);
 
     fSlopeSpectrum  = f[14];
@@ -110 +123 @@
     fEnergyMax      = f[16];
 
+    fEnergyCutoffHadrons   = f[17];
+    fEnergyCutoffMuons     = f[18];
+    fEnergyCutoffElectrons = f[19];
+    fEnergyCutoffPhotons   = f[20];
+
     // Implemented in CORSIKA Version >= 6.822
     fImpactMax = -1;
@@ -129 +147 @@
 
     // Implemented in CORSIKA Version >= 6.822
-    memcpy(fAtmosphericLayers, f+248, 5*4);
-
-    memcpy(fAtmosphericCoeffA, f+253, 5*4);
-    memcpy(fAtmosphericCoeffB, f+258, 5*4);
-    memcpy(fAtmosphericCoeffC, f+263, 5*4);
+    memcpy(fAtmosphericLayers, f.data()+248, 5*4);
+
+    memcpy(fAtmosphericCoeffA, f.data()+253, 5*4);
+    memcpy(fAtmosphericCoeffB, f.data()+258, 5*4);
+    memcpy(fAtmosphericCoeffC, f.data()+263, 5*4);
 
     return kTRUE;
@@ -178 +196 @@
     // WITH rounding: unbelievable!
     fCerenkovFlag = TMath::Nint(g[75]);
+    fCerenkovFileOption = TMath::Nint(g[90]);
 
     fZdMin = g[79];                // lower edge of theta in °
@@ -207 +226 @@
     fWavelengthMax = g[95];        // Cherenkov bandwidth upper end in nm
 
+    fThinningEnergyFractionH  = g[146]; // EFRCTHN
+    fThinningEnergyFractionEM = g[147]; // EFRCTHN*THINRAT
+    fThinningWeightLimitH     = g[148]; // WMAX
+    fThinningWeightLimitEM    = g[149]; // WMAX*WEITRAT
+    fThinningMaxRadius        = g[150]; // Max radial radius for thinning
+
     fViewConeInnerAngle = g[151];  // inner angle of view cone (°)
     fViewConeOuterAngle = g[152];  // outer angle of view cone (°)
@@ -285 +310 @@
     if (fImpactMax>0)
         *fLog << "Max.sim.Impact: " << fImpactMax << "cm" << endl;
+
+    *fLog << "Energy cutoff:  ";
+    *fLog << fEnergyCutoffHadrons   << "GeV (hadrons), ";
+    *fLog << fEnergyCutoffMuons     << "GeV (muons), ";
+    *fLog << fEnergyCutoffElectrons << "GeV (electrons), ";
+    *fLog << fEnergyCutoffPhotons   << "GeV (photons)";
+    *fLog << endl;
+
+    *fLog << "Thinning:       ";
+    if (fThinningWeightLimitH>0)
+    {
+        *fLog << "HADRONIC: E/Eth>" << fThinningEnergyFractionH  << " (w>" << fThinningWeightLimitH  << "), ";
+        *fLog << "EM: E/Eth>"       << fThinningEnergyFractionEM << " (w>" << fThinningWeightLimitEM << "), ";
+        *fLog << "R>" << fThinningMaxRadius << "cm";
+        *fLog << endl;
+    }
+    else
+        *fLog << "<off>" << endl;
 
     *fLog << "Options used:  ";
@@ -305 +348 @@
     *fLog << " [" << hex << fCerenkovFlag << "]" << dec << endl;
 
+    if (Has(kCerenkov))
+    {
+        *fLog << "File format:    ";
+        switch (fCerenkovFileOption)
+        {
+        case 0:
+            *fLog << "Cerenkov photons written to DAT-file.";
+            break;
+        case 1:
+            *fLog << "Cerenkov photons written to CER-file";
+            break;
+        case 2:
+            *fLog << "Cerenkov photons written to CER-file / Wavelength as 8th item in THIN option";
+            break;
+        default:
+            *fLog << "Cerenkov photons written to CER-file / Prod. height replaced by distance to array center.";
+            break;
+        }
+        *fLog << " [MCERFI=" << fCerenkovFileOption << "]" << endl;
+    }
+
     if (HasLayers())
     {

trunk/Mars/mcorsika/MCorsikaRunHeader.h

@@ -69 +69 @@
 
     UInt_t fCerenkovFlag;
+    UInt_t fCerenkovFileOption;       // MCERFI
+
+    Float_t fEnergyCutoffHadrons;     // [GeV]
+    Float_t fEnergyCutoffMuons;       // [GeV]
+    Float_t fEnergyCutoffElectrons;   // [GeV]
+    Float_t fEnergyCutoffPhotons;     // [GeV]
+
+    Float_t fThinningEnergyFractionH;  // Hadronic energy limit: EFRCTHN
+    Float_t fThinningEnergyFractionEM; // EM energy limit: EFRCTHN*THINRAT
+    Float_t fThinningWeightLimitH;     // Hadronic weight limit: WMAX
+    Float_t fThinningWeightLimitEM;    // EM weight limit: WMAX*WEITRAT
+    Float_t fThinningMaxRadius;        // [cm] Max radial raius for thinning
 
 public:
@@ -110 +122 @@
     Bool_t Has(CerenkovFlag_t opt) const { return fCerenkovFlag&opt ? 1 : 0; }
 
+    UInt_t GetCerenkovFileOption() const { return fCerenkovFileOption; }
+
     static Double_t EarthRadius() { return fgEarthRadius; }
 
@@ -123 +137 @@
 
     // I/O
-    Bool_t ReadEvt(MCorsikaFormat * fInFormat);
+    Bool_t ReadEvt(MCorsikaFormat * fInFormat, const uint32_t &blockLength);
     Bool_t ReadEventHeader(Float_t * g);
     Bool_t ReadEvtEnd(MCorsikaFormat * fInFormat, Bool_t runNumberVerify);
@@ -130 +144 @@
     void Print(Option_t *t=NULL) const;
 
-    ClassDef(MCorsikaRunHeader, 3)      // storage container for general info
+    ClassDef(MCorsikaRunHeader, 4)      // storage container for general info
 };
 #endif