]>
Commit | Line | Data |
---|---|---|
572b0139 | 1 | #ifndef ALIDIELECTRONSIGNALFUNC_H |
2 | #define ALIDIELECTRONSIGNALFUNC_H | |
3 | ||
4 | /* Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. * | |
5 | * See cxx source for full Copyright notice */ | |
6 | ||
7 | //############################################################# | |
8 | //# # | |
9 | //# Class AliDielectronSignalFunc # | |
10 | //# # | |
11 | //# Authors: # | |
12 | //# Anton Andronic, GSI / A.Andronic@gsi.de # | |
13 | //# Ionut C. Arsene, GSI / I.C.Arsene@gsi.de # | |
14 | //# Julian Book, Uni Ffm / Julian.Book@cern.ch # | |
15 | //# Frederick Kramer, Uni Ffm, / Frederick.Kramer@cern.ch # | |
16 | //# Magnus Mager, CERN / Magnus.Mager@cern.ch # | |
17 | //# WooJin J. Park, GSI / W.J.Park@gsi.de # | |
18 | //# Jens Wiechula, Uni HD / Jens.Wiechula@cern.ch # | |
19 | //# # | |
20 | //############################################################# | |
21 | ||
bc75eeb5 | 22 | /* |
23 | Class used for extracting the signal from an invariant mass spectrum. | |
24 | It implements the AliDielectronSignalBase class and it uses user provided | |
25 | functions to fit the unlike-sign spectrum (and the like-sign one). | |
26 | ||
27 | Example usage: | |
28 | ||
29 | AliDielectronSignalFunc *signalProcess = new AliDielectronSignalFunc(); | |
30 | TObjArray *histoArray = new TObjArray(); | |
31 | histoArray->Add(signalPP); // add the spectrum histograms to the array | |
32 | histoArray->Add(signalPM); // the order is important !!! | |
33 | histoArray->Add(signalMM); | |
34 | // set the extraction method | |
35 | // AliDielectronSignalBase::kFitted --> fit only the unlike-sign spectrum and extract everything from that | |
36 | // AliDielectronSignalBase::kLikeSign --> fit both the unlike- and like-sign spectra | |
37 | // AliDielectronSignalBase::kEventMixing --> fit both the unlike- and like-sign spectra from event mixing | |
38 | signalProcess->SetMethod(AliDielectronSignalBase::kLikeSign); | |
39 | // Initialize the functions to be used and pass them to the signal object | |
40 | // External preparation of the functions can(should) be done as this can be a 5 or more parameter fit | |
41 | TF1* gaus = new TF1("gaus", "gaus", 0., 4.); | |
42 | TF1* expo = new TF1("expo", "[0]*exp([1]*x)", 0., 4.); | |
43 | TF1* combined = new TF1("combined", "gaus + [3]*exp([4]*x)", 0.,4.); | |
44 | combined->SetParameter(1, 3.1); | |
45 | combined->SetParameter(1, 0.1); | |
46 | signalPP->Fit(expo, "SME", "", 2.4, 4.0); | |
47 | signalPM->Fit(gaus, "SME", "", 3.0, 3.15); | |
48 | Double_t pars[5]; | |
49 | gaus->GetParameters(&pars[0]); | |
50 | expo->GetParameters(&pars[3]); | |
51 | combined->SetParameters(pars); | |
52 | combined->SetParLimits(1, 3.05, 3.15); | |
53 | combined->SetParLimits(2, 0.03, 0.1); | |
54 | signalProcess->SetFunctions(combined, gaus, expo, 1, 2); | |
55 | ||
56 | signalProcess->SetFitRange(2.4,4.0); | |
57 | // Use the integral of the fit function to estimate the signal or not | |
58 | // The background will always be estimated from the fit | |
59 | // signalProcess->SetUseIntegral(kTRUE); | |
60 | signalProcess->SetFitOption("SME"); | |
61 | // Give the range where the signal is calculated | |
62 | signalProcess->SetIntegralRange(3.0,3.15); | |
63 | signalProcess->SetRebin(2); | |
64 | signalProcess->Process(histoArray); | |
65 | signalProcess->Draw("stat"); | |
66 | signalProcess->Print(); | |
67 | */ | |
68 | ||
69 | ||
572b0139 | 70 | #include <TVectorT.h> |
71 | #include <TString.h> | |
bc75eeb5 | 72 | #include <TH1F.h> |
572b0139 | 73 | |
74 | #include "AliDielectronSignalBase.h" | |
75 | ||
76 | class AliDielectronSignalFunc : public AliDielectronSignalBase { | |
77 | public: | |
78 | AliDielectronSignalFunc(); | |
79 | AliDielectronSignalFunc(const char*name, const char* title); | |
bc75eeb5 | 80 | AliDielectronSignalFunc(const AliDielectronSignalFunc &c); |
81 | AliDielectronSignalFunc &operator=(const AliDielectronSignalFunc &c); | |
572b0139 | 82 | |
83 | virtual ~AliDielectronSignalFunc(); | |
84 | ||
572b0139 | 85 | virtual void Process(TObjArray * const arrhist); |
bc75eeb5 | 86 | void ProcessFit(TObjArray * const arrhist); // fit the SE +- distribution |
87 | void ProcessLS(TObjArray * const arrhist); // substract the fitted SE like-sign background | |
3505bfad | 88 | void ProcessEM(TObjArray * const arrhist); // substract the fitted SE+ME like-sign background |
572b0139 | 89 | |
bc75eeb5 | 90 | void SetUseIntegral(Bool_t flag=kTRUE) {fUseIntegral = flag;}; |
91 | void SetFunctions(TF1 * const combined, TF1 * const sig=0, TF1 * const back=0, Int_t parM=1, Int_t parMres=2); | |
92 | void SetFitOption(const char* opt) { | |
93 | fFitOpt=opt; | |
94 | fFitOpt.ToLower(); | |
95 | if(!fFitOpt.Contains("s")) fFitOpt += "s"; | |
96 | } | |
572b0139 | 97 | void SetDefaults(Int_t type); |
bc75eeb5 | 98 | |
3505bfad | 99 | TF1* GetSignalFunction() const { return fFuncSignal; } |
100 | TF1* GetBackgroundFunction() const { return fFuncBackground; } | |
101 | TF1* GetCombinedFunction() const { return fFuncSigBack; } | |
572b0139 | 102 | |
103 | virtual void Draw(const Option_t* option = ""); | |
104 | ||
105 | private: | |
8df8e382 | 106 | |
bc75eeb5 | 107 | TF1 *fFuncSignal; // Function for the signal description |
108 | TF1 *fFuncBackground; // Function for the background description | |
109 | TF1 *fFuncSigBack; // Combined function signal plus background | |
110 | Int_t fParMass; // the index of the parameter corresponding to the resonance mass | |
111 | Int_t fParMassWidth; // the index of the parameter corresponding to the resonance mass width | |
572b0139 | 112 | |
bc75eeb5 | 113 | TString fFitOpt; // fit option used |
114 | Bool_t fUseIntegral; // use the integral of the fitted functions to extract signal and background | |
572b0139 | 115 | |
bc75eeb5 | 116 | ClassDef(AliDielectronSignalFunc,2) // Dielectron SignalFunc |
572b0139 | 117 | }; |
118 | ||
572b0139 | 119 | #endif |