source: tags/Mars-V1.0/mhflux/MHDisp.cc

Last change on this file was 7616, checked in by tbretz, 19 years ago
*** empty log message ***
File size: 21.9 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, 5/2005 <mailto:tbretz@astro.uni-wuerzburg.de>
19!
20! Copyright: MAGIC Software Development, 2000-2005
21!
22!
23\* ======================================================================== */
24
25//////////////////////////////////////////////////////////////////////////////
26//
27// MHDisp
28//
29// Create a false source plot using disp.
30//
31// Currently the use of this class requires to be after MFMagicCuts
32// in the tasklist. Switching of the M3L cut in MFMagicCuts is recommended.
33//
34// Class Version 3:
35// ----------------
36// + Double_t fScaleMin; // [deg] Minimum circle for integration of off-data for scaling
37// + Double_t fScaleMax; // [deg] Maximum circle for integration of off-data for scaling
38//
39//////////////////////////////////////////////////////////////////////////////
40#include "MHDisp.h"
41
42#include <TStyle.h>
43#include <TCanvas.h>
44
45#include <TF1.h>
46#include <TF2.h>
47#include <TProfile.h>
48
49#include "MLog.h"
50#include "MLogManip.h"
51
52#include "MMath.h"
53#include "MBinning.h"
54
55#include "MParList.h"
56#include "MParameters.h"
57
58#include "MHillas.h"
59#include "MSrcPosCam.h"
60#include "MPointingPos.h"
61#include "MPointingDev.h"
62
63ClassImp(MHDisp);
64
65using namespace std;
66
67// --------------------------------------------------------------------------
68//
69// Default Constructor
70//
71MHDisp::MHDisp(const char *name, const char *title)
72 : fDisp(0), fDeviation(0), fSmearing(-1), fWobble(kFALSE),
73 fScaleMin(0.325), fScaleMax(0.475)
74{
75 //
76 // set the name and title of this object
77 //
78 fName = name ? name : "MHDisp";
79 fTitle = title ? title : "3D-plot using Disp vs x, y";
80
81 fHist.SetName("Alpha");
82 fHist.SetTitle("3D-plot of ThetaSq vs x, y");
83 fHist.SetXTitle("x [\\circ]");
84 fHist.SetYTitle("y [\\circ]");
85 fHist.SetZTitle("Eq.cts");
86
87 fHist.SetDirectory(NULL);
88 fHistBg.SetDirectory(NULL);
89 fHistBg1.SetDirectory(NULL);
90 fHistBg2.SetDirectory(NULL);
91
92 fHist.SetBit(TH1::kNoStats);
93}
94
95// --------------------------------------------------------------------------
96//
97// Set binnings (takes BinningFalseSource) and prepare filling of the
98// histogram.
99//
100// Also search for MTime, MObservatory and MPointingPos
101//
102Bool_t MHDisp::SetupFill(const MParList *plist)
103{
104 if (!MHFalseSource::SetupFill(plist))
105 return kFALSE;
106
107 fDisp = (MParameterD*)plist->FindObject("Disp", "MParameterD");
108 if (!fDisp)
109 {
110 *fLog << err << "Disp [MParameterD] not found... abort." << endl;
111 return kFALSE;
112 }
113
114 fDeviation = (MPointingDev*)plist->FindObject("MPointingDev");
115 if (!fDeviation)
116 *fLog << warn << "MPointingDev not found... ignored." << endl;
117
118 MBinning binsx, binsy;
119 binsx.SetEdges(fHist, 'x');
120 binsy.SetEdges(fHist, 'y');
121 MH::SetBinning(&fHistBg, &binsx, &binsy);
122
123 if (!fHistOff)
124 {
125 MH::SetBinning(&fHistBg1, &binsx, &binsy);
126 MH::SetBinning(&fHistBg2, &binsx, &binsy);
127 }
128
129 return kTRUE;
130}
131
132// --------------------------------------------------------------------------
133//
134// Calculate the delta angle between fSrcPos->GetXY() and v.
135// Return result in deg.
136//
137Double_t MHDisp::DeltaPhiSrc(const TVector2 &v) const
138{
139 return TMath::Abs(fSrcPos->GetXY().DeltaPhi(v))*TMath::RadToDeg();
140}
141
142// --------------------------------------------------------------------------
143//
144// Fill the histogram. For details see the code or the class description
145//
146Bool_t MHDisp::Fill(const MParContainer *par, const Stat_t w)
147{
148 const MHillas *hil = dynamic_cast<const MHillas*>(par);
149 if (!hil)
150 {
151 *fLog << err << "MHDisp::Fill: No container specified!" << endl;
152 return kFALSE;
153 }
154
155 // Get camera rotation angle
156 Double_t rho = 0;
157 if (fTime && fObservatory && fPointPos)
158 rho = fPointPos->RotationAngle(*fObservatory, *fTime, fDeviation);
159
160 // Vector of length disp in direction of shower
161 // Move origin of vector to center-of-gravity of shower and derotate
162 TVector2 pos1 = hil->GetMean()*fMm2Deg + hil->GetNormAxis()*fDisp->GetVal();
163
164 Double_t w0 = 1;
165 if (fWobble)
166 {
167 const Double_t delta = DeltaPhiSrc(pos1);
168
169 // Skip off-data not in the same half than the source (here: anti-source)
170 // Could be replaced by some kind of theta cut...
171 if (!fHistOff)
172 {
173 if (delta>180-25)
174 return kTRUE;
175
176 // Because afterwards the plots are normalized with the number
177 // of entries the non-overlap region corresponds to half the
178 // contents, the overlap region to full contents. By taking
179 // the mean of both distributions we get the same result than
180 // we would have gotten using full off-events with a slightly
181 // increased uncertainty
182 // FIXME: The delta stuff could be replaced by a 2*antitheta cut...
183 w0 = delta>25 ? 1 : 2;
184 }
185
186 // Define by the source position which histogram to fill
187 if (DeltaPhiSrc(fFormerSrc)>90)
188 fHalf = !fHalf;
189 fFormerSrc = fSrcPos->GetXY();
190 }
191
192 // If on-data: Derotate source and P. Translate P to center.
193 TVector2 m;
194 if (fHistOff)
195 {
196 // Derotate all position around camera center by -rho
197 // Move origin of vector to center-of-gravity of shower and derotate
198 pos1=pos1.Rotate(-rho);
199
200 // Shift the source position to 0/0
201 if (fSrcPos)
202 {
203 // m: Position of the camera center in the FS plot
204 m = fSrcPos->GetXY().Rotate(-rho)*fMm2Deg;
205 pos1 -= m;
206 }
207 }
208
209 // -------------------------------------------------
210 // The following algorithm may look complicated...
211 // It is optimized for speed
212 // -------------------------------------------------
213
214 // Define a vector used to calculate rotated x and y component
215 const TVector2 rot(TMath::Sin(rho), TMath::Cos(rho));
216
217 // Fold event position with the PSF and fill it into histogram
218 const TAxis &axex = *fHist.GetXaxis();
219 const TAxis &axey = *fHist.GetYaxis();
220
221 const Int_t nx = axex.GetNbins();
222 const Int_t ny = axey.GetNbins();
223
224 // normg: Normalization for Gauss [sqrt(2*pi)*sigma]^2
225 const Double_t weight = axex.GetBinWidth(1)*axey.GetBinWidth(1)*w*w0;
226 const Double_t psf = 2*fSmearing*fSmearing;
227 const Double_t pi2 = fSmearing * TMath::Pi()/2;
228 const Double_t normg = pi2*pi2 * TMath::Sqrt(TMath::TwoPi()) / weight;
229 const Int_t bz = fHist.GetZaxis()->FindFixBin(0);
230
231 TH2 &bg = fHalf ? fHistBg1 : fHistBg2;
232
233 const Bool_t smear = fSmearing>0;
234 if (!smear)
235 {
236 if (!fHistOff)
237 bg.Fill(pos1.X(), pos1.Y(), w*w0);
238 else
239 fHist.Fill(pos1.X(), pos1.Y(), 0.0, w*w0);
240 }
241
242 // To calculate significance map smear with 2*theta-cut and
243 // do not normalize gauss, for event map use theta-cut/2 instead
244 if (smear || fHistOff)
245 {
246 for (int x=1; x<=nx; x++)
247 {
248 const Double_t cx = axex.GetBinCenter(x);
249 const Double_t px = cx-pos1.X();
250
251 for (int y=1; y<=ny; y++)
252 {
253 const Double_t cy = axey.GetBinCenter(y);
254 const Double_t sp = Sq(px, cy-pos1.Y());
255 if (smear)
256 {
257 const Double_t dp = sp/psf;
258
259 // Values below 1e-3 (>3.5sigma) are not filled into the histogram
260 if (dp<4)
261 {
262 const Double_t rc = TMath::Exp(-dp)/normg;
263 if (!fHistOff)
264 bg.AddBinContent(bg.GetBin(x, y), rc);
265 else
266 fHist.AddBinContent(fHist.GetBin(x, y, bz), rc);
267 }
268 }
269
270 if (!fHistOff)
271 continue;
272
273 // If we are filling the signal already (fHistOff!=NULL)
274 // we also fill the background by projecting the
275 // background in the camera into the sky plot.
276 const TVector2 v = TVector2(cx+m.X(), cy+m.Y());
277
278 // Speed up further: Xmax-fWobble
279 if (v.Mod()>axex.GetXmax()) // Gains 10% speed
280 continue;
281
282 const Int_t bx = axex.FindFixBin(v^rot);
283 const Int_t by = axey.FindFixBin(v*rot);
284 const Double_t bg = fHistOff->GetBinContent(bx, by, bz);
285
286 fHistBg.AddBinContent(fHistBg.GetBin(x, y), bg);
287 }
288 }
289 }
290
291 if (fHistOff)
292 fHistBg.SetEntries(fHistBg.GetEntries()+1);
293
294 if (!smear)
295 return kTRUE;
296
297 if (!fHistOff)
298 bg.SetEntries(bg.GetEntries()+1);
299 else
300 fHist.SetEntries(fHist.GetEntries()+1);
301
302 return kTRUE;
303}
304
305// --------------------------------------------------------------------------
306//
307// Compile the background in camera coordinates from the two background
308// histograms
309//
310Bool_t MHDisp::Finalize()
311{
312 if (fHistOff)
313 return kTRUE;
314
315 const Int_t bz = fHist.GetZaxis()->FindFixBin(0);
316
317 const Double_t n1 = fHistBg1.GetEntries();
318 const Double_t n2 = fHistBg2.GetEntries();
319
320 for (int x=1; x<=fHist.GetNbinsX(); x++)
321 for (int y=1; y<=fHist.GetNbinsY(); y++)
322 {
323 const Int_t bin1 = fHistBg1.GetBin(x, y);
324
325 const Double_t rc =
326 (n1==0?0:fHistBg1.GetBinContent(bin1)/n1)+
327 (n2==0?0:fHistBg2.GetBinContent(bin1)/n2);
328
329 fHist.SetBinContent(x, y, bz, rc/2);
330 }
331
332 fHist.SetEntries(n1+n2);
333
334 // Result corresponds to one smeared background event
335
336 return kTRUE;
337}
338
339// --------------------------------------------------------------------------
340//
341// Make sure that if the scale is changed by the context menu all subpads
342// are redrawn.
343//
344void MHDisp::Update()
345{
346 TVirtualPad *pad = gPad;
347 for (int i=1; i<=6; i++)
348 {
349 if (pad->GetPad(i))
350 pad->GetPad(i)->Modified();
351 }
352}
353
354// --------------------------------------------------------------------------
355//
356// Return the mean signal in h around (0,0) in a distance between
357// 0.325 and 0.475deg
358//
359Double_t MHDisp::GetOffSignal(TH1 &h) const
360{
361 const TAxis &axex = *h.GetXaxis();
362 const TAxis &axey = *h.GetYaxis();
363
364 Int_t cnt = 0;
365 Double_t sum = 0;
366 for (int x=0; x<h.GetNbinsX(); x++)
367 for (int y=0; y<h.GetNbinsY(); y++)
368 {
369 const Double_t d = TMath::Hypot(axex.GetBinCenter(x+1), axey.GetBinCenter(y+1));
370 if (d>fScaleMin && d<fScaleMax)
371 {
372 sum += h.GetBinContent(x+1,y+1);
373 cnt++;
374 }
375 }
376
377 return sum/cnt;
378}
379
380// --------------------------------------------------------------------------
381//
382// Fill h2 with the radial profile of h1 around (x0, y0)
383//
384void MHDisp::Profile(TH1 &h2, const TH2 &h1, Axis_t x0, Axis_t y0) const
385{
386 const TAxis &axex = *h1.GetXaxis();
387 const TAxis &axey = *h1.GetYaxis();
388
389 h2.Reset();
390
391 for (int x=1; x<=axex.GetNbins(); x++)
392 for (int y=1; y<=axey.GetNbins(); y++)
393 {
394 const Double_t dx = axex.GetBinCenter(x)-x0;
395 const Double_t dy = axey.GetBinCenter(y)-y0;
396
397 const Double_t r = TMath::Hypot(dx, dy);
398
399 h2.Fill(r, h1.GetBinContent(x, y));
400 }
401}
402
403// --------------------------------------------------------------------------
404//
405// Remove contents of histogram around a circle.
406//
407void MHDisp::MakeDot(TH2 &h2) const
408{
409 const TAxis &axex = *h2.GetXaxis();
410 const TAxis &axey = *h2.GetYaxis();
411
412 const Double_t rmax = (fWobble ? axex.GetXmax()-0.4 : axex.GetXmax()) - axex.GetBinWidth(1);
413
414 for (int x=1; x<=axex.GetNbins(); x++)
415 for (int y=1; y<=axey.GetNbins(); y++)
416 {
417 const Int_t bin = h2.GetBin(x,y);
418
419 const Double_t px = h2.GetBinCenter(x);
420 const Double_t py = h2.GetBinCenter(y);
421
422 if (rmax<TMath::Hypot(px, py))
423 h2.SetBinContent(bin, 0);
424 //else
425 // if (h2.GetBinContent(bin)==0)
426 // h2.SetBinContent(bin, 1e-10);
427 }
428}
429
430// --------------------------------------------------------------------------
431//
432// Calculate from signal and background the significance map
433//
434void MHDisp::MakeSignificance(TH2 &s, const TH2 &h1, const TH2 &h2, const Double_t scale) const
435{
436 const TAxis &axex = *s.GetXaxis();
437 const TAxis &axey = *s.GetYaxis();
438
439 const Int_t nx = axex.GetNbins();
440 const Int_t ny = axey.GetNbins();
441
442 const Int_t n = TMath::Nint(0.2/axex.GetBinWidth(1));
443
444 const Double_t sc = h2.GetEntries()/fHistOff->GetEntries();
445
446 for (int x=1; x<=nx; x++)
447 for (int y=1; y<=ny; y++)
448 {
449 Double_t S=0;
450
451 // Only calculate significances for pixels far enough
452 // from the border to get all expected pixels.
453 if (TMath::Hypot((float)x-0.5*nx, (float)y-0.5*ny)<0.5*axex.GetNbins()-n)
454 {
455 Double_t sig=0;
456 Double_t bg =0;
457
458 // Integral a region of n pixels around x/y
459 for (int dx=-n; dx<=n; dx++)
460 for (int dy=-n; dy<=n; dy++)
461 {
462 if (TMath::Hypot((float)dx, (float)dy)>n)
463 continue;
464
465 const Int_t bin = s.GetBin(x+dx,y+dy);
466
467 sig += h1.GetBinContent(bin);
468 bg += h2.GetBinContent(bin);
469 }
470
471 // Scale such, that the statistical error corresponds to
472 // the amount of off-data used to determin the background
473 // model, not to the background itself. This gives
474 // significances as calculated from the theta-sq plot.
475 S = sig>0 ? MMath::SignificanceLiMaSigned(sig, bg*scale/sc, sc) : 0;
476 }
477
478 s.SetBinContent(x, y, S);
479 }
480}
481
482void MHDisp::Profile1D(const char *name, const TH2 &h) const
483{
484 if (!gPad)
485 return;
486
487 TH1D *hrc = dynamic_cast<TH1D*>(gPad->FindObject(name));
488 if (!hrc)
489 return;
490
491 hrc->Reset();
492
493 //const_cast<TH2&>(h).SetMaximum();
494 const Double_t max = h.GetMaximum();
495
496 MBinning(51, -max*1.1, max*1.1).Apply(*hrc);
497
498 for (int x=1; x<=h.GetXaxis()->GetNbins(); x++)
499 for (int y=1; y<=h.GetXaxis()->GetNbins(); y++)
500 {
501 const Int_t bin = h.GetBin(x,y);
502 const Double_t sig = h.GetBinContent(bin);
503 if (sig!=0)
504 hrc->Fill(sig);
505 }
506
507 gPad->SetLogy(hrc->GetMaximum()>0);
508
509 if (!fHistOff || hrc->GetRMS()<0.1)
510 return;
511
512 // ---------- Fix mean ----------
513 // TF1 *g = (TF1*)gROOT->GetFunction("gaus");
514 // g->FixParameter(1, 0);
515 // hrc->Fit("gaus", "BQI");
516
517 hrc->Fit("gaus", "QI");
518
519 TF1 *f = hrc->GetFunction("gaus");
520 if (f)
521 {
522 f->SetLineWidth(1);
523 f->SetLineColor(kBlue);
524 }
525}
526
527void MHDisp::Paint(Option_t *o)
528{
529 // Compile Background if necessary
530 Finalize();
531
532 // Paint result
533 TVirtualPad *pad = gPad;
534
535 pad->cd(1);
536
537 // Project on data onto yx-plane
538 fHist.GetZaxis()->SetRange(0,9999);
539 TH2 *h1=(TH2*)fHist.Project3D("yx_on");
540
541 // Set Glow-palette (PaintSurface doesn't allow more than 99 colors)
542 MH::SetPalette(fHistOff?"glowsym":"glow1", 99);
543 h1->SetContour(99);
544
545 TH2 *hx=0;
546
547 Double_t scale = 1;
548 if (fHistOff)
549 {
550 // Project off data onto yx-plane and subtract it from on-data
551 fHistOff->GetZaxis()->SetRange(0,9999);
552 TH2 *h=(TH2*)fHistOff->Project3D("yx_off");
553
554 const Double_t h1off = GetOffSignal(*h1);
555 const Double_t hoff = GetOffSignal(fHistBg);
556
557 scale = hoff==0 ? -1 : -h1off/hoff;
558
559 hx = (TH2*)pad->GetPad(4)->FindObject("Alpha_yx_sig");
560 if (hx)
561 {
562 hx->SetContour(99);
563 MakeSignificance(*hx, *h1, fHistBg, TMath::Abs(scale));
564 MakeDot(*hx);
565 MakeSymmetric(hx);
566 }
567
568 h1->Add(&fHistBg, scale);
569
570 MakeDot(*h1);
571 MakeSymmetric(h1);
572
573 delete h;
574 }
575
576 pad->cd(3);
577 TH1 *h2 = (TH1*)gPad->FindObject("RadProf");
578
579 TString opt(o);
580 opt.ToLower();
581
582 if (h1 && h2)
583 {
584 Int_t ix, iy, iz;
585 h1->GetMaximumBin(ix, iy, iz);
586
587 const Double_t x0 = h1->GetXaxis()->GetBinCenter(ix);
588 const Double_t y0 = h1->GetYaxis()->GetBinCenter(iy);
589 //const Double_t w0 = h1->GetXaxis()->GetBinWidth(1);
590
591 Profile(*h2, *h1, 0, 0);
592 //Profile(*h2, *h1, x0, y0);
593
594 if (h2->GetEntries()>0)
595 {
596 // Replace with MAlphaFitter?
597 TF1 func("fcn", "gaus + [3]*x*x + [4]");
598 func.SetLineWidth(1);
599 func.SetLineColor(kBlue);
600
601 func.SetParLimits(2, h2->GetBinWidth(1), 1.0);
602
603 func.SetParameter(0, h2->GetBinContent(1));
604 func.FixParameter(1, 0);
605 func.SetParameter(2, 0.12);
606 if (fHistOff)
607 func.FixParameter(3, 0);
608 func.SetParameter(4, 0);//h2->GetBinContent(20));
609 h2->Fit(&func, "IQ", "", 0, 1.0);
610
611 h2->SetTitle(Form("P=(%.2f\\circ/%.2f\\circ) \\omega=%.2f\\circ f=%.2f", x0, y0, func.GetParameter(2), TMath::Abs(scale)));
612 }
613 }
614
615 pad->cd(5);
616 if (h1)
617 Profile1D("Exc", *h1);
618
619 pad->cd(6);
620 if (hx)
621 Profile1D("Sig", *hx);
622}
623
624void MHDisp::Draw(Option_t *o)
625{
626 TVirtualPad *pad = gPad ? gPad : MakeDefCanvas(this);
627 const Int_t col = pad->GetFillColor();
628 pad->SetBorderMode(0);
629
630 AppendPad(o);
631
632 // ----- Pad number 4 -----
633 TString name = Form("%s_4", pad->GetName());
634 TPad *p = new TPad(name,name, 0.525/*0.5025*/, 0.3355, 0.995, 0.995, col, 0, 0);
635 p->SetNumber(4);
636 p->Draw();
637 p->cd();
638
639 TH1 *h3 = fHist.Project3D("yx_sig");
640 h3->SetTitle("Significance Map");
641 h3->SetDirectory(NULL);
642 h3->SetXTitle(fHist.GetXaxis()->GetTitle());
643 h3->SetYTitle(fHist.GetYaxis()->GetTitle());
644 h3->SetMinimum(0);
645 h3->Draw("colz");
646 h3->SetBit(kCanDelete);
647
648 if (fHistOff)
649 GetCatalog()->Draw("mirror same *");
650
651 // ----- Pad number 1 -----
652 pad->cd();
653 name = Form("%s_1", pad->GetName());
654 p = new TPad(name,name, 0.005, 0.3355, 0.475/*0.4975*/, 0.995, col, 0, 0);
655 p->SetNumber(1);
656 p->Draw();
657 p->cd();
658
659 h3 = fHist.Project3D("yx_on");
660 h3->SetTitle("Distribution of equivalent events");
661 h3->SetDirectory(NULL);
662 h3->SetXTitle(fHist.GetXaxis()->GetTitle());
663 h3->SetYTitle(fHist.GetYaxis()->GetTitle());
664 h3->SetMinimum(0);
665 h3->Draw("colz");
666 h3->SetBit(kCanDelete);
667
668 if (fHistOff)
669 GetCatalog()->Draw("mirror same *");
670
671 // ----- Pad number 2 -----
672 pad->cd();
673 name = Form("%s_2", pad->GetName());
674 p = new TPad(name,name, 0.005, 0.005, 0.2485, 0.3315, col, 0, 0);
675 p->SetNumber(2);
676 p->Draw();
677 p->cd();
678 h3->Draw("surf3");
679
680 // ----- Pad number 3 -----
681 pad->cd();
682 name = Form("%s_3", pad->GetName());
683 p = new TPad(name,name, 0.2525, 0.005, 0.4985, 0.3315, col, 0, 0);
684 p->SetNumber(3);
685 p->SetGrid();
686 p->Draw();
687 p->cd();
688
689 const Double_t maxr = TMath::Hypot(h3->GetXaxis()->GetXmax(), h3->GetYaxis()->GetXmax());
690 const Int_t nbin = (h3->GetNbinsX()+h3->GetNbinsY())/2;
691 TProfile *h = new TProfile("RadProf", "Radial Profile", nbin, 0, maxr);
692 h->SetDirectory(0);
693 //TH1F *h = new TH1F("RadProf", "Radial Profile", nbin, 0, maxr);
694 //h->Sumw2();
695 h->SetXTitle("\\vartheta [\\circ]");
696 h->SetYTitle("<cts>/\\Delta R");
697 h->SetBit(kCanDelete);
698 h->Draw();
699
700 // ----- Pad number 5 -----
701 pad->cd();
702 name = Form("%s_5", pad->GetName());
703 p = new TPad(name,name, 0.5025, 0.005, 0.7485, 0.3315, col, 0, 0);
704 p->SetNumber(5);
705 p->SetGrid();
706 p->Draw();
707 p->cd();
708
709 h3 = new TH1D("Exc", "Distribution of excess events", 10, -1, 1);
710 h3->SetDirectory(NULL);
711 h3->SetXTitle("N");
712 h3->SetYTitle("Counts");
713 h3->Draw();
714 h3->SetBit(kCanDelete);
715
716 // ----- Pad number 6 -----
717 pad->cd();
718 name = Form("%s_6", pad->GetName());
719 p = new TPad(name,name, 0.7525, 0.005, 0.9995, 0.3315, col, 0, 0);
720 p->SetNumber(6);
721 p->SetGrid();
722 p->Draw();
723 p->cd();
724
725 h3 = new TH1D("Sig", "Distribution of significances", 10, -1, 1);
726 h3->SetDirectory(NULL);
727 h3->SetXTitle("N");
728 h3->SetYTitle("Counts");
729 h3->Draw();
730 h3->SetBit(kCanDelete);
731}
732
733// --------------------------------------------------------------------------
734//
735// The following resources are available:
736//
737// MHDisp.Smearing: 0.1
738// MHDisp.Wobble: on/off
739//
740Int_t MHDisp::ReadEnv(const TEnv &env, TString prefix, Bool_t print)
741{
742 Int_t rc = MHFalseSource::ReadEnv(env, prefix, print);
743 if (rc==kERROR)
744 return kERROR;
745
746 if (IsEnvDefined(env, prefix, "Smearing", print))
747 {
748 rc = kTRUE;
749 SetSmearing(GetEnvValue(env, prefix, "Smearing", fSmearing));
750 }
751 if (IsEnvDefined(env, prefix, "Wobble", print))
752 {
753 rc = kTRUE;
754 SetWobble(GetEnvValue(env, prefix, "Wobble", fWobble));
755 }
756 if (IsEnvDefined(env, prefix, "ScaleMin", print))
757 {
758 rc = kTRUE;
759 SetScaleMin(GetEnvValue(env, prefix, "ScaleMin", fScaleMin));
760 }
761 if (IsEnvDefined(env, prefix, "ScaleMax", print))
762 {
763 rc = kTRUE;
764 SetScaleMax(GetEnvValue(env, prefix, "ScaleMax", fScaleMax));
765 }
766
767 return rc;
768}
769
Note: See TracBrowser for help on using the repository browser.