Changeset 7188 for trunk/MagicSoft/Mars/mpedestal

Timestamp:

07/13/05 19:06:26 (20 years ago)

Author:

tbretz

Message:

*** empty log message ***

Location:

trunk/MagicSoft/Mars/mpedestal

Files:

• MExtractPedestal.cc (modified) (14 diffs)
• MExtractPedestal.h (modified) (4 diffs)
• MMcPedestalCopy.cc (modified) (2 diffs)
• MPedCalcPedRun.cc (modified) (7 diffs)
• MPedestalCam.cc (modified) (1 diff)
• MPedestalCam.h (modified) (1 diff)

trunk/MagicSoft/Mars/mpedestal/MExtractPedestal.cc

@@ -175 +175 @@
 
 const TString MExtractPedestal::fgNamePedestalCam = "MPedestalCam";
+const TString MExtractPedestal::fgNameRawEvtData  = "MRawEvtData";
 const UInt_t  MExtractPedestal::fgNumDump = 500;
 
@@ -191 +192 @@
 MExtractPedestal::MExtractPedestal(const char *name, const char *title)
     : fGeom(NULL), fPedestalsIn(NULL), fPedestalsInter(NULL), fPedestalsOut(NULL), 
-      fExtractor(NULL), fExtractWinFirst(0), fExtractWinSize(0)
+      fExtractor(NULL), fExtractWinFirst(0), fExtractWinSize(0), fUseSpecialPixels(kFALSE)
 {
     fName  = name  ? name  : "MExtractPedestal";
@@ -207 +208 @@
     SetNamePedestalCamIn();
     SetNamePedestalCamOut();
+    SetNamePedestalCamInter();
+    SetNameRawEvtData();
     SetNumDump();
 
@@ -311 +314 @@
 
   Clear();
-  
+
+  fRawEvt = (MRawEvtData*)pList->FindObject(AddSerialNumber(fNameRawEvtData));
+  if (!fRawEvt)
+      if (!fUseSpecialPixels)
+      {
+          *fLog << err << AddSerialNumber(fNameRawEvtData) << " not found... aborting." << endl;
+          return kFALSE;
+      }
+
+
   fRawEvt = (MRawEvtData*)pList->FindObject(AddSerialNumber("MRawEvtData"));
   if (!fRawEvt)
@@ -372 +384 @@
 Int_t MExtractPedestal::Process()
 {
+    //
+    // Necessary check for extraction of special pixels
+    // together with data which does not yet have them
+    //
+    if (fSumx.GetSize()==0)
+        return kTRUE;
+
     if (fExtractor)
         fExtractor->SetNoiseCalculation(fRandomCalculation);
@@ -405 +424 @@
 Bool_t MExtractPedestal::ReInit(MParList *pList)
 {
-  // If the size is not yet set, set the size
-  if (fSumx.GetSize()==0)
-  {
-      const Int_t npixels  = fPedestalsOut->GetSize();
-      const Int_t areas    = fPedestalsOut->GetNumAverageArea();
-      const Int_t sectors  = fPedestalsOut->GetNumAverageSector();
-
-      fSumx.  Set(npixels);
-      fSumx2. Set(npixels);
-      fSumAB0.Set(npixels);
-      fSumAB1.Set(npixels);
-
-      fAreaSumx.  Set(areas);
-      fAreaSumx2. Set(areas);
-      fAreaSumAB0.Set(areas);
-      fAreaSumAB1.Set(areas);
-      fAreaFilled.Set(areas);
-      fAreaValid .Set(areas);
-
-      fSectorSumx.  Set(sectors);
-      fSectorSumx2. Set(sectors);
-      fSectorSumAB0.Set(sectors);
-      fSectorSumAB1.Set(sectors);
-      fSectorFilled.Set(sectors);
-      fSectorValid .Set(sectors);
-
-      for (Int_t i=0; i<npixels; i++)
-      {
-          const UInt_t aidx   = (*fGeom)[i].GetAidx();
-          const UInt_t sector = (*fGeom)[i].GetSector();
-
-          fAreaValid  [aidx]  ++;
-          fSectorValid[sector]++;
-      }
-  }
-
-  if (fExtractor)
-  {
-      if (!((MTask*)fExtractor)->ReInit(pList))
-          return kFALSE;
-
-      SetExtractWindow(fExtractor->GetHiGainFirst(), (Int_t)TMath::Nint(fExtractor->GetNumHiGainSamples()));
-  }
-
-  return kTRUE;
+    // Necessary check for special pixels which might not yet have existed
+    if (!fRawEvt)
+    {
+        if (fRunHeader->GetFormatVersion() > 3)
+            return kTRUE;
+
+        *fLog << err << "ERROR - " << fNameRawEvtData << " [MRawEvtData] has not ";
+        *fLog << "been found and format version > 3... abort." << endl;
+        return kFALSE;
+    }
+
+    // If the size is not yet set, set the size
+    if (fSumx.GetSize()==0)
+    {
+        // Initialize the normal pixels (size of MPedestalCam already set by MGeomApply)
+        const Int_t npixels  = fPedestalsOut->GetSize();
+
+        fSumx.  Set(npixels);
+        fSumx2. Set(npixels);
+        fSumAB0.Set(npixels);
+        fSumAB1.Set(npixels);
+
+        if (fUseSpecialPixels)
+        {
+            // Initialize size of MPedestalCam in case of special pixels (not done by MGeomApply)
+            const UShort_t nspecial = fRunHeader->GetNumSpecialPixels();
+            if (nspecial == 0)
+            {
+                *fLog << warn << "WARNING - Number of special pixels is 0." << endl;
+                return kTRUE;
+            }
+
+            fPedestalsOut->InitSize((UInt_t)nspecial);
+        }
+        else
+        {
+            // Initialize the averaged areas and sectors (do not exist for special pixels)
+            const Int_t areas    = fPedestalsOut->GetNumAverageArea();
+            const Int_t sectors  = fPedestalsOut->GetNumAverageSector();
+
+            fAreaSumx.  Set(areas);
+            fAreaSumx2. Set(areas);
+            fAreaSumAB0.Set(areas);
+            fAreaSumAB1.Set(areas);
+            fAreaFilled.Set(areas);
+            fAreaValid .Set(areas);
+
+            fSectorSumx.  Set(sectors);
+            fSectorSumx2. Set(sectors);
+            fSectorSumAB0.Set(sectors);
+            fSectorSumAB1.Set(sectors);
+            fSectorFilled.Set(sectors);
+            fSectorValid .Set(sectors);
+
+            for (Int_t i=0; i<npixels; i++)
+            {
+                const UInt_t aidx   = (*fGeom)[i].GetAidx();
+                const UInt_t sector = (*fGeom)[i].GetSector();
+
+                fAreaValid  [aidx]  ++;
+                fSectorValid[sector]++;
+            }
+        }
+    }
+
+    if (fExtractor)
+    {
+        if (!((MTask*)fExtractor)->ReInit(pList))
+            return kFALSE;
+
+        SetExtractWindow(fExtractor->GetHiGainFirst(), (Int_t)TMath::Nint(fExtractor->GetNumHiGainSamples()));
+    }
+
+    return kTRUE;
 }
 
 Int_t MExtractPedestal::PostProcess()
 {
-  fPedestalsIn = NULL;
-  return fExtractor ? fExtractor->CallPostProcess() : kTRUE;
+    fPedestalsIn = NULL;
+    return fExtractor ? fExtractor->CallPostProcess() : kTRUE;
 }
 
@@ -486 +534 @@
     }
 
+    // find resource for fUseSpecialPixels
+    if (IsEnvDefined(env, prefix, "UseSpecialPixels", print))
+    {
+        SetUseSpecialPixels(GetEnvValue(env, prefix, "UseSpecialPixels", fUseSpecialPixels));
+        rc = kTRUE;
+    }
+
     // find resource for numeventsdump
     if (IsEnvDefined(env, prefix, "NumEventsDump", print))
@@ -577 +632 @@
 void MExtractPedestal::CalcPixResults(const UInt_t nevts, const UInt_t pixid)
 {
-    const Float_t sum  = fSumx[pixid];
-    const Float_t sum2 = fSumx2[pixid];
+    const Double_t sum  = fSumx[pixid];
+    const Double_t sum2 = fSumx2[pixid];
 
     // 1. Calculate the mean of the sums:
-    Float_t ped        = sum/nevts;
+    Double_t ped = sum/nevts;
 
     // 2. Calculate the Variance of the sums:
-    Float_t var = (sum2-sum*sum/nevts)/(nevts-1.);
+    Double_t var = (sum2-sum*sum/nevts)/(nevts-1.);
 
     // 3. Calculate the amplitude of the 150MHz "AB" noise
-    Float_t abOffs = (fSumAB0[pixid] - fSumAB1[pixid]) / nevts;
+    Double_t abOffs = (fSumAB0[pixid] - fSumAB1[pixid]) / nevts;
 
     // 4. Scale the mean, variance and AB-noise to the number of slices:
@@ -595 +650 @@
 
     // 5. Calculate the RMS from the Variance:
-    const Float_t rms = var<0 ? 0 : TMath::Sqrt(var);
+    const Double_t rms = var<0 ? 0 : TMath::Sqrt(var);
 
     (*fPedestalsOut)[pixid].Set(ped, rms, abOffs, nevts);
@@ -612 +667 @@
 void MExtractPedestal::CalcAreaResults(const UInt_t nevts, const UInt_t napix, const UInt_t aidx)
 {
-    const Float_t sum  = fAreaSumx[aidx];
-    const Float_t sum2 = fAreaSumx2[aidx];
+    const Double_t sum  = fAreaSumx[aidx];
+    const Double_t sum2 = fAreaSumx2[aidx];
 
     // 1. Calculate the mean of the sums:
-    Float_t ped        = sum/nevts;
+    Double_t ped = sum/nevts;
 
     // 2. Calculate the Variance of the sums:
-    Float_t var = (sum2-sum*sum/nevts)/(nevts-1.);
+    Double_t var = (sum2/napix-sum*sum/nevts)/(nevts-1.);
 
     // 3. Calculate the amplitude of the 150MHz "AB" noise
-    Float_t abOffs = (fAreaSumAB0[aidx] - fAreaSumAB1[aidx]) / nevts;
+    Double_t abOffs = (fAreaSumAB0[aidx] - fAreaSumAB1[aidx]) / nevts;
 
     // 4. Scale the mean, variance and AB-noise to the number of slices:
@@ -635 +690 @@
 
     // 6. Calculate the RMS from the Variance:
-    const Float_t rms = var<0 ? 0 : TMath::Sqrt(var);
+    const Double_t rms = var<0 ? 0 : TMath::Sqrt(var);
 
     fPedestalsOut->GetAverageArea(aidx).Set(ped, rms, abOffs, nevts);
@@ -652 +707 @@
 void MExtractPedestal::CalcSectorResults(const UInt_t nevts, const UInt_t nspix, const UInt_t sector)
 {
-    const Float_t sum  = fSectorSumx[sector];
-    const Float_t sum2 = fSectorSumx2[sector];
+    const Double_t sum  = fSectorSumx[sector];
+    const Double_t sum2 = fSectorSumx2[sector];
 
     // 1. Calculate the mean of the sums:
-    Float_t ped        = sum/nevts;
+    Double_t ped        = sum/nevts;
 
     // 2. Calculate the Variance of the sums:
-    Float_t var = (sum2-sum*sum/nevts)/(nevts-1.);
+    Double_t var = (sum2/nspix-sum*sum/nevts)/(nevts-1.);
 
     // 3. Calculate the amplitude of the 150MHz "AB" noise
-    Float_t abOffs = (fSectorSumAB0[sector] - fSectorSumAB1[sector]) / nevts;
+    Double_t abOffs = (fSectorSumAB0[sector] - fSectorSumAB1[sector]) / nevts;
 
     // 4. Scale the mean, variance and AB-noise to the number of slices:
@@ -675 +730 @@
 
     // 6. Calculate the RMS from the Variance:
-    const Float_t rms = var<0 ? 0 : TMath::Sqrt(var);
+    const Double_t rms = var<0 ? 0 : TMath::Sqrt(var);
 
     fPedestalsOut->GetAverageSector(sector).Set(ped, rms, abOffs, nevts);
@@ -688 +743 @@
     *fLog << "Intermediate Storage is       " << (fIntermediateStorage?"on":"off") << endl;
     *fLog << "Pedestal Update is            " << (fPedestalUpdate?"on":"off") << endl;
+    *fLog << "Special pixel mode            " << (fUseSpecialPixels?"on":"off") << endl;
     if (fPedestalUpdate)
     {

trunk/MagicSoft/Mars/mpedestal/MExtractPedestal.h

@@ -25 +25 @@
 private:
   static const TString  fgNamePedestalCam;  //! "MPedestalCam"
+  static const TString  fgNameRawEvtData;   //! "MRawEvtData"
   static const UInt_t   fgNumDump;          //!
 
   TString fNamePedestalCamIn;        // Name of the incoming 'MPedestalCam' container
   TString fNamePedestalCamOut;       // Name of the outgoing 'MPedestalCam' container
-  TString fNamePedestalCamInter;     // Name of the intermediate 'MPedestalCam' container  
+  TString fNamePedestalCamInter;     // Name of the intermediate 'MPedestalCam' container
+  TString fNameRawEvtData;           // Name of MRawEvtData
 
   Bool_t  fRandomCalculation;        // Is pedestalextraction by extractor random?
@@ -55 +57 @@
  
   Bool_t  fPedestalUpdate;           // Flag if the pedestal shall be updated after every fNumEventsDump
+  Bool_t  fUseSpecialPixels;         // Flag if the special pixels shall be treated
 
   MArrayD fSumx;                     // sum of values
@@ -103 +106 @@
 
   // names
-  void SetNamePedestalCamIn   (const char *name=fgNamePedestalCam.Data()) { fNamePedestalCamIn    = name; }
-  void SetNamePedestalCamInter(const char *name=fgNamePedestalCam.Data()) { fNamePedestalCamInter = name; }  
-  void SetNamePedestalCamOut  (const char *name=fgNamePedestalCam.Data()) { fNamePedestalCamOut   = name; }
+  void SetNamePedestalCamIn   (const char *name=fgNamePedestalCam) { fNamePedestalCamIn    = name; }
+  void SetNamePedestalCamInter(const char *name=fgNamePedestalCam) { fNamePedestalCamInter = name; }
+  void SetNamePedestalCamOut  (const char *name=fgNamePedestalCam) { fNamePedestalCamOut   = name; }
+  void SetNameRawEvtData      (const char *name=fgNameRawEvtData)  { fNameRawEvtData       = name; }
 
   // pointers
@@ -115 +119 @@
   // flags
   void SetIntermediateStorage (Bool_t b=kTRUE) { fIntermediateStorage = b; }
+  void SetUseSpecialPixels    (Bool_t b=kTRUE) { fUseSpecialPixels    = b; }
   void SetPedestalUpdate      (Bool_t b=kTRUE) { fPedestalUpdate      = b; }
   void SetRandomCalculation   (Bool_t b=kTRUE) { fRandomCalculation   = b; }

trunk/MagicSoft/Mars/mpedestal/MMcPedestalCopy.cc

@@ -79 +79 @@
 Bool_t MMcPedestalCopy::CheckRunType(MParList *pList) const
 {
-    const MRawRunHeader *run = (MRawRunHeader*)pList->FindObject("MRawRunHeader");
+    const MRawRunHeader *run = (MRawRunHeader*)pList->FindObject(AddSerialNumber("MRawRunHeader"));
     if (!run)
     {
-        *fLog << warn << dbginf << "Warning - cannot check file type, MRawRunHeader not found." << endl;
+      *fLog << warn << dbginf << "Warning - cannot check file type, " 
+            << AddSerialNumber("MRawRunHeader") << " not found." << endl;
         return kTRUE;
     }
@@ -138 +139 @@
       }
 
-    MMcRunHeader *mcrun = (MMcRunHeader*)pList->FindObject("MMcRunHeader");
+    MMcRunHeader *mcrun = (MMcRunHeader*)pList->FindObject(AddSerialNumber("MMcRunHeader"));
     if (!mcrun)
-        *fLog << warn << dbginf << "MMcRunHeader not found... assuming camera<0.7" << endl;
+      *fLog << warn << dbginf << AddSerialNumber("MMcRunHeader") 
+            << " not found... assuming camera<0.7" << endl;
 
     const int num = pedcam->GetSize();

trunk/MagicSoft/Mars/mpedestal/MPedCalcPedRun.cc

@@ -276 +276 @@
   
   //
-  // In the other case, produce the MPedestalCam to be use the subsequent (non-pedestal) events
+  // In the other case, produce the MPedestalCam to be use the subsequent
+  // (non-pedestal) events
   //
   *fLog << inf << "Finalizing pedestal calculations..." << flush;
@@ -299 +300 @@
 Int_t MPedCalcPedRun::Calc()
 {
-
   if (fIsNotPedRun && !IsPedBitSet())
     return kTRUE;
@@ -317 +317 @@
   while (pixel.Next())
   {
-      const UInt_t idx    = pixel.GetPixelId();
-      const UInt_t aidx   = (*fGeom)[idx].GetAidx();
-      const UInt_t sector = (*fGeom)[idx].GetSector();      
+      const UInt_t idx = pixel.GetPixelId();
 
       Float_t sum = 0.;
@@ -330 +328 @@
         CalcSums(pixel, sum, ab0, ab1);
 
-      fSumx[idx]           += sum;
+      fSumx[idx] += sum;
+
+      if (fIntermediateStorage)
+        (*fPedestalsInter)[idx].Set(sum,0.,0.,fUsedEvents);
+
+      const Float_t sqrsum = sum*sum;
+
+      fSumx2[idx]  += sqrsum;
+      fSumAB0[idx] += ab0;
+      fSumAB1[idx] += ab1;
+
+      if (fUseSpecialPixels)
+        continue;
+
+      const UInt_t aidx   = (*fGeom)[idx].GetAidx();
+      const UInt_t sector = (*fGeom)[idx].GetSector();      
+
       fAreaSumx[aidx]      += sum;
       fSectorSumx[sector]  += sum;
 
-      if (fIntermediateStorage)
-        (*fPedestalsInter)[idx].Set(sum,0.,0.,fUsedEvents);
-
-      const Float_t sqrsum   = sum*sum;
-      fSumx2[idx]           += sqrsum;
-      fAreaSumx2[aidx]      += sqrsum;
-      fSectorSumx2[sector]  += sqrsum;
-
-      fSumAB0[idx]         += ab0;
-      fSumAB1[idx]         += ab1;
+      fAreaSumx2[aidx]     += sqrsum;
+      fSectorSumx2[sector] += sqrsum;
 
       fAreaSumAB0[aidx]    += ab0;
@@ -422 +428 @@
 }
 
-
 // --------------------------------------------------------------------------
 //
@@ -429 +434 @@
 Int_t MPedCalcPedRun::Finalize()
 {
-
-  if (fUsedEvents == 0)
+    if (fUsedEvents == 0)
+        return kTRUE;
+
+    //
+    // Necessary check for extraction of special pixels
+    // together with data which does not yet have them
+    //
+    if (fSumx.GetSize()==0)
+        return kTRUE;
+
+    MRawEvtPixelIter pixel(fRawEvt);
+    while (pixel.Next())
+        CalcPixResults(fUsedEvents, pixel.GetPixelId());
+
+    if (!fUseSpecialPixels)
+    {
+
+        //
+        // Loop over the (two) area indices to get the averaged pedestal per aidx
+        //
+        for (UInt_t aidx=0; aidx<fAreaValid.GetSize(); aidx++)
+            if (fAreaValid[aidx]>0)
+                CalcAreaResults(fUsedEvents, fAreaValid[aidx], aidx);
+
+        //
+        // Loop over the (six) sector indices to get the averaged pedestal per sector
+        //
+        for (UInt_t sector=0; sector<fSectorValid.GetSize(); sector++)
+            if (fSectorValid[sector]>0)
+                CalcSectorResults(fUsedEvents, fSectorValid[sector], sector);
+    }
+
+    fPedestalsOut->SetTotalEntries(fUsedEvents*fExtractWinSize);
+    fPedestalsOut->SetReadyToSave();
+
     return kTRUE;
-  
-  MRawEvtPixelIter pixel(fRawEvt);
-  while (pixel.Next())
-    CalcPixResults(fUsedEvents, pixel.GetPixelId());
-  
-  //
-  // Loop over the (two) area indices to get the averaged pedestal per aidx
-  //
-  for (UInt_t aidx=0; aidx<fAreaValid.GetSize(); aidx++)
-    if (fAreaValid[aidx]>0)
-      CalcAreaResults(fUsedEvents, fAreaValid[aidx], aidx);
-  
-  //
-  // Loop over the (six) sector indices to get the averaged pedestal per sector
-  //
-  for (UInt_t sector=0; sector<fSectorValid.GetSize(); sector++)
-    if (fSectorValid[sector]>0)
-      CalcSectorResults(fUsedEvents, fSectorValid[sector], sector);
-  
-  fPedestalsOut->SetTotalEntries(fUsedEvents*fExtractWinSize);
-  fPedestalsOut->SetReadyToSave();
-  
-  return kTRUE;
 }
 
 Int_t MPedCalcPedRun::PostProcess()
 {
-  if (!Finalize())
-    return kFALSE;
-  
-  return MExtractPedestal::PostProcess();
+    if (!fRawEvt)
+        return kTRUE;
+
+    if (!Finalize())
+        return kFALSE;
+
+    return MExtractPedestal::PostProcess();
 }
 
@@ -485 +503 @@
 void MPedCalcPedRun::Print(Option_t *o) const
 {
-
     MExtractPedestal::Print(o);
 

trunk/MagicSoft/Mars/mpedestal/MPedestalCam.cc

@@ -276 +276 @@
 }
 
+void MPedestalCam::PrintArr(const TCollection &list) const
+{
+    Int_t id = 0;
+
+    TIter Next(&list);
+    MPedestalPix *pix = 0;
+
+    while ((pix=(MPedestalPix*)Next()))
+        *fLog << Form("Nr.: %4i  Pedestal: %5.2f    RMS: %5.2f    ABOffset: %5.2f ",
+                      id++, pix->GetPedestal(), pix->GetPedestalRms(), pix->GetPedestalABoffset()) << endl;
+}
+
 void MPedestalCam::Print(Option_t *o) const
 {
     *fLog << all << GetDescriptor() << ":" << endl;
-    int id = 0;
-
-    TIter Next(fArray);
-    MPedestalPix *pix;
-    while ((pix=(MPedestalPix*)Next()))
-    {
-        id++;
-
-        if (!pix->IsValid())
-            continue;
-
-        *fLog << id-1 << ": ";
-        *fLog << pix->GetPedestal() << " " << pix->GetPedestalRms() << endl;
-    }
+    *fLog << "Pixels:" << endl;    
+    *fLog << endl;
+
+    PrintArr(*fArray);
+
+    *fLog << endl;
+    *fLog << "Event-by-event averaged areas:" << endl;    
+    *fLog << endl;
+
+    PrintArr(*fAverageAreas);
+
+    *fLog << endl;
+    *fLog << "Event-by-event averaged sectors:" << endl;
+    *fLog << endl;
+
+    PrintArr(*fAverageSectors);
 }
 

trunk/MagicSoft/Mars/mpedestal/MPedestalCam.h

@@ -24 +24 @@
 
   UInt_t fTotalEntries;  // Total number of times, the Process was executed (to estimate the error of pedestal)
+
+  void PrintArr(const TCollection &list) const;
 
 public: