PWGDQ/dielectron/AliDielectronBtoJPSItoEleCDFfitHandler.cxx

   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
  29 void CDFFunction(Int_t &npar, Double_t *gin, Double_t &f, Double_t *x, Int_t iflag);
  30
  31 ClassImp(AliDielectronBtoJPSItoEleCDFfitHandler)
  32
  33         //______________________________________________________________________________
  34 void 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 //_________________________________________________________________________________________________
  45 AliDielectronBtoJPSItoEleCDFfitHandler::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         //
  57         for (Int_t i=0; i<20; ++i) fParamStartValues[i]=0x0;
  58 }
  59 //_________________________________________________________________________________________________
  60 AliDielectronBtoJPSItoEleCDFfitHandler::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         //
  73   for (Int_t i=0; i<20; ++i) fParamStartValues[i]=0x0;
  74
  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 //___________________________________________________________________________
 102 AliDielectronBtoJPSItoEleCDFfitHandler& AliDielectronBtoJPSItoEleCDFfitHandler::operator=(const AliDielectronBtoJPSItoEleCDFfitHandler& c)
 103 {
 104         //
 105         // Assignment operator
 106         //
 107         if (this!=&c) {
 108     for (Int_t i=0; i<20; ++i) fParamStartValues[i]=c.fParamStartValues[i];
 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 //_______________________________________________________________________________________
 121 AliDielectronBtoJPSItoEleCDFfitHandler::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         //
 134   for (Int_t i=0; i<20; ++i) fParamStartValues[i]=c.fParamStartValues[i];
 135 }
 136 //_______________________________________________________________________________________
 137 AliDielectronBtoJPSItoEleCDFfitHandler::~AliDielectronBtoJPSItoEleCDFfitHandler()
 138 {
 139         //
 140         //destructor
 141         //
 142         delete fLikely;
 143 }
 144 //_______________________________________________________________________________________
 145 Int_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 //_______________________________________________________________________________________
 187 void 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 //_______________________________________________________________________________________
 210 void AliDielectronBtoJPSItoEleCDFfitHandler::SetParamStartValues(Double_t inputparamvalues[20])
 211 {
 212         for(Int_t index=0; index < 20; index++) fParamStartValues[index] = inputparamvalues[index];
 213 }
 214 //_______________________________________________________________________________________
 215 void AliDielectronBtoJPSItoEleCDFfitHandler::SetResolutionConstants(Double_t* resolutionConst)
 216 {
 217         //
 218         // Sets constants for the resolution function
 219         //
 220         fLikely->SetResolutionConstants(resolutionConst);
 221
 222 }
 223 //_______________________________________________________________________________________
 224 void 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 //_______________________________________________________________________________________
 233 void 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 //_______________________________________________________________________________________
 241 void 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