PWG3/vertexingHF/AliBtoJPSItoEleCDFfitHandler.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 "AliBtoJPSItoEleCDFfitHandler.h"
  19 #include "AliBtoJPSItoEleCDFfitFCN.h"
  20
  21 //-------------------------------------------------------------------------
  22 //                      Class AliBtoJPSItoEleCDFfitHandler
  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(AliBtoJPSItoEleCDFfitHandler)
  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   AliBtoJPSItoEleCDFfitHandler* dummy = (AliBtoJPSItoEleCDFfitHandler *)TVirtualFitter::GetFitter()->GetObjectFit();
  40   dummy->CdfFCN(npar, gin, f, par, iflag);
  41 }
  42
  43
  44 //_________________________________________________________________________________________________
  45 AliBtoJPSItoEleCDFfitHandler::AliBtoJPSItoEleCDFfitHandler():
  46 fIsParamFixed(16),
  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 }
  58 //_________________________________________________________________________________________________
  59 AliBtoJPSItoEleCDFfitHandler::AliBtoJPSItoEleCDFfitHandler(Double_t* decaytime,
  60   Double_t* invariantmass, Int_t ncand) :
  61 fIsParamFixed(16),
  62 fPrintStatus(kFALSE),
  63 fUp(0),
  64 fX(decaytime),
  65 fM(invariantmass),
  66 fLikely(0x0),
  67 fNcand(ncand)
  68 {
  69   //
  70   // constructor
  71   //
  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));
  92 }
  93 //___________________________________________________________________________
  94 AliBtoJPSItoEleCDFfitHandler& AliBtoJPSItoEleCDFfitHandler::operator=(const AliBtoJPSItoEleCDFfitHandler& c)
  95 {
  96   //
  97   // Assignment operator
  98   //
  99   if (this!=&c) {
 100       fIsParamFixed = c.fIsParamFixed;
 101       fPrintStatus  = c.fPrintStatus;
 102       fUp           = c.fUp;
 103       fX            = c.fX;
 104       fM            = c.fM;
 105       fLikely       = c.fLikely;
 106       fNcand        = c.fNcand;
 107      }
 108   return *this;
 109 }
 110
 111 //_______________________________________________________________________________________
 112 AliBtoJPSItoEleCDFfitHandler::AliBtoJPSItoEleCDFfitHandler(const AliBtoJPSItoEleCDFfitHandler& c) :
 113 TNamed(c),
 114 fIsParamFixed(c.fIsParamFixed),
 115 fPrintStatus(c.fPrintStatus),
 116 fUp(c.fUp),
 117 fX(c.fX),
 118 fM(c.fM),
 119 fLikely(c.fLikely),
 120 fNcand(c.fNcand)
 121 {
 122   //
 123   // Copy Constructor
 124   //
 125 }
 126 //_______________________________________________________________________________________
 127 AliBtoJPSItoEleCDFfitHandler::~AliBtoJPSItoEleCDFfitHandler()
 128 {
 129   //
 130   //destructor
 131   //
 132   delete fLikely;
 133 }
 134 //_______________________________________________________________________________________
 135 Int_t AliBtoJPSItoEleCDFfitHandler::DoMinimization()
 136 {
 137   //
 138   // performs the minimization
 139   //
 140   static TVirtualFitter *fitter = TVirtualFitter::Fitter(this,16);
 141   fitter->SetFCN(CDFFunction);
 142
 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);
 162
 163   for(UInt_t indexparam = 0; indexparam < 16; indexparam++){
 164      if(IsParamFixed(indexparam))fitter->FixParameter((Int_t)indexparam);
 165   }
 166
 167   Double_t arglist[2]={10000,1.0};
 168   Int_t iret=fitter->ExecuteCommand("MIGRAD", arglist ,2);
 169   fitter->PrintResults(4,0);
 170
 171   AliInfo("Minimization procedure finished\n");
 172
 173   return iret;
 174 }
 175 //_______________________________________________________________________________________
 176 void AliBtoJPSItoEleCDFfitHandler::CdfFCN(Int_t & /* npar */,
 177               Double_t * /* gin */,Double_t &f,Double_t *par,Int_t /* iflag */)
 178 {
 179 //
 180 // Definition of the FCN to be used by minuit
 181 //
 182   fLikely->SetAllParameters(par);
 183   fLikely->ConvertFromSpherical();
 184   fLikely->ComputeIntegral();
 185   if(fPrintStatus)fLikely->PrintStatus();
 186
 187   TStopwatch t;
 188   t.Start();
 189
 190   f = fLikely->EvaluateLikelihood(fX,fM,fNcand);
 191
 192   t.Stop();
 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));
 196
 197   return;
 198 }
 199 //_______________________________________________________________________________________
 200 void AliBtoJPSItoEleCDFfitHandler::SetParamStartValues(Double_t inputparamvalues[16])
 201 {
 202   for(Int_t index=0; index < 16; index++) fParamStartValues[index] = inputparamvalues[index];
 203 }
 204 //_______________________________________________________________________________________
 205 void AliBtoJPSItoEleCDFfitHandler::SetResolutionConstants()
 206 {
 207   //
 208   // Sets constants for the resolution function
 209   //
 210   fLikely->SetResolutionConstants();
 211
 212 }
 213 //_______________________________________________________________________________________
 214 void AliBtoJPSItoEleCDFfitHandler::SetCrystalBallFunction(Bool_t okCB)
 215 {
 216   //
 217   // Sets the CB as the parametrization for the signal invariant mass spectrum
 218   // (otherwise Landau is chosen)
 219   //
 220   fLikely->SetCrystalBallFunction(okCB);
 221 }
 222 //_______________________________________________________________________________________
 223 void AliBtoJPSItoEleCDFfitHandler::SetMassWndHigh(Double_t limit)
 224 {
 225   //
 226   // Sets upper limit for the invariant mass window (under J/PSI mass peak)
 227   //
 228   fLikely->SetMassWndHigh(limit);
 229 }
 230 //_______________________________________________________________________________________
 231 void AliBtoJPSItoEleCDFfitHandler::SetMassWndLow(Double_t limit)
 232 {
 233   //
 234   // Sets lower limit for the invariant mass window (under J/PSI mass peak)
 235   //
 236   fLikely->SetMassWndLow(limit);
 237 }
 238