Index: trunk/MagicSoft/Mars/mpedestal/MPedCalcFromLoGain.cc =================================================================== --- trunk/MagicSoft/Mars/mpedestal/MPedCalcFromLoGain.cc (revision 4601) +++ trunk/MagicSoft/Mars/mpedestal/MPedCalcFromLoGain.cc (revision 4609) @@ -185,15 +185,17 @@ fGeom(NULL), fPedContainerName("MPedestalCam") { - fName = name ? name : "MPedCalcFromLoGain"; - fTitle = title ? title : "Task to calculate pedestals from pedestal runs raw data"; - - AddToBranchList("fHiGainPixId"); - AddToBranchList("fHiGainFadcSamples"); - - SetRange(fgHiGainFirst, fgHiGainLast, fgLoGainFirst, fgLoGainLast); - - SetMaxHiGainVar(fgMaxHiGainVar); - SetPedestalUpdate(kTRUE); - Clear(); + fName = name ? name : "MPedCalcFromLoGain"; + fTitle = title ? title : "Task to calculate pedestals from pedestal runs raw data"; + + AddToBranchList("fHiGainPixId"); + AddToBranchList("fLoGainPixId"); + AddToBranchList("fHiGainFadcSamples"); + AddToBranchList("fLoGainFadcSamples"); + + SetRange(); + SetMaxHiGainVar(); + SetPedestalUpdate(kTRUE); + + Clear(); } @@ -208,19 +210,19 @@ void MPedCalcFromLoGain::Clear(const Option_t *o) { - fRawEvt = NULL; - 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(); - } + fRawEvt = NULL; + 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(); + } } @@ -235,10 +237,10 @@ void MPedCalcFromLoGain::SetRange(Byte_t hifirst, Byte_t hilast, Byte_t lofirst, Byte_t lolast) { - MExtractor::SetRange(hifirst,hilast,lofirst,lolast); - - // - // Redo the checks if the window is still inside the ranges - // - SetWindowSize(fWindowSizeHiGain,fWindowSizeLoGain); + MExtractor::SetRange(hifirst, hilast, lofirst, lolast); + + // + // Redo the checks if the window is still inside the ranges + // + SetWindowSize(fWindowSizeHiGain,fWindowSizeLoGain); } @@ -247,5 +249,5 @@ void MPedCalcFromLoGain::SetMaxHiGainVar(Byte_t maxvar) { - fMaxHiGainVar = maxvar; + fMaxHiGainVar = maxvar; } @@ -265,43 +267,50 @@ void MPedCalcFromLoGain::SetWindowSize(Byte_t windowh, Byte_t windowl) { - fWindowSizeHiGain = windowh & ~1; - fWindowSizeLoGain = windowl & ~1; - - if (fWindowSizeHiGain != windowh) - *fLog << warn << GetDescriptor() << ": Hi Gain window size has to be even, set to: " - << int(fWindowSizeHiGain) << " samples " << endl; + fWindowSizeHiGain = windowh & ~1; + fWindowSizeLoGain = windowl & ~1; - if (fWindowSizeLoGain != windowl) - *fLog << warn << GetDescriptor() << ": Lo Gain window size has to be even, set to: " - << int(fWindowSizeLoGain) << " samples " << endl; + if (fWindowSizeHiGain != windowh) + { + *fLog << warn; + *fLog << GetDescriptor() << ": HiGain window has to be even, set to: "; + *fLog << int(fWindowSizeHiGain) << " samples " << endl; + } - const Byte_t availhirange = (fHiGainLast-fHiGainFirst+1) & ~1; - const Byte_t availlorange = (fLoGainLast-fLoGainFirst+1) & ~1; - - if (fWindowSizeHiGain > availhirange) - { - *fLog << warn << GetDescriptor() - << Form(": Hi Gain window size: %2i is bigger than available range: [%2i,%2i]", - (int)fWindowSizeHiGain, (int)fHiGainFirst, (int)fHiGainLast) << endl; - *fLog << warn << GetDescriptor() - << ": Will set window size to: " << (int)availhirange << endl; - fWindowSizeHiGain = availhirange; - } + if (fWindowSizeLoGain != windowl) + { + *fLog << warn; + *fLog << GetDescriptor() << ": Lo Gain window has to be even, set to: "; + *fLog << int(fWindowSizeLoGain) << " samples " << endl; + } + + const Byte_t availhirange = (fHiGainLast-fHiGainFirst+1) & ~1; + const Byte_t availlorange = (fLoGainLast-fLoGainFirst+1) & ~1; - if (fWindowSizeLoGain > availlorange) - { - *fLog << warn << GetDescriptor() - << Form(": Lo Gain window size: %2i is bigger than available range: [%2i,%2i]", - (int)fWindowSizeLoGain, (int)fLoGainFirst, (int)fLoGainLast) << endl; - *fLog << warn << GetDescriptor() - << ": Will set window size to: " << (int)availlorange << endl; - fWindowSizeLoGain = availlorange; - } - - fNumHiGainSamples = (Float_t)fWindowSizeHiGain; - fNumLoGainSamples = (Float_t)fWindowSizeLoGain; - - fSqrtHiGainSamples = TMath::Sqrt(fNumHiGainSamples); - fSqrtLoGainSamples = TMath::Sqrt(fNumLoGainSamples); + if (fWindowSizeHiGain > availhirange) + { + *fLog << warn; + *fLog << GetDescriptor() << ": HiGain window " << (int)fWindowSizeHiGain; + *fLog << " out of range [" << (int)fHiGainFirst; + *fLog << "," << (int)fHiGainLast << "]" << endl; + *fLog << "Will set window size to " << (int)availhirange << endl; + fWindowSizeHiGain = availhirange; + } + + if (fWindowSizeLoGain > availlorange) + { + *fLog << warn; + *fLog << GetDescriptor() << ": LoGain window " << (int)fWindowSizeLoGain; + *fLog << " out of range [" << (int)fLoGainFirst; + *fLog << "," << (int)fLoGainLast << "]" << endl; + *fLog << "Will set window size to " << (int)availlorange << endl; + fWindowSizeLoGain = availlorange; + } + /* + fNumHiGainSamples = (Float_t)fWindowSizeHiGain; + fNumLoGainSamples = (Float_t)fWindowSizeLoGain; + + fSqrtHiGainSamples = TMath::Sqrt(fNumHiGainSamples); + fSqrtLoGainSamples = TMath::Sqrt(fNumLoGainSamples); + */ } @@ -369,98 +378,63 @@ Bool_t MPedCalcFromLoGain::ReInit(MParList *pList) { - Int_t lastdesired = (Int_t)fLoGainLast; - Int_t lastavailable = (Int_t)fRunHeader->GetNumSamplesLoGain()-1; - - if (lastdesired > lastavailable) - { - const Int_t diff = lastdesired - lastavailable; - *fLog << endl; - *fLog << warn << GetDescriptor() - << Form(": Selected Lo Gain FADC Window [%2i,%2i] ranges out of the available limits: [0,%2i].", - (int)fLoGainFirst, lastdesired, lastavailable) << endl; - *fLog << GetDescriptor() << ": Will reduce the upper edge to " << (int)(fLoGainLast - diff) << endl; - SetRange(fHiGainFirst, fHiGainLast, fLoGainFirst, fLoGainLast-diff); - } - - lastdesired = (Int_t)fHiGainLast; - lastavailable = (Int_t)fRunHeader->GetNumSamplesHiGain()-1; - - if (lastdesired > lastavailable) - { - const Int_t diff = lastdesired - lastavailable; - *fLog << endl; - *fLog << warn << GetDescriptor() - << Form(": Selected Hi Gain Range [%2i,%2i] ranges out of the available limits: [0,%2i].", - (int)fHiGainFirst, lastdesired, lastavailable) << endl; - *fLog << warn << GetDescriptor() - << Form(": Will possibly use %2i samples from the Low-Gain for the High-Gain range", diff) - << endl; - fHiGainLast -= diff; - fHiLoLast = diff; - } - - lastdesired = (Int_t)fHiGainFirst+fWindowSizeHiGain-1; - lastavailable = (Int_t)fRunHeader->GetNumSamplesHiGain()-1; - - if (lastdesired > lastavailable) - { - const Int_t diff = lastdesired - lastavailable; - *fLog << endl; - *fLog << warn << GetDescriptor() - << Form(": Selected Hi Gain FADC Window size %2i ranges out of the available limits: [0,%2i].", - (int)fWindowSizeHiGain, lastavailable) << endl; - *fLog << warn << GetDescriptor() - << Form(": Will use %2i samples from the Low-Gain for the High-Gain extraction", diff) - << endl; - - if ((Int_t)fWindowSizeHiGain > diff) - { - fWindowSizeHiGain -= diff; - fWindowSizeLoGain += diff; - } - else - { - fWindowSizeLoGain += fWindowSizeHiGain; - fLoGainFirst = diff-fWindowSizeHiGain; - fWindowSizeHiGain = 0; - } - } - - - // 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); - fSumAB0.Set(npixels); - fSumAB1.Set(npixels); - fNumEventsUsed.Set(npixels); - fTotalCounter.Set(npixels); - - // Reset contents of arrays. - fSumx.Reset(); - fSumx2.Reset(); - fSumAB0.Reset(); - fSumAB1.Reset(); - fNumEventsUsed.Reset(); - fTotalCounter.Reset(); - } - - 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 - << " HiGain FADC samples starting with slice: " << (int)fHiGainFirst << endl; - *fLog << inf << GetDescriptor() << ": Taking " << (int)fWindowSizeLoGain - << " LoGain FADC samples starting with slice: " << (int)fLoGainFirst << endl; - - return kTRUE; + if (fRunHeader->GetNumSamplesHiGain()<1) + { + *fLog << err << "ERROR - Number of available HiGainSamples<1... abort." << endl; + return kFALSE; + } + if (fRunHeader->GetNumSamplesLoGain()<1) + { + *fLog << err << "ERROR - Number of available LoGainSamples<1... abort." << endl; + return kFALSE; + } + + fHiGainFirst = TMath::Max(fHiGainFirst, 0); + fHiGainLast = TMath::Min(fHiGainLast, fRunHeader->GetNumSamplesHiGain()-1); + + fLoGainFirst = TMath::Max(fLoGainFirst, 0); + fLoGainLast = TMath::Min(fLoGainLast, fRunHeader->GetNumSamplesLoGain()-1); + + fWindowSizeHiGain = TMath::Min(fWindowSizeHiGain, fHiGainLast-fHiGainFirst+1); + fWindowSizeLoGain = TMath::Min(fWindowSizeLoGain, fLoGainLast-fLoGainFirst+1); + + SetRange(fHiGainFirst, fHiGainLast, fLoGainFirst, fLoGainLast); + + *fLog << inf << endl; + *fLog << "Taking " << Form("%2d", (int)fWindowSizeHiGain) << " HiGain from " << (int)fHiGainFirst << endl; + *fLog << "Taking " << Form("%2d", (int)fWindowSizeLoGain) << " LoGain from " << (int)fLoGainFirst << endl; + + if (fWindowSizeHiGain==0) + { + *fLog << err << "ERROR - HiGain windows size == 0... abort." << endl; + return kFALSE; + } + if (fWindowSizeLoGain==0) + { + *fLog << err << "ERROR - HiGain windows size == 0... abort." << endl; + return kFALSE; + } + + // 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); + fSumAB0.Set(npixels); + fSumAB1.Set(npixels); + fNumEventsUsed.Set(npixels); + fTotalCounter.Set(npixels); + + // Reset contents of arrays. + fSumx.Reset(); + fSumx2.Reset(); + fSumAB0.Reset(); + fSumAB1.Reset(); + fNumEventsUsed.Reset(); + fTotalCounter.Reset(); + } + + return kTRUE; } @@ -515,10 +489,8 @@ // Find the maximum and minimum signal per slice in the high gain window do { - if (*ptr > max) { + if (*ptr > max) max = *ptr; - } - if (*ptr < min) { + if (*ptr < min) min = *ptr; - } } while (++ptr != end); @@ -527,12 +499,13 @@ continue; - Byte_t *firstSlice = pixel.GetLoGainSamples() + fLoGainFirst; - Byte_t *lastSlice = firstSlice + fWindowSizeLoGain; - - Byte_t *slice = firstSlice; + ptr = pixel.GetLoGainSamples() + fLoGainFirst; + end = ptr + fWindowSizeLoGain; + + Byte_t *firstSlice = ptr; + do { - sum += *slice; - sqr += *slice * *slice; - } while (++slice != lastSlice); + sum += *ptr; + sqr += *ptr * *ptr; + } while (++ptr != end); const Float_t msum = (Float_t)sum; @@ -569,5 +542,8 @@ if (fPedestalUpdate) + { + fPedestals->ReCalc(*fGeom); fPedestals->SetReadyToSave(); + } return kTRUE; @@ -581,8 +557,10 @@ { // Compute pedestals and rms from the whole run - if (fPedestalUpdate) + if (fPedestalUpdate || GetNumExecutions()<1) return kTRUE; - *fLog << flush << inf << "Calculating Pedestals..." << flush; + *fLog << flush << inf << "Calculating pedestals..." << flush; + + Double_t sum = 0; const Int_t npix = fGeom->GetNumPixels(); @@ -592,6 +570,11 @@ if (n>1) Calc(n, idx); - } - + sum += n; + } + + *fLog << flush << inf << "Calculating means..." << flush; + + fPedestals->SetTotalEntries((UInt_t)(sum/npix*(fWindowSizeLoGain+fWindowSizeHiGain))); + fPedestals->ReCalc(*fGeom); fPedestals->SetReadyToSave(); return kTRUE; @@ -621,27 +604,21 @@ SetWindowSize(hw, lw); - Int_t num = fNumEventsDump; if (IsEnvDefined(env, prefix, "NumEventsDump", print)) { - num = GetEnvValue(env, prefix, "NumEventsDump", num); + SetDumpEvents(GetEnvValue(env, prefix, "NumEventsDump", fNumEventsDump)); rc = kTRUE; } - SetDumpEvents(num); - - Byte_t max = fMaxHiGainVar; + if (IsEnvDefined(env, prefix, "MaxHiGainVar", print)) { - max = GetEnvValue(env, prefix, "MaxHiGainVar", max); + SetMaxHiGainVar(GetEnvValue(env, prefix, "MaxHiGainVar", fMaxHiGainVar)); rc = kTRUE; } - SetMaxHiGainVar(max); - - Bool_t upd = fPedestalUpdate; + if (IsEnvDefined(env, prefix, "PedestalUpdate", print)) { - upd = GetEnvValue(env, prefix, "PedestalUpdate", upd); + SetPedestalUpdate(GetEnvValue(env, prefix, "PedestalUpdate", fPedestalUpdate)); rc = kTRUE; } - SetPedestalUpdate(upd); return rc; Index: trunk/MagicSoft/Mars/mpedestal/MPedCalcFromLoGain.h =================================================================== --- trunk/MagicSoft/Mars/mpedestal/MPedCalcFromLoGain.h (revision 4601) +++ trunk/MagicSoft/Mars/mpedestal/MPedCalcFromLoGain.h (revision 4609) @@ -62,7 +62,7 @@ // 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 SetRange(Byte_t hifirst=fgHiGainFirst, Byte_t hilast=fgHiGainLast, Byte_t lofirst=fgLoGainFirst, Byte_t lolast=fgLoGainLast); + void SetWindowSize(Byte_t windowh=fgHiGainWindowSize, Byte_t windowl=fgLoGainWindowSize); + void SetMaxHiGainVar(Byte_t maxvar=fgMaxHiGainVar); void SetDumpEvents(UInt_t dumpevents = 0) {fNumEventsDump = dumpevents;} void SetPedestalUpdate(Bool_t pedupdate) {fPedestalUpdate = pedupdate;} Index: trunk/MagicSoft/Mars/mpedestal/MPedCalcPedRun.cc =================================================================== --- trunk/MagicSoft/Mars/mpedestal/MPedCalcPedRun.cc (revision 4601) +++ trunk/MagicSoft/Mars/mpedestal/MPedCalcPedRun.cc (revision 4609) @@ -394,6 +394,6 @@ const Int_t npixels = fPedestals->GetSize(); - const Int_t areas = fPedestals->GetAverageAreas(); - const Int_t sectors = fPedestals->GetAverageSectors(); +// const Int_t areas = fPedestals->GetAverageAreas(); +// const Int_t sectors = fPedestals->GetAverageSectors(); if (fSumx.GetSize()==0) @@ -402,11 +402,11 @@ fSumx2.Set(npixels); - fAreaSumx. Set(areas); - fAreaSumx2.Set(areas); - fAreaValid.Set(areas); - - fSectorSumx. Set(sectors); - fSectorSumx2.Set(sectors); - fSectorValid.Set(sectors); +// fAreaSumx. Set(areas); +// fAreaSumx2.Set(areas); +// fAreaValid.Set(areas); + +// fSectorSumx. Set(sectors); +// fSectorSumx2.Set(sectors); +// fSectorValid.Set(sectors); fSumx.Reset(); @@ -439,4 +439,66 @@ Int_t MPedCalcPedRun::Process() { + 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; + + if (fWindowSizeHiGain != 0) + { + do + { + sum += *ptr; + sqr += *ptr * *ptr; + } + while (++ptr != end); + } + + if (fWindowSizeLoGain != 0) + { + ptr = pixel.GetLoGainSamples() + fLoGainFirst; + end = ptr + fWindowSizeLoGain; + + do + { + sum += *ptr; + sqr += *ptr * *ptr; + } + while (++ptr != end); + } + + const Float_t msum = (Float_t)sum; + + // + // These three lines have been uncommented by Markus Gaug + // If anybody needs them, please contact me!! + // + // const Float_t higainped = msum/fNumHiGainSlices; + // const Float_t higainrms = TMath::Sqrt((msqr-msum*msum/fNumHiGainSlices)/(fNumHiGainSlices-1.)); + // (*fPedestals)[idx].Set(higainped, higainrms); + + fSumx[idx] += msum; + // + // The old version: + // + // const Float_t msqr = (Float_t)sqr; + // fSumx2[idx] += msqr; + // + // The new version: + // + const Float_t sqrsum = msum*msum; + fSumx2[idx] += sqrsum; + } + + fNumSamplesTot += fWindowSizeHiGain + fWindowSizeLoGain; + + return kTRUE; +/* MRawEvtPixelIter pixel(fRawEvt); @@ -445,5 +507,5 @@ const UInt_t idx = pixel.GetPixelId(); const UInt_t aidx = (*fGeom)[idx].GetAidx(); - const UInt_t sector = (*fGeom)[idx].GetSector(); + const UInt_t sector = (*fGeom)[idx].GetSector(); Byte_t *ptr = pixel.GetHiGainSamples() + fHiGainFirst; @@ -490,5 +552,5 @@ fSumx[idx] += msum; fAreaSumx[aidx] += msum; - fSectorSumx[sector] += msum; + fSectorSumx[sector] += msum; // // The old version: @@ -502,10 +564,10 @@ fSumx2[idx] += sqrsum; fAreaSumx2[aidx] += sqrsum; - fSectorSumx2[sector] += sqrsum; + fSectorSumx2[sector] += sqrsum; } fPedestals->SetReadyToSave(); fNumSamplesTot += fWindowSizeHiGain + fWindowSizeLoGain; - +*/ return kTRUE; } @@ -517,42 +579,82 @@ Int_t MPedCalcPedRun::PostProcess() { - // Compute pedestals and rms from the whole run - const ULong_t n = fNumSamplesTot; - const ULong_t nevts = GetNumExecutions(); - - MRawEvtPixelIter pixel(fRawEvt); - - while (pixel.Next()) - { - - const Int_t pixid = pixel.GetPixelId(); - const UInt_t aidx = (*fGeom)[pixid].GetAidx(); - const UInt_t sector = (*fGeom)[pixid].GetSector(); - - fAreaValid [aidx]++; - fSectorValid[sector]++; - - const Float_t sum = fSumx.At(pixid); - const Float_t sum2 = fSumx2.At(pixid); - const Float_t higainped = sum/n; - // - // The old version: - // - // const Float_t higainrms = TMath::Sqrt((sum2-sum*sum/n)/(n-1.)); - // - // The new version: - // - // 1. Calculate the Variance of the sums: - Float_t higainVar = (sum2-sum*sum/nevts)/(nevts-1.); - // 2. Scale the variance to the number of slices: - higainVar /= (Float_t)(fWindowSizeHiGain+fWindowSizeLoGain); - // 3. Calculate the RMS from the Variance: - (*fPedestals)[pixid].Set(higainped, higainVar < 0 ? 0. : TMath::Sqrt(higainVar)); - - } + const ULong_t nevts = GetNumExecutions(); + if (nevts<1) + { + *fLog << err << "ERROR - Not enough events recorded...abort." << endl; + return kFALSE; + } + + *fLog << flush << inf << "Calculating pedestals..." << flush; + + // Compute pedestals and rms from the whole run + const ULong_t n = fNumSamplesTot; + const ULong_t npix = fGeom->GetNumPixels(); + + for (UInt_t pixidx=0; pixidxSetTotalEntries(fNumSamplesTot); + fPedestals->ReCalc(*fGeom); + fPedestals->SetReadyToSave(); + + return kTRUE; // // Loop over the (two) area indices to get the averaged pedestal per aidx // + /* + + // Compute pedestals and rms from the whole run + const ULong_t n = fNumSamplesTot; + const ULong_t nevts = GetNumExecutions(); + + MRawEvtPixelIter pixel(fRawEvt); + + while (pixel.Next()) + { + const Int_t pixid = pixel.GetPixelId(); + const UInt_t aidx = (*fGeom)[pixid].GetAidx(); + const UInt_t sector = (*fGeom)[pixid].GetSector(); + + fAreaValid [aidx]++; + fSectorValid[sector]++; + + const Float_t sum = fSumx.At(pixid); + const Float_t sum2 = fSumx2.At(pixid); + const Float_t higainped = sum/n; + // + // The old version: + // + // const Float_t higainrms = TMath::Sqrt((sum2-sum*sum/n)/(n-1.)); + // + // The new version: + // + // 1. Calculate the Variance of the sums: + Float_t higainVar = (sum2-sum*sum/nevts)/(nevts-1.); + // 2. Scale the variance to the number of slices: + higainVar /= (Float_t)(fWindowSizeHiGain+fWindowSizeLoGain); + // 3. Calculate the RMS from the Variance: + (*fPedestals)[pixid].Set(higainped, higainVar < 0 ? 0. : TMath::Sqrt(higainVar), 0, nevts); + } for (Int_t aidx=0; aidxGetAverageArea(aidx).Set(higainped, higainrms); } - + */ // // Loop over the (six) sector indices to get the averaged pedestal per sector // + /* for (Int_t sector=0; sectorGetAverageSector(sector).Set(higainped, higainrms); } - fPedestals->SetTotalEntries(fNumSamplesTot); fPedestals->SetReadyToSave(); return kTRUE; + */ } Index: trunk/MagicSoft/Mars/mpedestal/MPedCalcPedRun.h =================================================================== --- trunk/MagicSoft/Mars/mpedestal/MPedCalcPedRun.h (revision 4601) +++ trunk/MagicSoft/Mars/mpedestal/MPedCalcPedRun.h (revision 4609) @@ -33,4 +33,5 @@ TArrayD fSumx; // sum of values TArrayD fSumx2; // sum of squared values + /* TArrayD fAreaSumx; // averaged sum of values per area idx TArrayD fAreaSumx2; // averaged sum of squared values per area idx @@ -39,4 +40,5 @@ TArrayD fSectorSumx2; // averaged sum of squared values per sector TArrayI fSectorValid; // number of valid pixel with sector idx + */ Int_t PreProcess (MParList *pList); Index: trunk/MagicSoft/Mars/mpedestal/MPedPhotCam.cc =================================================================== --- trunk/MagicSoft/Mars/mpedestal/MPedPhotCam.cc (revision 4601) +++ trunk/MagicSoft/Mars/mpedestal/MPedPhotCam.cc (revision 4609) @@ -229,10 +229,7 @@ void MPedPhotCam::ReCalc(const MGeomCam &geom, MBadPixelsCam *bad) { - const Int_t np = GetSize(); const Int_t ns = GetNumSectors(); const Int_t na = GetNumAreas(); - - TArrayI acnt(na); @@ -243,11 +240,8 @@ TArrayD ssumr(ns); - for (int i=0; i0) GetSector(i).Set(ssumx[i]/scnt[i], ssumr[i]/scnt[i], scnt[i]); - else GetSector(i).Set(-1., -1., 0); - } - - for (int i=0; i0) GetArea(i).Set(asumx[i]/acnt[i], asumr[i]/acnt[i], acnt[i]); - else GetArea(i).Set(-1., -1., 0); - } + } + + for (int i=0; i0) + GetSector(i).Set(ssumx[i]/scnt[i], ssumr[i]/scnt[i], scnt[i]); + else + GetSector(i).Clear(); + + for (int i=0; i0) + GetArea(i).Set(asumx[i]/acnt[i], asumr[i]/acnt[i], acnt[i]); + else + GetArea(i).Clear(); } Index: trunk/MagicSoft/Mars/mpedestal/MPedPhotCam.h =================================================================== --- trunk/MagicSoft/Mars/mpedestal/MPedPhotCam.h (revision 4601) +++ trunk/MagicSoft/Mars/mpedestal/MPedPhotCam.h (revision 4609) @@ -51,5 +51,5 @@ void Print(Option_t *o="") const; - void ReCalc(const MGeomCam &geom, MBadPixelsCam *bad); + void ReCalc(const MGeomCam &geom, MBadPixelsCam *bad=NULL); Bool_t GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type=0) const; Index: trunk/MagicSoft/Mars/mpedestal/MPedestalCam.cc =================================================================== --- trunk/MagicSoft/Mars/mpedestal/MPedestalCam.cc (revision 4601) +++ trunk/MagicSoft/Mars/mpedestal/MPedestalCam.cc (revision 4609) @@ -35,4 +35,7 @@ #include "MPedestalPix.h" +#include +#include + #include @@ -43,4 +46,8 @@ #include "MGeomCam.h" +#include "MGeomPix.h" + +#include "MBadPixelsCam.h" +#include "MBadPixelsPix.h" ClassImp(MPedestalCam); @@ -311,37 +318,96 @@ } -Bool_t MPedestalCam::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const { - - if (GetSize() <= idx) - return kFALSE; - - if (!(*this)[idx].IsValid()) - return kFALSE; - - switch (type) { +// -------------------------------------------------------------------------- +// +// Calculates the avarage pedestal and pedestal rms for all sectors +// and pixel sizes. The geometry container is used to get the necessary +// geometry information (sector number, size index) If the bad pixel +// container is given all pixels which have the flag 'bad' are ignored +// in the calculation of the sector and size average. +// +void MPedestalCam::ReCalc(const MGeomCam &geom, MBadPixelsCam *bad) +{ + const Int_t np = GetSize(); + const Int_t ns = geom.GetNumSectors(); + const Int_t na = geom.GetNumAreas(); + + TArrayI acnt(na); + TArrayI scnt(ns); + TArrayD asumx(na); + TArrayD ssumx(ns); + TArrayD asumr(na); + TArrayD ssumr(ns); + + for (int i=0; i0) + GetAverageSector(i).Set(ssumx[i]/scnt[i], ssumr[i]/scnt[i], 0, scnt[i]); + else + GetAverageSector(i).Clear(); + + for (int i=0; i0) + GetAverageArea(i).Set(asumx[i]/acnt[i], asumr[i]/acnt[i], 0, acnt[i]); + else + GetAverageArea(i).Clear(); +} + +Bool_t MPedestalCam::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const +{ + if (GetSize() <= idx) + return kFALSE; + + if (!(*this)[idx].IsValid()) + return kFALSE; + + switch (type) + { case 0: - val = (*this)[idx].GetPedestal(); - break; + val = (*this)[idx].GetPedestal(); + break; case 1: - val = fTotalEntries > 0 ? - (*this)[idx].GetPedestalRms()/TMath::Sqrt((Float_t)fTotalEntries) - : (*this)[idx].GetPedestalError(); - break; + val = fTotalEntries > 0 ? + (*this)[idx].GetPedestalRms()/TMath::Sqrt((Float_t)fTotalEntries) + : (*this)[idx].GetPedestalError(); + break; case 2: - val = (*this)[idx].GetPedestalRms(); - break; + val = (*this)[idx].GetPedestalRms(); + break; case 3: - val = fTotalEntries > 0 ? - (*this)[idx].GetPedestalRms()/TMath::Sqrt((Float_t)fTotalEntries)/2. - : (*this)[idx].GetPedestalRmsError(); - break; + val = fTotalEntries > 0 ? + (*this)[idx].GetPedestalRms()/TMath::Sqrt((Float_t)fTotalEntries)/2. + : (*this)[idx].GetPedestalRmsError(); + break; default: - return kFALSE; + return kFALSE; } - return kTRUE; + return kTRUE; } void MPedestalCam::DrawPixelContent(Int_t idx) const { - *fLog << warn << "MPedestalCam::DrawPixelContent - not available." << endl; -} + *fLog << warn << "MPedestalCam::DrawPixelContent - not available." << endl; +} Index: trunk/MagicSoft/Mars/mpedestal/MPedestalCam.h =================================================================== --- trunk/MagicSoft/Mars/mpedestal/MPedestalCam.h (revision 4601) +++ trunk/MagicSoft/Mars/mpedestal/MPedestalCam.h (revision 4609) @@ -13,4 +13,5 @@ class MGeomCam; class MPedestalPix; +class MBadPixelsCam; class MPedestalCam : public MParContainer, public MCamEvent @@ -51,4 +52,6 @@ void InitAverageSectors ( const UInt_t i ); + void ReCalc(const MGeomCam &geom, MBadPixelsCam *bad=NULL); + void Print(Option_t *o="") const; Index: trunk/MagicSoft/Mars/mpedestal/MPedestalPix.h =================================================================== --- trunk/MagicSoft/Mars/mpedestal/MPedestalPix.h (revision 4601) +++ trunk/MagicSoft/Mars/mpedestal/MPedestalPix.h (revision 4609) @@ -37,4 +37,6 @@ Float_t GetPedestalRmsError() const { return fNumEvents>0 ? fPedestalRms/TMath::Sqrt((Float_t)fNumEvents/2) : 0; } + UInt_t GetNumEvents() const { return fNumEvents; } + Bool_t IsValid() const;