Index: trunk/MagicSoft/Mars/mcalib/CalibLinkDef.h =================================================================== --- trunk/MagicSoft/Mars/mcalib/CalibLinkDef.h (revision 3264) +++ trunk/MagicSoft/Mars/mcalib/CalibLinkDef.h (revision 3265) @@ -8,5 +8,5 @@ #pragma link C++ class MCalibrateData+; -#pragma link C++ class MCalibrationPix+; +#pragma link C++ class MCalibrationChargePix+; #pragma link C++ class MCalibrationBlindPix+; #pragma link C++ class MCalibrationChargeCam+; @@ -15,5 +15,4 @@ #pragma link C++ class MHCalibrationBlindPixel+; -#pragma link C++ class MHCalibrationPixel+; #pragma link C++ class MHCalibrationChargePix+; #pragma link C++ class MHCalibrationChargeHiGainPix+; Index: trunk/MagicSoft/Mars/mcalib/MCalibrate.cc =================================================================== --- trunk/MagicSoft/Mars/mcalib/MCalibrate.cc (revision 3264) +++ trunk/MagicSoft/Mars/mcalib/MCalibrate.cc (revision 3265) @@ -52,5 +52,5 @@ #include "MCalibrationChargeCam.h" -#include "MCalibrationPix.h" +#include "MCalibrationChargePix.h" #include "MExtractedSignalCam.h" @@ -193,5 +193,5 @@ if(fCalibrationMode!=kNone) { - MCalibrationPix &pix = (*fCalibrations)[pixidx]; + MCalibrationChargePix &pix = (*fCalibrations)[pixidx]; if (!pix.IsChargeValid()) Index: trunk/MagicSoft/Mars/mcalib/MCalibrateData.cc =================================================================== --- trunk/MagicSoft/Mars/mcalib/MCalibrateData.cc (revision 3264) +++ trunk/MagicSoft/Mars/mcalib/MCalibrateData.cc (revision 3265) @@ -61,5 +61,5 @@ #include "MCalibrationChargeCam.h" -#include "MCalibrationPix.h" +#include "MCalibrationChargePix.h" #include "MExtractedSignalCam.h" @@ -234,5 +234,5 @@ if(fCalibrationMode!=kNone) { - MCalibrationPix &pix = (*fCalibrations)[pixidx]; + MCalibrationChargePix &pix = (*fCalibrations)[pixidx]; if (!pix.IsChargeValid()) Index: trunk/MagicSoft/Mars/mcalib/MCalibrationPix.cc =================================================================== --- trunk/MagicSoft/Mars/mcalib/MCalibrationPix.cc (revision 3264) +++ (revision ) @@ -1,885 +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-2001 -! -! -\* ======================================================================== */ - -///////////////////////////////////////////////////////////////////////////// -// // -// MCalibrationPix // -// // -// 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 "MCalibrationPix.h" - -#include "MLog.h" -#include "MLogManip.h" - -ClassImp(MCalibrationPix); - -using namespace std; - -const Float_t MCalibrationPix::gkElectronicPedRms = 1.5; -const Float_t MCalibrationPix::gkElectronicPedRmsErr = 0.3; -const Float_t MCalibrationPix::gkFFactor = 1.15; -const Float_t MCalibrationPix::gkFFactorErr = 0.02; -const Float_t MCalibrationPix::gkChargeLimit = 3.; -const Float_t MCalibrationPix::gkChargeErrLimit = 0.; -const Float_t MCalibrationPix::gkChargeRelErrLimit = 1.; -const Float_t MCalibrationPix::gkTimeLimit = 1.5; -const Float_t MCalibrationPix::gkTimeErrLimit = 3.; -const Float_t MCalibrationPix::gkConvFFactorRelErrLimit = 0.1; - -const Float_t MCalibrationPix::gkAverageQE = 0.18; -const Float_t MCalibrationPix::gkAverageQEErr = 0.02; -const Float_t MCalibrationPix::gkConversionHiLo = 10.; -const Float_t MCalibrationPix::gkConversionHiLoErr = 2.5; -// -------------------------------------------------------------------------- -// -// Default Constructor: -// -MCalibrationPix::MCalibrationPix(const char *name, const char *title) - : fPixId(-1), - fFlags(0) -{ - - fName = name ? name : "MCalibrationPixel"; - fTitle = title ? title : "Container of the MHCalibrationPixels and the fit results"; - - fHist = new MHCalibrationPixel("MHCalibrationPixel","Calibration Histograms Pixel "); - - if (!fHist) - *fLog << warn << dbginf << " Could not create MHCalibrationPixel " << endl; - - Clear(); - - // - // At the moment, we don't have a database, yet, - // so we get it from the configuration file - // - SetConversionHiLo(); - SetConversionHiLoErr(); - - SetAverageQE(); -} - -MCalibrationPix::~MCalibrationPix() -{ - delete fHist; -} - - -// ------------------------------------------------------------------------ -// -// Invalidate values -// -void MCalibrationPix::Clear(Option_t *o) -{ - - fHist->Reset(); - - CLRBIT(fFlags, kHiGainSaturation); - CLRBIT(fFlags, kExcluded); - CLRBIT(fFlags, kExcludeQualityCheck); - CLRBIT(fFlags, kChargeValid); - CLRBIT(fFlags, kFitted); - CLRBIT(fFlags, kOscillating); - CLRBIT(fFlags, kBlindPixelMethodValid); - CLRBIT(fFlags, kFFactorMethodValid); - CLRBIT(fFlags, kPINDiodeMethodValid); - CLRBIT(fFlags, kCombinedMethodValid); - - fCharge = -1.; - fChargeErr = -1.; - fSigmaCharge = -1.; - fSigmaChargeErr = -1.; - fRSigmaCharge = -1.; - fRSigmaChargeErr = -1.; - - fChargeProb = -1.; - fPed = -1.; - fPedRms = -1.; - fPedRmsErr = -1.; - - fNumHiGainSamples = -1.; - fNumLoGainSamples = -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 MCalibrationPix::DefinePixId(Int_t i) -{ - - fPixId = i; - fHist->ChangeHistId(i); - -} - - -// -------------------------------------------------------------------------- -// -// Set the pedestals from outside -// -void MCalibrationPix::SetPedestal(const Float_t ped, const Float_t pedrms, - const Float_t higainsamp, const Float_t logainsamp ) -{ - - fPed = ped; - fPedRms = pedrms; - - fNumHiGainSamples = higainsamp; - fNumLoGainSamples = logainsamp; - -} - -// -------------------------------------------------------------------------- -// -// Set the conversion factors from outside (only for MC) -// -void MCalibrationPix::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 MCalibrationPix::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 MCalibrationPix::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 MCalibrationPix::SetConversionPINDiodeMethod(Float_t c, Float_t err, Float_t sig) -{ - fMeanConversionPINDiodeMethod = c ; - fConversionPINDiodeMethodErr = err; - fSigmaConversionPINDiodeMethod = sig; -} - -// -------------------------------------------------------------------------- -// -// Set the Hi Gain Saturation Bit from outside (only for MC) -// -void MCalibrationPix::SetHiGainSaturation(Bool_t b) -{ - - if (b) - { - SETBIT(fFlags, kHiGainSaturation); - fHist->SetUseLoGain(1); - } - else - { - CLRBIT(fFlags, kHiGainSaturation); - fHist->SetUseLoGain(0); - } -} - -// -------------------------------------------------------------------------- -// -// Set the Excluded Bit from outside -// -void MCalibrationPix::SetExcluded(Bool_t b ) -{ - b ? SETBIT(fFlags, kExcluded) : CLRBIT(fFlags, kExcluded); -} - - -// -------------------------------------------------------------------------- -// -// Set the Excluded Bit from outside -// -void MCalibrationPix::SetExcludeQualityCheck(Bool_t b ) -{ - b ? SETBIT(fFlags, kExcludeQualityCheck) : CLRBIT(fFlags, kExcludeQualityCheck); -} - -// -------------------------------------------------------------------------- -// -// Set the Excluded Bit from outside -// -void MCalibrationPix::SetChargeValid(Bool_t b ) -{ - b ? SETBIT(fFlags, kChargeValid) : CLRBIT(fFlags, kChargeValid); -} - -// -------------------------------------------------------------------------- -// -// Set the Excluded Bit from outside -// -void MCalibrationPix::SetFitted(Bool_t b ) -{ - b ? SETBIT(fFlags, kFitted) : CLRBIT(fFlags, kFitted); -} - -// -------------------------------------------------------------------------- -// -// Set the Excluded Bit from outside -// -void MCalibrationPix::SetOscillating(Bool_t b ) -{ - b ? SETBIT(fFlags, kOscillating) : CLRBIT(fFlags, kOscillating); -} - -// -------------------------------------------------------------------------- -// -// Set the Excluded Bit from outside -// -void MCalibrationPix::SetBlindPixelMethodValid(Bool_t b ) -{ - b ? SETBIT(fFlags, kBlindPixelMethodValid) : CLRBIT(fFlags, kBlindPixelMethodValid); -} - -// -------------------------------------------------------------------------- -// -// Set the Excluded Bit from outside -// -void MCalibrationPix::SetFFactorMethodValid(Bool_t b ) -{ - b ? SETBIT(fFlags, kFFactorMethodValid) : CLRBIT(fFlags, kFFactorMethodValid); -} - -// -------------------------------------------------------------------------- -// -// Set the Excluded Bit from outside -// -void MCalibrationPix::SetPINDiodeMethodValid(Bool_t b ) -{ - b ? SETBIT(fFlags, kPINDiodeMethodValid) : CLRBIT(fFlags, kPINDiodeMethodValid); -} - -void MCalibrationPix::SetAbsTimeBordersHiGain(Byte_t f, Byte_t l) -{ - - fTimeFirstHiGain = f; - fTimeLastHiGain = l; - -} - -void MCalibrationPix::SetAbsTimeBordersLoGain(Byte_t f, Byte_t l) -{ - - fTimeFirstLoGain = f; - fTimeLastLoGain = l; - -} - -Bool_t MCalibrationPix::IsExcluded() const -{ - return TESTBIT(fFlags,kExcluded); -} - -Bool_t MCalibrationPix::IsExcludeQualityCheck() const -{ - return TESTBIT(fFlags,kExcludeQualityCheck); -} - -Bool_t MCalibrationPix::IsHiGainSaturation() const -{ - return TESTBIT(fFlags,kHiGainSaturation); -} - -Bool_t MCalibrationPix::IsChargeValid() const -{ - return TESTBIT(fFlags, kChargeValid); -} - -Bool_t MCalibrationPix::IsFitted() const -{ - return TESTBIT(fFlags, kFitted); -} - -Bool_t MCalibrationPix::IsOscillating() const -{ - return TESTBIT(fFlags, kOscillating); -} - -Bool_t MCalibrationPix::IsBlindPixelMethodValid() const -{ - return TESTBIT(fFlags, kBlindPixelMethodValid); -} - -Bool_t MCalibrationPix::IsFFactorMethodValid() const -{ - return TESTBIT(fFlags, kFFactorMethodValid); -} - -Bool_t MCalibrationPix::IsPINDiodeMethodValid() const -{ - return TESTBIT(fFlags, kPINDiodeMethodValid); -} - -Bool_t MCalibrationPix::IsCombinedMethodValid() const -{ - return TESTBIT(fFlags, kCombinedMethodValid); -} - - -// -------------------------------------------------------------------------- -// -// 1) Return if the charge distribution is already succesfully fitted -// or if the histogram is empty -// 2) Set a lower Fit range according to 1.5 Pedestal RMS in order to avoid -// possible remaining cosmics to spoil the fit. -// 3) Decide if the LoGain Histogram is fitted or the HiGain Histogram -// 4) Fit the histograms with a Gaussian -// 5) In case of failure set the bit kFitted to false -// 6) Retrieve the results and store them in this class -// 7) Calculate the number of photo-electrons after the F-Factor method -// 8) Calculate the errors of the F-Factor method -// -// The fits are declared valid (fFitValid = kTRUE), if: -// -// 1) Pixel has a fitted charge greater than 3*PedRMS -// 2) Pixel has a fit error greater than 0. -// 3) Pixel has a fit Probability greater than 0.0001 -// 4) Pixel has a charge sigma bigger than its Pedestal RMS -// 5) If FitTimes is used, -// the mean arrival time is at least 1.0 slices from the used edge slices -// (this stage is only performed in the times fit) -// -// If the histogram is empty, all values are set to -1. -// -// The conversion factor after the F-Factor method is declared valid, if: -// -// 1) fFitValid is kTRUE -// 2) Conversion Factor is bigger than 0. -// 3) The error of the conversion factor is smaller than 10% -// -Bool_t MCalibrationPix::FitCharge() -{ - - // - // 1) Return if the charge distribution is already succesfully fitted - // or if the histogram is empty - // - if (fHist->IsChargeFitOK() || fHist->IsEmpty()) - return kTRUE; - - // - // 2) Set a lower Fit range according to 1.5 Pedestal RMS in order to avoid - // possible remaining cosmics to spoil the fit. - // - // if (fPed && fPedRms) - // fHist->SetLowerFitRange(1.5*fPedRms); - // else - // *fLog << warn << "WARNING: Cannot set lower fit range: Pedestals not available" << endl; - - // - // 3) Decide if the LoGain Histogram is fitted or the HiGain Histogram - // - if (fHist->UseLoGain()) - SetHiGainSaturation(); - - // - // 4) Fit the Lo Gain histograms with a Gaussian - // - if (fHist->FitCharge()) - SETBIT(fFlags,kFitted); - else - { - *fLog << warn << "WARNING: Could not fit charges of pixel " << fPixId << endl; - // - // 5) In case of failure set the bit kFitted to false - // - CLRBIT(fFlags,kFitted); - } - - // - // 6) Retrieve the results and store them in this class - // If fFitted is false, we get the eans and RMS of the histogram!! - // - fCharge = fHist->GetChargeMean(); - fChargeErr = fHist->GetChargeMeanErr(); - fSigmaCharge = fHist->GetChargeSigma(); - fSigmaChargeErr = fHist->GetChargeSigmaErr(); - fChargeProb = fHist->GetChargeProb(); - - - fAbsTimeMean = fHist->GetAbsTimeMean(); - fAbsTimeMeanErr = fHist->GetAbsTimeMeanErr(); - fAbsTimeRms = fHist->GetAbsTimeRms(); - - if (CheckTimeFitValidity()) - SETBIT(fFlags,kTimeFitValid); - else - CLRBIT(fFlags,kTimeFitValid); - - // - //Calculate the conversion factors - // - if (IsHiGainSaturation()) - ApplyLoGainConversion(); - - if (CheckChargeValidity()) - SETBIT(fFlags,kChargeValid); - else - { - CLRBIT(fFlags,kChargeValid); - return kFALSE; - } - - return kTRUE; - -} - -// -// Calculate the number of photo-electrons after the F-Factor method -// Calculate the errors of the F-Factor method -// -Bool_t MCalibrationPix::CalcFFactorMethod() -{ - - if ( (fCharge == -1.) - || (fChargeErr < 0.) - || (fSigmaCharge < 0.) - || (fPedRms < 0.) ) - { - *fLog << warn << GetDescriptor() << "Cannot calculate the FFactor Method! " - << "Some of the needed parameters are not available "; - CLRBIT(fFlags,kFFactorMethodValid); - return kFALSE; - } - - // - // Square all variables in order to avoid applications of square root - // - // First the relative error squares - // - const Float_t chargeSquare = fCharge * fCharge; - const Float_t chargeSquareRelErrSquare = 4.* fChargeErr * fChargeErr / chargeSquare; - - const Float_t chargeRelErrSquare = fChargeErr * fChargeErr - / (fCharge * fCharge); - - 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) ; - // - // Now the absolute error squares - // - const Float_t sigmaSquare = fSigmaCharge * fSigmaCharge; - const Float_t sigmaSquareErrSquare = 4.*fSigmaChargeErr* fSigmaChargeErr * sigmaSquare; - - Float_t pedRmsSquare = fPedRms * fPedRms; - Float_t pedRmsSquareErrSquare = 4.*fPedRmsErr * fPedRmsErr * pedRmsSquare; - - if (!IsHiGainSaturation()) - { /* HiGain */ - - pedRmsSquare *= fNumHiGainSamples; - pedRmsSquareErrSquare *= fNumHiGainSamples*fNumHiGainSamples; - } - else - { /* LoGain */ - - // - // 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 = gkElectronicPedRms * gkElectronicPedRms; - const Float_t elecRmsSquareErrSquare = 4.*gkElectronicPedRmsErr * gkElectronicPedRmsErr * 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; - pedRmsSquareErrSquare = convertedNsbSquareErrSquare + elecRmsSquareErrSquare; - - // - // Now, correct for the number of used FADC slices in the LoGain: - // - pedRmsSquare *= fNumLoGainSamples; - pedRmsSquareErrSquare *= fNumLoGainSamples*fNumLoGainSamples; - // - // Correct also for the conversion to Hi-Gain: - // - pedRmsSquare *= fConversionHiLo*fConversionHiLo; - pedRmsSquareErrSquare *= fConversionHiLo*fConversionHiLo*fConversionHiLo*fConversionHiLo; - - } /* if (HiGainSaturation) */ - - // - // Calculate the reduced sigmas - // - const Float_t rsigmachargesquare = sigmaSquare - pedRmsSquare; - if (rsigmachargesquare <= 0.) - { - *fLog << warn - << "WARNING: Cannot apply F-Factor calibration: Reduced Sigma smaller than 0 in pixel " - << fPixId << endl; - CLRBIT(fFlags,kFFactorMethodValid); - return kFALSE; - } - - const Float_t rSigmaSquareRelErrSquare = (sigmaSquareErrSquare + pedRmsSquareErrSquare) - / (rsigmachargesquare * rsigmachargesquare) ; - - fRSigmaCharge = TMath::Sqrt(rsigmachargesquare); - fRSigmaChargeErr = TMath::Sqrt(sigmaSquareErrSquare + pedRmsSquareErrSquare); - - - // - // Calculate the number of phe's from the F-Factor method - // (independent on Hi Gain or Lo Gain) - // - fPheFFactorMethod = ffactorsquare * chargeSquare / rsigmachargesquare; - // - // 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 / fCharge ; - fConversionFFactorMethodErr = ( pheFFactorRelErrSquare + chargeRelErrSquare ) - * fMeanConversionFFactorMethod * fMeanConversionFFactorMethod; - - const Float_t convrelerror = fConversionFFactorMethodErr - / fMeanConversionFFactorMethod; - - if ( (fMeanConversionFFactorMethod > 0.) && (convrelerror < gkConvFFactorRelErrLimit)) - SETBIT(fFlags,kFFactorMethodValid); - - fSigmaConversionFFactorMethod = GetTotalFFactorFFactorMethod()*TMath::Sqrt(fMeanConversionFFactorMethod); - - // - // Calculate the Total F-Factor of the camera ( in photons ) - // - if (fPheFFactorMethod > 0) - { - fTotalFFactorFFactorMethod = (fRSigmaCharge/fCharge)*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; - - return kTRUE; -} - - -// -// The check returns kTRUE if: -// -// 0) Pixel has BIT fitted set: -// This means: -// a) No result is a nan -// b) The NDF is not smaller than fNDFLimit (5) -// c) The Probability is greater than gkProbLimit (default 0.001 == 99.9%) -// 1) Pixel has a fitted charge greater than 3*PedRMS -// 2) Pixel has a fit error greater than 0. -// 3) Pixel has a fitted charge greater its charge error -// 4) Pixel has a fit Probability greater than 0.0001 -// 5) Pixel has a charge sigma bigger than its Pedestal RMS -// -Bool_t MCalibrationPix::CheckChargeValidity() -{ - - if (!IsFitted()) - return kFALSE; - - if (IsExcludeQualityCheck()) - return kTRUE; - - Float_t pedestal; - - if (!IsHiGainSaturation()) /* higain */ - pedestal = GetPedRms()*TMath::Sqrt(fNumHiGainSamples); - else /* logain */ - pedestal = GetPedRms()*TMath::Sqrt(fNumLoGainSamples); - - - if (fCharge < gkChargeLimit*pedestal) - { - *fLog << warn << "WARNING: Fitted Charge is smaller than " - << gkChargeLimit << " Pedestal RMS in Pixel " << fPixId << endl; - return kFALSE; - } - - if (fChargeErr < gkChargeErrLimit) - { - *fLog << warn << "WARNING: Err of Fitted Charge is smaller than " - << gkChargeErrLimit << " in Pixel " << fPixId << endl; - return kFALSE; - } - - if (fCharge < gkChargeRelErrLimit*fChargeErr) - { - *fLog << warn << "WARNING: Fitted Charge is smaller than " - << gkChargeRelErrLimit << "* its error in Pixel " << fPixId << endl; - return kFALSE; - } - - if (!fHist->IsChargeFitOK()) - { - *fLog << warn << "WARNING: Probability of Fitted Charge too low in Pixel " - << fPixId << endl; - return kFALSE; - } - - if (fSigmaCharge < pedestal) - { - *fLog << warn << "WARNING: Sigma of Fitted Charge smaller than Pedestal RMS in Pixel " - << fPixId << endl; - return kFALSE; - } - return kTRUE; -} - -// -// The check return kTRUE if: -// -// 0) No value is nan -// 1) Pixel has a fitted rel. time smaller than 3*FADC slices -// 2) Pixel has a fit error greater than 0. -// 4) Pixel has a fit Probability greater than 0.001 -// 5) The absolute arrival time is at least 1.0 slices from the used edge slices -// -Bool_t MCalibrationPix::CheckTimeFitValidity() -{ - - - if (IsExcludeQualityCheck()) - return kTRUE; - - if (IsHiGainSaturation()) - { - - if (fAbsTimeMean < (Float_t)fTimeFirstLoGain+1) - { - *fLog << warn - << "WARNING: Some absolute times smaller than limit in Pixel " - << fPixId << " time: " << fAbsTimeMean - << " Limit: " << (Float_t)fTimeFirstLoGain+1. << endl; - return kFALSE; - } - - if (fAbsTimeMean > (Float_t)fTimeLastLoGain-1) - { - *fLog << warn - << "WARNING: Some absolute times bigger than limit in Pixel " - << fPixId << " time: " << fAbsTimeMean - << " Limit: " << (Float_t)fTimeLastLoGain-1. << endl; - return kFALSE; - } - - } - else - { - - if (fAbsTimeMean < (Float_t)fTimeFirstHiGain+1.) - { - *fLog << warn - << "WARNING: Some absolute times smaller than limit in Pixel " - << fPixId << " time: " << fAbsTimeMean - << " Limit: " << (Float_t)fTimeFirstHiGain+1. << endl; - // return kFALSE; - } - - if (fAbsTimeMean > (Float_t)fTimeLastHiGain-1.) - { - *fLog << warn - << "WARNING: Some absolute times bigger than limit in Pixel " - << fPixId << " time: " << fAbsTimeMean - << " Limit: " << (Float_t)fTimeLastHiGain-1. << endl; - // return kFALSE; - } - - } - - - - return kTRUE; -} - - -void MCalibrationPix::CheckOscillations() -{ - fHist->CheckOscillations(); -} - -void MCalibrationPix::ApplyLoGainConversion() -{ - - const Float_t chargeRelErrSquare = fChargeErr * fChargeErr - /( fCharge * fCharge ); - const Float_t sigmaRelErrSquare = fSigmaChargeErr * fSigmaChargeErr - /( fSigmaCharge * fSigmaCharge ); - const Float_t conversionRelErrSquare = fConversionHiLoErr * fConversionHiLoErr - /( fConversionHiLo * fConversionHiLo ); - - fCharge *= fConversionHiLo; - fChargeErr = TMath::Sqrt(chargeRelErrSquare + conversionRelErrSquare) * fCharge; - - fSigmaCharge *= fConversionHiLo; - fSigmaChargeErr = TMath::Sqrt(sigmaRelErrSquare + conversionRelErrSquare) * fSigmaCharge; - -} Index: trunk/MagicSoft/Mars/mcalib/MCalibrationPix.h =================================================================== --- trunk/MagicSoft/Mars/mcalib/MCalibrationPix.h (revision 3264) +++ (revision ) @@ -1,240 +1,0 @@ -#ifndef MARS_MCalibrationPix -#define MARS_MCalibrationPix - -#ifndef MARS_MParContainer -#include "MParContainer.h" -#endif - -#include "MHCalibrationPixel.h" - -class MCalibrationPix : 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 gkChargeLimit; // The limit (in units of PedRMS) for acceptance of the fitted mean charge - static const Float_t gkChargeErrLimit; // The limit (in units of PedRMS) for acceptance of the fitted charge sigma - static const Float_t gkChargeRelErrLimit; // The limit (in units of Error of fitted charge) for acceptance of the fitted mean - static const Float_t gkTimeLimit; // The limit (in units of FADC slices) for acceptance of the fitted time - static const Float_t gkTimeErrLimit; // The limit (in units of FADC slices) for acceptance of the fitted time sigma - static const Float_t gkConvFFactorRelErrLimit; // The limit for acceptance of the rel. error of the conversion factor with the FFactor method - - 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 fCharge; // The mean reduced charge after the fit - Float_t fChargeErr; // The error of reduced mean charge after the fit - Float_t fSigmaCharge; // The sigma of the mean charge after the fit - Float_t fSigmaChargeErr; // The error of the sigma of the mean charge after the fit - Float_t fRSigmaCharge; // The reduced squares of sigmas after the fit - Float_t fRSigmaChargeErr; // The reduced squares of sigmas after the fit - Float_t fChargeProb; // The probability of the fit function - - Float_t fPed; // The mean pedestal (from MPedestalPix) - Float_t fPedRms; // The pedestal RMS (from MPedestalPix) - Float_t fPedRmsErr; // The error of the pedestal RMS (from MPedestalPix) - - Float_t fAbsTimeMean; // The mean absolute arrival time - Float_t fAbsTimeMeanErr; // The error of the absolute mean 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 fNumHiGainSamples; // The number of hi-gain samples used for the signal extraction - Float_t fNumLoGainSamples; // The number of hi-gain samples used for the signal extraction - - enum { kHiGainSaturation, - kExcluded, kExcludeQualityCheck, - kChargeValid, kTimeFitValid, - kFitted, kOscillating, - kBlindPixelMethodValid, kFFactorMethodValid, - kPINDiodeMethodValid, kCombinedMethodValid }; - - MHCalibrationPixel *fHist; // Pointer to the histograms performing the fits, etc. - - Bool_t CheckChargeValidity(); - Bool_t CheckTimeFitValidity(); - -public: - - MCalibrationPix(const char *name=NULL, const char *title=NULL); - ~MCalibrationPix(); - - void Clear(Option_t *o=""); - - // Getter - MHCalibrationPixel *GetHist() const { return fHist; } - - // Setter - void SetPedestal(const Float_t ped, const Float_t pedrms, - const Float_t higainsamp, const Float_t logainsamp); - 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; } - - // 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( Bool_t b = kTRUE ); - void SetExcluded( Bool_t b = kTRUE ); - void SetExcludeQualityCheck( Bool_t b = kTRUE ); - void SetChargeValid( Bool_t b = kTRUE ); - void SetFitted( Bool_t b = kTRUE ); - void SetOscillating( Bool_t b = kTRUE ); - void SetBlindPixelMethodValid( Bool_t b = kTRUE ); - void SetFFactorMethodValid( Bool_t b = kTRUE ); - void SetPINDiodeMethodValid( Bool_t b = kTRUE ); - - void SetAbsTimeBordersHiGain( Byte_t f, Byte_t l ); - void SetAbsTimeBordersLoGain( Byte_t f, Byte_t l ); - - // Charges - Float_t GetCharge() const { return fCharge; } - Float_t GetChargeErr() const { return fChargeErr; } - Float_t GetChargeProb() const { return fChargeProb; } - Float_t GetSigmaCharge() const { return fSigmaCharge; } - Float_t GetSigmaChargeErr() const { return fSigmaChargeErr; } - Float_t GetRSigmaCharge() const { return fRSigmaCharge; } - Float_t GetRSigmaChargeErr() const { return fRSigmaChargeErr; } - - - Float_t GetAbsTimeMean() const { return fAbsTimeMean; } - Float_t GetAbsTimeMeanErr() const { return fAbsTimeMeanErr; } - 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 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 IsExcludeQualityCheck() const; - Bool_t IsHiGainSaturation() const; - Bool_t IsChargeValid() const; - Bool_t IsFitted() const; - Bool_t IsOscillating() const; - Bool_t IsBlindPixelMethodValid() const; - Bool_t IsPINDiodeMethodValid() const; - Bool_t IsFFactorMethodValid() const; - Bool_t IsCombinedMethodValid() const; - - // Fill histos - Bool_t FillChargeHiGain(Float_t q) const { return fHist->FillChargeHiGain(q); } - Bool_t FillAbsTimeHiGain(Float_t t) const { return fHist->FillAbsTimeHiGain(t); } - - Bool_t FillChargeLoGain(Float_t q) const { return fHist->FillChargeLoGain(q); } - Bool_t FillAbsTimeLoGain(Float_t t) const { return fHist->FillAbsTimeLoGain(t); } - - Bool_t FillGraphs(Float_t qhi,Float_t qlo) const { return fHist->FillGraphs(qhi,qlo); } - - void DefinePixId(Int_t i); - - // Fits - Bool_t FitCharge(); - - // Draws - void Draw(Option_t *opt="") { fHist->Draw(opt); } - TObject *DrawClone(Option_t *opt="") const { return fHist->DrawClone(opt); } - - // Miscellaneous - void ApplyLoGainConversion(); - void CheckOscillations(); - Bool_t CalcFFactorMethod(); - - ClassDef(MCalibrationPix, 1) // Container for Calibration of one pixel -}; - -#endif - Index: trunk/MagicSoft/Mars/mcalib/MMcCalibrationCalc.cc =================================================================== --- trunk/MagicSoft/Mars/mcalib/MMcCalibrationCalc.cc (revision 3264) +++ trunk/MagicSoft/Mars/mcalib/MMcCalibrationCalc.cc (revision 3265) @@ -47,5 +47,5 @@ #include "MParList.h" -#include "MCalibrationPix.h" +#include "MCalibrationChargePix.h" #include "MCalibrationChargeCam.h" @@ -246,5 +246,5 @@ for (int i=0; i