Index: /trunk/MagicSoft/Mars/Changelog =================================================================== --- /trunk/MagicSoft/Mars/Changelog (revision 3478) +++ /trunk/MagicSoft/Mars/Changelog (revision 3479) @@ -58,5 +58,10 @@ - added Print function - + * mcalib/MCalibrationChargeCam.[h,cc] + * mcalib/MCalibrationChargePix.h + * mcalib/MCalibrationChargeCalc.cc + - removed pointers to MGeomCam and MBadPixelsCam in MCalibrationChargeCam + - use two loops over pixels in order to determine mean flux of phe;s to + discard pixels far outside the normal distribution 2004/03/10: Abelardo Moralejo Index: /trunk/MagicSoft/Mars/mbadpixels/MBadPixelsCam.cc =================================================================== --- /trunk/MagicSoft/Mars/mbadpixels/MBadPixelsCam.cc (revision 3478) +++ /trunk/MagicSoft/Mars/mbadpixels/MBadPixelsCam.cc (revision 3479) @@ -17,4 +17,5 @@ ! ! Author(s): Thomas Bretz 1/2004 +! Markus Gaug 3/2004 ! ! Copyright: MAGIC Software Development, 2000-2004 @@ -225,6 +226,7 @@ if (!(count % 25)) - *fLog << endl; - } + *fLog << endl; + } + *fLog << endl; *fLog << count << " unreliable pixels :-(" << endl; Index: /trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCalc.cc =================================================================== --- /trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCalc.cc (revision 3478) +++ /trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCalc.cc (revision 3479) @@ -331,7 +331,4 @@ } - fCam->SetGeomCam(fGeom); - fCam->SetBadPixelsCam(fBadPixels); - fNumHiGainSamples = fSignals->GetNumUsedHiGainFADCSlices(); fNumLoGainSamples = fSignals->GetNumUsedLoGainFADCSlices(); @@ -437,21 +434,23 @@ } - pix.CheckChargeValidity(&bad); - pix.CheckTimeValidity(&bad); + pix.CheckChargeValidity (&bad); + pix.CheckTimeValidity (&bad); - if (bad.IsUnsuitableRun()) + if (bad.IsUnsuitable(MBadPixelsPix::kUnsuitableRun)) continue; - + nvalid++; if (!pix.CalcReducedSigma()) { - bad.SetUnsuitableRun(); + bad.SetUnsuitable(MBadPixelsPix::kUnsuitableRun); continue; } - pix.CalcFFactorMethod(); - } - + if (!pix.CalcFFactorMethod()) + bad.SetUnsuitable(MBadPixelsPix::kUnsuitableRun); + + } + @@ -525,7 +524,7 @@ // F-Factor calibration // - if (fCam->CalcMeanFluxPhotonsFFactorMethod()) - { - fCam->ApplyFFactorCalibration(); + if (fCam->CalcMeanFluxPhotonsFFactorMethod(*fGeom, *fBadPixels)) + { + fCam->ApplyFFactorCalibration(*fGeom,*fBadPixels); fCam->SetFFactorMethodValid(kTRUE); } @@ -536,6 +535,4 @@ } - - // // Blind Pixel calibration @@ -558,5 +555,5 @@ { fCam->SetBlindPixelMethodValid(kTRUE); - fCam->ApplyBlindPixelCalibration(); + fCam->ApplyBlindPixelCalibration(*fGeom,*fBadPixels); } } @@ -582,5 +579,5 @@ { fCam->SetPINDiodeMethodValid(kTRUE); - fCam->ApplyPINDiodeCalibration(); + fCam->ApplyPINDiodeCalibration(*fGeom,*fBadPixels); } } Index: /trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCam.cc =================================================================== --- /trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCam.cc (revision 3478) +++ /trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCam.cc (revision 3479) @@ -86,35 +86,32 @@ // 27: Excluded Pixels // 28: Pixels where the fit did not succeed --> results obtained only from the histograms -// 29: Pixels with apparently wrong results -// 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 +// 29: Number of probable pickup events in the Hi Gain +// 30: Number of probable pickup events in the Lo Gain // // Other classifications of pixels: // ================================ // -// 34: Pixels with saturated Hi-Gain +// 31: Pixels with saturated Hi-Gain // // Classification of validity of the calibrations: // =============================================== // -// 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 +// 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 // // Used Pedestals: // =============== // -// 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 +// 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 // // Calculated absolute arrival times (very low precision!): // ======================================================== // -// 42: Absolute Arrival time of the signal -// 43: RMS of the Absolute Arrival time of the signal +// 39: Absolute Arrival time of the signal +// 40: RMS of the Absolute Arrival time of the signal // ///////////////////////////////////////////////////////////////////////////// @@ -145,4 +142,5 @@ const Float_t MCalibrationChargeCam::gkAverageQEErr = 0.02; const Float_t MCalibrationChargeCam::fgConvFFactorRelErrLimit = 0.25; +const Float_t MCalibrationChargeCam::fgPheFFactorRelLimit = 5.; // -------------------------------------------------------------------------- // @@ -155,9 +153,5 @@ // MCalibrationChargeCam::MCalibrationChargeCam(const char *name, const char *title) - : fBlindPixel(NULL), - fPINDiode(NULL), - fGeomCam(NULL), - fBadPixels(NULL), - fOffsets(NULL), + : fOffsets(NULL), fSlopes(NULL), fOffvsSlope(NULL) @@ -176,4 +170,5 @@ SetAverageQE(); SetConvFFactorRelErrLimit(); + SetPheFFactorRelLimit(); } @@ -316,5 +311,5 @@ int id = 0; - *fLog << all << "Succesfully calibrated pixels:" << endl; + *fLog << all << "Calibrated pixels:" << endl; *fLog << all << endl; @@ -324,9 +319,9 @@ { - if ((*fBadPixels)[pix->GetPixId()].IsCalibrationResultOK() && !pix->IsExcluded()) + if (!pix->IsExcluded()) { - *fLog << all << "Pix " << pix->GetPixId() - << ": Ped. Rms: " << pix->GetPedRms() << " +- " << pix->GetPedRmsErr() + *fLog << all << "Pix " << pix->GetPixId() + << ": Ped. Rms: " << pix->GetPedRms() << " +- " << pix->GetPedRmsErr() << " Mean signal: " << pix->GetMeanCharge() << " +- " << pix->GetSigmaCharge() << " Reduced Sigma: " << pix->GetRSigmaCharge() @@ -338,57 +333,8 @@ } - *fLog << all << id << " succesful pixels :-))" << endl; + *fLog << all << id << " pixels" << endl; id = 0; - *fLog << all << endl; - *fLog << all << "Pixels with errors:" << endl; - *fLog << all << endl; - - TIter Next2(fPixels); - while ((pix=(MCalibrationChargePix*)Next2())) - { - - if (!pix->IsExcluded() && !(*fBadPixels)[pix->GetPixId()].IsCalibrationResultOK()) - { - - - *fLog << all << "Pix " << pix->GetPixId() - << ": Ped. Rms: " << pix->GetPedRms() << " +- " << pix->GetPedRmsErr() - << " Mean signal: " << pix->GetMeanCharge() << " +- " << pix->GetSigmaCharge() - << " Reduced Sigma: " << pix->GetRSigmaCharge() - << " Nr Phe's: " << pix->GetPheFFactorMethod() - << " Saturated? :" << pix->IsHiGainSaturation() - << endl; - id++; - } - } - *fLog << all << id << " pixels with errors :-((" << endl; - - *fLog << all << endl; - *fLog << all << "Pixels with oscillations:" << endl; - *fLog << all << endl; - - id = 0; - - TIter Next3(fPixels); - while ((pix=(MCalibrationChargePix*)Next3())) - { - - if ( (*fBadPixels)[pix->GetPixId()].IsCalibrationOscillating() && !pix->IsExcluded()) - { - - *fLog << all << "Pix " << pix->GetPixId() - << ": Ped. Rms: " << pix->GetPedRms() << " +- " << pix->GetPedRmsErr() - << " Mean signal: " << pix->GetMeanCharge() << " +- " << pix->GetSigmaCharge() - << " Reduced Sigma: " << pix->GetRSigmaCharge() - << " Nr Phe's: " << pix->GetPheFFactorMethod() - << " Saturated? :" << pix->IsHiGainSaturation() - << endl; - id++; - } - } - *fLog << all << id << " Oscillating pixels :-((" << endl; - - + *fLog << all << endl; *fLog << all << "Excluded pixels:" << endl; @@ -408,4 +354,5 @@ *fLog << all << id << " Excluded pixels " << endl; *fLog << endl; + *fLog << all << "Average Inner Pix:" << " Ped. Rms: " << fAverageInnerPix->GetPedRms() << " +- " << fAverageInnerPix->GetPedRmsErr() @@ -474,35 +421,32 @@ // 27: Excluded Pixels // 28: Pixels where the fit did not succeed --> results obtained only from the histograms -// 29: Pixels with apparently wrong results -// 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 +// 29: Number of probable pickup events in the Hi Gain +// 30: Number of probable pickup events in the Lo Gain // // Other classifications of pixels: // ================================ // -// 34: Pixels with saturated Hi-Gain +// 31: Pixels with saturated Hi-Gain // // Classification of validity of the calibrations: // =============================================== // -// 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 +// 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 // // Used Pedestals: // =============== // -// 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 +// 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 // // Calculated absolute arrival times (very low precision!): // ======================================================== // -// 42: Absolute Arrival time of the signal -// 43: RMS of the Absolute Arrival time of the signal +// 39: Absolute Arrival time of the signal +// 40: 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 @@ -540,10 +484,10 @@ break; case 4: - if ((*this)[idx].IsExcluded() || !(*fBadPixels)[idx].IsCalibrationSignalOK()) + if ((*this)[idx].IsExcluded()) return kFALSE; val = (*this)[idx].GetChargeProb(); break; case 5: - if ((*this)[idx].IsExcluded() || !(*fBadPixels)[idx].IsCalibrationSignalOK()) + if ((*this)[idx].IsExcluded()) return kFALSE; if ((*this)[idx].GetRSigmaCharge() == -1.) @@ -552,5 +496,5 @@ break; case 6: - if ((*this)[idx].IsExcluded() || !(*fBadPixels)[idx].IsCalibrationSignalOK()) + if ((*this)[idx].IsExcluded()) return kFALSE; if ((*this)[idx].GetRSigmaCharge() == -1.) @@ -559,5 +503,5 @@ break; case 7: - if ((*this)[idx].IsExcluded() || !(*fBadPixels)[idx].IsCalibrationSignalOK()) + if ((*this)[idx].IsExcluded()) return kFALSE; if ((*this)[idx].GetRSigmaCharge() == -1.) @@ -566,5 +510,5 @@ break; case 8: - if ((*this)[idx].IsExcluded() || !(*fBadPixels)[idx].IsCalibrationSignalOK()) + if ((*this)[idx].IsExcluded()) return kFALSE; if ((*this)[idx].GetRSigmaCharge() == -1.) @@ -582,126 +526,90 @@ break; case 9: - if ((*this)[idx].IsExcluded() - || !(*fBadPixels)[idx].IsCalibrationSignalOK() - || !(*this)[idx].IsFFactorMethodValid()) + if ((*this)[idx].IsExcluded() || !(*this)[idx].IsFFactorMethodValid()) return kFALSE; val = (*this)[idx].GetPheFFactorMethod(); break; case 10: - if ((*this)[idx].IsExcluded() - || !(*fBadPixels)[idx].IsCalibrationSignalOK() - || !(*this)[idx].IsFFactorMethodValid()) + if ((*this)[idx].IsExcluded() || !(*this)[idx].IsFFactorMethodValid()) return kFALSE; val = (*this)[idx].GetPheFFactorMethodErr(); break; case 11: - if ((*this)[idx].IsExcluded() - || !(*fBadPixels)[idx].IsCalibrationSignalOK() - || !(*this)[idx].IsFFactorMethodValid()) + if ((*this)[idx].IsExcluded() || !(*this)[idx].IsFFactorMethodValid()) return kFALSE; val = (*this)[idx].GetMeanConversionFFactorMethod(); break; case 12: - if ((*this)[idx].IsExcluded() - || !(*fBadPixels)[idx].IsCalibrationSignalOK() - || !(*this)[idx].IsFFactorMethodValid()) + if ((*this)[idx].IsExcluded() || !(*this)[idx].IsFFactorMethodValid()) return kFALSE; val = (*this)[idx].GetConversionFFactorMethodErr(); break; case 13: - if ((*this)[idx].IsExcluded() - || !(*fBadPixels)[idx].IsCalibrationSignalOK() - || !(*this)[idx].IsFFactorMethodValid()) + if ((*this)[idx].IsExcluded() || !(*this)[idx].IsFFactorMethodValid()) return kFALSE; val = (*this)[idx].GetTotalFFactorFFactorMethod(); break; case 14: - if ((*this)[idx].IsExcluded() - || !(*fBadPixels)[idx].IsCalibrationSignalOK() - || !(*this)[idx].IsFFactorMethodValid()) + if ((*this)[idx].IsExcluded() || !(*this)[idx].IsFFactorMethodValid()) return kFALSE; val = (*this)[idx].GetTotalFFactorErrFFactorMethod(); break; case 15: - if ((*this)[idx].IsExcluded() - || !(*fBadPixels)[idx].IsCalibrationSignalOK() - || !(*this)[idx].IsBlindPixelMethodValid()) + if ((*this)[idx].IsExcluded() || !(*this)[idx].IsBlindPixelMethodValid()) return kFALSE; val = fBlindPixel->GetMeanFluxInsidePlexiglass()*area; break; case 16: - if ((*this)[idx].IsExcluded() - || !(*fBadPixels)[idx].IsCalibrationSignalOK() - || !(*this)[idx].IsBlindPixelMethodValid()) + if ((*this)[idx].IsExcluded() || !(*this)[idx].IsBlindPixelMethodValid()) return kFALSE; val = fBlindPixel->GetMeanFluxErrInsidePlexiglass()*area; break; case 17: - if ((*this)[idx].IsExcluded() - || !(*fBadPixels)[idx].IsCalibrationSignalOK() - || !(*this)[idx].IsBlindPixelMethodValid()) + if ((*this)[idx].IsExcluded() || !(*this)[idx].IsBlindPixelMethodValid()) return kFALSE; val = (*this)[idx].GetMeanConversionBlindPixelMethod(); break; case 18: - if ((*this)[idx].IsExcluded() - || !(*fBadPixels)[idx].IsCalibrationSignalOK() - || !(*this)[idx].IsBlindPixelMethodValid()) + if ((*this)[idx].IsExcluded() || !(*this)[idx].IsBlindPixelMethodValid()) return kFALSE; val = (*this)[idx].GetConversionBlindPixelMethodErr(); break; case 19: - if ((*this)[idx].IsExcluded() - || !(*fBadPixels)[idx].IsCalibrationSignalOK() - || !(*this)[idx].IsBlindPixelMethodValid()) + if ((*this)[idx].IsExcluded() || !(*this)[idx].IsBlindPixelMethodValid()) return kFALSE; val = (*this)[idx].GetTotalFFactorBlindPixelMethod(); break; case 20: - if ((*this)[idx].IsExcluded() - || !(*fBadPixels)[idx].IsCalibrationSignalOK() - || !(*this)[idx].IsBlindPixelMethodValid()) + if ((*this)[idx].IsExcluded() || !(*this)[idx].IsBlindPixelMethodValid()) return kFALSE; val = (*this)[idx].GetTotalFFactorErrBlindPixelMethod(); break; case 21: - if ((*this)[idx].IsExcluded() - || !(*fBadPixels)[idx].IsCalibrationSignalOK() - || !(*this)[idx].IsPINDiodeMethodValid()) + if ((*this)[idx].IsExcluded() || !(*this)[idx].IsPINDiodeMethodValid()) return kFALSE; val = fPINDiode->GetMeanFluxOutsidePlexiglass()*area; break; case 22: - if ((*this)[idx].IsExcluded() - || !(*fBadPixels)[idx].IsCalibrationSignalOK() - || !(*this)[idx].IsPINDiodeMethodValid()) + if ((*this)[idx].IsExcluded() || !(*this)[idx].IsPINDiodeMethodValid()) return kFALSE; val = fPINDiode->GetMeanFluxErrOutsidePlexiglass()*area; break; case 23: - if ((*this)[idx].IsExcluded() - || !(*fBadPixels)[idx].IsCalibrationSignalOK() - || !(*this)[idx].IsPINDiodeMethodValid()) + if ((*this)[idx].IsExcluded() || !(*this)[idx].IsPINDiodeMethodValid()) return kFALSE; val = (*this)[idx].GetMeanConversionPINDiodeMethod(); break; case 24: - if ((*this)[idx].IsExcluded() - || !(*fBadPixels)[idx].IsCalibrationSignalOK() - || !(*this)[idx].IsPINDiodeMethodValid()) + if ((*this)[idx].IsExcluded() || !(*this)[idx].IsPINDiodeMethodValid()) return kFALSE; val = (*this)[idx].GetConversionPINDiodeMethodErr(); break; case 25: - if ((*this)[idx].IsExcluded() - || !(*fBadPixels)[idx].IsCalibrationSignalOK() - || !(*this)[idx].IsPINDiodeMethodValid()) + if ((*this)[idx].IsExcluded() || !(*this)[idx].IsPINDiodeMethodValid()) return kFALSE; val = (*this)[idx].GetTotalFFactorPINDiodeMethod(); break; case 26: - if ((*this)[idx].IsExcluded() - || !(*fBadPixels)[idx].IsCalibrationSignalOK() - || !(*this)[idx].IsPINDiodeMethodValid()) + if ((*this)[idx].IsExcluded() || !(*this)[idx].IsPINDiodeMethodValid()) return kFALSE; val = (*this)[idx].GetTotalFFactorErrPINDiodeMethod(); @@ -724,41 +632,17 @@ if ((*this)[idx].IsExcluded()) return kFALSE; - if (!(*fBadPixels)[idx].IsCalibrationSignalOK()) - val = 1; - else - return kFALSE; + val = (*this)[idx].GetHiGainNumPickup(); break; case 30: if ((*this)[idx].IsExcluded()) return kFALSE; - if ((*fBadPixels)[idx].IsCalibrationOscillating()) - val = 1; - else - return kFALSE; + val = (*this)[idx].GetLoGainNumPickup(); break; case 31: if ((*this)[idx].IsExcluded()) return kFALSE; - if ((*fBadPixels)[idx].IsCalibrationOscillating()) - val = 1; - else - return kFALSE; + val = (*this)[idx].IsHiGainSaturation(); 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; - val = (*this)[idx].IsHiGainSaturation(); - break; - case 35: if ((*this)[idx].IsExcluded()) return kFALSE; @@ -768,5 +652,5 @@ return kFALSE; break; - case 36: + case 33: if ((*this)[idx].IsExcluded()) return kFALSE; @@ -776,5 +660,5 @@ return kFALSE; break; - case 37: + case 34: if ((*this)[idx].IsExcluded()) return kFALSE; @@ -784,31 +668,31 @@ return kFALSE; break; + case 35: + if ((*this)[idx].IsExcluded()) + return kFALSE; + val = (*this)[idx].GetPed(); + break; + case 36: + if ((*this)[idx].IsExcluded()) + return kFALSE; + val = (*this)[idx].GetPedErr(); + break; + case 37: + if ((*this)[idx].IsExcluded()) + return kFALSE; + val = (*this)[idx].GetPedRms(); + break; case 38: if ((*this)[idx].IsExcluded()) return kFALSE; - val = (*this)[idx].GetPed(); + val = (*this)[idx].GetPedErr()/2.; break; case 39: if ((*this)[idx].IsExcluded()) return kFALSE; - val = (*this)[idx].GetPedErr(); + val = (*this)[idx].GetAbsTimeMean(); break; case 40: if ((*this)[idx].IsExcluded()) - return kFALSE; - val = (*this)[idx].GetPedRms(); - break; - case 41: - if ((*this)[idx].IsExcluded()) - return kFALSE; - val = (*this)[idx].GetPedErr()/2.; - break; - case 42: - if ((*this)[idx].IsExcluded() || !(*fBadPixels)[idx].IsCalibrationSignalOK()) - return kFALSE; - val = (*this)[idx].GetAbsTimeMean(); - break; - case 43: - if ((*this)[idx].IsExcluded() || !(*fBadPixels)[idx].IsCalibrationSignalOK()) return kFALSE; val = (*this)[idx].GetAbsTimeRms(); @@ -817,5 +701,7 @@ return kFALSE; } + return val!=-1.; + } @@ -836,7 +722,8 @@ // // Calculate the weighted mean of the phe's of all inner and outer pixels, respectively. -// Bad pixels are excluded from the calculation. -// -Bool_t MCalibrationChargeCam::CalcMeanFluxPhotonsFFactorMethod() +// Bad pixels are excluded from the calculation. Two loops are performed to exclude pixels +// which are fPheFFactorRelLimit sigmas from the mean. +// +Bool_t MCalibrationChargeCam::CalcMeanFluxPhotonsFFactorMethod(const MGeomCam &geom, const MBadPixelsCam &bad) { @@ -848,4 +735,7 @@ Float_t sumweightsouter = 0.; Float_t sumphesouter = 0.; + Int_t validinner = 0; + Int_t validouter = 0; + TIter Next(fPixels); @@ -859,10 +749,10 @@ const Int_t idx = pix->GetPixId(); - if(!(*fBadPixels)[idx].IsCalibrationResultOK()) + if(!bad[idx].IsCalibrationResultOK()) continue; const Float_t nphe = pix->GetPheFFactorMethod(); const Float_t npheerr = pix->GetPheFFactorMethodErr(); - const Float_t ratio = fGeomCam->GetPixRatio(idx); + const Float_t ratio = geom.GetPixRatio(idx); if (npheerr > 0.) @@ -876,4 +766,5 @@ sumweightsinner += weight; sumphesinner += weight*nphe; + validinner++; } else @@ -885,8 +776,9 @@ sumweightsouter += weight; sumphesouter += weight*nphe; + validouter++; } } /* if npheerr != 0 */ } /* while ((pix=(MCalibrationChargePix*)Next())) */ - + if (sumweightsinner <= 0. || sumphesinner <= 0.) { @@ -914,7 +806,111 @@ } - - const Float_t meanFluxPhotonsRelErrSquare = fMeanFluxPhesInnerPixelErr * fMeanFluxPhesInnerPixelErr - / (fMeanFluxPhesInnerPixel * fMeanFluxPhesInnerPixel); + Float_t meanFluxPhotonsRelErrSquare = fMeanFluxPhesInnerPixelErr * fMeanFluxPhesInnerPixelErr + / (fMeanFluxPhesInnerPixel * fMeanFluxPhesInnerPixel); + + fMeanFluxPhotonsInnerPixel = fMeanFluxPhesInnerPixel/fAverageQE; + fMeanFluxPhotonsInnerPixelErr = TMath::Sqrt(meanFluxPhotonsRelErrSquare + avQERelErrSquare) + * fMeanFluxPhotonsInnerPixel; + + fMeanFluxPhotonsOuterPixel = 4.*fMeanFluxPhotonsInnerPixel; + fMeanFluxPhotonsOuterPixelErr = 4.*fMeanFluxPhotonsInnerPixelErr; + + // + // Here starts the second loop discarting pixels out of the range: + // + const Float_t innererr = TMath::Sqrt((Float_t)validinner)*fPheFFactorRelLimit*fMeanFluxPhesInnerPixelErr; + const Float_t outererr = TMath::Sqrt((Float_t)validouter)*fPheFFactorRelLimit*fMeanFluxPhesOuterPixelErr; + + const Float_t lowerpheinnerlimit = fMeanFluxPhesInnerPixel - innererr; + const Float_t upperpheinnerlimit = fMeanFluxPhesInnerPixel + innererr; + + const Float_t lowerpheouterlimit = fMeanFluxPhesOuterPixel - outererr; + const Float_t upperpheouterlimit = fMeanFluxPhesOuterPixel + outererr; + + sumweightsinner = 0.; + sumphesinner = 0.; + sumweightsouter = 0.; + sumphesouter = 0.; + + TIter Next2(fPixels); + MCalibrationChargePix *pix2; + while ((pix2=(MCalibrationChargePix*)Next2())) + { + + if (!pix2->IsFFactorMethodValid()) + continue; + + const Int_t idx = pix2->GetPixId(); + + if(!bad[idx].IsCalibrationResultOK()) + continue; + + const Float_t nphe = pix2->GetPheFFactorMethod(); + const Float_t npheerr = pix2->GetPheFFactorMethodErr(); + const Float_t ratio = geom.GetPixRatio(idx); + + if (npheerr > 0.) + { + // + // first the inner pixels: + // + if (ratio == 1.) + { + + if (nphe < lowerpheinnerlimit || nphe > upperpheinnerlimit) + { + pix2->SetFFactorMethodValid(kFALSE); + continue; + } + + const Float_t weight = 1./npheerr/npheerr; + sumweightsinner += weight; + sumphesinner += weight*nphe; + } + else + { + // + // now the outers + // + if (nphe < lowerpheouterlimit || nphe > upperpheouterlimit) + { + pix2->SetFFactorMethodValid(kFALSE); + continue; + } + + const Float_t weight = 1./npheerr/npheerr; + sumweightsouter += weight; + sumphesouter += weight*nphe; + } + } /* if npheerr != 0 */ + } /* while ((pix2=(MCalibrationChargePix*)Next2())) */ + + if (sumweightsinner <= 0. || sumphesinner <= 0.) + { + *fLog << warn << " Mean number of phe's from inner pixels cannot be calculated: " + << " Sum of weights: " << sumweightsinner + << " Sum of weighted phes: " << sumphesinner << endl; + return kFALSE; + } + else + { + fMeanFluxPhesInnerPixel = sumphesinner/sumweightsinner; + fMeanFluxPhesInnerPixelErr = TMath::Sqrt(1./sumweightsinner); + + } + + if (sumweightsouter <= 0. || sumphesouter <= 0.) + { + *fLog << warn << " Mean number of phe's from outer pixels cannot be calculated: " + << " Sum of weights or sum of weighted phes is 0. " << endl; + } + else + { + fMeanFluxPhesOuterPixel = sumphesouter/sumweightsouter; + fMeanFluxPhesOuterPixelErr = TMath::Sqrt(1./sumweightsouter); + } + + meanFluxPhotonsRelErrSquare = fMeanFluxPhesInnerPixelErr * fMeanFluxPhesInnerPixelErr + / (fMeanFluxPhesInnerPixel * fMeanFluxPhesInnerPixel); fMeanFluxPhotonsInnerPixel = fMeanFluxPhesInnerPixel/fAverageQE; @@ -931,10 +927,8 @@ << fMeanFluxPhotonsInnerPixel << " +- " << fMeanFluxPhotonsInnerPixelErr << endl; - - return kTRUE; } -void MCalibrationChargeCam::ApplyFFactorCalibration() +void MCalibrationChargeCam::ApplyFFactorCalibration(const MGeomCam &geom, const MBadPixelsCam &bad) { @@ -947,13 +941,16 @@ { - if (!pix->IsFFactorMethodValid()) + const Int_t idx = pix->GetPixId(); + + if (!(bad[idx].IsCalibrationResultOK())) + { + pix->SetFFactorMethodValid(kFALSE); + continue; + } + + if (!(pix->IsFFactorMethodValid())) continue; - const Int_t idx = pix->GetPixId(); - - if(!(*fBadPixels)[idx].IsCalibrationResultOK()) - continue; - - const Float_t ratio = fGeomCam->GetPixRatio(idx); + const Float_t ratio = geom.GetPixRatio(idx); // // Calculate the conversion factor between PHOTONS and FADC counts @@ -1016,5 +1013,5 @@ -void MCalibrationChargeCam::ApplyBlindPixelCalibration() +void MCalibrationChargeCam::ApplyBlindPixelCalibration(const MGeomCam &geom, const MBadPixelsCam &bad) { @@ -1027,37 +1024,38 @@ { - if((*fBadPixels)[pix->GetPixId()].IsCalibrationResultOK()) + const Int_t idx = pix->GetPixId(); + + if (!(bad[idx].IsCalibrationResultOK())) { - - const Int_t idx = pix->GetPixId(); - - const Float_t charge = pix->GetMeanCharge(); - const Float_t area = (*fGeomCam)[idx].GetA(); - const Float_t chargeerr = pix->GetMeanChargeErr(); - - const Float_t nphot = flux * area; - const Float_t nphoterr = fluxerr * area; - const Float_t conversion = nphot/charge; - Float_t conversionerr; - - conversionerr = nphoterr/charge - * nphoterr/charge ; - conversionerr += chargeerr/charge - * chargeerr/charge - * conversion*conversion; - conversionerr = TMath::Sqrt(conversionerr); - - const Float_t conversionsigma = 0.; - - pix->SetConversionBlindPixelMethod(conversion, conversionerr, conversionsigma); - - if (conversionerr/conversion < 0.1) - pix->SetBlindPixelMethodValid(); + pix->SetBlindPixelMethodValid(kFALSE); + continue; } - } -} - - -void MCalibrationChargeCam::ApplyPINDiodeCalibration() + + const Float_t charge = pix->GetMeanCharge(); + const Float_t area = geom[idx].GetA(); + const Float_t chargeerr = pix->GetMeanChargeErr(); + + const Float_t nphot = flux * area; + const Float_t nphoterr = fluxerr * area; + const Float_t conversion = nphot/charge; + Float_t conversionerr; + + conversionerr = nphoterr/charge + * nphoterr/charge ; + conversionerr += chargeerr/charge + * chargeerr/charge + * conversion*conversion; + conversionerr = TMath::Sqrt(conversionerr); + + const Float_t conversionsigma = 0.; + + pix->SetConversionBlindPixelMethod(conversion, conversionerr, conversionsigma); + + if (conversionerr/conversion < 0.1) + pix->SetBlindPixelMethodValid(); + } +} + +void MCalibrationChargeCam::ApplyPINDiodeCalibration(const MGeomCam &geom, const MBadPixelsCam &bad) { @@ -1070,38 +1068,39 @@ { - if((*fBadPixels)[pix->GetPixId()].IsCalibrationResultOK()) + const Int_t idx = pix->GetPixId(); + + if (!(bad[idx].IsCalibrationResultOK())) { - - const Int_t idx = pix->GetPixId(); - - const Float_t charge = pix->GetMeanCharge(); - const Float_t area = (*fGeomCam)[idx].GetA(); - const Float_t chargeerr = pix->GetMeanChargeErr(); - - const Float_t nphot = flux * area; - const Float_t nphoterr = fluxerr * area; - const Float_t conversion = nphot/charge; - - Float_t conversionerr; - - conversionerr = nphoterr/charge - * nphoterr/charge ; - conversionerr += chargeerr/charge - * chargeerr/charge - * conversion*conversion; - if (conversionerr > 0.) - conversionerr = TMath::Sqrt(conversionerr); - - const Float_t conversionsigma = 0.; - - pix->SetConversionPINDiodeMethod(conversion, conversionerr, conversionsigma); - - if (conversionerr/conversion < 0.1) - pix->SetPINDiodeMethodValid(); - + pix->SetPINDiodeMethodValid(kFALSE); + continue; } - } -} - + + const Float_t charge = pix->GetMeanCharge(); + const Float_t area = geom[idx].GetA(); + const Float_t chargeerr = pix->GetMeanChargeErr(); + + const Float_t nphot = flux * area; + const Float_t nphoterr = fluxerr * area; + const Float_t conversion = nphot/charge; + + Float_t conversionerr; + + conversionerr = nphoterr/charge + * nphoterr/charge ; + conversionerr += chargeerr/charge + * chargeerr/charge + * conversion*conversion; + if (conversionerr > 0.) + conversionerr = TMath::Sqrt(conversionerr); + + const Float_t conversionsigma = 0.; + + pix->SetConversionPINDiodeMethod(conversion, conversionerr, conversionsigma); + + if (conversionerr/conversion < 0.1) + pix->SetPINDiodeMethodValid(); + + } +} Index: /trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCam.h =================================================================== --- /trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCam.h (revision 3478) +++ /trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCam.h (revision 3479) @@ -26,12 +26,15 @@ static const Float_t fgConvFFactorRelErrLimit; // The limit for acceptance of the rel. error of the conversion factor with the FFactor method - - Float_t fAverageQE; // The average quantum efficieny (see Class description) - Float_t fAverageQEErr; // The error of the average quantum efficieny (see Class description) + static const Float_t fgPheFFactorRelLimit; // The rel. limit for acceptance of a calculated number of phe's w.r.t the mean number (in sigma of the error) - Float_t fConvFFactorRelErrLimit; // The limit for acceptance of the rel. error of the conversion factor with the FFactor method - + 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 fConvFFactorRelErrLimit; // The limit for acceptance of the rel. error of the conversion factor with the FFactor method + Float_t fPheFFactorRelLimit; // The rel. limit for acceptance of a calculated number of phe's w.r.t the mean number (in sigma of the error). + Int_t fNumPixels; TClonesArray *fPixels; //-> Array of MCalibrationPix with fit results + MCalibrationChargePix *fAverageInnerPix; //-> Average Pixel of all events MCalibrationChargePix *fAverageOuterPix; //-> Average Pixel of all events @@ -42,7 +45,4 @@ const MCalibrationChargeBlindPix *fBlindPixel; //! Pointer to the Blind Pixel with fit results const MCalibrationChargePINDiode *fPINDiode; //! Pointer to the PIN Diode with fit results - - const MGeomCam *fGeomCam; //! Need geom cam to know which pixel in inner or outer - const MBadPixelsCam *fBadPixels; //! Bad Pixels cam, only for reading TH1D* fOffsets; //! @@ -79,8 +79,7 @@ const Float_t err=gkAverageQEErr) { fAverageQE = qe; fAverageQEErr = err; } + void SetNumPixelsExcluded( const UInt_t n ) { fNumExcludedPixels = n; } void SetConvFFactorRelErrLimit( const Float_t f=fgConvFFactorRelErrLimit ) { fConvFFactorRelErrLimit = f; } - void SetNumPixelsExcluded( const UInt_t n ) { fNumExcludedPixels = n; } - void SetGeomCam( const MGeomCam *geom) { fGeomCam = geom; } - void SetBadPixelsCam( const MBadPixelsCam *bad) { fBadPixels = bad; } + void SetPheFFactorRelLimit ( const Float_t f=fgPheFFactorRelLimit ) { fPheFFactorRelLimit = f; } void SetPINDiode ( const MCalibrationChargePINDiode *d ) { fPINDiode = d; } @@ -139,9 +138,9 @@ Bool_t GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type=0) const; - Bool_t CalcMeanFluxPhotonsFFactorMethod(); + Bool_t CalcMeanFluxPhotonsFFactorMethod(const MGeomCam &geom, const MBadPixelsCam &bad); - void ApplyPINDiodeCalibration(); - void ApplyBlindPixelCalibration(); - void ApplyFFactorCalibration(); + void ApplyPINDiodeCalibration(const MGeomCam &geom, const MBadPixelsCam &bad); + void ApplyBlindPixelCalibration(const MGeomCam &geom, const MBadPixelsCam &bad); + void ApplyFFactorCalibration(const MGeomCam &geom, const MBadPixelsCam &bad); ClassDef(MCalibrationChargeCam, 1) // Container for calibration information of the camera Index: /trunk/MagicSoft/Mars/mcalib/MCalibrationChargePix.cc =================================================================== --- /trunk/MagicSoft/Mars/mcalib/MCalibrationChargePix.cc (revision 3478) +++ /trunk/MagicSoft/Mars/mcalib/MCalibrationChargePix.cc (revision 3479) @@ -668,7 +668,4 @@ 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; @@ -703,5 +700,4 @@ return kTRUE; } -