PWGDQ/dielectron/AliDielectronBtoJPSItoEleCDFfitFCNfitter.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 "AliDielectronBtoJPSItoEleCDFfitFCNfitter.h"
  16 #include "TArrayD.h"
  17 #include "TFormula.h"
  18 #include "TF1.h"
  19 #include "TF2.h"
  20 #include "TH1F.h"
  21 #include "TH2F.h"
  22 #include "TMath.h"
  23
  24 AliDielectronBtoJPSItoEleCDFfitFCNfitter::AliDielectronBtoJPSItoEleCDFfitFCNfitter():
  25  fInvMass(0x0),
  26  fX(0x0),
  27  fX2D(0x0),
  28  fFCN(0x0),
  29  fFitOpt("")
  30 {
  31  //
  32  // ctor
  33  //
  34  fFCN = new AliDielectronBtoJPSItoEleCDFfitFCN();
  35  fFCN->SetCrystalBallFunction(kTRUE);
  36  for(Int_t i=0; i<kInvMassTotal; i++) fParameterInvMassToFix[i]=kFALSE;
  37  for(Int_t i=0; i<kPseudo; i++)  fParameterXToFix[i]=kFALSE;
  38  for(Int_t i=0; i<kPseudoBkg; i++)  fParameterXbkgToFix[i]=kFALSE;
  39 }
  40 //__________________________________________________________________________________________
  41 AliDielectronBtoJPSItoEleCDFfitFCNfitter::~AliDielectronBtoJPSItoEleCDFfitFCNfitter()
  42 {
  43  //
  44  // dtor
  45  //
  46  if(fInvMass) delete fInvMass;
  47  if(fX) delete fX;
  48  if(fX2D) delete fX2D;
  49  if(fFCN) delete fFCN;
  50 }
  51 //__________________________________________________________________________________________
  52 void AliDielectronBtoJPSItoEleCDFfitFCNfitter::SetInvMassSignalParameters(const Double_t massPar[5])
  53 {
  54  //
  55  // Setter of the inv mass signal distribution parameters
  56  //
  57 fFCN->SetCrystalBallMmean(massPar[0]); //9
  58 fFCN->SetCrystalBallNexp(massPar[1]); //10
  59 fFCN->SetCrystalBallSigma(massPar[2]); //11
  60 fFCN->SetCrystalBallAlpha(massPar[3]); //12
  61 fFCN->SetCrystalBallNorm(massPar[4]); //13
  62  //for(Int_t i=0; i<kInvMassSignal; i++) printf(" par%i %f \n",i,massPar[i]);
  63 }
  64 //__________________________________________________________________________________________
  65 void AliDielectronBtoJPSItoEleCDFfitFCNfitter::SetInvMassParameters(const Double_t massPar[9])
  66 {
  67  //
  68  // Setter of the inv mass total distribution parameters
  69  //
  70  Double_t massParSig[5];
  71  Double_t massParBkg[4];
  72
  73  for(Int_t iPar=0; iPar<9; iPar++){
  74   if(iPar<5) massParSig[iPar]=massPar[iPar];
  75   else massParBkg[iPar-5]=massPar[iPar];
  76  }
  77  SetInvMassSignalParameters(massParSig);
  78  SetInvMassBkgParameters(massParBkg);
  79 }
  80 //__________________________________________________________________________________________
  81 void AliDielectronBtoJPSItoEleCDFfitFCNfitter::SetInvMassBkgParameters(const Double_t massPar[4])
  82 {
  83  //
  84  // Setter of the Inv Mass background distribution parameters
  85  //
  86 fFCN->SetBkgInvMassNorm(massPar[0]);//14
  87 fFCN->SetBkgInvMassMean(massPar[1]);//15
  88 fFCN->SetBkgInvMassSlope(massPar[2]);//16
  89 fFCN->SetBkgInvMassConst(massPar[3]);//17
  90 }
  91 //__________________________________________________________________________________________
  92 void AliDielectronBtoJPSItoEleCDFfitFCNfitter::GetInvMassParameters(TArrayD &massPar)
  93 {
  94  //
  95  // Getter of the inv mass total distribution parameters
  96  //
  97  massPar.Reset();
  98  massPar.Set(kInvMassTotal);
  99  massPar.AddAt(fFCN->GetCrystalBallMmean(),0); //9
 100  massPar.AddAt(fFCN->GetCrystalBallNexp(),1); //10
 101  massPar.AddAt(fFCN->GetCrystalBallSigma(),2); //11
 102  massPar.AddAt(fFCN->GetCrystalBallAlpha(),3); //12
 103  massPar.AddAt(fFCN->GetCrystalBallNorm(),4); //13
 104  massPar.AddAt(fFCN->GetBkgInvMassNorm(),5);//14
 105  massPar.AddAt(fFCN->GetBkgInvMassMean(),6);//15
 106  massPar.AddAt(fFCN->GetBkgInvMassSlope(),7); //16
 107  massPar.AddAt(fFCN->GetBkgInvMassConst(),8); //17
 108 }
 109 //__________________________________________________________________________________________
 110 void AliDielectronBtoJPSItoEleCDFfitFCNfitter::GetInvMassSignalParameters(TArrayD &massPar)
 111 {
 112  //
 113  // Getter of the Inv Mass signal distribution parameters
 114  //
 115  massPar.Reset();
 116  massPar.Set(kInvMassSignal);
 117  massPar.AddAt(fFCN->GetCrystalBallMmean(),0); //9
 118  massPar.AddAt(fFCN->GetCrystalBallNexp(),1); //10
 119  massPar.AddAt(fFCN->GetCrystalBallSigma(),2); //11
 120  massPar.AddAt(fFCN->GetCrystalBallAlpha(),3); //12
 121  massPar.AddAt(fFCN->GetCrystalBallNorm(),4); //13
 122
 123 }
 124 //__________________________________________________________________________________________
 125 void AliDielectronBtoJPSItoEleCDFfitFCNfitter::GetInvMassBkgParameters(TArrayD &massPar)
 126 {
 127  //
 128  // Getter of the Inv Mass background distribution parameters
 129  //
 130  massPar.Reset();
 131  massPar.Set(kInvMassBkg);
 132  massPar.AddAt(fFCN->GetBkgInvMassNorm(),0); //14
 133  massPar.AddAt(fFCN->GetBkgInvMassMean(),1); //15
 134  massPar.AddAt(fFCN->GetCrystalBallSigma(),2); //11
 135  massPar.AddAt(fFCN->GetCrystalBallAlpha(),3); //12
 136 }
 137
 138 //__________________________________________________________________________________________
 139 void AliDielectronBtoJPSItoEleCDFfitFCNfitter::SetPseudoProperBkgParameters(Double_t xBkgPar[7])
 140 {
 141  //
 142  // Setter of the psudoproper background distribution parameters
 143  //
 144  fFCN->SetResWeight(xBkgPar[0]);//0
 145  fFCN->SetFPlus(xBkgPar[1]);//1
 146  fFCN->SetFMinus(xBkgPar[2]);//2
 147  fFCN->SetFSym(xBkgPar[3]);//3
 148  fFCN->SetLamPlus(xBkgPar[4]); //4
 149  fFCN->SetLamMinus(xBkgPar[5]); // 5
 150  fFCN->SetLamSym(xBkgPar[6]); // 6
 151
 152 }
 153 //__________________________________________________________________________________________
 154 void AliDielectronBtoJPSItoEleCDFfitFCNfitter::GetPseudoProperBkgParameters(TArrayD &xBkgPar)
 155 {
 156  //
 157  // Getter of the psudoproper background distribution parameters
 158  //
 159  xBkgPar.Reset();
 160  xBkgPar.Set(kPseudoBkg);
 161  xBkgPar.AddAt(fFCN->GetResWeight(),0);
 162  xBkgPar.AddAt(fFCN->GetFPlus(),1);
 163  xBkgPar.AddAt(fFCN->GetFMinus(),2);
 164  xBkgPar.AddAt(fFCN->GetFSym(),3);
 165  xBkgPar.AddAt(fFCN->GetLamPlus(),4);
 166  xBkgPar.AddAt(fFCN->GetLamMinus(),5);
 167  xBkgPar.AddAt(fFCN->GetLamSym(),6);
 168 }
 169 //__________________________________________________________________________________________
 170 TH1F* AliDielectronBtoJPSItoEleCDFfitFCNfitter::FitInvMassSignal(Double_t norm,Double_t mMin, Double_t mMax)
 171 {
 172  //
 173  // Fit on the invariant mass signal distribution ( Crystal ball params )
 174  //
 175  if(!fInvMass){
 176   printf("no histogram...exiting\n");
 177   return 0x0;
 178  }
 179
 180  TArrayD pars;
 181  GetInvMassSignalParameters(pars);
 182
 183  TF1 * invMass = new TF1("invMassSignal",this,&AliDielectronBtoJPSItoEleCDFfitFCNfitter::CDFInvMassSignal,mMin,mMax,kInvMassSignal+1);
 184  for(Int_t i=0; i<kInvMassSignal ; i++) {
 185   if(fParameterInvMassToFix[i]) invMass->FixParameter(i,pars.At(i));
 186   else invMass->SetParameter(i,pars.At(i));
 187  }
 188  invMass->SetParameter(kInvMassSignal,norm);
 189  fInvMass->Fit("invMassSignal",Form("R%s",fFitOpt.Data()),"",mMin,mMax);
 190  return fInvMass;
 191 }
 192 //__________________________________________________________________________________________
 193 TH1F * AliDielectronBtoJPSItoEleCDFfitFCNfitter::FitInvMass(Double_t norm[2],Double_t mMin, Double_t mMax)
 194 {
 195  //
 196  // Fit on the invariant mass total distribution
 197  //
 198  if(!fInvMass){
 199   printf("no histogram...exiting\n");
 200   return 0x0;
 201  }
 202
 203  TArrayD pars;
 204  GetInvMassParameters(pars);
 205
 206  TF1 * invMassTot = new TF1("invMassTot",this,&AliDielectronBtoJPSItoEleCDFfitFCNfitter::CDFInvMassTotal,mMin,mMax,11);
 207
 208  for(Int_t i=0; i<kInvMassTotal ; i++) {
 209   if(fParameterInvMassToFix[i]) invMassTot->FixParameter(i,pars.At(i));
 210   else invMassTot->SetParameter(i,pars.At(i));
 211  }
 212
 213  if(!fParameterInvMassToFix[5]) invMassTot->SetParLimits(5,0,999999999);
 214
 215
 216  invMassTot->SetParameter(9,norm[0]);
 217  invMassTot->SetParLimits(9,0,999999999);
 218  invMassTot->SetParameter(10,norm[1]);
 219  invMassTot->SetParLimits(10,0,999999999);
 220
 221  fInvMass->Fit("invMassTot",Form("R%s",fFitOpt.Data()),"",mMin,mMax);
 222
 223  //printf("Inv Mass Params : ");
 224  //for(Int_t i=0; i< 11; i++) printf("p%i %f ",i,invMassTot->GetParameter(i));
 225  //printf("\n");
 226
 227  return fInvMass;
 228 }
 229 //__________________________________________________________________________________________
 230 void AliDielectronBtoJPSItoEleCDFfitFCNfitter::GetPseudoProperParameters(TArrayD &xPar, Int_t type)
 231 {
 232  //
 233  // Getter of the psudoproper background distribution parameters per dielectron type (FF, FS, SS)
 234  //
 235  xPar.Reset();
 236  xPar.Set(kPseudo);
 237  xPar.AddAt(fFCN->GetNormGaus1ResFunc(type),0);
 238  xPar.AddAt(fFCN->GetNormGaus2ResFunc(type),3);
 239  xPar.AddAt(fFCN->GetResMean1(type),1);
 240  xPar.AddAt(fFCN->GetResSigma1(type),2);
 241  xPar.AddAt(fFCN->GetResMean2(type),4);
 242  xPar.AddAt(fFCN->GetResSigma2(type),5);
 243  xPar.AddAt(fFCN->GetResAlfa(type),6);
 244  xPar.AddAt(fFCN->GetResLambda(type),7);
 245  xPar.AddAt(fFCN->GetResNormExp(type),8);
 246 }
 247 //__________________________________________________________________________________________
 248 TH1F * AliDielectronBtoJPSItoEleCDFfitFCNfitter::FitResolutionFunction(Double_t xMin, Double_t xMax, Int_t type, Double_t norm)
 249 {
 250  //
 251  // Fit on the resolution function distribution per dielectron type
 252  //
 253  if(!fX){
 254   printf("no Pseudoproper decayhistogram...exiting\n");
 255   return 0x0;
 256  }
 257
 258  TArrayD pars;
 259  GetPseudoProperParameters(pars,type);
 260
 261  TF1 * xPrompt = new TF1("xPrompt",this,&AliDielectronBtoJPSItoEleCDFfitFCNfitter::CDFResolutionFunction,xMin,xMax,11);
 262  for(Int_t i=0; i<kPseudo ; i++) {
 263   if(fParameterXToFix[i]) xPrompt->FixParameter(i,pars.At(i));
 264   else xPrompt->SetParameter(i,pars.At(i));
 265   printf("%f ",pars.At(i));
 266  }
 267  printf("\n");
 268  if(!fParameterXToFix[3]) xPrompt->SetParLimits(3,0,999999999);
 269  xPrompt->FixParameter(9,type);
 270  xPrompt->SetParameter(10,norm);
 271  fX->Fit("xPrompt",Form("R%s",fFitOpt.Data()),"",xMin,xMax);
 272
 273 // printf("Res Func type %i : ",type);
 274 // for(Int_t i=0; i< 11; i++) printf("p%i %6.3f ",i,xPrompt->GetParameter(i));
 275 // printf("\n");
 276  return fX;
 277 }
 278
 279 //__________________________________________________________________________________________
 280 TH2F * AliDielectronBtoJPSItoEleCDFfitFCNfitter::FitBkgPsudoProper(Double_t xMin, Double_t xMax, Double_t norm)
 281 {
 282  //
 283  // Fit on the resolution function distribution of the background ( in 2 D, whose one D is the dielectron
 284  // type)
 285  //
 286
 287  if(!fX2D){
 288   printf("no BKG Pseudoproper decay 2D histogram ...exiting\n");
 289   return 0x0;
 290  }
 291
 292  TArrayD bkgParam;
 293  GetPseudoProperBkgParameters(bkgParam);
 294
 295  TF2 *psProperBkgnd = new TF2("bkgPseudoProper",this,
 296    &AliDielectronBtoJPSItoEleCDFfitFCNfitter::PsProperBackFitFunc,xMin,xMax,-0.5,2.5,8,"AliDielectronBtoJPSItoEleCDFfitFCNfitter","PsProperBackFitFunc");
 297
 298  for(Int_t i=0; i< kPseudoBkg ; i++) {
 299   if(fParameterXbkgToFix[i]) psProperBkgnd->FixParameter(i,bkgParam.At(i));
 300   else psProperBkgnd->SetParameter(i,bkgParam.At(i));
 301  }
 302  psProperBkgnd->SetParameter(kPseudoBkg,norm);
 303
 304  fX2D->Fit("bkgPseudoProper",Form("R%s",fFitOpt.Data()),"LEGO2",xMin,xMax);
 305
 306  return fX2D;
 307 }
 308
 309 //__________________________________________________________________________________________
 310 Double_t AliDielectronBtoJPSItoEleCDFfitFCNfitter::CDFInvMassSignal(const Double_t *x,const Double_t *par)
 311 {
 312  //
 313  // function
 314  //
 315  SetInvMassSignalParameters(par);
 316  return par[5]*fFCN->EvaluateCDFInvMassSigDistr(x[0]);
 317 }
 318 //__________________________________________________________________________________________
 319 Double_t AliDielectronBtoJPSItoEleCDFfitFCNfitter::CDFInvMassBkg(const Double_t *x,const Double_t *par)
 320 {
 321  //
 322  // function
 323  //
 324  SetInvMassBkgParameters(par);
 325  return fFCN->EvaluateCDFInvMassBkgDistr(x[0]);
 326 }
 327 //__________________________________________________________________________________________
 328 Double_t AliDielectronBtoJPSItoEleCDFfitFCNfitter::CDFInvMassTotal(const Double_t *x,const Double_t *par)
 329 {
 330  //
 331  // function (it doesn't use the function integral->pay attention to it)
 332  //
 333  SetInvMassParameters(par);
 334  //return ((EvaluateCDFInvMassSigDistr(x[0])))*par[9] + (EvaluateCDFInvMassBkgDistr(x[0]))*par[10];
 335  return (((fFCN->EvaluateCDFInvMassSigDistr(x[0])))*par[9] + (1-par[9])*(fFCN->EvaluateCDFInvMassBkgDistr(x[0])))*par[10];
 336 }
 337 //__________________________________________________________________________________________
 338 Double_t AliDielectronBtoJPSItoEleCDFfitFCNfitter::CDFResolutionFunction(const Double_t *x,const Double_t *par)
 339 {
 340  //
 341  // function
 342  //
 343 //  Double_t resParam[kPseudo];
 344 //  for(Int_t i=0; i<kPseudo; i++) {
 345 //   resParam[i] = par[i];
 346 //  }
 347  fFCN->SetResolutionConstants(par,(Int_t)par[9]);
 348  return par[10]*fFCN->ResolutionFunc(x[0],-1.,(Int_t)(par[9]));
 349 }
 350 //__________________________________________________________________________________________
 351 Double_t AliDielectronBtoJPSItoEleCDFfitFCNfitter::PsProperBackFitFunc(const Double_t* x, const Double_t* par)
 352 {
 353  //
 354  // function
 355  //
 356  Int_t type = TMath::Nint(x[1]);
 357  Double_t params[7];
 358  for(Int_t i=0; i<7; i++) params[i]=par[i];
 359  SetPseudoProperBkgParameters(params);
 360  return fFCN->EvaluateCDFDecayTimeBkgDistr(x[0],type)*par[7]*fFCN->GetResWeight(type);
 361 }