Changeset 9347 for trunk/MagicSoft/Mars/mcorsika

Timestamp:

02/18/09 12:07:41 (16 years ago)

Author:

tbretz

Message:

*** empty log message ***

Location:

trunk/MagicSoft/Mars/mcorsika

Files:

• MCorsikaRunHeader.cc (modified) (8 diffs)
• MCorsikaRunHeader.h (modified) (3 diffs)

trunk/MagicSoft/Mars/mcorsika/MCorsikaRunHeader.cc

@@ -37 +37 @@
 //  + Float_t fMagneticFieldZ
 //  + Float_t fMagneticFieldAz
+//  + Float_t fAtmosphericCoeffA[5]
+//  + Float_t fAtmosphericCoeffB[5]
+//  + Float_t fAtmosphericCoeffC[5]
+//  + UInt_t  fCerenkovFlag
 //
 ////////////////////////////////////////////////////////////////////////////
@@ -53 +57 @@
 
 using namespace std;
+
+const Double_t MCorsikaRunHeader::fgEarthRadius = 637131500; // [cm] Earth radius as defined in CORSIKA
 
 // --------------------------------------------------------------------------
@@ -81 +87 @@
     }
 
-    Float_t f[68];
-    fin.read((char*)f, 68);
+    Float_t f[272*4];
+    fin.read((char*)f, 272*4);
 
     fRunNumber = TMath::Nint(f[0]);
@@ -105 +111 @@
     fEnergyMax      = f[16];
 
-    fin.seekg(1020, ios::cur);     // skip the remaining data of this block
+    // 0/10 not supported?
+
+    memcpy(fAtmosphericCoeffA, f+253, 5*4);
+    memcpy(fAtmosphericCoeffB, f+258, 5*4);
+    memcpy(fAtmosphericCoeffC, f+263, 5*4);
 
     // -------------------- Read first event header -------------------
@@ -157 +167 @@
     fMagneticFieldAz =  g[91];  // Azimuth angle of magnetic north expressed in telescope coordinates
 
+    // WITH rounding: unbelievable!
+    fCerenkovFlag = TMath::Nint(g[75]);
+
     fZdMin = g[79];                // lower edge of theta in °
     fZdMax = g[80];                // upper edge of theta in °
@@ -173 +186 @@
     fImpactMax = g[86];
 
-    fWavelengthMin = g[95];        // Cherenkov bandwidth lower end in nm
-    fWavelengthMax = g[96];        // Cherenkov bandwidth upper end in nm
+    fWavelengthMin = g[94];        // Cherenkov bandwidth lower end in nm
+    fWavelengthMax = g[95];        // Cherenkov bandwidth upper end in nm
 
     fViewConeInnerAngle = g[151];  // inner angle of view cone (°)
@@ -251 +264 @@
 
     *fLog << "Spectrum:      Slope=" << fSlopeSpectrum << "  (" << fEnergyMin << "GeV-" << fEnergyMax << "GeV)" <<  endl;
+    *fLog << "Wavelength:    " << fWavelengthMin << "nm - " << fWavelengthMax << "nm" << endl;
 
     if (fViewConeOuterAngle>0)
@@ -269 +283 @@
         *fLog << endl;
     }
-}
-
+
+    *fLog << "Atm.Coeff A:  ";
+    for (int i=0; i<5; i++)
+        *fLog << " " << fAtmosphericCoeffA[i];
+    *fLog << endl;
+    *fLog << "Atm.Coeff B:  ";
+    for (int i=0; i<5; i++)
+        *fLog << " " << fAtmosphericCoeffB[i];
+    *fLog << endl;
+    *fLog << "Atm.Coeff C:  ";
+    for (int i=0; i<5; i++)
+        *fLog << " " << fAtmosphericCoeffC[i];
+    *fLog << endl;
+}
+

trunk/MagicSoft/Mars/mcorsika/MCorsikaRunHeader.h

@@ -14 +14 @@
 {
     friend class MCorsikaEvtHeader;
+public:
+    enum CerenkovFlag_t
+    {
+        kCerenkov   = BIT(0),
+        kIact       = BIT(1),
+        kCeffic     = BIT(2),
+        kAtmext     = BIT(3),
+        kRefraction = BIT(4),
+        kVolumedet  = BIT(5),
+        kCurved     = BIT(6),
+        kSlant      = BIT(8)
+    };
 
 private:
+    static const Double_t fgEarthRadius; // Take same Earth radius as in CORSIKA (cm)
+
     UInt_t  fRunNumber;               // Run number
     UInt_t  fParticleID;              // Particle ID (see MMcEvtBasic or CORSIKA manual)
@@ -45 +59 @@
     Float_t fViewConeInnerAngle;      // [deg]
     Float_t fViewConeOuterAngle;      // [deg]
+
+    Float_t fAtmosphericCoeffA[5];    // [g/cm²] AATM (see Corsika Manual for details)
+    Float_t fAtmosphericCoeffB[5];    // [g/cm²] BATM (see Corsika Manual for details)
+    Float_t fAtmosphericCoeffC[5];    // [cm]    CATM (see Corsika Manual for details)
+
+    UInt_t fCerenkovFlag;
 
 public:
@@ -80 +100 @@
     Bool_t HasViewCone() const { return fViewConeOuterAngle>0; }
 
+    Float_t GetObsLevel(UInt_t i=0) const { return i>9 ? -1 : fObsLevel[i]; }
+
+    Bool_t Has(CerenkovFlag_t opt) { return fCerenkovFlag&opt ? 1 : 0; }
+
+    static Double_t EarthRadius() { return fgEarthRadius; }
+
+    // Preliminary!
+    const Float_t *GetAtmosphericCoeffA() const { return fAtmosphericCoeffA; }
+    const Float_t *GetAtmosphericCoeffB() const { return fAtmosphericCoeffB; }
+    const Float_t *GetAtmosphericCoeffC() const { return fAtmosphericCoeffC; }
+
+    UInt_t GetNumAtmosphericModel() const { return (fCerenkovFlag>>10)&0x3ff; }
+
     // I/O
     Bool_t ReadEvt(istream& fin);