Index: /trunk/MagicSoft/Mars/Changelog =================================================================== --- /trunk/MagicSoft/Mars/Changelog (revision 3263) +++ /trunk/MagicSoft/Mars/Changelog (revision 3264) @@ -41,4 +41,5 @@ * mcalib/MCalibrationChargeCalc.[h,cc] * mcalib/MCalibrationChargeCam.[h,cc] + * mcalib/MCalibrationChargePix.[h,cc] * mcalib/MCalibrationChargePINDiode.[h,cc] * mcalib/MCalibrationCalc.[h,cc] Index: /trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCalc.cc =================================================================== --- /trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCalc.cc (revision 3263) +++ /trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCalc.cc (revision 3264) @@ -16,7 +16,7 @@ ! ! -! Author(s): Markus Gaug 09/2003 +! Author(s): Markus Gaug 02/2004 ! -! Copyright: MAGIC Software Development, 2000-2001 +! Copyright: MAGIC Software Development, 2000-2004 ! ! @@ -31,5 +31,5 @@ // MExtractedSignalCam. The pedestals by an MPedestalCam. // -// The output container MCalibrationCam holds one entry of type MCalibrationPix +// The output container MCalibrationCam holds one entry of type MCalibrationChargePix // for every pixel. It is filled in the following way: // @@ -39,9 +39,9 @@ // // ReInit: MCalibrationCam::InitSize(NumPixels) is called which allocates -// memory in a TClonesArray of type MCalibrationPix +// memory in a TClonesArray of type MCalibrationChargePix // Initialize number of used FADC slices // Optionally exclude pixels from calibration // -// Process: Every MCalibrationPix holds a histogram class, +// Process: Every MCalibrationChargePix holds a histogram class, // MHCalibrationPixel which itself hold histograms of type: // HCharge(npix) (distribution of summed FADC time slice @@ -97,5 +97,5 @@ #include "MCalibrationChargeCam.h" #include "MCalibrationChargePINDiode.h" -#include "MCalibrationPix.h" +#include "MCalibrationChargePix.h" #include "MExtractedSignalCam.h" @@ -142,7 +142,4 @@ SETBIT(fFlags, kUseQualityChecks); SETBIT(fFlags, kHiLoGainCalibration); - - CLRBIT(fFlags, kHiGainOverFlow); - CLRBIT(fFlags, kLoGainOverFlow); fNumBlindPixelSinglePhe = 0; @@ -270,5 +267,5 @@ { - MCalibrationPix &pix = (*fCam)[i]; + MCalibrationChargePix &pix = (*fCam)[i]; pix.DefinePixId(i); @@ -278,6 +275,4 @@ fSignals->GetLastUsedSliceLoGain()); - if (!TESTBIT(fFlags,kUseQualityChecks)) - pix.SetExcludeQualityCheck(); // Exclude the blind pixel and the PIN Diode from normal pixel calibration: @@ -330,5 +325,5 @@ // Exclude pixel // - MCalibrationPix &pix = (*fCam)[pixel]; + MCalibrationChargePix &pix = (*fCam)[pixel]; pix.SetExcluded(); @@ -412,53 +407,4 @@ hist->Fill(i,*ptr++); - pixel.Reset(); - - while (pixel.Next()) - { - - const UInt_t pixid = pixel.GetPixelId(); - - MCalibrationPix &pix = (*fCam)[pixid]; - MExtractedSignalPix &sig = (*fSignals) [pixid]; - - const Float_t sumhi = sig.GetExtractedSignalHiGain(); - const Float_t sumlo = sig.GetExtractedSignalLoGain(); - - Float_t abstime = 0.; - - if (sig.IsLoGainUsed()) - abstime = (Float_t)pixel.GetIdxMaxLoGainSample(); - else - abstime = (Float_t)pixel.GetIdxMaxHiGainSample(); - - if (pix.IsExcluded()) - continue; - - pix.FillGraphs(sumhi,sumlo); - - if (sig.IsLoGainUsed()) - { - - if (!pix.FillChargeLoGain(sumlo)) - *fLog << warn << "Could not fill Lo Gain Charge of pixel: " << pixid - << " signal = " << sumlo << endl; - - if (!pix.FillAbsTimeLoGain(abstime)) - *fLog << warn << "Could not fill Lo Gain Abs. Time of pixel: " - << pixid << " time = " << abstime << endl; - } - else - { - if (!pix.FillChargeHiGain(sumhi)) - *fLog << warn << "Could not fill Hi Gain Charge of pixel: " << pixid - << " signal = " << sumhi << endl; - - if (!pix.FillAbsTimeHiGain(abstime)) - *fLog << warn << "Could not fill Hi Gain Abs. Time of pixel: " - << pixid << " time = " << abstime << endl; - } - - } /* while (pixel.Next()) */ - return kTRUE; } @@ -467,9 +413,4 @@ { *fLog << inf << GetDescriptor() << ": Cut Histogram Edges" << endl; - - // - // Cut edges to make fits and viewing of the hists easier - // - fCam->CutEdges(); // @@ -550,5 +491,5 @@ { - MCalibrationPix &pix = (*fCam)[pixid]; + MCalibrationChargePix &pix = (*fCam)[pixid]; // @@ -563,24 +504,32 @@ const Float_t ped = (*fPedestals)[pixid].GetPedestal(); const Float_t prms = (*fPedestals)[pixid].GetPedestalRms(); + const Int_t num = fPedestals->GetTotalEntries(); // // set them in the calibration camera // - pix.SetPedestal(ped,prms,(Float_t)fNumHiGainSamples,(Float_t)fNumLoGainSamples); - - // - // perform the Gauss fits to the charges - // - pix.FitCharge(); - - // - // check also for oscillations - // - pix.CheckOscillations(); - - // - // calculate the F-Factor method - // + if (pix.IsHiGainSaturation()) + { + pix.SetPedestal(ped * (Float_t)fNumLoGainSamples, + prms * TMath::Sqrt((Float_t)fNumLoGainSamples), + prms * (Float_t)fNumLoGainSamples / num); + pix.SetNumLoGainSamples((Float_t)fNumLoGainSamples); + pix.ApplyLoGainConversion(); + } + else + { + pix.SetPedestal(ped * (Float_t)fNumHiGainSamples, + prms * TMath::Sqrt((Float_t)fNumHiGainSamples), + prms * (Float_t)fNumHiGainSamples / num); + } + + if (!pix.CheckChargeValidity() || !pix.CheckTimeValidity()) + continue; + + if (!pix.CalcReducedSigma()) + continue; + pix.CalcFFactorMethod(); + } @@ -597,14 +546,21 @@ *fLog << inf << GetDescriptor() << ": Skipping Blind Pixel Calibration! " << endl; - if (!fPINDiode->CalcFluxOutsidePlexiglass()) + if (!fPINDiode->CheckChargeFitValidity() || !fPINDiode->CheckTimeFitValidity()) { - *fLog << warn << "Could not calculate the flux of photons from the PIN Diode, will skip PIN Diode Calibration " << endl; + *fLog << warn << "Could not calculate the flux of photons from the PIN Diode, charge fit not valid " << endl; fCam->SetPINDiodeMethodValid(kFALSE); } else - { - fCam->SetPINDiodeMethodValid(kTRUE); + { + if (!fPINDiode->CalcFluxOutsidePlexiglass()) + { + *fLog << warn << "Could not calculate the flux of photons from the PIN Diode, will skip PIN Diode Calibration " << endl; + fCam->SetPINDiodeMethodValid(kFALSE); + } + else + { + fCam->SetPINDiodeMethodValid(kTRUE); + } } - fCam->SetReadyToSave(); @@ -613,2 +569,4 @@ + + Index: /trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCalc.h =================================================================== --- /trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCalc.h (revision 3263) +++ /trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCalc.h (revision 3264) @@ -67,12 +67,10 @@ enum { kUseBlindPixelFit, kUseQualityChecks, - kHiLoGainCalibration, - kHiGainOverFlow, - kLoGainOverFlow }; + kHiLoGainCalibration }; + Int_t PreProcess(MParList *pList); Bool_t ReInit(MParList *pList); - Int_t PreProcess(MParList *pList); - Int_t Process(); - Int_t PostProcess(); + Int_t Process(); + Int_t PostProcess(); public: @@ -91,5 +89,4 @@ // Setters - void SetConversionHiLo ( const Float_t conv ) { fConversionHiLo = conv; } void SetBlindPixelSinglePheCut ( const Int_t cut=fgBlindPixelSinglePheCut) { fBlindPixelSinglePheCut = cut; } Index: /trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCam.cc =================================================================== --- /trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCam.cc (revision 3263) +++ /trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCam.cc (revision 3264) @@ -126,9 +126,8 @@ #include "MGeomPix.h" -#include "MCalibrationPix.h" +#include "MCalibrationChargePix.h" #include "MCalibrationBlindPix.h" #include "MCalibrationChargePINDiode.h" -#include "MHCalibrationPixel.h" ClassImp(MCalibrationChargeCam); @@ -174,5 +173,5 @@ fTitle = title ? title : "Storage container for the Calibration Information in the camera"; - fPixels = new TClonesArray("MCalibrationPix",1); + fPixels = new TClonesArray("MCalibrationChargePix",1); fBlindPixel = new MCalibrationBlindPix(); @@ -252,7 +251,7 @@ // Get i-th pixel (pixel number) // -MCalibrationPix &MCalibrationChargeCam::operator[](UInt_t i) -{ - return *static_cast(fPixels->UncheckedAt(i)); +MCalibrationChargePix &MCalibrationChargeCam::operator[](UInt_t i) +{ + return *static_cast(fPixels->UncheckedAt(i)); } @@ -261,7 +260,7 @@ // Get i-th pixel (pixel number) // -const MCalibrationPix &MCalibrationChargeCam::operator[](UInt_t i) const -{ - return *static_cast(fPixels->UncheckedAt(i)); +const MCalibrationChargePix &MCalibrationChargeCam::operator[](UInt_t i) const +{ + return *static_cast(fPixels->UncheckedAt(i)); } @@ -329,13 +328,13 @@ TIter Next(fPixels); - MCalibrationPix *pix; - while ((pix=(MCalibrationPix*)Next())) + MCalibrationChargePix *pix; + while ((pix=(MCalibrationChargePix*)Next())) { - if (pix->IsChargeValid() && !pix->IsExcluded() && !pix->IsOscillating() && pix->IsFitted()) + if (pix->IsChargeValid() && !pix->IsExcluded() && !pix->IsOscillating()) { *fLog << all << pix->GetPixId() << " Pedestals: " << pix->GetPed() << " +- " - << pix->GetPedRms() << " Reduced Charge: " << pix->GetCharge() << " +- " + << pix->GetPedRms() << " Reduced Charge: " << pix->GetMeanCharge() << " +- " << pix->GetSigmaCharge() << " Reduced Sigma: " << pix->GetRSigmaCharge() << " Nr Phe's: " << pix->GetPheFFactorMethod() << endl; @@ -352,12 +351,12 @@ TIter Next2(fPixels); - while ((pix=(MCalibrationPix*)Next2())) + while ((pix=(MCalibrationChargePix*)Next2())) { - if (!pix->IsChargeValid() && !pix->IsExcluded() && !pix->IsFitted()) + if (!pix->IsExcluded() && !pix->IsChargeValid()) { *fLog << all << pix->GetPixId() << " Pedestals: " << pix->GetPed() << " +- " - << pix->GetPedRms() << " Reduced Charge: " << pix->GetCharge() << " +- " + << pix->GetPedRms() << " Reduced Charge: " << pix->GetMeanCharge() << " +- " << pix->GetSigmaCharge() << " Reduced Sigma: " << pix->GetRSigmaCharge() << endl; id++; @@ -373,12 +372,12 @@ TIter Next3(fPixels); - while ((pix=(MCalibrationPix*)Next3())) + while ((pix=(MCalibrationChargePix*)Next3())) { - - if (pix->IsOscillating() && !pix->IsExcluded()) + + if ( pix->IsOscillating() && !pix->IsExcluded()) { *fLog << all << pix->GetPixId() << " Pedestals: " << pix->GetPed() << " +- " - << pix->GetPedRms() << " Reduced Charge: " << pix->GetCharge() << " +- " + << pix->GetPedRms() << " Reduced Charge: " << pix->GetMeanCharge() << " +- " << pix->GetSigmaCharge() << " Reduced Sigma: " << pix->GetRSigmaCharge() << endl; id++; @@ -392,36 +391,18 @@ *fLog << all << endl; + id = 0; + TIter Next4(fPixels); - while ((pix=(MCalibrationPix*)Next4())) - if (pix->IsExcluded()) - *fLog << all << pix->GetPixId() << endl; - - *fLog << all << fNumExcludedPixels << " excluded pixels " << endl; -} - -// -------------------------------------------------------------------------- -// -// Return true if pixel is inside bounds of the TClonesArray fPixels -// -Bool_t MCalibrationChargeCam::IsPixelUsed(Int_t idx) const -{ - if (!CheckBounds(idx)) - return kFALSE; - - return kTRUE; -} - -// -------------------------------------------------------------------------- -// -// Return true if pixel has already been fitted once (independent of the result) -// -Bool_t MCalibrationChargeCam::IsPixelFitted(Int_t idx) const -{ - - if (!CheckBounds(idx)) - return kFALSE; - - return (*this)[idx].IsFitted(); -} + while ((pix=(MCalibrationChargePix*)Next4())) + { + if (pix->IsExcluded()) + { + *fLog << all << pix->GetPixId() << endl; + id++; + } + } + *fLog << all << id << " Excluded pixels " << endl; +} + // -------------------------------------------------------------------------- @@ -435,11 +416,4 @@ fBlindPixel->GetHist()->CutAllEdges(); - - TIter Next(fPixels); - MCalibrationPix *pix; - while ((pix=(MCalibrationPix*)Next())) - { - pix->GetHist()->CutAllEdges(); - } return; @@ -495,33 +469,34 @@ // 28: Pixels where the fit did not succeed --> results obtained only from the histograms // 29: Pixels with succeeded fit, but apparently wrong results -// 30: Pixels with un-expected behavior in the fourier spectrum (e.g. oscillations) +// 30: Pixels with un-expected behavior in the Hi Gain fourier spectrum (e.g. oscillations) +// 31: Pixels with un-expected behavior in the Lo Gain fourier spectrum (e.g. oscillations)a +// 32: Number of probable pickup events in the Hi Gain +// 33: Number of probable pickup events in the Lo Gain // // Other classifications of pixels: // ================================ // -// 31: Pixels with saturated Hi-Gain +// 34: Pixels with saturated Hi-Gain // // Classification of validity of the calibrations: // =============================================== // -// 32: Pixels with valid calibration by the F-Factor-Method -// 33: Pixels with valid calibration by the Blind Pixel-Method -// 34: Pixels with valid calibration by the PIN Diode-Method +// 35: Pixels with valid calibration by the F-Factor-Method +// 36: Pixels with valid calibration by the Blind Pixel-Method +// 37: Pixels with valid calibration by the PIN Diode-Method // // Used Pedestals: // =============== // -// 35: Mean Pedestal over the entire range of signal extraction -// 36: Error on the Mean Pedestal over the entire range of signal extraction -// 37: Pedestal RMS over the entire range of signal extraction -// 38: Error on the Pedestal RMS over the entire range of signal extraction +// 38: Mean Pedestal over the entire range of signal extraction +// 39: Error on the Mean Pedestal over the entire range of signal extraction +// 40: Pedestal RMS over the entire range of signal extraction +// 41: Error on the Pedestal RMS over the entire range of signal extraction // // Calculated absolute arrival times (very low precision!): // ======================================================== // -// 39: Absolute Arrival time of the signal -// 40: Error on the Absolute Arrival time of the signal -// 41: RMS of the Absolute Arrival time of the signal -// 42: Error on the RMS of the Absolute Arrival time of the signal +// 42: Absolute Arrival time of the signal +// 43: RMS of the Absolute Arrival time of the signal // Bool_t MCalibrationChargeCam::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const @@ -541,10 +516,10 @@ if ((*this)[idx].IsExcluded()) return kFALSE; - val = (*this)[idx].GetCharge(); + val = (*this)[idx].GetMeanCharge(); break; case 1: if ((*this)[idx].IsExcluded()) return kFALSE; - val = (*this)[idx].GetChargeErr(); + val = (*this)[idx].GetMeanChargeErr(); break; case 2: @@ -576,5 +551,5 @@ if ((*this)[idx].IsExcluded()) return kFALSE; - val = (*this)[idx].GetRSigmaCharge() / (*this)[idx].GetCharge(); + val = (*this)[idx].GetRSigmaCharge() / (*this)[idx].GetMeanCharge(); break; case 8: @@ -585,10 +560,10 @@ / ((*this)[idx].GetRSigmaCharge() * (*this)[idx].GetRSigmaCharge() ); // relative error Charge square - val += (*this)[idx].GetChargeErr() * (*this)[idx].GetChargeErr() - / ((*this)[idx].GetCharge() * (*this)[idx].GetCharge() ); + val += (*this)[idx].GetMeanChargeErr() * (*this)[idx].GetMeanChargeErr() + / ((*this)[idx].GetMeanCharge() * (*this)[idx].GetMeanCharge() ); // calculate relative error out of squares val = TMath::Sqrt(val) ; // multiply with value to get absolute error - val *= (*this)[idx].GetRSigmaCharge() / (*this)[idx].GetCharge(); + val *= (*this)[idx].GetRSigmaCharge() / (*this)[idx].GetMeanCharge(); break; case 9: @@ -709,5 +684,5 @@ if ((*this)[idx].IsExcluded()) return kFALSE; - if ((*this)[idx].IsOscillating()) + if ((*this)[idx].IsHiGainOscillating()) val = 1; else @@ -715,4 +690,22 @@ break; case 31: + if ((*this)[idx].IsExcluded()) + return kFALSE; + if ((*this)[idx].IsLoGainOscillating()) + val = 1; + else + return kFALSE; + break; + case 32: + if ((*this)[idx].IsExcluded()) + return kFALSE; + val = (*this)[idx].GetHiGainNumPickup(); + break; + case 33: + if ((*this)[idx].IsExcluded()) + return kFALSE; + val = (*this)[idx].GetLoGainNumPickup(); + break; + case 34: if ((*this)[idx].IsExcluded()) return kFALSE; @@ -722,5 +715,5 @@ return kFALSE; break; - case 32: + case 35: if ((*this)[idx].IsExcluded()) return kFALSE; @@ -730,5 +723,5 @@ return kFALSE; break; - case 33: + case 36: if ((*this)[idx].IsExcluded()) return kFALSE; @@ -738,5 +731,5 @@ return kFALSE; break; - case 34: + case 37: if ((*this)[idx].IsExcluded()) return kFALSE; @@ -746,45 +739,33 @@ return kFALSE; break; - case 35: + case 38: if ((*this)[idx].IsExcluded()) return kFALSE; val = (*this)[idx].GetPed(); break; - case 36: - if ((*this)[idx].IsExcluded()) - return kFALSE; - val = 1.; - // val = (*this)[idx].GetPedError(); - break; - case 37: + case 39: + if ((*this)[idx].IsExcluded()) + return kFALSE; + val = (*this)[idx].GetPedErr(); + break; + case 40: if ((*this)[idx].IsExcluded()) return kFALSE; val = (*this)[idx].GetPedRms(); break; - case 38: - if ((*this)[idx].IsExcluded()) - return kFALSE; - val = 1.; - // val = (*this)[idx].GetPedRmsError(); - break; - case 39: + case 41: + if ((*this)[idx].IsExcluded()) + return kFALSE; + val = (*this)[idx].GetPedErr()/2.; + break; + case 42: if ((*this)[idx].IsExcluded()) return kFALSE; val = (*this)[idx].GetAbsTimeMean(); break; - case 40: - if ((*this)[idx].IsExcluded()) - return kFALSE; - val = (*this)[idx].GetAbsTimeMeanErr(); - break; - case 41: + case 43: if ((*this)[idx].IsExcluded()) return kFALSE; val = (*this)[idx].GetAbsTimeRms(); - break; - case 42: - if ((*this)[idx].IsExcluded()) - return kFALSE; - val = (*this)[idx].GetAbsTimeMeanErr()/TMath::Sqrt(2.); break; default: @@ -893,6 +874,6 @@ TIter Next(fPixels); - MCalibrationPix *pix; - while ((pix=(MCalibrationPix*)Next())) + MCalibrationChargePix *pix; + while ((pix=(MCalibrationChargePix*)Next())) { @@ -902,7 +883,7 @@ const Int_t idx = pix->GetPixId(); - const Float_t charge = pix->GetCharge(); + const Float_t charge = pix->GetMeanCharge(); const Float_t area = (*fGeomCam)[idx].GetA(); - const Float_t chargeerr = pix->GetChargeErr(); + const Float_t chargeerr = pix->GetMeanChargeErr(); const Float_t nphot = fMeanFluxInsidePlexiglass*area; @@ -937,6 +918,6 @@ TIter Next(fPixels); - MCalibrationPix *pix; - while ((pix=(MCalibrationPix*)Next())) + MCalibrationChargePix *pix; + while ((pix=(MCalibrationChargePix*)Next())) { @@ -946,7 +927,7 @@ const Int_t idx = pix->GetPixId(); - const Float_t charge = pix->GetCharge(); + const Float_t charge = pix->GetMeanCharge(); const Float_t area = (*fGeomCam)[idx].GetA(); - const Float_t chargeerr = pix->GetChargeErr(); + const Float_t chargeerr = pix->GetMeanChargeErr(); const Float_t nphot = flux * area; @@ -980,6 +961,4 @@ { - if (ipx < 0 || !IsPixelFitted(ipx)) - return kFALSE; if (!IsFluxInsidePlexiglassAvailable()) @@ -997,7 +976,4 @@ Bool_t MCalibrationChargeCam::GetConversionFactorFFactor(Int_t ipx, Float_t &mean, Float_t &err, Float_t &sigma) { - - if (ipx < 0 || !IsPixelFitted(ipx)) - return kFALSE; Float_t conv = (*this)[ipx].GetMeanConversionFFactorMethod(); @@ -1024,7 +1000,4 @@ Bool_t MCalibrationChargeCam::GetConversionFactorPINDiode(Int_t ipx, Float_t &mean, Float_t &err, Float_t &sigma) { - - if (ipx < 0 || !IsPixelFitted(ipx)) - return kFALSE; mean = (*this)[ipx].GetMeanConversionPINDiodeMethod(); @@ -1048,8 +1021,4 @@ Bool_t MCalibrationChargeCam::GetConversionFactorCombined(Int_t ipx, Float_t &mean, Float_t &err, Float_t &sigma) { - - if (ipx < 0 || !IsPixelFitted(ipx)) - return kFALSE; - return kFALSE; Index: /trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCam.h =================================================================== --- /trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCam.h (revision 3263) +++ /trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCam.h (revision 3264) @@ -13,7 +13,7 @@ class TClonesArray; -class MCalibrationPix; class MCalibrationBlindPix; class MCalibrationChargePINDiode; +class MCalibrationChargePix; class MCalibrationChargeCam : public MParContainer, public MCamEvent { @@ -103,7 +103,4 @@ Bool_t GetConversionFactorCombined( Int_t ipx, Float_t &mean, Float_t &err, Float_t &sigma ); - Bool_t IsPixelUsed(Int_t idx) const; - Bool_t IsPixelFitted(Int_t idx) const; - Bool_t IsBlindPixelMethodValid() const; Bool_t IsPINDiodeMethodValid() const; @@ -112,6 +109,6 @@ // Others - MCalibrationPix &operator[](UInt_t i); - const MCalibrationPix &operator[](UInt_t i) const; + MCalibrationChargePix &operator[](UInt_t i); + const MCalibrationChargePix &operator[](UInt_t i) const; void CutEdges(); Index: /trunk/MagicSoft/Mars/mcalib/MCalibrationChargePINDiode.h =================================================================== --- /trunk/MagicSoft/Mars/mcalib/MCalibrationChargePINDiode.h (revision 3263) +++ /trunk/MagicSoft/Mars/mcalib/MCalibrationChargePINDiode.h (revision 3264) @@ -80,7 +80,4 @@ kFluxOutsidePlexiglassAvailable }; - Bool_t CheckChargeFitValidity(); - Bool_t CheckTimeFitValidity(); - public: @@ -91,5 +88,5 @@ // Setters - void SetColor( const PulserColor_t color ) { fColor = color; } + void SetColor( const PulserColor_t color ) { fColor = color; } void SetMeanCharge ( const Float_t f ) { fMeanCharge = f; } @@ -146,4 +143,6 @@ Bool_t IsMeanTimeInLastBin() const; + Bool_t CheckChargeFitValidity(); + Bool_t CheckTimeFitValidity(); Bool_t CalcFluxOutsidePlexiglass(); @@ -153,2 +152,10 @@ #endif /* MARS_MCalibrationChargePINDiode */ + + + + + + + + Index: /trunk/MagicSoft/Mars/mcalib/MCalibrationChargePix.cc =================================================================== --- /trunk/MagicSoft/Mars/mcalib/MCalibrationChargePix.cc (revision 3264) +++ /trunk/MagicSoft/Mars/mcalib/MCalibrationChargePix.cc (revision 3264) @@ -0,0 +1,806 @@ +/* ======================================================================== *\ +! +! * +! * This file is part of MARS, the MAGIC Analysis and Reconstruction +! * Software. It is distributed to you in the hope that it can be a useful +! * and timesaving tool in analysing Data of imaging Cerenkov telescopes. +! * It is distributed WITHOUT ANY WARRANTY. +! * +! * Permission to use, copy, modify and distribute this software and its +! * documentation for any purpose is hereby granted without fee, +! * provided that the above copyright notice appear in all copies and +! * that both that copyright notice and this permission notice appear +! * in supporting documentation. It is provided "as is" without express +! * or implied warranty. +! * +! +! +! Author(s): Markus Gaug 02/2004 +! +! Copyright: MAGIC Software Development, 2000-2004 +! +! +\* ======================================================================== */ + +///////////////////////////////////////////////////////////////////////////// +// // +// MCalibrationChargePix // +// // +// This is the storage container to hold informations about the pedestal // +// (offset) value of one Pixel (PMT). // +// // +// The following values are initialized to meaningful values: +// +// - The Electronic Rms to 1.5 per FADC slice +// - The uncertainty about the Electronic RMS to 0.3 per slice +// - The F-Factor is assumed to have been measured in Munich to 1.13 - 1.17. +// with the Munich definition of the F-Factor, thus: +// F = Sigma(Out)/Mean(Out) * Mean(In)/Sigma(In) +// Mean F-Factor = 1.15 +// Error F-Factor = 0.02 +// +// - Average QE: (email David Paneque, 14.2.04): +// +// The conversion factor that comes purely from QE folded to a Cherenkov +// spectrum has to be multiplied by: +// * Plexiglass window -->> 0.96 X 0.96 +// * PMT photoelectron collection efficiency -->> 0.9 +// * Light guides efficiency -->> 0.94 +// +// Concerning the light guides effiency estimation... Daniel Ferenc +// is preparing some work (simulations) to estimate it. Yet so far, he has +// been busy with other stuff, and this work is still UNfinished. +// +// The estimation I did comes from: +// 1) Reflectivity of light guide walls is 85 % (aluminum) +// 2) At ZERO degree light incidence, 37% of the light hits such walls +// (0.15X37%= 5.6% of light lost) +// 3) When increasing the light incidence angle, more and more light hits +// the walls. +// +// However, the loses due to larger amount of photons hitting the walls is more +// or less counteracted by the fact that more and more photon trajectories cross +// the PMT photocathode twice, increasing the effective sensitivity of the PMT. +// +// Jurgen Gebauer did some quick measurements about this issue. I attach a +// plot. You can see that the angular dependence is (more or less) in agreement +// with a CosTheta function (below 20-25 degrees), +// which is the variation of teh entrance window cross section. So, in +// first approximation, no loses when increasing light incidence angle; +// and therefore, the factor 0.94. +// +// So, summarizing... I would propose the following conversion factors +// (while working with CT1 cal box) in order to get the final number of photons +// from the detected measured size in ADC counts. +// +// Nph = ADC * FmethodConversionFactor / ConvPhe-PhFactor +// +// FmethodConversionFactor ; measured for individual pmts +// +// ConvPhe-PhFactor = 0.98 * 0.23 * 0.90 * 0.94 * 0.96 * 0.96 = 0.18 +// +// I would not apply any smearing of this factor (which we have in nature), +// since we might be applying it to PMTs in the totally wrong direction. +// +// +///////////////////////////////////////////////////////////////////////////// +#include "MCalibrationChargePix.h" + +#include "MLog.h" +#include "MLogManip.h" + +ClassImp(MCalibrationChargePix); + +using namespace std; + +const Float_t MCalibrationChargePix::gkElectronicPedRms = 1.5; +const Float_t MCalibrationChargePix::gkElectronicPedRmsErr = 0.3; +const Float_t MCalibrationChargePix::gkFFactor = 1.15; +const Float_t MCalibrationChargePix::gkFFactorErr = 0.02; + +const Float_t MCalibrationChargePix::gkAverageQE = 0.18; +const Float_t MCalibrationChargePix::gkAverageQEErr = 0.02; +const Float_t MCalibrationChargePix::gkConversionHiLo = 10.; +const Float_t MCalibrationChargePix::gkConversionHiLoErr = 2.5; + +const Float_t MCalibrationChargePix::fgChargeLimit = 3.; +const Float_t MCalibrationChargePix::fgChargeErrLimit = 0.; +const Float_t MCalibrationChargePix::fgChargeRelErrLimit = 1.; +const Float_t MCalibrationChargePix::fgTimeLimit = 1.5; +const Float_t MCalibrationChargePix::fgTimeErrLimit = 3.; +const Float_t MCalibrationChargePix::fgConvFFactorRelErrLimit = 0.1; +// -------------------------------------------------------------------------- +// +// Default Constructor: +// +MCalibrationChargePix::MCalibrationChargePix(const char *name, const char *title) + : fPixId(-1), + fFlags(0) +{ + + fName = name ? name : "MCalibrationChargePix"; + fTitle = title ? title : "Container of the fit results of MHCalibrationChargePixs "; + + Clear(); + + // + // At the moment, we don't have a database, yet, + // so we get it from the configuration file + // + SetConversionHiLo(); + SetConversionHiLoErr(); + + SetAverageQE(); + SetChargeLimit(); + SetChargeErrLimit(); + + SetChargeRelErrLimit(); + SetTimeLimit(); + SetTimeErrLimit(); + SetConvFFactorRelErrLimit(); +} + +// ------------------------------------------------------------------------ +// +// Invalidate values +// +void MCalibrationChargePix::Clear(Option_t *o) +{ + + SetHiGainSaturation ( kFALSE ); + SetLoGainSaturation ( kFALSE ); + SetHiGainFitted ( kFALSE ); + SetLoGainFitted ( kFALSE ); + SetHiGainOscillating ( kFALSE ); + SetLoGainOscillating ( kFALSE ); + SetExcluded ( kFALSE ); + SetChargeValid ( kFALSE ); + SetBlindPixelMethodValid ( kFALSE ); + SetFFactorMethodValid ( kFALSE ); + SetPINDiodeMethodValid ( kFALSE ); + SetCombinedMethodValid ( kFALSE ); + + SetMeanTimeInFirstBin ( kFALSE ); + SetMeanTimeInLastBin ( kFALSE ); + + fHiGainMeanCharge = -1.; + fHiGainMeanChargeErr = -1.; + fHiGainSigmaCharge = -1.; + fHiGainSigmaChargeErr = -1.; + fHiGainChargeProb = -1.; + + fLoGainMeanCharge = -1.; + fLoGainMeanChargeErr = -1.; + fLoGainSigmaCharge = -1.; + fLoGainSigmaChargeErr = -1.; + fLoGainChargeProb = -1.; + + fRSigmaCharge = -1.; + fRSigmaChargeErr = -1.; + + fHiGainNumPickup = -1; + fLoGainNumPickup = -1; + + fNumLoGainSamples = -1.; + + fPed = -1.; + fPedRms = -1.; + fPedErr = -1.; + + fLoGainPedRms = -1.; + fLoGainPedRmsErr = -1.; + + fTimeFirstHiGain = 0 ; + fTimeLastHiGain = 0 ; + fTimeFirstLoGain = 0 ; + fTimeLastLoGain = 0 ; + + fAbsTimeMean = -1.; + fAbsTimeRms = -1.; + + fPheFFactorMethod = -1.; + fPheFFactorMethodErr = -1.; + + fMeanConversionFFactorMethod = -1.; + fMeanConversionBlindPixelMethod = -1.; + fMeanConversionPINDiodeMethod = -1.; + fMeanConversionCombinedMethod = -1.; + + fConversionFFactorMethodErr = -1.; + fConversionBlindPixelMethodErr = -1.; + fConversionPINDiodeMethodErr = -1.; + fConversionCombinedMethodErr = -1.; + + fSigmaConversionFFactorMethod = -1.; + fSigmaConversionBlindPixelMethod = -1.; + fSigmaConversionPINDiodeMethod = -1.; + fSigmaConversionCombinedMethod = -1.; + + fTotalFFactorFFactorMethod = -1.; + fTotalFFactorBlindPixelMethod = -1.; + fTotalFFactorPINDiodeMethod = -1.; + fTotalFFactorCombinedMethod = -1.; + +} + + +void MCalibrationChargePix::DefinePixId(Int_t i) +{ + fPixId = i; +} + + +// -------------------------------------------------------------------------- +// +// Set the pedestals from outside +// +void MCalibrationChargePix::SetPedestal(const Float_t ped, const Float_t pedrms, const Float_t pederr) +{ + + fPed = ped; + fPedRms = pedrms; + fPedErr = pederr; +} + + +// -------------------------------------------------------------------------- +// +// Set the conversion factors from outside (only for MC) +// +void MCalibrationChargePix::SetConversionFFactorMethod(Float_t c, Float_t err, Float_t sig) +{ + fMeanConversionFFactorMethod = c; + fConversionFFactorMethodErr = err; + fSigmaConversionFFactorMethod = sig; +} + +// -------------------------------------------------------------------------- +// +// Set the conversion factors from outside (only for MC) +// +void MCalibrationChargePix::SetConversionCombinedMethod(Float_t c, Float_t err, Float_t sig) +{ + fMeanConversionCombinedMethod = c; + fConversionCombinedMethodErr = err; + fSigmaConversionCombinedMethod = sig; +} + + +// -------------------------------------------------------------------------- +// +// Set the conversion factors from outside (only for MC) +// +void MCalibrationChargePix::SetConversionBlindPixelMethod(Float_t c, Float_t err, Float_t sig) +{ + fMeanConversionBlindPixelMethod = c; + fConversionBlindPixelMethodErr = err; + fSigmaConversionBlindPixelMethod = sig; +} + +// -------------------------------------------------------------------------- +// +// Set the conversion factors from outside (only for MC) +// +void MCalibrationChargePix::SetConversionPINDiodeMethod(Float_t c, Float_t err, Float_t sig) +{ + fMeanConversionPINDiodeMethod = c ; + fConversionPINDiodeMethodErr = err; + fSigmaConversionPINDiodeMethod = sig; +} + +// -------------------------------------------------------------------------- +// +// Set the Hi Gain Saturation Bit from outside +// +void MCalibrationChargePix::SetHiGainSaturation(Bool_t b) +{ + b ? SETBIT(fFlags, kHiGainSaturation) : CLRBIT(fFlags, kHiGainSaturation); +} + +// -------------------------------------------------------------------------- +// +// Set the Lo Gain Saturation Bit from outside +// +void MCalibrationChargePix::SetLoGainSaturation(Bool_t b) +{ + b ? SETBIT(fFlags, kLoGainSaturation) : CLRBIT(fFlags, kLoGainSaturation); +} + +// -------------------------------------------------------------------------- +// +// Set the Excluded Bit from outside +// +void MCalibrationChargePix::SetExcluded(Bool_t b ) +{ + b ? SETBIT(fFlags, kExcluded) : CLRBIT(fFlags, kExcluded); +} + +// -------------------------------------------------------------------------- +// +// Set the Excluded Bit from outside +// +void MCalibrationChargePix::SetChargeValid(Bool_t b ) +{ + b ? SETBIT(fFlags, kChargeValid) : CLRBIT(fFlags, kChargeValid); +} + +// -------------------------------------------------------------------------- +// +// Set the Fitted Bit from outside +// +void MCalibrationChargePix::SetHiGainFitted(Bool_t b ) +{ + b ? SETBIT(fFlags, kHiGainFitted) : CLRBIT(fFlags, kHiGainFitted); +} + +// -------------------------------------------------------------------------- +// +// Set the Fitted Bit from outside +// +void MCalibrationChargePix::SetLoGainFitted(const Bool_t b ) +{ + b ? SETBIT(fFlags, kLoGainFitted) : CLRBIT(fFlags, kLoGainFitted); +} + +// -------------------------------------------------------------------------- +// +// Set the Ocillating Bit from outside +// +void MCalibrationChargePix::SetHiGainOscillating(const Bool_t b ) +{ + b ? SETBIT(fFlags, kHiGainOscillating) : CLRBIT(fFlags, kHiGainOscillating); +} + +// -------------------------------------------------------------------------- +// +// Set the Ocillating Bit from outside +// +void MCalibrationChargePix::SetLoGainOscillating(const Bool_t b ) +{ + b ? SETBIT(fFlags, kLoGainOscillating) : CLRBIT(fFlags, kLoGainOscillating); +} + +// -------------------------------------------------------------------------- +// +// Set the Excluded Bit from outside +// +void MCalibrationChargePix::SetBlindPixelMethodValid(const Bool_t b ) +{ + b ? SETBIT(fFlags, kBlindPixelMethodValid) : CLRBIT(fFlags, kBlindPixelMethodValid); +} + +// -------------------------------------------------------------------------- +// +// Set the Excluded Bit from outside +// +void MCalibrationChargePix::SetFFactorMethodValid(const Bool_t b ) +{ + b ? SETBIT(fFlags, kFFactorMethodValid) : CLRBIT(fFlags, kFFactorMethodValid); +} + +// -------------------------------------------------------------------------- +// +// Set the Excluded Bit from outside +// +void MCalibrationChargePix::SetPINDiodeMethodValid(const Bool_t b ) +{ + b ? SETBIT(fFlags, kPINDiodeMethodValid) : CLRBIT(fFlags, kPINDiodeMethodValid); +} + +// -------------------------------------------------------------------------- +// +// Set the Excluded Bit from outside +// +void MCalibrationChargePix::SetCombinedMethodValid(const Bool_t b ) +{ + b ? SETBIT(fFlags, kCombinedMethodValid) : CLRBIT(fFlags, kCombinedMethodValid); +} + +void MCalibrationChargePix::SetAbsTimeBordersHiGain(const Byte_t f, const Byte_t l) +{ + + fTimeFirstHiGain = f; + fTimeLastHiGain = l; + +} + +void MCalibrationChargePix::SetAbsTimeBordersLoGain(const Byte_t f, const Byte_t l) +{ + + fTimeFirstLoGain = f; + fTimeLastLoGain = l; + +} + +void MCalibrationChargePix::SetMeanTimeInFirstBin(const Bool_t b) +{ + b ? SETBIT(fFlags,kMeanTimeInFirstBin) : CLRBIT(fFlags,kMeanTimeInFirstBin); +} + +void MCalibrationChargePix::SetMeanTimeInLastBin(const Bool_t b) +{ + b ? SETBIT(fFlags,kMeanTimeInLastBin) : CLRBIT(fFlags,kMeanTimeInLastBin); +} + +Float_t MCalibrationChargePix::GetMeanCharge() const +{ + return IsHiGainSaturation() ? fLoGainMeanCharge : fHiGainMeanCharge ; +} + +Float_t MCalibrationChargePix::GetMeanChargeErr() const +{ + return IsHiGainSaturation() ? fLoGainMeanChargeErr : fHiGainMeanChargeErr ; +} + +Float_t MCalibrationChargePix::GetChargeProb() const +{ + return IsHiGainSaturation() ? fLoGainChargeProb : fHiGainChargeProb ; +} + +Float_t MCalibrationChargePix::GetSigmaCharge() const +{ + return IsHiGainSaturation() ? fLoGainSigmaCharge : fHiGainSigmaCharge ; +} + +Float_t MCalibrationChargePix::GetSigmaChargeErr() const +{ + return IsHiGainSaturation() ? fLoGainSigmaChargeErr : fHiGainSigmaChargeErr ; +} + +Bool_t MCalibrationChargePix::IsFitted() const +{ + return IsHiGainSaturation() ? IsLoGainFitted() : IsHiGainFitted(); +} + +Bool_t MCalibrationChargePix::IsOscillating() const +{ + return IsHiGainSaturation() ? IsLoGainOscillating() : IsHiGainOscillating(); +} + +Bool_t MCalibrationChargePix::IsMeanTimeInFirstBin() const +{ + return TESTBIT(fFlags,kMeanTimeInFirstBin); +} + +Bool_t MCalibrationChargePix::IsMeanTimeInLastBin() const +{ + return TESTBIT(fFlags,kMeanTimeInLastBin); +} + + + +Bool_t MCalibrationChargePix::IsExcluded() const +{ + return TESTBIT(fFlags,kExcluded); +} + +Bool_t MCalibrationChargePix::IsHiGainSaturation() const +{ + return TESTBIT(fFlags,kHiGainSaturation); +} + +Bool_t MCalibrationChargePix::IsLoGainSaturation() const +{ + return TESTBIT(fFlags,kLoGainSaturation); +} + +Bool_t MCalibrationChargePix::IsChargeValid() const +{ + return TESTBIT(fFlags, kChargeValid); +} + +Bool_t MCalibrationChargePix::IsHiGainFitted() const +{ + return TESTBIT(fFlags, kHiGainFitted); +} + +Bool_t MCalibrationChargePix::IsLoGainFitted() const +{ + return TESTBIT(fFlags, kLoGainFitted); +} + +Bool_t MCalibrationChargePix::IsHiGainOscillating() const +{ + return TESTBIT(fFlags, kHiGainOscillating); +} + +Bool_t MCalibrationChargePix::IsLoGainOscillating() const +{ + return TESTBIT(fFlags, kLoGainOscillating); +} + +Bool_t MCalibrationChargePix::IsBlindPixelMethodValid() const +{ + return TESTBIT(fFlags, kBlindPixelMethodValid); +} + +Bool_t MCalibrationChargePix::IsFFactorMethodValid() const +{ + return TESTBIT(fFlags, kFFactorMethodValid); +} + +Bool_t MCalibrationChargePix::IsPINDiodeMethodValid() const +{ + return TESTBIT(fFlags, kPINDiodeMethodValid); +} + +Bool_t MCalibrationChargePix::IsCombinedMethodValid() const +{ + return TESTBIT(fFlags, kCombinedMethodValid); +} + + +// +// The check return kTRUE if: +// +// 1) Pixel has a fitted charge greater than fChargeLimit*PedRMS +// 2) Pixel has a fit error greater than fChargeErrLimit +// 3) Pixel has a fitted charge greater its fChargeRelErrLimit times its charge error +// 4) Pixel has a charge sigma bigger than its Pedestal RMS +// +Bool_t MCalibrationChargePix::CheckChargeValidity() +{ + + if (!IsFitted()) + return kFALSE; + + if (GetMeanCharge() < fChargeLimit*GetPedRms()) + { + *fLog << warn << "WARNING: Fitted Charge is smaller than " + << fChargeLimit << " Pedestal RMS in Pixel " << fPixId << endl; + return kFALSE; + } + + if (GetMeanChargeErr() < fChargeErrLimit) + { + *fLog << warn << "WARNING: Error of Fitted Charge is smaller than " + << fChargeErrLimit << " in Pixel " << fPixId << endl; + return kFALSE; + } + + if (GetMeanCharge() < fChargeRelErrLimit*GetMeanChargeErr()) + { + *fLog << warn << "WARNING: Fitted Charge is smaller than " + << fChargeRelErrLimit << "* its error in Pixel " << fPixId << endl; + return kFALSE; + } + + if (GetSigmaCharge() < GetPedRms()) + { + *fLog << warn << "WARNING: Sigma of Fitted Charge smaller than Pedestal RMS in Pixel " << fPixId << endl; + return kFALSE; + } + + SetChargeValid(); + return kTRUE; +} + +// +// The check returns kTRUE if: +// +// The mean arrival time is at least 1.0 slices from the used edge slices +// +Bool_t MCalibrationChargePix::CheckTimeValidity() +{ + + const Float_t loweredge = IsHiGainSaturation() ? fTimeFirstHiGain : fTimeFirstLoGain; + const Float_t upperedge = IsHiGainSaturation() ? fTimeLastHiGain : fTimeLastLoGain; + + if ( fAbsTimeMean < loweredge+1.) + { + *fLog << warn << "WARNING: Mean ArrivalTime in first extraction bin of the Pixel " << fPixId << endl; + SetMeanTimeInFirstBin(); + return kFALSE; + } + + if ( fAbsTimeMean > upperedge-1.) + { + *fLog << warn << "WARNING: Mean ArrivalTime in last extraction bin of the Pixel " << fPixId << endl; + SetMeanTimeInLastBin(); + return kFALSE; + } + + return kTRUE; +} + +void MCalibrationChargePix::CalcLoGainPed() +{ + + Float_t pedRmsSquare = fPedRms * fPedRms; + Float_t pedRmsSquareErrSquare = fPedErr * fPedErr * pedRmsSquare; // fPedRmsErr = fPedErr/2. + Float_t pedRmsSquareErr; + + // + // We do not know the Lo Gain Pedestal RMS, so we have to retrieve it + // from the HI GAIN (all calculation per slice up to now): + // + // We extract the pure NSB contribution: + // + const Float_t elecRmsSquare = fElectronicPedRms * fElectronicPedRms; + const Float_t elecRmsSquareErrSquare = 4.*fElectronicPedRmsErr * fElectronicPedRmsErr * elecRmsSquare; + + Float_t nsbSquare = pedRmsSquare - elecRmsSquare; + Float_t nsbSquareRelErrSquare = (pedRmsSquareErrSquare + elecRmsSquareErrSquare) + / (nsbSquare * nsbSquare) ; + + if (nsbSquare < 0.) + nsbSquare = 0.; + + // + // Now, we divide the NSB by the conversion factor and + // add it quadratically to the electronic noise + // + const Float_t conversionSquare = fConversionHiLo * fConversionHiLo; + const Float_t conversionSquareRelErrSquare = 4.*fConversionHiLoErr * fConversionHiLoErr / conversionSquare; + + const Float_t convertedNsbSquare = nsbSquare / conversionSquare; + const Float_t convertedNsbSquareErrSquare = (nsbSquareRelErrSquare + conversionSquareRelErrSquare) + * convertedNsbSquare * convertedNsbSquare; + + pedRmsSquare = convertedNsbSquare + elecRmsSquare; + pedRmsSquareErr = TMath::Sqrt( convertedNsbSquareErrSquare + elecRmsSquareErrSquare ); + + // + // Correct also for the conversion to Hi-Gain: + // + fLoGainPedRms = TMath::Sqrt(pedRmsSquare) * fConversionHiLo; + fLoGainPedRmsErr = 0.5 * pedRmsSquareErr / fLoGainPedRms * fConversionHiLo; +} + +Bool_t MCalibrationChargePix::CalcReducedSigma() +{ + + const Float_t sigmaSquare = GetSigmaCharge() * GetSigmaCharge(); + const Float_t sigmaSquareErrSquare = 4.*GetSigmaChargeErr()* GetSigmaChargeErr() * sigmaSquare; + + Float_t pedRmsSquare; + Float_t pedRmsSquareErrSquare; + + if (IsHiGainSaturation()) + { + CalcLoGainPed(); + + pedRmsSquare = fLoGainPedRms * fLoGainPedRms; + pedRmsSquareErrSquare = 4.* fLoGainPedRmsErr * fLoGainPedRmsErr * pedRmsSquare; + } /* if (HiGainSaturation) */ + else + { + + pedRmsSquare = fPedRms * fPedRms; + pedRmsSquareErrSquare = fPedErr * fPedErr * pedRmsSquare; // fPedRmsErr = fPedErr/2. + } + // + // Calculate the reduced sigmas + // + const Float_t rsigmachargesquare = sigmaSquare - pedRmsSquare; + + if (rsigmachargesquare <= 0.) + { + *fLog << warn + << "WARNING: Cannot calculate the reduced sigma: smaller than 0 in pixel " + << fPixId << endl; + return kFALSE; + } + + fRSigmaCharge = TMath::Sqrt(rsigmachargesquare); + fRSigmaChargeErr = TMath::Sqrt(sigmaSquareErrSquare + pedRmsSquareErrSquare); + + return kTRUE; +} + +// +// Calculate the number of photo-electrons after the F-Factor method +// Calculate the errors of the F-Factor method +// +Bool_t MCalibrationChargePix::CalcFFactorMethod() +{ + + + if (fRSigmaCharge < 0.) + return kFALSE; + + // + // Square all variables in order to avoid applications of square root + // + // First the relative error squares + // + const Float_t chargeSquare = GetMeanCharge() * GetMeanCharge(); + const Float_t chargeSquareRelErrSquare = 4.* GetMeanChargeErr() * GetMeanChargeErr() / chargeSquare; + + const Float_t chargeRelErrSquare = GetMeanChargeErr() * GetMeanChargeErr() + / (GetMeanCharge() * GetMeanCharge()); + + const Float_t ffactorsquare = gkFFactor * gkFFactor; + const Float_t ffactorsquareRelErrSquare = 4.*gkFFactorErr * gkFFactorErr / ffactorsquare; + + const Float_t avQERelErrSquare = fAverageQEErr * fAverageQEErr / fAverageQE / fAverageQE; + + const Float_t avQEFFactor = TMath::Sqrt( ( 1. - fAverageQE ) / fAverageQE ); + const Float_t avQEFFactorErr = 1./ ( 2. * avQEFFactor ) * fAverageQEErr + / ( fAverageQE * fAverageQE ); + const Float_t avQEFFactorRelErrSquare = avQEFFactorErr * avQEFFactorErr + / ( avQEFFactor * avQEFFactor) ; + + const Float_t rsigmaSquare = fRSigmaCharge * fRSigmaCharge; + const Float_t rsigmaSquareRelErrSquare = 4.* fRSigmaChargeErr * fRSigmaChargeErr / rsigmaSquare; + + // + // Calculate the number of phe's from the F-Factor method + // (independent on Hi Gain or Lo Gain) + // + fPheFFactorMethod = ffactorsquare * chargeSquare / rsigmaSquare; + + // + // Calculate the number of photons from the F-Factor method + // FIXME: This is a preliminary solution, the qe shall be + // calibrated itself! + // + fPheFFactorMethod /= fAverageQE; + + const Float_t pheFFactorRelErrSquare = ffactorsquareRelErrSquare + + chargeSquareRelErrSquare + + rsigmaSquareRelErrSquare + + avQERelErrSquare; + + fPheFFactorMethodErr = TMath::Sqrt(pheFFactorRelErrSquare) * fPheFFactorMethod; + + fMeanConversionFFactorMethod = fPheFFactorMethod / GetMeanCharge() ; + fConversionFFactorMethodErr = ( pheFFactorRelErrSquare + chargeRelErrSquare ) + * fMeanConversionFFactorMethod * fMeanConversionFFactorMethod; + + const Float_t convrelerror = fConversionFFactorMethodErr + / fMeanConversionFFactorMethod; + + if ( (fMeanConversionFFactorMethod > 0.) && (convrelerror < fConvFFactorRelErrLimit)) + SETBIT(fFlags,kFFactorMethodValid); + + fSigmaConversionFFactorMethod = GetTotalFFactorFFactorMethod()*TMath::Sqrt(fMeanConversionFFactorMethod); + + // + // Calculate the Total F-Factor of the camera ( in photons ) + // + if (fPheFFactorMethod > 0) + { + fTotalFFactorFFactorMethod = (fRSigmaCharge/GetMeanCharge())*TMath::Sqrt(fPheFFactorMethod); + fTotalFFactorFFactorMethod *= avQEFFactor; + } + + // + // Calculate the error of the Total F-Factor of the camera ( in photons ) + // + const Float_t rsigmaChargeRelErrSquare = fRSigmaChargeErr * fRSigmaChargeErr + / (fRSigmaCharge * fRSigmaCharge) ; + + fTotalFFactorErrFFactorMethod = TMath::Sqrt( rsigmaChargeRelErrSquare + + chargeRelErrSquare + + pheFFactorRelErrSquare + + avQEFFactorRelErrSquare ); + + fTotalFFactorErrFFactorMethod *= fTotalFFactorFFactorMethod; + + SetFFactorMethodValid(); + return kTRUE; +} + + + +void MCalibrationChargePix::ApplyLoGainConversion() +{ + + const Float_t chargeRelErrSquare = fLoGainMeanChargeErr * fLoGainMeanChargeErr + /( fLoGainMeanCharge * fLoGainMeanCharge ); + const Float_t sigmaRelErrSquare = fLoGainSigmaChargeErr * fLoGainSigmaChargeErr + /( fLoGainSigmaCharge * fLoGainSigmaCharge ); + const Float_t conversionRelErrSquare = fConversionHiLoErr * fConversionHiLoErr + /( fConversionHiLo * fConversionHiLo ); + + fLoGainMeanCharge *= fConversionHiLo; + fLoGainMeanChargeErr = TMath::Sqrt(chargeRelErrSquare + conversionRelErrSquare) * fLoGainMeanCharge; + + fLoGainSigmaCharge *= fConversionHiLo; + fLoGainSigmaChargeErr = TMath::Sqrt(sigmaRelErrSquare + conversionRelErrSquare) * fLoGainSigmaCharge; + + fElectronicPedRms = gkElectronicPedRms * TMath::Sqrt(fNumLoGainSamples); + fElectronicPedRmsErr = gkElectronicPedRmsErr * TMath::Sqrt(fNumLoGainSamples); + +} Index: /trunk/MagicSoft/Mars/mcalib/MCalibrationChargePix.h =================================================================== --- /trunk/MagicSoft/Mars/mcalib/MCalibrationChargePix.h (revision 3264) +++ /trunk/MagicSoft/Mars/mcalib/MCalibrationChargePix.h (revision 3264) @@ -0,0 +1,302 @@ +#ifndef MARS_MCalibrationChargePix +#define MARS_MCalibrationChargePix + +#ifndef MARS_MParContainer +#include "MParContainer.h" +#endif + +class MCalibrationChargePix : public MParContainer +{ +private: + + static const Float_t gkAverageQE; // The average quantum efficieny agreed on for the first analysis + static const Float_t gkAverageQEErr; // The error of average quantum efficieny + + static const Float_t gkConversionHiLo; // The default conversion factor HI gain - Lo Gain + static const Float_t gkConversionHiLoErr; // The error of the default conversion factor + + static const Float_t gkElectronicPedRms; // The pure electronic component of the RMS + static const Float_t gkElectronicPedRmsErr; // The error of the pure electronic component of the RMS + static const Float_t gkFFactor; // The laboratory F-factor of the PMTs + static const Float_t gkFFactorErr; // The laboratory F-factor Error of the PMTs + + static const Float_t fgChargeLimit; // The default limit (in units of PedRMS) for acceptance of the fitted mean charge + static const Float_t fgChargeErrLimit; // The default limit (in units of PedRMS) for acceptance of the fitted charge sigma + static const Float_t fgChargeRelErrLimit; // The default limit (in units of Error of fitted charge) for acceptance of the fitted mean + + static const Float_t fgTimeLimit; // The limit (in units of FADC slices) for acceptance of the fitted time + static const Float_t fgTimeErrLimit; // The limit (in units of FADC slices) for acceptance of the fitted time sigma + static const Float_t fgConvFFactorRelErrLimit; // The limit for acceptance of the rel. error of the conversion factor with the FFactor method + + Float_t fChargeLimit; // The limit (in units of PedRMS) for acceptance of the fitted mean charge + Float_t fChargeErrLimit; // The limit (in units of PedRMS) for acceptance of the fitted charge sigma + Float_t fChargeRelErrLimit; // The limit (in units of Error of fitted charge) for acceptance of the fitted mean + + Float_t fTimeLimit; // The limit (in units of FADC slices) for acceptance of the fitted time + Float_t fTimeErrLimit; // The limit (in units of FADC slices) for acceptance of the fitted time sigma + Float_t fConvFFactorRelErrLimit; // The limit for acceptance of the rel. error of the conversion factor with the FFactor method + + Float_t fElectronicPedRms; // The pure electronic component of the RMS + Float_t fElectronicPedRmsErr; // The error of the pure electronic component of the RMS + + Int_t fPixId; // the pixel Id + + UInt_t fFlags; // Flag for the set bits + + Float_t fAverageQE; // The average quantum efficieny (see Class description) + Float_t fAverageQEErr; // The error of the average quantum efficieny (see Class description) + + Float_t fHiGainMeanCharge; // The mean reduced charge after the fit + Float_t fHiGainMeanChargeErr; // The error of reduced mean charge after the fit + Float_t fHiGainSigmaCharge; // The sigma of the mean charge after the fit + Float_t fHiGainSigmaChargeErr; // The error of the sigma of the mean charge after the fit + Float_t fHiGainChargeProb; // The probability of the fit function + + Float_t fLoGainMeanCharge; // The mean reduced charge after the fit + Float_t fLoGainMeanChargeErr; // The error of reduced mean charge after the fit + Float_t fLoGainSigmaCharge; // The sigma of the mean charge after the fit + Float_t fLoGainSigmaChargeErr; // The error of the sigma of the mean charge after the fit + Float_t fLoGainChargeProb; // The probability of the fit function + + Float_t fRSigmaCharge; // The reduced squares of sigmas after the fit + Float_t fRSigmaChargeErr; // The reduced squares of sigmas after the fit + + Float_t fPed; // The mean pedestal (from MPedestalPix) times number of FADC slices + Float_t fPedErr; // The error of the pedestal + Float_t fPedRms; // The pedestal RMS (from MPedestalPix) times sqrt of number of FADC slices + + Float_t fLoGainPedRms; // The pedestal RMS of the low gain + Float_t fLoGainPedRmsErr; // The pedestal RMS Error of the low gain + + Float_t fAbsTimeMean; // The mean absolute arrival time + Float_t fAbsTimeRms; // The rms of the mean absolute arrival time + + Byte_t fTimeFirstHiGain; // The first used FADC slice + Byte_t fTimeLastHiGain; // The last used FADC slice + + Byte_t fTimeFirstLoGain; // The first used FADC slice + Byte_t fTimeLastLoGain; // The last used FADC slice + + Float_t fPheFFactorMethod; // The number of Phe's calculated (F-factor method) + Float_t fPheFFactorMethodErr; // The error on the number of Phe's calculated (F-factor method) + + Float_t fMeanConversionFFactorMethod; // The conversion factor to Phe's (F-factor method) + Float_t fMeanConversionBlindPixelMethod; // The conversion factor to Ph's (Blind Pixel method) + Float_t fMeanConversionPINDiodeMethod; // The conversion factor to Ph's (PIN Diode method) + Float_t fMeanConversionCombinedMethod; // The conversion factor to Ph's (all methods combined) + + Float_t fConversionFFactorMethodErr; // The error of the conversion factor to Phe's (F-factor method) + Float_t fConversionBlindPixelMethodErr; // The error of the conversion factor to Ph's (Blind Pixel method) + Float_t fConversionPINDiodeMethodErr; // The error of the conversion factor to Ph's (PIN Diode method) + Float_t fConversionCombinedMethodErr; // The error of the conversion factor to Ph's (all methods combined) + + Float_t fSigmaConversionFFactorMethod; // The sigma of conversion factor to Phe's (F-factor method) + Float_t fSigmaConversionBlindPixelMethod; // The conversion factor to Ph's (Blind Pixel method) + Float_t fSigmaConversionPINDiodeMethod; // The conversion factor to Ph's (PIN Diode method) + Float_t fSigmaConversionCombinedMethod; // The conversion factor to Ph's (all methods combined) + + Float_t fTotalFFactorFFactorMethod; // The total F-Factor to Ph's (F-factor method) + Float_t fTotalFFactorBlindPixelMethod; // The total F-Factor to Ph's (Blind Pixel method) + Float_t fTotalFFactorPINDiodeMethod; // The total F-Factor to Ph's (PIN Diode method) + Float_t fTotalFFactorCombinedMethod; // The total F-Factor to Ph's (all methods combined) + + Float_t fTotalFFactorErrFFactorMethod; // The error of the total F-Factor to Ph's (F-factor method) + Float_t fTotalFFactorErrBlindPixelMethod; // The error of the total F-Factor to Ph's (Blind Pixel method) + Float_t fTotalFFactorErrPINDiodeMethod; // The error of the total F-Factor to Ph's (PIN Diode method) + Float_t fTotalFFactorErrCombinedMethod; // The error of the total F-Factor to Ph's (all methods combined) + + Float_t fTotalFFactor; // The F-Factor of the total readout system (Sigma(out)/mean(out)*Mean(in)/sigma(in) + Float_t fTotalFFactorErr; // The error on the F-Factor of the total readout system + + Float_t fConversionHiLo; // The conversion factor between Hi Gain and Lo Gain + Float_t fConversionHiLoErr; // The error of the conversion factor between Hi Gain and Lo Gain + + Float_t fNumLoGainSamples; + + Int_t fHiGainNumPickup; + Int_t fLoGainNumPickup; + + enum { kHiGainSaturation, kLoGainSaturation, + kExcluded, + kChargeValid, kTimeFitValid, + kHiGainFitted, kLoGainFitted, + kHiGainOscillating, kLoGainOscillating, + kMeanTimeInFirstBin, kMeanTimeInLastBin, + kBlindPixelMethodValid, kFFactorMethodValid, + kPINDiodeMethodValid, kCombinedMethodValid }; + + void CalcLoGainPed(); + +public: + + MCalibrationChargePix(const char *name=NULL, const char *title=NULL); + ~MCalibrationChargePix() {} + + void Clear(Option_t *o=""); + + // Setter + void SetPedestal(const Float_t ped, const Float_t pedrms, const Float_t pederr); + + void SetConversionHiLo( const Float_t c = gkConversionHiLo) { fConversionHiLo = c; } + void SetConversionHiLoErr( const Float_t e = gkConversionHiLoErr) { fConversionHiLoErr = e; } + void SetAverageQE( const Float_t qe= gkAverageQE, + const Float_t err=gkAverageQEErr) { fAverageQE = qe; + fAverageQEErr = err; } + void SetChargeLimit ( const Float_t f=fgChargeLimit ) { fChargeLimit = f; } + void SetChargeErrLimit ( const Float_t f=fgChargeErrLimit ) { fChargeErrLimit = f; } + void SetChargeRelErrLimit( const Float_t f=fgChargeRelErrLimit ) { fChargeRelErrLimit = f; } + + void SetTimeLimit ( const Float_t f=fgTimeLimit ) { fTimeLimit = f; } + void SetTimeErrLimit ( const Float_t f=fgTimeErrLimit ) { fTimeErrLimit = f; } + void SetConvFFactorRelErrLimit( const Float_t f=fgConvFFactorRelErrLimit ) { fConvFFactorRelErrLimit = f; } + + void SetMeanTimeInFirstBin( const Bool_t b = kTRUE ); + void SetMeanTimeInLastBin ( const Bool_t b = kTRUE ); + + // Charges + void SetHiGainMeanCharge ( const Float_t f ) { fHiGainMeanCharge = f; } + void SetHiGainMeanChargeErr ( const Float_t f ) { fHiGainMeanChargeErr = f; } + void SetHiGainChargeProb ( const Float_t f ) { fHiGainChargeProb = f; } + void SetHiGainSigmaCharge ( const Float_t f ) { fHiGainSigmaCharge = f; } + void SetHiGainSigmaChargeErr ( const Float_t f ) { fHiGainSigmaChargeErr = f; } + + void SetLoGainMeanCharge ( const Float_t f ) { fLoGainMeanCharge = f; } + void SetLoGainMeanChargeErr ( const Float_t f ) { fLoGainMeanChargeErr = f; } + void SetLoGainChargeProb ( const Float_t f ) { fLoGainChargeProb = f; } + void SetLoGainSigmaCharge ( const Float_t f ) { fLoGainSigmaCharge = f; } + void SetLoGainSigmaChargeErr ( const Float_t f ) { fLoGainSigmaChargeErr = f; } + + void SetHiGainNumPickup ( const Int_t i ) { fHiGainNumPickup = i; } + void SetLoGainNumPickup ( const Int_t i ) { fLoGainNumPickup = i; } + + // Times + void SetAbsTimeMean ( const Float_t f ) { fAbsTimeMean = f; } + void SetAbsTimeRms ( const Float_t f ) { fAbsTimeRms = f; } + + void SetNumLoGainSamples ( const Float_t f ) { fNumLoGainSamples = f; } + + // Setters for MC + void SetConversionFFactorMethod ( Float_t c, Float_t err, Float_t sig ); + void SetConversionBlindPixelMethod( Float_t c, Float_t err, Float_t sig ); + void SetConversionPINDiodeMethod ( Float_t c, Float_t err, Float_t sig ); + void SetConversionCombinedMethod ( Float_t c, Float_t err, Float_t sig ); + + // Bit Setters + void SetHiGainSaturation ( const Bool_t b = kTRUE ); + void SetLoGainSaturation ( const Bool_t b = kTRUE ); + void SetExcluded ( const Bool_t b = kTRUE ); + void SetChargeValid ( const Bool_t b = kTRUE ); + void SetHiGainFitted ( const Bool_t b = kTRUE ); + void SetLoGainFitted ( const Bool_t b = kTRUE ); + void SetHiGainOscillating ( const Bool_t b = kTRUE ); + void SetLoGainOscillating ( const Bool_t b = kTRUE ); + void SetBlindPixelMethodValid( const Bool_t b = kTRUE ); + void SetFFactorMethodValid ( const Bool_t b = kTRUE ); + void SetPINDiodeMethodValid ( const Bool_t b = kTRUE ); + void SetCombinedMethodValid ( const Bool_t b = kTRUE ); + + void SetAbsTimeBordersHiGain( Byte_t f, Byte_t l ); + void SetAbsTimeBordersLoGain( Byte_t f, Byte_t l ); + + // Charges + Float_t GetHiGainMeanCharge() const { return fHiGainMeanCharge; } + Float_t GetHiGainMeanChargeErr() const { return fHiGainMeanChargeErr; } + Float_t GetHiGainChargeProb() const { return fHiGainChargeProb; } + Float_t GetHiGainSigmaCharge() const { return fHiGainSigmaCharge; } + Float_t GetHiGainSigmaChargeErr() const { return fHiGainSigmaChargeErr; } + + Float_t GetLoGainMeanCharge() const { return fLoGainMeanCharge; } + Float_t GetLoGainMeanChargeErr() const { return fLoGainMeanChargeErr; } + Float_t GetLoGainChargeProb() const { return fLoGainChargeProb; } + Float_t GetLoGainSigmaCharge() const { return fLoGainSigmaCharge; } + Float_t GetLoGainSigmaChargeErr() const { return fLoGainSigmaChargeErr; } + + Float_t GetMeanCharge() const; + Float_t GetMeanChargeErr() const; + Float_t GetChargeProb() const; + Float_t GetSigmaCharge() const; + Float_t GetSigmaChargeErr() const; + + Float_t GetRSigmaCharge() const { return fRSigmaCharge; } + Float_t GetRSigmaChargeErr() const { return fRSigmaChargeErr; } + + Int_t GetHiGainNumPickup() const { return fHiGainNumPickup; } + Int_t GetLoGainNumPickup() const { return fLoGainNumPickup; } + + Float_t GetAbsTimeMean() const { return fAbsTimeMean; } + Float_t GetAbsTimeRms() const { return fAbsTimeRms; } + + // Conversion Factors + Float_t GetConversionHiLo() const { return fConversionHiLo; } + Float_t GetConversionHiLoErr() const { return fConversionHiLoErr; } + + Float_t GetMeanConversionBlindPixelMethod() const { return fMeanConversionBlindPixelMethod ; } + Float_t GetConversionBlindPixelMethodErr() const { return fConversionBlindPixelMethodErr ; } + Float_t GetSigmaConversionBlindPixelMethod() const { return fSigmaConversionBlindPixelMethod ; } + + Float_t GetMeanConversionFFactorMethod() const { return fMeanConversionFFactorMethod; } + Float_t GetConversionFFactorMethodErr() const { return fConversionFFactorMethodErr; } + Float_t GetSigmaConversionFFactorMethod() const { return fSigmaConversionFFactorMethod; } + + Float_t GetMeanConversionPINDiodeMethod() const { return fMeanConversionPINDiodeMethod ; } + Float_t GetConversionPINDiodeMethodErr() const { return fConversionPINDiodeMethodErr ; } + Float_t GetSigmaConversionPINDiodeMethod() const { return fSigmaConversionPINDiodeMethod ; } + + Float_t GetMeanConversionCombinedMethod() const { return fMeanConversionCombinedMethod ; } + Float_t GetConversionCombinedMethodErr() const { return fConversionCombinedMethodErr ; } + Float_t GetSigmaConversionCombinedMethod() const { return fSigmaConversionCombinedMethod ; } + + Float_t GetPheFFactorMethod() const { return fPheFFactorMethod; } + Float_t GetPheFFactorMethodErr() const { return fPheFFactorMethodErr; } + + Int_t GetPixId() const { return fPixId; } + + Float_t GetPed() const { return fPed; } + Float_t GetPedErr() const { return fPedErr; } + Float_t GetPedRms() const { return fPedRms; } + + Float_t GetTotalFFactorFFactorMethod() const { return fTotalFFactorFFactorMethod; } + Float_t GetTotalFFactorErrFFactorMethod() const { return fTotalFFactorErrFFactorMethod; } + + Float_t GetTotalFFactorBlindPixelMethod() const { return fTotalFFactorBlindPixelMethod; } + Float_t GetTotalFFactorErrBlindPixelMethod() const { return fTotalFFactorErrBlindPixelMethod; } + + Float_t GetTotalFFactorPINDiodeMethod() const { return fTotalFFactorPINDiodeMethod; } + Float_t GetTotalFFactorErrPINDiodeMethod() const { return fTotalFFactorErrPINDiodeMethod; } + + Float_t GetTotalFFactorCombinedMethod() const { return fTotalFFactorCombinedMethod; } + Float_t GetTotalFFactorErrCombinedMethod() const { return fTotalFFactorErrCombinedMethod; } + + Bool_t IsExcluded() const; + Bool_t IsHiGainSaturation() const; + Bool_t IsLoGainSaturation() const; + Bool_t IsChargeValid() const; + Bool_t IsHiGainFitted() const; + Bool_t IsLoGainFitted() const; + Bool_t IsFitted() const; + Bool_t IsHiGainOscillating() const; + Bool_t IsLoGainOscillating() const; + Bool_t IsOscillating() const; + Bool_t IsBlindPixelMethodValid() const; + Bool_t IsPINDiodeMethodValid() const; + Bool_t IsFFactorMethodValid() const; + Bool_t IsCombinedMethodValid() const; + + Bool_t IsMeanTimeInFirstBin() const; + Bool_t IsMeanTimeInLastBin() const; + + void DefinePixId(Int_t i); + + // Miscellaneous + void ApplyLoGainConversion(); + + Bool_t CheckChargeValidity(); + Bool_t CheckTimeValidity(); + Bool_t CalcReducedSigma(); + Bool_t CalcFFactorMethod(); + + ClassDef(MCalibrationChargePix, 1) // Container for Calibration of one pixel +}; + +#endif + Index: /trunk/MagicSoft/Mars/mcalib/MHCalibrationChargeHiGainPix.cc =================================================================== --- /trunk/MagicSoft/Mars/mcalib/MHCalibrationChargeHiGainPix.cc (revision 3263) +++ /trunk/MagicSoft/Mars/mcalib/MHCalibrationChargeHiGainPix.cc (revision 3264) @@ -49,10 +49,10 @@ using namespace std; -const Int_t MHCalibrationChargeHiGainPix::fgChargeNbins = 200; +const Int_t MHCalibrationChargeHiGainPix::fgChargeNbins = 2000; const Axis_t MHCalibrationChargeHiGainPix::fgChargeFirst = -0.5; -const Axis_t MHCalibrationChargeHiGainPix::fgChargeLast = 199.5; -const Int_t MHCalibrationChargeHiGainPix::fgAbsTimeNbins = 30; +const Axis_t MHCalibrationChargeHiGainPix::fgChargeLast = 1999.5; +const Int_t MHCalibrationChargeHiGainPix::fgAbsTimeNbins = 20; const Axis_t MHCalibrationChargeHiGainPix::fgAbsTimeFirst = -0.5; -const Axis_t MHCalibrationChargeHiGainPix::fgAbsTimeLast = 29.5; +const Axis_t MHCalibrationChargeHiGainPix::fgAbsTimeLast = 19.5; // -------------------------------------------------------------------------- // @@ -60,9 +60,8 @@ // MHCalibrationChargeHiGainPix::MHCalibrationChargeHiGainPix(const char *name, const char *title) - : fPix(NULL) { fName = name ? name : "MHCalibrationChargeHiGainPix"; - fTitle = title ? title : "Fill the FADC sums of the HiGainPix events and perform the fits"; + fTitle = title ? title : "Fill the FADC sums of the HiGainPix events and perform the fits Pixel "; SetChargeNbins(); @@ -73,5 +72,4 @@ SetAbsTimeFirst(); SetAbsTimeLast(); - } @@ -97,14 +95,44 @@ +void MHCalibrationChargeHiGainPix::ChangeHistId(Int_t id) +{ + + fPixId = id; + + fHGausHist.SetName(Form("%s%d", fHGausHist.GetName(), id)); + fHGausHist.SetTitle(Form("%s%d", fHGausHist.GetTitle(), id)); + + fHAbsTime.SetName(Form("%s%d", fHAbsTime.GetName(), id)); + fHAbsTime.SetTitle(Form("%s%d", fHAbsTime.GetTitle(), id)); + + fName = Form("%s%d", fName.Data(), id); + fTitle = Form("%s%d", fTitle.Data(), id); +} + + +void MHCalibrationChargeHiGainPix::Draw(const Option_t *opt) +{ + + TString option(opt); + option.ToLower(); + + Int_t win = 1; + + TVirtualPad *pad = gPad ? gPad : MH::MakeDefCanvas(this,600, 600); + + pad->SetTicks(); + pad->SetBorderMode(0); + + + pad->Divide(1,win); + pad->cd(1); + + if (!IsEmpty()) + pad->SetLogy(); + + MHGausEvents::Draw(opt); + +} - - - - - - - - - Index: /trunk/MagicSoft/Mars/mcalib/MHCalibrationChargeHiGainPix.h =================================================================== --- /trunk/MagicSoft/Mars/mcalib/MHCalibrationChargeHiGainPix.h (revision 3263) +++ /trunk/MagicSoft/Mars/mcalib/MHCalibrationChargeHiGainPix.h (revision 3264) @@ -7,10 +7,8 @@ #endif -class MCalibrationChargePix; +class TH1F; class MHCalibrationChargeHiGainPix : public MHCalibrationChargePix { private: - - MCalibrationChargePix *fPix; //! Storage container of the results static const Int_t fgChargeNbins; @@ -38,5 +36,8 @@ void SetAbsTimeLast (const Axis_t last =fgAbsTimeLast) { fAbsTimeLast = last; } - void SetPix (MCalibrationChargePix *pix ) { fPix = pix; } + void ChangeHistId(Int_t i); + + // Draw + void Draw(Option_t *opt=""); ClassDef(MHCalibrationChargeHiGainPix, 1) // Histogram class for a Calibration Pixel with extracted Hi Gain @@ -47,2 +48,5 @@ + + + Index: /trunk/MagicSoft/Mars/mcalib/MHCalibrationChargeLoGainPix.cc =================================================================== --- /trunk/MagicSoft/Mars/mcalib/MHCalibrationChargeLoGainPix.cc (revision 3263) +++ /trunk/MagicSoft/Mars/mcalib/MHCalibrationChargeLoGainPix.cc (revision 3264) @@ -52,7 +52,7 @@ const Axis_t MHCalibrationChargeLoGainPix::fgChargeFirst = -0.5; const Axis_t MHCalibrationChargeLoGainPix::fgChargeLast = 199.5; -const Int_t MHCalibrationChargeLoGainPix::fgAbsTimeNbins = 30; +const Int_t MHCalibrationChargeLoGainPix::fgAbsTimeNbins = 15; const Axis_t MHCalibrationChargeLoGainPix::fgAbsTimeFirst = -0.5; -const Axis_t MHCalibrationChargeLoGainPix::fgAbsTimeLast = 29.5; +const Axis_t MHCalibrationChargeLoGainPix::fgAbsTimeLast = 14.5; // -------------------------------------------------------------------------- // @@ -60,9 +60,7 @@ // MHCalibrationChargeLoGainPix::MHCalibrationChargeLoGainPix(const char *name, const char *title) - : fPix(NULL) { - fName = name ? name : "MHCalibrationChargeLoGainPix"; - fTitle = title ? title : "Fill the FADC sums of the Low Gain events and perform the fits"; + fTitle = title ? title : "Fill the FADC sums of the Low Gain events and perform the fits Pixel "; SetChargeNbins(); @@ -73,5 +71,4 @@ SetAbsTimeFirst(); SetAbsTimeLast(); - } @@ -98,4 +95,20 @@ +void MHCalibrationChargeLoGainPix::ChangeHistId(Int_t id) +{ + + fPixId = id; + + fHGausHist.SetName(Form("%s%d", fHGausHist.GetName(), id)); + fHGausHist.SetTitle(Form("%s%d", fHGausHist.GetTitle(), id)); + + fHAbsTime.SetName(Form("%s%d", fHAbsTime.GetName(), id)); + fHAbsTime.SetTitle(Form("%s%d", fHAbsTime.GetTitle(), id)); + + fName = Form("%s%d", fName.Data(), id); + fTitle = Form("%s%d", fTitle.Data(), id); +} + + Index: /trunk/MagicSoft/Mars/mcalib/MHCalibrationChargeLoGainPix.h =================================================================== --- /trunk/MagicSoft/Mars/mcalib/MHCalibrationChargeLoGainPix.h (revision 3263) +++ /trunk/MagicSoft/Mars/mcalib/MHCalibrationChargeLoGainPix.h (revision 3264) @@ -7,10 +7,7 @@ #endif -class MCalibrationChargePix; class MHCalibrationChargeLoGainPix : public MHCalibrationChargePix { private: - - MCalibrationChargePix *fPix; //! Storage container of the results static const Int_t fgChargeNbins; @@ -38,5 +35,5 @@ void SetAbsTimeLast (const Axis_t last =fgAbsTimeLast) { fAbsTimeLast = last; } - void SetPix (MCalibrationChargePix *pix ) { fPix = pix; } + void ChangeHistId(Int_t i); ClassDef(MHCalibrationChargeLoGainPix, 1) // Histogram class for a Calibration Pixel with extracted Low Gain @@ -47,2 +44,7 @@ + + + + + Index: /trunk/MagicSoft/Mars/mcalib/MHCalibrationChargePix.cc =================================================================== --- /trunk/MagicSoft/Mars/mcalib/MHCalibrationChargePix.cc (revision 3263) +++ /trunk/MagicSoft/Mars/mcalib/MHCalibrationChargePix.cc (revision 3264) @@ -22,5 +22,4 @@ ! \* ======================================================================== */ - ////////////////////////////////////////////////////////////////////////////// // @@ -58,4 +57,5 @@ const Axis_t MHCalibrationChargePix::fgAbsTimeLast = 14.5; +const Float_t MHCalibrationChargePix::fgPickupLimit = 5.; const Int_t MHCalibrationChargePix::fgPulserFrequency = 200; // -------------------------------------------------------------------------- @@ -77,4 +77,6 @@ SetAbsTimeFirst(); SetAbsTimeLast(); + + SetPickupLimit(); fHAbsTime.UseCurrentStyle(); @@ -104,10 +106,11 @@ } - void MHCalibrationChargePix::Clear(Option_t *o) { - fPixId = -1; + fPixId = -1; + fSaturated = 0; + fPickup = 0.; MHGausEvents::Clear(); @@ -174,5 +177,8 @@ if (option.Contains("time")) - win++; + { + option.ReplaceAll("events",""); + win++; + } pad->SetTicks(); @@ -193,4 +199,43 @@ } +Bool_t MHCalibrationChargePix::RepeatFit(const Option_t *option) +{ + + // + // Get new fitting ranges + // + Axis_t rmin = GetMean() - fPickupLimit * GetSigma(); + Axis_t rmax = GetMean() + fPickupLimit * GetSigma(); + + GetFGausFit()->SetRange(rmin,rmax); + + GetHGausHist()->Fit(GetFGausFit(),option); + + SetMean ( GetFGausFit()->GetParameter(1) ); + SetMeanErr ( GetFGausFit()->GetParameter(2) ); + SetSigma ( GetFGausFit()->GetParError(1) ); + SetSigmaErr ( GetFGausFit()->GetParError(2) ); + SetProb ( GetFGausFit()->GetProb() ); + + // + // The fit result is accepted under condition: + // 1) The results are not nan's + // 2) The NDF is not smaller than fNDFLimit (5) + // 3) The Probability is greater than fProbLimit (default 0.001 == 99.9%) + // + if ( TMath::IsNaN ( GetMean() ) + || TMath::IsNaN ( GetMeanErr() ) + || TMath::IsNaN ( GetProb() ) + || TMath::IsNaN ( GetSigma() ) + || TMath::IsNaN ( GetSigmaErr() ) + || GetFGausFit()->GetNDF() < fNDFLimit + || GetProb() < fProbLimit ) + return kFALSE; + + SetGausFitOK(kTRUE); + return kTRUE; + +} + void MHCalibrationChargePix::BypassFit() { @@ -200,8 +245,20 @@ // with the fit, we take the histogram values // - SetMean(fHGausHist.GetMean()); - SetMeanErr(fHGausHist.GetRMS()/fHGausHist.GetEntries()); - SetSigma(fHGausHist.GetRMS()); - SetSigmaErr(fHGausHist.GetRMS()/fHGausHist.GetEntries()/2.); -} - + SetMean ( fHGausHist.GetMean() ); + SetMeanErr ( fHGausHist.GetRMS() / fHGausHist.GetEntries() ); + SetSigma ( fHGausHist.GetRMS() ); + SetSigmaErr ( fHGausHist.GetRMS() / fHGausHist.GetEntries() / 2. ); +} + +void MHCalibrationChargePix::CountPickup() +{ + fPickup = GetHGausHist()->Integral(GetHGausHist()->GetXaxis()->FindBin(GetMean()+fPickupLimit*GetSigma()), + GetHGausHist()->GetXaxis()->GetLast(), + "width"); +} + + + + + + Index: /trunk/MagicSoft/Mars/mcalib/MHCalibrationChargePix.h =================================================================== --- /trunk/MagicSoft/Mars/mcalib/MHCalibrationChargePix.h (revision 3263) +++ /trunk/MagicSoft/Mars/mcalib/MHCalibrationChargePix.h (revision 3264) @@ -20,4 +20,6 @@ static const Axis_t fgAbsTimeLast; + static const Float_t fgPickupLimit; + protected: @@ -36,5 +38,8 @@ Axis_t fAbsTimeLast; - Bool_t fSaturated; + Float_t fPickupLimit; // The limit in number of sigmas from the fitted mean above which events are considered as pickup + + Int_t fSaturated; + Float_t fPickup; public: @@ -46,4 +51,5 @@ virtual void Reset(); virtual void Init(); + virtual void ChangeHistId(Int_t i); // Setters @@ -56,7 +62,9 @@ virtual void SetAbsTimeLast( const Axis_t last =fgAbsTimeLast) { fAbsTimeLast = last; } + virtual void SetPickupLimit( const Float_t lim =fgPickupLimit) { fPickupLimit = lim; } + void SetPulserFrequency(Float_t f=fgPulserFrequency); - void SetSaturated( const Bool_t b) { fSaturated = b; } + void SetSaturated (const Int_t i) { fSaturated += i; } // Getters @@ -69,17 +77,20 @@ const Float_t GetAbsTimeRms() const; - const Bool_t IsSaturated() const { return fSaturated; } + const Int_t GetSaturated() const { return fSaturated; } + const Float_t GetPickup() const { return fPickup; } - void ChangeHistId(Int_t i); - // Fill histos Bool_t FillAbsTime(const Float_t t); // Fits - void BypassFit(); + void BypassFit(); + Bool_t RepeatFit(const Option_t *option="RQ0"); // Draws virtual void Draw(Option_t *opt=""); - + + // Miscelleaneous + void CountPickup(); + ClassDef(MHCalibrationChargePix, 1) // Base class for charge calibrated pixel }; Index: unk/MagicSoft/Mars/mcalib/MHCalibrationPixel.cc =================================================================== --- /trunk/MagicSoft/Mars/mcalib/MHCalibrationPixel.cc (revision 3263) +++ (revision ) @@ -1,937 +1,0 @@ -/* ======================================================================== *\ -! -! * -! * This file is part of MARS, the MAGIC Analysis and Reconstruction -! * Software. It is distributed to you in the hope that it can be a useful -! * and timesaving tool in analysing Data of imaging Cerenkov telescopes. -! * It is distributed WITHOUT ANY WARRANTY. -! * -! * Permission to use, copy, modify and distribute this software and its -! * documentation for any purpose is hereby granted without fee, -! * provided that the above copyright notice appear in all copies and -! * that both that copyright notice and this permission notice appear -! * in supporting documentation. It is provided "as is" without express -! * or implied warranty. -! * -! -! -! Author(s): Markus Gaug 11/2003 -! -! Copyright: MAGIC Software Development, 2000-2002 -! -! -\* ======================================================================== */ - -////////////////////////////////////////////////////////////////////////////// -// -// MHCalibrationPixel -// -// Performs all the necessary fits to extract the mean number of photons -// out of the derived light flux -// -////////////////////////////////////////////////////////////////////////////// -#include "MHCalibrationPixel.h" -#include "MHCalibrationConfig.h" - -#include -#include -#include - -#include -#include -#include -#include -#include - -#include "MFFT.h" - -#include "MLog.h" -#include "MLogManip.h" - -ClassImp(MHCalibrationPixel); - -using namespace std; - -const Int_t MHCalibrationPixel::fChargeNbinsHiGain = 2100; -const Int_t MHCalibrationPixel::fChargeNbinsLoGain = 1010; -const Int_t MHCalibrationPixel::fAbsTimeNbins = 16; -const Axis_t MHCalibrationPixel::fAbsTimeFirst = - 0.25; -const Axis_t MHCalibrationPixel::fAbsTimeLast = 15.75; -const Float_t MHCalibrationPixel::fProbLimit = 0.001; -const Int_t MHCalibrationPixel::fNDFLimit = 5; - -const Int_t MHCalibrationPixel::fPSDNbins = 30; -const Int_t MHCalibrationPixel::fPulserFrequency = 500; - -// -------------------------------------------------------------------------- -// -// Default Constructor. -// -MHCalibrationPixel::MHCalibrationPixel(const char *name, const char *title) - : fPixId(-1), - fHivsLoGain(NULL), - fPSDHiGain(NULL), - fPSDLoGain(NULL), - fChargeGausFit(NULL), - fHPSD(NULL), - fPSDExpFit(NULL), - fChargeXaxis(NULL), - fPSDXaxis(NULL), - fFitLegend(NULL), - fCurrentSize(1024) -{ - - fName = name ? name : "MHCalibrationPixel"; - fTitle = title ? title : "Fill the accumulated charges of all events and perform fits"; - - fChargeFirstHiGain = -100.5; - fChargeLastHiGain = 1999.5; - fChargeFirstLoGain = -100.5; - fChargeLastLoGain = 9999.5; - - // Create a large number of bins, later we will rebin - fHChargeHiGain = new TH1F("HChargeHiGain","Distribution of Summed FADC Hi Gain Slices Pixel ", - fChargeNbinsHiGain,fChargeFirstHiGain,fChargeLastHiGain); - fHChargeLoGain = new TH1F("HChargeLoGain","Distribution of Summed FADC Lo Gain Slices Pixel ", - fChargeNbinsLoGain,fChargeFirstLoGain,fChargeLastLoGain); - - fHChargeLoGain->SetXTitle("Sum FADC Slices (Lo Gain)"); - fHChargeHiGain->SetXTitle("Sum FADC Slices (Hi Gain)"); - - fHChargeLoGain->SetYTitle("Nr. of events"); - fHChargeHiGain->SetYTitle("Nr. of events"); - - fHChargeHiGain->Sumw2(); - fHChargeLoGain->Sumw2(); - - // Absolute Times - fHAbsTimeHiGain = new TH1F("HAbsTimeHiGain","Distribution of Absolute Arrival Hi Gain Times Pixel ", - fAbsTimeNbins,fAbsTimeFirst,fAbsTimeLast); - fHAbsTimeLoGain = new TH1F("HAbsTimeLoGain","Distribution of Absolute Arrival Lo Gain Times Pixel ", - fAbsTimeNbins,fAbsTimeFirst,fAbsTimeLast); - - fHAbsTimeHiGain->SetXTitle("Absolute Arrival Time [Hi Gain FADC slice nr]"); - fHAbsTimeLoGain->SetXTitle("Absolute Arrival Time [Lo Gain FADC slice nr]"); - - fHAbsTimeHiGain->SetYTitle("Nr. of events"); - fHAbsTimeLoGain->SetYTitle("Nr. of events"); - - - fHChargeHiGain->SetDirectory(NULL); - fHChargeLoGain->SetDirectory(NULL); - fHAbsTimeHiGain->SetDirectory(NULL); - fHAbsTimeLoGain->SetDirectory(NULL); - - fHiGains = new TArrayF(fCurrentSize); - fLoGains = new TArrayF(fCurrentSize); - - Clear(); - -} - -MHCalibrationPixel::~MHCalibrationPixel() -{ - - delete fHChargeHiGain; - delete fHAbsTimeHiGain; - - delete fHChargeLoGain; - delete fHAbsTimeLoGain; - - delete fHiGains; - delete fLoGains; - - if (fChargeGausFit) - delete fChargeGausFit; - if (fPSDExpFit) - delete fPSDExpFit; - if (fHPSD) - delete fHPSD; - if (fFitLegend) - delete fFitLegend; - if (fHivsLoGain) - delete fHivsLoGain; - if (fChargeXaxis) - delete fChargeXaxis; - if (fPSDXaxis) - delete fPSDXaxis; -} - - -void MHCalibrationPixel::Clear(Option_t *o) -{ - fTotalEntries = 0; - fCurrentSize = 1024; - - fChargeFirstHiGain = -100.5; - fChargeLastHiGain = 1999.5; - fChargeFirstLoGain = -100.5; - fChargeLastLoGain = 9999.5; - - fChargeChisquare = -1; - fChargeProb = -1; - fChargeNdf = -1; - - fAbsTimeFirstHiGain = -1; - fAbsTimeFirstLoGain = -1; - fAbsTimeLastHiGain = -1; - fAbsTimeLastLoGain = -1; - - fAbsTimeMean = -1; - fAbsTimeMeanErr = -1; - fAbsTimeRms = -1; - - fOffset = 0.; - fSlope = 0.; - - if (fChargeGausFit) - { - delete fChargeGausFit; - fChargeGausFit=0; - } - if (fPSDExpFit) - { - delete fPSDExpFit; - fPSDExpFit=0; - } - if (fHPSD) - { - delete fHPSD; - fHPSD=0; - } - if (fFitLegend) - { - delete fFitLegend; - fFitLegend=0; - } - if (fHivsLoGain) - { - delete fHivsLoGain; - fHivsLoGain=0; - } - if (fChargeXaxis) - { - delete fChargeXaxis; - fChargeXaxis=0; - } - if (fPSDXaxis) - { - delete fPSDXaxis; - fPSDXaxis=0; - } - if (fPSDHiGain) - { - delete fPSDHiGain; - fPSDHiGain=0; - } - if (fPSDLoGain) - { - delete fPSDLoGain; - fPSDLoGain=0; - } - - CLRBIT(fFlags, kUseLoGain); - CLRBIT(fFlags, kChargeFitOK); - CLRBIT(fFlags, kOscillating); -} - - -void MHCalibrationPixel::Reset() -{ - Clear(); - - fHChargeHiGain->Reset(); - fHChargeLoGain->Reset(); - - fHAbsTimeHiGain->Reset(); - fHAbsTimeLoGain->Reset(); - - fHiGains->Set(1024); - fLoGains->Set(1024); - - fHiGains->Reset(0.); - fLoGains->Reset(0.); -} - -void MHCalibrationPixel::SetUseLoGain(Bool_t b) -{ - if (b) - SETBIT(fFlags, kUseLoGain) ; - else - CLRBIT(fFlags, kUseLoGain); -} - - -Bool_t MHCalibrationPixel::CheckOscillations() -{ - - if (fPSDExpFit) - return IsOscillating(); - - // - // The number of entries HAS to be a potence of 2, - // so we can only cut out from the last potence of 2 to the rest. - // Another possibility would be to fill everything with - // zeros, but that gives a low frequency peak, which we would - // have to cut out later again. - // - // So, we have to live with the possibility that at the end - // of the calibration run, something has happened without noticing - // it... - // - - // This cuts only the non-used zero's, but MFFT will later cut the rest - CutArrayBorder(fHiGains); - CutArrayBorder(fLoGains); - - MFFT fourier; - - fPSDLoGain = fourier.PowerSpectrumDensity(fLoGains); - fPSDHiGain = fourier.PowerSpectrumDensity(fHiGains); - - if (IsUseLoGain()) - fHPSD = ProjectArray(*fPSDLoGain, fPSDNbins, - Form("%s%d","PSDProjection",fPixId), - Form("%s%d","Power Spectrum Density Projection LoGain ",fPixId)); - else - fHPSD = ProjectArray(*fPSDHiGain, fPSDNbins, - Form("%s%d","PSDProjection",fPixId), - Form("%s%d","Power Spectrum Density Projection HiGain ",fPixId)); - - // - // First guesses for the fit (should be as close to reality as possible, - // - const Double_t xmax = fHPSD->GetXaxis()->GetXmax(); - - fPSDExpFit = new TF1(Form("%s%d","PSDExpFit",fPixId),"exp([0]-[1]*x)",0.,xmax); - - const Double_t slope_guess = (TMath::Log(fHPSD->GetEntries())+1.)/xmax; - const Double_t offset_guess = slope_guess*xmax; - - fPSDExpFit->SetParameters(offset_guess, slope_guess); - fPSDExpFit->SetParNames("Offset","Slope"); - fPSDExpFit->SetParLimits(0,offset_guess/2.,2.*offset_guess); - fPSDExpFit->SetParLimits(1,slope_guess/1.5,1.5*slope_guess); - fPSDExpFit->SetRange(0.,xmax); - - fHPSD->Fit(fPSDExpFit,"RQL0"); - - fPSDProb = fPSDExpFit->GetProb(); - - if (fPSDProb < gkProbLimit) - { - SETBIT(fFlags,kOscillating); - return kTRUE; - } - - CLRBIT(fFlags,kOscillating); - - return kFALSE; -} - -void MHCalibrationPixel::CreatePSDXaxis(Int_t n) -{ - - if (fPSDXaxis) - return; - - fPSDXaxis = new TArrayF(n); - - for (Int_t i=0;iAddAt((Float_t)(fPulserFrequency*i)/2./n,i); -} - -void MHCalibrationPixel::CreateChargeXaxis(Int_t n) -{ - - if (!fChargeXaxis) - { - fChargeXaxis = new TArrayF(n); - for (Int_t i=0;iAddAt((Float_t)i/fPulserFrequency,i); - return; - } - - if (fChargeXaxis->GetSize() == n) - return; - - const Int_t diff = fChargeXaxis->GetSize()-n; - - fChargeXaxis->Set(n); - if (diff < 0) - for (Int_t i=n-1;i>n-diff-1;i--) - fChargeXaxis->AddAt((Float_t)i,i); -} - -void MHCalibrationPixel::CutArrayBorder(TArrayF *array) -{ - - Int_t i; - - for (i=array->GetSize()-1;i>=0;i--) - if (array->At(i) != 0) - { - array->Set(i+1); - break; - } -} - - - -Bool_t MHCalibrationPixel::IsEmpty() const -{ - return !(fHChargeHiGain->GetEntries() || fHChargeLoGain->GetEntries()); -} - -Bool_t MHCalibrationPixel::IsUseLoGain() const -{ - return TESTBIT(fFlags,kUseLoGain); -} - -Bool_t MHCalibrationPixel::IsChargeFitOK() const -{ - return TESTBIT(fFlags,kChargeFitOK); -} - -Bool_t MHCalibrationPixel::IsOscillating() -{ - - if (fPSDExpFit) - return TESTBIT(fFlags,kOscillating); - - return CheckOscillations(); - -} - -Bool_t MHCalibrationPixel::FillChargeLoGain(Float_t q) -{ - return (fHChargeLoGain->Fill(q) > -1); -} - -Bool_t MHCalibrationPixel::FillAbsTimeLoGain(Float_t t) -{ - return (fHAbsTimeLoGain->Fill(t)> -1); -} - -Bool_t MHCalibrationPixel::FillChargeHiGain(Float_t q) -{ - return (fHChargeHiGain->Fill(q) > -1); -} - -Bool_t MHCalibrationPixel::FillAbsTimeHiGain(Float_t t) -{ - return (fHAbsTimeHiGain->Fill(t) > -1); -} - -void MHCalibrationPixel::ChangeHistId(Int_t id) -{ - - // Change only if the names have not yet been changed - if (fPixId == -1) - { - - // - // Names Hi gain Histograms - // - fHChargeHiGain->SetName( Form("%s%d",fHChargeHiGain->GetName(), id)); - fHAbsTimeHiGain->SetName( Form("%s%d",fHAbsTimeHiGain->GetName(), id)); - - // - // Title Hi gain Histograms - // - fHChargeHiGain->SetTitle( Form("%s%d",fHChargeHiGain->GetTitle(), id)); - fHAbsTimeHiGain->SetTitle( Form("%s%d",fHAbsTimeHiGain->GetTitle(), id)); - - // - // Names Low Gain Histograms - // - fHChargeLoGain->SetName( Form("%s%d",fHChargeLoGain->GetName(),id)); - fHAbsTimeLoGain->SetName( Form("%s%d",fHAbsTimeLoGain->GetName(),id)); - - // - // Titles Low Gain Histograms - // - fHChargeLoGain->SetTitle( Form("%s%d",fHChargeLoGain->GetTitle(),id)); - fHAbsTimeLoGain->SetTitle( Form("%s%d",fHAbsTimeLoGain->GetTitle(),id)); - - fPixId = id; - } - -} - - -Bool_t MHCalibrationPixel::UseLoGain() -{ - - if (fHChargeHiGain->Integral() > fHChargeLoGain->Integral()) - { - CLRBIT(fFlags,kUseLoGain); - return kFALSE; - } - else - { - SETBIT(fFlags,kUseLoGain); - return kTRUE; - } -} - -Bool_t MHCalibrationPixel::FillGraphs(Float_t qhi,Float_t qlo) -{ - - if (fTotalEntries >= fCurrentSize) - { - fCurrentSize *= 2; - - fHiGains->Set(fCurrentSize); - fLoGains->Set(fCurrentSize); - } - - fHiGains->AddAt(qhi,fTotalEntries); - fLoGains->AddAt(qlo,fTotalEntries); - - fTotalEntries++; - - return kTRUE; - -} - - -// ------------------------------------------------------------------------- -// -// Draw a legend with the fit results -// -void MHCalibrationPixel::DrawLegend() -{ - - fFitLegend = new TPaveText(0.05,0.05,0.95,0.95); - - if (IsChargeFitOK()) - fFitLegend->SetFillColor(80); - else - fFitLegend->SetFillColor(2); - - fFitLegend->SetLabel("Results of the Gauss Fit:"); - fFitLegend->SetTextSize(0.05); - - const TString line1 = - Form("Mean: Q_{#mu} = %2.2f #pm %2.2f",fChargeMean,fChargeMeanErr); - TText *t1 = fFitLegend->AddText(line1); - t1->SetBit(kCanDelete); - - const TString line4 = - Form("Sigma: #sigma_{Q} = %2.2f #pm %2.2f",fChargeSigma,fChargeSigmaErr); - TText *t2 = fFitLegend->AddText(line4); - t2->SetBit(kCanDelete); - - const TString line7 = - Form("#chi^{2} / N_{dof}: %4.2f / %3i",fChargeChisquare,fChargeNdf); - TText *t3 = fFitLegend->AddText(line7); - t3->SetBit(kCanDelete); - - const TString line8 = - Form("Probability: %4.3f ",fChargeProb); - TText *t4 = fFitLegend->AddText(line8); - t4->SetBit(kCanDelete); - - if (IsChargeFitOK()) - { - TText *t5 = fFitLegend->AddText("Result of the Fit: OK"); - t5->SetBit(kCanDelete); - } - else - { - TText *t6 = fFitLegend->AddText("Result of the Fit: NOT OK"); - t6->SetBit(kCanDelete); - } - - fFitLegend->SetBit(kCanDelete); - fFitLegend->Draw(); - -} - - -TObject *MHCalibrationPixel::DrawClone(Option_t *option) const -{ - - gROOT->SetSelectedPad(NULL); - - MHCalibrationPixel *newobj = (MHCalibrationPixel*)Clone(); - - if (!newobj) - return 0; - newobj->SetBit(kCanDelete); - - if (strlen(option)) - newobj->Draw(option); - else - newobj->Draw(GetDrawOption()); - - return newobj; -} - - - -// ------------------------------------------------------------------------- -// -// Draw the histogram -// -void MHCalibrationPixel::Draw(Option_t *opt) -{ - - if (!fHivsLoGain) - FitHiGainvsLoGain(); - - gStyle->SetOptFit(0); - gStyle->SetOptStat(111111); - - gROOT->SetSelectedPad(NULL); - - TCanvas *c = MH::MakeDefCanvas(this,600,900); - c->SetBit(kCanDelete); - - c->Divide(2,5); - - c->cd(1); - gPad->SetBorderMode(0); - gPad->SetTicks(); - - if (fHChargeHiGain->Integral() > 0) - gPad->SetLogy(); - - fHChargeHiGain->Draw(opt); - - if (IsUseLoGain()) - { - - c->cd(2); - gPad->SetTicks(); - - if (fHChargeLoGain->Integral() > 0) - gPad->SetLogy(); - - fHChargeLoGain->Draw(opt); - - if (fChargeGausFit) - { - fChargeGausFit->SetLineColor(IsChargeFitOK() ? kGreen : kRed); - fChargeGausFit->Draw("same"); - } - - c->cd(3); - gROOT->SetSelectedPad(NULL); - gStyle->SetOptFit(); - if (fHivsLoGain) - fHivsLoGain->Draw("prof"); - - c->cd(4); - DrawLegend(); - } - else - { - if (fChargeGausFit) - { - fChargeGausFit->SetLineColor(IsChargeFitOK() ? kGreen : kRed); - fChargeGausFit->Draw("same"); - } - - c->cd(2); - gPad->SetTicks(); - - if (fHChargeLoGain->Integral() > 0) - gPad->SetLogy(1); - - fHChargeLoGain->Draw(opt); - - c->cd(3); - DrawLegend(); - - c->cd(4); - - gROOT->SetSelectedPad(NULL); - gStyle->SetOptFit(); - if (fHivsLoGain) - fHivsLoGain->Draw("prof"); - } - - c->cd(5); - gPad->SetTicks(); - fHAbsTimeHiGain->Draw(opt); - - c->cd(6); - gPad->SetTicks(); - fHAbsTimeLoGain->Draw(opt); - - CreateChargeXaxis(fHiGains->GetSize()); - - c->cd(7); - gPad->SetTicks(); - TGraph *gr1 = new TGraph(fChargeXaxis->GetSize(), - fChargeXaxis->GetArray(), - fHiGains->GetArray()); - gr1->SetTitle("Evolution of HiGain charges with time"); - gr1->GetXaxis()->SetTitle("Time [s]"); - gr1->GetYaxis()->SetTitle("Sum FADC slices"); - gr1->SetBit(kCanDelete); - gr1->Draw("AL"); - - CreateChargeXaxis(fLoGains->GetSize()); - - c->cd(8); - gPad->SetTicks(); - TGraph *gr2 = new TGraph(fChargeXaxis->GetSize(), - fChargeXaxis->GetArray(), - fLoGains->GetArray()); - gr2->SetTitle("Evolution of LoGain charges with time"); - gr2->GetXaxis()->SetTitle("Time [s]"); - gr2->GetYaxis()->SetTitle("Sum FADC slices"); - gr2->SetBit(kCanDelete); - gr2->Draw("AL"); - - c->Modified(); - c->Update(); - - c->cd(9); - - TArrayF *array; - TString title = "Power Spectrum Density "; - - if(IsUseLoGain()) - { - if (!fPSDLoGain) - return; - array = fPSDLoGain; - title += "LoGain"; - } - else - { - if (!fPSDHiGain) - return; - array = fPSDHiGain; - title += "HiGain"; - } - - if (!fPSDXaxis) - CreatePSDXaxis(array->GetSize()); - - TGraph *gr3 = new TGraph(fPSDXaxis->GetSize(),fPSDXaxis->GetArray(),array->GetArray()); - gr3->SetTitle(title.Data()); - gr3->GetXaxis()->SetTitle("Freq. [Hz]"); - gr3->GetYaxis()->SetTitle("P(f)"); - gr3->SetBit(kCanDelete); - gr3->Draw("AL"); - - c->Modified(); - c->Update(); - - c->cd(10); - - gStyle->SetOptStat(111111); - gStyle->SetOptFit(); - gPad->SetTicks(); - - if (fHPSD->Integral() > 0) - gPad->SetLogy(); - - fHPSD->Draw(opt); - fHPSD->GetXaxis()->SetTitle("P(f)"); - fHPSD->GetYaxis()->SetTitle("Counts"); - - if (fPSDExpFit) - { - fPSDExpFit->SetLineColor(IsOscillating() ? kRed : kGreen); - fPSDExpFit->Draw("same"); - } - - c->Modified(); - c->Update(); - - return; -} - -void MHCalibrationPixel::FitHiGainvsLoGain() -{ - - *fLog << inf << "Fit Hi Gain vs Lo Gain " << endl; - - if (fHivsLoGain) - return; - - if ((!fHiGains) || (!fLoGains) || (fHiGains->GetSize() == 0) || (fLoGains->GetSize() == 0)) - return; - - gStyle->SetOptFit(); - - fHivsLoGain = new TProfile("HiGainvsLoGain","Plot the High Gain vs. Low Gain",100,0.,1000.,0.,1000.); - fHivsLoGain->GetXaxis()->SetTitle("Sum of Charges High Gain"); - fHivsLoGain->GetYaxis()->SetTitle("Sum of Charges Low Gain"); - fHivsLoGain->SetName(Form("%s%d",fHivsLoGain->GetName(),fPixId)); - - for (Int_t i=0;iFill(fHiGains->At(i),fLoGains->At(i),1); - - fHivsLoGain->Fit("pol1","rq","",fHivsLoGain->GetMean()-2.5*fHivsLoGain->GetRMS(), - fHivsLoGain->GetMean()+2.5*fHivsLoGain->GetRMS()); - - fOffset = fHivsLoGain->GetFunction("pol1")->GetParameter(0); - fSlope = fHivsLoGain->GetFunction("pol1")->GetParameter(1); - -} - - - -Bool_t MHCalibrationPixel::FitCharge(Option_t *option) -{ - - if (fChargeGausFit) - return kFALSE; - - // - // Get the fitting ranges - // - Axis_t rmin = fChargeFirstHiGain; - Axis_t rmax = fChargeLastHiGain; - TH1F *hist = fHChargeHiGain; - - if (TESTBIT(fFlags,kUseLoGain)) - { - rmin = fChargeFirstLoGain; - rmax = fChargeLastLoGain; - hist = fHChargeLoGain; - } - - // - // First guesses for the fit (should be as close to reality as possible, - // - const Stat_t entries = hist->Integral("width"); - const Double_t mu_guess = hist->GetBinCenter(hist->GetMaximumBin()); - const Double_t sigma_guess = (rmax-rmin)/2.; - const Double_t area_guess = entries/gkSq2Pi/sigma_guess; - - fChargeGausFit = new TF1(Form("%s%d","ChargeGausFit",fPixId),"gaus",rmin,rmax); - - if (!fChargeGausFit) - { - *fLog << warn << dbginf << "WARNING: Could not create fit function for Gauss fit" << endl; - return kFALSE; - } - - fChargeGausFit->SetParameters(area_guess,mu_guess,sigma_guess); - fChargeGausFit->SetParNames("Area","#mu","#sigma"); - fChargeGausFit->SetParLimits(0,0.,entries); - fChargeGausFit->SetParLimits(1,rmin,rmax); - fChargeGausFit->SetParLimits(2,0.,rmax-rmin); - fChargeGausFit->SetRange(rmin,rmax); - - hist->Fit(fChargeGausFit,option); - - // - // If we are not able to fit, try once again - // - if (fChargeGausFit->GetProb() < fProbLimit) - { - - Axis_t rtry = fChargeGausFit->GetParameter(1) - 2.5*fChargeGausFit->GetParameter(2); - rmin = (rtry < rmin ? rmin : rtry); - rmax = fChargeGausFit->GetParameter(1) + 2.5*fChargeGausFit->GetParameter(2); - - fChargeGausFit->SetRange(rmin,rmax); - hist->Fit(fChargeGausFit,option); - } - - fChargeChisquare = fChargeGausFit->GetChisquare(); - fChargeNdf = fChargeGausFit->GetNDF(); - fChargeProb = fChargeGausFit->GetProb(); - fChargeMean = fChargeGausFit->GetParameter(1); - fChargeMeanErr = fChargeGausFit->GetParError(1); - fChargeSigma = fChargeGausFit->GetParameter(2); - fChargeSigmaErr = fChargeGausFit->GetParError(2); - - // - // From the absolute time, we only take the mean and RMS - // - fAbsTimeMean = (Float_t)fHAbsTimeHiGain->GetMean(); - fAbsTimeRms = (Float_t)fHAbsTimeHiGain->GetRMS(); - fAbsTimeMeanErr = (Float_t)fAbsTimeRms / TMath::Sqrt(fHAbsTimeHiGain->GetEntries()); - - if (TESTBIT(fFlags,kUseLoGain)) - { - fAbsTimeMean = fHAbsTimeLoGain->GetMean(); - fAbsTimeRms = fHAbsTimeLoGain->GetRMS(); - fAbsTimeMeanErr = fAbsTimeRms / TMath::Sqrt(fHAbsTimeLoGain->GetEntries()); - } - - // - // The fit result is accepted under condition: - // 1) The results are not nan's - // 2) The NDF is not smaller than fNDFLimit (5) - // 3) The Probability is greater than fProbLimit (default 0.001 == 99.9%) - // - // Otherwise take means and RMS of the histograms - // - if ( TMath::IsNaN(fChargeMean) - || TMath::IsNaN(fChargeMeanErr) - || TMath::IsNaN(fChargeProb) - || TMath::IsNaN(fChargeSigma) - || TMath::IsNaN(fChargeSigmaErr) - || (fChargeNdf < fNDFLimit) - || (fChargeProb < fProbLimit) ) - { - - fChargeMean = hist->GetMean(); - fChargeMeanErr = hist->GetRMS()/hist->GetEntries(); - fChargeSigma = hist->GetRMS(); - fChargeSigmaErr = fChargeMeanErr/2.; - - CLRBIT(fFlags,kChargeFitOK); - return kFALSE; - } - - SETBIT(fFlags,kChargeFitOK); - return kTRUE; -} - - -void MHCalibrationPixel::CutAllEdges() -{ - - Int_t nbins = 30; - - StripZeros(fHChargeHiGain,nbins); - - fChargeFirstHiGain = fHChargeHiGain->GetBinLowEdge(fHChargeHiGain->GetXaxis()->GetFirst()); - fChargeLastHiGain = fHChargeHiGain->GetBinLowEdge(fHChargeHiGain->GetXaxis()->GetLast()) - +fHChargeHiGain->GetBinWidth(0); - - StripZeros(fHChargeLoGain,nbins); - - fChargeFirstLoGain = fHChargeLoGain->GetBinLowEdge(fHChargeLoGain->GetXaxis()->GetFirst()); - fChargeLastLoGain = fHChargeLoGain->GetBinLowEdge(fHChargeLoGain->GetXaxis()->GetLast()) - +fHChargeLoGain->GetBinWidth(0); - - - StripZeros(fHAbsTimeHiGain,0); - - fAbsTimeFirstHiGain = fHAbsTimeHiGain->GetBinLowEdge(fHAbsTimeHiGain->GetXaxis()->GetFirst()); - fAbsTimeLastHiGain = fHAbsTimeHiGain->GetBinLowEdge(fHAbsTimeHiGain->GetXaxis()->GetLast()) - +fHAbsTimeHiGain->GetBinWidth(0); - - StripZeros(fHAbsTimeLoGain,0); - - fAbsTimeFirstLoGain = fHAbsTimeLoGain->GetBinLowEdge(fHAbsTimeLoGain->GetXaxis()->GetFirst()); - fAbsTimeLastLoGain = fHAbsTimeLoGain->GetBinLowEdge(fHAbsTimeLoGain->GetXaxis()->GetLast()) - +fHAbsTimeLoGain->GetBinWidth(0); - - -} - -void MHCalibrationPixel::PrintChargeFitResult() -{ - - *fLog << all << "Results of the Summed Charges Fit: " << endl; - *fLog << all << "Chisquare: " << fChargeChisquare << endl; - *fLog << all << "DoF: " << fChargeNdf << endl; - *fLog << all << "Probability: " << fChargeProb << endl; - *fLog << all << endl; - -} - Index: unk/MagicSoft/Mars/mcalib/MHCalibrationPixel.h =================================================================== --- /trunk/MagicSoft/Mars/mcalib/MHCalibrationPixel.h (revision 3263) +++ (revision ) @@ -1,173 +1,0 @@ -#ifndef MARS_MHCalibrationPixel -#define MARS_MHCalibrationPixel - -#ifndef MARS_MH -#include "MH.h" -#endif - -class TArrayF; -class TH1F; -class TH1I; -class TF1; -class TProfile; -class TPaveText; - -class MHCalibrationPixel : public MH -{ - -private: - - static const Int_t fChargeNbinsHiGain; - static const Int_t fChargeNbinsLoGain; - static const Int_t fAbsTimeNbins; - static const Axis_t fAbsTimeFirst; - static const Axis_t fAbsTimeLast; - static const Float_t fProbLimit; - static const Int_t fNDFLimit; - - static const Int_t fPSDNbins; - static const Int_t fPulserFrequency; - - Int_t fPixId; // Pixel Nr - - TProfile* fHivsLoGain; //-> - - Double_t fOffset; - Double_t fSlope; - -protected: - - TH1F* fHChargeHiGain; //-> Summed FADC slices High Gain - TH1F* fHAbsTimeHiGain; //-> Mean arrival time in number of FADC sice - - TH1F* fHChargeLoGain; //-> Summed FADC slices Low Gain - TH1F* fHAbsTimeLoGain; //-> Mean arrival time in number of FADC sice - - TArrayF* fPSDHiGain; //-> Power spectrum density of fHiGains - TArrayF* fPSDLoGain; //-> Power spectrum density of fLoGains - - TF1* fChargeGausFit; //-> - - TH1I* fHPSD; //-> - TF1* fPSDExpFit; //-> - - TArrayF *fHiGains; //-> - TArrayF *fLoGains; //-> - TArrayF *fChargeXaxis; // - TArrayF *fPSDXaxis; // - - TPaveText *fFitLegend; //-> - - Int_t fTotalEntries; // Number of entries - Int_t fCurrentSize; - - Axis_t fChargeFirstHiGain; - Axis_t fChargeLastHiGain; - Axis_t fChargeFirstLoGain; - Axis_t fChargeLastLoGain; - - Double_t fChargeChisquare; - Double_t fChargeProb; - Int_t fChargeNdf; - - Double_t fChargeMean; - Double_t fChargeMeanErr; - Double_t fChargeSigma; - Double_t fChargeSigmaErr; - - Float_t fAbsTimeMean; - Float_t fAbsTimeMeanErr; - Float_t fAbsTimeRms; - - Float_t fAbsTimeFirstHiGain; - Float_t fAbsTimeFirstLoGain; - Float_t fAbsTimeLastHiGain; - Float_t fAbsTimeLastLoGain; - - Float_t fPSDProb; - - Byte_t fFlags; - - enum { kUseLoGain, kChargeFitOK, kOscillating }; - - virtual void DrawLegend(); - virtual void CreateChargeXaxis(Int_t n); - virtual void CreatePSDXaxis(Int_t n); - virtual void CutArrayBorder(TArrayF *array); - -public: - - MHCalibrationPixel(const char *name=NULL, const char *title=NULL); - ~MHCalibrationPixel(); - - void Clear(Option_t *o=""); - void Reset(); - - void ChangeHistId(Int_t i); - - // Setters - void SetUseLoGain(Bool_t b = kTRUE); - - // Getters - const Double_t GetChargeMean() const { return fChargeMean; } - const Double_t GetChargeMeanErr() const { return fChargeMeanErr; } - const Double_t GetChargeSigma() const { return fChargeSigma; } - const Double_t GetChargeSigmaErr() const { return fChargeSigmaErr; } - const Double_t GetChargeChiSquare() const { return fChargeChisquare; } - const Double_t GetChargeProb() const { return fChargeProb; } - const Int_t GetChargeNdf() const { return fChargeNdf; } - - const Float_t GetAbsTimeFirstHiGain() const { return fAbsTimeFirstHiGain; } - const Float_t GetAbsTimeFirstLoGain() const { return fAbsTimeFirstLoGain; } - const Float_t GetAbsTimeLastHiGain() const { return fAbsTimeLastHiGain; } - const Float_t GetAbsTimeLastLoGain() const { return fAbsTimeLastLoGain; } - - const Float_t GetAbsTimeMean() const { return fAbsTimeMean; } - const Float_t GetAbsTimeMeanErr() const { return fAbsTimeMeanErr; } - const Float_t GetAbsTimeRms() const { return fAbsTimeRms; } - - const TH1F *GetHCharge() { return fHChargeHiGain; } - const TH1F *GetHCharge() const { return fHChargeHiGain; } - - const TH1F *GetHAbsTime() { return fHAbsTimeHiGain; } - const TH1F *GetHAbsTime() const { return fHAbsTimeHiGain; } - - Double_t GetOffset() const { return fOffset; } - Double_t GetSlope() const { return fSlope; } - - Bool_t UseLoGain(); - Bool_t CheckOscillations(); - - Bool_t IsChargeFitOK() const; - Bool_t IsOscillating(); - Bool_t IsUseLoGain() const; - Bool_t IsEmpty() const; - - // Fill histos - Bool_t FillChargeLoGain(Float_t q); - Bool_t FillAbsTimeLoGain(Float_t t); - - Bool_t FillChargeHiGain(Float_t q); - Bool_t FillAbsTimeHiGain(Float_t t); - - Bool_t FillGraphs(Float_t qhi, Float_t qlo); - - // Fits - Bool_t FitCharge(Option_t *option="RQ0"); - - void FitHiGainvsLoGain(); - - // Draws - virtual void Draw(Option_t *option=""); - TObject *DrawClone(Option_t *option="") const; - - // Prints - void PrintChargeFitResult(); - - // Others - virtual void CutAllEdges(); - - ClassDef(MHCalibrationPixel, 1) // Histograms for each calibrated pixel -}; - -#endif Index: /trunk/MagicSoft/Mars/mcalib/MHGausEvents.cc =================================================================== --- /trunk/MagicSoft/Mars/mcalib/MHGausEvents.cc (revision 3263) +++ /trunk/MagicSoft/Mars/mcalib/MHGausEvents.cc (revision 3264) @@ -91,8 +91,8 @@ using namespace std; -const Float_t MHGausEvents::fgProbLimit = 0.01; +const Float_t MHGausEvents::fgProbLimit = 0.005; const Int_t MHGausEvents::fgNDFLimit = 2; -const Int_t MHGausEvents::fgPowerProbabilityBins = 25; -const Int_t MHGausEvents::fgBinsAfterStripping = 25; +const Int_t MHGausEvents::fgPowerProbabilityBins = 20; +const Int_t MHGausEvents::fgBinsAfterStripping = 40; // -------------------------------------------------------------------------- // @@ -423,5 +423,5 @@ const Stat_t entries = fHGausHist.Integral("width"); const Double_t mu_guess = fHGausHist.GetBinCenter(fHGausHist.GetMaximumBin()); - const Double_t sigma_guess = (rmax-rmin)/2.; + const Double_t sigma_guess = fHGausHist.GetRMS(); const Double_t area_guess = entries/TMath::Sqrt(TMath::TwoPi())/sigma_guess; @@ -436,5 +436,5 @@ fFGausFit->SetParameters(area_guess,mu_guess,sigma_guess); fFGausFit->SetParNames("Area","#mu","#sigma"); - fFGausFit->SetParLimits(0,0.,entries); + fFGausFit->SetParLimits(0,0.,area_guess*1.5); fFGausFit->SetParLimits(1,rmin,rmax); fFGausFit->SetParLimits(2,0.,rmax-rmin); @@ -442,10 +442,10 @@ fHGausHist.Fit(fFGausFit,option); - + + fMean = fFGausFit->GetParameter(1); fSigma = fFGausFit->GetParameter(2); fMeanErr = fFGausFit->GetParError(1); fSigmaErr = fFGausFit->GetParError(2); - fProb = fFGausFit->GetProb(); // @@ -455,5 +455,5 @@ // 3) The Probability is greater than fProbLimit (default 0.001 == 99.9%) // - if ( TMath::IsNaN(fMean) + if ( TMath::IsNaN(fMean) || TMath::IsNaN(fMeanErr) || TMath::IsNaN(fProb) @@ -578,5 +578,5 @@ if (!IsEmpty()) - pad->SetLogy(); + gPad->SetLogy(); fHGausHist.Draw(opt); @@ -587,5 +587,4 @@ fFGausFit->Draw("same"); } - switch (win) { Index: /trunk/MagicSoft/Mars/mcalib/MHGausEvents.h =================================================================== --- /trunk/MagicSoft/Mars/mcalib/MHGausEvents.h (revision 3263) +++ /trunk/MagicSoft/Mars/mcalib/MHGausEvents.h (revision 3264) @@ -25,6 +25,4 @@ const static Int_t fgBinsAfterStripping; // Default number of bins for the Gauss Histogram after stripping off the zeros at both end - Float_t fProbLimit; // Probability limit for judgement if fit is OK - Int_t fNDFLimit; // NDF limit for judgement if fit is OK Int_t fPowerProbabilityBins; // number of bins for the projected power spectrum Int_t fBinsAfterStripping; // number of bins for the Gauss Histogram after stripping off the zeros at both end @@ -59,7 +57,10 @@ TArrayF fEvents; // Array which holds the entries of GausHist + Float_t fProbLimit; // Probability limit for judgement if fit is OK + Int_t fNDFLimit; // NDF limit for judgement if fit is OK + // Setters - void SetPowerProbabilityBins(const Int_t nbins=fgPowerProbabilityBins) { fPowerProbabilityBins = nbins; } - void SetBinsAfterStripping(const Int_t nbins=fgBinsAfterStripping) { fBinsAfterStripping = nbins; } + void SetPowerProbabilityBins ( const Int_t nbins=fgPowerProbabilityBins ) { fPowerProbabilityBins = nbins; } + void SetBinsAfterStripping ( const Int_t nbins=fgBinsAfterStripping ) { fBinsAfterStripping = nbins; } void DrawEvents(); // Draw a graph of the array fEvents @@ -77,8 +78,9 @@ void SetEventFrequency(const Float_t f=0) { fEventFrequency = f; } - void SetMean( const Double_t d ) { fMean = d; } - void SetMeanErr( const Double_t d ) { fMeanErr = d; } - void SetSigma( const Double_t d ) { fSigma = d; } + void SetMean ( const Double_t d ) { fMean = d; } + void SetMeanErr ( const Double_t d ) { fMeanErr = d; } + void SetSigma ( const Double_t d ) { fSigma = d; } void SetSigmaErr( const Double_t d ) { fSigmaErr = d; } + void SetProb ( const Double_t d ) { fProb = d; } void SetProbLimit( const Float_t lim=fgProbLimit ) { fProbLimit = lim; }