Index: trunk/MagicSoft/Mars/Changelog =================================================================== --- trunk/MagicSoft/Mars/Changelog (revision 2879) +++ trunk/MagicSoft/Mars/Changelog (revision 2880) @@ -20,4 +20,5 @@ plexiglass is not intended to be calculated any more (previous bug) + - exclude excluded pixels in eventloop Index: trunk/MagicSoft/Mars/mcalib/MCalibrationCalc.cc =================================================================== --- trunk/MagicSoft/Mars/mcalib/MCalibrationCalc.cc (revision 2879) +++ trunk/MagicSoft/Mars/mcalib/MCalibrationCalc.cc (revision 2880) @@ -261,6 +261,4 @@ { - *fLog << "HALLO: " << fExcludedPixelsFile.Data() << endl; - fExcludedPixelsFile = gSystem->ExpandPathName(fExcludedPixelsFile.Data()); @@ -388,4 +386,9 @@ const UInt_t pixid = pixel.GetPixelId(); + MCalibrationPix &pix = (*fCalibrations)[pixid]; + + if (pix.IsExcluded()) + continue; + MExtractedSignalPix &sig = (*fSignals)[pixid]; @@ -394,5 +397,4 @@ Float_t mtime = sig.GetMeanArrivalTime(); - MCalibrationPix &pix = (*fCalibrations)[pixid]; switch(pixid) Index: trunk/MagicSoft/Mars/mcalib/MCalibrationCam.cc =================================================================== --- trunk/MagicSoft/Mars/mcalib/MCalibrationCam.cc (revision 2879) +++ trunk/MagicSoft/Mars/mcalib/MCalibrationCam.cc (revision 2880) @@ -303,4 +303,9 @@ *fLog << all << endl; + TIter Next3(fPixels); + while ((pix=(MCalibrationPix*)Next3())) + if (pix->IsExcluded()) + *fLog << all << pix->GetPixId() << endl; + *fLog << all << fNumExcludedPixels << " excluded pixels " << endl; @@ -425,5 +430,5 @@ // -------------------------------------------------------------------------- // - +// // Bool_t MCalibrationCam::CalcNumPhotInsidePlexiglass() @@ -442,4 +447,6 @@ *TMath::Power(10,gkCalibrationBlindPixelAttGreen) // correct for absorption * gkCalibrationInnerPixelArea; // correct for area + + break; case kECBlue: Index: trunk/MagicSoft/Mars/mcalib/MHCalibrationBlindPixel.cc =================================================================== --- trunk/MagicSoft/Mars/mcalib/MHCalibrationBlindPixel.cc (revision 2879) +++ trunk/MagicSoft/Mars/mcalib/MHCalibrationBlindPixel.cc (revision 2880) @@ -67,5 +67,15 @@ // MHCalibrationBlindPixel::MHCalibrationBlindPixel(const char *name, const char *title) - : fFitLegend(NULL), fFitOK(kFALSE), + : fHBlindPixelCharge(NULL), + fHBlindPixelTime(NULL), + fHBlindPixelChargevsN(NULL), + fHBlindPixelPSD(NULL), + fSinglePheFit(NULL), + fTimeGausFit(NULL), + fSinglePhePedFit(NULL), + fBlindPixelChargefirst(0.), + fBlindPixelChargelast(0.), + fBlindPixelChargenbins(0), + fFitLegend(NULL), fFitOK(kFALSE), fLambda(0.), fMu0(0.), fMu1(0.), fSigma0(0.), fSigma1(0.), fLambdaErr(0.), fMu0Err(0.), fMu1Err(0.), fSigma0Err(0.), fSigma1Err(0.), @@ -122,5 +132,15 @@ delete fHBlindPixelTime; delete fHBlindPixelChargevsN; - + + if (fHBlindPixelPSD) + delete fHBlindPixelPSD; + + if (fSinglePheFit) + delete fSinglePheFit; + if (fTimeGausFit) + delete fTimeGausFit; + if(fSinglePhePedFit) + delete fSinglePhePedFit; + } @@ -208,22 +228,33 @@ } - -TObject *MHCalibrationBlindPixel::DrawClone(Option_t *opt) const +TObject *MHCalibrationBlindPixel::DrawClone(Option_t *option) const { gROOT->SetSelectedPad(NULL); - - TObject *newobj = Clone(); - + + // TVirtualPad *padsav = gPad; + + // TObject *newobj = (TObject *)IsA()->New(); + // Copy(*newobj); + // TObject *newobj = Clone(); + MHCalibrationBlindPixel *newobj = (MHCalibrationBlindPixel*)Clone(); + + if (!newobj) + return 0; newobj->SetBit(kCanDelete); - if (!newobj) - return 0; - - // ((MHCalibrationBlindPixel*)newobj)->Draw(opt); - newobj->Draw(opt); + + if (strlen(option)) + newobj->Draw(option); + // ((MHCalibrationBlindPixel*)newobj)->Draw(option); + // newobj->Draw(option); + else + newobj->Draw(GetDrawOption()); + // ((MHCalibrationBlindPixel*)newobj)->Draw(GetDrawOption()); + // newobj->Draw(GetDrawOption()); + // if (padsav) padsav->cd(); + return newobj; } - // ------------------------------------------------------------------------- @@ -239,25 +270,26 @@ TCanvas *c = MakeDefCanvas(this,550,700); - c->Divide(2,2); + c->Divide(2,3); gROOT->SetSelectedPad(NULL); c->cd(1); + gPad->SetLogx(0); gPad->SetLogy(1); gPad->SetTicks(); - fHBlindPixelCharge->DrawCopy(opt); + fHBlindPixelCharge->Draw(opt); if (fFitOK) - fSinglePheFit.SetLineColor(kGreen); + fSinglePheFit->SetLineColor(kGreen); else - fSinglePheFit.SetLineColor(kRed); + fSinglePheFit->SetLineColor(kRed); - fSinglePheFit.DrawCopy("same"); + fSinglePheFit->Draw("same"); c->Modified(); c->Update(); - fSinglePhePedFit.SetLineColor(kBlue); - fSinglePhePedFit.DrawCopy("same"); + fSinglePhePedFit->SetLineColor(kBlue); + fSinglePhePedFit->Draw("same"); c->cd(2); @@ -269,25 +301,22 @@ gPad->SetLogy(1); gPad->SetBorderMode(0); - fHBlindPixelTime->DrawCopy(opt); - - if (fHBlindPixelTime->GetFunction("GausTime")) - { - TF1 *tfit = fHBlindPixelTime->GetFunction("GausTime"); - if (tfit->GetProb() < 0.01) - tfit->SetLineColor(kRed); - else - tfit->SetLineColor(kGreen); - - tfit->DrawCopy("same"); - c->Modified(); - c->Update(); - } + fHBlindPixelTime->Draw(opt); c->cd(4); - fHBlindPixelChargevsN->DrawCopy(opt); + fHBlindPixelChargevsN->Draw(opt); c->Modified(); c->Update(); + + c->cd(5); + + // fHBlindPixelPSD = (TH1F*)MFFT::PowerSpectrumDensity(fHBlindPixelChargevsN); + // TH1F *hist = MFFT::PowerSpectrumDensity((*this)->fHBlindPixelChargevsN); + + // fHBlindPixelPSD->Draw(opt); + c->Modified(); + c->Update(); + } @@ -305,15 +334,14 @@ } - TF1 simulateSinglePhe; - - InitFit(simulateSinglePhe,fBlindPixelChargefirst,fBlindPixelChargelast); - + if (!InitFit(fBlindPixelChargefirst,fBlindPixelChargelast)) + return kFALSE; + for (Int_t i=0;i<10000; i++) - fHBlindPixelCharge->Fill(simulateSinglePhe.GetRandom()); + fHBlindPixelCharge->Fill(fSinglePheFit->GetRandom()); return kTRUE; } -void MHCalibrationBlindPixel::InitFit(TF1& f, Axis_t min, Axis_t max) +Bool_t MHCalibrationBlindPixel::InitFit(Axis_t min, Axis_t max) { @@ -337,37 +365,37 @@ 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); + fSinglePheFit = new TF1("SinglePheFit",&fPoissonKto4,min,max,6); + 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","Area"); + fSinglePheFit->SetParLimits(0,0.,1.); + fSinglePheFit->SetParLimits(1,min,(max-min)/1.5); + fSinglePheFit->SetParLimits(2,(max-min)/2.,(max-0.05*(max-min))); + fSinglePheFit->SetParLimits(3,1.0,(max-min)/2.0); + fSinglePheFit->SetParLimits(4,1.0,(max-min)/2.5); + fSinglePheFit->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); + fSinglePheFit = new TF1("SinglePheFit",&fPoissonKto5,min,max,6); + 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","Area"); + fSinglePheFit->SetParLimits(0,0.,1.); + fSinglePheFit->SetParLimits(1,min,(max-min)/1.5); + fSinglePheFit->SetParLimits(2,(max-min)/2.,(max-0.05*(max-min))); + fSinglePheFit->SetParLimits(3,1.0,(max-min)/2.0); + fSinglePheFit->SetParLimits(4,1.0,(max-min)/2.5); + fSinglePheFit->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); + fSinglePheFit = new TF1("SinglePheFit",&fPoissonKto6,min,max,6); + 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","Area"); + fSinglePheFit->SetParLimits(0,0.,1.); + fSinglePheFit->SetParLimits(1,min,(max-min)/1.5); + fSinglePheFit->SetParLimits(2,(max-min)/2.,(max-0.05*(max-min))); + fSinglePheFit->SetParLimits(3,1.0,(max-min)/2.0); + fSinglePheFit->SetParLimits(4,1.0,(max-min)/2.5); + fSinglePheFit->SetParLimits(5,norm-0.1,norm+0.1); break; @@ -380,9 +408,15 @@ case kEMichele: break; - + + default: + *fLog << warn << "WARNING: Could not find Fit Function for Blind Pixel " << endl; + return kFALSE; + break; } -} - -void MHCalibrationBlindPixel::ExitFit(TF1 &f) + + return kTRUE; +} + +void MHCalibrationBlindPixel::ExitFit(TF1 *f) { @@ -398,15 +432,15 @@ case kEPoisson6: case kEPoisson7: - fLambda = fSinglePheFit.GetParameter(0); - fMu0 = fSinglePheFit.GetParameter(1); - fMu1 = fSinglePheFit.GetParameter(2); - fSigma0 = fSinglePheFit.GetParameter(3); - fSigma1 = fSinglePheFit.GetParameter(4); + 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); + fLambdaErr = fSinglePheFit->GetParError(0); + fMu0Err = fSinglePheFit->GetParError(1); + fMu1Err = fSinglePheFit->GetParError(2); + fSigma0Err = fSinglePheFit->GetParError(3); + fSigma1Err = fSinglePheFit->GetParError(4); break; default: @@ -426,24 +460,26 @@ rmax = (rmax != 0.) ? rmax : fBlindPixelChargelast; - InitFit(fSinglePheFit,rmin,rmax); - - fHBlindPixelCharge->Fit(&fSinglePheFit,opt); + if (!InitFit(rmin,rmax)) + return kFALSE; + + fHBlindPixelCharge->Fit(fSinglePheFit,opt); + ExitFit(fSinglePheFit); - fProb = fSinglePheFit.GetProb(); - fChisquare = fSinglePheFit.GetChisquare(); - fNdf = fSinglePheFit.GetNDF(); + fProb = fSinglePheFit->GetProb(); + fChisquare = fSinglePheFit->GetChisquare(); + fNdf = fSinglePheFit->GetNDF(); // Perform the cross-check fitting only the pedestal: - fSinglePhePedFit = TF1("GausPed","gaus",rmin,fMu0); - fHBlindPixelCharge->Fit(&fSinglePhePedFit,opt); + fSinglePhePedFit = new 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); + 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; @@ -512,26 +548,26 @@ const Double_t area_guess = entries/gkSq2Pi; - 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); + 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); *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 2879) +++ trunk/MagicSoft/Mars/mcalib/MHCalibrationBlindPixel.h (revision 2880) @@ -10,4 +10,8 @@ #endif +#ifndef MARS_MFFT +#include "MFFT.h" +#endif + #ifndef ROOT_TH1 #include "TH1.h" @@ -25,4 +29,6 @@ #include "TPaveText.h" #endif + + const Double_t NoWay = 10000000000.0; @@ -34,14 +40,15 @@ private: - TH1F* fHBlindPixelCharge; //-> Histogram with the single Phe spectrum - TH1I* fHBlindPixelTime; //-> Variance of summed FADC slices - TH1I* fHBlindPixelChargevsN; //-> Summed Charge vs. Event Nr. - - TF1 fSinglePheFit; - TF1 fTimeGausFit; - TF1 fSinglePhePedFit; + TH1F* fHBlindPixelCharge; // Histogram with the single Phe spectrum + TH1I* fHBlindPixelTime; // Variance of summed FADC slices + TH1I* fHBlindPixelChargevsN; // Summed Charge vs. Event Nr. + TH1F* fHBlindPixelPSD; // Power spectrum density of fHBlindPixelChargevsN + + TF1 *fSinglePheFit; + TF1 *fTimeGausFit; + TF1 *fSinglePhePedFit; Axis_t fBlindPixelChargefirst; - Axis_t fBlindPixelChargelast; + Axis_t fBlindPixelChargelast; Int_t fBlindPixelChargenbins; @@ -49,5 +56,5 @@ void DrawLegend(); - TPaveText *fFitLegend; //! + TPaveText *fFitLegend; Bool_t fFitOK; @@ -73,6 +80,6 @@ Double_t fSigmaTimeErr; - Double_t fLambdaCheck; - Double_t fLambdaCheckErr; + Double_t fLambdaCheck; + Double_t fLambdaCheckErr; public: @@ -112,4 +119,7 @@ const Bool_t IsFitOK() { return fFitOK; } + const TH1I *GetHBlindPixelChargevsN() const { return fHBlindPixelChargevsN; } + const TH1F *GetHBlindPixelPSD() const { return fHBlindPixelPSD; } + // Draws TObject *DrawClone(Option_t *option="") const; @@ -137,6 +147,6 @@ private: - void InitFit(TF1& f, Axis_t min, Axis_t max); - void ExitFit(TF1& f); + Bool_t InitFit(Axis_t min, Axis_t max); + void ExitFit(TF1 *f); inline static Double_t fFitFuncMichele(Double_t *x, Double_t *par)