Add fast merging option (Diego)
[u/mrichter/AliRoot.git] / PWG3 / dielectron / AliDielectronBtoJPSItoEleCDFfitHandler.cxx
CommitLineData
ba15fdfb 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#include <TVirtualFitter.h>
16#include <TStopwatch.h>
17#include "AliLog.h"
18#include "AliDielectronBtoJPSItoEleCDFfitHandler.h"
19#include "AliDielectronBtoJPSItoEleCDFfitFCN.h"
20
21//-------------------------------------------------------------------------
22// Class AliDielectronBtoJPSItoEleCDFfitHandler
23// Class to perform unbinned log-likelihood fit
24//
25// Origin: C. Di Giglio
26// Contact: carmelo.digiglio@ba.infn.it , giuseppe.bruno@ba.infn.it
27//-------------------------------------------------------------------------
28
29void CDFFunction(Int_t &npar, Double_t *gin, Double_t &f, Double_t *x, Int_t iflag);
30
31ClassImp(AliDielectronBtoJPSItoEleCDFfitHandler)
32
33 //______________________________________________________________________________
34void CDFFunction(Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t iflag)
35{
36 // This function is called by minuit
37 // The corresponding member method is called
38 // using SetObjectFit/GetObjectFit methods of TMinuit
39 AliDielectronBtoJPSItoEleCDFfitHandler* dummy = (AliDielectronBtoJPSItoEleCDFfitHandler *)TVirtualFitter::GetFitter()->GetObjectFit();
40 dummy->CdfFCN(npar, gin, f, par, iflag);
41}
42
43
44//_________________________________________________________________________________________________
45AliDielectronBtoJPSItoEleCDFfitHandler::AliDielectronBtoJPSItoEleCDFfitHandler():
46 fIsParamFixed(20),
47 fPrintStatus(kFALSE),
48 fUp(0),
49 fX(0x0),
50 fM(0x0),
51 fLikely(0x0),
52 fNcand(0)
53{
54 //
55 // default constructor
56 //
5cc8c825 57 for (Int_t i=0; i<20; ++i) fParamStartValues[i]=0x0;
ba15fdfb 58}
59//_________________________________________________________________________________________________
60AliDielectronBtoJPSItoEleCDFfitHandler::AliDielectronBtoJPSItoEleCDFfitHandler(Double_t* decaytime,
61 Double_t* invariantmass, Int_t ncand) :
62 fIsParamFixed(20),
63 fPrintStatus(kFALSE),
64 fUp(0),
65 fX(decaytime),
66 fM(invariantmass),
67 fLikely(0x0),
68 fNcand(ncand)
69{
70 //
71 // constructor
72 //
5cc8c825 73 for (Int_t i=0; i<20; ++i) fParamStartValues[i]=0x0;
74
ba15fdfb 75 AliInfo("\n+++\n+++ Minimization object AliDielectronBtoJPSItoEleCDFfitHandler created\n+++\n");
76 fLikely = new AliDielectronBtoJPSItoEleCDFfitFCN();
77 AliInfo("\n+++\n+++ CDF fit function object AliDielectronBtoJPSItoEleCDFfitFCN created\n+++\n");
78 AliInfo("Parameter 0 ----> fWeightRes");
79 AliInfo("Parameter 1 ----> fPos");
80 AliInfo("Parameter 2 ----> fNeg");
81 AliInfo("Parameter 3 ----> fSym");
82 AliInfo("Parameter 4 ----> fOneOvLamPlus");
83 AliInfo("Parameter 5 ----> fOneOvLamMinus");
84 AliInfo("Parameter 6 ----> fOneOvLamSym");
85 AliInfo("Parameter 7 ----> fB");
86 AliInfo("Parameter 8 ----> fFsig");
87 AliInfo("Parameter 9 ----> fMmean");
88 AliInfo("Parameter 10 ----> fNexp");
89 AliInfo("Parameter 11 ----> fSigma");
90 AliInfo("Parameter 12 ----> fAlpha");
91 AliInfo("Parameter 13 ----> fNorm");
92 AliInfo("Parameter 14 ----> fBkgNorm");
93 AliInfo("Parameter 15 ----> fBkgMean");
94 AliInfo("Parameter 16 ----> fBkgSlope");
95 AliInfo("Parameter 17 ----> fBkgConst");
96 AliInfo("Parameter 18 ----> fNormGaus1");
97 AliInfo("Parameter 19 ----> fNormGaus2");
98
99 AliInfo(Form("\n+++\n+++ Number of candidates ---> %d\n+++\n ", ncand));
100}
101//___________________________________________________________________________
102AliDielectronBtoJPSItoEleCDFfitHandler& AliDielectronBtoJPSItoEleCDFfitHandler::operator=(const AliDielectronBtoJPSItoEleCDFfitHandler& c)
103{
104 //
105 // Assignment operator
106 //
107 if (this!=&c) {
5cc8c825 108 for (Int_t i=0; i<20; ++i) fParamStartValues[i]=c.fParamStartValues[i];
ba15fdfb 109 fIsParamFixed = c.fIsParamFixed;
110 fPrintStatus = c.fPrintStatus;
111 fUp = c.fUp;
112 fX = c.fX;
113 fM = c.fM;
114 fLikely = c.fLikely;
115 fNcand = c.fNcand;
116 }
117 return *this;
118}
119
120//_______________________________________________________________________________________
121AliDielectronBtoJPSItoEleCDFfitHandler::AliDielectronBtoJPSItoEleCDFfitHandler(const AliDielectronBtoJPSItoEleCDFfitHandler& c) :
122 TNamed(c),
123 fIsParamFixed(c.fIsParamFixed),
124 fPrintStatus(c.fPrintStatus),
125 fUp(c.fUp),
126 fX(c.fX),
127 fM(c.fM),
128 fLikely(c.fLikely),
129 fNcand(c.fNcand)
130{
131 //
132 // Copy Constructor
133 //
5cc8c825 134 for (Int_t i=0; i<20; ++i) fParamStartValues[i]=c.fParamStartValues[i];
ba15fdfb 135}
136//_______________________________________________________________________________________
137AliDielectronBtoJPSItoEleCDFfitHandler::~AliDielectronBtoJPSItoEleCDFfitHandler()
138{
139 //
140 //destructor
141 //
142 delete fLikely;
143}
144//_______________________________________________________________________________________
145Int_t AliDielectronBtoJPSItoEleCDFfitHandler::DoMinimization()
146{
147 //
148 // performs the minimization
149 //
150 static TVirtualFitter *fitter = TVirtualFitter::Fitter(this,20);
151 fitter->SetFCN(CDFFunction);
152
153 fitter->SetParameter(0,"fWeightRes",fParamStartValues[0], 1.e-08, 0., 1.e+06);
154 fitter->SetParameter(1,"fPos",fParamStartValues[1], 1.e-08, 0.,1.e+06);
155 fitter->SetParameter(2,"fNeg",fParamStartValues[2], 1.e-08, 0.,1.e+06);
156 fitter->SetParameter(3,"fSym",fParamStartValues[3], 1.e-08, 0.,1.e+06);
157 fitter->SetParameter(4,"fOneOvLamPlus",fParamStartValues[4], 1.e-10, 0.0000001, 5.e+01);
158 fitter->SetParameter(5,"fOneOvLamMinus",fParamStartValues[5], 1.e-10, 0.00000001, 5.e+01);
159 fitter->SetParameter(6,"fOneOvLamSym",fParamStartValues[6], 1.e-10, 0.00000001, 5.e+01);
160 fitter->SetParameter(7,"fB",fParamStartValues[7], 1.e-10, 0., 1.);
161 fitter->SetParameter(8,"fFsig",fParamStartValues[8], 1.e-10, 0., 1.);
162 fitter->SetParameter(9,"fMmean",fParamStartValues[9], 1.e-08, 0., 1.e+04);
163 fitter->SetParameter(10,"fNexp",fParamStartValues[10], 1.e-08, 0., 1.e+02);
164 fitter->SetParameter(11,"fSigma",fParamStartValues[11], 1.e-08, 0., 1.e+04);
165 fitter->SetParameter(12,"fAlpha",fParamStartValues[12], 1.e-08, 0., 1.e+04);
166 fitter->SetParameter(13,"fNorm",fParamStartValues[13], 1.e-08, 0., 1.e+04);
167 fitter->SetParameter(14,"fBkgNorm",fParamStartValues[14], 1.e-08, 0., 1.e+04);
168 fitter->SetParameter(15,"fBkgMean",fParamStartValues[15], 1.e-08, 0., 1.e+04);
169 fitter->SetParameter(16,"fBkgSlope",fParamStartValues[16], 1.e-08, 0., 1.e+04);
170 fitter->SetParameter(17,"fBkgSlope",fParamStartValues[17], 1.e-08, 0., 1.e+04);
171 fitter->SetParameter(18,"fNormGaus1",fParamStartValues[18], 1.e-08, 0., 1.e+05);
172 fitter->SetParameter(19,"fNormGaus2",fParamStartValues[19], 1.e-08, 0., 1.e+05);
173
174 for(UInt_t indexparam = 0; indexparam < 20; indexparam++){
175 if(IsParamFixed(indexparam))fitter->FixParameter((Int_t)indexparam);
176 }
177
178 Double_t arglist[2]={10000,1.0};
179 Int_t iret=fitter->ExecuteCommand("MIGRAD", arglist ,2);
180 fitter->PrintResults(4,0);
181
182 AliInfo("Minimization procedure finished\n");
183
184 return iret;
185}
186//_______________________________________________________________________________________
187void AliDielectronBtoJPSItoEleCDFfitHandler::CdfFCN(Int_t & /* npar */,
188 Double_t * /* gin */,Double_t &f,Double_t *par,Int_t /* iflag */)
189{
190 //
191 // Definition of the FCN to be used by minuit
192 //
193 fLikely->SetAllParameters(par);
194 fLikely->ComputeMassIntegral();
195 if(fPrintStatus)fLikely->PrintStatus();
196
197 TStopwatch t;
198 t.Start();
199
200 f = fLikely->EvaluateLikelihood(fX,fM,fNcand);
201
202 t.Stop();
203 AliDebug(2,Form("Real time spent to calculate function == %f \n", t.RealTime()));
204 AliDebug(2,Form("CPU time spent to calculate function == %f \n", t.CpuTime()));
205 AliDebug(2,Form("Actual value of the AliDielectronBtoJPSItoEleCDFfitFCN function == %f \n", f));
206
207 return;
208}
209//_______________________________________________________________________________________
210void AliDielectronBtoJPSItoEleCDFfitHandler::SetParamStartValues(Double_t inputparamvalues[20])
211{
212 for(Int_t index=0; index < 20; index++) fParamStartValues[index] = inputparamvalues[index];
213}
214//_______________________________________________________________________________________
215void AliDielectronBtoJPSItoEleCDFfitHandler::SetResolutionConstants(Double_t* resolutionConst)
216{
217 //
218 // Sets constants for the resolution function
219 //
220 fLikely->SetResolutionConstants(resolutionConst);
221
222}
223//_______________________________________________________________________________________
224void AliDielectronBtoJPSItoEleCDFfitHandler::SetCrystalBallFunction(Bool_t okCB)
225{
226 //
227 // Sets the CB as the parametrization for the signal invariant mass spectrum
228 // (otherwise Landau is chosen)
229 //
230 fLikely->SetCrystalBallFunction(okCB);
231}
232//_______________________________________________________________________________________
233void AliDielectronBtoJPSItoEleCDFfitHandler::SetMassWndHigh(Double_t limit)
234{
235 //
236 // Sets upper limit for the invariant mass window (under J/PSI mass peak)
237 //
238 fLikely->SetMassWndHigh(limit);
239}
240//_______________________________________________________________________________________
241void AliDielectronBtoJPSItoEleCDFfitHandler::SetMassWndLow(Double_t limit)
242{
243 //
244 // Sets lower limit for the invariant mass window (under J/PSI mass peak)
245 //
246 fLikely->SetMassWndLow(limit);
247}
248