Index: trunk/MagicSoft/Mars/mpedestal/MPedCalcFromLoGain.cc =================================================================== --- trunk/MagicSoft/Mars/mpedestal/MPedCalcFromLoGain.cc (revision 4565) +++ trunk/MagicSoft/Mars/mpedestal/MPedCalcFromLoGain.cc (revision 4601) @@ -114,4 +114,14 @@ // fLoGainWindowSize = fgLoGainWindowSize = 0 // +// Variables: +// fgHiGainFirst; First FADC slice Hi-Gain (currently set to: 3) +// fgHiGainLast: Last FADC slice Hi-Gain (currently set to: 14) +// fgLoGainFirst: First FADC slice Lo-Gain (currently set to: 3) +// fgLoGainLast: Last FADC slice Lo-Gain (currently set to: 14) +// fgHiGainWindowSize: The extraction window Hi-Gain +// fgLoGainWindowSize: The extraction window Lo-Gain +// fgMaxHiGainVar: The maximum difference between the highest and lowest slice +// in the high gain window allowed in order to use low gain +// // Input Containers: // MRawEvtData @@ -201,4 +211,16 @@ fRunHeader = NULL; fPedestals = NULL; + + // If the size is yet set, set the size + if (fSumx.GetSize()>0) + { + // Reset contents of arrays. + fSumx.Reset(); + fSumx2.Reset(); + fSumAB0.Reset(); + fSumAB1.Reset(); + fNumEventsUsed.Reset(); + fTotalCounter.Reset(); + } } @@ -297,40 +319,34 @@ // - MPedestalCam // -Int_t MPedCalcFromLoGain::PreProcess( MParList *pList ) -{ - Clear(); - - fRawEvt = (MRawEvtData*)pList->FindObject("MRawEvtData"); - if (!fRawEvt) - { - *fLog << err << "MRawEvtData not found... aborting." << endl; - return kFALSE; - } - - fRunHeader = (MRawRunHeader*)pList->FindObject(AddSerialNumber("MRawRunHeader")); - if (!fRunHeader) - { - *fLog << err << AddSerialNumber("MRawRunHeader") << " not found... aborting." << endl; - return kFALSE; - } - - fGeom = (MGeomCam*)pList->FindObject("MGeomCam"); - if (!fGeom) - { - *fLog << err << "MGeomCam not found... aborting." << endl; - return kFALSE; - } - - *fLog << " searching for the container " << fPedContainerName << endl; - - fPedestals = (MPedestalCam*)pList->FindCreateObj( AddSerialNumber("MPedestalCam"),fPedContainerName); - if (!fPedestals) - { - *fLog << err << fPedContainerName << " can not be created" << endl; - return kFALSE; - } - - - return kTRUE; +Int_t MPedCalcFromLoGain::PreProcess(MParList *pList) +{ + Clear(); + + fRawEvt = (MRawEvtData*)pList->FindObject("MRawEvtData"); + if (!fRawEvt) + { + *fLog << err << "MRawEvtData not found... aborting." << endl; + return kFALSE; + } + + fRunHeader = (MRawRunHeader*)pList->FindObject(AddSerialNumber("MRawRunHeader")); + if (!fRunHeader) + { + *fLog << err << AddSerialNumber("MRawRunHeader") << " not found... aborting." << endl; + return kFALSE; + } + + fGeom = (MGeomCam*)pList->FindObject("MGeomCam"); + if (!fGeom) + { + *fLog << err << "MGeomCam not found... aborting." << endl; + return kFALSE; + } + + fPedestals = (MPedestalCam*)pList->FindCreateObj("MPedestalCam", AddSerialNumber(fPedContainerName)); + if (!fPedestals) + return kFALSE; + + return kTRUE; } @@ -353,5 +369,4 @@ Bool_t MPedCalcFromLoGain::ReInit(MParList *pList) { - Int_t lastdesired = (Int_t)fLoGainLast; Int_t lastavailable = (Int_t)fRunHeader->GetNumSamplesLoGain()-1; @@ -413,8 +428,9 @@ - Int_t npixels = fPedestals->GetSize(); - + // If the size is not yet set, set the size if (fSumx.GetSize()==0) - { + { + const Int_t npixels = fPedestals->GetSize(); + fSumx. Set(npixels); fSumx2.Set(npixels); @@ -423,5 +439,6 @@ fNumEventsUsed.Set(npixels); fTotalCounter.Set(npixels); - + + // Reset contents of arrays. fSumx.Reset(); fSumx2.Reset(); @@ -430,14 +447,13 @@ fNumEventsUsed.Reset(); fTotalCounter.Reset(); - } - - if (fWindowSizeHiGain == 0 && fWindowSizeLoGain == 0) - { - *fLog << err << GetDescriptor() + } + + if (fWindowSizeHiGain==0 && fWindowSizeLoGain==0) + { + *fLog << err << GetDescriptor() << ": Number of extracted Slices is 0, abort ... " << endl; return kFALSE; - } - - + } + *fLog << endl; *fLog << inf << GetDescriptor() << ": Taking " << (int)fWindowSizeHiGain @@ -449,4 +465,29 @@ } +void MPedCalcFromLoGain::Calc(ULong_t n, UInt_t idx) +{ + const ULong_t nsamplestot = n*fWindowSizeLoGain; + + const Float_t sum = fSumx.At(idx); + const Float_t sum2 = fSumx2.At(idx); + const Float_t ped = sum/nsamplestot; + + // 1. Calculate the Variance of the sums: + Float_t var = (sum2-sum*sum/n)/(n-1.); + + // 2. Scale the variance to the number of slices: + var /= (Float_t)(fWindowSizeLoGain); + + // 3. Calculate the RMS from the Variance: + const Float_t rms = var<0 ? 0 : TMath::Sqrt(var); + + // 4. Calculate the amplitude of the 150MHz "AB" noise + const Float_t abOffs = (fSumAB0[idx] - fSumAB1[idx]) / nsamplestot; + + (*fPedestals)[idx].Set(ped, rms, abOffs, n); + + fTotalCounter[idx]++; +} + // -------------------------------------------------------------------------- // @@ -457,97 +498,78 @@ Int_t MPedCalcFromLoGain::Process() { - MRawEvtPixelIter pixel(fRawEvt); - - while (pixel.Next()) { + MRawEvtPixelIter pixel(fRawEvt); + + while (pixel.Next()) + { + const UInt_t idx = pixel.GetPixelId(); + + Byte_t *ptr = pixel.GetHiGainSamples() + fHiGainFirst; + Byte_t *end = ptr + fWindowSizeHiGain; + + UInt_t sum = 0; + UInt_t sqr = 0; + + UInt_t max = 0; + UInt_t min = 255; - const UInt_t idx = pixel.GetPixelId(); - - Byte_t *ptr = pixel.GetHiGainSamples() + fHiGainFirst; - Byte_t *end = ptr + fWindowSizeHiGain; - - UInt_t sum = 0; - UInt_t sqr = 0; - - UInt_t max = 0; - UInt_t min = 255; - - // Find the maximum and minimum signal per slice in the high gain window - do { - if (*ptr > max) { - max = *ptr; - } - if (*ptr < min) { - min = *ptr; - } - } while (++ptr != end); - - // If the maximum in the high gain window is smaller than - if ((max-min < fMaxHiGainVar) && (max < 255)) { - - Byte_t *firstSlice = pixel.GetLoGainSamples() + fLoGainFirst; - Byte_t *lastSlice = firstSlice + fWindowSizeLoGain; - - Byte_t *slice = firstSlice; - do { - sum += *slice; - sqr += *slice * *slice; - } while (++slice != lastSlice); - - const Float_t msum = (Float_t)sum; - const Float_t sqrsum = msum*msum; - - fSumx[idx] += msum; - fSumx2[idx] += sqrsum; - fNumEventsUsed[idx]++; - - // Calculate the amplitude of the 150MHz "AB" noise - - Int_t abFlag = (fRunHeader->GetNumSamplesHiGain() - + fLoGainFirst - + pixel.HasABFlag()) & 0x1; - for (Int_t islice=0; islice= fNumEventsDump)) { - - const ULong_t n = fNumEventsDump; - - const ULong_t nsamplestot = n*fWindowSizeLoGain; - - const Float_t sum = fSumx.At(idx); - const Float_t sum2 = fSumx2.At(idx); - const Float_t ped = sum/(nsamplestot); - - // 1. Calculate the Variance of the sums: - Float_t var = (sum2-sum*sum/n)/(n-1.); - - // 2. Scale the variance to the number of slices: - var /= (Float_t)(fWindowSizeLoGain); - - // 3. Calculate the RMS from the Variance: - Float_t rms = var < 0 ? 0. : TMath::Sqrt(var); - - // 4. Calculate the amplitude of the 150MHz "AB" noise - Float_t abOffs = (fSumAB0[idx] - fSumAB1[idx]) / nsamplestot; - - (*fPedestals)[idx].Set(ped, rms, abOffs, n); - - fTotalCounter[idx]++; - fNumEventsUsed[idx]=0; - fSumx[idx]=0; - fSumx2[idx]=0; - fSumAB0[idx]=0; - fSumAB1[idx]=0; - } - } - } - - fPedestals->SetReadyToSave(); - - return kTRUE; + // Find the maximum and minimum signal per slice in the high gain window + do { + if (*ptr > max) { + max = *ptr; + } + if (*ptr < min) { + min = *ptr; + } + } while (++ptr != end); + + // If the maximum in the high gain window is smaller than + if (max-min>=fMaxHiGainVar || max>=255) + continue; + + Byte_t *firstSlice = pixel.GetLoGainSamples() + fLoGainFirst; + Byte_t *lastSlice = firstSlice + fWindowSizeLoGain; + + Byte_t *slice = firstSlice; + do { + sum += *slice; + sqr += *slice * *slice; + } while (++slice != lastSlice); + + const Float_t msum = (Float_t)sum; + const Float_t sqrsum = msum*msum; + + fSumx[idx] += msum; + fSumx2[idx] += sqrsum; + fNumEventsUsed[idx]++; + + // Calculate the amplitude of the 150MHz "AB" noise + + Int_t abFlag = (fRunHeader->GetNumSamplesHiGain() + + fLoGainFirst + + pixel.HasABFlag()) & 0x1; + for (Int_t islice=0; isliceSetReadyToSave(); + + return kTRUE; } @@ -558,53 +580,69 @@ Int_t MPedCalcFromLoGain::PostProcess() { - - // Compute pedestals and rms from the whole run - - if (!fPedestalUpdate) { - - MRawEvtPixelIter pixel(fRawEvt); - - while (pixel.Next()) { - - const Int_t idx = pixel.GetPixelId(); - - const ULong_t n = fNumEventsUsed[idx]; - - if (n < 2) - continue; - - const ULong_t nsamplestot = n*fWindowSizeLoGain; - - const Float_t sum = fSumx.At(idx); - const Float_t sum2 = fSumx2.At(idx); - const Float_t ped = sum/(nsamplestot); - - // 1. Calculate the Variance of the sums: - Float_t var = (sum2-sum*sum/n)/(n-1.); - - // 2. Scale the variance to the number of slices: - var /= (Float_t)(fWindowSizeLoGain); - - // 3. Calculate the RMS from the Variance: - Float_t rms = var < 0 ? 0. : TMath::Sqrt(var); - - // 4. Calculate the amplitude of the 150MHz "AB" noise - Float_t abOffs = (fSumAB0[idx] - fSumAB1[idx]) / nsamplestot; - - (*fPedestals)[idx].Set(ped, rms, abOffs, n); - - fTotalCounter[idx]++; + // Compute pedestals and rms from the whole run + if (fPedestalUpdate) + return kTRUE; + + *fLog << flush << inf << "Calculating Pedestals..." << flush; + + const Int_t npix = fGeom->GetNumPixels(); + for (Int_t idx=0; idx1) + Calc(n, idx); } fPedestals->SetReadyToSave(); - } - - fSumx.Reset(); - fSumx2.Reset(); - fSumAB0.Reset(); - fSumAB1.Reset(); - fNumEventsUsed.Reset(); - fTotalCounter.Reset(); - - return kTRUE; -} + return kTRUE; +} + +Int_t MPedCalcFromLoGain::ReadEnv(const TEnv &env, TString prefix, Bool_t print) +{ + if (MExtractor::ReadEnv(env, prefix, print)==kERROR) + return kERROR; + + Bool_t rc=kFALSE; + + Byte_t hw = fWindowSizeHiGain; + Byte_t lw = fWindowSizeLoGain; + if (IsEnvDefined(env, prefix, "WindowSizeHiGain", print)) + { + hw = GetEnvValue(env, prefix, "WindowSizeHiGain", hw); + rc = kTRUE; + } + if (IsEnvDefined(env, prefix, "WindowSizeLoGain", print)) + { + lw = GetEnvValue(env, prefix, "WindowSizeLoGain", lw); + rc = kTRUE; + } + + if (rc) + SetWindowSize(hw, lw); + + Int_t num = fNumEventsDump; + if (IsEnvDefined(env, prefix, "NumEventsDump", print)) + { + num = GetEnvValue(env, prefix, "NumEventsDump", num); + rc = kTRUE; + } + SetDumpEvents(num); + + Byte_t max = fMaxHiGainVar; + if (IsEnvDefined(env, prefix, "MaxHiGainVar", print)) + { + max = GetEnvValue(env, prefix, "MaxHiGainVar", max); + rc = kTRUE; + } + SetMaxHiGainVar(max); + + Bool_t upd = fPedestalUpdate; + if (IsEnvDefined(env, prefix, "PedestalUpdate", print)) + { + upd = GetEnvValue(env, prefix, "PedestalUpdate", upd); + rc = kTRUE; + } + SetPedestalUpdate(upd); + + return rc; +} Index: trunk/MagicSoft/Mars/mpedestal/MPedCalcFromLoGain.h =================================================================== --- trunk/MagicSoft/Mars/mpedestal/MPedCalcFromLoGain.h (revision 4565) +++ trunk/MagicSoft/Mars/mpedestal/MPedCalcFromLoGain.h (revision 4601) @@ -19,52 +19,59 @@ class MPedCalcFromLoGain : public MExtractor { - static const Byte_t fgHiGainFirst; // First FADC slice Hi-Gain (currently set to: 3) - static const Byte_t fgHiGainLast; // Last FADC slice Hi-Gain (currently set to: 14) - static const Byte_t fgLoGainFirst; // First FADC slice Lo-Gain (currently set to: 3) - static const Byte_t fgLoGainLast; // Last FADC slice Lo-Gain (currently set to: 14) - static const Byte_t fgHiGainWindowSize; // The extraction window Hi-Gain - static const Byte_t fgLoGainWindowSize; // The extraction window Lo-Gain - static const Byte_t fgMaxHiGainVar; // The maximum difference between the highest and lowest slice - // in the high gain window allowed in order to use low gain + static const Byte_t fgHiGainFirst; // First FADC slice Hi-Gain (currently set to: 3) + static const Byte_t fgHiGainLast; // Last FADC slice Hi-Gain (currently set to: 14) + static const Byte_t fgLoGainFirst; // First FADC slice Lo-Gain (currently set to: 3) + static const Byte_t fgLoGainLast; // Last FADC slice Lo-Gain (currently set to: 14) + static const Byte_t fgHiGainWindowSize; // The extraction window Hi-Gain + static const Byte_t fgLoGainWindowSize; // The extraction window Lo-Gain + static const Byte_t fgMaxHiGainVar; // The maximum difference between the highest and lowest slice - Int_t fNumEventsDump; // Number of event after which MPedestalCam gets updated + Int_t fNumEventsDump; // Number of event after which MPedestalCam gets updated - Byte_t fMaxHiGainVar; - Byte_t fWindowSizeHiGain; // Number of Hi Gain slices in window - Byte_t fWindowSizeLoGain; // Number of Lo Gain slices in window + Byte_t fMaxHiGainVar; + Byte_t fWindowSizeHiGain; // Number of Hi Gain slices in window + Byte_t fWindowSizeLoGain; // Number of Lo Gain slices in window - Bool_t fPedestalUpdate; + Bool_t fPedestalUpdate; - MGeomCam *fGeom; // Camera geometry - TString fPedContainerName; // name of the 'MPedestalCam' container - - TArrayI fNumEventsUsed; // Number of events used for pedestal calc for each pixel - TArrayI fTotalCounter; // Counter for dumping values to Pedestal Container - TArrayD fSumx; // sum of values - TArrayD fSumx2; // sum of squared values - TArrayD fSumAB0; // sum of ABFlag=0 slices - TArrayD fSumAB1; // sum of ABFlag=1 slices + MGeomCam *fGeom; // Camera geometry + TString fPedContainerName; // name of the 'MPedestalCam' container - - Int_t PreProcess (MParList *pList); - Bool_t ReInit (MParList *pList); - Int_t Process (); - Int_t PostProcess(); - + TArrayI fNumEventsUsed; // Number of events used for pedestal calc for each pixel + TArrayI fTotalCounter; // Counter for dumping values to Pedestal Container + TArrayD fSumx; // sum of values + TArrayD fSumx2; // sum of squared values + TArrayD fSumAB0; // sum of ABFlag=0 slices + TArrayD fSumAB1; // sum of ABFlag=1 slices + + // MParContainer + Int_t PreProcess (MParList *pList); + Bool_t ReInit (MParList *pList); + Int_t Process (); + Int_t PostProcess(); + Int_t ReadEnv(const TEnv &env, TString prefix, Bool_t print); + + // Calculation + void Calc(ULong_t n, UInt_t idx); + public: - MPedCalcFromLoGain(const char *name=NULL, const char *title=NULL); - - void Clear(const Option_t *o=""); - void SetRange(Byte_t hifirst=0, Byte_t hilast=0, Byte_t lofirst=0, Byte_t lolast=0); - void SetWindowSize(Byte_t windowh=0, Byte_t windowl=0); - void SetMaxHiGainVar(Byte_t maxvar=0); - void SetDumpEvents(UInt_t dumpevents = 0) {fNumEventsDump = dumpevents;} - void SetPedestalUpdate(Bool_t pedupdate) {fPedestalUpdate = pedupdate;} - - void SetPedContainerName(const char *name) { fPedContainerName = name; } + MPedCalcFromLoGain(const char *name=NULL, const char *title=NULL); - TArrayI *GetNumEventsUsed() {return &fNumEventsUsed;}; - - ClassDef(MPedCalcFromLoGain, 0) // Task to calculate pedestals from pedestal runs raw data + // TObject + void Clear(const Option_t *o=""); + + // Setter + void SetRange(Byte_t hifirst=0, Byte_t hilast=0, Byte_t lofirst=0, Byte_t lolast=0); + void SetWindowSize(Byte_t windowh=0, Byte_t windowl=0); + void SetMaxHiGainVar(Byte_t maxvar=0); + void SetDumpEvents(UInt_t dumpevents = 0) {fNumEventsDump = dumpevents;} + void SetPedestalUpdate(Bool_t pedupdate) {fPedestalUpdate = pedupdate;} + + void SetPedContainerName(const char *name) { fPedContainerName = name; } + + // Getter + TArrayI *GetNumEventsUsed() { return &fNumEventsUsed; } + + ClassDef(MPedCalcFromLoGain, 0) // Task to calculate pedestals from pedestal runs raw data }; Index: trunk/MagicSoft/Mars/mpedestal/MPedCalcPedRun.cc =================================================================== --- trunk/MagicSoft/Mars/mpedestal/MPedCalcPedRun.cc (revision 4565) +++ trunk/MagicSoft/Mars/mpedestal/MPedCalcPedRun.cc (revision 4601) @@ -443,5 +443,4 @@ while (pixel.Next()) { - const UInt_t idx = pixel.GetPixelId(); const UInt_t aidx = (*fGeom)[idx].GetAidx(); @@ -620,2 +619,30 @@ return kTRUE; } + +Int_t MPedCalcPedRun::ReadEnv(const TEnv &env, TString prefix, Bool_t print) +{ + if (MExtractor::ReadEnv(env, prefix, print)==kERROR) + return kERROR; + + Byte_t hw = fWindowSizeHiGain; + Byte_t lw = fWindowSizeLoGain; + + Bool_t rc = kFALSE; + + if (IsEnvDefined(env, prefix, "WindowSizeHiGain", print)) + { + hw = GetEnvValue(env, prefix, "WindowSizeHiGain", hw); + rc=kTRUE; + } + + if (IsEnvDefined(env, prefix, "WindowSizeLoGain", print)) + { + lw = GetEnvValue(env, prefix, "WindowSizeLoGain", lw); + rc=kTRUE; + } + + if (rc) + SetWindowSize(hw, lw); + + return rc; +} Index: trunk/MagicSoft/Mars/mpedestal/MPedCalcPedRun.h =================================================================== --- trunk/MagicSoft/Mars/mpedestal/MPedCalcPedRun.h (revision 4565) +++ trunk/MagicSoft/Mars/mpedestal/MPedCalcPedRun.h (revision 4601) @@ -44,4 +44,5 @@ Int_t Process (); Int_t PostProcess(); + Int_t ReadEnv(const TEnv &env, TString prefix, Bool_t print); public: