1 /**************************************************************************
2 * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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>
18 #include "AliBtoJPSItoEleCDFfitHandler.h"
19 #include "AliBtoJPSItoEleCDFfitFCN.h"
21 //-------------------------------------------------------------------------
22 // Class AliBtoJPSItoEleCDFfitHandler
23 // Class to perform unbinned log-likelihood fit
25 // Origin: C. Di Giglio
26 // Contact: carmelo.digiglio@ba.infn.it , giuseppe.bruno@ba.infn.it
27 //-------------------------------------------------------------------------
29 void CDFFunction(Int_t &npar, Double_t *gin, Double_t &f, Double_t *x, Int_t iflag);
31 ClassImp(AliBtoJPSItoEleCDFfitHandler)
33 //______________________________________________________________________________
34 void CDFFunction(Int_t &npar, Double_t *gin, Double_t &f, Double_t *par, Int_t iflag)
36 // This function is called by minuit
37 // The corresponding member method is called
38 // using SetObjectFit/GetObjectFit methods of TMinuit
39 AliBtoJPSItoEleCDFfitHandler* dummy = (AliBtoJPSItoEleCDFfitHandler *)TVirtualFitter::GetFitter()->GetObjectFit();
40 dummy->CdfFCN(npar, gin, f, par, iflag);
44 //_________________________________________________________________________________________________
45 AliBtoJPSItoEleCDFfitHandler::AliBtoJPSItoEleCDFfitHandler():
55 // default constructor
58 //_________________________________________________________________________________________________
59 AliBtoJPSItoEleCDFfitHandler::AliBtoJPSItoEleCDFfitHandler(Double_t* decaytime,
60 Double_t* invariantmass, Int_t ncand) :
72 AliInfo("\n+++\n+++ Minimization object AliBtoJPSItoEleCDFfitHandler created\n+++\n");
73 fLikely = new AliBtoJPSItoEleCDFfitFCN();
74 AliInfo("\n+++\n+++ CDF fit function object AliBtoJPSItoEleCDFfitFCN created\n+++\n");
75 AliInfo("Parameter 0 ----> fRadius");
76 AliInfo("Parameter 1 ----> fTheta");
77 AliInfo("Parameter 2 ----> fPhi");
78 AliInfo("Parameter 3 ----> fOneOvLamPlus");
79 AliInfo("Parameter 4 ----> fOneOvLamMinus");
80 AliInfo("Parameter 5 ----> fOneOvLamSym");
81 AliInfo("Parameter 6 ----> fMSlope");
82 AliInfo("Parameter 7 ----> fB");
83 AliInfo("Parameter 8 ----> fFsig");
84 AliInfo("Parameter 9 ----> fMmean");
85 AliInfo("Parameter 10 ----> fNexp");
86 AliInfo("Parameter 11 ----> fSigma");
87 AliInfo("Parameter 12 ----> fAlpha");
88 AliInfo("Parameter 13 ----> fNorm");
89 AliInfo("Parameter 14 ----> fSigmaResol");
90 AliInfo("Parameter 15 ----> fNResol");
91 AliInfo(Form("\n+++\n+++ Number of candidates ---> %d\n+++\n ", ncand));
93 //___________________________________________________________________________
94 AliBtoJPSItoEleCDFfitHandler& AliBtoJPSItoEleCDFfitHandler::operator=(const AliBtoJPSItoEleCDFfitHandler& c)
97 // Assignment operator
100 fIsParamFixed = c.fIsParamFixed;
101 fPrintStatus = c.fPrintStatus;
111 //_______________________________________________________________________________________
112 AliBtoJPSItoEleCDFfitHandler::AliBtoJPSItoEleCDFfitHandler(const AliBtoJPSItoEleCDFfitHandler& c) :
114 fIsParamFixed(c.fIsParamFixed),
115 fPrintStatus(c.fPrintStatus),
126 //_______________________________________________________________________________________
127 AliBtoJPSItoEleCDFfitHandler::~AliBtoJPSItoEleCDFfitHandler()
134 //_______________________________________________________________________________________
135 Int_t AliBtoJPSItoEleCDFfitHandler::DoMinimization()
138 // performs the minimization
140 static TVirtualFitter *fitter = TVirtualFitter::Fitter(this,16);
141 fitter->SetFCN(CDFFunction);
143 fitter->SetParameter(0,"fRadius",fParamStartValues[0], 1.e-06, 0., 1.);
144 fitter->SetParameter(1,"fTheta",fParamStartValues[1], 1.e-06, 0.,2*TMath::Pi());
145 fitter->SetParameter(2,"fPhi",fParamStartValues[2], 1.e-06, 0.,2*TMath::Pi());
146 // fitter->SetParameter(3,"fOneOvLamPlus",fParamStartValues[3], 1.e-10, 0., 5.e+01);
147 // fitter->SetParameter(4,"fOneOvLamMinus",fParamStartValues[4], 1.e-10, 0., 5.e+01);
148 // fitter->SetParameter(5,"fOneOvLamSym",fParamStartValues[5], 1.e-10, 0., 5.e+01);
149 fitter->SetParameter(3,"fOneOvLamPlus",fParamStartValues[3], 1.e-10, 0.0000001, 5.e+01);
150 fitter->SetParameter(4,"fOneOvLamMinus",fParamStartValues[4], 1.e-10, 0.00000001, 5.e+01);
151 fitter->SetParameter(5,"fOneOvLamSym",fParamStartValues[5], 1.e-10, 0.00000001, 5.e+01);
152 fitter->SetParameter(6,"fMSlope",fParamStartValues[6], 1.e-04, -2.5, 2.5);
153 fitter->SetParameter(7,"fB",fParamStartValues[7], 1.e-08, 0., 1.);
154 fitter->SetParameter(8,"fFsig",fParamStartValues[8], 1.e-08, 0., 1.);
155 fitter->SetParameter(9,"fMmean",fParamStartValues[9], 1.e-08, 0., 1.e+04);
156 fitter->SetParameter(10,"fNexp",fParamStartValues[10], 1.e-08, 0., 1.e+02);
157 fitter->SetParameter(11,"fSigma",fParamStartValues[11], 1.e-08, 0., 1.e+04);
158 fitter->SetParameter(12,"fAlpha",fParamStartValues[12], 1.e-08, 0., 1.e+04);
159 fitter->SetParameter(13,"fNorm",fParamStartValues[13], 1.e-08, 0., 1.e+01);
160 fitter->SetParameter(14,"fSigmaResol",fParamStartValues[14], 1.e-08, 0., 1.e+04);
161 fitter->SetParameter(15,"fNResol",fParamStartValues[15], 1.e-08, 0., 1.e+05);
163 for(UInt_t indexparam = 0; indexparam < 16; indexparam++){
164 if(IsParamFixed(indexparam))fitter->FixParameter((Int_t)indexparam);
167 Double_t arglist[2]={10000,1.0};
168 Int_t iret=fitter->ExecuteCommand("MIGRAD", arglist ,2);
169 fitter->PrintResults(4,0);
171 AliInfo("Minimization procedure finished\n");
175 //_______________________________________________________________________________________
176 void AliBtoJPSItoEleCDFfitHandler::CdfFCN(Int_t & /* npar */,
177 Double_t * /* gin */,Double_t &f,Double_t *par,Int_t /* iflag */)
180 // Definition of the FCN to be used by minuit
182 fLikely->SetAllParameters(par);
183 fLikely->ConvertFromSpherical();
184 fLikely->ComputeIntegral();
185 if(fPrintStatus)fLikely->PrintStatus();
190 f = fLikely->EvaluateLikelihood(fX,fM,fNcand);
193 AliDebug(2,Form("Real time spent to calculate function == %f \n", t.RealTime()));
194 AliDebug(2,Form("CPU time spent to calculate function == %f \n", t.CpuTime()));
195 AliDebug(2,Form("Actual value of the AliBtoJPSItoEleCDFfitFCN function == %f \n", f));
199 //_______________________________________________________________________________________
200 void AliBtoJPSItoEleCDFfitHandler::SetParamStartValues(Double_t inputparamvalues[16])
202 for(Int_t index=0; index < 16; index++) fParamStartValues[index] = inputparamvalues[index];
204 //_______________________________________________________________________________________
205 void AliBtoJPSItoEleCDFfitHandler::SetResolutionConstants()
208 // Sets constants for the resolution function
210 fLikely->SetResolutionConstants();
213 //_______________________________________________________________________________________
214 void AliBtoJPSItoEleCDFfitHandler::SetCrystalBallFunction(Bool_t okCB)
217 // Sets the CB as the parametrization for the signal invariant mass spectrum
218 // (otherwise Landau is chosen)
220 fLikely->SetCrystalBallFunction(okCB);
222 //_______________________________________________________________________________________
223 void AliBtoJPSItoEleCDFfitHandler::SetMassWndHigh(Double_t limit)
226 // Sets upper limit for the invariant mass window (under J/PSI mass peak)
228 fLikely->SetMassWndHigh(limit);
230 //_______________________________________________________________________________________
231 void AliBtoJPSItoEleCDFfitHandler::SetMassWndLow(Double_t limit)
234 // Sets lower limit for the invariant mass window (under J/PSI mass peak)
236 fLikely->SetMassWndLow(limit);