Changeset 2629 for trunk/MagicSoft

Timestamp:

12/10/03 08:24:42 (21 years ago)

Author:

gaug

Message:

*** empty log message ***

Location:

trunk/MagicSoft/Mars

Files:

• macros/calibration.C (modified) (1 diff)
• manalysis/MCalibrationPix.cc (modified) (1 diff)
• manalysis/MPedCalcPedRun.cc (modified) (8 diffs)
• manalysis/MPedCalcPedRun.h (modified) (3 diffs)
• mhist/MHCalibrationPixel.cc (modified) (1 diff)
• mhist/MHCalibrationPixel.h (modified) (1 diff)

trunk/MagicSoft/Mars/macros/calibration.C

@@ -218 +218 @@
     disp15.SetYTitle("Conversion Factor [Ph/FADC count]");
 
+    disp1.SetPrettyPalette();
+    disp3.SetPrettyPalette();
+    disp5.SetPrettyPalette();
+    disp6.SetPrettyPalette();
+    disp7.SetPrettyPalette();
+    disp8.SetPrettyPalette();
+    disp9.SetPrettyPalette();
+    disp10.SetPrettyPalette();
+    disp11.SetPrettyPalette();
+    disp12.SetPrettyPalette();
+    disp13.SetPrettyPalette();
+    disp14.SetPrettyPalette();
+    disp15.SetPrettyPalette();
+
+
+
     MStatusDisplay *d2 = new MStatusDisplay;
  

trunk/MagicSoft/Mars/manalysis/MCalibrationPix.cc

@@ -107 +107 @@
 
 
-  if (fHiGainSaturation)
+  if (fHist->UseLoGain())
     {
+      SetHiGainSaturation();
+
       if(!fHist->FitChargeLoGain())
         {

trunk/MagicSoft/Mars/manalysis/MPedCalcPedRun.cc

@@ -50 +50 @@
 #include "MPedestalCam.h"
 
+#include "MGeomCamMagic.h"
+
 ClassImp(MPedCalcPedRun);
 
@@ -61 +63 @@
     AddToBranchList("fHiGainPixId");
     AddToBranchList("fHiGainFadcSamples");
+    
+    fCounter=0;
 }
 
@@ -75 +79 @@
     if (!fPedestals)
         return kFALSE;
+
+    MGeomCamMagic magiccam;
+
+    fSumx.Set(magiccam.GetNumPixels());
+    fSumx2.Set(magiccam.GetNumPixels());
+    
+    for(UInt_t i=0;i<magiccam.GetNumPixels();i++)
+    {
+        fSumx.AddAt(0,i);
+        fSumx2.AddAt(0,i);
+    }
 
     return kTRUE;
@@ -110 +125 @@
     while (pixel.Next())
     {
-             Byte_t *ptr = pixel.GetHiGainSamples();
-        const Byte_t *end = ptr + fRawEvt->GetNumHiGainSamples();
+        Byte_t shift=(fNumHiGainSamples/2*2==fNumHiGainSamples) ? 0:1;
+              Byte_t *ptr = pixel.GetHiGainSamples();
+        const Byte_t *end = ptr + fRawEvt->GetNumHiGainSamples()-shift;
 
         const Float_t higainped = CalcHiGainMean(ptr, end);
@@ -118 +134 @@
         const UInt_t pixid = pixel.GetPixelId();
         MPedestalPix &pix = (*fPedestals)[pixid];
-
+        
+        // cumulate the sum of pedestals and of pedestal squares
+        fSumx.AddAt(higainped+fSumx.At(pixid),pixid);
+        fSumx2.AddAt(GetSumx2(ptr, end)+fSumx2.At(pixid),pixid);
+
+        // set the value of the pedestal and rms computed from the processed event
         pix.Set(higainped, higainrms);
-
-    }
+    }
+
+
+    fCounter++;
 
     fPedestals->SetReadyToSave();
@@ -128 +151 @@
 }
 
+Int_t MPedCalcPedRun::PostProcess()
+{
+// Compute pedestals and rms from the whole run
+
+    MRawEvtPixelIter pixel(fRawEvt);
+
+    while (pixel.Next())
+    {
+        const UInt_t pixid = pixel.GetPixelId();
+        MPedestalPix &pix = (*fPedestals)[pixid];
+
+        const Int_t N = fCounter;
+        const Float_t sum = fSumx.At(pixid);
+        const Float_t sum2 = fSumx2.At(pixid);
+        const Float_t higainped = sum/N;
+        const Float_t higainrms = sqrt(1./(N-1.)*(sum2-sum*sum/N));
+        pix.Set(higainped,higainrms);
+
+    }
+    return kTRUE;
+
+}
+
+
 Float_t MPedCalcPedRun::CalcHiGainMean(Byte_t *ptr, const Byte_t *end) const
 {
     Int_t sum=0;
+    Byte_t EvenNumSamples=(fNumHiGainSamples/2*2==fNumHiGainSamples) ? fNumHiGainSamples:fNumHiGainSamples-1;
 
     do sum += *ptr;
     while (++ptr != end);
-
-    return (Float_t)sum/fNumHiGainSamples;
+   
+    return (Float_t)sum/EvenNumSamples;
+}
+
+
+
+Float_t MPedCalcPedRun::GetSumx2(Byte_t *ptr, const Byte_t *end) const
+{
+    Float_t square = 0;
+
+    // Take an even number of time slices to avoid biases due to A/B effect
+    Byte_t EvenNumSamples=(fNumHiGainSamples/2*2==fNumHiGainSamples) ? fNumHiGainSamples:fNumHiGainSamples-1;
+
+    do
+    {
+        const Float_t val = (Float_t)(*ptr);
+
+        square += val*val;
+    } while (++ptr != end);
+
+    return square/EvenNumSamples;
 }
 
@@ -142 +209 @@
 {
     Float_t rms = 0;
+    Byte_t EvenNumSamples=(fNumHiGainSamples/2*2==fNumHiGainSamples) ? fNumHiGainSamples:fNumHiGainSamples-1;
 
     do
@@ -150 +218 @@
     } while (++ptr != end);
 
-    return sqrt(rms/fNumHiGainSamples);
-}
+    return sqrt(rms/EvenNumSamples);
+}
+
+
+
 /*
 Float_t MPedCalcPedRun::CalcHiGainMeanErr(Float_t higainrms) const

trunk/MagicSoft/Mars/manalysis/MPedCalcPedRun.h

@@ -14 +14 @@
 #endif
 
+#include <TArrayF.h>
+
 class MRawRunHeader;
 class MRawEvtData;
@@ -22 +24 @@
     Byte_t fNumHiGainSamples;
 
+    Int_t fCounter; // keep number of processed events
+
     MRawRunHeader *fRunheader; // raw event run header
     MRawEvtData  *fRawEvt;     // raw event data (time slices)
     MPedestalCam *fPedestals;  // Pedestals of all pixels in the camera
 
+    TArrayF fSumx;   // sum of values
+    TArrayF fSumx2;  // sum of squared values
+
     Float_t CalcHiGainMean(Byte_t *ptr, const Byte_t *end) const;
     Float_t CalcHiGainRms(Byte_t *ptr, const Byte_t *end, Float_t higainped) const;
-    //Float_t CalcHiGainMeanErr(Float_t higainrms) const;
+    Float_t GetSumx2(Byte_t* ptr, const Byte_t* end) const; 
+   //Float_t CalcHiGainMeanErr(Float_t higainrms) const;
     //Float_t CalcHiGainRmsErr(Float_t higainrms) const;
 
@@ -35 +43 @@
     Int_t PreProcess(MParList *pList);
     Int_t Process();
+    Int_t PostProcess();
 
 public:

trunk/MagicSoft/Mars/mhist/MHCalibrationPixel.cc

@@ -264 +264 @@
 }
 
+
+Bool_t MHCalibrationPixel::UseLoGain()
+{
+
+  if (fHChargeHiGain->GetEntries() > fHChargeLoGain->GetEntries())
+    {
+      fUseLoGain = kFALSE;
+      return kFALSE;
+    }
+  else
+    {
+      fUseLoGain = kTRUE;
+      return kTRUE;
+    }
+}
 
 

trunk/MagicSoft/Mars/mhist/MHCalibrationPixel.h

@@ -93 +93 @@
   Bool_t FillChargevsNHiGain(Float_t q, Int_t n) { return fHChargevsNHiGain->Fill(n,q) > -1; }
 
-  void SetUseLoGain()                      { fUseLoGain = kTRUE; }
+  void   SetUseLoGain()                          { fUseLoGain = kTRUE; }
+  Bool_t UseLoGain();
 
   const TH1F *GetHCharge()                 { return fHChargeHiGain;    }