PWGGA/GammaConv/AliConversionMesonCuts.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Authors: Svein Lindal, Daniel Lohner                                   *
   5  * Version 1.0                                                            *
   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
  16 ////////////////////////////////////////////////
  17 //---------------------------------------------
  18 // Class handling all kinds of selection cuts for
  19 // Gamma Conversion analysis
  20 //---------------------------------------------
  21 ////////////////////////////////////////////////
  22
  23 #include "AliConversionMesonCuts.h"
  24
  25 #include "AliKFVertex.h"
  26 #include "AliAODTrack.h"
  27 #include "AliESDtrack.h"
  28 #include "AliAnalysisManager.h"
  29 #include "AliInputEventHandler.h"
  30 #include "AliMCEventHandler.h"
  31 #include "AliAODHandler.h"
  32 #include "AliPIDResponse.h"
  33 #include "TH1.h"
  34 #include "TH2.h"
  35 #include "AliStack.h"
  36 #include "AliAODConversionMother.h"
  37 #include "TObjString.h"
  38 #include "AliAODEvent.h"
  39 #include "AliESDEvent.h"
  40 #include "AliCentrality.h"
  41 #include "TList.h"
  42 class iostream;
  43
  44 using namespace std;
  45
  46 ClassImp(AliConversionMesonCuts)
  47
  48
  49 const char* AliConversionMesonCuts::fgkCutNames[AliConversionMesonCuts::kNCuts] = {
  50    "MesonKind", //0
  51    "BackgroundScheme", //1
  52    "NumberOfBGEvents", //2
  53    "DegreesForRotationMethod", //3
  54    "RapidityMesonCut", //4
  55    "RCut", //5
  56    "AlphaMesonCut", //6
  57    "Chi2MesonCut", //7
  58    "SharedElectronCuts", //8
  59    "RejectToCloseV0s", //9
  60    "UseMCPSmearing", //10
  61 };
  62
  63
  64 //________________________________________________________________________
  65 AliConversionMesonCuts::AliConversionMesonCuts(const char *name,const char *title) :
  66    AliAnalysisCuts(name,title),
  67    fHistograms(NULL),
  68    fMesonKind(0),
  69    fMaxR(200),
  70    fChi2CutMeson(1000),
  71    fAlphaMinCutMeson(0),
  72    fAlphaCutMeson(1),
  73    fRapidityCutMeson(1),
  74    fUseRotationMethodInBG(kFALSE),
  75    fDoBG(kTRUE),
  76    fdoBGProbability(kFALSE),
  77    fUseTrackMultiplicityForBG(kFALSE),
  78    fnDegreeRotationPMForBG(0),
  79    fNumberOfBGEvents(0),
  80    fOpeningAngle(0.005),
  81    fDoToCloseV0sCut(kFALSE),
  82    fminV0Dist(200.),
  83    fDoSharedElecCut(kFALSE),
  84    fUseMCPSmearing(kFALSE),
  85    fPBremSmearing(0),
  86    fPSigSmearing(0),
  87    fPSigSmearingCte(0),
  88    fBrem(NULL),
  89    fRandom(0),
  90    fElectronLabelArraySize(500),
  91    fElectronLabelArray(NULL),
  92    fBackgroundHandler(0),
  93    fCutString(NULL),
  94    hMesonCuts(NULL),
  95    hMesonBGCuts(NULL)
  96
  97 {
  98    for(Int_t jj=0;jj<kNCuts;jj++){fCuts[jj]=0;}
  99    fCutString=new TObjString((GetCutNumber()).Data());
 100    fElectronLabelArray = new Int_t[fElectronLabelArraySize];
 101    if (fBrem == NULL){
 102       fBrem = new TF1("fBrem","pow(-log(x),[0]/log(2.0)-1.0)/TMath::Gamma([0]/log(2.0))",0.00001,0.999999999);
 103       // tests done with 1.0e-14
 104       fBrem->SetParameter(0,fPBremSmearing);
 105       fBrem->SetNpx(100000);
 106    }
 107
 108 }
 109
 110 //________________________________________________________________________
 111 AliConversionMesonCuts::AliConversionMesonCuts(const AliConversionMesonCuts &ref) :
 112    AliAnalysisCuts(ref),
 113    fHistograms(NULL),
 114    fMesonKind(ref.fMesonKind),
 115    fMaxR(ref.fMaxR),
 116    fChi2CutMeson(ref.fChi2CutMeson),
 117    fAlphaMinCutMeson(ref.fAlphaMinCutMeson),
 118    fAlphaCutMeson(ref.fAlphaCutMeson),
 119    fRapidityCutMeson(ref.fRapidityCutMeson),
 120    fUseRotationMethodInBG(ref.fUseRotationMethodInBG),
 121    fDoBG(ref.fDoBG),
 122    fdoBGProbability(ref.fdoBGProbability),
 123    fUseTrackMultiplicityForBG(ref.fUseTrackMultiplicityForBG),
 124    fnDegreeRotationPMForBG(ref.fnDegreeRotationPMForBG),
 125    fNumberOfBGEvents(ref. fNumberOfBGEvents),
 126    fOpeningAngle(ref.fOpeningAngle),
 127    fDoToCloseV0sCut(ref.fDoToCloseV0sCut),
 128    fminV0Dist(ref.fminV0Dist),
 129    fDoSharedElecCut(ref.fDoSharedElecCut),
 130    fUseMCPSmearing(ref.fUseMCPSmearing),
 131    fPBremSmearing(ref.fPBremSmearing),
 132    fPSigSmearing(ref.fPSigSmearing),
 133    fPSigSmearingCte(ref.fPSigSmearingCte),
 134    fBrem(NULL),
 135    fRandom(ref.fRandom),
 136    fElectronLabelArraySize(ref.fElectronLabelArraySize),
 137    fElectronLabelArray(NULL),
 138    fBackgroundHandler(ref.fBackgroundHandler),
 139    fCutString(NULL),
 140    hMesonCuts(NULL),
 141    hMesonBGCuts(NULL)
 142
 143 {
 144    // Copy Constructor
 145    for(Int_t jj=0;jj<kNCuts;jj++){fCuts[jj]=ref.fCuts[jj];}
 146    fCutString=new TObjString((GetCutNumber()).Data());
 147    fElectronLabelArray = new Int_t[fElectronLabelArraySize];
 148    if (fBrem == NULL)fBrem = (TF1*)ref.fBrem->Clone("fBrem");
 149    // Histograms are not copied, if you need them, call InitCutHistograms
 150 }
 151
 152
 153 //________________________________________________________________________
 154 AliConversionMesonCuts::~AliConversionMesonCuts() {
 155    // Destructor
 156    //Deleting fHistograms leads to seg fault it it's added to output collection of a task
 157    // if(fHistograms)
 158    //   delete fHistograms;
 159    // fHistograms = NULL;
 160    if(fCutString != NULL){
 161       delete fCutString;
 162       fCutString = NULL;
 163    }
 164    if(fElectronLabelArray){
 165       delete fElectronLabelArray;
 166       fElectronLabelArray = NULL;
 167    }
 168
 169 }
 170
 171 //________________________________________________________________________
 172 void AliConversionMesonCuts::InitCutHistograms(TString name){
 173
 174    // Initialize Cut Histograms for QA (only initialized and filled if function is called)
 175    TH1::AddDirectory(kFALSE);
 176
 177    if(fHistograms != NULL){
 178       delete fHistograms;
 179       fHistograms=NULL;
 180    }
 181    if(fHistograms==NULL){
 182       fHistograms=new TList();
 183       fHistograms->SetOwner(kTRUE);
 184       if(name=="")fHistograms->SetName(Form("ConvMesonCuts_%s",GetCutNumber().Data()));
 185       else fHistograms->SetName(Form("%s_%s",name.Data(),GetCutNumber().Data()));
 186    }
 187
 188    // Meson Cuts
 189    hMesonCuts=new TH1F(Form("MesonCuts %s",GetCutNumber().Data()),"MesonCuts",10,-0.5,9.5);
 190    hMesonCuts->GetXaxis()->SetBinLabel(1,"in");
 191    hMesonCuts->GetXaxis()->SetBinLabel(2,"undef rapidity");
 192    hMesonCuts->GetXaxis()->SetBinLabel(3,"rapidity cut");
 193    hMesonCuts->GetXaxis()->SetBinLabel(4,"opening angle");
 194    hMesonCuts->GetXaxis()->SetBinLabel(5,"alpha max");
 195    hMesonCuts->GetXaxis()->SetBinLabel(6,"alpha min");
 196    hMesonCuts->GetXaxis()->SetBinLabel(7,"out");
 197    fHistograms->Add(hMesonCuts);
 198
 199    hMesonBGCuts=new TH1F(Form("MesonBGCuts %s",GetCutNumber().Data()),"MesonBGCuts",10,-0.5,9.5);
 200    hMesonBGCuts->GetXaxis()->SetBinLabel(1,"in");
 201    hMesonBGCuts->GetXaxis()->SetBinLabel(2,"undef rapidity");
 202    hMesonBGCuts->GetXaxis()->SetBinLabel(3,"rapidity cut");
 203    hMesonBGCuts->GetXaxis()->SetBinLabel(4,"opening angle");
 204    hMesonBGCuts->GetXaxis()->SetBinLabel(5,"alpha max");
 205    hMesonBGCuts->GetXaxis()->SetBinLabel(6,"alpha min");
 206    hMesonBGCuts->GetXaxis()->SetBinLabel(7,"out");
 207    fHistograms->Add(hMesonBGCuts);
 208
 209    TH1::AddDirectory(kTRUE);
 210 }
 211
 212 //________________________________________________________________________
 213 Bool_t AliConversionMesonCuts::MesonIsSelectedMC(TParticle *fMCMother,AliStack *fMCStack, Double_t fRapidityShift){
 214    // Returns true for all pions within acceptance cuts for decay into 2 photons
 215    // If bMCDaughtersInAcceptance is selected, it requires in addition that both daughter photons are within acceptance cuts
 216
 217    if(!fMCStack)return kFALSE;
 218
 219    if(fMCMother->GetPdgCode()==111 || fMCMother->GetPdgCode()==221){
 220       if(fMCMother->R()>fMaxR)  return kFALSE; // cuts on distance from collision point
 221
 222       Double_t rapidity = 10.;
 223       if(fMCMother->Energy() - fMCMother->Pz() == 0 || fMCMother->Energy() + fMCMother->Pz() == 0){
 224          rapidity=8.-fRapidityShift;
 225       } else{
 226          rapidity = 0.5*(TMath::Log((fMCMother->Energy()+fMCMother->Pz()) / (fMCMother->Energy()-fMCMother->Pz())))-fRapidityShift;
 227       }
 228
 229       // Rapidity Cut
 230       if(abs(rapidity)>fRapidityCutMeson)return kFALSE;
 231
 232       // Select only -> 2y decay channel
 233       if(fMCMother->GetNDaughters()!=2)return kFALSE;
 234
 235       for(Int_t i=0;i<2;i++){
 236          TParticle *MDaughter=fMCStack->Particle(fMCMother->GetDaughter(i));
 237          // Is Daughter a Photon?
 238          if(MDaughter->GetPdgCode()!=22)return kFALSE;
 239          // Is Photon in Acceptance?
 240          //   if(bMCDaughtersInAcceptance){
 241          //     if(!PhotonIsSelectedMC(MDaughter,fMCStack)){return kFALSE;}
 242          //   }
 243       }
 244       return kTRUE;
 245    }
 246    return kFALSE;
 247 }
 248
 249 //________________________________________________________________________
 250 Bool_t AliConversionMesonCuts::MesonIsSelectedMCDalitz(TParticle *fMCMother,AliStack *fMCStack, Double_t fRapidityShift){
 251    // Returns true for all pions within acceptance cuts for decay into 2 photons
 252    // If bMCDaughtersInAcceptance is selected, it requires in addition that both daughter photons are within acceptance cuts
 253
 254    if(!fMCStack)return kFALSE;
 255
 256    if(fMCMother->GetPdgCode()==111 || fMCMother->GetPdgCode()==221){
 257
 258       if(fMCMother->R()>fMaxR)  return kFALSE; // cuts on distance from collision point
 259
 260       Double_t rapidity = 10.;
 261       if(fMCMother->Energy() - fMCMother->Pz() == 0 || fMCMother->Energy() + fMCMother->Pz() == 0){
 262          rapidity=8.-fRapidityShift;
 263       }
 264       else{
 265          rapidity = 0.5*(TMath::Log((fMCMother->Energy()+fMCMother->Pz()) / (fMCMother->Energy()-fMCMother->Pz())))-fRapidityShift;
 266       }
 267
 268       // Rapidity Cut
 269       if(abs(rapidity)>fRapidityCutMeson)return kFALSE;
 270
 271       // Select only -> Dalitz decay channel
 272       if(fMCMother->GetNDaughters()!=3)return kFALSE;
 273
 274       Int_t daughterPDGs[3] = {0,0,0};
 275       Int_t index = 0;
 276
 277       //                iParticle->GetFirstDaughter(); idxPi0 <= iParticle->GetLastDaughter()
 278
 279       for(Int_t i=fMCMother->GetFirstDaughter(); i<= fMCMother->GetLastDaughter();i++){
 280          TParticle *MDaughter=fMCStack->Particle(i);
 281          // Is Daughter a Photon or an electron?
 282          daughterPDGs[index]=MDaughter->GetPdgCode();
 283          index++;
 284       }
 285       for (Int_t j=0;j<2;j++){
 286
 287          for (Int_t i=0;i<2;i++){
 288             if (daughterPDGs[i] > daughterPDGs[i+1]){
 289                Int_t interMed = daughterPDGs[i] ;
 290                daughterPDGs[i] = daughterPDGs[i+1];
 291                daughterPDGs[i+1] = interMed;
 292             }
 293          }
 294       }
 295       if (daughterPDGs[0] == -11 && daughterPDGs[1] == 11 && daughterPDGs[2] == 22) return kTRUE;
 296    }
 297    return kFALSE;
 298 }
 299
 300
 301 ///________________________________________________________________________
 302 Bool_t AliConversionMesonCuts::MesonIsSelected(AliAODConversionMother *pi0,Bool_t IsSignal, Double_t fRapidityShift)
 303 {
 304    // Selection of reconstructed Meson candidates
 305    // Use flag IsSignal in order to fill Fill different
 306    // histograms for Signal and Background
 307    TH1 *hist=0x0;
 308
 309    if(IsSignal){hist=hMesonCuts;}
 310    else{hist=hMesonBGCuts;}
 311
 312    Int_t cutIndex=0;
 313
 314    if(hist)hist->Fill(cutIndex);
 315    cutIndex++;
 316
 317    // Undefined Rapidity -> Floating Point exception
 318    if((pi0->E()+pi0->Pz())/(pi0->E()-pi0->Pz())<=0){
 319       if(hist)hist->Fill(cutIndex);
 320       cutIndex++;
 321       return kFALSE;
 322    }
 323    else{
 324       // PseudoRapidity Cut --> But we cut on Rapidity !!!
 325       cutIndex++;
 326       if(abs(pi0->Rapidity()-fRapidityShift)>fRapidityCutMeson){
 327          if(hist)hist->Fill(cutIndex);
 328          return kFALSE;
 329       }
 330    }
 331    cutIndex++;
 332
 333    // Opening Angle Cut
 334    //fOpeningAngle=2*TMath::ATan(0.134/pi0->P());// physical minimum opening angle
 335    if( pi0->GetOpeningAngle() < fOpeningAngle){
 336       if(hist)hist->Fill(cutIndex);
 337       return kFALSE;
 338    }
 339    cutIndex++;
 340
 341    // Alpha Max Cut
 342    if(pi0->GetAlpha()>fAlphaCutMeson){
 343       if(hist)hist->Fill(cutIndex);
 344       return kFALSE;
 345    }
 346    cutIndex++;
 347
 348    // Alpha Min Cut
 349    if(pi0->GetAlpha()<fAlphaMinCutMeson){
 350       if(hist)hist->Fill(cutIndex);
 351       return kFALSE;
 352    }
 353    cutIndex++;
 354
 355    if(hist)hist->Fill(cutIndex);
 356    return kTRUE;
 357 }
 358
 359
 360
 361 ///________________________________________________________________________
 362 ///________________________________________________________________________
 363 Bool_t AliConversionMesonCuts::UpdateCutString() {
 364    ///Update the cut string (if it has been created yet)
 365
 366    if(fCutString && fCutString->GetString().Length() == kNCuts) {
 367       fCutString->SetString(GetCutNumber());
 368    } else {
 369       return kFALSE;
 370    }
 371    return kTRUE;
 372 }
 373
 374 ///________________________________________________________________________
 375 Bool_t AliConversionMesonCuts::InitializeCutsFromCutString(const TString analysisCutSelection ) {
 376    // Initialize Cuts from a given Cut string
 377    AliInfo(Form("Set Meson Cutnumber: %s",analysisCutSelection.Data()));
 378    if(analysisCutSelection.Length()!=kNCuts) {
 379       AliError(Form("Cut selection has the wrong length! size is %d, number of cuts is %d", analysisCutSelection.Length(), kNCuts));
 380       return kFALSE;
 381    }
 382    if(!analysisCutSelection.IsDigit()){
 383       AliError("Cut selection contains characters");
 384       return kFALSE;
 385    }
 386
 387    const char *cutSelection = analysisCutSelection.Data();
 388 #define ASSIGNARRAY(i)  fCuts[i] = cutSelection[i] - '0'
 389    for(Int_t ii=0;ii<kNCuts;ii++){
 390       ASSIGNARRAY(ii);
 391    }
 392
 393    // Set Individual Cuts
 394    for(Int_t ii=0;ii<kNCuts;ii++){
 395       if(!SetCut(cutIds(ii),fCuts[ii]))return kFALSE;
 396    }
 397
 398    //PrintCuts();
 399    return kTRUE;
 400 }
 401 ///________________________________________________________________________
 402 Bool_t AliConversionMesonCuts::SetCut(cutIds cutID, const Int_t value) {
 403    ///Set individual cut ID
 404
 405    //cout << "Updating cut  " << fgkCutNames[cutID] << " (" << cutID << ") to " << value << endl;
 406    switch (cutID) {
 407    case kMesonKind:
 408       if( SetMesonKind(value)) {
 409          fCuts[kMesonKind] = value;
 410          UpdateCutString();
 411          return kTRUE;
 412       } else return kFALSE;
 413    case kchi2MesonCut:
 414       if( SetChi2MesonCut(value)) {
 415          fCuts[kchi2MesonCut] = value;
 416          UpdateCutString();
 417          return kTRUE;
 418       } else return kFALSE;
 419    case kalphaMesonCut:
 420       if( SetAlphaMesonCut(value)) {
 421          fCuts[kalphaMesonCut] = value;
 422          UpdateCutString();
 423          return kTRUE;
 424       } else return kFALSE;
 425    case kRCut:
 426       if( SetRCut(value)) {
 427          fCuts[kRCut] = value;
 428          UpdateCutString();
 429          return kTRUE;
 430       } else return kFALSE;
 431
 432    case kRapidityMesonCut:
 433       if( SetRapidityMesonCut(value)) {
 434          fCuts[kRapidityMesonCut] = value;
 435          UpdateCutString();
 436          return kTRUE;
 437       } else return kFALSE;
 438
 439    case kBackgroundScheme:
 440       if( SetBackgroundScheme(value)) {
 441          fCuts[kBackgroundScheme] = value;
 442          UpdateCutString();
 443          return kTRUE;
 444       } else return kFALSE;
 445
 446    case kDegreesForRotationMethod:
 447       if( SetNDegreesForRotationMethod(value)) {
 448          fCuts[kDegreesForRotationMethod] = value;
 449          UpdateCutString();
 450          return kTRUE;
 451       } else return kFALSE;
 452
 453    case kNumberOfBGEvents:
 454       if( SetNumberOfBGEvents(value)) {
 455          fCuts[kNumberOfBGEvents] = value;
 456          UpdateCutString();
 457          return kTRUE;
 458       } else return kFALSE;
 459
 460    case kuseMCPSmearing:
 461       if( SetMCPSmearing(value)) {
 462          fCuts[kuseMCPSmearing] = value;
 463          UpdateCutString();
 464          return kTRUE;
 465       } else return kFALSE;
 466    case kElecShare:
 467       if( SetSharedElectronCut(value)) {
 468          fCuts[kElecShare] = value;
 469          UpdateCutString();
 470          return kTRUE;
 471       } else return kFALSE;
 472    case kToCloseV0s:
 473       if( SetToCloseV0sCut(value)) {
 474          fCuts[kToCloseV0s] = value;
 475          UpdateCutString();
 476          return kTRUE;
 477       } else return kFALSE;
 478
 479    case kNCuts:
 480       cout << "Error:: Cut id out of range"<< endl;
 481       return kFALSE;
 482    }
 483
 484    cout << "Error:: Cut id " << cutID << " not recognized "<< endl;
 485    return kFALSE;
 486
 487 }
 488
 489
 490 ///________________________________________________________________________
 491 void AliConversionMesonCuts::PrintCuts() {
 492    // Print out current Cut Selection
 493    for(Int_t ic = 0; ic < kNCuts; ic++) {
 494       printf("%-30s : %d \n", fgkCutNames[ic], fCuts[ic]);
 495    }
 496 }
 497
 498 ///________________________________________________________________________
 499 Bool_t AliConversionMesonCuts::SetMesonKind(Int_t mesonKind){
 500    // Set Cut
 501    switch(mesonKind){
 502    case 0:
 503       fMesonKind=0;
 504       break;
 505    case 1:
 506       fMesonKind=1;;
 507       break;
 508    default:
 509       cout<<"Warning: Meson kind not defined"<<mesonKind<<endl;
 510       return kFALSE;
 511    }
 512    return kTRUE;
 513 }
 514
 515 ///________________________________________________________________________
 516 Bool_t AliConversionMesonCuts::SetRCut(Int_t rCut){
 517    // Set Cut
 518    switch(rCut){
 519    case 0:
 520       fMaxR = 180.;
 521       break;
 522    case 1:
 523       fMaxR = 180.;
 524       break;
 525    case 2:
 526       fMaxR = 180.;
 527       break;
 528    case 3:
 529       fMaxR = 70.;
 530       break;
 531    case 4:
 532       fMaxR = 70.;
 533       break;
 534    case 5:
 535       fMaxR = 180.;
 536       break;
 537       // High purity cuts for PbPb
 538    case 9:
 539       fMaxR = 180.;
 540       break;
 541    default:
 542       cout<<"Warning: rCut not defined"<<rCut<<endl;
 543       return kFALSE;
 544    }
 545    return kTRUE;
 546 }
 547
 548 ///________________________________________________________________________
 549 Bool_t AliConversionMesonCuts::SetChi2MesonCut(Int_t chi2MesonCut){
 550    // Set Cut
 551    switch(chi2MesonCut){
 552    case 0:  // 100.
 553       fChi2CutMeson = 100.;
 554       break;
 555    case 1:  // 50.
 556       fChi2CutMeson = 50.;
 557       break;
 558    case 2:  // 30.
 559       fChi2CutMeson = 30.;
 560       break;
 561    case 3:
 562       fChi2CutMeson = 200.;
 563       break;
 564    case 4:
 565       fChi2CutMeson = 500.;
 566       break;
 567    case 5:
 568       fChi2CutMeson = 1000.;
 569       break;
 570    default:
 571       cout<<"Warning: Chi2MesonCut not defined "<<chi2MesonCut<<endl;
 572       return kFALSE;
 573    }
 574    return kTRUE;
 575 }
 576
 577
 578 ///________________________________________________________________________
 579 Bool_t AliConversionMesonCuts::SetAlphaMesonCut(Int_t alphaMesonCut)
 580 {   // Set Cut
 581    switch(alphaMesonCut){
 582    case 0:      // 0- 0.7
 583       fAlphaMinCutMeson  = 0.0;
 584       fAlphaCutMeson     = 0.7;
 585       break;
 586    case 1:      // 0-0.5
 587       fAlphaMinCutMeson  = 0.0;
 588       fAlphaCutMeson     = 0.5;
 589       break;
 590    case 2:      // 0.5-1
 591       fAlphaMinCutMeson  = 0.5;
 592       fAlphaCutMeson     = 1.;
 593       break;
 594    case 3:      // 0.0-1
 595       fAlphaMinCutMeson  = 0.0;
 596       fAlphaCutMeson     = 1.;
 597       break;
 598    case 4:      // 0-0.65
 599       fAlphaMinCutMeson  = 0.0;
 600       fAlphaCutMeson     = 0.65;
 601       break;
 602    case 5:      // 0-0.75
 603       fAlphaMinCutMeson  = 0.0;
 604       fAlphaCutMeson     = 0.75;
 605       break;
 606    case 6:      // 0-0.8
 607       fAlphaMinCutMeson  = 0.0;
 608       fAlphaCutMeson     = 0.8;
 609       break;
 610    case 7:      // 0.0-0.85
 611       fAlphaMinCutMeson  = 0.0;
 612       fAlphaCutMeson     = 0.85;
 613       break;
 614    case 8:      // 0.0-0.6
 615       fAlphaMinCutMeson  = 0.0;
 616       fAlphaCutMeson     = 0.6;
 617       break;
 618    case 9:      // 0.0-0.3
 619       fAlphaMinCutMeson  = 0.0;
 620       fAlphaCutMeson     = 0.3;
 621       break;
 622    default:
 623       cout<<"Warning: AlphaMesonCut not defined "<<alphaMesonCut<<endl;
 624       return kFALSE;
 625    }
 626    return kTRUE;
 627 }
 628
 629 ///________________________________________________________________________
 630 Bool_t AliConversionMesonCuts::SetRapidityMesonCut(Int_t RapidityMesonCut){
 631    // Set Cut
 632    switch(RapidityMesonCut){
 633    case 0:  //
 634       fRapidityCutMeson   = 0.9;
 635       break;
 636    case 1:  //
 637       fRapidityCutMeson   = 0.8;
 638       break;
 639    case 2:  //
 640       fRapidityCutMeson   = 0.7;
 641       break;
 642    case 3:  //
 643       fRapidityCutMeson   = 0.6;
 644       break;
 645    case 4:  //
 646       fRapidityCutMeson   = 0.5;
 647       break;
 648    case 5:  //
 649       fRapidityCutMeson   = 0.85;
 650       break;
 651    case 6:  //
 652       fRapidityCutMeson   = 0.75;
 653       break;
 654    case 7:  //
 655       fRapidityCutMeson   = 0.3;
 656       break;
 657    case 8:  //
 658       fRapidityCutMeson   = 0.35;
 659       break;
 660    case 9:  //
 661       fRapidityCutMeson   = 0.4;
 662       break;
 663    default:
 664       cout<<"Warning: RapidityMesonCut not defined "<<RapidityMesonCut<<endl;
 665       return kFALSE;
 666    }
 667    return kTRUE;
 668 }
 669
 670
 671 ///________________________________________________________________________
 672 Bool_t AliConversionMesonCuts::SetBackgroundScheme(Int_t BackgroundScheme){
 673    // Set Cut
 674    switch(BackgroundScheme){
 675    case 0: //Rotation
 676       fUseRotationMethodInBG=kTRUE;
 677       fdoBGProbability=kFALSE;
 678       break;
 679    case 1: // mixed event with V0 multiplicity
 680       fUseRotationMethodInBG=kFALSE;
 681       fUseTrackMultiplicityForBG=kFALSE;
 682       fdoBGProbability=kFALSE;
 683       break;
 684    case 2: // mixed event with track multiplicity
 685       fUseRotationMethodInBG=kFALSE;
 686       fUseTrackMultiplicityForBG=kTRUE;
 687       fdoBGProbability=kFALSE;
 688       break;
 689    case 3: //Rotation
 690       fUseRotationMethodInBG=kTRUE;
 691       fdoBGProbability=kTRUE;
 692       break;
 693    case 4: //No BG calculation
 694       cout << "no BG calculation should be done" << endl;
 695       fUseRotationMethodInBG=kFALSE;
 696       fdoBGProbability=kFALSE;
 697       fDoBG=kFALSE;
 698       break;
 699    case 5: //Rotation
 700       fUseRotationMethodInBG=kTRUE;
 701       fdoBGProbability=kFALSE;
 702       fBackgroundHandler = 1;
 703       break;
 704    case 6: // mixed event with V0 multiplicity
 705       fUseRotationMethodInBG=kFALSE;
 706       fUseTrackMultiplicityForBG=kFALSE;
 707       fdoBGProbability=kFALSE;
 708       fBackgroundHandler = 1;
 709       break;
 710    case 7: // mixed event with track multiplicity
 711       fUseRotationMethodInBG=kFALSE;
 712       fUseTrackMultiplicityForBG=kTRUE;
 713       fdoBGProbability=kFALSE;
 714       fBackgroundHandler = 1;
 715       break;
 716    case 8: //Rotation
 717       fUseRotationMethodInBG=kTRUE;
 718       fdoBGProbability=kTRUE;
 719       fBackgroundHandler = 1;
 720       break;
 721    default:
 722       cout<<"Warning: BackgroundScheme not defined "<<BackgroundScheme<<endl;
 723       return kFALSE;
 724    }
 725    return kTRUE;
 726 }
 727
 728
 729 ///________________________________________________________________________
 730 Bool_t AliConversionMesonCuts::SetNDegreesForRotationMethod(Int_t DegreesForRotationMethod){
 731    // Set Cut
 732    switch(DegreesForRotationMethod){
 733    case 0:
 734       fnDegreeRotationPMForBG = 5;
 735       break;
 736    case 1:
 737       fnDegreeRotationPMForBG = 10;
 738       break;
 739    case 2:
 740       fnDegreeRotationPMForBG = 15;
 741       break;
 742    case 3:
 743       fnDegreeRotationPMForBG = 20;
 744       break;
 745    default:
 746       cout<<"Warning: DegreesForRotationMethod not defined "<<DegreesForRotationMethod<<endl;
 747       return kFALSE;
 748    }
 749    fCuts[kDegreesForRotationMethod]=DegreesForRotationMethod;
 750    return kTRUE;
 751 }
 752
 753 ///________________________________________________________________________
 754 Bool_t AliConversionMesonCuts::SetNumberOfBGEvents(Int_t NumberOfBGEvents){
 755    // Set Cut
 756    switch(NumberOfBGEvents){
 757    case 0:
 758       fNumberOfBGEvents = 5;
 759       break;
 760    case 1:
 761       fNumberOfBGEvents = 10;
 762       break;
 763    case 2:
 764       fNumberOfBGEvents = 15;
 765       break;
 766    case 3:
 767       fNumberOfBGEvents = 20;
 768       break;
 769    case 4:
 770       fNumberOfBGEvents = 2;
 771       break;
 772    case 5:
 773       fNumberOfBGEvents = 50;
 774       break;
 775    case 6:
 776       fNumberOfBGEvents = 80;
 777       break;
 778    case 7:
 779       fNumberOfBGEvents = 100;
 780       break;
 781    default:
 782       cout<<"Warning: NumberOfBGEvents not defined "<<NumberOfBGEvents<<endl;
 783       return kFALSE;
 784    }
 785    return kTRUE;
 786 }
 787 ///________________________________________________________________________
 788 Bool_t AliConversionMesonCuts::SetSharedElectronCut(Int_t sharedElec) {
 789
 790    switch(sharedElec){
 791    case 0:
 792       fDoSharedElecCut = kFALSE;
 793       break;
 794    case 1:
 795       fDoSharedElecCut = kTRUE;
 796       break;
 797    default:
 798       cout<<"Warning: Shared Electron Cut not defined "<<sharedElec<<endl;
 799       return kFALSE;
 800    }
 801
 802    return kTRUE;
 803 }
 804
 805 ///________________________________________________________________________
 806 Bool_t AliConversionMesonCuts::SetToCloseV0sCut(Int_t toClose) {
 807
 808    switch(toClose){
 809    case 0:
 810       fDoToCloseV0sCut = kFALSE;
 811       fminV0Dist = 250;
 812       break;
 813    case 1:
 814       fDoToCloseV0sCut = kTRUE;
 815       fminV0Dist = 1;
 816       break;
 817    case 2:
 818       fDoToCloseV0sCut = kTRUE;
 819       fminV0Dist = 2;
 820       break;
 821    case 3:
 822       fDoToCloseV0sCut = kTRUE;
 823       fminV0Dist = 3;
 824       break;
 825    default:
 826       cout<<"Warning: Shared Electron Cut not defined "<<toClose<<endl;
 827       return kFALSE;
 828    }
 829    return kTRUE;
 830 }
 831
 832 ///________________________________________________________________________
 833 Bool_t AliConversionMesonCuts::SetMCPSmearing(Int_t useMCPSmearing)
 834 {// Set Cut
 835    switch(useMCPSmearing){
 836    case 0:
 837       fUseMCPSmearing=0;
 838       fPBremSmearing=1.;
 839       fPSigSmearing=0.;
 840       fPSigSmearingCte=0.;
 841       break;
 842    case 1:
 843       fUseMCPSmearing=1;
 844       fPBremSmearing=1.0e-14;
 845       fPSigSmearing=0.;
 846       fPSigSmearingCte=0.;
 847       break;
 848    case 2:
 849       fUseMCPSmearing=1;
 850       fPBremSmearing=1.0e-15;
 851       fPSigSmearing=0.0;
 852       fPSigSmearingCte=0.;
 853       break;
 854    case 3:
 855       fUseMCPSmearing=1;
 856       fPBremSmearing=1.;
 857       fPSigSmearing=0.003;
 858       fPSigSmearingCte=0.002;
 859       break;
 860    case 4:
 861       fUseMCPSmearing=1;
 862       fPBremSmearing=1.;
 863       fPSigSmearing=0.003;
 864       fPSigSmearingCte=0.007;
 865       break;
 866    case 5:
 867       fUseMCPSmearing=1;
 868       fPBremSmearing=1.;
 869       fPSigSmearing=0.003;
 870       fPSigSmearingCte=0.016;
 871       break;
 872    case 6:
 873       fUseMCPSmearing=1;
 874       fPBremSmearing=1.;
 875       fPSigSmearing=0.007;
 876       fPSigSmearingCte=0.016;
 877       break;
 878    case 7:
 879       fUseMCPSmearing=1;
 880       fPBremSmearing=1.0e-16;
 881       fPSigSmearing=0.0;
 882       fPSigSmearingCte=0.;
 883       break;
 884    case 8:
 885       fUseMCPSmearing=1;
 886       fPBremSmearing=1.;
 887       fPSigSmearing=0.007;
 888       fPSigSmearingCte=0.014;
 889       break;
 890    case 9:
 891       fUseMCPSmearing=1;
 892       fPBremSmearing=1.;
 893       fPSigSmearing=0.007;
 894       fPSigSmearingCte=0.011;
 895       break;
 896
 897    default:
 898       AliError("Warning: UseMCPSmearing not defined");
 899       return kFALSE;
 900    }
 901    return kTRUE;
 902 }
 903 ///________________________________________________________________________
 904 TString AliConversionMesonCuts::GetCutNumber(){
 905    // returns TString with current cut number
 906    TString a(kNCuts);
 907    for(Int_t ii=0;ii<kNCuts;ii++){
 908       a.Append(Form("%d",fCuts[ii]));
 909    }
 910    return a;
 911 }
 912
 913 ///________________________________________________________________________
 914 void AliConversionMesonCuts::FillElectonLabelArray(AliAODConversionPhoton* photon, Int_t nV0){
 915
 916    Int_t posLabel = photon->GetTrackLabelPositive();
 917    Int_t negLabel = photon->GetTrackLabelNegative();
 918
 919    fElectronLabelArray[nV0*2] = posLabel;
 920    fElectronLabelArray[(nV0*2)+1] = negLabel;
 921 }
 922
 923 ///________________________________________________________________________
 924 Bool_t AliConversionMesonCuts::RejectSharedElectronV0s(AliAODConversionPhoton* photon, Int_t nV0, Int_t nV0s){
 925
 926    Int_t posLabel = photon->GetTrackLabelPositive();
 927    Int_t negLabel = photon->GetTrackLabelNegative();
 928
 929    for(Int_t i = 0; i<nV0s*2;i++){
 930       if(i==nV0*2)     continue;
 931       if(i==(nV0*2)+1) continue;
 932       if(fElectronLabelArray[i] == posLabel){
 933          return kFALSE;}
 934       if(fElectronLabelArray[i] == negLabel){
 935          return kFALSE;}
 936    }
 937
 938    return kTRUE;
 939 }
 940 ///________________________________________________________________________
 941 Bool_t AliConversionMesonCuts::RejectToCloseV0s(AliAODConversionPhoton* photon, TList *photons, Int_t nV0){
 942    Double_t posX = photon->GetConversionX();
 943    Double_t posY = photon->GetConversionY();
 944    Double_t posZ = photon->GetConversionZ();
 945
 946    for(Int_t i = 0;i<photons->GetEntries();i++){
 947       if(nV0 == i) continue;
 948       AliAODConversionPhoton *photonComp = (AliAODConversionPhoton*) photons->At(i);
 949       Double_t posCompX = photonComp->GetConversionX();
 950       Double_t posCompY = photonComp->GetConversionY();
 951       Double_t posCompZ = photonComp->GetConversionZ();
 952
 953       Double_t dist = pow((posX - posCompX),2)+pow((posY - posCompY),2)+pow((posZ - posCompZ),2);
 954
 955       if(dist < fminV0Dist*fminV0Dist){
 956          if(photon->GetChi2perNDF() < photonComp->GetChi2perNDF()) return kTRUE;
 957          else {
 958             return kFALSE;}
 959       }
 960
 961    }
 962    return kTRUE;
 963 }
 964
 965 ///________________________________________________________________________
 966 void AliConversionMesonCuts::SmearParticle(AliAODConversionPhoton* photon)
 967 {
 968
 969    if (photon==NULL) return;
 970    Double_t facPBrem = 1.;
 971    Double_t facPSig = 0.;
 972
 973    Double_t phi=0.;
 974    Double_t theta=0.;
 975    Double_t P=0.;
 976
 977
 978    P=photon->P();
 979    phi=photon->Phi();
 980    if( photon->P()!=0){
 981       theta=acos( photon->Pz()/ photon->P());
 982    }
 983
 984    if( fPSigSmearing != 0. || fPSigSmearingCte!=0. ){
 985       facPSig = TMath::Sqrt(fPSigSmearingCte*fPSigSmearingCte+fPSigSmearing*fPSigSmearing*P*P)*fRandom.Gaus(0.,1.);
 986    }
 987
 988    if( fPBremSmearing != 1.){
 989       if(fBrem!=NULL){
 990          facPBrem = fBrem->GetRandom();
 991       }
 992    }
 993
 994    photon->SetPx(facPBrem* (1+facPSig)* P*sin(theta)*cos(phi)) ;
 995    photon->SetPy(facPBrem* (1+facPSig)* P*sin(theta)*sin(phi)) ;
 996    photon->SetPz(facPBrem* (1+facPSig)* P*cos(theta)) ;
 997    photon->SetE(photon->P());
 998 }