/* ======================================================================== *\ ! ! * ! * 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 ! ! \* ======================================================================== */ ///////////////////////////////////////////////////////////////////////////// // // MHCalibrationCam // // Base class for camera calibration classes. Incorporates the TOrdCollection's: // - fHiGainArray (for calibrated High Gains per pixel) // - fLoGainArray (for calibrated Low Gains per pixel) // - fAverageHiGainAreas (for averaged High Gains events per camera area index) // - fAverageLoGainAreas (for averaged High Gains events per camera area index) // - fAverageHiGainSectors (for averaged High Gains events per camera sector ) // - fAverageLoGainSectors (for averaged High Gains events per camera sector ) // These TOrdCollection's are called by their default constructors, thus no objects // are created, until the derived class does so. // // The corresponding operators: [],() and the operators GetAverageHiGainArea(), // GetAverageLoGainArea(), GetAverageHiGainSector() and GetAverageLoGainSector() // have to be cast to the corresponding class. It is assumed that all classes // dealing with calibration pixels derive from MHCalibrationPix. // // The following flag can be set: // - SetAverageing() for calculating the event-by-event averages. // - SetDebug() for debug output. // - SetLoGain() for the case that low-gain slices are available, but // MRawRunHeader::GetNumLoGainSlices() gives still 0. // - SetCheckSize() for testing the sizes of the f*Arrays in ReInit() if they // are coinciding with the equiv. sizes in MGeomCam // // The flag kLoGain steers if the low-gain signal is treated at all or not. // The flag kAverageing steers if the event-by-event averages are treated at all. // ///////////////////////////////////////////////////////////////////////////// #include "MHCalibrationCam.h" #include "MHCalibrationPix.h" #include #include #include #include #include #include #include #include "MLog.h" #include "MLogManip.h" #include "MCalibrationIntensityCam.h" #include "MCalibrationCam.h" #include "MCalibrationPix.h" #include "MBadPixelsIntensityCam.h" #include "MBadPixelsCam.h" #include "MBadPixelsPix.h" #include "MGeomCam.h" #include "MGeomPix.h" #include "MParList.h" #include "MRawRunHeader.h" ClassImp(MHCalibrationCam); using namespace std; const Int_t MHCalibrationCam::fgPulserFrequency = 500; const Float_t MHCalibrationCam::fgProbLimit = 0.0001; const Float_t MHCalibrationCam::fgOverflowLimit = 0.005; const TString MHCalibrationCam::gsHistName = "Hist"; const TString MHCalibrationCam::gsHistTitle = ""; const TString MHCalibrationCam::gsHistXTitle = ""; const TString MHCalibrationCam::gsHistYTitle = "Nr. events"; // -------------------------------------------------------------------------- // // Default Constructor. // // Sets: // - all pointers to NULL // // Initializes and sets owner of: // - fHiGainArray, fLoGainArray // - fAverageHiGainAreas, fAverageLoGainAreas // - fAverageHiGainSectors, fAverageLoGainSectors // // Initializes: // - fPulserFrequency to fgPulserFrequency // - fProbLimit to fgProbLimit // - fOverflowLimit to fgOverflowLimit // // - SetAveregeing (kTRUE); // - SetDebug (kFALSE); // - SetLoGain (kTRUE); //- SetOscillations(kTRUE); //- SetSizeCheck (kTRUE); //- SetInterlaced (kFALSE); // MHCalibrationCam::MHCalibrationCam(const char *name, const char *title) : fHistName(gsHistName),fHistTitle(gsHistTitle), fHistXTitle(gsHistXTitle),fHistYTitle(gsHistYTitle), fColor(MCalibrationCam::kNONE), fIntensBad(NULL), fBadPixels(NULL), fIntensCam(NULL), fCam(NULL), fGeom(NULL), fRunHeader(NULL) { fHiGainArray = new TOrdCollection; fHiGainArray->SetOwner(); fLoGainArray = new TOrdCollection; fLoGainArray->SetOwner(); fAverageHiGainAreas = new TOrdCollection; fAverageHiGainAreas->SetOwner(); fAverageLoGainAreas = new TOrdCollection; fAverageLoGainAreas->SetOwner(); fAverageHiGainSectors = new TOrdCollection; fAverageHiGainSectors->SetOwner(); fAverageLoGainSectors = new TOrdCollection; fAverageLoGainSectors->SetOwner(); SetPulserFrequency(); SetProbLimit(); SetOverflowLimit(); SetAverageing (kTRUE); SetDebug (kFALSE); SetLoGain (kTRUE); SetOscillations(kTRUE); SetSizeCheck (kTRUE); } // -------------------------------------------------------------------------- // // Deletes the TOrdCollection of: // - fHiGainArray, fLoGainArray // - fAverageHiGainAreas, fAverageLoGainAreas // - fAverageHiGainSectors, fAverageLoGainSectors // MHCalibrationCam::~MHCalibrationCam() { delete fHiGainArray; delete fLoGainArray; delete fAverageHiGainAreas; delete fAverageLoGainAreas; delete fAverageHiGainSectors; delete fAverageLoGainSectors; /* Remove ( fHiGainArray ); Remove ( fLoGainArray ); Remove ( fAverageHiGainAreas ); Remove ( fAverageLoGainAreas ); Remove ( fAverageHiGainSectors ); Remove ( fAverageLoGainSectors ); */ } void MHCalibrationCam::Remove(TOrdCollection *col) { if (!col) return; TOrdCollectionIter Next(col); Int_t count = 0; while (MHCalibrationPix *obj = (MHCalibrationPix*)Next()) { *fLog << ++count << " " << obj << flush; if (obj && obj->IsOnHeap()) { obj->Draw(); delete obj; } } delete col; } // -------------------------------------------------------------------------- // // Returns size of fHiGainArray // const Int_t MHCalibrationCam::GetSize() const { return fHiGainArray->GetSize(); } // -------------------------------------------------------------------------- // // Get i-th High Gain pixel (pixel number) // MHCalibrationPix &MHCalibrationCam::operator[](UInt_t i) { return *static_cast(fHiGainArray->At(i)); } // -------------------------------------------------------------------------- // // Get i-th High Gain pixel (pixel number) // const MHCalibrationPix &MHCalibrationCam::operator[](UInt_t i) const { return *static_cast(fHiGainArray->At(i)); } // -------------------------------------------------------------------------- // // Get i-th Low Gain pixel (pixel number) // MHCalibrationPix &MHCalibrationCam::operator()(UInt_t i) { return *static_cast(fLoGainArray->At(i)); } // -------------------------------------------------------------------------- // // Get i-th Low Gain pixel (pixel number) // const MHCalibrationPix &MHCalibrationCam::operator()(UInt_t i) const { return *static_cast(fLoGainArray->At(i)); } // -------------------------------------------------------------------------- // // Returns the current size of the TOrdCollection fAverageHiGainAreas // independently if the MHCalibrationPix is filled with values or not. // const Int_t MHCalibrationCam::GetAverageAreas() const { return fAverageHiGainAreas->GetSize(); } // -------------------------------------------------------------------------- // // Get i-th High Gain pixel Area (area number) // MHCalibrationPix &MHCalibrationCam::GetAverageHiGainArea(UInt_t i) { return *static_cast(fAverageHiGainAreas->At(i)); } // -------------------------------------------------------------------------- // // Get i-th High Gain pixel Area (area number) // const MHCalibrationPix &MHCalibrationCam::GetAverageHiGainArea(UInt_t i) const { return *static_cast(fAverageHiGainAreas->At(i)); } // -------------------------------------------------------------------------- // // Get i-th Low Gain pixel Area (area number) // MHCalibrationPix &MHCalibrationCam::GetAverageLoGainArea(UInt_t i) { return *static_cast(fAverageLoGainAreas->At(i)); } // -------------------------------------------------------------------------- // // Get i-th Low Gain pixel Area (area number) // const MHCalibrationPix &MHCalibrationCam::GetAverageLoGainArea(UInt_t i) const { return *static_cast(fAverageLoGainAreas->At(i)); } // -------------------------------------------------------------------------- // // Returns the current size of the TOrdCollection fAverageHiGainSectors // independently if the MHCalibrationPix is filled with values or not. // const Int_t MHCalibrationCam::GetAverageSectors() const { return fAverageHiGainSectors->GetSize(); } // -------------------------------------------------------------------------- // // Get i-th High Gain Sector (sector number) // MHCalibrationPix &MHCalibrationCam::GetAverageHiGainSector(UInt_t i) { return *static_cast(fAverageHiGainSectors->At(i)); } // -------------------------------------------------------------------------- // // Get i-th High Gain Sector (sector number) // const MHCalibrationPix &MHCalibrationCam::GetAverageHiGainSector(UInt_t i) const { return *static_cast(fAverageHiGainSectors->At(i)); } // -------------------------------------------------------------------------- // // Get i-th Low Gain Sector (sector number) // MHCalibrationPix &MHCalibrationCam::GetAverageLoGainSector(UInt_t i) { return *static_cast(fAverageLoGainSectors->At(i)); } // -------------------------------------------------------------------------- // // Get i-th Low Gain Sector (sector number) // const MHCalibrationPix &MHCalibrationCam::GetAverageLoGainSector(UInt_t i) const { return *static_cast(fAverageLoGainSectors->At(i)); } // -------------------------------------------------------------------------- // // Calls ResetHistTitles() // // Calls Reset() for each entry in: // - fHiGainArray, fLoGainArray // - fAverageHiGainAreas, fAverageLoGainAreas // - fAverageHiGainSectors, fAverageLoGainSectors // void MHCalibrationCam::ResetHists() { ResetHistTitles(); if (fHiGainArray) { fHiGainArray->ForEach(MHCalibrationPix,Reset)(); } if (IsAverageing()) { if (fAverageHiGainAreas) { fAverageHiGainAreas->ForEach(MHCalibrationPix,Reset)(); } if (fAverageHiGainSectors) { fAverageHiGainSectors->ForEach(MHCalibrationPix,Reset)(); } } if (!IsLoGain()) return; if (fLoGainArray) { fLoGainArray->ForEach(MHCalibrationPix,Reset)(); } if (IsAverageing()) { if (fAverageLoGainAreas) { fAverageLoGainAreas->ForEach(MHCalibrationPix,Reset)(); } if (fAverageLoGainSectors) { fAverageLoGainSectors->ForEach(MHCalibrationPix,Reset)(); } } } // -------------------------------------------------------------------------- // // Resets the histogram titles for each entry in: // - fHiGainArray, fLoGainArray // - fAverageHiGainAreas, fAverageLoGainAreas // - fAverageHiGainSectors, fAverageLoGainSectors // void MHCalibrationCam::ResetHistTitles() { TH1F *h; if (fHiGainArray) for (Int_t i=0;iGetSize(); i++) { MHCalibrationPix &pix = (*this)[i]; h = pix.GetHGausHist(); h->SetName (Form("%s%s%s%4i","H",fHistName.Data(),"HiGainPix",i)); h->SetTitle(Form("%s%s%4i%s",fHistTitle.Data()," High Gain Pixel ",i," Runs: ")); h->SetXTitle(fHistXTitle.Data()); h->SetYTitle(fHistYTitle.Data()); } if (IsAverageing()) { if (fAverageHiGainAreas) for (Int_t j=0; jGetSize(); j++) { MHCalibrationPix &pix = GetAverageHiGainArea(j); h = pix.GetHGausHist(); h->SetName (Form("%s%s%s%d","H",fHistName.Data(),"HiGainArea",j)); h->SetXTitle(fHistXTitle.Data()); h->SetYTitle(fHistYTitle.Data()); if (fGeom->InheritsFrom("MGeomCamMagic")) h->SetTitle(Form("%s%s%s%s",fHistTitle.Data()," averaged on event-by-event basis ", j==0 ? "Inner Pixels " : "Outer Pixels ","High Gain Runs: ")); else h->SetTitle(Form("%s%s%d%s",fHistTitle.Data(), " averaged on event-by-event basis High Gain Area Idx ",j," Runs: ")); } if (fAverageHiGainSectors) for (Int_t j=0; jGetSize(); j++) { MHCalibrationPix &pix = GetAverageHiGainSector(j); h = pix.GetHGausHist(); h->SetName (Form("%s%s%s%2i","H",fHistName.Data(),"HiGainSector",j)); h->SetTitle(Form("%s%s%2i%s",fHistTitle.Data(), " averaged on event-by-event basis High Gain Sector ",j," Runs: ")); h->SetXTitle(fHistXTitle.Data()); h->SetYTitle(fHistYTitle.Data()); } } if (!IsLoGain()) return; if (fLoGainArray) for (Int_t i=0;iGetSize(); i++) { MHCalibrationPix &pix = (*this)(i); h = pix.GetHGausHist(); h->SetName (Form("%s%s%s%4i","H",fHistName.Data(),"LoGainPix",i)); h->SetTitle(Form("%s%s%4i%s",fHistTitle.Data()," Low Gain Pixel ",i," Runs: ")); h->SetXTitle(fHistXTitle.Data()); h->SetYTitle(fHistYTitle.Data()); } if (IsAverageing()) { if (fAverageLoGainAreas) for (Int_t j=0; jGetSize(); j++) { MHCalibrationPix &pix = GetAverageLoGainArea(j); h = pix.GetHGausHist(); h->SetName (Form("%s%s%s%d","H",fHistName.Data(),"LoGainArea",j)); h->SetXTitle(fHistXTitle.Data()); h->SetYTitle(fHistYTitle.Data()); if (fGeom->InheritsFrom("MGeomCamMagic")) h->SetTitle(Form("%s%s%s%s",fHistTitle.Data()," averaged on event-by-event basis ", j==0 ? "Inner Pixels " : "Outer Pixels ","Low Gain Runs: ")); else h->SetTitle(Form("%s%s%d%s",fHistTitle.Data(), " averaged on event-by-event basis Low Gain Area Idx ",j," Runs: ")); } if (fAverageLoGainSectors) for (Int_t j=0; jGetSize(); j++) { MHCalibrationPix &pix = GetAverageLoGainSector(j); h = pix.GetHGausHist(); h->SetName (Form("%s%s%s%2i","H",fHistName.Data(),"LoGainSector",j)); h->SetTitle(Form("%s%s%2i%s",fHistTitle.Data(), " averaged on event-by-event basis Low Gain Sector ",j," Runs: ")); h->SetXTitle(fHistXTitle.Data()); h->SetYTitle(fHistYTitle.Data()); } } } // -------------------------------------------------------------------------- // // Gets the pointers to: // - MGeomCam // // Calls SetupHists(const MParList *pList) // // Calls Delete-Function of: // - MHCalibrationCam::fHiGainArray, MHCalibrationCam::fLoGainArray // - MHCalibrationCam::fAverageHiGainAreas, MHCalibrationCam::fAverageLoGainAreas // - MHCalibrationCam::fAverageHiGainSectors, MHCalibrationCam::fAverageLoGainSectors // Bool_t MHCalibrationCam::SetupFill(const MParList *pList) { fGeom = (MGeomCam*)pList->FindObject("MGeomCam"); if (!fGeom) { *fLog << err << GetDescriptor() << ": MGeomCam not found... aborting." << endl; return kFALSE; } fRunHeader = (MRawRunHeader*)pList->FindObject("MRawRunHeader"); if (!fRunHeader) { *fLog << warn << GetDescriptor() << ": MRawRunHeader not found... will not store run numbers." << endl; } return SetupHists(pList); } // -------------------------------------------------------------------------- // // Searches MRawEvtHeader to find the correct pulser colour // // Gets or creates the pointers to: // - MBadPixelsIntensityCam // - MBadPixelsCam // // Searches pointer to: // - MArrivalTimeCam // // Initializes, if empty to MArrivalTimeCam::GetSize() for: // - MHCalibrationCam::fHiGainArray, MHCalibrationCam::fLoGainArray // // Initializes, if empty to MGeomCam::GetNumAreas() for: // - MHCalibrationCam::fAverageHiGainAreas, MHCalibrationCam::fAverageLoGainAreas // // Initializes, if empty to MGeomCam::GetNumSectors() for: // - MHCalibrationCam::fAverageHiGainSectors, MHCalibrationCam::fAverageLoGainSectors // // Initializes TArray's to MGeomCam::GetNumAreas and MGeomCam::GetNumSectors, respectively // Fills with number of valid pixels (if !MBadPixelsPix::IsBad()): // - MHCalibrationCam::fAverageAreaNum[area index] // - MHCalibrationCam::fAverageSectorNum[area index] // // Calls InitializeHists() for every entry in: // - MHCalibrationCam::fHiGainArray // - MHCalibrationCam::fAverageHiGainAreas // - MHCalibrationCam::fAverageHiGainSectors // // Sets Titles and Names for the Histograms // - MHCalibrationCam::fAverageHiGainAreas // - MHCalibrationCam::fAverageHiGainSectors // // Retrieves the run numbers from MRawRunHeader and stores them in fRunNumbers // Bool_t MHCalibrationCam::ReInit(MParList *pList) { const Int_t npixels = fGeom->GetNumPixels(); const Int_t nsectors = fGeom->GetNumSectors(); const Int_t nareas = fGeom->GetNumAreas(); fIntensBad = (MBadPixelsIntensityCam*)pList->FindObject(AddSerialNumber("MBadPixelsIntensityCam")); if (fIntensBad) *fLog << inf << "Found MBadPixelsIntensityCam ... " << endl; else { fBadPixels = (MBadPixelsCam*)pList->FindObject("MBadPixelsCam"); if (!fBadPixels) { fBadPixels = (MBadPixelsCam*)pList->FindCreateObj(AddSerialNumber("MBadPixelsCam")); if (!fBadPixels) { *fLog << err << "Cannot find nor create MBadPixelsCam ... abort." << endl; return kFALSE; } else fBadPixels->InitSize(npixels); } } if (IsAverageing()) { // // The function TArrayF::Set() already sets all entries to 0. // fAverageAreaNum. Set(nareas); fAverageAreaSat. Set(nareas); fAverageAreaSigma. Set(nareas); fAverageAreaSigmaVar. Set(nareas); fAverageAreaRelSigma. Set(nareas); fAverageAreaRelSigmaVar.Set(nareas); fAverageSectorNum. Set(nsectors); for (Int_t aidx=0; aidxGetRunNumber(); if (IsLoGain()) SetLoGain(fRunHeader->GetNumSamplesLoGain()); } if (!ReInitHists(pList)) return kFALSE; ResetHistTitles(); if (!fRunHeader) return kTRUE; for (Int_t i=0; iGetSize(); i++) { TH1F *h = (*this)[i].GetHGausHist(); h->SetTitle( Form("%s%i%s", h->GetTitle(),fRunNumbers[fRunNumbers.GetSize()-1]," ")); } if (IsLoGain()) for (Int_t i=0; iGetSize(); i++) { TH1F *h = (*this)(i).GetHGausHist(); h->SetTitle( Form("%s%i%s", h->GetTitle(),fRunNumbers[fRunNumbers.GetSize()-1]," ")); } if (!IsAverageing()) return kTRUE; for (Int_t j=0; jSetTitle( Form("%s%i%s", h->GetTitle(),fRunNumbers[fRunNumbers.GetSize()-1]," ")); } if (IsLoGain()) for (Int_t j=0; jSetTitle( Form("%s%i%s", h->GetTitle(),fRunNumbers[fRunNumbers.GetSize()-1]," ")); } for (Int_t j=0; jSetTitle( Form("%s%i%s", h->GetTitle(),fRunNumbers[fRunNumbers.GetSize()-1]," ")); } if (IsLoGain()) for (Int_t j=0; jSetTitle( Form("%s%i%s", h->GetTitle(),fRunNumbers[fRunNumbers.GetSize()-1]," ")); } return kTRUE; } //-------------------------------------------------------------------------------------- // // Initializes the High Gain Arrays: // // - For every entry in the expanded arrays: // * Initialize an MHCalibrationPix // * Set Binning from fNbins, fFirst and fLast // * Set Histgram names and titles from fHistName and fHistTitle // * Set X-axis and Y-axis titles with fHistXTitle and fHistYTitle // * Call InitHists // void MHCalibrationCam::InitHiGainArrays(const Int_t npixels, const Int_t nareas, const Int_t nsectors) { if (fHiGainArray->GetSize()==0) { for (Int_t i=0; iAddAt(new MHCalibrationPix(Form("%sHiGainPix%04d",fHistName.Data(),i), Form("%s High Gain Pixel %4d",fHistTitle.Data(),i)),i); MHCalibrationPix &pix = (*this)[i]; pix.SetNbins(fNbins); pix.SetFirst(fFirst); pix.SetLast (fLast); pix.SetProbLimit(fProbLimit); MBadPixelsPix &bad = fIntensBad ? (*fIntensBad)[i] : (*fBadPixels)[i]; InitHists(pix,bad,i); } } if (!IsAverageing()) return; if (fAverageHiGainAreas->GetSize()==0) { for (Int_t j=0; jAddAt(new MHCalibrationPix(Form("%sHiGainArea%d",fHistName.Data(),j), Form("%s High Gain Area Idx %d",fHistTitle.Data(),j)),j); MHCalibrationPix &pix = GetAverageHiGainArea(j); pix.SetNbins(fNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nareas)); pix.SetFirst(fFirst); pix.SetLast (fLast); if (fGeom && fGeom->InheritsFrom("MGeomCamMagic")) { pix.InitBins(); pix.SetEventFrequency(fPulserFrequency); } else InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadArea(j) : fCam->GetAverageBadArea(j),j); } } if (fAverageHiGainSectors->GetSize()==0) { for (Int_t j=0; jAddAt(new MHCalibrationPix(Form("%sHiGainSector%02d",fHistName.Data(),j), Form("%s High Gain Sector %02d",fHistTitle.Data(),j)),j); MHCalibrationPix &pix = GetAverageHiGainSector(j); pix.SetNbins(fNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nsectors)); pix.SetFirst(fFirst); pix.SetLast (fLast); InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadSector(j) : fCam->GetAverageBadSector(j),j); } } } //-------------------------------------------------------------------------------------- // // Return, if IsLoGain() is kFALSE // // Initializes the Low Gain Arrays: // // - For every entry in the expanded arrays: // * Initialize an MHCalibrationPix // * Set Binning from fNbins, fFirst and fLast // * Set Histgram names and titles from fHistName and fHistTitle // * Set X-axis and Y-axis titles with fHistXTitle and fHistYTitle // * Call InitHists // void MHCalibrationCam::InitLoGainArrays(const Int_t npixels, const Int_t nareas, const Int_t nsectors) { if (!IsLoGain()) return; if (fLoGainArray->GetSize()==0) { for (Int_t i=0; iAddAt(new MHCalibrationPix(Form("%sLoGainPix%04d",fHistName.Data(),i), Form("%s Low Gain Pixel%04d",fHistTitle.Data(),i)),i); MHCalibrationPix &pix = (*this)(i); pix.SetNbins(fNbins); pix.SetFirst(fFirst); pix.SetLast (fLast); pix.SetProbLimit(fProbLimit); MBadPixelsPix &bad = fIntensBad ? (*fIntensBad)[i] : (*fBadPixels)[i]; InitHists(pix,bad,i); } } if (!IsAverageing()) return; if (fAverageLoGainAreas->GetSize()==0) { for (Int_t j=0; jAddAt(new MHCalibrationPix(Form("%sLoGainArea%d",fHistName.Data(),j), Form("%s Low Gain Area Idx %d",fHistTitle.Data(),j)),j); MHCalibrationPix &pix = GetAverageLoGainArea(j); pix.SetNbins(fNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nareas)); pix.SetFirst(fFirst); pix.SetLast (fLast); if (fGeom && fGeom->InheritsFrom("MGeomCamMagic")) { pix.InitBins(); pix.SetEventFrequency(fPulserFrequency); } else InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadArea(j) : fCam->GetAverageBadArea(j),j); } } if (fAverageLoGainSectors->GetSize()==0) { for (Int_t j=0; jAddAt(new MHCalibrationPix(Form("%sLoGainSector%02d",fHistName.Data(),j), Form("%s Low Gain Sector %02d",fHistTitle.Data(),j)),j); MHCalibrationPix &pix = GetAverageLoGainSector(j); pix.SetNbins(fNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nsectors)); pix.SetFirst(fFirst); pix.SetLast (fLast); InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadSector(j) : fCam->GetAverageBadSector(j),j); } } } //-------------------------------------------------------------------------------- // // Retrieves from MGeomCam: // - number of pixels // - number of pixel areas // - number of sectors // // Return kFALSE, if sizes of the TOrdCollections do not match npixels, nareas or nsectors // // Call FillHists() // Bool_t MHCalibrationCam::Fill(const MParContainer *par, const Stat_t w) { if (!IsSizeCheck()) return FillHists(par,w); const Int_t npixels = fGeom->GetNumPixels(); const Int_t nareas = fGeom->GetNumAreas(); const Int_t nsectors = fGeom->GetNumSectors(); // // Hi-Gain OrdCollections // if (fHiGainArray->GetSize() != npixels) { *fLog << err << "ERROR - Size mismatch in number of pixels... abort." << endl; return kFALSE; } if (IsLoGain()) { if (fLoGainArray->GetSize() != npixels) { *fLog << err << "ERROR - Size mismatch in number of pixels... abort." << endl; return kFALSE; } } if (!IsAverageing()) return FillHists(par,w); if (fAverageHiGainAreas->GetSize() != nareas) { *fLog << err << "ERROR - Size mismatch in number of areas ... abort." << endl; return kFALSE; } if (fAverageHiGainSectors->GetSize() != nsectors) { *fLog << err << "ERROR - Size mismatch in number of sectors ... abort." << endl; return kFALSE; } if (IsLoGain()) { if (fAverageLoGainAreas->GetSize() != nareas) { *fLog << err << "ERROR - Size mismatch in number of areas ... abort." << endl; return kFALSE; } if (fAverageLoGainSectors->GetSize() != nsectors) { *fLog << err << "ERROR - Size mismatch in number of sectors ... abort." << endl; return kFALSE; } } return FillHists(par, w); } // -------------------------------------------------------------------------- // // 0) Ask if fHiGainArray and fLoGainArray have been initialized, // otherwise return kFALSE. // 1) FinalizeHists() // 2) FinalizeBadPixels() // 3) CalcAverageSigma() // Bool_t MHCalibrationCam::Finalize() { if (fHiGainArray->GetSize() == 0 && fLoGainArray->GetSize() == 0) { *fLog << err << GetDescriptor() << ": ERROR: Both (HiGain and LoGain) histogram arrays have not been initialized... abort." << endl; return kFALSE; } for (Int_t i=0; iGetSize(); i++) { TH1F *h = GetAverageHiGainArea(i).GetHGausHist(); switch (fColor) { case MCalibrationCam::kNONE: break; case MCalibrationCam::kBLUE: h->SetTitle( Form("%s%s", h->GetTitle(),"BLUE ")); break; case MCalibrationCam::kGREEN: h->SetTitle( Form("%s%s", h->GetTitle(),"GREEN ")); break; case MCalibrationCam::kUV: h->SetTitle( Form("%s%s", h->GetTitle(),"UV ")); break; case MCalibrationCam::kCT1: h->SetTitle( Form("%s%s", h->GetTitle(),"CT1-Pulser ")); break; } } for (Int_t i=0; iGetSize(); i++) { TH1F *h = GetAverageLoGainArea(i).GetHGausHist(); switch (fColor) { case MCalibrationCam::kNONE: break; case MCalibrationCam::kBLUE: h->SetTitle( Form("%s%s", h->GetTitle(),"BLUE ")); break; case MCalibrationCam::kGREEN: h->SetTitle( Form("%s%s", h->GetTitle(),"GREEN ")); break; case MCalibrationCam::kUV: h->SetTitle( Form("%s%s", h->GetTitle(),"UV ")); break; case MCalibrationCam::kCT1: h->SetTitle( Form("%s%s", h->GetTitle(),"CT1-Pulser ")); break; } } if (!FinalizeHists()) return kFALSE; FinalizeBadPixels(); CalcAverageSigma(); return kTRUE; } // ------------------------------------------------------------- // // If MBadPixelsPix::IsUnsuitable(MBadPixelsPix::kUnsuitableRun): // - calls MHCalibrationPix::SetExcluded() // // Calls: // - MHGausEvents::InitBins() // - MHCalibrationPix::ChangeHistId(i) // - MHCalibrationPix::SetEventFrequency(fPulserFrequency) // void MHCalibrationCam::InitHists(MHCalibrationPix &hist, MBadPixelsPix &bad, const Int_t i) { if (bad.IsUnsuitable(MBadPixelsPix::kUnsuitableRun)) hist.SetExcluded(); hist.InitBins(); hist.SetEventFrequency(fPulserFrequency); } // ------------------------------------------------------------- // // - Searches for the CalibrationIntensity*Cam corresponding to 'name'. // - In case, it does not exist in the parameter list, it searches // for the corresponding MCalibration*Cam. // - Initializes the MCalibration*Cam, if not yet done. // Bool_t MHCalibrationCam::InitCams( MParList *plist, const TString name ) { TString intensname = "MCalibrationIntensity"; intensname += name; intensname += "Cam"; TString ordname = "MCalibration"; ordname += name; ordname += "Cam"; fIntensCam = (MCalibrationIntensityCam*)plist->FindObject(AddSerialNumber(intensname)); if (fIntensCam) *fLog << inf << "Found " << intensname << "... " << endl; else { fCam = (MCalibrationCam*)plist->FindObject(AddSerialNumber(ordname)); if (!fCam) { fCam = (MCalibrationCam*)plist->FindCreateObj(AddSerialNumber(ordname)); if (!fCam) { *fLog << err << "Cannot find nor create " << ordname << " ... abort." << endl; return kFALSE; } fCam->Init(*fGeom); } } return kTRUE; } // -------------------------------------------------------------------------- // // Calls FitHiGainHists for every entry in: // - fHiGainArray // - fAverageHiGainAreas // - fAverageHiGainSectors // void MHCalibrationCam::FitHiGainArrays(MCalibrationCam &calcam, MBadPixelsCam &badcam, MBadPixelsPix::UncalibratedType_t fittyp, MBadPixelsPix::UncalibratedType_t osctyp) { for (Int_t i=0; iGetSize(); i++) { MHCalibrationPix &hist = (*this)[i]; if (hist.IsExcluded()) continue; MCalibrationPix &pix = calcam[i]; MBadPixelsPix &bad = badcam[i]; FitHiGainHists(hist,pix,bad,fittyp,osctyp); } if (!IsAverageing()) return; for (Int_t j=0; jGetSize(); j++) { MHCalibrationPix &hist = GetAverageHiGainArea(j); MCalibrationPix &pix = calcam.GetAverageArea(j); MBadPixelsPix &bad = calcam.GetAverageBadArea(j); FitHiGainHists(hist,pix,bad,fittyp,osctyp); } for (Int_t j=0; jGetSize(); j++) { MHCalibrationPix &hist = GetAverageHiGainSector(j); MCalibrationPix &pix = calcam.GetAverageSector(j); MBadPixelsPix &bad = calcam.GetAverageBadSector(j); FitHiGainHists(hist,pix,bad,fittyp,osctyp); } } // -------------------------------------------------------------------------- // // Calls FitLoGainHists for every entry in: // - fLoGainArray // - fAverageLoGainAreas // - fAverageLoGainSectors // void MHCalibrationCam::FitLoGainArrays(MCalibrationCam &calcam, MBadPixelsCam &badcam, MBadPixelsPix::UncalibratedType_t fittyp, MBadPixelsPix::UncalibratedType_t osctyp) { if (!IsLoGain()) return; for (Int_t i=0; iGetSize(); i++) { MHCalibrationPix &hist = (*this)(i); if (hist.IsExcluded()) continue; MCalibrationPix &pix = calcam[i]; MBadPixelsPix &bad = badcam[i]; FitLoGainHists(hist,pix,bad,fittyp,osctyp); } if (!IsAverageing()) return; for (Int_t j=0; jGetSize(); j++) { MHCalibrationPix &hist = GetAverageLoGainArea(j); MCalibrationPix &pix = calcam.GetAverageArea(j); MBadPixelsPix &bad = calcam.GetAverageBadArea(j); FitLoGainHists(hist,pix,bad,fittyp,osctyp); } for (Int_t j=0; jGetSize(); j++) { MHCalibrationPix &hist = GetAverageLoGainSector(j); MCalibrationPix &pix = calcam.GetAverageSector(j); MBadPixelsPix &bad = calcam.GetAverageBadSector(j); FitLoGainHists(hist,pix,bad,fittyp,osctyp); } } //------------------------------------------------------------ // // For all averaged areas, the fitted sigma is multiplied with the square root of // the number involved pixels // void MHCalibrationCam::CalcAverageSigma() { if (!fCam) return; if (!IsAverageing()) return; for (UInt_t j=0; jGetNumAreas(); j++) { MCalibrationPix &pix = fCam->GetAverageArea(j); const Float_t numsqr = TMath::Sqrt((Float_t)fAverageAreaNum[j]); fAverageAreaSigma[j] = pix.GetSigma () * numsqr; fAverageAreaSigmaVar[j] = pix.GetSigmaErr () * pix.GetSigmaErr() * numsqr; pix.SetSigma (fAverageAreaSigma[j]); pix.SetSigmaVar(fAverageAreaSigmaVar[j]); fAverageAreaRelSigma [j] = fAverageAreaSigma[j] / pix.GetMean(); fAverageAreaRelSigmaVar[j] = fAverageAreaSigmaVar[j] / (fAverageAreaSigma[j]*fAverageAreaSigma[j]); fAverageAreaRelSigmaVar[j] += pix.GetMeanRelVar(); fAverageAreaRelSigmaVar[j] *= fAverageAreaRelSigma[j]; } } // --------------------------------------------------------------------------- // // Returns if the histogram is empty and sets the following flag: // - MBadPixelsPix::SetUnsuitable(MBadPixelsPix::kUnsuitableRun) // // Fits the histograms with a Gaussian, in case of failure // calls MHCalibrationPix::RepeatFit(), in case of repeated failure // calls MHCalibrationPix::BypassFit() and sets the following flags: // - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::UncalibratedType_t fittyp ) // - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun ) // // Creates the fourier spectrum and tests MHGausEvents::IsFourierSpectrumOK(). // In case no, sets the following flags: // - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::UncalibratedType_t osctyp ) // - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun ) // // Retrieves the results and store them in MCalibrationPix // void MHCalibrationCam::FitHiGainHists(MHCalibrationPix &hist, MCalibrationPix &pix, MBadPixelsPix &bad, MBadPixelsPix::UncalibratedType_t fittyp, MBadPixelsPix::UncalibratedType_t osctyp) { if (hist.IsEmpty() || hist.IsOnlyOverflow() || hist.IsOnlyUnderflow()) { *fLog << err << GetDescriptor() << ": Only overflow or underflow in high-gain pixel: " << pix.GetPixId() << endl; return; } // // 2) Fit the Hi Gain histograms with a Gaussian // if (!hist.FitGaus()) // // 3) In case of failure set the bit Fitted to false and take histogram means and RMS // if (!hist.RepeatFit()) { hist.BypassFit(); bad.SetUncalibrated( fittyp ); } // // 4) Check for oscillations // if (IsOscillations()) { hist.CreateFourierSpectrum(); if (!hist.IsFourierSpectrumOK()) bad.SetUncalibrated( osctyp ); } // // 5) Retrieve the results and store them in this class // pix.SetHiGainMean ( hist.GetMean() ); pix.SetHiGainMeanVar ( hist.GetMeanErr() * hist.GetMeanErr() ); pix.SetHiGainRms ( hist.GetHistRms() ); pix.SetHiGainSigma ( hist.GetSigma() ); pix.SetHiGainSigmaVar ( hist.GetSigmaErr()* hist.GetSigmaErr() ); pix.SetHiGainProb ( hist.GetProb() ); pix.SetHiGainNumBlackout( hist.GetBlackout() ); pix.SetHiGainNumPickup ( hist.GetPickup() ); if (IsDebug()) { *fLog << dbginf << GetDescriptor() << ": ID " << GetName() << " HiGainSaturation: " << pix.IsHiGainSaturation() << " HiGainMean: " << hist.GetMean () << " HiGainMeanErr: " << hist.GetMeanErr () << " HiGainMeanSigma: " << hist.GetSigma () << " HiGainMeanSigmaErr: " << hist.GetSigmaErr() << " HiGainMeanProb: " << hist.GetProb () << " HiGainNumBlackout: " << hist.GetBlackout() << " HiGainNumPickup : " << hist.GetPickup () << endl; } } // --------------------------------------------------------------------------- // // Returns if the histogram is empty and sets the following flag: // - MBadPixelsPix::SetUnsuitable(MBadPixelsPix::kUnsuitableRun) // // Fits the histograms with a Gaussian, in case of failure // calls MHCalibrationPix::RepeatFit(), in case of repeated failure // calls MHCalibrationPix::BypassFit() and sets the following flags: // - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::UncalibratedType_t fittyp ) // - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun ) // // Creates the fourier spectrum and tests MHGausEvents::IsFourierSpectrumOK(). // In case no, sets the following flags: // - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::UncalibratedType_t osctyp ) // - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun ) // // Retrieves the results and store them in MCalibrationPix // void MHCalibrationCam::FitLoGainHists(MHCalibrationPix &hist, MCalibrationPix &pix, MBadPixelsPix &bad, MBadPixelsPix::UncalibratedType_t fittyp, MBadPixelsPix::UncalibratedType_t osctyp) { if (hist.IsEmpty() || hist.IsOnlyOverflow() || hist.IsOnlyUnderflow()) { *fLog << err << GetDescriptor() << ": Only overflow or underflow in low-gain pixel: " << pix.GetPixId() << endl; return; } // // 2) Fit the Hi Gain histograms with a Gaussian // if (!hist.FitGaus()) // // 3) In case of failure set the bit Fitted to false and take histogram means and RMS // if (!hist.RepeatFit()) { hist.BypassFit(); if (pix.IsHiGainSaturation()) bad.SetUncalibrated( fittyp ); } // // 4) Check for oscillations // if (IsOscillations()) { hist.CreateFourierSpectrum(); if (!hist.IsFourierSpectrumOK()) bad.SetUncalibrated( osctyp ); } // // 5) Retrieve the results and store them in this class // pix.SetLoGainMean ( hist.GetMean() ); pix.SetLoGainMeanVar ( hist.GetMeanErr() * hist.GetMeanErr() ); pix.SetLoGainRms ( hist.GetHistRms() ); pix.SetLoGainSigma ( hist.GetSigma() ); pix.SetLoGainSigmaVar ( hist.GetSigmaErr() * hist.GetSigmaErr() ); pix.SetLoGainProb ( hist.GetProb() ); pix.SetLoGainNumBlackout( hist.GetBlackout() ); pix.SetLoGainNumPickup ( hist.GetPickup() ); if (IsDebug()) { *fLog << dbginf << GetDescriptor() << "ID: " << hist.GetName() << " HiGainSaturation: " << pix.IsHiGainSaturation() << " LoGainMean: " << hist.GetMean () << " LoGainMeanErr: " << hist.GetMeanErr () << " LoGainMeanSigma: " << hist.GetSigma () << " LoGainMeanSigmaErr: " << hist.GetSigmaErr() << " LoGainMeanProb: " << hist.GetProb () << " LoGainNumBlackout: " << hist.GetBlackout() << " LoGainNumPickup : " << hist.GetPickup () << endl; } } // ----------------------------------------------------------------------------- // // Default draw: // // Displays the averaged areas, both High Gain and Low Gain // // Calls the Draw of the fAverageHiGainAreas and fAverageLoGainAreas objects with options // void MHCalibrationCam::Draw(const Option_t *opt) { if (!IsAverageing()) return; const Int_t nareas = fAverageHiGainAreas->GetSize(); if (nareas == 0) return; TVirtualPad *pad = gPad ? gPad : MH::MakeDefCanvas(this); pad->SetBorderMode(0); pad->Divide(IsLoGain() ? 2 : 1,nareas); for (Int_t i=0; icd(IsLoGain() ? 2*(i+1)-1 : i+1); GetAverageHiGainArea(i).Draw(opt); if (!fAverageAreaSat[i]) DrawAverageSigma(fAverageAreaSat[i], i, fAverageAreaSigma[i], fAverageAreaSigmaVar[i], fAverageAreaRelSigma[i], fAverageAreaRelSigmaVar[i]); if (IsLoGain()) { pad->cd(2*(i+1)); GetAverageLoGainArea(i).Draw(opt); } if (fAverageAreaSat[i]) DrawAverageSigma(fAverageAreaSat[i], i, fAverageAreaSigma[i], fAverageAreaSigmaVar[i], fAverageAreaRelSigma[i], fAverageAreaRelSigmaVar[i]); } } // ----------------------------------------------------------------------------- // // Default draw: // // Displays a TPaveText with the re-normalized sigmas of the average area // void MHCalibrationCam::DrawAverageSigma(Bool_t sat, Bool_t inner, Float_t sigma, Float_t sigmavar, Float_t relsigma, Float_t relsigmavar) const { if (sigma != 0 && sigmavar >= 0 && relsigmavar >= 0.) { TPad *newpad = new TPad("newpad","transparent",0,0,1,1); newpad->SetFillStyle(4000); newpad->Draw(); newpad->cd(); TPaveText *text = new TPaveText(sat? 0.1 : 0.35,0.7,sat ? 0.4 : 0.7,1.0); text->SetTextSize(0.07); const TString line1 = Form("%s%s%s",inner ? "Outer" : "Inner", " Pixels ", sat ? "Low Gain" : "High Gain"); TText *txt1 = text->AddText(line1.Data()); const TString line2 = Form("#sigma per pix: %2.2f #pm %2.2f",sigma,TMath::Sqrt(sigmavar)); TText *txt2 = text->AddText(line2.Data()); const TString line3 = Form("Rel. #sigma per pix: %2.2f #pm %2.2f",relsigma,TMath::Sqrt(relsigmavar)); TText *txt3 = text->AddText(line3.Data()); text->Draw(""); text->SetBit(kCanDelete); txt1->SetBit(kCanDelete); txt2->SetBit(kCanDelete); txt3->SetBit(kCanDelete); newpad->SetBit(kCanDelete); } } Int_t MHCalibrationCam::ReadEnv(const TEnv &env, TString prefix, Bool_t print) { Bool_t rc = kFALSE; if (IsEnvDefined(env, prefix, "Debug", print)) { SetDebug(GetEnvValue(env, prefix, "Debug", IsDebug())); rc = kTRUE; } if (IsEnvDefined(env, prefix, "LoGain", print)) { SetDebug(GetEnvValue(env, prefix, "LoGain", IsLoGain())); rc = kTRUE; } if (IsEnvDefined(env, prefix, "Oscillations", print)) { SetOscillations(GetEnvValue(env, prefix, "Oscillations", IsOscillations())); rc = kTRUE; } if (IsEnvDefined(env, prefix, "SizeCheck", print)) { SetSizeCheck(GetEnvValue(env, prefix, "SizeCheck", IsSizeCheck())); rc = kTRUE; } if (IsEnvDefined(env, prefix, "Averageing", print)) { SetAverageing(GetEnvValue(env, prefix, "Averageing", IsAverageing())); rc = kTRUE; } if (IsEnvDefined(env, prefix, "Nbins", print)) { SetNbins(GetEnvValue(env, prefix, "Nbins", fNbins)); rc = kTRUE; } if (IsEnvDefined(env, prefix, "First", print)) { SetFirst(GetEnvValue(env, prefix, "First", fFirst)); rc = kTRUE; } if (IsEnvDefined(env, prefix, "Last", print)) { SetLast(GetEnvValue(env, prefix, "Last", fLast)); rc = kTRUE; } if (IsEnvDefined(env, prefix, "ProbLimit", print)) { SetProbLimit(GetEnvValue(env, prefix, "ProbLimit", fProbLimit)); rc = kTRUE; } if (IsEnvDefined(env, prefix, "OverflowLimit", print)) { SetOverflowLimit(GetEnvValue(env, prefix, "OverflowLimit", fOverflowLimit)); rc = kTRUE; } if (IsEnvDefined(env, prefix, "PulserFrequency", print)) { SetPulserFrequency(GetEnvValue(env, prefix, "PulserFrequency", fPulserFrequency)); rc = kTRUE; } return rc; }