CbmRoot
PairAnalysisSignalFit.cxx
Go to the documentation of this file.
1 // Dielectron SignalFit //
3 // //
4 // //
5 /*
6 
7  Class used for extracting the signal from an invariant mass spectrum.
8  It implements the PairAnalysisSignalBase and -Ext classes and it uses user provided
9  functions to fit the spectrum with a combined signa+background fit.
10  Used invariant mass spectra are provided via an array of histograms. There are serveral method
11  to estimate the background and to extract the raw yield from the background subtracted spectra.
12 
13  Example usage:
14 
15  PairAnalysisSignalFit *sig = new PairAnalysisSignalFit();
16 
17 
18  1) invariant mass input spectra
19 
20  1.1) Assuming a PairAnalysisCF container as data format (check class for more details)
21  PairAnalysisCFdraw *cf = new PairAnalysisCFdraw("path/to/the/output/file.root");
22  TObjArray *arrHists = cf->CollectMinvProj(cf->FindStep("Config"));
23 
24  1.2) Assuming a PairAnalysisHF grid as data format (check class for more details)
25  PairAnalysisHFhelper *hf = new PairAnalysisHFhelper("path/to/the/output/file.root", "ConfigName");
26  TObjArray *arrHists = hf->CollectHistos(PairAnalysisVarManager::kM);
27 
28  1.3) Assuming a single histograms
29  TObjArray *histoArray = new TObjArray();
30  arrHists->Add(signalPP); // add the spectrum histograms to the array
31  arrHists->Add(signalPM); // the order is important !!!
32  arrHists->Add(signalMM);
33 
34 
35  2) background estimation
36 
37  2.1) set the method for the background estimation (methods can be found in PairAnalysisSignalBase)
38  sig->SetMethod(PairAnalysisSignalBase::kFitted);
39  2.2) rebin the spectras if needed
40  // sig->SetRebin(2);
41  2.3) add any background function you like
42  TF1 *fB = new TF1("fitBgrd","pol3",minFit,maxFit);
43 
44 
45  3) configure the signal extraction
46 
47  3.1) chose one of the signal functions (MCshape, CrystalBall, Gauss)
48  TF1 *fS = new TF1("fitSign",PairAnalysisFunction::PeakFunCB,minFit,maxFit,5); // has 5 parameters
49  // TF1 *fS = new TF1("fitSign",PairAnalysisFunction::PeakFunGaus,minFit,maxFit,3); // has 3 parameters
50  // sig->SetMCSignalShape(hMCsign);
51  // TF1 *fS = new TF1("fitSign",PairAnalysisFunction::PeakFunMC,minFit,maxFit,1); // requires a MC shape
52  3.2) set the method for the signal extraction (methods can be found in PairAnalysisSignalBase)
53  depending on the method serveral inputs are needed (e.g. MC shape, PDG code of the particle of interest)
54  // sig->SetParticleOfInterest(443); //default is jpsi
55  // sig->SetMCSignalShape(signalMC);
56  // sig->SetIntegralRange(minInt, maxInt);
57  sig->SetExtractionMethod(PairAnalysisSignal::BinCounting); // this is the default
58 
59 
60  4) combined fit of bgrd+signal
61 
62  4.1) combine the two functions
63  sig->CombineFunc(fS,fB);
64  4.2) apply fitting ranges and the fit options
65  sig->SetFitRange(minFit, maxFit);
66  sig->SetFitOption("NR");
67 
68 
69  5) start the processing
70 
71  sig->Process(arrHists);
72  sig->Print(""); // print values and errors extracted
73 
74 
75  6) access the spectra and values created
76 
77  6.1) standard spectra as provided filled in PairAnalysisSignalExt
78  TH1F *hsign = (TH1F*) sig->GetUnlikeSignHistogram(); // same as the input (rebinned)
79  TH1F *hbgrd = (TH1F*) sig->GetBackgroundHistogram(); // filled histogram with fitBgrd
80  TH1F *hextr = (TH1F*) sig->GetSignalHistogram(); // after backgound extraction (rebinned)
81  TObject *oPeak = (TObject*) sig->GetPeakShape(); // can be a TF1 or TH1 depending on the method
82  6.2) flow spectra
83  TF1 *fFitSign = sig->GetCombinedFunction(); // combined fit function
84  TF1 *fFitExtr = sig->GetSignalFunction(); // signal function
85  TF1 *fFitBgrd = sig->GetBackgroundFunction(); // background function
86  6.3) access the extracted values and errors
87  sig->GetValues(); or GetErrors(); // yield extraction
88 
89 */
90 // //
92 
93 #include <TF1.h>
94 #include <TFitResult.h>
95 #include <TGraph.h>
96 #include <TH1.h>
97 #include <TList.h>
98 #include <TMath.h>
99 #include <TPaveText.h>
100 #include <TString.h>
101 //#include <../hist/hist/src/TF1Helper.h>
102 
103 #include "PairAnalysisSignalFit.h"
104 
106 
107  //TH1F* PairAnalysisSignalFunc::fgHistSimPM=0x0;
108  // // TF1* PairAnalysisSignalFunc::fFuncSignal=0x0;
109  // // TF1* PairAnalysisSignalFunc::fFuncBackground=0x0;
110  // // Int_t PairAnalysisSignalFunc::fNparPeak=0;
111  // // Int_t PairAnalysisSignalFunc::fNparBgnd=0;
112 
113 
116 // ,PairAnalysisFunction()
117 {
118  //
119  // Default Constructor
120  //
121 }
122 
123 //______________________________________________
125  const char* title)
126  : PairAnalysisSignalExt(name, title) //,
127 // PairAnalysisFunction()
128 {
129  //
130  // Named Constructor
131  //
132 }
133 
134 //______________________________________________
136  //
137  // Default Destructor
138  //
139 }
140 
141 
142 //______________________________________________
143 void PairAnalysisSignalFit::Process(TObjArray* const arrhist) {
144  //
145  // Fit the invariant mass histograms and retrieve the signal and background
146  //
147  switch (fMethod) {
148  case kFitted: ProcessFit(arrhist); break;
149 
150  case kLikeSignFit: ProcessFitLS(arrhist); break;
151 
152  case kEventMixingFit: ProcessFitEM(arrhist); break;
153 
154  case kLikeSign:
155  case kLikeSignArithm:
156  case kLikeSignRcorr:
158  case kEventMixing:
159  case kRotation: PairAnalysisSignalExt::Process(arrhist); break;
160 
161  default: Error("Process", "Background substraction method not known!");
162  }
163 }
164 
165 //______________________________________________
166 void PairAnalysisSignalFit::ProcessFit(TObjArray* const arrhist) {
167  //
168  // Fit the +- invariant mass distribution only
169  // Here we assume that the combined fit function is a sum of the signal and background functions
170  // and that the signal function is always the first term of this sum
171  //
172 
173  fHistDataPM = (TH1F*) (arrhist->At(1))->Clone("histPM"); // +- SE
174  fHistDataPM->Sumw2();
175  if (fRebin > 1) fHistDataPM->Rebin(fRebin);
176 
177  fHistSignal = new TH1F("HistSignal",
178  "Fit substracted signal",
179  fHistDataPM->GetXaxis()->GetNbins(),
180  fHistDataPM->GetXaxis()->GetXmin(),
181  fHistDataPM->GetXaxis()->GetXmax());
182  fHistBackground = new TH1F("HistBackground",
183  "Fit contribution",
184  fHistDataPM->GetXaxis()->GetNbins(),
185  fHistDataPM->GetXaxis()->GetXmin(),
186  fHistDataPM->GetXaxis()->GetXmax());
187 
188  // the starting parameters of the fit function and their limits can be tuned
189  // by the user in its macro
190  fHistDataPM->Fit(fFuncSigBack, fFitOpt.Data(), "", fFitMin, fFitMax);
191  TFitResultPtr pmFitPtr =
192  fHistDataPM->Fit(fFuncSigBack, fFitOpt.Data(), "", fFitMin, fFitMax);
193  //TFitResult *pmFitResult = pmFitPtr.Get(); // used only with TF1Helper
194  //fFuncBackground->SetParameters(fFuncSigBack->GetParameters());
195  fFuncSignal->SetParameters(fFuncSigBack->GetParameters());
196  fFuncBackground->SetParameters(fFuncSigBack->GetParameters()
197  + fFuncSignal->GetNpar());
198 
199  // fill the background spectrum
201  // set the error for the background histogram
202  fHistBackground->Fit(fFuncBackground, "0qR", "", fFitMin, fFitMax);
203  Double_t inte = fFuncBackground->IntegralError(fIntMin, fIntMax)
204  / fHistDataPM->GetBinWidth(1);
205  Double_t binte =
206  inte
207  / TMath::Sqrt(
208  (fHistDataPM->FindBin(fIntMax) - fHistDataPM->FindBin(fIntMin)) + 1);
209  for (Int_t iBin = fHistDataPM->FindBin(fIntMin);
210  iBin <= fHistDataPM->FindBin(fIntMax);
211  iBin++) {
212  fHistBackground->SetBinError(iBin, binte);
213  }
214 
215  for (Int_t iBin = 1; iBin <= fHistDataPM->GetXaxis()->GetNbins(); iBin++) {
216  // Double_t m = fHistDataPM->GetBinCenter(iBin);
217  Double_t pm = fHistDataPM->GetBinContent(iBin);
218  Double_t epm = fHistDataPM->GetBinError(iBin);
219  Double_t bknd = fHistBackground->GetBinContent(iBin);
220  Double_t ebknd = fHistBackground->GetBinError(iBin);
221  Double_t signal = pm - bknd;
222  Double_t error = TMath::Sqrt(epm * epm + ebknd);
223  // theres no signal extraction outside the fit region
224  if (fHistDataPM->GetBinLowEdge(iBin) > fFitMax
225  || fHistDataPM->GetBinLowEdge(iBin) + fHistDataPM->GetBinWidth(iBin)
226  < fFitMin) {
227  signal = 0.0;
228  error = 0.0;
229  }
230  fHistSignal->SetBinContent(iBin, signal);
231  fHistSignal->SetBinError(iBin, error);
232  }
233 
234  if (fUseIntegral) {
235  // signal
236  fValues(0) =
237  fFuncSignal->Integral(fIntMin, fIntMax) / fHistDataPM->GetBinWidth(1);
238  fErrors(0) = fFuncSignal->IntegralError(fIntMin, fIntMax)
239  / fHistDataPM->GetBinWidth(1);
240  // fErrors(0) = 0;
241  // background
242  fValues(1) =
243  fFuncBackground->Integral(fIntMin, fIntMax) / fHistDataPM->GetBinWidth(1);
244  fErrors(1) = fFuncBackground->IntegralError(fIntMin, fIntMax)
245  / fHistDataPM->GetBinWidth(1);
246  } else {
247  // signal
248  fValues(0) = fHistSignal->IntegralAndError(
249  fHistSignal->FindBin(fIntMin), fHistSignal->FindBin(fIntMax), fErrors(0));
250  // background
251  fValues(1) =
252  fHistBackground->IntegralAndError(fHistBackground->FindBin(fIntMin),
253  fHistBackground->FindBin(fIntMax),
254  fErrors(1));
255  }
256  // S/B and significance
258  fValues(4) = fFuncSigBack->GetParameter(fParMass);
259  fErrors(4) = fFuncSigBack->GetParError(fParMass);
260  fValues(5) = fFuncSigBack->GetParameter(fParMassWidth);
261  fErrors(5) = fFuncSigBack->GetParError(fParMassWidth);
262 
263  // flag
264  fProcessed = kTRUE;
265 
266  fHistBackground->GetListOfFunctions()->Add(fFuncBackground);
267 }
268 
269 //______________________________________________
270 void PairAnalysisSignalFit::ProcessFitLS(TObjArray* const arrhist) {
271  //
272  // Substract background using the like-sign spectrum
273  //
274  fHistDataPP = (TH1F*) (arrhist->At(0))->Clone("histPP"); // ++ SE
275  fHistDataPM = (TH1F*) (arrhist->At(1))->Clone("histPM"); // +- SE
276  fHistDataMM = (TH1F*) (arrhist->At(2))->Clone("histMM"); // -- SE
277  if (fRebin > 1) {
278  fHistDataPP->Rebin(fRebin);
279  fHistDataPM->Rebin(fRebin);
280  fHistDataMM->Rebin(fRebin);
281  }
282  fHistDataPP->Sumw2();
283  fHistDataPM->Sumw2();
284  fHistDataMM->Sumw2();
285 
286  fHistSignal = new TH1F("HistSignal",
287  "Like-Sign substracted signal",
288  fHistDataPM->GetXaxis()->GetNbins(),
289  fHistDataPM->GetXaxis()->GetXmin(),
290  fHistDataPM->GetXaxis()->GetXmax());
291  fHistBackground = new TH1F("HistBackground",
292  "Like-sign contribution",
293  fHistDataPM->GetXaxis()->GetNbins(),
294  fHistDataPM->GetXaxis()->GetXmin(),
295  fHistDataPM->GetXaxis()->GetXmax());
296 
297  // fit the +- mass distribution
298  fHistDataPM->Fit(fFuncSigBack, fFitOpt.Data(), "", fFitMin, fFitMax);
299  fHistDataPM->Fit(fFuncSigBack, fFitOpt.Data(), "", fFitMin, fFitMax);
300  // declare the variables where the like-sign fit results will be stored
301  // TFitResult *ppFitResult = 0x0;
302  // TFitResult *mmFitResult = 0x0;
303  // fit the like sign background
304  TF1* funcClonePP = (TF1*) fFuncBackground->Clone("funcClonePP");
305  TF1* funcCloneMM = (TF1*) fFuncBackground->Clone("funcCloneMM");
306  fHistDataPP->Fit(funcClonePP, fFitOpt.Data(), "", fFitMin, fFitMax);
307  // TFitResultPtr ppFitPtr = fHistDataPP->Fit(funcClonePP, fFitOpt.Data(), "", fFitMin, fFitMax);
308  // ppFitResult = ppFitPtr.Get();
309  fHistDataMM->Fit(funcCloneMM, fFitOpt.Data(), "", fFitMin, fFitMax);
310  // TFitResultPtr mmFitPtr = fHistDataMM->Fit(funcCloneMM, fFitOpt.Data(), "", fFitMin, fFitMax);
311  // mmFitResult = mmFitPtr.Get();
312 
313  for (Int_t iBin = 1; iBin <= fHistDataPM->GetXaxis()->GetNbins(); iBin++) {
314  Double_t m = fHistDataPM->GetBinCenter(iBin);
315  Double_t pm = fHistDataPM->GetBinContent(iBin);
316  Double_t pp = funcClonePP->Eval(m);
317  Double_t mm = funcCloneMM->Eval(m);
318  Double_t epm = fHistDataPM->GetBinError(iBin);
319  // TODO: use TFitResults for errors?
320  Double_t epp = 0;
321  Double_t emm = 0;
322 
323  Double_t signal = pm - 2.0 * TMath::Sqrt(pp * mm);
324  Double_t background = 2.0 * TMath::Sqrt(pp * mm);
325 
326  // error propagation on the signal calculation above
327  Double_t esignal =
328  TMath::Sqrt(epm * epm + (mm / pp) * epp + (pp / mm) * emm);
329  Double_t ebackground = TMath::Sqrt((mm / pp) * epp + (pp / mm) * emm);
330  fHistSignal->SetBinContent(iBin, signal);
331  fHistSignal->SetBinError(iBin, esignal);
332  fHistBackground->SetBinContent(iBin, background);
333  fHistBackground->SetBinError(iBin, ebackground);
334  }
335 
336  // signal
337  fValues(0) = fHistSignal->IntegralAndError(
338  fHistSignal->FindBin(fIntMin), fHistSignal->FindBin(fIntMax), fErrors(0));
339  // background
340  fValues(1) =
341  fHistBackground->IntegralAndError(fHistBackground->FindBin(fIntMin),
342  fHistBackground->FindBin(fIntMax),
343  fErrors(1));
344  // S/B and significance
346  fValues(4) = fFuncSigBack->GetParameter(fParMass);
347  fErrors(4) = fFuncSigBack->GetParError(fParMass);
348  fValues(5) = fFuncSigBack->GetParameter(fParMassWidth);
349  fErrors(5) = fFuncSigBack->GetParError(fParMassWidth);
350 
351  // flag
352  fProcessed = kTRUE;
353 }
354 
355 //______________________________________________
356 void PairAnalysisSignalFit::ProcessFitEM(TObjArray* const arrhist) {
357  //
358  // Substract background with the event mixing technique
359  //
360  arrhist->GetEntries(); // just to avoid the unused parameter warning
361  Error("ProcessFitEM",
362  "Event mixing for background substraction method not yet implemented!");
363 }
364 
365 //______________________________________________
366 void PairAnalysisSignalFit::Draw(const Option_t* option) {
367  //
368  // Draw the fitted function
369  //
370  // TODO: update
371  TString drawOpt(option);
372  drawOpt.ToLower();
373 
374  Bool_t optStat = drawOpt.Contains("stat");
375 
376  fFuncSigBack->SetNpx(200);
377  fFuncSigBack->SetRange(fIntMin, fIntMax);
378  fFuncBackground->SetNpx(200);
379  fFuncBackground->SetRange(fIntMin, fIntMax);
380 
381  TGraph* grSig = new TGraph(fFuncSigBack);
382  grSig->SetFillColor(kGreen);
383  grSig->SetFillStyle(3001);
384 
385  TGraph* grBack = new TGraph(fFuncBackground);
386  grBack->SetFillColor(kRed);
387  grBack->SetFillStyle(3001);
388 
389  grSig->SetPoint(0, grBack->GetX()[0], grBack->GetY()[0]);
390  grSig->SetPoint(grSig->GetN() - 1,
391  grBack->GetX()[grBack->GetN() - 1],
392  grBack->GetY()[grBack->GetN() - 1]);
393 
394  grBack->SetPoint(0, grBack->GetX()[0], 0.);
395  grBack->SetPoint(grBack->GetN() - 1, grBack->GetX()[grBack->GetN() - 1], 0.);
396 
397  fFuncSigBack->SetRange(fFitMin, fFitMax);
398  fFuncBackground->SetRange(fFitMin, fFitMax);
399 
400  if (!drawOpt.Contains("same")) {
401  if (fHistDataPM) {
402  fHistDataPM->Draw();
403  grSig->Draw("f");
404  } else {
405  grSig->Draw("af");
406  }
407  } else {
408  grSig->Draw("f");
409  }
410  if (fMethod == kFitted || fMethod == kFittedMC) grBack->Draw("f");
411  fFuncSigBack->Draw("same");
412  fFuncSigBack->SetLineWidth(2);
413  if (fMethod == kLikeSign) {
414  fHistDataPP->SetLineWidth(2);
415  fHistDataPP->SetLineColor(6);
416  fHistDataPP->Draw("same");
417  fHistDataMM->SetLineWidth(2);
418  fHistDataMM->SetLineColor(8);
419  fHistDataMM->Draw("same");
420  }
421 
422  if (fMethod == kFitted || fMethod == kFittedMC) fFuncBackground->Draw("same");
423 
424  if (optStat) DrawStats();
425 }
426 
427 //______________________________________________
428 void PairAnalysisSignalFit::Print(Option_t* /*option*/) const {
429  //
430  // Print the statistics
431  //
432  PairAnalysisSignalBase::Print();
433  printf("Fit int. : %.5g - %.5g \n", fFitMin, fFitMax);
434  printf("Fit chi/%d: %.5g \n", fDof, fChi2Dof);
435 }
PairAnalysisSignalExt::fIntMax
Double_t fIntMax
Definition: PairAnalysisSignalExt.h:236
PairAnalysisSignalExt::kFittedMC
@ kFittedMC
Definition: PairAnalysisSignalExt.h:29
PairAnalysisSignalExt::kLikeSignRcorr
@ kLikeSignRcorr
Definition: PairAnalysisSignalExt.h:33
PairAnalysisSignalExt::kLikeSignFit
@ kLikeSignFit
Definition: PairAnalysisSignalExt.h:35
PairAnalysisFunction::fFitOpt
TString fFitOpt
Definition: PairAnalysisFunction.h:116
PairAnalysisSignalExt::fErrors
TVectorD fErrors
Definition: PairAnalysisSignalExt.h:233
PairAnalysisSignalFit.h
PairAnalysisSignalFit::ProcessFitEM
void ProcessFitEM(TObjArray *const arrhist)
Definition: PairAnalysisSignalFit.cxx:356
PairAnalysisFunction::fDof
Int_t fDof
Definition: PairAnalysisFunction.h:120
PairAnalysisSignalFit
Definition: PairAnalysisSignalFit.h:20
PairAnalysisSignalExt
Definition: PairAnalysisSignalExt.h:25
PairAnalysisSignalExt::kRotation
@ kRotation
Definition: PairAnalysisSignalExt.h:38
PairAnalysisSignalExt::fHistDataPP
TH1 * fHistDataPP
Definition: PairAnalysisSignalExt.h:219
PairAnalysisSignalExt::fHistDataPM
TH1 * fHistDataPM
Definition: PairAnalysisSignalExt.h:218
PairAnalysisFunction::fFitMax
Double_t fFitMax
Definition: PairAnalysisFunction.h:107
PairAnalysisSignalExt::fMethod
EBackgroundMethod fMethod
Definition: PairAnalysisSignalExt.h:245
PairAnalysisSignalExt::kEventMixingFit
@ kEventMixingFit
Definition: PairAnalysisSignalExt.h:37
PairAnalysisFunction::fChi2Dof
Double_t fChi2Dof
Definition: PairAnalysisFunction.h:121
PairAnalysisSignalExt::kLikeSign
@ kLikeSign
Definition: PairAnalysisSignalExt.h:31
PairAnalysisFunction::fParMass
Int_t fParMass
Definition: PairAnalysisFunction.h:111
PairAnalysisSignalExt::fHistSignal
TH1 * fHistSignal
Definition: PairAnalysisSignalExt.h:212
PairAnalysisFunction::fFitMin
Double_t fFitMin
Definition: PairAnalysisFunction.h:106
PairAnalysisFunction::fFuncBackground
TF1 * fFuncBackground
Definition: PairAnalysisFunction.h:103
PairAnalysisFunction::fFuncSigBack
TF1 * fFuncSigBack
Definition: PairAnalysisFunction.h:104
PairAnalysisSignalFit::Draw
virtual void Draw(const Option_t *option="")
Definition: PairAnalysisSignalFit.cxx:366
PairAnalysisSignalExt::kEventMixing
@ kEventMixing
Definition: PairAnalysisSignalExt.h:36
PairAnalysisSignalExt::fHistDataMM
TH1 * fHistDataMM
Definition: PairAnalysisSignalExt.h:220
PairAnalysisSignalExt::kFitted
@ kFitted
Definition: PairAnalysisSignalExt.h:30
PairAnalysisSignalExt::fRebin
Int_t fRebin
Definition: PairAnalysisSignalExt.h:240
PairAnalysisFunction::fFuncSignal
TF1 * fFuncSignal
Definition: PairAnalysisFunction.h:102
PairAnalysisSignalExt::fProcessed
Bool_t fProcessed
Definition: PairAnalysisSignalExt.h:265
PairAnalysisSignalExt::fIntMin
Double_t fIntMin
Definition: PairAnalysisSignalExt.h:235
PairAnalysisSignalExt::fValues
TVectorD fValues
Definition: PairAnalysisSignalExt.h:232
PairAnalysisSignalExt::Process
void Process(TObjArray *const arrhist)
Definition: PairAnalysisSignalExt.cxx:596
PairAnalysisFunction::fUseIntegral
Bool_t fUseIntegral
Definition: PairAnalysisFunction.h:117
PairAnalysisSignalFit::Process
virtual void Process(TObjArray *const arrhist)
Definition: PairAnalysisSignalFit.cxx:143
PairAnalysisSignalExt::kLikeSignArithmRcorr
@ kLikeSignArithmRcorr
Definition: PairAnalysisSignalExt.h:34
PairAnalysisSignalFit::~PairAnalysisSignalFit
virtual ~PairAnalysisSignalFit()
Definition: PairAnalysisSignalFit.cxx:135
PairAnalysisSignalFit::Print
virtual void Print(Option_t *option="") const
Definition: PairAnalysisSignalFit.cxx:428
m
__m128 m
Definition: L1/vectors/P4_F32vec4.h:26
PairAnalysisSignalExt::kLikeSignArithm
@ kLikeSignArithm
Definition: PairAnalysisSignalExt.h:32
ClassImp
ClassImp(PairAnalysisSignalFit) PairAnalysisSignalFit
Definition: PairAnalysisSignalFit.cxx:105
PairAnalysisSignalFit::ProcessFitLS
void ProcessFitLS(TObjArray *const arrhist)
Definition: PairAnalysisSignalFit.cxx:270
PairAnalysisSignalExt::fHistBackground
TH1 * fHistBackground
Definition: PairAnalysisSignalExt.h:215
PairAnalysisFunction::fParMassWidth
Int_t fParMassWidth
Definition: PairAnalysisFunction.h:113
PairAnalysisSignalExt::SetSignificanceAndSOB
void SetSignificanceAndSOB()
Definition: PairAnalysisSignalExt.h:314
PairAnalysisSignalExt::DrawStats
TPaveText * DrawStats(Double_t x1=0., Double_t y1=0., Double_t x2=0., Double_t y2=0., TString opt="pRnbsSmrc")
Definition: PairAnalysisSignalExt.cxx:140
PairAnalysisSignalFit::ProcessFit
void ProcessFit(TObjArray *const arrhist)
Definition: PairAnalysisSignalFit.cxx:166
PairAnalysisSignalFit::PairAnalysisSignalFit
PairAnalysisSignalFit()