]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGDQ/dielectron/AliDielectronSignalFunc.cxx
Fix in LHC magnets initialization convention: for pA beams the provided energy
[u/mrichter/AliRoot.git] / PWGDQ / dielectron / AliDielectronSignalFunc.cxx
CommitLineData
572b0139 1/*************************************************************************
2* Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
3* *
4* Author: The ALICE Off-line Project. *
5* Contributors are mentioned in the code where appropriate. *
6* *
7* Permission to use, copy, modify and distribute this software and its *
8* documentation strictly for non-commercial purposes is hereby granted *
9* without fee, provided that the above copyright notice appears in all *
10* copies and that both the copyright notice and this permission notice *
11* appear in the supporting documentation. The authors make no claims *
12* about the suitability of this software for any purpose. It is *
13* provided "as is" without express or implied warranty. *
14**************************************************************************/
15
16///////////////////////////////////////////////////////////////////////////
17// Dielectron SignalFunc //
18// //
19// //
20/*
21Dielectron signal extraction class using functions as input.
22
23A function to describe the signal as well as one to describe the background
24has to be deployed by the user. Alternatively on of the default implementaions
25can be used.
26
27*/
28// //
29///////////////////////////////////////////////////////////////////////////
30
31#include <TF1.h>
32#include <TH1.h>
33#include <TGraph.h>
34#include <TMath.h>
35#include <TString.h>
36#include <TPaveText.h>
bc75eeb5 37#include <TList.h>
38#include <TFitResult.h>
cafdf10b 39//#include <../hist/hist/src/TF1Helper.h>
572b0139 40
41#include <AliLog.h>
42
43#include "AliDielectronSignalFunc.h"
44
45ClassImp(AliDielectronSignalFunc)
46
5720c765 47TH1F* AliDielectronSignalFunc::fgHistSimPM=0x0;
e43f22fb 48TF1* AliDielectronSignalFunc::fFuncSignal=0x0;
49TF1* AliDielectronSignalFunc::fFuncBackground=0x0;
50Int_t AliDielectronSignalFunc::fNparPeak=0;
51Int_t AliDielectronSignalFunc::fNparBgnd=0;
52
5720c765 53
572b0139 54AliDielectronSignalFunc::AliDielectronSignalFunc() :
3505bfad 55AliDielectronSignalBase(),
3505bfad 56fFuncSigBack(0x0),
57fParMass(1),
58fParMassWidth(2),
59fFitOpt("SMNQE"),
5720c765 60fUseIntegral(kFALSE),
61fPolDeg(0),
e43f22fb 62fDof(0),
5720c765 63fChi2Dof(0.0)
572b0139 64{
65 //
66 // Default Constructor
67 //
68
69}
70
71//______________________________________________
72AliDielectronSignalFunc::AliDielectronSignalFunc(const char* name, const char* title) :
3505bfad 73AliDielectronSignalBase(name, title),
3505bfad 74fFuncSigBack(0x0),
75fParMass(1),
76fParMassWidth(2),
77fFitOpt("SMNQE"),
5720c765 78fUseIntegral(kFALSE),
79fPolDeg(0),
e43f22fb 80fDof(0),
5720c765 81fChi2Dof(0.0)
572b0139 82{
83 //
84 // Named Constructor
85 //
86}
87
88//______________________________________________
89AliDielectronSignalFunc::~AliDielectronSignalFunc()
90{
91 //
92 // Default Destructor
93 //
bc75eeb5 94 if(fFuncSigBack) delete fFuncSigBack;
572b0139 95}
96
97
98//______________________________________________
bc75eeb5 99void AliDielectronSignalFunc::Process(TObjArray * const arrhist)
572b0139 100{
101 //
bc75eeb5 102 // Fit the invariant mass histograms and retrieve the signal and background
572b0139 103 //
bc75eeb5 104 switch(fMethod) {
105 case kFitted :
106 ProcessFit(arrhist);
107 break;
3505bfad 108
bc75eeb5 109 case kLikeSign :
110 ProcessLS(arrhist);
111 break;
3505bfad 112
bc75eeb5 113 case kEventMixing :
3505bfad 114 ProcessEM(arrhist);
bc75eeb5 115 break;
3505bfad 116
bc75eeb5 117 default :
118 AliError("Background substraction method not known!");
572b0139 119 }
bc75eeb5 120}
5720c765 121//______________________________________________________________________________
e43f22fb 122Double_t AliDielectronSignalFunc::PeakFunMC(const Double_t *x, const Double_t *par) {
5720c765 123 // Fit MC signal shape
124 // parameters
125 // [0]: scale for simpeak
126
127 Double_t xx = x[0];
128
129 TH1F *hPeak = fgHistSimPM;
130 if (!hPeak) {
ac390e40 131 printf("E-AliDielectronSignalFunc::PeakFun: No histogram for peak fit defined!\n");
132 return 0.0;
5720c765 133 }
134
135 Int_t idx = hPeak->FindBin(xx);
136 if ((idx <= 1) ||
137 (idx >= hPeak->GetNbinsX())) {
138 return 0.0;
139 }
140
e43f22fb 141 return (par[0+fNparBgnd] * hPeak->GetBinContent(idx));
5720c765 142
5720c765 143}
144
145//______________________________________________________________________________
146Double_t AliDielectronSignalFunc::PeakFunCB(const Double_t *x, const Double_t *par) {
e43f22fb 147 // Crystal Ball fit function
5720c765 148
e43f22fb 149 Double_t alpha = par[0+fNparBgnd];
150 Double_t n = par[1+fNparBgnd];
151 Double_t meanx = par[2+fNparBgnd];
152 Double_t sigma = par[3+fNparBgnd];
153 Double_t nn = par[4+fNparBgnd];
5720c765 154
155 Double_t a = TMath::Power((n/TMath::Abs(alpha)), n) * TMath::Exp(-.5*alpha*alpha);
156 Double_t b = n/TMath::Abs(alpha) - TMath::Abs(alpha);
157
158 Double_t arg = (x[0] - meanx)/sigma;
159 Double_t fitval = 0;
160
161 if (arg > -1.*alpha) {
162 fitval = nn * TMath::Exp(-.5*arg*arg);
163 } else {
164 fitval = nn * a * TMath::Power((b-arg), (-1*n));
165 }
166
167 return fitval;
168}
169
e43f22fb 170//______________________________________________________________________________
171Double_t AliDielectronSignalFunc::PeakFunGaus(const Double_t *x, const Double_t *par) {
172 // Gaussian fit function
173
174 Double_t n = par[0+fNparBgnd];
175 Double_t mean = par[1+fNparBgnd];
176 Double_t sigma = par[2+fNparBgnd];
177 Double_t xx = x[0];
178
179 return ( n*exp(-0.5*TMath::Power((xx-mean)/sigma,2)) );
180}
572b0139 181
bc75eeb5 182//______________________________________________
183void AliDielectronSignalFunc::ProcessFit(TObjArray * const arrhist) {
184 //
185 // Fit the +- invariant mass distribution only
186 // Here we assume that the combined fit function is a sum of the signal and background functions
187 // and that the signal function is always the first term of this sum
188 //
3505bfad 189
bc75eeb5 190 fHistDataPM = (TH1F*)(arrhist->At(1))->Clone("histPM"); // +- SE
191 fHistDataPM->Sumw2();
192 if(fRebin>1)
193 fHistDataPM->Rebin(fRebin);
572b0139 194
5720c765 195 fHistSignal = new TH1F("HistSignal", "Fit substracted signal",
3505bfad 196 fHistDataPM->GetXaxis()->GetNbins(),
197 fHistDataPM->GetXaxis()->GetXmin(), fHistDataPM->GetXaxis()->GetXmax());
5720c765 198 fHistBackground = new TH1F("HistBackground", "Fit contribution",
3505bfad 199 fHistDataPM->GetXaxis()->GetNbins(),
200 fHistDataPM->GetXaxis()->GetXmin(), fHistDataPM->GetXaxis()->GetXmax());
201
202 // the starting parameters of the fit function and their limits can be tuned
bc75eeb5 203 // by the user in its macro
204 fHistDataPM->Fit(fFuncSigBack, fFitOpt.Data(), "", fFitMin, fFitMax);
205 TFitResultPtr pmFitPtr = fHistDataPM->Fit(fFuncSigBack, fFitOpt.Data(), "", fFitMin, fFitMax);
cafdf10b 206 //TFitResult *pmFitResult = pmFitPtr.Get(); // used only with TF1Helper
e43f22fb 207 fFuncBackground->SetParameters(fFuncSigBack->GetParameters());
208 fFuncSignal->SetParameters(fFuncSigBack->GetParameters()+fFuncBackground->GetNpar());
209 // fFuncBackground->SetParameters(fFuncSigBack->GetParameters()+fFuncSignal->GetNpar());
572b0139 210
bc75eeb5 211 for(Int_t iBin=1; iBin<=fHistDataPM->GetXaxis()->GetNbins(); iBin++) {
212 Double_t m = fHistDataPM->GetBinCenter(iBin);
213 Double_t pm = fHistDataPM->GetBinContent(iBin);
214 Double_t epm = fHistDataPM->GetBinError(iBin);
215 Double_t bknd = fFuncBackground->Eval(m);
216 Double_t ebknd = 0;
217 for(Int_t iPar=fFuncSignal->GetNpar(); iPar<fFuncSigBack->GetNpar(); iPar++) {
cafdf10b 218/* TF1Helper problem on alien compilation
bc75eeb5 219 for(Int_t jPar=iPar; jPar<fFuncSigBack->GetNpar(); jPar++) {
220 TF1 gradientIpar("gradientIpar",
3505bfad 221 ROOT::TF1Helper::TGradientParFunction(iPar-fFuncSignal->GetNpar(),fFuncBackground),0,0,0);
bc75eeb5 222 TF1 gradientJpar("gradientJpar",
3505bfad 223 ROOT::TF1Helper::TGradientParFunction(jPar-fFuncSignal->GetNpar(),fFuncBackground),0,0,0);
224 ebknd += pmFitResult->CovMatrix(iPar,jPar)*
225 gradientIpar.Eval(m)*gradientJpar.Eval(m)*
226 (iPar==jPar ? 1.0 : 2.0);
bc75eeb5 227 }
cafdf10b 228*/
bc75eeb5 229 }
230 Double_t signal = pm-bknd;
231 Double_t error = TMath::Sqrt(epm*epm+ebknd);
232 fHistSignal->SetBinContent(iBin, signal);
233 fHistSignal->SetBinError(iBin, error);
234 fHistBackground->SetBinContent(iBin, bknd);
235 fHistBackground->SetBinError(iBin, TMath::Sqrt(ebknd));
236 }
3505bfad 237
bc75eeb5 238 if(fUseIntegral) {
239 // signal
240 fValues(0) = fFuncSignal->Integral(fIntMin, fIntMax)/fHistDataPM->GetBinWidth(1);
241 fErrors(0) = 0;
242 for(Int_t iPar=0; iPar<fFuncSignal->GetNpar(); iPar++) {
cafdf10b 243/* TF1Helper problem on alien compilation
bc75eeb5 244 for(Int_t jPar=iPar; jPar<fFuncSignal->GetNpar(); jPar++) {
245 TF1 gradientIpar("gradientIpar",
3505bfad 246 ROOT::TF1Helper::TGradientParFunction(iPar,fFuncSignal),0,0,0);
bc75eeb5 247 TF1 gradientJpar("gradientJpar",
3505bfad 248 ROOT::TF1Helper::TGradientParFunction(jPar,fFuncSignal),0,0,0);
249 fErrors(0) += pmFitResult->CovMatrix(iPar,jPar)*
250 gradientIpar.Integral(fIntMin,fIntMax)*gradientJpar.Integral(fIntMin,fIntMax)*
251 (iPar==jPar ? 1.0 : 2.0);
bc75eeb5 252 }
cafdf10b 253*/
3505bfad 254 }
bc75eeb5 255 // background
256 fValues(1) = fFuncBackground->Integral(fIntMin, fIntMax)/fHistDataPM->GetBinWidth(1);
257 fErrors(1) = 0;
258 for(Int_t iPar=fFuncSignal->GetNpar(); iPar<fFuncSigBack->GetNpar(); iPar++) {
cafdf10b 259/* TF1Helper problem on alien compilation
bc75eeb5 260 for(Int_t jPar=iPar; jPar<fFuncSigBack->GetNpar(); jPar++) {
261 TF1 gradientIpar("gradientIpar",
3505bfad 262 ROOT::TF1Helper::TGradientParFunction(iPar-fFuncSignal->GetNpar(),fFuncBackground),0,0,0);
bc75eeb5 263 TF1 gradientJpar("gradientJpar",
3505bfad 264 ROOT::TF1Helper::TGradientParFunction(jPar-fFuncSignal->GetNpar(),fFuncBackground),0,0,0);
265 fErrors(1) += pmFitResult->CovMatrix(iPar,jPar)*
266 gradientIpar.Integral(fIntMin, fIntMax)*gradientJpar.Integral(fIntMin, fIntMax)*
267 (iPar==jPar ? 1.0 : 2.0);
bc75eeb5 268 }
cafdf10b 269*/
bc75eeb5 270 }
572b0139 271 }
bc75eeb5 272 else {
273 // signal
274 fValues(0) = fHistSignal->IntegralAndError(fHistSignal->FindBin(fIntMin),
3505bfad 275 fHistSignal->FindBin(fIntMax), fErrors(0));
bc75eeb5 276 // background
277 fValues(1) = fHistBackground->IntegralAndError(fHistBackground->FindBin(fIntMin),
3505bfad 278 fHistBackground->FindBin(fIntMax),
279 fErrors(1));
8df8e382 280 }
bc75eeb5 281 // S/B and significance
282 SetSignificanceAndSOB();
283 fValues(4) = fFuncSigBack->GetParameter(fParMass);
284 fErrors(4) = fFuncSigBack->GetParError(fParMass);
285 fValues(5) = fFuncSigBack->GetParameter(fParMassWidth);
286 fErrors(5) = fFuncSigBack->GetParError(fParMassWidth);
3505bfad 287
bc75eeb5 288 fProcessed = kTRUE;
e43f22fb 289
290 fHistBackground->GetListOfFunctions()->Add(fFuncBackground);
291
572b0139 292}
293
294//______________________________________________
bc75eeb5 295void AliDielectronSignalFunc::ProcessLS(TObjArray * const arrhist) {
572b0139 296 //
bc75eeb5 297 // Substract background using the like-sign spectrum
572b0139 298 //
bc75eeb5 299 fHistDataPP = (TH1F*)(arrhist->At(0))->Clone("histPP"); // ++ SE
300 fHistDataPM = (TH1F*)(arrhist->At(1))->Clone("histPM"); // +- SE
3505bfad 301 fHistDataMM = (TH1F*)(arrhist->At(2))->Clone("histMM"); // -- SE
302 if (fRebin>1) {
bc75eeb5 303 fHistDataPP->Rebin(fRebin);
304 fHistDataPM->Rebin(fRebin);
305 fHistDataMM->Rebin(fRebin);
3505bfad 306 }
bc75eeb5 307 fHistDataPP->Sumw2();
308 fHistDataPM->Sumw2();
309 fHistDataMM->Sumw2();
3505bfad 310
311 fHistSignal = new TH1F("HistSignal", "Like-Sign substracted signal",
312 fHistDataPM->GetXaxis()->GetNbins(),
313 fHistDataPM->GetXaxis()->GetXmin(), fHistDataPM->GetXaxis()->GetXmax());
314 fHistBackground = new TH1F("HistBackground", "Like-sign contribution",
315 fHistDataPM->GetXaxis()->GetNbins(),
316 fHistDataPM->GetXaxis()->GetXmin(), fHistDataPM->GetXaxis()->GetXmax());
317
bc75eeb5 318 // fit the +- mass distribution
319 fHistDataPM->Fit(fFuncSigBack, fFitOpt.Data(), "", fFitMin, fFitMax);
320 fHistDataPM->Fit(fFuncSigBack, fFitOpt.Data(), "", fFitMin, fFitMax);
321 // declare the variables where the like-sign fit results will be stored
45b2b1b8 322// TFitResult *ppFitResult = 0x0;
323// TFitResult *mmFitResult = 0x0;
bc75eeb5 324 // fit the like sign background
325 TF1 *funcClonePP = (TF1*)fFuncBackground->Clone("funcClonePP");
326 TF1 *funcCloneMM = (TF1*)fFuncBackground->Clone("funcCloneMM");
327 fHistDataPP->Fit(funcClonePP, fFitOpt.Data(), "", fFitMin, fFitMax);
45b2b1b8 328 fHistDataPP->Fit(funcClonePP, fFitOpt.Data(), "", fFitMin, fFitMax);
329// TFitResultPtr ppFitPtr = fHistDataPP->Fit(funcClonePP, fFitOpt.Data(), "", fFitMin, fFitMax);
330// ppFitResult = ppFitPtr.Get();
331 fHistDataMM->Fit(funcCloneMM, fFitOpt.Data(), "", fFitMin, fFitMax);
bc75eeb5 332 fHistDataMM->Fit(funcCloneMM, fFitOpt.Data(), "", fFitMin, fFitMax);
45b2b1b8 333// TFitResultPtr mmFitPtr = fHistDataMM->Fit(funcCloneMM, fFitOpt.Data(), "", fFitMin, fFitMax);
334// mmFitResult = mmFitPtr.Get();
3505bfad 335
bc75eeb5 336 for(Int_t iBin=1; iBin<=fHistDataPM->GetXaxis()->GetNbins(); iBin++) {
337 Double_t m = fHistDataPM->GetBinCenter(iBin);
338 Double_t pm = fHistDataPM->GetBinContent(iBin);
339 Double_t pp = funcClonePP->Eval(m);
340 Double_t mm = funcCloneMM->Eval(m);
341 Double_t epm = fHistDataPM->GetBinError(iBin);
342 Double_t epp = 0;
343 for(Int_t iPar=0; iPar<funcClonePP->GetNpar(); iPar++) {
cafdf10b 344/* TF1Helper problem on alien compilation
bc75eeb5 345 for(Int_t jPar=iPar; jPar<funcClonePP->GetNpar(); jPar++) {
346 TF1 gradientIpar("gradientIpar",
3505bfad 347 ROOT::TF1Helper::TGradientParFunction(iPar,funcClonePP),0,0,0);
bc75eeb5 348 TF1 gradientJpar("gradientJpar",
3505bfad 349 ROOT::TF1Helper::TGradientParFunction(jPar,funcClonePP),0,0,0);
350 epp += ppFitResult->CovMatrix(iPar,jPar)*
351 gradientIpar.Eval(m)*gradientJpar.Eval(m)*
352 (iPar==jPar ? 1.0 : 2.0);
bc75eeb5 353 }
cafdf10b 354*/
bc75eeb5 355 }
356 Double_t emm = 0;
357 for(Int_t iPar=0; iPar<funcCloneMM->GetNpar(); iPar++) {
cafdf10b 358/* TF1Helper problem on alien compilation
bc75eeb5 359 for(Int_t jPar=iPar; jPar<funcCloneMM->GetNpar(); jPar++) {
360 TF1 gradientIpar("gradientIpar",
3505bfad 361 ROOT::TF1Helper::TGradientParFunction(iPar,funcCloneMM),0,0,0);
bc75eeb5 362 TF1 gradientJpar("gradientJpar",
3505bfad 363 ROOT::TF1Helper::TGradientParFunction(jPar,funcCloneMM),0,0,0);
364 emm += mmFitResult->CovMatrix(iPar,jPar)*
365 gradientIpar.Eval(m)*gradientJpar.Eval(m)*
366 (iPar==jPar ? 1.0 : 2.0);
bc75eeb5 367 }
cafdf10b 368*/
bc75eeb5 369 }
3505bfad 370
bc75eeb5 371 Double_t signal = pm-2.0*TMath::Sqrt(pp*mm);
372 Double_t background = 2.0*TMath::Sqrt(pp*mm);
373 // error propagation on the signal calculation above
374 Double_t esignal = TMath::Sqrt(epm*epm+(mm/pp)*epp+(pp/mm)*emm);
375 Double_t ebackground = TMath::Sqrt((mm/pp)*epp+(pp/mm)*emm);
376 fHistSignal->SetBinContent(iBin, signal);
377 fHistSignal->SetBinError(iBin, esignal);
378 fHistBackground->SetBinContent(iBin, background);
379 fHistBackground->SetBinError(iBin, ebackground);
380 }
2a14a7b1 381
bc75eeb5 382 // signal
383 fValues(0) = fHistSignal->IntegralAndError(fHistSignal->FindBin(fIntMin),
3505bfad 384 fHistSignal->FindBin(fIntMax), fErrors(0));
bc75eeb5 385 // background
386 fValues(1) = fHistBackground->IntegralAndError(fHistBackground->FindBin(fIntMin),
3505bfad 387 fHistBackground->FindBin(fIntMax),
388 fErrors(1));
bc75eeb5 389 // S/B and significance
390 SetSignificanceAndSOB();
391 fValues(4) = fFuncSigBack->GetParameter(fParMass);
392 fErrors(4) = fFuncSigBack->GetParError(fParMass);
393 fValues(5) = fFuncSigBack->GetParameter(fParMassWidth);
394 fErrors(5) = fFuncSigBack->GetParError(fParMassWidth);
3505bfad 395
bc75eeb5 396 fProcessed = kTRUE;
572b0139 397}
bc75eeb5 398
399//______________________________________________
3505bfad 400void AliDielectronSignalFunc::ProcessEM(TObjArray * const arrhist) {
bc75eeb5 401 //
402 // Substract background with the event mixing technique
403 //
3505bfad 404 arrhist->GetEntries(); // just to avoid the unused parameter warning
405 AliError("Event mixing for background substraction method not implemented!");
406}
bc75eeb5 407
572b0139 408//______________________________________________
bc75eeb5 409void AliDielectronSignalFunc::SetFunctions(TF1 * const combined, TF1 * const sig, TF1 * const back,
572b0139 410 Int_t parM, Int_t parMres)
411{
412 //
413 // Set the signal, background functions and combined fit function
414 // Note: The process method assumes that the first n parameters in the
415 // combined fit function correspond to the n parameters of the signal function
416 // and the n+1 to n+m parameters to the m parameters of the background function!!!
bc75eeb5 417
572b0139 418 if (!sig||!back||!combined) {
419 AliError("Both, signal and background function need to be set!");
61d106d3 420 return;
572b0139 421 }
bc75eeb5 422 fFuncSignal=sig;
423 fFuncBackground=back;
424 fFuncSigBack=combined;
425 fParMass=parM;
426 fParMassWidth=parMres;
e43f22fb 427
572b0139 428}
429
430//______________________________________________
431void AliDielectronSignalFunc::SetDefaults(Int_t type)
432{
433 //
434 // Setup some default functions:
435 // type = 0: gaus signal + linear background in 2.5 - 4 GeV inv. mass
436 // type = 1: gaus signal + exponential background in 2.5 - 4 GeV inv. mass
437 // type = 2: half gaussian, half exponential signal function
438 // type = 3: Crystal-Ball function
439 // type = 4: Crystal-Ball signal + exponential background
440 //
572b0139 441
442 if (type==0){
bc75eeb5 443 fFuncSignal=new TF1("DieleSignal","gaus",2.5,4);
444 fFuncBackground=new TF1("DieleBackground","pol1",2.5,4);
445 fFuncSigBack=new TF1("DieleCombined","gaus+pol1(3)",2.5,4);
572b0139 446
bc75eeb5 447 fFuncSigBack->SetParameters(1,3.1,.05,2.5,1);
448 fFuncSigBack->SetParLimits(0,0,10000000);
449 fFuncSigBack->SetParLimits(1,3.05,3.15);
450 fFuncSigBack->SetParLimits(2,.02,.1);
3505bfad 451 }
bc75eeb5 452 else if (type==1){
453 fFuncSignal=new TF1("DieleSignal","gaus",2.5,4);
454 fFuncBackground=new TF1("DieleBackground","[0]*exp(-(x-[1])/[2])",2.5,4);
455 fFuncSigBack=new TF1("DieleCombined","gaus+[3]*exp(-(x-[4])/[5])",2.5,4);
572b0139 456
bc75eeb5 457 fFuncSigBack->SetParameters(1,3.1,.05,1,2.5,1);
458 fFuncSigBack->SetParLimits(0,0,10000000);
459 fFuncSigBack->SetParLimits(1,3.05,3.15);
460 fFuncSigBack->SetParLimits(2,.02,.1);
3505bfad 461 }
bc75eeb5 462 else if (type==2){
572b0139 463 // half gaussian, half exponential signal function
464 // exponential background
bc75eeb5 465 fFuncSignal = new TF1("DieleSignal","(x<[1])*([0]*(exp(-0.5*((x-[1])/[2])^2)+exp((x-[1])/[3])*(1-exp(-0.5*((x-[1])/[2])^2))))+(x>=[1])*([0]*exp(-0.5*((x-[1])/[2])^2))",2.5,4);
466 fFuncBackground = new TF1("DieleBackground","[0]*exp(-(x-[1])/[2])+[3]",2.5,4);
467 fFuncSigBack = new TF1("DieleCombined","(x<[1])*([0]*(exp(-0.5*((x-[1])/[2])^2)+exp((x-[1])/[3])*(1-exp(-0.5*((x-[1])/[2])^2))))+(x>=[1])*([0]*exp(-0.5*((x-[1])/[2])^2))+[4]*exp(-(x-[5])/[6])+[7]",2.5,4);
468 fFuncSigBack->SetParameters(1.,3.1,.05,.1,1,2.5,1,0);
572b0139 469
bc75eeb5 470 fFuncSigBack->SetParLimits(0,0,10000000);
471 fFuncSigBack->SetParLimits(1,3.05,3.15);
472 fFuncSigBack->SetParLimits(2,.02,.1);
473 fFuncSigBack->FixParameter(6,2.5);
474 fFuncSigBack->FixParameter(7,0);
572b0139 475 }
476}
477
478
479//______________________________________________
480void AliDielectronSignalFunc::Draw(const Option_t* option)
481{
482 //
483 // Draw the fitted function
484 //
bc75eeb5 485
572b0139 486 TString drawOpt(option);
487 drawOpt.ToLower();
3505bfad 488
572b0139 489 Bool_t optStat=drawOpt.Contains("stat");
490
bc75eeb5 491 fFuncSigBack->SetNpx(200);
492 fFuncSigBack->SetRange(fIntMin,fIntMax);
493 fFuncBackground->SetNpx(200);
494 fFuncBackground->SetRange(fIntMin,fIntMax);
572b0139 495
bc75eeb5 496 TGraph *grSig=new TGraph(fFuncSigBack);
572b0139 497 grSig->SetFillColor(kGreen);
498 grSig->SetFillStyle(3001);
3505bfad 499
bc75eeb5 500 TGraph *grBack=new TGraph(fFuncBackground);
572b0139 501 grBack->SetFillColor(kRed);
502 grBack->SetFillStyle(3001);
3505bfad 503
572b0139 504 grSig->SetPoint(0,grBack->GetX()[0],grBack->GetY()[0]);
505 grSig->SetPoint(grSig->GetN()-1,grBack->GetX()[grBack->GetN()-1],grBack->GetY()[grBack->GetN()-1]);
506
507 grBack->SetPoint(0,grBack->GetX()[0],0.);
508 grBack->SetPoint(grBack->GetN()-1,grBack->GetX()[grBack->GetN()-1],0.);
509
bc75eeb5 510 fFuncSigBack->SetRange(fFitMin,fFitMax);
511 fFuncBackground->SetRange(fFitMin,fFitMax);
572b0139 512
513 if (!drawOpt.Contains("same")){
bc75eeb5 514 if (fHistDataPM){
515 fHistDataPM->Draw();
8df8e382 516 grSig->Draw("f");
517 } else {
518 grSig->Draw("af");
519 }
572b0139 520 } else {
521 grSig->Draw("f");
522 }
5720c765 523 if(fMethod==kFitted || fMethod==kFittedMC) grBack->Draw("f");
bc75eeb5 524 fFuncSigBack->Draw("same");
525 fFuncSigBack->SetLineWidth(2);
526 if(fMethod==kLikeSign) {
527 fHistDataPP->SetLineWidth(2);
528 fHistDataPP->SetLineColor(6);
529 fHistDataPP->Draw("same");
530 fHistDataMM->SetLineWidth(2);
3505bfad 531 fHistDataMM->SetLineColor(8);
bc75eeb5 532 fHistDataMM->Draw("same");
533 }
3505bfad 534
5720c765 535 if(fMethod==kFitted || fMethod==kFittedMC)
bc75eeb5 536 fFuncBackground->Draw("same");
3505bfad 537
8df8e382 538 if (optStat) DrawStats();
bc75eeb5 539
572b0139 540}
e43f22fb 541
542
543//______________________________________________________________________________
544void AliDielectronSignalFunc::CombineFunc(TF1 * const peak, TF1 * const bgnd) {
545 //
546 // combine the bgnd and the peak function
547 //
548
549 if (!peak||!bgnd) {
550 AliError("Both, signal and background function need to be set!");
551 return;
552 }
553 fFuncSignal=peak;
554 fFuncBackground=bgnd;
555
556 fNparPeak = fFuncSignal->GetNpar();
557 fNparBgnd = fFuncBackground->GetNpar();
558
559 fFuncSigBack = new TF1("BgndPeak",AliDielectronSignalFunc::PeakBgndFun, 0.,5.,fNparPeak+fNparBgnd);
560 return;
561}
562
563//______________________________________________________________________________
564Double_t AliDielectronSignalFunc::PeakBgndFun(const Double_t *x, const Double_t *par) {
565 //
566 // merge peak and bgnd functions
567 //
568 return (fFuncBackground->EvalPar(x,par) + fFuncSignal->EvalPar(x,par));
569}
570