source: trunk/MagicSoft/Mars/mhist/MHCalibrationPixel.cc@ 2632

Last change on this file since 2632 was 2629, checked in by gaug, 21 years ago
*** empty log message ***
File size: 21.3 KB
Line 
1/* ======================================================================== *\
2!
3! *
4! * This file is part of MARS, the MAGIC Analysis and Reconstruction
5! * Software. It is distributed to you in the hope that it can be a useful
6! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
7! * It is distributed WITHOUT ANY WARRANTY.
8! *
9! * Permission to use, copy, modify and distribute this software and its
10! * documentation for any purpose is hereby granted without fee,
11! * provided that the above copyright notice appear in all copies and
12! * that both that copyright notice and this permission notice appear
13! * in supporting documentation. It is provided "as is" without express
14! * or implied warranty.
15! *
16!
17!
18! Author(s): Markus Gaug 11/2003 <mailto:markus@ifae.es>
19!
20! Copyright: MAGIC Software Development, 2000-2002
21!
22!
23\* ======================================================================== */
24
25//////////////////////////////////////////////////////////////////////////////
26// //
27// MHCalibrationPixel //
28// //
29// Performs all the necessary fits to extract the mean number of photons //
30// out of the derived light flux //
31// //
32//////////////////////////////////////////////////////////////////////////////
33
34#include "MHCalibrationPixel.h"
35#include "MHCalibrationConfig.h"
36#include "MCalibrationFits.h"
37
38#include <TStyle.h>
39#include <TMath.h>
40
41#include <TFitter.h>
42
43#include <TF1.h>
44#include <TH2.h>
45#include <TCanvas.h>
46#include <TPad.h>
47#include <TPaveText.h>
48
49#include "MParList.h"
50
51#include "MLog.h"
52#include "MLogManip.h"
53
54ClassImp(MHCalibrationPixel);
55
56using namespace std;
57
58// --------------------------------------------------------------------------
59//
60// Default Constructor.
61//
62MHCalibrationPixel::MHCalibrationPixel(const char *name, const char *title)
63 : fPixId(-1),
64 fChargeGausFit(NULL),
65 fTimeGausFit(NULL),
66 fFitLegend(NULL),
67 fLowerFitRange(-2000.),
68 fChargeFirstHiGain(-2000.5),
69 fChargeLastHiGain(9999.5),
70 fChargeNbinsHiGain(12000),
71 fChargeFirstLoGain(-2000.5),
72 fChargeLastLoGain(9999.5),
73 fChargeNbinsLoGain(1200),
74 fFitOK(kFALSE),
75 fChargeChisquare(-1.),
76 fChargeProb(-1.),
77 fChargeNdf(-1),
78 fTimeChisquare(-1.),
79 fTimeProb(-1.),
80 fTimeNdf(-1),
81 fUseLoGain(kFALSE)
82{
83
84 fName = name ? name : "MHCalibrationPixel";
85 fTitle = title ? title : "Fill the accumulated charges and times of all events and perform fits";
86
87 // Create a large number of bins, later we will rebin
88 fHChargeHiGain = new TH1F("HChargeHiGain","Distribution of Summed FADC Hi Gain Slices Pixel ",
89 fChargeNbinsHiGain,fChargeFirstHiGain,fChargeLastHiGain);
90 fHChargeHiGain->SetXTitle("Sum FADC Slices (Hi Gain)");
91 fHChargeHiGain->SetYTitle("Nr. of events");
92 fHChargeHiGain->Sumw2();
93
94 fHChargeHiGain->SetDirectory(NULL);
95
96 fHChargeLoGain = new TH1F("HChargeLoGain","Distribution of Summed FADC Lo Gain Slices Pixel ",
97 fChargeNbinsLoGain,fChargeFirstLoGain,fChargeLastLoGain);
98 fHChargeLoGain->SetXTitle("Sum FADC Slices (Lo Gain)");
99 fHChargeLoGain->SetYTitle("Nr. of events");
100 fHChargeLoGain->Sumw2();
101
102 fHChargeLoGain->SetDirectory(NULL);
103
104 Axis_t tfirst = -0.5;
105 Axis_t tlast = 15.5;
106 Int_t ntbins = 16;
107
108 fHTimeHiGain = new TH1I("HTimeHiGain","Distribution of Mean Arrival Hi Gain Times Pixel ",
109 ntbins,tfirst,tlast);
110 fHTimeHiGain->SetXTitle("Mean Arrival Times [Hi Gain FADC slice nr]");
111 fHTimeHiGain->SetYTitle("Nr. of events");
112 fHTimeHiGain->Sumw2();
113
114 fHTimeHiGain->SetDirectory(NULL);
115
116 fHTimeLoGain = new TH1I("HTimeLoGain","Distribution of Mean Arrival Lo Gain Times Pixel ",
117 ntbins,tfirst,tlast);
118 fHTimeLoGain->SetXTitle("Mean Arrival Times [Lo Gain FADC slice nr]");
119 fHTimeLoGain->SetYTitle("Nr. of events");
120 fHTimeLoGain->Sumw2();
121
122 fHTimeLoGain->SetDirectory(NULL);
123
124 // We define a reasonable number and later enlarge it if necessary
125 Int_t nqbins = 20000;
126 Axis_t nfirst = -0.5;
127 Axis_t nlast = (Axis_t)nqbins - 0.5;
128
129 fHChargevsNHiGain = new TH1I("HChargevsNHiGain","Sum of Hi Gain Charges vs. Event Number Pixel ",
130 nqbins,nfirst,nlast);
131 fHChargevsNHiGain->SetXTitle("Event Nr.");
132 fHChargevsNHiGain->SetYTitle("Sum of Hi Gain FADC slices");
133
134 fHChargevsNHiGain->SetDirectory(NULL);
135
136 fHChargevsNLoGain = new TH1I("HChargevsNLoGain","Sum of Lo Gain Charges vs. Event Number Pixel ",
137 nqbins,nfirst,nlast);
138 fHChargevsNLoGain->SetXTitle("Event Nr.");
139 fHChargevsNLoGain->SetYTitle("Sum of Lo Gain FADC slices");
140
141 fHChargevsNLoGain->SetDirectory(NULL);
142
143}
144
145MHCalibrationPixel::~MHCalibrationPixel()
146{
147
148 delete fHChargeHiGain;
149 delete fHTimeHiGain;
150 delete fHChargevsNHiGain;
151
152 delete fHChargeLoGain;
153 delete fHTimeLoGain;
154 delete fHChargevsNLoGain;
155
156 if (fChargeGausFit)
157 delete fChargeGausFit;
158 if (fTimeGausFit)
159 delete fTimeGausFit;
160 if (fFitLegend)
161 delete fFitLegend;
162
163}
164
165
166void MHCalibrationPixel::ChangeHistId(Int_t id)
167{
168
169 fPixId = id;
170
171 TString nameQHiGain = TString(fHChargeHiGain->GetName());
172 nameQHiGain += id;
173 fHChargeHiGain->SetName(nameQHiGain.Data());
174
175 TString nameTHiGain = TString(fHTimeHiGain->GetName());
176 nameTHiGain += id;
177 fHTimeHiGain->SetName(nameTHiGain.Data());
178
179 TString nameQvsNHiGain = TString(fHChargevsNHiGain->GetName());
180 nameQvsNHiGain += id;
181 fHChargevsNHiGain->SetName(nameQvsNHiGain.Data());
182
183 TString titleQHiGain = TString(fHChargeHiGain->GetTitle());
184 titleQHiGain += id;
185 fHChargeHiGain->SetTitle(titleQHiGain.Data());
186
187 TString titleTHiGain = TString(fHTimeHiGain->GetTitle());
188 titleTHiGain += id;
189 fHTimeHiGain->SetTitle(titleTHiGain.Data());
190
191 TString titleQvsNHiGain = TString(fHChargevsNHiGain->GetTitle());
192 titleQvsNHiGain += id;
193 fHChargevsNHiGain->SetTitle(titleQvsNHiGain.Data());
194
195
196 TString nameQLoGain = TString(fHChargeLoGain->GetName());
197 nameQLoGain += id;
198 fHChargeLoGain->SetName(nameQLoGain.Data());
199
200 TString nameTLoGain = TString(fHTimeLoGain->GetName());
201 nameTLoGain += id;
202 fHTimeLoGain->SetName(nameTLoGain.Data());
203
204 TString nameQvsNLoGain = TString(fHChargevsNLoGain->GetName());
205 nameQvsNLoGain += id;
206 fHChargevsNLoGain->SetName(nameQvsNLoGain.Data());
207
208 TString titleQLoGain = TString(fHChargeLoGain->GetTitle());
209 titleQLoGain += id;
210 fHChargeLoGain->SetTitle(titleQLoGain.Data());
211
212 TString titleTLoGain = TString(fHTimeLoGain->GetTitle());
213 titleTLoGain += id;
214 fHTimeLoGain->SetTitle(titleTLoGain.Data());
215
216 TString titleQvsNLoGain = TString(fHChargevsNLoGain->GetTitle());
217 titleQvsNLoGain += id;
218 fHChargevsNLoGain->SetTitle(titleQvsNLoGain.Data());
219}
220
221
222void MHCalibrationPixel::Reset()
223{
224
225 for (Int_t i = fHChargeHiGain->FindBin(fChargeFirstHiGain);
226 i <= fHChargeHiGain->FindBin(fChargeLastHiGain); i++)
227 fHChargeHiGain->SetBinContent(i, 1.e-20);
228
229 for (Int_t i = 0; i < 16; i++)
230 fHTimeHiGain->SetBinContent(i, 1.e-20);
231
232 fChargeLastHiGain = 9999.5;
233
234 fHChargeLoGain->GetXaxis()->SetRangeUser(0.,fChargeLastLoGain);
235
236 for (Int_t i = fHChargeLoGain->FindBin(fChargeFirstLoGain);
237 i <= fHChargeLoGain->FindBin(fChargeLastLoGain); i++)
238 fHChargeLoGain->SetBinContent(i, 1.e-20);
239
240 for (Int_t i = 0; i < 16; i++)
241 fHTimeLoGain->SetBinContent(i, 1.e-20);
242
243 fChargeLastLoGain = 9999.5;
244
245 fHChargeLoGain->GetXaxis()->SetRangeUser(0.,fChargeLastLoGain);
246
247 return;
248}
249
250
251// -------------------------------------------------------------------------
252//
253// Set the binnings and prepare the filling of the histograms
254//
255Bool_t MHCalibrationPixel::SetupFill(const MParList *plist)
256{
257
258 fHChargeHiGain->Reset();
259 fHTimeHiGain->Reset();
260 fHChargeLoGain->Reset();
261 fHTimeLoGain->Reset();
262
263 return kTRUE;
264}
265
266
267Bool_t MHCalibrationPixel::UseLoGain()
268{
269
270 if (fHChargeHiGain->GetEntries() > fHChargeLoGain->GetEntries())
271 {
272 fUseLoGain = kFALSE;
273 return kFALSE;
274 }
275 else
276 {
277 fUseLoGain = kTRUE;
278 return kTRUE;
279 }
280}
281
282
283// -------------------------------------------------------------------------
284//
285// Draw a legend with the fit results
286//
287void MHCalibrationPixel::DrawLegend()
288{
289
290 fFitLegend = new TPaveText(0.05,0.05,0.95,0.95);
291
292 if (fFitOK)
293 fFitLegend->SetFillColor(80);
294 else
295 fFitLegend->SetFillColor(2);
296
297 fFitLegend->SetLabel("Results of the Gauss Fit:");
298 fFitLegend->SetTextSize(0.05);
299
300 const TString line1 =
301 Form("Mean: Q_{#mu} = %2.2f #pm %2.2f",fChargeMean,fChargeMeanErr);
302
303 fFitLegend->AddText(line1);
304
305
306 const TString line4 =
307 Form("Sigma: #sigma_{Q} = %2.2f #pm %2.2f",fChargeSigma,fChargeSigmaErr);
308
309 fFitLegend->AddText(line4);
310
311
312 const TString line7 =
313 Form("#chi^{2} / N_{dof}: %4.2f / %3i",fChargeChisquare,fChargeNdf);
314
315 fFitLegend->AddText(line7);
316
317
318 const TString line8 =
319 Form("Probability: %4.3f ",fChargeProb);
320
321 fFitLegend->AddText(line8);
322
323 if (fFitOK)
324 fFitLegend->AddText("Result of the Fit: OK");
325 else
326 fFitLegend->AddText("Result of the Fit: NOT OK");
327
328 fFitLegend->SetBit(kCanDelete);
329 fFitLegend->Draw();
330
331}
332
333
334// -------------------------------------------------------------------------
335//
336// Draw the histogram
337//
338void MHCalibrationPixel::Draw(Option_t *opt)
339{
340
341 gStyle->SetOptFit(0);
342 gStyle->SetOptStat(1111111);
343
344 TCanvas *c = MakeDefCanvas(this,600,900);
345
346 gROOT->SetSelectedPad(NULL);
347
348 c->Divide(2,4);
349
350 c->cd(1);
351 gPad->SetBorderMode(0);
352 gPad->SetLogy(1);
353 gPad->SetTicks();
354
355 fHChargeHiGain->DrawCopy(opt);
356
357 c->Modified();
358 c->Update();
359
360 if (fUseLoGain)
361 {
362 c->cd(2);
363 gPad->SetLogy(1);
364 gPad->SetTicks();
365 fHChargeLoGain->DrawCopy(opt);
366
367 if (fChargeGausFit)
368 {
369 if (fFitOK)
370 fChargeGausFit->SetLineColor(kGreen);
371 else
372 fChargeGausFit->SetLineColor(kRed);
373
374 fChargeGausFit->DrawCopy("same");
375 }
376
377 c->cd(4);
378 DrawLegend();
379
380 }
381 else
382 {
383 if (fChargeGausFit)
384 {
385 if (fFitOK)
386 fChargeGausFit->SetLineColor(kGreen);
387 else
388 fChargeGausFit->SetLineColor(kRed);
389
390 fChargeGausFit->DrawCopy("same");
391 }
392 c->cd(2);
393 gPad->SetLogy(1);
394 gPad->SetTicks();
395 fHChargeLoGain->DrawCopy(opt);
396 c->cd(3);
397 DrawLegend();
398 }
399
400 c->Modified();
401 c->Update();
402
403 c->cd(5);
404 gStyle->SetOptStat(1111111);
405
406 gPad->SetLogy(1);
407 fHTimeHiGain->DrawCopy(opt);
408 c->Modified();
409 c->Update();
410
411 if (fUseLoGain)
412 {
413
414 c->cd(6);
415 gPad->SetLogy(1);
416 fHTimeLoGain->DrawCopy(opt);
417 c->Modified();
418 c->Update();
419
420 if (fTimeGausFit)
421 {
422 if (fTimeChisquare > 1.)
423 fTimeGausFit->SetLineColor(kRed);
424 else
425 fTimeGausFit->SetLineColor(kGreen);
426
427 fTimeGausFit->DrawCopy("same");
428 }
429 }
430 else
431 {
432 if (fTimeGausFit)
433 {
434 if (fTimeChisquare > 1.)
435 fTimeGausFit->SetLineColor(kRed);
436 else
437 fTimeGausFit->SetLineColor(kGreen);
438
439 fTimeGausFit->DrawCopy("same");
440 c->Modified();
441 c->Update();
442 }
443
444 c->cd(6);
445 gPad->SetLogy(1);
446 fHTimeLoGain->DrawCopy(opt);
447 }
448 c->Modified();
449 c->Update();
450
451 c->cd(7);
452 fHChargevsNHiGain->DrawCopy(opt);
453 c->Modified();
454 c->Update();
455
456 c->cd(8);
457 fHChargevsNLoGain->DrawCopy(opt);
458 c->Modified();
459 c->Update();
460
461}
462
463
464
465Bool_t MHCalibrationPixel::FitTimeHiGain(Axis_t rmin, Axis_t rmax, Option_t *option)
466{
467
468 if (fTimeGausFit)
469 return kFALSE;
470
471 rmin = (rmin != 0.) ? rmin : 4.;
472 rmax = (rmax != 0.) ? rmax : 9.;
473
474 const Stat_t entries = fHTimeHiGain->GetEntries();
475 const Double_t mu_guess = fHTimeHiGain->GetBinCenter(fHTimeHiGain->GetMaximumBin());
476 const Double_t sigma_guess = (rmax - rmin)/2.;
477 const Double_t area_guess = entries/gkSq2Pi;
478
479 TString name = TString("GausTime");
480 name += fPixId;
481 fTimeGausFit = new TF1(name.Data(),"gaus",rmin,rmax);
482
483 if (!fTimeGausFit)
484 {
485 *fLog << err << dbginf << "Could not create fit function for Gauss fit" << endl;
486 return kFALSE;
487 }
488
489 fTimeGausFit->SetParameters(area_guess,mu_guess,sigma_guess);
490 fTimeGausFit->SetParNames("Area","#mu","#sigma");
491 fTimeGausFit->SetParLimits(0,0.,entries);
492 fTimeGausFit->SetParLimits(1,rmin,rmax);
493 fTimeGausFit->SetParLimits(2,0.,(rmax-rmin));
494 fTimeGausFit->SetRange(rmin,rmax);
495
496 fHTimeHiGain->Fit(fTimeGausFit,option);
497
498 rmin = fTimeGausFit->GetParameter(1) - 3.*fTimeGausFit->GetParameter(2);
499 rmax = fTimeGausFit->GetParameter(1) + 3.*fTimeGausFit->GetParameter(2);
500 fTimeGausFit->SetRange(rmin,rmax);
501
502 fHTimeHiGain->Fit(fTimeGausFit,option);
503
504 fTimeChisquare = fTimeGausFit->GetChisquare();
505 fTimeNdf = fTimeGausFit->GetNDF();
506 fTimeProb = fTimeGausFit->GetProb();
507
508 fTimeMean = fTimeGausFit->GetParameter(1);
509 fTimeSigma = fTimeGausFit->GetParameter(2);
510
511 if (fTimeChisquare > 1.)
512 {
513 *fLog << warn << "Fit of the Arrival times failed ! " << endl;
514 return kFALSE;
515 }
516
517 return kTRUE;
518
519}
520
521Bool_t MHCalibrationPixel::FitTimeLoGain(Axis_t rmin, Axis_t rmax, Option_t *option)
522{
523
524 if (fTimeGausFit)
525 return kFALSE;
526
527 rmin = (rmin != 0.) ? rmin : 4.;
528 rmax = (rmax != 0.) ? rmax : 9.;
529
530 const Stat_t entries = fHTimeLoGain->GetEntries();
531 const Double_t mu_guess = fHTimeLoGain->GetBinCenter(fHTimeLoGain->GetMaximumBin());
532 const Double_t sigma_guess = (rmax - rmin)/2.;
533 const Double_t area_guess = entries/gkSq2Pi;
534
535 TString name = TString("GausTime");
536 name += fPixId;
537 fTimeGausFit = new TF1(name.Data(),"gaus",rmin,rmax);
538
539 if (!fTimeGausFit)
540 {
541 *fLog << err << dbginf << "Could not create fit function for Gauss fit" << endl;
542 return kFALSE;
543 }
544
545 fTimeGausFit->SetParameters(area_guess,mu_guess,sigma_guess);
546 fTimeGausFit->SetParNames("Area","#mu","#sigma");
547 fTimeGausFit->SetParLimits(0,0.,entries);
548 fTimeGausFit->SetParLimits(1,rmin,rmax);
549 fTimeGausFit->SetParLimits(2,0.,(rmax-rmin));
550 fTimeGausFit->SetRange(rmin,rmax);
551
552 fHTimeLoGain->Fit(fTimeGausFit,option);
553
554 rmin = fTimeGausFit->GetParameter(1) - 3.*fTimeGausFit->GetParameter(2);
555 rmax = fTimeGausFit->GetParameter(1) + 3.*fTimeGausFit->GetParameter(2);
556 fTimeGausFit->SetRange(rmin,rmax);
557
558 fHTimeLoGain->Fit(fTimeGausFit,option);
559
560 fTimeChisquare = fTimeGausFit->GetChisquare();
561 fTimeNdf = fTimeGausFit->GetNDF();
562 fTimeProb = fTimeGausFit->GetProb();
563
564 fTimeMean = fTimeGausFit->GetParameter(1);
565 fTimeSigma = fTimeGausFit->GetParameter(2);
566
567 if (fTimeChisquare > 1.)
568 {
569 *fLog << warn << "Fit of the Arrival times failed ! " << endl;
570 return kFALSE;
571 }
572
573 return kTRUE;
574
575}
576
577Bool_t MHCalibrationPixel::FitChargeHiGain(Option_t *option)
578{
579
580 if (fChargeGausFit)
581 return kFALSE;
582
583 //
584 // Get the fitting ranges
585 //
586 Axis_t rmin = (fLowerFitRange != 0.) ? fLowerFitRange : fChargeFirstHiGain;
587 Axis_t rmax = 0.;
588
589 //
590 // First guesses for the fit (should be as close to reality as possible,
591 // otherwise the fit goes gaga because of high number of dimensions ...
592 //
593 const Stat_t entries = fHChargeHiGain->GetEntries();
594 const Double_t area_guess = entries/gkSq2Pi;
595 const Double_t mu_guess = fHChargeHiGain->GetBinCenter(fHChargeHiGain->GetMaximumBin());
596 const Double_t sigma_guess = mu_guess/15.;
597
598 TString name = TString("ChargeGausFit");
599 name += fPixId;
600
601 fChargeGausFit = new TF1(name.Data(),"gaus",rmin,fChargeLastHiGain);
602
603 if (!fChargeGausFit)
604 {
605 *fLog << err << dbginf << "Could not create fit function for Gauss fit" << endl;
606 return kFALSE;
607 }
608
609 fChargeGausFit->SetParameters(area_guess,mu_guess,sigma_guess);
610 fChargeGausFit->SetParNames("Area","#mu","#sigma");
611 fChargeGausFit->SetParLimits(0,0.,entries);
612 fChargeGausFit->SetParLimits(1,rmin,fChargeLastHiGain);
613 fChargeGausFit->SetParLimits(2,0.,fChargeLastHiGain-rmin);
614 fChargeGausFit->SetRange(rmin,fChargeLastHiGain);
615
616 fHChargeHiGain->Fit(fChargeGausFit,option);
617
618 Axis_t rtry = fChargeGausFit->GetParameter(1) - 2.5*fChargeGausFit->GetParameter(2);
619
620 rmin = (rtry < rmin ? rmin : rtry);
621 rmax = fChargeGausFit->GetParameter(1) + 2.5*fChargeGausFit->GetParameter(2);
622 fChargeGausFit->SetRange(rmin,rmax);
623
624 fHChargeHiGain->Fit(fChargeGausFit,option);
625
626 // rmin = fChargeGausFit->GetParameter(1) - 2.5*fChargeGausFit->GetParameter(2);
627 // rmax = fChargeGausFit->GetParameter(1) + 2.5*fChargeGausFit->GetParameter(2);
628 // fChargeGausFit->SetRange(rmin,rmax);
629
630 // fHChargeHiGain->Fit(fChargeGausFit,option);
631
632 fChargeChisquare = fChargeGausFit->GetChisquare();
633 fChargeNdf = fChargeGausFit->GetNDF();
634 fChargeProb = fChargeGausFit->GetProb();
635 fChargeMean = fChargeGausFit->GetParameter(1);
636 fChargeMeanErr = fChargeGausFit->GetParError(1);
637 fChargeSigma = fChargeGausFit->GetParameter(2);
638 fChargeSigmaErr = fChargeGausFit->GetParError(2);
639
640 //
641 // The fit result is accepted under condition
642 // The Probability is greater than gkProbLimit (default 0.01 == 99%)
643 //
644 if (fChargeProb < gkProbLimit)
645 {
646 *fLog << warn << "Prob: " << fChargeProb << " is smaller than the allowed value: " << gkProbLimit << endl;
647 fFitOK = kFALSE;
648 return kFALSE;
649 }
650
651
652 fFitOK = kTRUE;
653
654 return kTRUE;
655}
656
657
658Bool_t MHCalibrationPixel::FitChargeLoGain(Option_t *option)
659{
660
661 if (fChargeGausFit)
662 return kFALSE;
663
664 //
665 // Get the fitting ranges
666 //
667 Axis_t rmin = (fLowerFitRange != 0.) ? fLowerFitRange : fChargeFirstLoGain;
668 Axis_t rmax = 0.;
669
670 //
671 // First guesses for the fit (should be as close to reality as possible,
672 // otherwise the fit goes gaga because of high number of dimensions ...
673 //
674 const Stat_t entries = fHChargeLoGain->GetEntries();
675 const Double_t area_guess = entries/gkSq2Pi;
676 const Double_t mu_guess = fHChargeLoGain->GetBinCenter(fHChargeLoGain->GetMaximumBin());
677 const Double_t sigma_guess = mu_guess/15.;
678
679 TString name = TString("ChargeGausFit");
680 name += fPixId;
681
682 fChargeGausFit = new TF1(name.Data(),"gaus",rmin,fChargeLastLoGain);
683
684 if (!fChargeGausFit)
685 {
686 *fLog << err << dbginf << "Could not create fit function for Gauss fit" << endl;
687 return kFALSE;
688 }
689
690 fChargeGausFit->SetParameters(area_guess,mu_guess,sigma_guess);
691 fChargeGausFit->SetParNames("Area","#mu","#sigma");
692 fChargeGausFit->SetParLimits(0,0.,entries);
693 fChargeGausFit->SetParLimits(1,rmin,fChargeLastLoGain);
694 fChargeGausFit->SetParLimits(2,0.,fChargeLastLoGain-rmin);
695 fChargeGausFit->SetRange(rmin,fChargeLastLoGain);
696
697 fHChargeLoGain->Fit(fChargeGausFit,option);
698
699 Axis_t rtry = fChargeGausFit->GetParameter(1) - 2.5*fChargeGausFit->GetParameter(2);
700
701 rmin = (rtry < rmin ? rmin : rtry);
702 rmax = fChargeGausFit->GetParameter(1) + 2.5*fChargeGausFit->GetParameter(2);
703 fChargeGausFit->SetRange(rmin,rmax);
704
705 fHChargeLoGain->Fit(fChargeGausFit,option);
706
707 // rmin = fChargeGausFit->GetParameter(1) - 2.5*fChargeGausFit->GetParameter(2);
708 // rmax = fChargeGausFit->GetParameter(1) + 2.5*fChargeGausFit->GetParameter(2);
709 // fChargeGausFit->SetRange(rmin,rmax);
710
711 // fHChargeLoGain->Fit(fChargeGausFit,option);
712
713 fChargeChisquare = fChargeGausFit->GetChisquare();
714 fChargeNdf = fChargeGausFit->GetNDF();
715 fChargeProb = fChargeGausFit->GetProb();
716 fChargeMean = fChargeGausFit->GetParameter(1);
717 fChargeMeanErr = fChargeGausFit->GetParError(1);
718 fChargeSigma = fChargeGausFit->GetParameter(2);
719 fChargeSigmaErr = fChargeGausFit->GetParError(2);
720
721 //
722 // The fit result is accepted under condition
723 // The Probability is greater than gkProbLimit (default 0.01 == 99%)
724 //
725 if (fChargeProb < gkProbLimit)
726 {
727 *fLog << warn << "Prob: " << fChargeProb << " is smaller than the allowed value: " << gkProbLimit << endl;
728 fFitOK = kFALSE;
729 return kFALSE;
730 }
731
732
733 fFitOK = kTRUE;
734
735 return kTRUE;
736}
737
738
739void MHCalibrationPixel::CutAllEdges()
740{
741
742 Int_t nbins = 30;
743
744 CutEdges(fHChargeHiGain,nbins);
745
746 fChargeFirstHiGain = fHChargeHiGain->GetBinLowEdge(fHChargeHiGain->GetXaxis()->GetFirst());
747 fChargeLastHiGain = fHChargeHiGain->GetBinLowEdge(fHChargeHiGain->GetXaxis()->GetLast())
748 +fHChargeHiGain->GetBinWidth(0);
749 fChargeNbinsHiGain = nbins;
750
751 CutEdges(fHChargeLoGain,nbins);
752
753 fChargeFirstLoGain = fHChargeLoGain->GetBinLowEdge(fHChargeLoGain->GetXaxis()->GetFirst());
754 fChargeLastLoGain = fHChargeLoGain->GetBinLowEdge(fHChargeLoGain->GetXaxis()->GetLast())
755 +fHChargeLoGain->GetBinWidth(0);
756 fChargeNbinsLoGain = nbins;
757
758 CutEdges(fHChargevsNHiGain,0);
759 CutEdges(fHChargevsNLoGain,0);
760
761}
762
763void MHCalibrationPixel::PrintChargeFitResult()
764{
765
766 *fLog << "Results of the Summed Charges Fit: " << endl;
767 *fLog << "Chisquare: " << fChargeChisquare << endl;
768 *fLog << "DoF: " << fChargeNdf << endl;
769 *fLog << "Probability: " << fChargeProb << endl;
770 *fLog << endl;
771
772}
773
774void MHCalibrationPixel::PrintTimeFitResult()
775{
776
777 *fLog << "Results of the Time Slices Fit: " << endl;
778 *fLog << "Chisquare: " << fTimeChisquare << endl;
779 *fLog << "Ndf: " << fTimeNdf << endl;
780 *fLog << "Probability: " << fTimeProb << endl;
781 *fLog << endl;
782
783}
Note: See TracBrowser for help on using the repository browser.