/* ======================================================================== *\ ! ! * ! * 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 ! ! \* ======================================================================== */ ///////////////////////////////////////////////////////////////////////////// // // MHCalibrationChargeCam // // Fills the extracted signals of MExtractedSignalCam into the MHCalibrationPix-classes // MHCalibrationChargeHiGainPix and MHCalibrationChargeLoGainPix for every: // // - Pixel, stored in the TObjArray's MHCalibrationCam::fHiGainArray and // MHCalibrationCam::fLoGainArray // // - Average pixel per AREA index (e.g. inner and outer for the MAGIC camera), // stored in the TObjArray's MHCalibrationCam::fAverageHiGainAreas and // MHCalibrationCam::fAverageLoGainAreas // // - Average pixel per camera SECTOR (e.g. sectors 1-6 for the MAGIC camera), // stored in the TObjArray's MHCalibrationCam::fAverageHiGainSectors and // MHCalibrationCam::fAverageLoGainSectors // // Every signal is taken from MExtractedSignalCam and filled into a histogram and // an array, in order to perform a Fourier analysis (see MHGausEvents). // The signals are moreover averaged on an event-by-event basis and written into // the corresponding average pixels. // // Additionally, the (FADC slice) position of the maximum is stored in an Absolute // Arrival Time histogram. This histogram serves for a rough cross-check if the // signal does not lie at or outside the edges of the extraction window. // // The Charge histograms are fitted to a Gaussian, mean and sigma with its errors // and the fit probability are extracted. If none of these values are NaN's and // if the probability is bigger than MHGausEvents::fProbLimit (default: 0.5%), // the fit is declared valid. // Otherwise, the fit is repeated within ranges of the previous mean // +- MHCalibrationPix::fPickupLimit (default: 5) sigma (see MHCalibrationPix::RepeatFit()) // In case this does not make the fit valid, the histogram means and RMS's are // taken directly (see MHCalibrationPix::BypassFit()) and the following flags are set: // - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kHiGainNotFitted ) or // - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainNotFitted ) and // - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun ) // // Outliers of more than MHCalibrationPix::fPickupLimit (default: 5) sigmas // from the mean are counted as Pickup events (stored in MHCalibrationPix::fPickup) // // Unless more than fNumHiGainSaturationLimit (default: 1%) of the overall FADC // slices show saturation, the following flag is set: // - MCalibrationChargePix::SetHiGainSaturation(); // In that case, the calibration constants are derived from the low-gain results. // // If more than fNumLoGainSaturationLimit (default: 1%) of the overall // low-gain FADC slices saturate, the following flags are set: // - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainSaturation ) and // - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnsuitableRun ) // // The class also fills arrays with the signal vs. event number, creates a fourier // spectrum and investigates if the projected fourier components follow an exponential // distribution. In case that the probability of the exponential fit is less than // MHGausEvents::fProbLimit (default: 0.5%), the following flags are set: // - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kHiGainOscillating ) or // - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainOscillating ) and // - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun ) // // This same procedure is performed for the average pixels. // // The following results are written into MCalibrationChargeCam: // // - MCalibrationPix::SetHiGainSaturation() // - MCalibrationPix::SetHiGainMean() // - MCalibrationPix::SetHiGainMeanErr() // - MCalibrationPix::SetHiGainSigma() // - MCalibrationPix::SetHiGainSigmaErr() // - MCalibrationPix::SetHiGainProb() // - MCalibrationPix::SetHiGainNumPickup() // // - MCalibrationPix::SetLoGainMean() // - MCalibrationPix::SetLoGainMeanErr() // - MCalibrationPix::SetLoGainSigma() // - MCalibrationPix::SetLoGainSigmaErr() // - MCalibrationPix::SetLoGainProb() // - MCalibrationPix::SetLoGainNumPickup() // // - MCalibrationChargePix::SetAbsTimeMean() // - MCalibrationChargePix::SetAbsTimeRms() // // For all averaged areas, the fitted sigma is multiplied with the square root of // the number involved pixels in order to be able to compare it to the average of // sigmas in the camera. // ///////////////////////////////////////////////////////////////////////////// #include "MHCalibrationChargeCam.h" #include "MHCalibrationCam.h" #include "MLog.h" #include "MLogManip.h" #include "MParList.h" #include "MHCalibrationChargePix.h" #include "MHCalibrationPix.h" #include "MCalibrationIntensityCam.h" #include "MCalibrationChargeCam.h" #include "MCalibrationChargePix.h" #include "MGeomCam.h" #include "MGeomPix.h" #include "MBadPixelsCam.h" #include "MBadPixelsPix.h" #include "MRawEvtData.h" #include "MRawRunHeader.h" #include "MRawEvtPixelIter.h" #include "MExtractedSignalCam.h" #include "MExtractedSignalPix.h" #include "MArrayI.h" #include "MArrayD.h" #include #include #include #include #include #include #include #include #include ClassImp(MHCalibrationChargeCam); using namespace std; const Axis_t MHCalibrationChargeCam::fgChargeHiGainFirst = -100.5; const Axis_t MHCalibrationChargeCam::fgChargeHiGainLast = 999.5; const Axis_t MHCalibrationChargeCam::fgChargeLoGainFirst = -150.5; const Axis_t MHCalibrationChargeCam::fgChargeLoGainLast = 499.5; const Float_t MHCalibrationChargeCam::fgNumHiGainSaturationLimit = 0.01; const Float_t MHCalibrationChargeCam::fgNumLoGainSaturationLimit = 0.005; const Float_t MHCalibrationChargeCam::fgTimeLowerLimit = 1.; const Float_t MHCalibrationChargeCam::fgTimeUpperLimit = 2.; // 1Led Green, 1 LED blue, 5 LEDs blue, 10 LEDs blue, 10 LEDs UV, CT1, 5Leds Green const Float_t MHCalibrationChargeCam::gkHiGainInnerRefLines[7] = { 245., 323. , 1065., 1467., 180., 211. , 533.5}; const Float_t MHCalibrationChargeCam::gkHiGainOuterRefLines[7] = { 217., 307.5, 932. , 1405., 167., 183.5, 405.5}; const Float_t MHCalibrationChargeCam::gkLoGainInnerRefLines[7] = { 20.8, 28.0 , 121. , 200.2, 16.5, 13.5 , 41.7 }; const Float_t MHCalibrationChargeCam::gkLoGainOuterRefLines[7] = { 18.9, 26.0 , 108.3, 198. , 14.0, 11. , 42. }; // -------------------------------------------------------------------------- // // Default Constructor. // // Sets: // - all pointers to NULL // // Initializes: // - fNumHiGainSaturationLimit to fgNumHiGainSaturationLimit // - fNumLoGainSaturationLimit to fgNumLoGainSaturationLimit // - fTimeLowerLimit to fgTimeLowerLimit // - fTimeUpperLimit to fgTimeUpperLimit // MHCalibrationChargeCam::MHCalibrationChargeCam(const char *name, const char *title) : fRawEvt(NULL) { fName = name ? name : "MHCalibrationChargeCam"; fTitle = title ? title : "Class to fill the calibration histograms "; SetNumHiGainSaturationLimit(fgNumHiGainSaturationLimit); SetNumLoGainSaturationLimit(fgNumLoGainSaturationLimit); SetTimeLowerLimit(); SetTimeUpperLimit(); } // -------------------------------------------------------------------------- // // Gets the pointers to: // - MRawEvtData // Bool_t MHCalibrationChargeCam::SetupHists(const MParList *pList) { fRawEvt = (MRawEvtData*)pList->FindObject("MRawEvtData"); if (!fRawEvt) { *fLog << err << dbginf << "MRawEvtData not found... aborting." << endl; return kFALSE; } return kTRUE; } // -------------------------------------------------------------------------- // // Gets or creates the pointers to: // - MExtractedSignalCam // - MCalibrationChargeCam // - MBadPixelsCam // // Initializes the number of used FADC slices from MExtractedSignalCam // into MCalibrationChargeCam and test for changes in that variable // // Initializes, if empty to MGeomCam::GetNumPixels(): // - MHCalibrationCam::fHiGainArray, MHCalibrationCam::fLoGainArray // - MHCalibrationCam::fHiGainOverflow, MHCalibrationCam::fLoGainOverflow // // Initializes, if empty to MGeomCam::GetNumAreas() for: // - MHCalibrationCam::fAverageHiGainAreas, MHCalibrationCam::fAverageLoGainAreas // // Initializes, if empty to MGeomCam::GetNumSectors() for: // - MHCalibrationCam::fAverageHiGainSectors, MHCalibrationCam::fAverageLoGainSectors // // Calls MHCalibrationCam::InitHists() for every entry in: // - MHCalibrationCam::fHiGainArray, MHCalibrationCam::fLoGainArray // - MHCalibrationCam::fAverageHiGainAreas, MHCalibrationCam::fAverageLoGainAreas // - MHCalibrationCam::fAverageHiGainSectors, MHCalibrationCam::fAverageLoGainSectors // // Sets Titles and Names for the Charge Histograms: // - MHCalibrationCam::fAverageHiGainAreas // - MHCalibrationCam::fAverageHiGainSectors // Bool_t MHCalibrationChargeCam::ReInitHists(MParList *pList) { MExtractedSignalCam *signal = (MExtractedSignalCam*)pList->FindObject(AddSerialNumber("MExtractedSignalCam")); if (!signal) { *fLog << err << "MExtractedSignalCam not found... abort." << endl; return kFALSE; } fIntensCam = (MCalibrationIntensityCam*)pList->FindObject(AddSerialNumber("MCalibrationIntensityChargeCam")); if (fIntensCam) *fLog << inf << "Found MCalibrationIntensityChargeCam ... " << endl; else { fCam = (MCalibrationCam*)pList->FindObject(AddSerialNumber("MCalibrationChargeCam")); if (!fCam) { fCam = (MCalibrationCam*)pList->FindCreateObj(AddSerialNumber("MCalibrationChargeCam")); if (!fCam) { *fLog << err << "Cannot find nor create MCalibrationChargeCam ... abort." << endl; return kFALSE; } fCam->Init(*fGeom); } } fFirstHiGain = signal->GetFirstUsedSliceHiGain(); fLastHiGain = signal->GetLastUsedSliceHiGain(); fFirstLoGain = signal->GetFirstUsedSliceLoGain(); fLastLoGain = signal->GetLastUsedSliceLoGain(); /* const Float_t numhigain = signal->GetNumUsedHiGainFADCSlices(); const Float_t numlogain = signal->GetNumUsedLoGainFADCSlices(); if (fCam) { if (fCam->GetNumHiGainFADCSlices() == 0.) fCam->SetNumHiGainFADCSlices ( numhigain ); else if (fCam->GetNumHiGainFADCSlices() != numhigain) { *fLog << err << GetDescriptor() << ": Number of High Gain FADC extraction slices has changed, abort..." << endl; return kFALSE; } if (fCam->GetNumLoGainFADCSlices() == 0.) fCam->SetNumLoGainFADCSlices ( numlogain ); else if (fCam->GetNumLoGainFADCSlices() != numlogain) { *fLog << err << GetDescriptor() << ": Number of Low Gain FADC extraction slices has changes, abort..." << endl; return kFALSE; } } */ const Int_t npixels = fGeom->GetNumPixels(); const Int_t nsectors = fGeom->GetNumSectors(); const Int_t nareas = fGeom->GetNumAreas(); const Int_t higainsamples = fRunHeader->GetNumSamplesHiGain(); const Int_t logainsamples = fRunHeader->GetNumSamplesLoGain(); if (fHiGainArray->GetEntries()==0) { fHiGainArray->Expand(npixels); for (Int_t i=0; iSetName ("HChargeHiGainPix"); pix.GetHGausHist()->SetTitle("Signals High Gain Pixel "); pix.GetHAbsTime()->SetName ("HAbsTimeHiGainPix"); pix.GetHAbsTime()->SetTitle("Absolute Arrival Times High Gain Pixel "); InitHists((*this)[i],(*fBadPixels)[i],i); } } if (fLoGainArray->GetEntries()==0 && IsLoGain()) { fLoGainArray->Expand(npixels); for (Int_t i=0; iSetName ("HChargeLoGainPix"); pix.GetHGausHist()->SetTitle("Signals Low Gain Pixel "); pix.GetHAbsTime()->SetName ("HAbsTimeLoGainPix"); pix.GetHAbsTime()->SetTitle("Absolute Arrival Times Low Gain Pixel "); // // Adapt the range for the case, the intense blue is used: // FIXME: this is a nasty workaround, but for the moment necessary // in order to avoid default memory space. // if (fGeom->InheritsFrom("MGeomCamMagic")) { if ( fColor == MCalibrationCam::kBLUE) { pix.SetLast(999.5); pix.SetNbins((Int_t)(999.5-fgChargeLoGainFirst)); } } InitHists(pix,(*fBadPixels)[i],i); } } if (fAverageHiGainAreas->GetEntries()==0) { fAverageHiGainAreas->Expand(nareas); for (Int_t j=0; jSetName ("HChargeHiGainAreaIdx"); pix.GetHGausHist()->SetTitle("Signals averaged on event-by-event basis High Gain Area Idx "); pix.GetHAbsTime()->SetName ("HAbsTimeHiGainAreaIdx"); pix.GetHAbsTime()->SetTitle("Absolute Arrival Times average Hi Gain Area Idx "); if (fGeom->InheritsFrom("MGeomCamMagic")) { pix.GetHGausHist()->SetTitle(Form("%s%s%s","Signals averaged on event-by-event basis ", j==0 ? "Inner Pixels " : "Outer Pixels ","High Gain Runs: ")); pix.InitBins(); pix.SetEventFrequency(fPulserFrequency); } else { InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadArea(j) : fCam->GetAverageBadArea(j),j); } } } if (fAverageLoGainAreas->GetEntries()==0 && IsLoGain()) { fAverageLoGainAreas->Expand(nareas); for (Int_t j=0; jSetName ("HChargeLoGainAreaIdx"); pix.GetHGausHist()->SetTitle("Signals averaged on event-by-event basis Low Gain Area Idx "); pix.GetHAbsTime()->SetName ("HAbsTimeLoGainAreaIdx"); pix.GetHAbsTime()->SetTitle("Absolute Arrival Times average Low Gain Area Idx "); // // Adapt the range for the case, the intense blue is used: // FIXME: this is a nasty workaround, but for the moment necessary // in order to avoid default memory space. // if (fGeom->InheritsFrom("MGeomCamMagic")) { if ( fColor == MCalibrationCam::kBLUE) { pix.SetFirst(-10.5); pix.SetLast(999.5); pix.SetNbins(3030); } } if (fGeom->InheritsFrom("MGeomCamMagic")) { pix.GetHGausHist()->SetTitle(Form("%s%s%s","Signals averaged on event-by-event basis ", j==0 ? "Inner Pixels " : "Outer Pixels ","Low Gain Runs: ")); pix.InitBins(); pix.SetEventFrequency(fPulserFrequency); } else { InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadArea(j) : fCam->GetAverageBadArea(j),j); } } } if (fAverageHiGainSectors->GetEntries()==0) { fAverageHiGainSectors->Expand(nsectors); for (Int_t j=0; jSetName ("HChargeHiGainPix"); pix.GetHGausHist()->SetTitle("Signals averaged on event-by-event basis HiGain Sector "); pix.GetHAbsTime()->SetName ("HAbsTimeHiGainPix"); pix.GetHAbsTime()->SetTitle("Absolute Arrival Time average HiGain Sector "); InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadSector(j) : fCam->GetAverageBadSector(j),j); } } if (fAverageLoGainSectors->GetEntries()==0 && IsLoGain()) { fAverageLoGainSectors->Expand(nsectors); for (Int_t j=0; jSetName ("HChargeLoGainPix"); pix.GetHGausHist()->SetTitle("Signals averaged on event-by-event basis LoGain Sector "); pix.GetHAbsTime()->SetName ("HAbsTimeLoGainPix"); pix.GetHAbsTime()->SetTitle("Absolute Arrival Time average LoGain Sector "); // // Adapt the range for the case, the intense blue is used: // FIXME: this is a nasty workaround, but for the moment necessary // in order to avoid default memory space. // if (fGeom->InheritsFrom("MGeomCamMagic")) { if ( fColor == MCalibrationCam::kBLUE) { pix.SetFirst(-10.5); pix.SetLast(999.5); pix.SetNbins(3030); } } InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadSector(j) : fCam->GetAverageBadSector(j),j); } } fSumhiarea .Set(nareas); fSumloarea .Set(nareas); fTimehiarea .Set(nareas); fTimeloarea .Set(nareas); fSumhisector.Set(nsectors); fSumlosector.Set(nsectors); fTimehisector.Set(nsectors); fTimelosector.Set(nsectors); fSathiarea .Set(nareas); fSatloarea .Set(nareas); fSathisector.Set(nsectors); fSatlosector.Set(nsectors); return kTRUE; } // -------------------------------------------------------------------------- // // Retrieves from MExtractedSignalCam: // - first used LoGain FADC slice // // Retrieves from MGeomCam: // - number of pixels // - number of pixel areas // - number of sectors // // For all TObjArray's (including the averaged ones), the following steps are performed: // // 1) Fill Charges histograms (MHGausEvents::FillHistAndArray()) with: // - MExtractedSignalPix::GetExtractedSignalHiGain(); // - MExtractedSignalPix::GetExtractedSignalLoGain(); // // 2) Set number of saturated slices (MHCalibrationChargePix::SetSaturated()) with: // - MExtractedSignalPix::GetNumHiGainSaturated(); // - MExtractedSignalPix::GetNumLoGainSaturated(); // // 3) Fill AbsTime histograms (MHCalibrationChargePix::FillAbsTime()) with: // - MRawEvtPixelIter::GetIdxMaxHiGainSample(); // - MRawEvtPixelIter::GetIdxMaxLoGainSample(first slice); // Bool_t MHCalibrationChargeCam::FillHists(const MParContainer *par, const Stat_t w) { MExtractedSignalCam *signal = (MExtractedSignalCam*)par; if (!signal) { *fLog << err << "No argument in MExtractedSignalCam::Fill... abort." << endl; return kFALSE; } const UInt_t npixels = fGeom->GetNumPixels(); const UInt_t nareas = fGeom->GetNumAreas(); const UInt_t nsectors = fGeom->GetNumSectors(); const UInt_t lofirst = signal->GetFirstUsedSliceLoGain(); fSumhiarea .Reset(); fSumloarea .Reset(); fTimehiarea .Reset(); fTimeloarea .Reset(); fSumhisector.Reset(); fSumlosector.Reset(); fTimehisector.Reset(); fTimelosector.Reset(); fSathiarea .Reset(); fSatloarea .Reset(); fSathisector.Reset(); fSatlosector.Reset(); for (UInt_t i=0; i 0.5 ? 1 : 0); hipix.FillAbsTime (fTimehiarea[j]/npix); if (IsLoGain()) { MHCalibrationChargePix &lopix = (MHCalibrationChargePix&)GetAverageLoGainArea(j); if (IsOscillations()) lopix.FillHistAndArray(fSumloarea [j]/npix); else lopix.FillHist(fSumloarea [j]/npix); lopix.SetSaturated ((Float_t)fSatloarea [j]/npix > 0.5 ? 1 : 0); lopix.FillAbsTime (fTimeloarea[j]/npix); } } for (UInt_t j=0; j 0.5 ? 1 : 0); hipix.FillAbsTime (fTimehisector[j]/npix); if (IsLoGain()) { MHCalibrationChargePix &lopix = (MHCalibrationChargePix&)GetAverageLoGainSector(j); if (IsOscillations()) lopix.FillHistAndArray(fSumlosector [j]/npix); else lopix.FillHist(fSumlosector [j]/npix); lopix.SetSaturated ((Float_t)fSatlosector[j]/npix > 0.5 ? 1 : 0); lopix.FillAbsTime (fTimelosector[j]/npix); } } return kTRUE; } // -------------------------------------------------------------------------- // // For all TObjArray's (including the averaged ones), the following steps are performed: // // 1) Returns if the pixel is excluded. // 2) Tests saturation. In case yes, set the flag: MCalibrationPix::SetHiGainSaturation() // or the flag: MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainSaturated ) // 3) Store the absolute arrival times in the MCalibrationChargePix's. If flag // MCalibrationPix::IsHiGainSaturation() is set, the Low-Gain arrival times are stored, // otherwise the Hi-Gain ones. // 4) Calls to MHCalibrationCam::FitHiGainArrays() and MCalibrationCam::FitLoGainArrays() // with the flags: // - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kHiGainNotFitted ) // - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainNotFitted ) // - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kHiGainOscillating ) // - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainOscillating ) // Bool_t MHCalibrationChargeCam::FinalizeHists() { *fLog << endl; for (Int_t i=0; iGetSize(); i++) { MHCalibrationChargePix &histhi = (MHCalibrationChargePix&)(*this)[i]; if (histhi.IsExcluded()) continue; MCalibrationChargePix &pix = fIntensCam ? (MCalibrationChargePix&)(*fIntensCam)[i] : (MCalibrationChargePix&)(*fCam)[i]; if (histhi.GetSaturated() > fNumHiGainSaturationLimit*histhi.GetHGausHist()->GetEntries()) { pix.SetHiGainSaturation(); if (IsOscillations()) histhi.CreateFourierSpectrum(); continue; } MBadPixelsPix &bad = (*fBadPixels)[i]; Stat_t overflow = histhi.GetHGausHist()->GetBinContent(histhi.GetHGausHist()->GetNbinsX()+1); if (overflow > 0.1) { *fLog << warn << GetDescriptor() << ": HiGain Histogram Overflow occurred " << overflow << " times in pixel: " << i << " (without saturation!) " << endl; bad.SetUncalibrated( MBadPixelsPix::kHiGainOverFlow ); } overflow = histhi.GetHGausHist()->GetBinContent(0); if (overflow > 0.1) { *fLog << warn << GetDescriptor() << ": HiGain Histogram Underflow occurred " << overflow << " times in pixel: " << i << " (without saturation!) " << endl; bad.SetUncalibrated( MBadPixelsPix::kHiGainOverFlow ); } FinalizeAbsTimes(histhi, pix, bad, fFirstHiGain, fLastHiGain); } if (IsLoGain()) for (Int_t i=0; iGetSize(); i++) { MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)(*this)(i); MBadPixelsPix &bad = (*fBadPixels)[i]; if (histlo.IsExcluded()) continue; if (histlo.GetSaturated() > fNumLoGainSaturationLimit*histlo.GetHGausHist()->GetEntries()) { *fLog << warn << "Saturated Lo Gain histogram in pixel: " << i << endl; bad.SetUncalibrated( MBadPixelsPix::kLoGainSaturation ); if (IsOscillations()) histlo.CreateFourierSpectrum(); continue; } Stat_t overflow = histlo.GetHGausHist()->GetBinContent(histlo.GetHGausHist()->GetNbinsX()+1); if (overflow > 0.1) { *fLog << warn << GetDescriptor() << ": Lo-Gain Histogram Overflow occurred " << overflow << " times in pixel: " << i << " (without saturation!) " << endl; bad.SetUncalibrated( MBadPixelsPix::kLoGainOverFlow ); } overflow = histlo.GetHGausHist()->GetBinContent(0); if (overflow > 0.1) { *fLog << warn << GetDescriptor() << ": Lo-Gain Histogram Underflow occurred " << overflow << " times in pixel: " << i << " (without saturation!) " << endl; bad.SetUncalibrated( MBadPixelsPix::kLoGainOverFlow ); } MCalibrationChargePix &pix = fIntensCam ? (MCalibrationChargePix&)(*fIntensCam)[i] : (MCalibrationChargePix&)(*fCam)[i]; if (pix.IsHiGainSaturation()) FinalizeAbsTimes(histlo, pix, bad, fFirstLoGain, fLastLoGain); } for (Int_t j=0; jGetSize(); j++) { MHCalibrationChargePix &histhi = (MHCalibrationChargePix&)GetAverageHiGainArea(j); MCalibrationChargePix &pix = fIntensCam ? (MCalibrationChargePix&)fIntensCam->GetAverageArea(j) : (MCalibrationChargePix&)fCam->GetAverageArea(j); if (histhi.GetSaturated() > fNumHiGainSaturationLimit*histhi.GetHGausHist()->GetEntries()) { pix.SetHiGainSaturation(); if (IsOscillations()) histhi.CreateFourierSpectrum(); continue; } MBadPixelsPix &bad = fIntensCam ? fIntensCam->GetAverageBadArea(j) : fCam->GetAverageBadArea(j); FinalizeAbsTimes(histhi, pix, bad, fFirstHiGain, fLastHiGain); } if (IsLoGain()) for (Int_t j=0; jGetSize(); j++) { MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)GetAverageLoGainArea(j); if (histlo.GetSaturated() > fNumLoGainSaturationLimit*histlo.GetHGausHist()->GetEntries()) { *fLog << warn << "Saturated Lo Gain histogram in area idx: " << j << endl; histlo.CreateFourierSpectrum(); continue; } MCalibrationChargePix &pix = fIntensCam ? (MCalibrationChargePix&)fIntensCam->GetAverageArea(j) : (MCalibrationChargePix&)fCam->GetAverageArea(j) ; if (pix.IsHiGainSaturation()) { MBadPixelsPix &bad = fIntensCam ? fIntensCam->GetAverageBadArea(j) : fCam->GetAverageBadArea(j); FinalizeAbsTimes(histlo, pix, bad, fFirstLoGain, fLastLoGain); } } for (Int_t j=0; jGetSize(); j++) { MHCalibrationChargePix &histhi = (MHCalibrationChargePix&)GetAverageHiGainSector(j); MCalibrationChargePix &pix = fIntensCam ? (MCalibrationChargePix&)fIntensCam->GetAverageSector(j) : (MCalibrationChargePix&)fCam->GetAverageSector(j); if (histhi.GetSaturated() > fNumHiGainSaturationLimit*histhi.GetHGausHist()->GetEntries()) { pix.SetHiGainSaturation(); if (IsOscillations()) histhi.CreateFourierSpectrum(); continue; } MBadPixelsPix &bad = fIntensCam ? fIntensCam->GetAverageBadSector(j) : fCam->GetAverageBadSector(j); FinalizeAbsTimes(histhi, pix, bad, fFirstHiGain, fLastHiGain); } if (IsLoGain()) for (Int_t j=0; jGetSize(); j++) { MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)GetAverageLoGainSector(j); MBadPixelsPix &bad = fIntensCam ? fIntensCam->GetAverageBadSector(j) : fCam->GetAverageBadSector(j); if (histlo.GetSaturated() > fNumLoGainSaturationLimit*histlo.GetHGausHist()->GetEntries()) { *fLog << warn << "Saturated Lo Gain histogram in sector: " << j << endl; bad.SetUncalibrated( MBadPixelsPix::kLoGainSaturation ); if (IsOscillations()) histlo.CreateFourierSpectrum(); continue; } MCalibrationChargePix &pix = fIntensCam ? (MCalibrationChargePix&)fIntensCam->GetAverageSector(j) : (MCalibrationChargePix&)fCam->GetAverageSector(j); if (pix.IsHiGainSaturation()) FinalizeAbsTimes(histlo, pix, bad, fFirstLoGain, fLastLoGain); } // // Perform the fitting for the High Gain (done in MHCalibrationCam) // FitHiGainArrays(fIntensCam ? (MCalibrationCam&)(*fIntensCam->GetCam()) : (MCalibrationCam&)(*fCam), *fBadPixels, MBadPixelsPix::kHiGainNotFitted, MBadPixelsPix::kHiGainOscillating); // // Perform the fitting for the Low Gain (done in MHCalibrationCam) // if (IsLoGain()) FitLoGainArrays(fIntensCam ? (MCalibrationCam&)(*fIntensCam->GetCam()) : (MCalibrationCam&)(*fCam), *fBadPixels, MBadPixelsPix::kLoGainNotFitted, MBadPixelsPix::kLoGainOscillating); return kTRUE; } // -------------------------------------------------------------------------------- // // Fill the absolute time results into MCalibrationChargePix // // Check absolute time validity: // - Mean arrival time is at least fTimeLowerLimit slices from the lower edge // - Mean arrival time is at least fUpperLimit slices from the upper edge // void MHCalibrationChargeCam::FinalizeAbsTimes(MHCalibrationChargePix &hist, MCalibrationChargePix &pix, MBadPixelsPix &bad, Byte_t first, Byte_t last) { const Float_t mean = hist.GetAbsTimeMean(); const Float_t rms = hist.GetAbsTimeRms(); pix.SetAbsTimeMean ( mean ); pix.SetAbsTimeRms ( rms ); const Float_t lowerlimit = (Float_t)first + fTimeLowerLimit; const Float_t upperlimit = (Float_t)last + fTimeUpperLimit; if ( mean < lowerlimit) { *fLog << warn << GetDescriptor() << Form("%s%3.1f%s%2.1f%s%4i",": Mean ArrivalTime: ",mean," smaller than ",fTimeLowerLimit, " FADC slices from lower edge in pixel ",hist.GetPixId()) << endl; bad.SetUncalibrated( MBadPixelsPix::kMeanTimeInFirstBin ); } if ( mean > upperlimit ) { *fLog << warn << GetDescriptor() << Form("%s%3.1f%s%2.1f%s%4i",": Mean ArrivalTime: ",mean," greater than ",fTimeUpperLimit, " FADC slices from upper edge in pixel ",hist.GetPixId()) << endl; bad.SetUncalibrated( MBadPixelsPix::kMeanTimeInLast2Bins ); } } // -------------------------------------------------------------------------- // // Sets all pixels to MBadPixelsPix::kUnsuitableRun, if following flags are set: // - MBadPixelsPix::kLoGainSaturation // // Sets all pixels to MBadPixelsPix::kUnreliableRun, if following flags are set: // - if MBadPixelsPix::kHiGainNotFitted and !MCalibrationPix::IsHiGainSaturation() // - if MBadPixelsPix::kHiGainOscillating and !MCalibrationPix::IsHiGainSaturation() // - if MBadPixelsPix::kLoGainNotFitted and MCalibrationPix::IsLoGainSaturation() // - if MBadPixelsPix::kLoGainOscillating and MCalibrationPix::IsLoGainSaturation() // void MHCalibrationChargeCam::FinalizeBadPixels() { for (Int_t i=0; iGetSize(); i++) { MBadPixelsPix &bad = (*fBadPixels)[i]; MCalibrationPix &pix = fIntensCam ? (*fIntensCam)[i] : (*fCam)[i]; if (bad.IsUncalibrated( MBadPixelsPix::kHiGainNotFitted )) if (!pix.IsHiGainSaturation()) bad.SetUnsuitable( MBadPixelsPix::kUnreliableRun ); if (bad.IsUncalibrated( MBadPixelsPix::kLoGainNotFitted )) if (pix.IsHiGainSaturation()) bad.SetUnsuitable( MBadPixelsPix::kUnreliableRun ); if (bad.IsUncalibrated( MBadPixelsPix::kLoGainSaturation )) bad.SetUnsuitable( MBadPixelsPix::kUnsuitableRun ); if (IsOscillations()) { if (bad.IsUncalibrated( MBadPixelsPix::kHiGainOscillating )) bad.SetUnsuitable( MBadPixelsPix::kUnreliableRun ); if (bad.IsUncalibrated( MBadPixelsPix::kLoGainOscillating )) if (pix.IsHiGainSaturation()) bad.SetUnsuitable( MBadPixelsPix::kUnreliableRun ); } } } // -------------------------------------------------------------------------- // // Dummy, needed by MCamEvent // Bool_t MHCalibrationChargeCam::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const { return kTRUE; } // -------------------------------------------------------------------------- // // Calls MHCalibrationPix::DrawClone() for pixel idx // void MHCalibrationChargeCam::DrawPixelContent(Int_t idx) const { (*this)[idx].DrawClone(); } // ----------------------------------------------------------------------------- // // Default draw: // // Displays the averaged areas, both High Gain and Low Gain // // Calls the Draw of the fAverageHiGainAreas and fAverageLoGainAreas objects with options // void MHCalibrationChargeCam::Draw(const Option_t *opt) { const Int_t nareas = fAverageHiGainAreas->GetEntries(); if (nareas == 0) return; TString option(opt); option.ToLower(); if (!option.Contains("datacheck")) { MHCalibrationCam::Draw(opt); return; } // // From here on , the datacheck - Draw // TVirtualPad *pad = gPad ? gPad : MH::MakeDefCanvas(this); pad->SetBorderMode(0); pad->Divide(1,nareas); // // Loop over inner and outer pixels // for (Int_t i=0; icd(i+1); MHCalibrationChargePix &hipix = (MHCalibrationChargePix&)GetAverageHiGainArea(i); // // Ask for Hi-Gain saturation // if (hipix.GetSaturated() > fNumHiGainSaturationLimit*hipix.GetHGausHist()->GetEntries() && IsLoGain()) { MHCalibrationChargePix &lopix = (MHCalibrationChargePix&)GetAverageLoGainArea(i); DrawDataCheckPixel(lopix,i ? gkLoGainOuterRefLines : gkLoGainInnerRefLines); } else DrawDataCheckPixel(hipix,i ? gkHiGainOuterRefLines : gkHiGainInnerRefLines); } } // -------------------------------------------------------------------------- // // Our own clone function is necessary since root 3.01/06 or Mars 0.4 // I don't know the reason. // // Creates new MHCalibrationCam // TObject *MHCalibrationChargeCam::Clone(const char *name) const { const Int_t navhi = fAverageHiGainAreas->GetEntries(); const Int_t navlo = fAverageLoGainAreas->GetEntries(); const Int_t nsehi = fAverageHiGainSectors->GetEntries(); const Int_t nselo = fAverageLoGainSectors->GetEntries(); // // FIXME, this might be done faster and more elegant, by direct copy. // MHCalibrationChargeCam *cam = new MHCalibrationChargeCam(); cam->fAverageHiGainAreas->Expand(navhi); cam->fAverageLoGainAreas->Expand(navlo); cam->fAverageHiGainSectors->Expand(nsehi); cam->fAverageLoGainSectors->Expand(nselo); cam->fAverageHiGainAreas->Expand(navhi); cam->fAverageLoGainAreas->Expand(navlo); cam->fAverageHiGainSectors->Expand(nsehi); cam->fAverageLoGainSectors->Expand(nselo); for (int i=0; ifAverageHiGainAreas) [i] = (*fAverageHiGainAreas) [i]->Clone(); for (int i=0; ifAverageLoGainAreas) [i] = (*fAverageLoGainAreas) [i]->Clone(); for (int i=0; ifAverageHiGainSectors)[i] = (*fAverageHiGainSectors)[i]->Clone(); for (int i=0; ifAverageLoGainSectors)[i] = (*fAverageLoGainSectors)[i]->Clone(); cam->fAverageAreaNum = fAverageAreaNum; cam->fAverageAreaSat = fAverageAreaSat; cam->fAverageAreaSigma = fAverageAreaSigma; cam->fAverageAreaSigmaVar = fAverageAreaSigmaVar; cam->fAverageAreaRelSigma = fAverageAreaRelSigma; cam->fAverageAreaRelSigmaVar = fAverageAreaRelSigmaVar; cam->fAverageSectorNum = fAverageSectorNum; cam->fRunNumbers = fRunNumbers; cam->fColor = fColor; cam->fPulserFrequency = fPulserFrequency; return cam; } void MHCalibrationChargeCam::DrawDataCheckPixel(MHCalibrationChargePix &pix, const Float_t refline[]) { TVirtualPad *newpad = gPad; newpad->Divide(1,2); newpad->cd(1); gPad->SetTicks(); if (!pix.IsEmpty() && !pix.IsOnlyOverflow() && !pix.IsOnlyUnderflow()) gPad->SetLogy(); TH1F *hist = pix.GetHGausHist(); TH2D *null = new TH2D("Null",hist->GetTitle(),100,pix.GetFirst() > -1. ? 0. : 100.,pix.GetLast()/2., 100,0.,hist->GetEntries()/10.); null->SetDirectory(NULL); null->SetBit(kCanDelete); null->SetStats(kFALSE); // // set the labels bigger // TAxis *xaxe = null->GetXaxis(); TAxis *yaxe = null->GetYaxis(); xaxe->CenterTitle(); yaxe->CenterTitle(); xaxe->SetTitleSize(0.07); yaxe->SetTitleSize(0.07); xaxe->SetTitleOffset(0.7); yaxe->SetTitleOffset(0.55); xaxe->SetLabelSize(0.06); yaxe->SetLabelSize(0.06); xaxe->SetTitle(hist->GetXaxis()->GetTitle()); yaxe->SetTitle(hist->GetYaxis()->GetTitle()); null->Draw(); hist->Draw("same"); gStyle->SetOptFit(); if (pix.GetFGausFit()) { switch ( fColor ) { case MCalibrationCam::kGREEN: pix.GetFGausFit()->SetLineColor(kGreen); break; case MCalibrationCam::kBLUE: pix.GetFGausFit()->SetLineColor(kBlue); break; case MCalibrationCam::kUV: pix.GetFGausFit()->SetLineColor(106); break; case MCalibrationCam::kCT1: pix.GetFGausFit()->SetLineColor(006); break; default: pix.GetFGausFit()->SetLineColor(kRed); } pix.GetFGausFit()->Draw("same"); } DisplayRefLines(null,refline); newpad->cd(2); gPad->SetTicks(); pix.DrawEvents(); return; } void MHCalibrationChargeCam::DisplayRefLines(const TH2D *hist, const Float_t refline[]) const { TLine *green1 = new TLine(refline[0],0.,refline[0],hist->GetYaxis()->GetXmax()); green1->SetBit(kCanDelete); green1->SetLineColor(kGreen); green1->SetLineStyle(2); green1->SetLineWidth(3); green1->Draw(); TLine *green5 = new TLine(refline[6],0.,refline[6],hist->GetYaxis()->GetXmax()); green5->SetBit(kCanDelete); green5->SetLineColor(8); green5->SetLineStyle(2); green5->SetLineWidth(3); green5->Draw(); TLine *blue1 = new TLine(refline[1],0.,refline[1],hist->GetYaxis()->GetXmax()); blue1->SetBit(kCanDelete); blue1->SetLineColor(227); blue1->SetLineStyle(2); blue1->SetLineWidth(3); blue1->Draw(); TLine *blue5 = new TLine(refline[2],0.,refline[2],hist->GetYaxis()->GetXmax()); blue5->SetBit(kCanDelete); blue5->SetLineColor(68); blue5->SetLineStyle(2); blue5->SetLineWidth(3); blue5->Draw(); TLine *blue10 = new TLine(refline[3],0.,refline[3],hist->GetYaxis()->GetXmax()); blue10->SetBit(kCanDelete); blue10->SetLineColor(4); blue10->SetLineStyle(2); blue10->SetLineWidth(3); blue10->Draw(); TLine *uv10 = new TLine(refline[4],0.,refline[4],hist->GetYaxis()->GetXmax()); uv10->SetBit(kCanDelete); uv10->SetLineColor(106); uv10->SetLineStyle(2); uv10->SetLineWidth(3); uv10->Draw(); TLine *ct1 = new TLine(refline[5],0.,refline[5],hist->GetYaxis()->GetXmax()); ct1->SetBit(kCanDelete); ct1->SetLineColor(6); ct1->SetLineStyle(2); ct1->SetLineWidth(3); ct1->Draw(); TLegend *leg = new TLegend(0.7,0.35,0.9,0.99); leg->SetBit(kCanDelete); leg->AddEntry(green1,"1 Led GREEN","l"); leg->AddEntry(green5,"5 Leds GREEN","l"); leg->AddEntry(blue1,"1 Led BLUE","l"); leg->AddEntry(blue5,"5 Leds BLUE","l"); leg->AddEntry(blue10,"10 Leds BLUE","l"); leg->AddEntry(uv10,"10 Leds UV","l"); leg->AddEntry(ct1,"CT1-Pulser","l"); leg->Draw(); }