Index: trunk/MagicSoft/Mars/Changelog =================================================================== --- trunk/MagicSoft/Mars/Changelog (revision 2851) +++ trunk/MagicSoft/Mars/Changelog (revision 2852) @@ -4,4 +4,33 @@ -*-*- END OF LINE -*-*- + + 2004/01/19: Markus Gaug + + * mcalib/MHCalibrationBlindPixel.[h,cc] + - fixed a bug in fFitLegend due to which program crashed by calls + to DrawClone + - Modified way to change the fit function + + * mcalib/MHCalibrationPixel.[h,cc] + - reordered function calls + - removed SetupFill + + * mcalib/MHCalibrationPINDiode.h + - reordered function calls + + * mcalib/MCalibrationPix.[h,cc] + - limits to define fFitValid now as variables in class + + * mcalib/MCalibrationCam.[h,cc] + - reordered function calls + - incorporate option to exclude pixels + + * mcalib/MCalibrationBlindPix.h + - Modified way to change the fit function + + * mcalib/MCalibrationCalc.[h,cc] + - Modified way to change the fit function + - incorporate option to exclude pixels from configuration file + 2004/01/19: Javier Rico Index: trunk/MagicSoft/Mars/mcalib/MCalibrationBlindPix.h =================================================================== --- trunk/MagicSoft/Mars/mcalib/MCalibrationBlindPix.h (revision 2851) +++ trunk/MagicSoft/Mars/mcalib/MCalibrationBlindPix.h (revision 2852) @@ -37,4 +37,5 @@ void Clear(Option_t *o=""); + // Getters Float_t GetLambda() const { return fLambda; } Float_t GetMu0() const { return fMu0; } @@ -52,17 +53,22 @@ Float_t GetErrTime() const { return fErrTime; } + MHCalibrationBlindPixel *GetHist() const { return fHist; } + + Bool_t IsFitOK() { return fHist->IsFitOK(); } + + // Fill histos Bool_t FillCharge(Float_t q) { return fHist->FillBlindPixelCharge(q); } Bool_t FillTime(Int_t t) { return fHist->FillBlindPixelTime(t); } Bool_t FillRChargevsTime(Float_t rq, Int_t t) { return fHist->FillBlindPixelChargevsN(rq,t); } - Bool_t IsFitOK() { return fHist->IsFitOK(); } - + // Fits Bool_t FitCharge(); Bool_t FitTime(); + void ChangeFitFunc(MHCalibrationBlindPixel::FitFunc_t f) { fHist->ChangeFitFunc(f); } + // Draws void Draw(Option_t *opt="") { fHist->Draw(opt); } TObject *DrawClone(Option_t *opt="") const { return fHist->DrawClone(opt); } - MHCalibrationBlindPixel *GetHist() const { return fHist; } ClassDef(MCalibrationBlindPix, 1) // Storage Container for Calibration information of one pixel Index: trunk/MagicSoft/Mars/mcalib/MCalibrationCalc.cc =================================================================== --- trunk/MagicSoft/Mars/mcalib/MCalibrationCalc.cc (revision 2851) +++ trunk/MagicSoft/Mars/mcalib/MCalibrationCalc.cc (revision 2852) @@ -98,4 +98,7 @@ #include "MTime.h" #include "TMath.h" +#include "TSystem.h" + +#include ClassImp(MCalibrationCalc); @@ -111,4 +114,5 @@ fEvents(0), fHistOverFlow(0), fCosmics(0), fNumHiGainSamples(0), fNumLoGainSamples(0), fConversionHiLo(0.), + fNumExcludedPixels(0), fColor(kEBlue) { @@ -231,26 +235,87 @@ fNumHiGainSamples = fSignals->GetNumUsedHiGainFADCSlices(); fNumLoGainSamples = fSignals->GetNumUsedLoGainFADCSlices(); - fSqrtHiGainSamples = TMath::Sqrt((Float_t)fNumHiGainSamples); - fCalibrations->InitSize(cam->GetNumPixels()); - - for (UInt_t i=0;iGetNumPixels();i++) - { - + UInt_t npixels = cam->GetNumPixels(); + + fCalibrations->InitSize(npixels); + + for (UInt_t i=0; iSetTimeFitRangesHiGain(fSignals->GetFirstUsedSliceHiGain(), - fSignals->GetLastUsedSliceHiGain()); + fSignals->GetLastUsedSliceHiGain()); hist->SetTimeFitRangesLoGain(fSignals->GetFirstUsedSliceLoGain(), - fSignals->GetLastUsedSliceLoGain()); - + fSignals->GetLastUsedSliceLoGain()); + + } + + // + // Look for file to exclude pixels from analysis + // + if (!fExcludedPixelsFile.IsNull()) + { + + *fLog << "HALLO: " << fExcludedPixelsFile.Data() << endl; + + fExcludedPixelsFile = gSystem->ExpandPathName(fExcludedPixelsFile.Data()); + + // + // Initialize reading the file + // + ifstream in(fExcludedPixelsFile.Data(),ios::in); + + if (in) + { + *fLog << inf << "Use excluded pixels from file: '" << fExcludedPixelsFile.Data() << "'" << endl; + // + // Read the file and count the number of entries + // + UInt_t pixel = 0; + + while (++fNumExcludedPixels) + { + + in >> pixel; + + if (!in.good()) + break; + // + // Check for out of range + // + if (pixel > npixels) + { + *fLog << warn << "WARNING: To be excluded pixel: " << pixel + << " is out of range " << endl; + continue; + } + // + // Exclude pixel + // + MCalibrationPix &pix = (*fCalibrations)[pixel]; + pix.SetExcluded(); + + *fLog << GetDescriptor() << inf << ": Exclude Pixel: " << pixel << endl; + + } + + if (--fNumExcludedPixels == 0) + *fLog << warn << "WARNING: File '" << fExcludedPixelsFile.Data() + << "'" << " is empty " << endl; + else + fCalibrations->SetNumPixelsExcluded(fNumExcludedPixels); + + } + else + *fLog << warn << dbginf << "Cannot open file '" << fExcludedPixelsFile.Data() << "'" << endl; } return kTRUE; - } + // -------------------------------------------------------------------------- Index: trunk/MagicSoft/Mars/mcalib/MCalibrationCalc.h =================================================================== --- trunk/MagicSoft/Mars/mcalib/MCalibrationCalc.h (revision 2851) +++ trunk/MagicSoft/Mars/mcalib/MCalibrationCalc.h (revision 2852) @@ -14,4 +14,14 @@ #include "MTask.h" #endif + +#ifndef ROOT_TArrayI +#include "TArrayI.h" +#endif + +#ifndef MARS_MCalibrationCam +#include "MCalibrationCam.h" +#endif + +#include "TString.h" class MRawEvtData; @@ -45,10 +55,10 @@ Float_t fSqrtHiGainSamples; + Float_t fConversionHiLo; Byte_t fFlags; // Flag for the fits used - - Float_t fConversionHiLo; + + TString fExcludedPixelsFile; + UInt_t fNumExcludedPixels; - enum { kUseTimeFits, kUseBlindPixelFit, kUsePinDiodeFit }; - public: @@ -68,14 +78,29 @@ MCalibrationCalc(const char *name=NULL, const char *title=NULL); - void SetSkipTimeFits(Bool_t b=kTRUE) + +private: + + enum { kUseTimeFits, kUseBlindPixelFit, kUsePinDiodeFit }; + +public: + + // Skipping fits + void SkipTimeFits(Bool_t b=kTRUE) {b ? CLRBIT(fFlags, kUseTimeFits) : SETBIT(fFlags, kUseTimeFits);} - void SetSkipBlindPixelFit(Bool_t b=kTRUE) + void SkipBlindPixelFit(Bool_t b=kTRUE) {b ? CLRBIT(fFlags, kUseBlindPixelFit) : SETBIT(fFlags, kUseBlindPixelFit);} - void SetSkipPinDiodeFit(Bool_t b=kTRUE) + void SkipPinDiodeFit(Bool_t b=kTRUE) {b ? CLRBIT(fFlags, kUsePinDiodeFit) : SETBIT(fFlags, kUsePinDiodeFit);} + // Setters void SetPulserColor(PulserColor_t color) { fColor = color; } + void SetConversionHiLo(Float_t conv) { fConversionHiLo = conv; } - void SetConversionHiLo(Float_t conv) { fConversionHiLo = conv; } + // Getters + MCalibrationBlindPix *GetBlindPixel() const { return fCalibrations->GetBlindPixel(); } + MCalibrationPINDiode *GetPINDiode() const { return fCalibrations->GetPINDiode(); } + + // Exclude pixels from configuration file + void ExcludePixelsFromAsciiFile(const char *file) { fExcludedPixelsFile = file; } ClassDef(MCalibrationCalc, 1) // Task to fill the Calibration Containers from raw data Index: trunk/MagicSoft/Mars/mcalib/MCalibrationCam.cc =================================================================== --- trunk/MagicSoft/Mars/mcalib/MCalibrationCam.cc (revision 2851) +++ trunk/MagicSoft/Mars/mcalib/MCalibrationCam.cc (revision 2852) @@ -75,5 +75,7 @@ fOffsets(NULL), fSlopes(NULL), - fOffvsSlope(NULL) + fOffvsSlope(NULL), + fNumExcludedPixels(0) + { fName = name ? name : "MCalibrationCam"; @@ -256,5 +258,5 @@ { - if (pix->IsFitValid()) + if (pix->IsFitValid() && !pix->IsExcluded()) { @@ -282,5 +284,5 @@ { - if (!pix->IsFitValid()) + if (!pix->IsFitValid() && !pix->IsExcluded()) { @@ -297,4 +299,10 @@ *fLog << all << id << " pixels with errors :-((" << endl; + *fLog << all << endl; + *fLog << all << "Excluded pixels:" << endl; + *fLog << all << endl; + + *fLog << all << fNumExcludedPixels << " excluded pixels " << endl; + } Index: trunk/MagicSoft/Mars/mcalib/MCalibrationCam.h =================================================================== --- trunk/MagicSoft/Mars/mcalib/MCalibrationCam.h (revision 2851) +++ trunk/MagicSoft/Mars/mcalib/MCalibrationCam.h (revision 2852) @@ -62,4 +62,6 @@ TH2D* fOffvsSlope; //! + UInt_t fNumExcludedPixels; + public: @@ -78,27 +80,15 @@ void InitSize(const Int_t i); + + // Setters + void SetColor(const CalibrationColor_t color) { fColor = color; } + void SetNumPixelsExcluded(const UInt_t n) { fNumExcludedPixels = n; } + + // Getters Int_t GetSize() const; - UInt_t GetNumPixels() const { return fNumPixels; } - Bool_t IsPixelUsed(Int_t idx) const; - Bool_t IsPixelFitted(Int_t idx) const; - - MCalibrationPix &operator[](Int_t i); - MCalibrationPix &operator[](Int_t i) const; - - Bool_t CheckBounds(Int_t i) const; - - void Print(Option_t *o="") const; - - void CutEdges(); - - Bool_t GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type=0) const; - void DrawPixelContent(Int_t num) const; - MCalibrationBlindPix *GetBlindPixel() const { return fBlindPixel; } MCalibrationPINDiode *GetPINDiode() const { return fPINDiode; } - - void SetColor(CalibrationColor_t color) { fColor = color; } Bool_t GetConversionFactorFFactor(Int_t ipx, Float_t &mean, Float_t &err, Float_t &sigma); @@ -107,8 +97,26 @@ Bool_t GetConversionFactorCombined(Int_t ipx, Float_t &mean, Float_t &err, Float_t &sigma); + Bool_t IsPixelUsed(Int_t idx) const; + Bool_t IsPixelFitted(Int_t idx) const; + + // Others + MCalibrationPix &operator[](Int_t i); + MCalibrationPix &operator[](Int_t i) const; + + void CutEdges(); + Bool_t CheckBounds(Int_t i) const; + + // Prints + void Print(Option_t *o="") const; + + // Draws + void DrawPixelContent(Int_t num) const; + void DrawHiLoFits(); + + // Others + Bool_t GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type=0) const; + Bool_t CalcNumPhotInsidePlexiglass(); Bool_t CalcNumPhotOutsidePlexiglass(); - - void DrawHiLoFits(); ClassDef(MCalibrationCam, 1) // Storage Container for all calibration information of the camera Index: trunk/MagicSoft/Mars/mcalib/MCalibrationPix.cc =================================================================== --- trunk/MagicSoft/Mars/mcalib/MCalibrationPix.cc (revision 2851) +++ trunk/MagicSoft/Mars/mcalib/MCalibrationPix.cc (revision 2852) @@ -89,5 +89,7 @@ fBlindPixelMethodValid(kFALSE), fFFactorMethodValid(kFALSE), - fPINDiodeMethodValid(kFALSE) + fPINDiodeMethodValid(kFALSE), + fChargeLimit(5.), + fChargeErrLimit(0.) { @@ -375,13 +377,15 @@ equivpedestal /= fConversionHiLo; - if (fCharge < 5.*equivpedestal) - { - *fLog << warn << "WARNING: Fitted Charge is smaller than 5 Pedestal RMS in Pixel " << fPixId << endl; - return kFALSE; - } - - if (fErrCharge < 0.) - { - *fLog << warn << "WARNING: Error of Fitted Charge is smaller than 0 in Pixel " << fPixId << endl; + if (fCharge < fChargeLimit*equivpedestal) + { + *fLog << warn << "WARNING: Fitted Charge is smaller than " + << fChargeLimit << " Pedestal RMS in Pixel " << fPixId << endl; + return kFALSE; + } + + if (fErrCharge < fChargeErrLimit) + { + *fLog << warn << "WARNING: Error of Fitted Charge is smaller than " + << fChargeErrLimit << " in Pixel " << fPixId << endl; return kFALSE; } @@ -511,5 +515,5 @@ { *fLog << warn << "WARNING: Could not fit Lo Gain times of pixel " << fPixId << endl; - fHist->PrintTimeFitResult(); + // fHist->PrintTimeFitResult(); return kFALSE; } @@ -524,5 +528,5 @@ { *fLog << warn << "WARNING: Could not fit Hi Gain times of pixel " << fPixId << endl; - fHist->PrintTimeFitResult(); + // fHist->PrintTimeFitResult(); return kFALSE; } Index: trunk/MagicSoft/Mars/mcalib/MCalibrationPix.h =================================================================== --- trunk/MagicSoft/Mars/mcalib/MCalibrationPix.h (revision 2851) +++ trunk/MagicSoft/Mars/mcalib/MCalibrationPix.h (revision 2852) @@ -63,4 +63,7 @@ MHCalibrationPixel *fHist; //! Pointer to the histograms performing the fits, etc. + Float_t fChargeLimit; + Float_t fChargeErrLimit; + Bool_t CheckChargeFitValidity(); Bool_t CheckTimeFitValidity(); Index: trunk/MagicSoft/Mars/mcalib/MHCalibrationBlindPixel.cc =================================================================== --- trunk/MagicSoft/Mars/mcalib/MHCalibrationBlindPixel.cc (revision 2851) +++ trunk/MagicSoft/Mars/mcalib/MHCalibrationBlindPixel.cc (revision 2852) @@ -51,4 +51,6 @@ #include +#include + #include "MBinning.h" #include "MParList.h" @@ -65,12 +67,11 @@ // MHCalibrationBlindPixel::MHCalibrationBlindPixel(const char *name, const char *title) - : fSinglePheFit(NULL), fTimeGausFit(NULL), fSinglePhePedFit(NULL), - fLambda(0.), fMu0(0.), fMu1(0.), fSigma0(0.), fSigma1(0.), - fLambdaErr(0.), fMu0Err(0.), fMu1Err(0.), fSigma0Err(0.), fSigma1Err(0.), - fChisquare(0.), fProb(0.), fNdf(0), - fMeanTime(0.), fMeanTimeErr(0.), fSigmaTime(0.), fSigmaTimeErr(0.), - fLambdaCheck(0.), fLambdaCheckErr(0.), - fgSinglePheFitFunc(&fKto4),fgSinglePheFitNPar(6) - + : fFitLegend(NULL), fFitOK(kFALSE), + fLambda(0.), fMu0(0.), fMu1(0.), fSigma0(0.), fSigma1(0.), + fLambdaErr(0.), fMu0Err(0.), fMu1Err(0.), fSigma0Err(0.), fSigma1Err(0.), + fChisquare(0.), fProb(0.), fNdf(0), + fMeanTime(0.), fMeanTimeErr(0.), fSigmaTime(0.), fSigmaTimeErr(0.), + fLambdaCheck(0.), fLambdaCheckErr(0.), + fFitFunc(kEPoisson4) { @@ -115,13 +116,11 @@ { + if (fFitLegend) + delete fFitLegend; + delete fHBlindPixelCharge; delete fHBlindPixelTime; - - if (fSinglePheFit) - delete fSinglePheFit; - if (fSinglePhePedFit) - delete fSinglePhePedFit; - if (fTimeGausFit) - delete fTimeGausFit; + delete fHBlindPixelChargevsN; + } @@ -154,42 +153,57 @@ const TString line1 = Form("Mean: #lambda = %2.2f #pm %2.2f",GetLambda(),GetLambdaErr()); - fFitLegend->AddText(line1); + TText *t1 = fFitLegend->AddText(line1.Data()); + t1->SetBit(kCanDelete); const TString line6 = Form("Mean #lambda (check) = %2.2f #pm %2.2f",GetLambdaCheck(),GetLambdaCheckErr()); - fFitLegend->AddText(line6); + TText *t2 = fFitLegend->AddText(line6.Data()); + t2->SetBit(kCanDelete); const TString line2 = Form("Pedestal: #mu_{0} = %2.2f #pm %2.2f",GetMu0(),GetMu0Err()); - fFitLegend->AddText(line2); + TText *t3 = fFitLegend->AddText(line2.Data()); + t3->SetBit(kCanDelete); const TString line3 = Form("Width Pedestal: #sigma_{0} = %2.2f #pm %2.2f",GetSigma0(),GetSigma0Err()); - fFitLegend->AddText(line3); + TText *t4 = fFitLegend->AddText(line3.Data()); + t4->SetBit(kCanDelete); const TString line4 = Form("1^{st} Phe-peak: #mu_{1} = %2.2f #pm %2.2f",GetMu1(),GetMu1Err()); - fFitLegend->AddText(line4); + TText *t5 = fFitLegend->AddText(line4.Data()); + t5->SetBit(kCanDelete); const TString line5 = Form("Width 1^{st} Phe-peak: #sigma_{1} = %2.2f #pm %2.2f",GetSigma1(),GetSigma1Err()); - fFitLegend->AddText(line5); + TText *t6 = fFitLegend->AddText(line5.Data()); + t6->SetBit(kCanDelete); const TString line7 = Form("#chi^{2} / N_{dof}: %4.2f / %3i",GetChiSquare(),GetNdf()); - fFitLegend->AddText(line7); + TText *t7 = fFitLegend->AddText(line7.Data()); + t7->SetBit(kCanDelete); const TString line8 = Form("Probability: %4.2f ",GetProb()); - fFitLegend->AddText(line8); + TText *t8 = fFitLegend->AddText(line8.Data()); + t8->SetBit(kCanDelete); if (fFitOK) - fFitLegend->AddText("Result of the Fit: OK"); - else - fFitLegend->AddText("Result of the Fit: NOT OK"); - + { + TText *t = fFitLegend->AddText(0.,0.,"Result of the Fit: OK"); + t->SetBit(kCanDelete); + } + else + { + TText *t = fFitLegend->AddText("Result of the Fit: NOT OK"); + t->SetBit(kCanDelete); + } + fFitLegend->SetBit(kCanDelete); fFitLegend->Draw(); - + + return; } @@ -201,9 +215,10 @@ TObject *newobj = Clone(); + newobj->SetBit(kCanDelete); - if (!newobj) return 0; + // ((MHCalibrationBlindPixel*)newobj)->Draw(opt); newobj->Draw(opt); return newobj; @@ -220,5 +235,5 @@ gStyle->SetOptFit(0); - gStyle->SetOptStat(1111111); + gStyle->SetOptStat(111111); TCanvas *c = MakeDefCanvas(this,550,700); @@ -234,22 +249,15 @@ fHBlindPixelCharge->DrawCopy(opt); - if (fSinglePheFit) - { - if (fFitOK) - fSinglePheFit->SetLineColor(kGreen); - else - fSinglePheFit->SetLineColor(kRed); - - fSinglePheFit->DrawCopy("same"); - c->Modified(); - c->Update(); - - if (fSinglePhePedFit) - { - fSinglePhePedFit->SetLineColor(kBlue); - fSinglePhePedFit->DrawCopy("same"); - } - } - + if (fFitOK) + fSinglePheFit.SetLineColor(kGreen); + else + fSinglePheFit.SetLineColor(kRed); + + fSinglePheFit.DrawCopy("same"); + c->Modified(); + c->Update(); + + fSinglePhePedFit.SetLineColor(kBlue); + fSinglePhePedFit.DrawCopy("same"); c->cd(2); @@ -257,10 +265,10 @@ c->Modified(); c->Update(); - + c->cd(3); gPad->SetLogy(1); gPad->SetBorderMode(0); fHBlindPixelTime->DrawCopy(opt); - + if (fHBlindPixelTime->GetFunction("GausTime")) { @@ -270,5 +278,5 @@ else tfit->SetLineColor(kGreen); - + tfit->DrawCopy("same"); c->Modified(); @@ -296,44 +304,18 @@ return kFALSE; } - - TF1 *simulateSinglePhe = new TF1("simulateSinglePhe",fgSinglePheFitFunc, - fBlindPixelChargefirst,fBlindPixelChargelast,fgSinglePheFitNPar); - - simulateSinglePhe->SetParameters(lambda,mu0,mu1,sigma0,sigma1); - simulateSinglePhe->SetParNames("#lambda","#mu_0","#mu_1","#sigma_0","#sigma_1"); - simulateSinglePhe->SetNpx(fBlindPixelChargenbins); - + + TF1 simulateSinglePhe; + + InitFit(simulateSinglePhe,fBlindPixelChargefirst,fBlindPixelChargelast); + for (Int_t i=0;i<10000; i++) - { - fHBlindPixelCharge->Fill(simulateSinglePhe->GetRandom()); - } + fHBlindPixelCharge->Fill(simulateSinglePhe.GetRandom()); return kTRUE; } - -void MHCalibrationBlindPixel::ChangeFitFunc(BlindPixelFitFunc fitfunc, Int_t par) -{ - - fgSinglePheFitFunc = fitfunc; - fgSinglePheFitNPar = par; - -} - - - -Bool_t MHCalibrationBlindPixel::FitSinglePhe(Axis_t rmin, Axis_t rmax, Option_t *opt) -{ - - if (fSinglePheFit) - return kFALSE; - - - // - // Get the fitting ranges - // - rmin = (rmin != 0.) ? rmin : fBlindPixelChargefirst; - rmax = (rmax != 0.) ? rmax : fBlindPixelChargelast; - +void MHCalibrationBlindPixel::InitFit(TF1& f, Axis_t min, Axis_t max) +{ + // // First guesses for the fit (should be as close to reality as possible, @@ -348,53 +330,120 @@ const Double_t norm = entries/gkSq2Pi; - fSinglePheFit = new TF1("SinglePheFit",fgSinglePheFitFunc,rmin,rmax,fgSinglePheFitNPar); - // fSinglePheFit = new TF1("SinglePheFit",fgSinglePheFitFunc,rmin,rmax,fgSinglePheFitNPar+1); - // fSinglePheFit->SetParameters(lambda_guess,mu_0_guess,mu_1_guess,si_0_guess,si_1_guess); - fSinglePheFit->SetParameters(lambda_guess,mu_0_guess,mu_1_guess,si_0_guess,si_1_guess,norm); - // fSinglePheFit->SetParNames("#lambda","#mu_0","#mu_1","#sigma_0","#sigma_1"); - fSinglePheFit->SetParNames("#lambda","#mu_0","#mu_1","#sigma_0","#sigma_1","Area"); - fSinglePheFit->SetParLimits(0,0.,1.); - fSinglePheFit->SetParLimits(1,rmin,(rmax-rmin)/1.5); - fSinglePheFit->SetParLimits(2,(rmax-rmin)/2.,(rmax-0.05*(rmax-rmin))); - fSinglePheFit->SetParLimits(3,1.0,(rmax-rmin)/2.0); - fSinglePheFit->SetParLimits(4,1.0,(rmax-rmin)/2.5); - fSinglePheFit->SetParLimits(5,norm-0.1,norm+0.1); - // fSinglePheFit->SetParLimits(5,entries/gkSq2Pi,entries/gkSq2Pi); - // fSinglePheFit->SetParLimits(5,0.,1.5*entries); - // - // Normalize the histogram to facilitate faster fitting of the area - // For speed reasons, fKto4 is normalized to Sqrt(2 pi). - // Therefore also normalize the histogram to that value - // - // fHBlindPixelCharge->Scale(gkSq2Pi*(float)bins/npx/entries); - // fHBlindPixelCharge->Scale(gkSq2Pi/entries); - // norm *= (Float_t)fSinglePheFit->GetNpx()/(Float_t)fBlindPixelChargenbins; - - fHBlindPixelCharge->Fit(fSinglePheFit,opt); - - fLambda = fSinglePheFit->GetParameter(0); - fMu0 = fSinglePheFit->GetParameter(1); - fMu1 = fSinglePheFit->GetParameter(2); - fSigma0 = fSinglePheFit->GetParameter(3); - fSigma1 = fSinglePheFit->GetParameter(4); - - fLambdaErr = fSinglePheFit->GetParError(0); - fMu0Err = fSinglePheFit->GetParError(1); - fMu1Err = fSinglePheFit->GetParError(2); - fSigma0Err = fSinglePheFit->GetParError(3); - fSigma1Err = fSinglePheFit->GetParError(4); - - fProb = fSinglePheFit->GetProb(); - fChisquare = fSinglePheFit->GetChisquare(); - fNdf = fSinglePheFit->GetNDF(); + // + // Initialize + // + switch (fFitFunc) + { + + case kEPoisson4: + f = TF1("SinglePheFit",&fPoissonKto4,min,max,6); + f.SetParameters(lambda_guess,mu_0_guess,mu_1_guess,si_0_guess,si_1_guess,norm); + f.SetParNames("#lambda","#mu_0","#mu_1","#sigma_0","#sigma_1","Area"); + f.SetParLimits(0,0.,1.); + f.SetParLimits(1,min,(max-min)/1.5); + f.SetParLimits(2,(max-min)/2.,(max-0.05*(max-min))); + f.SetParLimits(3,1.0,(max-min)/2.0); + f.SetParLimits(4,1.0,(max-min)/2.5); + f.SetParLimits(5,norm-0.1,norm+0.1); + break; + + case kEPoisson5: + f = TF1("SinglePheFit",&fPoissonKto5,min,max,6); + f.SetParameters(lambda_guess,mu_0_guess,mu_1_guess,si_0_guess,si_1_guess,norm); + f.SetParNames("#lambda","#mu_0","#mu_1","#sigma_0","#sigma_1","Area"); + f.SetParLimits(0,0.,1.); + f.SetParLimits(1,min,(max-min)/1.5); + f.SetParLimits(2,(max-min)/2.,(max-0.05*(max-min))); + f.SetParLimits(3,1.0,(max-min)/2.0); + f.SetParLimits(4,1.0,(max-min)/2.5); + f.SetParLimits(5,norm-0.1,norm+0.1); + break; + + case kEPoisson6: + f = TF1("SinglePheFit",&fPoissonKto6,min,max,6); + f.SetParameters(lambda_guess,mu_0_guess,mu_1_guess,si_0_guess,si_1_guess,norm); + f.SetParNames("#lambda","#mu_0","#mu_1","#sigma_0","#sigma_1","Area"); + f.SetParLimits(0,0.,1.); + f.SetParLimits(1,min,(max-min)/1.5); + f.SetParLimits(2,(max-min)/2.,(max-0.05*(max-min))); + f.SetParLimits(3,1.0,(max-min)/2.0); + f.SetParLimits(4,1.0,(max-min)/2.5); + f.SetParLimits(5,norm-0.1,norm+0.1); + break; + + case kEPoisson7: + break; + + case kEPolya: + break; + + case kEMichele: + break; + + } +} + +void MHCalibrationBlindPixel::ExitFit(TF1 &f) +{ + + + // + // Finalize + // + switch (fFitFunc) + { + + case kEPoisson4: + case kEPoisson5: + case kEPoisson6: + case kEPoisson7: + fLambda = fSinglePheFit.GetParameter(0); + fMu0 = fSinglePheFit.GetParameter(1); + fMu1 = fSinglePheFit.GetParameter(2); + fSigma0 = fSinglePheFit.GetParameter(3); + fSigma1 = fSinglePheFit.GetParameter(4); + + fLambdaErr = fSinglePheFit.GetParError(0); + fMu0Err = fSinglePheFit.GetParError(1); + fMu1Err = fSinglePheFit.GetParError(2); + fSigma0Err = fSinglePheFit.GetParError(3); + fSigma1Err = fSinglePheFit.GetParError(4); + break; + default: + break; + } + +} + + +Bool_t MHCalibrationBlindPixel::FitSinglePhe(Axis_t rmin, Axis_t rmax, Option_t *opt) +{ + + // + // Get the fitting ranges + // + rmin = (rmin != 0.) ? rmin : fBlindPixelChargefirst; + rmax = (rmax != 0.) ? rmax : fBlindPixelChargelast; + + InitFit(fSinglePheFit,rmin,rmax); + + fHBlindPixelCharge->Fit(&fSinglePheFit,opt); + + ExitFit(fSinglePheFit); + + fProb = fSinglePheFit.GetProb(); + fChisquare = fSinglePheFit.GetChisquare(); + fNdf = fSinglePheFit.GetNDF(); // Perform the cross-check fitting only the pedestal: - fSinglePhePedFit = new TF1("GausPed","gaus",rmin,fMu0); - fHBlindPixelCharge->Fit(fSinglePhePedFit,opt); - - Double_t pedarea = fSinglePhePedFit->GetParameter(0)*gkSq2Pi*fSinglePhePedFit->GetParameter(2); + fSinglePhePedFit = TF1("GausPed","gaus",rmin,fMu0); + fHBlindPixelCharge->Fit(&fSinglePhePedFit,opt); + + const Stat_t entries = fHBlindPixelCharge->Integral("width"); + + Double_t pedarea = fSinglePhePedFit.GetParameter(0)*gkSq2Pi*fSinglePhePedFit.GetParameter(2); fLambdaCheck = TMath::Log(entries/pedarea); - fLambdaCheckErr = fSinglePhePedFit->GetParError(0)/fSinglePhePedFit->GetParameter(0) - + fSinglePhePedFit->GetParError(2)/fSinglePhePedFit->GetParameter(2); + fLambdaCheckErr = fSinglePhePedFit.GetParError(0)/fSinglePhePedFit.GetParameter(0) + + fSinglePhePedFit.GetParError(2)/fSinglePhePedFit.GetParameter(2); *fLog << inf << "Results of the Blind Pixel Fit: " << endl; @@ -429,6 +478,5 @@ } else - *fLog << inf << GetDescriptor() << ": " << contSinglePhe - << " in Single Photo-Electron peak " << endl; + *fLog << inf << contSinglePhe << " in Single Photo-Electron peak " << endl; fFitOK = kTRUE; @@ -456,7 +504,4 @@ { - if (fTimeGausFit) - return kFALSE; - rmin = (rmin != 0.) ? rmin : 4.; rmax = (rmax != 0.) ? rmax : 9.; @@ -467,27 +512,26 @@ const Double_t area_guess = entries/gkSq2Pi; - fTimeGausFit = new TF1("GausTime","gaus",rmin,rmax); - fTimeGausFit->SetParameters(area_guess,mu_guess,sigma_guess); - fTimeGausFit->SetParNames("Area","#mu","#sigma"); - fTimeGausFit->SetParLimits(0,0.,entries); - fTimeGausFit->SetParLimits(1,rmin,rmax); - fTimeGausFit->SetParLimits(2,0.,rmax-rmin); - - fHBlindPixelTime->Fit(fTimeGausFit,opt); - - rmin = fTimeGausFit->GetParameter(1) - 2.*fTimeGausFit->GetParameter(2); - rmax = fTimeGausFit->GetParameter(1) + 2.*fTimeGausFit->GetParameter(2); - fTimeGausFit->SetRange(rmin,rmax); - - fHBlindPixelTime->Fit(fTimeGausFit,opt); - - fMeanTime = fTimeGausFit->GetParameter(2); - fSigmaTime = fTimeGausFit->GetParameter(3); - fMeanTimeErr = fTimeGausFit->GetParError(2); - fSigmaTimeErr = fTimeGausFit->GetParError(3); + fTimeGausFit = TF1("GausTime","gaus",rmin,rmax); + fTimeGausFit.SetParameters(area_guess,mu_guess,sigma_guess); + fTimeGausFit.SetParNames("Area","#mu","#sigma"); + fTimeGausFit.SetParLimits(0,0.,entries); + fTimeGausFit.SetParLimits(1,rmin,rmax); + fTimeGausFit.SetParLimits(2,0.,rmax-rmin); + + fHBlindPixelTime->Fit(&fTimeGausFit,opt); + rmin = fTimeGausFit.GetParameter(1) - 2.*fTimeGausFit.GetParameter(2); + rmax = fTimeGausFit.GetParameter(1) + 2.*fTimeGausFit.GetParameter(2); + fTimeGausFit.SetRange(rmin,rmax); + + fHBlindPixelTime->Fit(&fTimeGausFit,opt); + + fMeanTime = fTimeGausFit.GetParameter(2); + fSigmaTime = fTimeGausFit.GetParameter(3); + fMeanTimeErr = fTimeGausFit.GetParError(2); + fSigmaTimeErr = fTimeGausFit.GetParError(3); *fLog << inf << "Results of the Times Fit: " << endl; - *fLog << inf << "Chisquare: " << fTimeGausFit->GetChisquare() << endl; - *fLog << inf << "Ndf: " << fTimeGausFit->GetNDF() << endl; + *fLog << inf << "Chisquare: " << fTimeGausFit.GetChisquare() << endl; + *fLog << inf << "Ndf: " << fTimeGausFit.GetNDF() << endl; return kTRUE; Index: trunk/MagicSoft/Mars/mcalib/MHCalibrationBlindPixel.h =================================================================== --- trunk/MagicSoft/Mars/mcalib/MHCalibrationBlindPixel.h (revision 2851) +++ trunk/MagicSoft/Mars/mcalib/MHCalibrationBlindPixel.h (revision 2852) @@ -38,7 +38,7 @@ TH1I* fHBlindPixelChargevsN; //-> Summed Charge vs. Event Nr. - TF1 *fSinglePheFit; - TF1 *fTimeGausFit; - TF1 *fSinglePhePedFit; + TF1 fSinglePheFit; + TF1 fTimeGausFit; + TF1 fSinglePhePedFit; Axis_t fBlindPixelChargefirst; @@ -49,5 +49,5 @@ void DrawLegend(); - TPaveText *fFitLegend; + TPaveText *fFitLegend; //! Bool_t fFitOK; @@ -81,12 +81,12 @@ ~MHCalibrationBlindPixel(); - typedef Double_t (*BlindPixelFitFunc)(Double_t *, Double_t *); - Bool_t FillBlindPixelCharge(Float_t q) { return fHBlindPixelCharge->Fill(q) > -1; } Bool_t FillBlindPixelTime(Int_t t) { return fHBlindPixelTime->Fill(t) > -1; } Bool_t FillBlindPixelChargevsN(Stat_t rq, Int_t t) { return fHBlindPixelChargevsN->Fill(t,rq) > -1; } + + //Getters const Double_t GetLambda() const { return fLambda; } - const Double_t GetLambdaCheck() const { return fLambdaCheck; } + const Double_t GetLambdaCheck() const { return fLambdaCheck; } const Double_t GetMu0() const { return fMu0; } const Double_t GetMu1() const { return fMu1; } @@ -95,5 +95,5 @@ const Double_t GetLambdaErr() const { return fLambdaErr; } - const Double_t GetLambdaCheckErr() const { return fLambdaCheckErr; } + const Double_t GetLambdaCheckErr() const { return fLambdaCheckErr; } const Double_t GetMu0Err() const { return fMu0Err; } const Double_t GetMu1Err() const { return fMu1Err; } @@ -102,6 +102,6 @@ const Double_t GetChiSquare() const { return fChisquare; } - const Double_t GetProb() const { return fProb; } - const Int_t GetNdf() const { return fNdf; } + const Double_t GetProb() const { return fProb; } + const Int_t GetNdf() const { return fNdf; } const Double_t GetMeanTime() const { return fMeanTime; } @@ -110,28 +110,35 @@ const Double_t GetSigmaTimeErr() const { return fSigmaTimeErr; } - Bool_t SimulateSinglePhe(Double_t lambda, - Double_t mu0, - Double_t mu1, - Double_t sigma0, - Double_t sigma1); - + const Bool_t IsFitOK() { return fFitOK; } + + // Draws + TObject *DrawClone(Option_t *option="") const; + void Draw(Option_t *option=""); + + // Fits + enum FitFunc_t { kEPoisson4, kEPoisson5, kEPoisson6, kEPoisson7, kEPolya, kEMichele }; + +private: + FitFunc_t fFitFunc; + +public: Bool_t FitSinglePhe(Axis_t rmin=0, Axis_t rmax=0, Option_t *opt="RL0+Q"); Bool_t FitTime(Axis_t rmin=0., Axis_t rmax=0.,Option_t *opt="R0+Q"); - - void ChangeFitFunc(BlindPixelFitFunc fitfunc, Int_t par=6); - + void ChangeFitFunc(FitFunc_t func) { fFitFunc = func; } + + // Simulation + Bool_t SimulateSinglePhe(Double_t lambda, + Double_t mu0,Double_t mu1, + Double_t sigma0,Double_t sigma1); + + // Others void CutAllEdges(); - TObject *DrawClone(Option_t *option="") const; - void Draw(Option_t *option=""); - - Bool_t IsFitOK() { return fFitOK; } - private: - - BlindPixelFitFunc fgSinglePheFitFunc; //! In the beginning, - Int_t fgSinglePheFitNPar; //! we want to be flexible using different functions - - inline static Double_t fAna(Double_t *x, Double_t *par) + + void InitFit(TF1& f, Axis_t min, Axis_t max); + void ExitFit(TF1& f); + + inline static Double_t fFitFuncMichele(Double_t *x, Double_t *par) { @@ -189,5 +196,5 @@ } - inline static Double_t fKto4(Double_t *x, Double_t *par) + inline static Double_t fPoissonKto4(Double_t *x, Double_t *par) { @@ -246,5 +253,5 @@ - inline static Double_t fKto5(Double_t *x, Double_t *par) + inline static Double_t fPoissonKto5(Double_t *x, Double_t *par) { @@ -311,5 +318,5 @@ - inline static Double_t fKto6(Double_t *x, Double_t *par) + inline static Double_t fPoissonKto6(Double_t *x, Double_t *par) { @@ -382,5 +389,5 @@ } - ClassDef(MHCalibrationBlindPixel, 1) + ClassDef(MHCalibrationBlindPixel, 1) // Histograms from the Calibration Blind Pixel }; Index: trunk/MagicSoft/Mars/mcalib/MHCalibrationPINDiode.h =================================================================== --- trunk/MagicSoft/Mars/mcalib/MHCalibrationPINDiode.h (revision 2851) +++ trunk/MagicSoft/Mars/mcalib/MHCalibrationPINDiode.h (revision 2852) @@ -26,6 +26,6 @@ private: - TH1I* fHPCharge; //-> Histogram containing the summed 32 PINDiode slices - TH1F* fHErrCharge; //-> Variance of summed FADC slices + TH1I* fHPCharge; //-> Histogram containing the summed 32 PINDiode slices + TH1F* fHErrCharge; //-> Variance of summed FADC slices TH1I* fHPTime; //-> Histogram with time evolution of summed charges @@ -45,5 +45,5 @@ const Double_t GetErrTime() const { return fVarGausFit->GetParameter(3); } - ClassDef(MHCalibrationPINDiode, 0) + ClassDef(MHCalibrationPINDiode, 0) // Histograms from the Calibration PIN Diode }; Index: trunk/MagicSoft/Mars/mcalib/MHCalibrationPixel.cc =================================================================== --- trunk/MagicSoft/Mars/mcalib/MHCalibrationPixel.cc (revision 2851) +++ trunk/MagicSoft/Mars/mcalib/MHCalibrationPixel.cc (revision 2852) @@ -278,20 +278,4 @@ -// ------------------------------------------------------------------------- -// -// Set the binnings and prepare the filling of the histograms -// -Bool_t MHCalibrationPixel::SetupFill(const MParList *plist) -{ - - fHChargeHiGain->Reset(); - fHTimeHiGain->Reset(); - fHChargeLoGain->Reset(); - fHTimeLoGain->Reset(); - - return kTRUE; -} - - Bool_t MHCalibrationPixel::UseLoGain() { @@ -627,9 +611,6 @@ if (fTimeChisquare > 20.) // Cannot use Probability because Ndf is sometimes < 1 - { - *fLog << warn << "WARNING: Fit of the Arrival times failed ! " << endl; - return kFALSE; - } - + return kFALSE; + return kTRUE; Index: trunk/MagicSoft/Mars/mcalib/MHCalibrationPixel.h =================================================================== --- trunk/MagicSoft/Mars/mcalib/MHCalibrationPixel.h (revision 2851) +++ trunk/MagicSoft/Mars/mcalib/MHCalibrationPixel.h (revision 2852) @@ -103,19 +103,7 @@ void ChangeHistId(Int_t i); - Bool_t SetupFill(const MParList *pList); - Bool_t Fill(const MParContainer *, const Stat_t w=1) { return kTRUE; } - + // Setters void SetPointInGraph(Float_t qhi, Float_t qlo); - - Bool_t FillChargeLoGain(Float_t q) { return (fHChargeLoGain->Fill(q) > -1); } - Bool_t FillTimeLoGain(Int_t t) { return (fHTimeLoGain->Fill(t) > -1); } - Bool_t FillChargevsNLoGain(Float_t q, Int_t n) { return (fHChargevsNLoGain->Fill(n,q) > -1); } - - Bool_t FillChargeHiGain(Float_t q) { return (fHChargeHiGain->Fill(q) > -1); } - Bool_t FillTimeHiGain(Int_t t) { return (fHTimeHiGain->Fill(t) > -1); } - Bool_t FillChargevsNHiGain(Float_t q, Int_t n) { return (fHChargevsNHiGain->Fill(n,q) > -1); } - - void SetUseLoGain() { fUseLoGain = kTRUE; } - Bool_t UseLoGain(); + void SetUseLoGain(Bool_t b = kTRUE) { fUseLoGain = b; } void SetTimeFitRangesHiGain(Byte_t low, Byte_t up) { fTimeLowerFitRangeHiGain = low, @@ -123,5 +111,8 @@ void SetTimeFitRangesLoGain(Byte_t low, Byte_t up) { fTimeLowerFitRangeLoGain = low, fTimeUpperFitRangeLoGain = up ; } - + + void SetLowerFitRange(Axis_t min) { fLowerFitRange = min; } + + // Getters const TH1F *GetHCharge() { return fHChargeHiGain; } const TH1F *GetHCharge() const { return fHChargeHiGain; } @@ -158,5 +149,21 @@ Double_t GetOffset() { return fOffset; } Double_t GetSlope() { return fSlope; } + + Bool_t UseLoGain(); + + Bool_t IsFitOK() { return fFitOK; } + Bool_t IsEmpty() { return !( (fHChargeHiGain->GetEntries()) + || (fHChargeLoGain->GetEntries()) ); } + // Fill histos + Bool_t FillChargeLoGain(Float_t q) { return (fHChargeLoGain->Fill(q) > -1); } + Bool_t FillTimeLoGain(Int_t t) { return (fHTimeLoGain->Fill(t) > -1); } + Bool_t FillChargevsNLoGain(Float_t q, Int_t n) { return (fHChargevsNLoGain->Fill(n,q) > -1); } + + Bool_t FillChargeHiGain(Float_t q) { return (fHChargeHiGain->Fill(q) > -1); } + Bool_t FillTimeHiGain(Int_t t) { return (fHTimeHiGain->Fill(t) > -1); } + Bool_t FillChargevsNHiGain(Float_t q, Int_t n) { return (fHChargevsNHiGain->Fill(n,q) > -1); } + + // Fits Bool_t FitChargeHiGain(Option_t *option="RQ0"); Bool_t FitTimeHiGain(Axis_t rmin=0, Axis_t rmax=0, Option_t *option="RQ0"); @@ -167,18 +174,16 @@ void FitHiGainvsLoGain(); + // Draws virtual void Draw(Option_t *option=""); + + // Prints + void PrintChargeFitResult(); + void PrintTimeFitResult(); + + // Others virtual void CutAllEdges(); virtual void Reset(); - void SetLowerFitRange(Axis_t min) { fLowerFitRange = min; } - - void PrintChargeFitResult(); - void PrintTimeFitResult(); - - Bool_t IsFitOK() { return fFitOK; } - Bool_t IsEmpty() { return !( (fHChargeHiGain->GetEntries()) - || (fHChargeLoGain->GetEntries()) ); } - - ClassDef(MHCalibrationPixel, 1) + ClassDef(MHCalibrationPixel, 1) // Histograms for each calibrated pixel };