source: trunk/MagicSoft/Mars/mhist/MHEffectiveOnTime.cc@ 4934

Last change on this file since 4934 was 4927, checked in by tbretz, 20 years ago
*** empty log message ***
File size: 19.7 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): Thomas Bretz, 8/2002 <mailto:tbretz@astro.uni-wuerzburg.de>
19! Author(s): Wolfgang Wittek, 1/2002 <mailto:wittek@mppmu.mpg.de>
20!
21! Copyright: MAGIC Software Development, 2000-2004
22!
23!
24\* ======================================================================== */
25
26//////////////////////////////////////////////////////////////////////////////
27//
28// MHEffectiveOnTime
29//
30// Filling this you will get the effective on-time versus theta and
31// observation time.
32//
33// From this histogram the effective on-time is determined by a fit.
34// The result of the fit (see Fit()) and the fit-parameters (like chi^2)
35// are stored in corresponding histograms
36//
37// To determin the efective on time a poisson fit is done. For more details
38// please have a look into the source code of FitH() it should be simple
39// to understand. In this function a Delta-T distribution is fitted, while
40// Delta-T is the time between two consecutive events.
41//
42// The fit is done for projections of a 2D histogram in Theta and Delta T.
43// So you get the effective on time versus theta.
44//
45// To get the effective on-time versus time a histogram is filled with
46// the Delta-T distribution of a number of events set by SetNumEvents().
47// The default is 12000 (roughly 1min at 200Hz)
48//
49// For each "time-bin" the histogram is fitted and the resulting effective
50// on-time is stored in the fHEffOnTime histogram. Each entry in this
51// histogram is the effective observation time between the upper and
52// lower edges of the bins.
53//
54// In addition the calculated effective on time is stored in a
55// "MEffectiveOnTime [MParameterDerr]" and the corresponding time-stamp
56// (the upper edge of the bin) "MTimeEffectiveOnTime [MTime]"
57//
58// The class takes two binnings from the Parameter list; if these binnings
59// are not available the defaultbinning is used:
60// MBinning("BinningDeltaT"); // Units of seconds
61// MBinning("BinningTheta"); // Units of degrees
62//
63//
64// Usage:
65// ------
66// MFillH fill("MHEffectiveOnTime", "MTime");
67// tlist.AddToList(&fill);
68//
69//
70// Input Container:
71// MPointingPos
72//
73// Output Container:
74// MEffectiveOnTime [MParameterDerr]
75// MTimeEffectiveOnTime [MTime]
76//
77//////////////////////////////////////////////////////////////////////////////
78#include "MHEffectiveOnTime.h"
79
80#include <TF1.h>
81#include <TMinuit.h>
82#include <TRandom.h>
83
84#include <TLatex.h>
85#include <TCanvas.h>
86#include <TPaveStats.h>
87
88#include "MTime.h"
89#include "MParameters.h"
90#include "MPointingPos.h"
91
92#include "MBinning.h"
93#include "MParList.h"
94
95#include "MLog.h"
96#include "MLogManip.h"
97
98ClassImp(MHEffectiveOnTime);
99
100using namespace std;
101
102// --------------------------------------------------------------------------
103//
104// Default Constructor. It initializes all histograms.
105//
106MHEffectiveOnTime::MHEffectiveOnTime(const char *name, const char *title)
107 : fPointPos(0), fTime(0), fParam(0), fIsFinalized(kFALSE),
108 fNumEvents(200*60), fNameProjDeltaT(Form("DeltaT_%p", this)),
109 fNameProjTheta(Form("Theta_%p", this))
110{
111 //
112 // set the name and title of this object
113 //
114 fName = name ? name : "MHEffectiveOnTime";
115 fTitle = title ? title : "Histogram to determin effective On-Time vs Time and Zenith Angle";
116
117 // Main histogram
118 fH2DeltaT.SetName("DeltaT");
119 fH2DeltaT.SetXTitle("\\Delta t [s]");
120 fH2DeltaT.SetYTitle("\\Theta [\\circ]");
121 fH2DeltaT.SetZTitle("Count");
122 fH2DeltaT.UseCurrentStyle();
123 fH2DeltaT.SetDirectory(NULL);
124
125 // Main histogram
126 fH1DeltaT.SetName("DeltaT");
127 fH1DeltaT.SetXTitle("\\Delta t [s]");
128 fH1DeltaT.SetYTitle("Counts");
129 fH1DeltaT.UseCurrentStyle();
130 fH1DeltaT.SetDirectory(NULL);
131
132 // effective on time versus theta
133 fHEffOnTheta.SetName("EffOnTheta");
134 fHEffOnTheta.SetTitle("Effective On Time T_{eff}");
135 fHEffOnTheta.SetXTitle("\\Theta [\\circ]");
136 fHEffOnTheta.SetYTitle("T_{eff} [s]");
137 fHEffOnTheta.UseCurrentStyle();
138 fHEffOnTheta.SetDirectory(NULL);
139 fHEffOnTheta.GetYaxis()->SetTitleOffset(1.2);
140 //fHEffOn.Sumw2();
141
142 // effective on time versus time
143 fHEffOnTime.SetName("EffOnTime");
144 fHEffOnTime.SetTitle("Effective On Time T_{eff}");
145 fHEffOnTime.SetXTitle("Time");
146 fHEffOnTime.SetYTitle("T_{eff} [s]");
147 fHEffOnTime.UseCurrentStyle();
148 fHEffOnTime.SetDirectory(NULL);
149 fHEffOnTime.GetYaxis()->SetTitleOffset(1.2);
150 fHEffOnTime.GetXaxis()->SetLabelSize(0.033);
151 fHEffOnTime.GetXaxis()->SetTimeFormat("%H:%M:%S %F1995-01-01 00:00:00 GMT");
152 fHEffOnTime.GetXaxis()->SetTimeDisplay(1);
153 fHEffOnTime.Sumw2();
154
155 // chi2 probability
156 fHProbTheta.SetName("ProbTheta");
157 fHProbTheta.SetTitle("\\chi^{2} Probability of Fit");
158 fHProbTheta.SetXTitle("\\Theta [\\circ]");
159 fHProbTheta.SetYTitle("p [%]");
160 fHProbTheta.UseCurrentStyle();
161 fHProbTheta.SetDirectory(NULL);
162 fHProbTheta.GetYaxis()->SetTitleOffset(1.2);
163 fHProbTheta.SetMaximum(101);
164
165 // chi2 probability
166 fHProbTime.SetName("ProbTime");
167 fHProbTime.SetTitle("\\chi^{2} Probability of Fit");
168 fHProbTime.SetXTitle("Time");
169 fHProbTime.SetYTitle("p [%]");
170 fHProbTime.UseCurrentStyle();
171 fHProbTime.SetDirectory(NULL);
172 fHProbTime.GetYaxis()->SetTitleOffset(1.2);
173 fHProbTime.GetXaxis()->SetLabelSize(0.033);
174 fHProbTime.GetXaxis()->SetTimeFormat("%H:%M:%S %F1995-01-01 00:00:00 GMT");
175 fHProbTime.GetXaxis()->SetTimeDisplay(1);
176 fHProbTime.SetMaximum(101);
177
178 // lambda versus theta
179 fHLambda.SetName("lambda");
180 fHLambda.SetTitle("lambda of Effective On Time Fit");
181 fHLambda.SetXTitle("\\Theta [\\circ]");
182 fHLambda.SetYTitle("\\lambda [Hz]");
183 fHLambda.UseCurrentStyle();
184 fHLambda.SetDirectory(NULL);
185 // fHLambda.Sumw2();
186
187 // N0 versus theta
188 fHN0.SetName("N0");
189 fHN0.SetTitle("Ideal number of events N_{0}");
190 fHN0.SetXTitle("\\Theta [\\circ]");
191 fHN0.SetYTitle("N_{0}");
192 fHN0.UseCurrentStyle();
193 fHN0.SetDirectory(NULL);
194 //fHN0del.Sumw2();
195
196 // chi2/NDF versus theta
197 /*
198 fHChi2.SetName("Chi2/NDF");
199 fHChi2.SetTitle("\\chi^{2}/NDF of Effective On Time Fit");
200 fHChi2.SetXTitle("\\Theta [\\circ]");
201 fHChi2.SetYTitle("\\chi^{2}/NDF");
202 fHChi2.UseCurrentStyle();
203 fHChi2.SetDirectory(NULL);
204 fHChi2.GetYaxis()->SetTitleOffset(1.2);
205 //fHChi2.Sumw2();
206 */
207
208 // setup binning
209 MBinning btheta("BinningTheta");
210 btheta.SetEdgesCos(101, 0, 60);
211
212 MBinning btime("BinningDeltaT");
213 btime.SetEdges(50, 0, 0.1);
214
215 MH::SetBinning(&fH2DeltaT, &btime, &btheta);
216
217 btime.Apply(fH1DeltaT);
218
219 btheta.Apply(fHEffOnTheta);
220 btheta.Apply(fHLambda);
221 btheta.Apply(fHN0);
222 btheta.Apply(fHProbTheta);
223 //btheta.Apply(fHChi2);
224}
225
226// --------------------------------------------------------------------------
227//
228// Set the binnings and prepare the filling of the histogram
229//
230Bool_t MHEffectiveOnTime::SetupFill(const MParList *plist)
231{
232 fPointPos = (MPointingPos*)plist->FindObject("MPointingPos");
233 if (!fPointPos)
234 {
235 *fLog << err << dbginf << "MPointingPos not found... aborting." << endl;
236 return kFALSE;
237 }
238
239 // FIXME: Remove const-qualifier from base-class!
240 fTime = (MTime*)const_cast<MParList*>(plist)->FindCreateObj("MTime", "MTimeEffectiveOnTime");
241 if (!fTime)
242 return kFALSE;
243 fParam = (MParameterDerr*)const_cast<MParList*>(plist)->FindCreateObj("MParameterDerr", "MEffectiveOnTime");
244 if (!fParam)
245 return kFALSE;
246
247 const MBinning* binsdtime = (MBinning*)plist->FindObject("BinningDeltaT");
248 const MBinning* binstheta = (MBinning*)plist->FindObject("BinningTheta");
249 if (binsdtime)
250 binsdtime->Apply(fH1DeltaT);
251 if (binstheta)
252 {
253 binstheta->Apply(fHEffOnTheta);
254 binstheta->Apply(fHLambda);
255 binstheta->Apply(fHN0);
256 binstheta->Apply(fHProbTheta);
257 //binstheta->Apply(fHChi2);
258 }
259 if (binstheta && binsdtime)
260 SetBinning(&fH2DeltaT, binsdtime, binstheta);
261
262 return kTRUE;
263}
264
265// --------------------------------------------------------------------------
266//
267// Fit a single Delta-T distribution. See source code for more details
268//
269Bool_t MHEffectiveOnTime::FitH(TH1D *h, Double_t *res, Bool_t paint) const
270{
271 const Double_t Nm = h->Integral();
272
273 // FIXME: Do fit only if contents of bin has changed
274 if (Nm<=0)
275 return kFALSE;
276
277 // determine range (yq[0], yq[1]) of time differences
278 // where fit should be performed;
279 // require a fraction >=xq[0] of all entries to lie below yq[0]
280 // and a fraction <=xq[1] of all entries to lie below yq[1];
281 // within the range (yq[0], yq[1]) there must be no empty bin;
282 // choose pedestrian approach as long as GetQuantiles is not available
283 Double_t xq[2] = { 0.05, 0.95 };
284 Double_t yq[2];
285 h->GetQuantiles(2, yq, xq);
286
287 // Nmdel = Nm * binwidth, with Nm = number of observed events
288 const Double_t Nmdel = h->Integral("width");
289
290 //
291 // Setup Poisson function for the fit:
292 // lambda [Hz], N0 = ideal no of evts, del = bin width of dt
293 //
294 // parameter 0 = lambda
295 // parameter 1 = N0*del
296 //
297 TF1 func("Poisson", " [1]*[2] * [0] * exp(-[0] *x)", yq[0], yq[1]);
298 func.SetParNames("lambda", "N0", "del");
299
300 func.SetParameter(0, 100); // Hz
301 func.SetParameter(1, Nm);
302 func.FixParameter(2, Nmdel/Nm);
303
304 // options : 0 do not plot the function
305 // I use integral of function in bin rather than value at bin center
306 // R use the range specified in the function range
307 // Q quiet mode
308 h->Fit(&func, "0IRQ");
309
310 const Double_t chi2 = func.GetChisquare();
311 const Int_t NDF = func.GetNDF();
312
313 // was fit successful ?
314 const Bool_t ok = NDF>0 && chi2<2.5*NDF;
315
316 if (paint)
317 {
318 func.SetLineWidth(2);
319 func.SetLineColor(ok ? kGreen : kRed);
320 func.Paint("same");
321 }
322
323 if (!ok)
324 return kFALSE;
325
326 const Double_t lambda = func.GetParameter(0);
327 const Double_t N0 = func.GetParameter(1);
328 const Double_t prob = func.GetProb();
329
330 /*
331 *fLog << all << "Nm/lambda=" << Nm/lambda << " chi2/NDF=";
332 *fLog << (NDF ? chi2/NDF : 0.0) << " lambda=";
333 *fLog << lambda << " N0=" << N0 << endl;
334 */
335
336 Double_t emat[2][2];
337 gMinuit->mnemat((Double_t*)emat, 2);
338
339 const Double_t dldl = emat[0][0];
340 //const Double_t dN0dN0 = emat[1][1];
341
342 const Double_t teff = Nm/lambda;
343 const Double_t dteff = teff * TMath::Sqrt(dldl/(lambda*lambda) + 1.0/Nm);
344 const Double_t dl = TMath::Sqrt(dldl);
345
346 //const Double_t kappa = Nm/N0;
347 //const Double_t Rdead = 1.0 - kappa;
348 //const Double_t dRdead = kappa * TMath::Sqrt(dN0dN0/(N0*N0) + 1.0/Nm);
349
350 // the effective on time is Nm/lambda
351 res[0] = teff;
352 res[1] = dteff;
353
354 // plot chi2-probability of fit
355 res[2] = prob*100;
356
357 // plot chi2/NDF of fit
358 //res[3] = NDF ? chi2/NDF : 0.0;
359
360 // lambda of fit
361 res[3] = lambda;
362 res[4] = dl;
363
364 // N0 of fit
365 res[5] = N0;
366
367 // Rdead (from fit) is the fraction from real time lost by the dead time
368 //fHRdead.SetBinContent(i, Rdead);
369 //fHRdead.SetBinError (i,dRdead);
370
371 return kTRUE;
372}
373
374// --------------------------------------------------------------------------
375//
376// Fit a all bins of the distribution in theta. Store the result in the
377// Theta-Histograms
378//
379void MHEffectiveOnTime::FitThetaBins()
380{
381 fHEffOnTheta.Reset();
382 fHProbTheta.Reset();
383 fHLambda.Reset();
384 fHN0.Reset();
385
386 const TString name = Form("CalcTheta%d", (UInt_t)gRandom->Uniform(999999999));
387
388 // nbins = number of Theta bins
389 const Int_t nbins = fH2DeltaT.GetNbinsY();
390
391 TH1D *h=0;
392 for (int i=1; i<=nbins; i++)
393 {
394 // TH1D &h = *hist->ProjectionX("Calc-theta", i, i);
395 h = fH2DeltaT.ProjectionX(name, i, i, "E");
396
397 Double_t res[6];
398 if (!FitH(h, res))
399 continue;
400
401 // the effective on time is Nm/lambda
402 fHEffOnTheta.SetBinContent(i, res[0]);
403 fHEffOnTheta.SetBinError (i, res[1]);
404
405 // plot chi2-probability of fit
406 fHProbTheta.SetBinContent(i, res[2]);
407
408 // plot chi2/NDF of fit
409 //fHChi2.SetBinContent(i, res[3]);
410
411 // lambda of fit
412 fHLambda.SetBinContent(i, res[3]);
413 fHLambda.SetBinError (i, res[4]);
414
415 // N0 of fit
416 fHN0.SetBinContent(i, res[5]);
417
418 // Rdead (from fit) is the fraction from real time lost by the dead time
419 //fHRdead.SetBinContent(i, Rdead);
420 //fHRdead.SetBinError (i,dRdead);
421 }
422
423 // Histogram is reused via gROOT->FindObject()
424 // Need to be deleted only once
425 if (h)
426 delete h;
427}
428
429// --------------------------------------------------------------------------
430//
431// Fit the single-time-bin histogram. Store the result in the
432// Time-Histograms
433//
434void MHEffectiveOnTime::FitTimeBin()
435{
436 //
437 // Fit histogram
438 //
439 Double_t res[6];
440 if (!FitH(&fH1DeltaT, res))
441 return;
442
443 // Reset Histogram
444 fH1DeltaT.Reset();
445
446 //
447 // Prepare Histogram
448 //
449
450 // Get x-axis
451 TAxis &x = *fHEffOnTime.GetXaxis();
452
453 // Get number of bins
454 const Int_t n = x.GetNbins();
455
456 // Enhance binning
457 MBinning bins;
458 bins.SetEdges(x);
459 bins.AddEdge(fLastTime.GetAxisTime());
460 bins.Apply(fHEffOnTime);
461 bins.Apply(fHProbTime);
462
463 //
464 // Fill histogram
465 //
466 fHEffOnTime.SetBinContent(n, res[0]);
467 fHEffOnTime.SetBinError(n, res[1]);
468
469 fHProbTime.SetBinContent(n, res[2]);
470
471 //
472 // Now prepare output
473 //
474 fParam->SetVal(res[0], res[1]);
475 fParam->SetReadyToSave();
476
477 *fTime = fLastTime;
478
479 // Include the current event
480 fTime->Plus1ns();
481
482 *fLog << fLastTime << ": Val=" << res[0] << " Err=" << res[1] << endl;
483}
484
485// --------------------------------------------------------------------------
486//
487// Fill the histogram
488//
489Bool_t MHEffectiveOnTime::Fill(const MParContainer *par, const Stat_t w)
490{
491 const MTime *time = dynamic_cast<const MTime*>(par);
492 if (!time)
493 {
494 *fLog << err << "ERROR - MHEffectiveOnTime::Fill without argument or container doesn't inherit from MTime... abort." << endl;
495 return kFALSE;
496 }
497
498 //
499 // If this is the first call we have to initialize the time-histogram
500 //
501 if (fLastTime==MTime())
502 {
503 MBinning bins;
504 bins.SetEdges(2, time->GetAxisTime()-fNumEvents/200, time->GetAxisTime());
505 bins.Apply(fHEffOnTime);
506 bins.Apply(fHProbTime);
507
508 fParam->SetVal(0, 0);
509 fParam->SetReadyToSave();
510
511 *fTime = *time;
512
513 // Make this 1ns before the first event!
514 fTime->Minus1ns();
515 }
516
517 //
518 // Fill time difference into the histograms
519 //
520 const Double_t dt = *time-fLastTime;
521
522 fH2DeltaT.Fill(dt, fPointPos->GetZd(), w);
523 fH1DeltaT.Fill(dt, w);
524
525 fLastTime = *time;
526
527 //
528 // If we reached the event number limit for the time-bins fit the histogram
529 //
530 if (fH1DeltaT.GetEntries()>=fNumEvents)
531 FitTimeBin();
532
533 return kTRUE;
534}
535
536// --------------------------------------------------------------------------
537//
538// Fit the theta projections of the 2D histogram and the 1D Delta-T
539// distribution
540//
541Bool_t MHEffectiveOnTime::Finalize()
542{
543 FitThetaBins();
544 FitTimeBin();
545
546 fIsFinalized = kTRUE;
547
548 return kTRUE;
549}
550
551// --------------------------------------------------------------------------
552//
553// Paint the integral and the error on top of the histogram
554//
555void MHEffectiveOnTime::PaintText(Double_t val, Double_t error) const
556{
557 TLatex text(0.45, 0.94, Form("T_{eff} = %.1fs \\pm %.1fs", val, error));
558 text.SetBit(TLatex::kTextNDC);
559 text.SetTextSize(0.04);
560 text.Paint();
561}
562
563// --------------------------------------------------------------------------
564//
565// Prepare painting the histograms
566//
567void MHEffectiveOnTime::Paint(Option_t *opt)
568{
569 TH1D *h=0;
570 TPaveStats *st=0;
571
572 TString o(opt);
573 if (o==(TString)"fit")
574 {
575 TVirtualPad *pad = gPad;
576
577 for (int x=0; x<2; x++)
578 for (int y=0; y<3; y++)
579 {
580 TVirtualPad *p=gPad->GetPad(x+1)->GetPad(y+1);
581 if ((st = (TPaveStats*)p->GetPrimitive("stats")))
582 {
583 if (st->GetOptStat()==11)
584 continue;
585
586 const Double_t y1 = st->GetY1NDC();
587 const Double_t y2 = st->GetY2NDC();
588 const Double_t x1 = st->GetX1NDC();
589 const Double_t x2 = st->GetX2NDC();
590
591 st->SetY1NDC((y2-y1)/3+y1);
592 st->SetX1NDC((x2-x1)/3+x1);
593 st->SetOptStat(11);
594 }
595 }
596
597 pad->GetPad(1)->cd(1);
598 if ((h = (TH1D*)gPad->FindObject(fNameProjDeltaT)))
599 {
600 h = fH2DeltaT.ProjectionX(fNameProjDeltaT, -1, 9999, "E");
601 if (h->GetEntries()>0)
602 gPad->SetLogy();
603 }
604
605 pad->GetPad(2)->cd(1);
606 if ((h = (TH1D*)gPad->FindObject(fNameProjTheta)))
607 fH2DeltaT.ProjectionY(fNameProjTheta, -1, 9999, "E");
608
609 if (!fIsFinalized)
610 FitThetaBins();
611 return;
612 }
613 if (o==(TString)"paint")
614 {
615 if ((h = (TH1D*)gPad->FindObject(fNameProjDeltaT)))
616 {
617 Double_t res[6];
618 if (FitH(h, res, kTRUE))
619 PaintText(res[0], res[1]);
620 }
621 return;
622 }
623
624 h=0;
625 if (o==(TString)"theta")
626 h = &fHEffOnTheta;
627 if (o==(TString)"time")
628 h = &fHEffOnTime;
629
630 if (!h)
631 return;
632
633 Double_t error = 0;
634 for (int i=0; i<h->GetXaxis()->GetNbins(); i++)
635 error += h->GetBinError(i);
636
637 PaintText(h->Integral(), error);
638}
639
640// --------------------------------------------------------------------------
641//
642// Draw the histogram
643//
644void MHEffectiveOnTime::Draw(Option_t *opt)
645{
646 TVirtualPad *pad = gPad ? gPad : MakeDefCanvas(this);
647 pad->SetBorderMode(0);
648
649 AppendPad("fit");
650
651 pad->Divide(2, 1, 0, 0);
652
653 TH1 *h;
654
655 pad->cd(1);
656 gPad->SetBorderMode(0);
657 gPad->Divide(1, 3, 0, 0);
658 pad->GetPad(1)->cd(1);
659 gPad->SetBorderMode(0);
660 h = fH2DeltaT.ProjectionX(fNameProjDeltaT, -1, 9999, "E");
661 h->SetTitle("Distribution of \\Delta t [s]");
662 h->SetXTitle("\\Delta t [s]");
663 h->SetYTitle("Counts");
664 h->SetDirectory(NULL);
665 h->SetMarkerStyle(kFullDotMedium);
666 h->SetBit(kCanDelete);
667 h->Draw();
668 AppendPad("paint");
669
670 pad->GetPad(1)->cd(2);
671 gPad->SetBorderMode(0);
672 fHProbTime.Draw();
673 AppendPad("paint");
674
675 pad->GetPad(1)->cd(3);
676 gPad->SetBorderMode(0);
677 fHEffOnTime.Draw();
678 AppendPad("time");
679
680 pad->cd(2);
681 gPad->SetBorderMode(0);
682 gPad->Divide(1, 3, 0, 0);
683
684 pad->GetPad(2)->cd(1);
685 gPad->SetBorderMode(0);
686 h = fH2DeltaT.ProjectionY(fNameProjTheta, -1, 9999, "E");
687 h->SetTitle("Distribution of \\Theta [\\circ]");
688 h->SetXTitle("\\Theta [\\circ]");
689 h->SetYTitle("Counts");
690 h->SetDirectory(NULL);
691 h->SetMarkerStyle(kFullDotMedium);
692 h->SetBit(kCanDelete);
693 h->GetYaxis()->SetTitleOffset(1.1);
694 h->Draw();
695
696 pad->GetPad(2)->cd(2);
697 gPad->SetBorderMode(0);
698 fHProbTheta.Draw();
699
700 pad->GetPad(2)->cd(3);
701 gPad->SetBorderMode(0);
702 fHEffOnTheta.Draw();
703 AppendPad("theta");
704}
Note: See TracBrowser for help on using the repository browser.