Index: trunk/MagicSoft/Mars/Changelog =================================================================== --- trunk/MagicSoft/Mars/Changelog (revision 5748) +++ trunk/MagicSoft/Mars/Changelog (revision 5749) @@ -21,4 +21,12 @@ -*-*- END OF LINE -*-*- 2005/01/08 Markus Gaug + + * mhcalib/MHCalibrationHiLoCam.[h,cc] + * mcalib/MCalibrationHiLoCam.[h,cc] + * mcalib/MCalibrationHiLoPix.[h,cc] + * mcalib/MIntensityCalibrationHiLoCam.[h,cc] + * mhcalib/Makefile, mcalib/Makefile, mhcalib/HCalibLinkDef.h, + mcalib/CalibLinkDef.h + - new classes for calibration of conversion high-gain vs. low-gain * msignal/MExtractedSignalPix.[h,cc] Index: trunk/MagicSoft/Mars/mcalib/CalibLinkDef.h =================================================================== --- trunk/MagicSoft/Mars/mcalib/CalibLinkDef.h (revision 5748) +++ trunk/MagicSoft/Mars/mcalib/CalibLinkDef.h (revision 5749) @@ -18,4 +18,5 @@ #pragma link C++ class MCalibrationIntensityRelTimeCam+; #pragma link C++ class MCalibrationIntensityTestCam+; +#pragma link C++ class MCalibrationIntensityHiLoCam+; #pragma link C++ class MCalibrationCam+; #pragma link C++ class MCalibrationPix+; @@ -24,4 +25,6 @@ #pragma link C++ class MCalibrationRelTimePix+; #pragma link C++ class MCalibrationPedCam+; +#pragma link C++ class MCalibrationHiLoCam+; +#pragma link C++ class MCalibrationHiLoPix+; #pragma link C++ class MCalibrationQECamMagic+; #pragma link C++ class MCalibrationQECam+; Index: trunk/MagicSoft/Mars/mcalib/MCalibrationHiLoCam.cc =================================================================== --- trunk/MagicSoft/Mars/mcalib/MCalibrationHiLoCam.cc (revision 5749) +++ trunk/MagicSoft/Mars/mcalib/MCalibrationHiLoCam.cc (revision 5749) @@ -0,0 +1,256 @@ +/* ======================================================================== *\ +! +! * +! * 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-2004 +! +! +\* ======================================================================== */ + +///////////////////////////////////////////////////////////////////////////// +// +// MCalibrationHiLoCam +// +// Storage container for ratio between higain and logain charge extraction +// calibration results of the whole camera. +// +// Individual pixels have to be cast when retrieved e.g.: +// MCalibrationHiLoPix &avpix = (MCalibrationHiLoPix&)(*fHiLoCam)[i] +// +// The following "calibration" constants can be retrieved from each pixel: +// - GetHiLoRatio(): The mean signal ratio between high-gain and low-gain +// - GetHiLoRatioErr(): The Gauss sigma of histogrammed signal ratios +// +// Averaged values over one whole area index (e.g. inner or outer pixels for +// the MAGIC camera), can be retrieved via: +// MCalibrationHiLoPix &avpix = (MCalibrationHiLoPix&)fRelCam->GetAverageArea(i) +// +// Averaged values over one whole camera sector can be retrieved via: +// MCalibrationHiLoPix &avpix = (MCalibrationHiLoPix&)fRelCam->GetAverageSector(i) +// +// Note the averageing has been done on an event-by-event basis. Resulting +// Sigma's of the Gauss fit have been multiplied with the square root of the number +// of involved pixels in order to make a direct comparison possible with the mean of +// sigmas. +// +// See also: MHCalibrationHiLoPix, MHCalibrationHiLoCam +// +// The calculated values (types of GetPixelContent) are: +// +// Fitted values: +// ============== +// +// 0: HiLo Ratio +// 1: Error of HiLo Ratio +// 2: Sigma of HiLo Ratio == Resolution +// 3: Error of Sigma of HiLo Ratio +// +// Useful variables derived from the fit results: +// ============================================= +// +// 4: Returned probability of Gauss fit to distribution +// +///////////////////////////////////////////////////////////////////////////// +#include "MCalibrationHiLoCam.h" +#include "MCalibrationCam.h" + +#include + +#include "MLog.h" +#include "MLogManip.h" + +#include "MGeomCam.h" +#include "MGeomPix.h" + +#include "MCalibrationHiLoPix.h" + +ClassImp(MCalibrationHiLoCam); + +using namespace std; + +// -------------------------------------------------------------------------- +// +// Default constructor. +// +MCalibrationHiLoCam::MCalibrationHiLoCam(const char *name, const char *title) +{ + + fName = name ? name : "MCalibrationHiLoCam"; + fTitle = title ? title : "Container for High Gain vs. Low Gain amplification ratio"; + +} + +void MCalibrationHiLoCam::Add(const UInt_t a, const UInt_t b) +{ + for (UInt_t i=a; iAddAt(new MCalibrationHiLoPix,i); +} + + +void MCalibrationHiLoCam::AddArea(const UInt_t a, const UInt_t b) +{ + for (UInt_t i=a; iAddAt(new MCalibrationHiLoPix,i); +} + +void MCalibrationHiLoCam::AddSector(const UInt_t a, const UInt_t b) +{ + for (UInt_t i=a; iAddAt(new MCalibrationHiLoPix,i); +} + +// -------------------------------------------------------------------------- +// +// Print first the well fitted pixels +// and then the ones which are not Valid +// +void MCalibrationHiLoCam::Print(Option_t *o) const +{ + + *fLog << all << GetDescriptor() << ":" << endl; + int id = 0; + + *fLog << all << "Calibrated pixels:" << endl; + *fLog << all << endl; + + TIter Next(fPixels); + MCalibrationHiLoPix *pix; + while ((pix=(MCalibrationHiLoPix*)Next())) + { + + if (!pix->IsExcluded()) + { + + *fLog << all + << Form("%s%4i%s%4.2f%s%4.2f%s%4.2f%s%4.2f","Pix ",pix->GetPixId(), + ": Ratio: ",pix->GetHiLoRatio()," +- ",pix->GetHiLoRatioErr(), + " Sigma: ",pix->GetHiLoSigma()," +- ",pix->GetHiLoSigmaErr()) + << endl; + id++; + } + } + + *fLog << all << id << " pixels" << endl; + id = 0; + + + *fLog << all << endl; + *fLog << all << "Excluded pixels:" << endl; + *fLog << all << endl; + + id = 0; + + TIter Next4(fPixels); + while ((pix=(MCalibrationHiLoPix*)Next4())) + { + if (pix->IsExcluded()) + { + *fLog << all << pix->GetPixId() << endl; + id++; + } + } + *fLog << all << id << " Excluded pixels " << endl; + *fLog << endl; + + TIter Next5(fAverageAreas); + while ((pix=(MCalibrationHiLoPix*)Next5())) + { + *fLog << all + << Form("%s%4i%s%4.2f%s%4.2f%s%4.2f%s%4.2f","Average Area ",pix->GetPixId(), + ": Ratio: ",pix->GetHiLoRatio()," +- ",pix->GetHiLoRatioErr(), + " Sigma: ",pix->GetHiLoSigma()," +- ",pix->GetHiLoSigmaErr()) + << endl; + } + + TIter Next6(fAverageSectors); + while ((pix=(MCalibrationHiLoPix*)Next5())) + { + *fLog << all + << Form("%s%4i%s%4.2f%s%4.2f%s%4.2f%s%4.2f","Average Sector ",pix->GetPixId(), + ": Ratio: ",pix->GetHiLoRatio()," +- ",pix->GetHiLoRatioErr(), + " Sigma: ",pix->GetHiLoSigma()," +- ",pix->GetHiLoSigmaErr()) + << endl; + } +} + + +// -------------------------------------------------------------------------- +// +// The types are as follows: +// +// Fitted values: +// ============== +// +// 0: Fitted HiLo +// 1: Error of fitted HiLo +// 2: Sigma of fitted HiLo +// 3: Error of Sigma of fitted HiLo +// +// Useful variables derived from the fit results: +// ============================================= +// +// 4: Returned probability of Gauss fit to HiLo distribution +// +Bool_t MCalibrationHiLoCam::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const +{ + + if (idx > GetSize()) + return kFALSE; + + Float_t area = cam[idx].GetA(); + + if (area == 0) + return kFALSE; + + MCalibrationHiLoPix &pix = (MCalibrationHiLoPix&)(*this)[idx]; + + switch (type) + { + case 0: + if (pix.IsExcluded()) + return kFALSE; + val = pix.GetMean(); + break; + case 1: + if (pix.IsExcluded()) + return kFALSE; + val = pix.GetMeanErr(); + break; + case 2: + if (pix.IsExcluded()) + return kFALSE; + val = pix.GetSigma(); + break; + case 3: + if (pix.IsExcluded()) + return kFALSE; + val = pix.GetSigmaErr(); + break; + case 4: + if (pix.IsExcluded()) + return kFALSE; + val = pix.GetProb(); + break; + default: + return kFALSE; + } + + return val!=-1.; + +} + Index: trunk/MagicSoft/Mars/mcalib/MCalibrationHiLoCam.h =================================================================== --- trunk/MagicSoft/Mars/mcalib/MCalibrationHiLoCam.h (revision 5749) +++ trunk/MagicSoft/Mars/mcalib/MCalibrationHiLoCam.h (revision 5749) @@ -0,0 +1,29 @@ +#ifndef MARS_MCalibrationHiLoCam +#define MARS_MCalibrationHiLoCam + +#ifndef MARS_MCalibrationCam +#include "MCalibrationCam.h" +#endif + +class MCalibrationHiLoCam : public MCalibrationCam +{ +private: + + void Add(const UInt_t a, const UInt_t b); + void AddArea(const UInt_t a, const UInt_t b); + void AddSector(const UInt_t a, const UInt_t b); + +public: + + MCalibrationHiLoCam(const char *name=NULL, const char *title=NULL); + + // Prints + void Print(Option_t *o="") const; + + // Others + Bool_t GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type=0) const; + + ClassDef(MCalibrationHiLoCam, 1) // Container Rel. Arrival Time Calibration Results Camera +}; + +#endif Index: trunk/MagicSoft/Mars/mcalib/MCalibrationHiLoPix.cc =================================================================== --- trunk/MagicSoft/Mars/mcalib/MCalibrationHiLoPix.cc (revision 5749) +++ trunk/MagicSoft/Mars/mcalib/MCalibrationHiLoPix.cc (revision 5749) @@ -0,0 +1,54 @@ +/* ======================================================================== *\ +! +! * +! * 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 +! +\* ======================================================================== */ +///////////////////////////////////////////////////////////////////////////// +// // +// MCalibrationHiLoPix // +// // +// Storage container for high-gain vs. low-gain charge calibration results // +// of one Pixel (PMT). // +// The following "calibration" constants can be retrieved: // +// - GetHiLoRatio(): The mean conversion High-gain vs. Low-gain +// with which the low-gain result has to be multiplied +// - GetHiLoSigma(): The Gauss sigma of histogrammed High-gain vs. Low-gain ratios +// +// See also: MHCalibrationHiLoPix, MHCalibrationHiLoCam // +// // +///////////////////////////////////////////////////////////////////////////// +#include "MCalibrationHiLoPix.h" + +ClassImp(MCalibrationHiLoPix); + +using namespace std; + +// -------------------------------------------------------------------------- +// +// Default Constructor +// +MCalibrationHiLoPix::MCalibrationHiLoPix(const char *name, const char *title) +{ + + fName = name ? name : "MCalibrationHiLoPix"; + fTitle = title ? title : "Results of MHCalibrationHiLoPix "; + +} + Index: trunk/MagicSoft/Mars/mcalib/MCalibrationHiLoPix.h =================================================================== --- trunk/MagicSoft/Mars/mcalib/MCalibrationHiLoPix.h (revision 5749) +++ trunk/MagicSoft/Mars/mcalib/MCalibrationHiLoPix.h (revision 5749) @@ -0,0 +1,26 @@ +#ifndef MARS_MCalibrationHiLoPix +#define MARS_MCalibrationHiLoPix + +#ifndef MARS_MCalibrationPix +#include "MCalibrationPix.h" +#endif + +class MCalibrationHiLoPix : public MCalibrationPix +{ +private: + +public: + + MCalibrationHiLoPix(const char *name=NULL, const char *title=NULL); + ~MCalibrationHiLoPix() {} + + Float_t GetHiLoRatio() const { return GetMean(); } + Float_t GetHiLoRatioErr() const { return GetMeanErr(); } + Float_t GetHiLoSigma() const { return GetSigma(); } + Float_t GetHiLoSigmaErr() const { return GetSigmaErr(); } + + ClassDef(MCalibrationHiLoPix, 1) // Container HiLo conversion Calibration Results Pixel +}; + +#endif + Index: trunk/MagicSoft/Mars/mcalib/MCalibrationIntensityHiLoCam.cc =================================================================== --- trunk/MagicSoft/Mars/mcalib/MCalibrationIntensityHiLoCam.cc (revision 5749) +++ trunk/MagicSoft/Mars/mcalib/MCalibrationIntensityHiLoCam.cc (revision 5749) @@ -0,0 +1,252 @@ +/* ======================================================================== *\ +! +! * +! * 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-2004 +! +! +\* ======================================================================== */ +///////////////////////////////////////////////////////////////////////////// +// +// MCalibrationIntensityHiLoCam +// +// Storage container for intensity charge calibration results. +// +// Individual MCalibrationHiLoCam's can be retrieved with: +// - GetCam() yielding the current cam. +// - GetCam("name") yielding the current camera with name "name". +// - GetCam(i) yielding the i-th camera. +// +// See also: MCalibrationIntensityCam, MCalibrationHiLoCam, +// MCalibrationHiLoPix, MCalibrationHiLoCalc, MCalibrationQECam +// MHCalibrationHiLoPix, MHCalibrationHiLoCam +// +///////////////////////////////////////////////////////////////////////////// +#include "MCalibrationIntensityHiLoCam.h" +#include "MCalibrationHiLoCam.h" +#include "MCalibrationChargeCam.h" +#include "MCalibrationHiLoPix.h" +#include "MCalibrationChargePix.h" + +#include "MGeomCam.h" +#include "MGeomPix.h" + +#include "MLogManip.h" + +#include +#include + +ClassImp(MCalibrationIntensityHiLoCam); + +using namespace std; + +// -------------------------------------------------------------------------- +// +// Default constructor. +// +MCalibrationIntensityHiLoCam::MCalibrationIntensityHiLoCam(const char *name, const char *title) +{ + + fName = name ? name : "MCalibrationIntensityHiLoCam"; + fTitle = title ? title : "Results of the Intensity Calibration"; + + InitSize(1); +} + +// ------------------------------------------------------------------- +// +// Add MCalibrationHiLoCam's in the ranges from - to. +// +void MCalibrationIntensityHiLoCam::Add(const UInt_t from, const UInt_t to) +{ + for (UInt_t i=from; iAddAt(new MCalibrationHiLoCam,i); +} + +TGraphErrors *MCalibrationIntensityHiLoCam::GetHiLoRatioVsCharge( const UInt_t pixid, const MCalibrationIntensityChargeCam &chargecam, + const MCalibrationCam::PulserColor_t col) +{ + + if (chargecam.GetSize() != GetSize()) + { + *fLog << err << GetDescriptor() << ": Size mismatch between MCalibrationIntensityHiLoCam " + << "and MCalibrationIntensityChargeCam. " << endl; + return NULL; + } + + Int_t size = CountNumEntries(col); + + if (size == 0) + return NULL; + + if (size != chargecam.CountNumEntries(col)) + { + *fLog << err << GetDescriptor() << ": Size mismatch in colour between MCalibrationIntensityHiLoCam " + << "and MCalibrationIntensityChargeCam. " << endl; + return NULL; + } + + const Int_t nvalid = chargecam.CountNumValidEntries(pixid,col); + + if (nvalid == 0) + { + *fLog << err << GetDescriptor() << ": Only un-calibrated events in pixel: " << pixid << endl; + return NULL; + } + + TArrayF ratio(nvalid); + TArrayF ratioerr(nvalid); + TArrayF sig(nvalid); + TArrayF sigerr(nvalid); + + Int_t cnt = 0; + + for (Int_t i=0;iGetPulserColor() != col) + continue; + // + // Get the calibration pix from the calibration cam + // + MCalibrationChargePix &pix = (MCalibrationChargePix&)(*cam)[pixid]; + MCalibrationHiLoPix &relpix = (MCalibrationHiLoPix&)(*hilocam)[pixid]; + // + // Don't use bad pixels + // + if (!pix.IsFFactorMethodValid()) + continue; + // + ratio[cnt] = relpix.GetHiLoRatio(); + ratioerr[cnt] = relpix.GetHiLoRatioErr(); + // + sig [cnt] = pix.GetPheFFactorMethod(); + sigerr[cnt] = pix.GetPheFFactorMethodErr(); + cnt++; + } + + TGraphErrors *gr = new TGraphErrors(nvalid, + sig.GetArray(),ratio.GetArray(), + sigerr.GetArray(),ratioerr.GetArray()); + gr->SetTitle(Form("%s%3i","Pixel ",pixid)); + gr->GetXaxis()->SetTitle(" [1]"); + gr->GetYaxis()->SetTitle("HiLo Ratio [1]"); + return gr; +} + + +TGraphErrors *MCalibrationIntensityHiLoCam::GetHiLoRatioVsChargePerArea( const Int_t aidx,const MCalibrationIntensityChargeCam &chargecam, const MGeomCam &geom, const MCalibrationCam::PulserColor_t col) +{ + + Int_t size = CountNumEntries(col); + + if (size == 0) + return NULL; + + if (size != chargecam.CountNumEntries(col)) + { + *fLog << err << GetDescriptor() << ": Size mismatch in colour between MCalibrationIntensityHiLoCam " + << "and MCalibrationIntensityChargeCam. " << endl; + return NULL; + } + + TArrayF ratio(size); + TArrayF ratioerr(size); + TArrayF sig(size); + TArrayF sigerr(size); + + Int_t cnt = 0; + + for (Int_t i=0;iGetPulserColor() != col) + continue; + + const MCalibrationChargePix &apix = (MCalibrationChargePix&)cam->GetAverageArea(aidx); + const Float_t phe = apix.GetPheFFactorMethod(); + const Float_t pheerr = apix.GetPheFFactorMethodErr(); + + sig[cnt] = phe; + sigerr[cnt] = pheerr; + + Double_t ratiool = 0.; + Double_t ratiool2 = 0.; + Double_t var = 0.; + Int_t num = 0; + // + // Get the area calibration pix from the calibration cam + // + for (Int_t i=0; iGetSize(); i++) + { + const MCalibrationChargePix &pix = (MCalibrationChargePix&)(*cam)[i]; + const MCalibrationHiLoPix &relpix = (MCalibrationHiLoPix&)(*hilocam)[i]; + // + // Don't use bad pixels + // + if (!pix.IsFFactorMethodValid()) + continue; + // + // + if (aidx != geom[i].GetAidx()) + continue; + + ratiool += relpix.GetHiLoRatio(); + ratiool2 += relpix.GetHiLoRatio()*relpix.GetHiLoRatio(); + num++; + } + + if (num > 1) + { + ratiool /= num; + var = (ratiool2 - ratiool*ratiool*num) / (num-1); + + ratio[cnt] = ratiool; + if (var > 0.) + ratioerr[cnt] = TMath::Sqrt(var); + else + ratioerr[cnt] = 0.; + } + else + { + ratio[cnt] = -1.; + ratioerr[cnt] = 0.; + } + cnt++; + } + + TGraphErrors *gr = new TGraphErrors(size, + sig.GetArray(),ratio.GetArray(), + sigerr.GetArray(),ratioerr.GetArray()); + gr->SetTitle(Form("%s%3i","Area Index ",aidx)); + gr->GetXaxis()->SetTitle(" [1]"); + gr->GetYaxis()->SetTitle("HiLo Ratio [1]"); + return gr; +} Index: trunk/MagicSoft/Mars/mcalib/MCalibrationIntensityHiLoCam.h =================================================================== --- trunk/MagicSoft/Mars/mcalib/MCalibrationIntensityHiLoCam.h (revision 5749) +++ trunk/MagicSoft/Mars/mcalib/MCalibrationIntensityHiLoCam.h (revision 5749) @@ -0,0 +1,32 @@ +#ifndef MARS_MCalibrationIntensityHiLoCam +#define MARS_MCalibrationIntensityHiLoCam + +#ifndef MARS_MCalibrationIntensityCam +#include "MCalibrationIntensityCam.h" +#endif + +#ifndef MARS_MCalibrationIntensityChargeCam +#include "MCalibrationIntensityChargeCam.h" +#endif + +class TGraphErrors; +class MCalibrationIntensityHiLoCam : public MCalibrationIntensityCam +{ +private: + + void Add(const UInt_t from, const UInt_t to); + +public: + + MCalibrationIntensityHiLoCam(const char *name=NULL, const char *title=NULL); + + TGraphErrors *GetHiLoRatioVsCharge( const UInt_t pixid, const MCalibrationIntensityChargeCam &chargecam, + const MCalibrationCam::PulserColor_t col=MCalibrationCam::kNONE); + TGraphErrors *GetHiLoRatioVsChargePerArea( const Int_t aidx, const MCalibrationIntensityChargeCam &chargecam, + const MGeomCam &geom, + const MCalibrationCam::PulserColor_t col=MCalibrationCam::kNONE); + + ClassDef(MCalibrationIntensityHiLoCam, 1) // Container Intensity HiLo Calibration Results Camera +}; + +#endif Index: trunk/MagicSoft/Mars/mcalib/Makefile =================================================================== --- trunk/MagicSoft/Mars/mcalib/Makefile (revision 5748) +++ trunk/MagicSoft/Mars/mcalib/Makefile (revision 5749) @@ -42,4 +42,5 @@ MCalibrationIntensityRelTimeCam.cc \ MCalibrationIntensityTestCam.cc \ + MCalibrationIntensityHiLoCam.cc \ MCalibrationCam.cc \ MCalibrationPix.cc \ @@ -54,4 +55,6 @@ MCalibrationTestCam.cc \ MCalibrationTestPix.cc \ + MCalibrationHiLoCam.cc \ + MCalibrationHiLoPix.cc \ MCalibrationChargeCalc.cc \ MCalibrationChargeCam.cc \ Index: trunk/MagicSoft/Mars/mhcalib/HCalibLinkDef.h =================================================================== --- trunk/MagicSoft/Mars/mhcalib/HCalibLinkDef.h (revision 5748) +++ trunk/MagicSoft/Mars/mhcalib/HCalibLinkDef.h (revision 5749) @@ -14,4 +14,5 @@ #pragma link C++ class MHCalibrationRelTimeCam+; #pragma link C++ class MHCalibrationTestCam+; +#pragma link C++ class MHCalibrationHiLoCam+; #pragma link C++ class MHPedestalCam+; #pragma link C++ class MHCalibrationTestTimeCam+; Index: trunk/MagicSoft/Mars/mhcalib/MHCalibrationHiLoCam.cc =================================================================== --- trunk/MagicSoft/Mars/mhcalib/MHCalibrationHiLoCam.cc (revision 5749) +++ trunk/MagicSoft/Mars/mhcalib/MHCalibrationHiLoCam.cc (revision 5749) @@ -0,0 +1,559 @@ +/* ======================================================================== *\ +! +! * +! * 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 +! +! +\* ======================================================================== */ +///////////////////////////////////////////////////////////////////////////// +// +// MHCalibrationHiLoCam +// +// Fills the extracted high-gain low-gain charge ratios of MArrivalTimeCam into +// the MHCalibrationPix-classes MHCalibrationPix for every: +// +// - Pixel, stored in the TObjArray's MHCalibrationCam::fHiGainArray +// or MHCalibrationCam::fHiGainArray, respectively, depending if +// MArrivalTimePix::IsLoGainUsed() is set. +// +// - Average pixel per AREA index (e.g. inner and outer for the MAGIC camera), +// stored in the TObjArray's MHCalibrationCam::fAverageHiGainAreas and +// MHCalibrationCam::fAverageHiGainAreas +// +// - Average pixel per camera SECTOR (e.g. sectors 1-6 for the MAGIC camera), +// stored in the TObjArray's MHCalibrationCam::fAverageHiGainSectors +// and MHCalibrationCam::fAverageHiGainSectors +// +// The 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::kHiLoNotFitted ) 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) +// +// The class also fills arrays with the signal vs. event number, creates a fourier +// spectrum (see MHGausEvents::CreateFourierSpectrum()) 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::kHiLoOscillating ) and +// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun ) +// +// This same procedure is performed for the average pixels. +// +// The following results are written into MCalibrationHiLoCam: +// +// - MCalibrationPix::SetMean() +// - MCalibrationPix::SetMeanErr() +// - MCalibrationPix::SetSigma() +// - MCalibrationPix::SetSigmaErr() +// - MCalibrationPix::SetProb() +// - MCalibrationPix::SetNumPickup() +// +// 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 "MHCalibrationHiLoCam.h" +#include "MHCalibrationPix.h" + +#include "MLog.h" +#include "MLogManip.h" + +#include "MParList.h" + +#include "MCalibrationIntensityHiLoCam.h" + +#include "MCalibrationHiLoCam.h" +#include "MCalibrationCam.h" +#include "MCalibrationPix.h" + +#include "MExtractedSignalCam.h" +#include "MExtractedSignalPix.h" + +#include "MRawEvtPixelIter.h" + +#include "MGeomCam.h" +#include "MGeomPix.h" + +#include "MBadPixelsIntensityCam.h" +#include "MBadPixelsCam.h" +#include "MBadPixelsPix.h" + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +ClassImp(MHCalibrationHiLoCam); + +using namespace std; + +const Int_t MHCalibrationHiLoCam::fgNbins = 200; +const Axis_t MHCalibrationHiLoCam::fgFirst = 0.; +const Axis_t MHCalibrationHiLoCam::fgLast = 20.; +const Float_t MHCalibrationHiLoCam::fgProbLimit = 0.; +const TString MHCalibrationHiLoCam::gsHistName = "HiLo"; +const TString MHCalibrationHiLoCam::gsHistTitle = "HiGain vs. LoGain"; +const TString MHCalibrationHiLoCam::gsHistXTitle = "Amplification Ratio [1]"; +const TString MHCalibrationHiLoCam::gsHistYTitle = "Nr. events"; +const Byte_t MHCalibrationHiLoCam::fgLowerLim = 200; +const Byte_t MHCalibrationHiLoCam::fgUpperLim = 252; +// -------------------------------------------------------------------------- +// +// Default Constructor. +// +// Sets: +// - fNbins to fgNbins +// - fFirst to fgFirst +// - fLast to fgLast +// +// - fHistName to gsHistName +// - fHistTitle to gsHistTitle +// - fHistXTitle to gsHistXTitle +// - fHistYTitle to gsHistYTitle +// +// - fLowerLimt to fgLowerLim +// - fUpperLimt to fgUpperLim +// +MHCalibrationHiLoCam::MHCalibrationHiLoCam(const char *name, const char *title) +{ + + fName = name ? name : "MHCalibrationHiLoCam"; + fTitle = title ? title : "Histogram class for the high-gain vs. low-gain amplification ratio calibration"; + + SetNbins(fgNbins); + SetFirst(fgFirst); + SetLast (fgLast ); + + SetProbLimit(fgProbLimit); + + SetHistName (gsHistName .Data()); + SetHistTitle (gsHistTitle .Data()); + SetHistXTitle(gsHistXTitle.Data()); + SetHistYTitle(gsHistYTitle.Data()); + + SetLowerLim(); + SetUpperLim(); +} + +// -------------------------------------------------------------------------- +// +// Creates new MHCalibrationHiLoCam only with the averaged areas: +// the rest has to be retrieved directly, e.g. via: +// MHCalibrationHiLoCam *cam = MParList::FindObject("MHCalibrationHiLoCam"); +// - cam->GetAverageSector(5).DrawClone(); +// - (*cam)[100].DrawClone() +// +TObject *MHCalibrationHiLoCam::Clone(const char *) const +{ + + MHCalibrationHiLoCam *cam = new MHCalibrationHiLoCam(); + + // + // Copy the data members + // + cam->fColor = fColor; + cam->fRunNumbers = fRunNumbers; + cam->fPulserFrequency = fPulserFrequency; + cam->fFlags = fFlags; + cam->fNbins = fNbins; + cam->fFirst = fFirst; + cam->fLast = fLast; + + // + // Copy the MArrays + // + cam->fAverageAreaRelSigma = fAverageAreaRelSigma; + cam->fAverageAreaRelSigmaVar = fAverageAreaRelSigmaVar; + cam->fAverageAreaSat = fAverageAreaSat; + cam->fAverageAreaSigma = fAverageAreaSigma; + cam->fAverageAreaSigmaVar = fAverageAreaSigmaVar; + cam->fAverageAreaNum = fAverageAreaNum; + cam->fAverageSectorNum = fAverageSectorNum; + + if (!IsAverageing()) + return cam; + + const Int_t navhi = fAverageHiGainAreas->GetSize(); + + for (int i=0; ifAverageHiGainAreas->AddAt(GetAverageHiGainArea(i).Clone(),i); + + return cam; +} + +// -------------------------------------------------------------------------- +// +// Gets the pointers to: +// - MRawEvtData +// +Bool_t MHCalibrationHiLoCam::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: +// - MCalibrationHiLoCam +// +// Searches pointer to: +// - MExtractedSignalCam +// +// Calls: +// - MHCalibrationCam::InitHiGainArrays() +// +// Sets: +// - SetLoGain(kFALSE); +// - fSumarea to nareas +// - fSumsector to nareas +// - fNumarea to nareas +// - fNumsector to nareas +// +Bool_t MHCalibrationHiLoCam::ReInitHists(MParList *pList) +{ + + fIntensCam = (MCalibrationIntensityCam*)pList->FindObject(AddSerialNumber("MCalibrationIntensityHiLoCam")); + if (fIntensCam) + *fLog << inf << "Found MCalibrationIntensityHiLoCam ... " << endl; + else + { + fCam = (MCalibrationCam*)pList->FindObject(AddSerialNumber("MCalibrationHiLoCam")); + if (!fCam) + { + fCam = (MCalibrationCam*)pList->FindCreateObj(AddSerialNumber("MCalibrationHiLoCam")); + if (!fCam) + { + *fLog << err << "Cannot find nor create MCalibrationHiLoCam ... abort." << endl; + return kFALSE; + } + fCam->Init(*fGeom); + } + } + + MExtractedSignalCam *signal = (MExtractedSignalCam*)pList->FindObject("MExtractedSignalCam"); + if (!signal) + { + *fLog << err << "MExtractedSignalCam not found... abort." << endl; + return kFALSE; + } + + const Int_t npixels = fGeom->GetNumPixels(); + const Int_t nsectors = fGeom->GetNumSectors(); + const Int_t nareas = fGeom->GetNumAreas(); + + InitHiGainArrays(npixels,nareas,nsectors); + + fSumarea .Set(nareas); + fSumsector.Set(nsectors); + fNumarea .Set(nareas); + fNumsector.Set(nsectors); + + SetLoGain(kFALSE); + + return kTRUE; +} + +// ------------------------------------------------------------------------------- +// +// Retrieves pointer to MExtractedSignalCam: +// +// Retrieves from MGeomCam: +// - number of pixels +// - number of pixel areas +// - number of sectors +// +// Fills histograms (MHGausEvents::FillHistAndArray()) with: +// - MExtractedSignalPix::GetExtractedSignalHiGain(pixid) / MExtractedSignalPix::GetExtractedSignalLoGain; +// if the high-gain signal lies in between the limits: fLowerLim and fUpperLim +// +Bool_t MHCalibrationHiLoCam::FillHists(const MParContainer *par, const Stat_t w) +{ + + MExtractedSignalCam *signal = (MExtractedSignalCam*)par; + if (!signal) + { + gLog << err << "No argument in MExtractedSignal::Fill... abort." << endl; + return kFALSE; + } + + const Int_t nareas = fGeom->GetNumAreas(); + const Int_t nsectors = fGeom->GetNumSectors(); + + fSumarea .Reset(); + fSumsector.Reset(); + fNumarea .Reset(); + fNumsector.Reset(); + + MRawEvtPixelIter pixel(fRawEvt); + + while (pixel.Next()) + { + + const Byte_t max = pixel.GetMaxHiGainSample(); + + if (max < fLowerLim || max > fUpperLim) + continue; + + const UInt_t pixidx = pixel.GetPixelId(); + + MHCalibrationPix &hist = (*this)[pixidx]; + + if (hist.IsExcluded()) + continue; + + const MExtractedSignalPix &pix = (*signal)[pixidx]; + const Int_t aidx = (*fGeom)[pixidx].GetAidx(); + const Int_t sector = (*fGeom)[pixidx].GetSector(); + + const Float_t ratio = pix.GetExtractedSignalHiGain() / pix.GetExtractedSignalLoGain(); + + hist.FillHistAndArray(ratio) ; + fSumarea [aidx] += ratio; + fNumarea [aidx] ++; + fSumsector[sector] += ratio; + fNumsector[sector] ++; + } + + for (Int_t j=0; jGetCam() : fCam; + MBadPixelsCam *badcam = fIntensBad ? fIntensBad->GetCam() : fBadPixels; + + const Int_t nareas = fAverageHiGainAreas->GetSize(); + const Int_t nsectors = fAverageHiGainSectors->GetSize(); + + TArrayI satarea(nareas); + TArrayI satsect(nsectors); + fNumarea .Reset(); + fNumsector.Reset(); + + for (Int_t i=0; iGetSize(); i++) + { + + MHCalibrationPix &hist = (*this)[i]; + + if (hist.IsExcluded()) + continue; + + const Int_t aidx = (*fGeom)[i].GetAidx(); + const Int_t sector = (*fGeom)[i].GetSector(); + + fNumarea[aidx]++; + fNumsector[sector]++; + // + // Check histogram overflow + // + CheckOverflow(hist); + } + + for (Int_t j=0; jGetSize(); j++) + { + + MHCalibrationPix &hist = GetAverageHiGainSector(j); + // + // Check histogram overflow + // + CheckOverflow(hist); + } + + FitHiGainArrays(*hilocam,*badcam, + MBadPixelsPix::kHiLoNotFitted, + MBadPixelsPix::kHiLoOscillating); + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Sets all pixels to MBadPixelsPix::kUnreliableRun, if following flags are set: +// - MBadPixelsPix::kHiLoNotFitted +// - MBadPixelsPix::kHiLoOscillating +// +void MHCalibrationHiLoCam::FinalizeBadPixels() +{ + + MBadPixelsCam *badcam = fIntensBad ? fIntensBad->GetCam() : fBadPixels; + + for (Int_t i=0; iGetSize(); i++) + { + MBadPixelsPix &bad = (*badcam)[i]; + + if (bad.IsUncalibrated( MBadPixelsPix::kHiLoNotFitted )) + bad.SetUnsuitable( MBadPixelsPix::kUnreliableRun ); + + if (bad.IsUncalibrated( MBadPixelsPix::kHiLoOscillating)) + bad.SetUnsuitable( MBadPixelsPix::kUnreliableRun ); + + } +} + +// -------------------------------------------------------------------------- +// +// The types are as follows: +// +// Fitted values: +// ============== +// +// 0: Fitted Mean High-Gain Low-Gain Charge Ratio in FADC slices (MHGausEvents::GetMean() +// 1: Error Mean High-Gain Low-Gain Charge Ratio in FADC slices (MHGausEvents::GetMeanErr() +// 2: Sigma fitted High-Gain Low-Gain Charge Ratio in FADC slices (MHGausEvents::GetSigma() +// 3: Error Sigma High-Gain Low-Gain Charge Ratio in FADC slices (MHGausEvents::GetSigmaErr() +// +// Useful variables derived from the fit results: +// ============================================= +// +// 4: Returned probability of Gauss fit (calls: MHGausEvents::GetProb()) +// +// Localized defects: +// ================== +// +// 5: Gaus fit not OK (calls: MHGausEvents::IsGausFitOK()) +// 6: Fourier spectrum not OK (calls: MHGausEvents::IsFourierSpectrumOK()) +// +Bool_t MHCalibrationHiLoCam::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const +{ + + if (fHiGainArray->GetSize() <= idx) + return kFALSE; + + const MHCalibrationPix &pix = (*this)[idx]; + + switch (type) + { + case 0: + val = pix.GetMean(); + break; + case 1: + val = pix.GetMeanErr(); + break; + case 2: + val = pix.GetSigma(); + break; + case 3: + val = pix.GetSigmaErr(); + break; + case 4: + val = pix.GetProb(); + break; + case 5: + if (!pix.IsGausFitOK()) + val = 1.; + break; + case 6: + if (!pix.IsFourierSpectrumOK()) + val = 1.; + break; + default: + return kFALSE; + } + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Calls MHCalibrationPix::DrawClone() for pixel idx +// +void MHCalibrationHiLoCam::DrawPixelContent(Int_t idx) const +{ + (*this)[idx].DrawClone(); +} + +void MHCalibrationHiLoCam::CheckOverflow( MHCalibrationPix &pix ) +{ + + if (pix.IsExcluded()) + return; + + TH1F *hist = pix.GetHGausHist(); + + Stat_t overflow = hist->GetBinContent(hist->GetNbinsX()+1); + if (overflow > 0.0005*hist->GetEntries()) + { + *fLog << warn << "Hist-overflow " << overflow + << " times in " << pix.GetName() << " (w/o saturation!) " << endl; + } + + overflow = hist->GetBinContent(0); + if (overflow > 0.0005*hist->GetEntries()) + { + *fLog << warn << "Hist-underflow " << overflow + << " times in " << pix.GetName() << " (w/o saturation!) " << endl; + } +} + Index: trunk/MagicSoft/Mars/mhcalib/MHCalibrationHiLoCam.h =================================================================== --- trunk/MagicSoft/Mars/mhcalib/MHCalibrationHiLoCam.h (revision 5749) +++ trunk/MagicSoft/Mars/mhcalib/MHCalibrationHiLoCam.h (revision 5749) @@ -0,0 +1,76 @@ +#ifndef MARS_MHCalibrationHiLoCam +#define MARS_MHCalibrationHiLoCam + +#ifndef MARS_MHCalibrationCam +#include "MHCalibrationCam.h" +#endif + +#ifndef MARS_MArrayI +#include "MArrayI.h" +#endif + +#ifndef MARS_MArrayD +#include "MArrayD.h" +#endif + +class TH1F; + +class MGeomCam; +class MRawEvtData; +class MHCalibrationHiLoCam : public MHCalibrationCam +{ + +private: + + static const Int_t fgNbins; //! Default number of bins (now set to: 900 ) + static const Axis_t fgFirst; //! Default lower histogram limit (now set to: -13.5 ) + static const Axis_t fgLast; //! Default upper histogram limit (now set to: 13.5 ) + static const Float_t fgProbLimit; //! The default for fProbLimit (now set to: 0.0) + + static const TString gsHistName; //! Default Histogram names + static const TString gsHistTitle; //! Default Histogram titles + static const TString gsHistXTitle; //! Default Histogram x-axis titles + static const TString gsHistYTitle; //! Default Histogram y-axis titles + + static const Byte_t fgLowerLim; //! Default for fLowerLimt (now set to: 200) + static const Byte_t fgUpperLim; //! Default for fUpperLimt (now set to: 252) + + Byte_t fLowerLim; // Lower allowed FADC value limit for max. slicey + Byte_t fUpperLim; // Upper allowed FADC value limit for max. slicey + + MArrayD fSumarea ; // + MArrayD fSumsector; // + MArrayI fNumarea ; // + MArrayI fNumsector; // + + MRawEvtData *fRawEvt; //! Raw event data + + Bool_t SetupHists(const MParList *pList); + Bool_t ReInitHists(MParList *pList); + Bool_t FillHists(const MParContainer *par, const Stat_t w=1); + Bool_t FinalizeHists(); + void FinalizeBadPixels(); + + void CheckOverflow( MHCalibrationPix &pix); + +public: + + MHCalibrationHiLoCam(const char *name=NULL, const char *title=NULL); + ~MHCalibrationHiLoCam() {} + + // Clone + TObject *Clone(const char *name="") const; + + Byte_t GetLowerLim() const { return fLowerLim; } + Byte_t GetUpperLim() const { return fUpperLim; } + + void SetLowerLim( const UInt_t i=fgLowerLim) { fLowerLim = i; } + void SetUpperLim( const UInt_t i=fgUpperLim) { fUpperLim = i; } + + Bool_t GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type=0) const; + void DrawPixelContent(Int_t idx) const; + + ClassDef(MHCalibrationHiLoCam, 1) // Histogram class for High Gain vs. Low Gain Calibration +}; + +#endif Index: trunk/MagicSoft/Mars/mhcalib/Makefile =================================================================== --- trunk/MagicSoft/Mars/mhcalib/Makefile (revision 5748) +++ trunk/MagicSoft/Mars/mhcalib/Makefile (revision 5749) @@ -40,4 +40,5 @@ MHCalibrationChargePINDiode.cc \ MHCalibrationRelTimeCam.cc \ + MHCalibrationHiLoCam.cc \ MHCalibrationTestCam.cc \ MHPedestalCam.cc \