Changeset 10060 for trunk/Mars/msim

Timestamp:

11/29/10 14:24:23 (14 years ago)

Author:

rohlfs

Message:

two new command line arguments of readcorsika: -A=arrayNo and -T=telescopeNo. New design of program flow in MCorsikaRead: It is now determined by the order of the data blocks in the input file.

Location:

trunk/Mars/msim

Files:

• MPhotonData.cc (modified) (2 diffs)
• MPhotonEvent.cc (modified) (1 diff)
• MPhotonEvent.h (modified) (1 diff)

trunk/Mars/msim/MPhotonData.cc

@@ -222 +222 @@
 {
     const UInt_t n = TMath::Nint(f[0]);
+
     if (n==0)
         // FIXME: Do we need to decode the rest anyway?
@@ -227 +228 @@
 
     // Check reuse
-    const Int_t reuse = (n/1000)%100; // Force this to be 1!
-    if (reuse!=i)
-        return kCONTINUE;
+    if (i >=0)
+      {
+       const Int_t reuse = (n/1000)%100; // Force this to be 1!
+       if (reuse!=i)
+           return kCONTINUE;
+      }
 
     // This seems to be special to mmcs

trunk/Mars/msim/MPhotonEvent.cc

@@ -507 +507 @@
         operator[](i).SimWavelength(wmin, wmax);
 }
-
-// --------------------------------------------------------------------------
-//
-// Read the Event section from the file
-//
-Int_t MPhotonEvent::ReadCorsikaEvt(MCorsikaFormat *fInFormat, Int_t id)
-{
-    Int_t n = 0;
-
-    // --- old I/O ---
-    // Read only + Reflector (no absorption)
-    // Muons:   1.06GB/115s =  9.2MB/s (100kEvs)
-    // Gammas:  1.57GB/275s =  5.7MB/s (  1kEvs)
-
-    // Read only:
-    // Gammas:  1.57GB/158s =  9.9MB/s (  1kEvs)
-    // Muons:   1.06GB/ 77s = 13.8MB/s (100kEvs)
-
-    // --- new I/O ---
-    // Read only (don't allocate storage space):
-    // Gammas:  1.57GB/143s = 11.0MB/s (  1kEvs)
-    // Muons:   1.06GB/ 77s = 13.8MB/s (100kEvs)
-
-    // Read only in blocks (with storage allocation):
-    // Gammas:  1.57GB/28s  =  56MB/s (  1kEvs)
-    // Muons:   1.06GB/5.2s = 204MB/s (100kEvs)
-
-    // Read only in blocks (without storage allocation):
-    // similar to just copy
-
-    // Copy with cp
-    // 1.57GB/ 5s   CPU
-    // 1.57GB/28s   REAL
-    // 1.06GB/ 3s   CPU
-    // 1.06GB/22s   REAL
-    Float_t *buffer = 0;
-
-    if (fInFormat->IsEventioFormat())
-    {
-        while (1)
-        {
-            const Int_t rc = fInFormat->GetNextEvent(&buffer, id);
-            if (rc==kERROR)
-                return kERROR;
-            if (rc==kFALSE)
-                break;
-
-            // Loop over number of photons in bunch
-            while (Add(n).FillEventIO(buffer))
-                n++;
-        }
-    }
-    else
-    {
-        while (1)
-        {
-            const Int_t rc1 = fInFormat->GetNextEvent(&buffer);
-            if (rc1==kERROR)
-                return kERROR;
-            if (rc1==kFALSE)
-                break;
-
-            const Int_t rc2 = Add(n).FillCorsika(buffer, id);
-            switch (rc2)
-            {
-            case kCONTINUE:  continue;        // No data in this bunch... skip it.
-            case kERROR:     return kERROR;   // Error occured
-            //case kFALSE:     return kFALSE;   // End of stream
-            }
-
-            // This is a photon we would like to keep later.
-            // Increase the counter by one
-            n++;
-        }
-    }
-
-/*
-
-    while (1)
-    {
-        Float_t buffer[273];
-        Float_t * ptr = buffer;
-
-
-        if (!fInFormat->ReadData(273, buffer))
-            return kFALSE;
-
-        if (!memcmp(ptr, "EVTE", 4))
-            {
-
-            fInFormat->UnreadLastData();
-            break;
-            }
-
-        Float_t *end = ptr + 273;
-
-        // Loop over all sub-blocks (photons) in the block and if
-        // they contain valid data add them to the array
-        while (ptr<end)
-        {
-            // Get/Add the n-th entry from the array and
-            // fill it with the current 7 floats
-            const Int_t rc = Add(n).FillCorsika(ptr);
-            ptr += 7;
-
-            switch (rc)
-            {
-            case kCONTINUE:  continue;        // No data in this bunch... skip it.
-            case kERROR:     return kERROR;   // Error occured
-            //case kFALSE:     return kFALSE;   // End of stream
-            }
-
-            // This is a photon we would like to keep later.
-            // Increase the counter by one
-            n++;
-        }
-    }
-
-*/
-    Resize(n);
-    fData.UnSort();
-
-    SetReadyToSave();
-
-    //*fLog << all << "Number of photon bunches: " << fData.GetEntriesFast() << endl;
-    return kTRUE;
-
-}
-
-Int_t MPhotonEvent::ReadCorsikaEvt(istream &fin, Int_t i)
-{
-    Int_t n = 0;
-
-    // --- old I/O ---
-    // Read only + Reflector (no absorption)
-    // Muons:   1.06GB/115s =  9.2MB/s (100kEvs)
-    // Gammas:  1.57GB/275s =  5.7MB/s (  1kEvs)
-
-    // Read only:
-    // Gammas:  1.57GB/158s =  9.9MB/s (  1kEvs)
-    // Muons:   1.06GB/ 77s = 13.8MB/s (100kEvs)
-
-    // --- new I/O ---
-    // Read only (don't allocate storage space):
-    // Gammas:  1.57GB/143s = 11.0MB/s (  1kEvs)
-    // Muons:   1.06GB/ 77s = 13.8MB/s (100kEvs)
-
-    // Read only in blocks (with storage allocation):
-    // Gammas:  1.57GB/28s  =  56MB/s (  1kEvs)
-    // Muons:   1.06GB/5.2s = 204MB/s (100kEvs)
-
-    // Read only in blocks (without storage allocation):
-    // similar to just copy
-
-    // Copy with cp
-    // 1.57GB/ 5s   CPU
-    // 1.57GB/28s   REAL
-    // 1.06GB/ 3s   CPU
-    // 1.06GB/22s   REAL
-
-    while (1)
-    {
-        // Stprage for one block
-        char c[273*4];
-
-        // Read the first four byte to check whether the next block
-        // doen't belong to the event anymore
-        fin.read(c, 4);
-        if (!fin)
-            return kFALSE;
-
-        // Check if the event is finished
-        if (!memcmp(c, "EVTE", 4))
-            break;
-
-        // Now read the rest of the data
-        fin.read(c+4, 272*4);
-
-        Float_t *ptr = reinterpret_cast<Float_t*>(c);
-        Float_t *end = ptr + 273;
-
-        // Loop over all sub-blocks (photons) in the block and if
-        // they contain valid data add them to the array
-        while (ptr<end)
-        {
-            // Get/Add the n-th entry from the array and
-            // fill it with the current 7 floats
-            const Int_t rc = Add(n).FillCorsika(ptr, i);
-            ptr += 7;
-
-            switch (rc)
-            {
-            case kCONTINUE:  continue;        // No data in this bunch... skip it.
-            case kERROR:     return kERROR;   // Error occured
-            //case kFALSE:     return kFALSE;   // End of stream
-            }
-
-            // This is a photon we would like to keep later.
-            // Increase the counter by one
-            n++;
-        }
-    }
-/*
-    while (1)
-    {
-        // Check the first four bytes
-        char c[4];
-        fin.read(c, 4);
-
-        // End of stream
-        if (!fin)
-            return kFALSE;
-
-        // Check if we found the end of the event
-        if (!memcmp(c, "EVTE", 4))
-            break;
-
-        // The first for byte contained data already --> go back
-        fin.seekg(-4, ios::cur);
-
-        // Do not modify this. It is optimized for execution
-        // speed and flexibility!
-        MPhotonData &ph = Add(n);
-        // It checks how many entries the lookup table has. If it has enough
-        // entries and the entry was already allocated, we can re-use it,
-        // otherwise we have to allocate it.
-
-        // Now we read a single cherenkov bunch. Note that for speed reason we have not
-        // called the constructor if the event was already constructed (virtual table
-        // set), consequently we must make sure that ReadCorsikaEvent does reset
-        // all data mebers no matter whether they are read or not.
-        const Int_t rc = ph.ReadCorsikaEvt(fin);
-
-        // Evaluate result from reading event
-        switch (rc)
-        {
-        case kCONTINUE:  continue;        // No data in this bunch... skip it.
-        case kFALSE:     return kFALSE;   // End of stream
-        case kERROR:     return kERROR;   // Error occured
-        }
-
-        // FIXME: If fNumPhotons!=1 add the photon more than once
-
-        // Now increase the number of entries which are kept,
-        // i.e. keep this photon(s)
-        n++;
-    }
-  */
-
-    Resize(n);
-    fData.UnSort();
-
-    SetReadyToSave();
-
-    //*fLog << all << "Number of photon bunches: " << fData.GetEntriesFast() << endl;
-    return kTRUE;
-}
-
-// --------------------------------------------------------------------------
-/*
-Int_t MPhotonEvent::ReadRflEvt(std::istream &fin, Int_t i)
-{
-    Int_t n = 0;
-
-    while (1)
-    {
-        // Check the first four bytes
-        char c[13];
-        fin.read(c, 13);
-
-        // End of stream
-        if (!fin)
-            return kFALSE;
-
-        // Check if we found the end of the event
-        if (!memcmp(c, "\nEND---EVENT\n", 13))
-            break;
-
-        // The first for byte contained data already --> go back
-        fin.seekg(-13, ios::cur);
-
-        // Do not modify this. It is optimized for execution
-        // speed and flexibility!
-        //TObject *o = n<fData.GetSize() && fData.UncheckedAt(n) ? fData.UncheckedAt(n) : fData.New(n);
-
-        // Now we read a single cherenkov bunch
-        //const Int_t rc = static_cast<MPhotonData*>(o)->ReadRflEvt(fin);
-        const Int_t rc = Add(n).ReadRflEvt(fin, i);
-
-        // Evaluate result from reading event
-        switch (rc)
-        {
-        case kCONTINUE:  continue;        // No data in this bunch... skip it.
-        case kFALSE:     return kFALSE;   // End of stream
-        case kERROR:     return kERROR;   // Error occured
-        }
-
-        // Now increase the number of entries which are kept,
-        // i.e. keep this photon(s)
-        n++;
-    }
-
-    Shrink(n);
-
-    SetReadyToSave();
-
-    // *fLog << all << "Number of photon bunches: " << fData.GetEntriesFast() << endl;
-    return kTRUE;
-}
-*/
+Int_t MPhotonEvent::ReadEventIoEvt(MCorsikaFormat *fInFormat)
+{
+   Int_t  bunchHeader[3];
+   fInFormat->Read(bunchHeader, 3 * sizeof(Int_t));
+
+   Int_t n = 0;
+
+   for (int bunch = 0; bunch < bunchHeader[2]; bunch++)
+      {
+      Float_t buffer[8];
+      fInFormat->Read(buffer, 8 * sizeof(Float_t));
+
+      if (Add(n).FillEventIO(buffer))
+         n++;
+      }
+
+   Resize(n);
+   fData.UnSort();
+
+   SetReadyToSave();
+
+   //*fLog << all << "Number of photon bunches: " << fData.GetEntriesFast() << endl;
+   return kTRUE;
+
+}
+
+Int_t MPhotonEvent::ReadCorsikaEvt(Float_t * data, Int_t numEvents, Int_t arrayIdx)
+{
+   Int_t n = 0;
+
+   for (Int_t event = 0; event < numEvents; event++)
+      {
+      const Int_t rc = Add(n).FillCorsika(data + 7 * event, arrayIdx);
+
+      switch (rc)
+      {
+      case kCONTINUE:  continue;        // No data in this bunch... skip it.
+      case kERROR:     return kERROR;   // Error occured
+      }
+
+      // This is a photon we would like to keep later.
+      // Increase the counter by one
+      n++;
+      }
+
+   Resize(n);
+   fData.UnSort();
+
+   SetReadyToSave();
+
+   //*fLog << all << "Number of photon bunches: " << fData.GetEntriesFast() << endl;
+   return kTRUE;
+
+}
 
 // --------------------------------------------------------------------------

trunk/Mars/msim/MPhotonEvent.h

@@ -60 +60 @@
 
     // I/O
-    Int_t ReadCorsikaEvt(MCorsikaFormat *fInFormat, Int_t i);
-    Int_t ReadCorsikaEvt(istream &fin, Int_t i);
-    //Int_t ReadRflEvt(istream &fin, Int_t i);
+    Int_t ReadEventIoEvt(MCorsikaFormat *fInFormat);
+    Int_t ReadCorsikaEvt(Float_t * data, Int_t numEvents, Int_t arrayIdx);
 
     // TObject