PWGHF/correlationHF/AliAnalysisTaskSEDplusCorrelations.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-2008, 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
  16 /* $Id$ */
  17
  18 //*************************************************************************
  19 // Class AliAnalysisTaskSEDplusCorrelations
  20 // AliAnalysisTaskSE for the D+ candidates Invariant Mass Histogram and
  21 //comparison of heavy-flavour decay candidates
  22 // to MC truth (kinematics stored in the AOD)
  23 // Authors: Sadhana Dash (correlation)
  24 /////////////////////////////////////////////////////////////
  25
  26 #include <TClonesArray.h>
  27 #include <TNtuple.h>
  28 #include <TCanvas.h>
  29 #include <TList.h>
  30 #include <TString.h>
  31 #include <TDatabasePDG.h>
  32
  33 #include <AliAnalysisDataSlot.h>
  34 #include <AliAnalysisDataContainer.h>
  35 #include "AliAnalysisManager.h"
  36 #include "AliRDHFCutsDplustoKpipi.h"
  37 #include "AliAODHandler.h"
  38 #include "AliAODEvent.h"
  39 #include "AliAODPidHF.h"
  40 #include "AliAODVertex.h"
  41 #include "AliAODTrack.h"
  42 #include "AliAODRecoDecayHF3Prong.h"
  43 #include "AliAnalysisVertexingHF.h"
  44 #include "AliAnalysisTaskSE.h"
  45 #include "AliAnalysisTaskSEDplusCorrelations.h"
  46 #include "AliNormalizationCounter.h"
  47 #include "AliVParticle.h"
  48 #include "AliHFAssociatedTrackCuts.h"
  49 #include "AliReducedParticle.h"
  50 #include "AliHFCorrelator.h"
  51
  52
  53 using std::cout;
  54 using std::endl;
  55
  56 ClassImp(AliAnalysisTaskSEDplusCorrelations)
  57
  58
  59 //________________________________________________________________________
  60 AliAnalysisTaskSEDplusCorrelations::AliAnalysisTaskSEDplusCorrelations():
  61 AliAnalysisTaskSE(),
  62   fOutput(0),
  63   fCorrelator(0),
  64   fSelect(0),
  65   fDisplacement(0),
  66   fHistNEvents(0),
  67   fMCSources(0),
  68   fK0Origin(0),
  69   fKaonOrigin(0),
  70   fInvMassK0S(0),
  71   fEventMix(0),
  72   fPtVsMass(0),
  73   fPtVsMassTC(0),
  74   fYVsPt(0),
  75   fYVsPtTC(0),
  76   fYVsPtSig(0),
  77   fYVsPtSigTC(0),
  78   fUpmasslimit(1.965),
  79   fLowmasslimit(1.765),
  80   fNPtBins(0),
  81   fBinWidth(0.002),
  82   fListCuts(0),
  83   fListCutsAsso(0),
  84   fRDCutsAnalysis(0),
  85   fCuts(0),
  86   fCounter(0),
  87   fReadMC(kFALSE),
  88   fUseBit(kTRUE),
  89   fMixing(kFALSE),
  90   fSystem(kFALSE)
  91 {
  92   // Default constructor
  93
  94   for(Int_t i=0;i<3*kMaxPtBins;i++){
  95     if(fMassVsdPhiHistHad[i])fMassVsdPhiHistHad[i]=0;
  96     if(fMassVsdEtaHistHad[i])fMassVsdEtaHistHad[i]=0;
  97     if(fMassVsdPhiHistKaon[i])fMassVsdPhiHistKaon[i]=0;
  98     if(fMassVsdEtaHistKaon[i])fMassVsdEtaHistKaon[i]=0;
  99     if(fMassVsdPhiHistKshort[i])fMassVsdPhiHistKshort[i]=0;
 100     if(fMassVsdEtaHistKshort[i])fMassVsdEtaHistKshort[i]=0;
 101     if(fMassVsdPhiHistLeadHad[i])fMassVsdPhiHistLeadHad[i]=0;
 102     if(fMassVsdEtaHistLeadHad[i])fMassVsdEtaHistLeadHad[i]=0;
 103     if(fMassHistK0S[i])fMassHistK0S[i]=0;
 104     if(fLeadPt[i])fLeadPt[i]=0;
 105     if(fMassHist[i])fMassHist[i]=0;
 106     if(fMassHistTCPlus[i])fMassHistTCPlus[i]=0;
 107     if(fMassHistTCMinus[i])fMassHistTCMinus[i]=0;
 108   }
 109
 110   for(Int_t i=0;i<kMaxPtBins+1;i++){
 111     fArrayBinLimits[i]=0;
 112   }
 113
 114 }
 115
 116 //________________________________________________________________________
 117 AliAnalysisTaskSEDplusCorrelations::AliAnalysisTaskSEDplusCorrelations(const char *name,AliRDHFCutsDplustoKpipi *dpluscutsana, AliHFAssociatedTrackCuts *asstrkcuts):
 118   AliAnalysisTaskSE(name),
 119   fOutput(0),
 120   fCorrelator(0),
 121   fSelect(0),
 122   fDisplacement(0),
 123   fHistNEvents(0),
 124   fMCSources(0),
 125   fK0Origin(0),
 126   fKaonOrigin(0),
 127   fInvMassK0S(0),
 128   fEventMix(0),
 129   fPtVsMass(0),
 130   fPtVsMassTC(0),
 131   fYVsPt(0),
 132   fYVsPtTC(0),
 133   fYVsPtSig(0),
 134   fYVsPtSigTC(0),
 135   fUpmasslimit(1.965),
 136   fLowmasslimit(1.765),
 137   fNPtBins(0),
 138   fBinWidth(0.002),
 139   fListCuts(0),
 140   fListCutsAsso(0),
 141   fRDCutsAnalysis(dpluscutsana),
 142   fCuts(asstrkcuts),
 143   fCounter(0),
 144   fReadMC(kFALSE),
 145   fUseBit(kTRUE),
 146   fMixing(kFALSE),
 147   fSystem(kFALSE)
 148 {
 149   //
 150   // Standrd constructor
 151   //
 152   fNPtBins=fRDCutsAnalysis->GetNPtBins();
 153
 154
 155   for(Int_t i=0;i<3*kMaxPtBins;i++){
 156     if(fMassVsdPhiHistHad[i])fMassVsdPhiHistHad[i]=0;
 157     if(fMassVsdEtaHistHad[i])fMassVsdEtaHistHad[i]=0;
 158     if(fMassVsdPhiHistKaon[i])fMassVsdPhiHistKaon[i]=0;
 159     if(fMassVsdEtaHistKaon[i])fMassVsdEtaHistKaon[i]=0;
 160     if(fMassVsdPhiHistKshort[i])fMassVsdPhiHistKshort[i]=0;
 161     if(fMassVsdEtaHistKshort[i])fMassVsdEtaHistKshort[i]=0;
 162     if(fMassVsdPhiHistLeadHad[i])fMassVsdPhiHistLeadHad[i]=0;
 163     if(fMassVsdEtaHistLeadHad[i])fMassVsdEtaHistLeadHad[i]=0;
 164     if(fMassHistK0S[i])fMassHistK0S[i]=0;
 165     if(fLeadPt[i])fLeadPt[i]=0;
 166     if(fMassHist[i])fMassHist[i]=0;
 167     if(fMassHistTC[i])fMassHistTC[i]=0;
 168     if(fMassHistTCPlus[i])fMassHistTCPlus[i]=0;
 169     if(fMassHistTCMinus[i])fMassHistTCMinus[i]=0;
 170
 171
 172   }
 173
 174   for(Int_t i=0;i<kMaxPtBins+1;i++){
 175     fArrayBinLimits[i]=0;
 176   }
 177
 178
 179   // Default constructor
 180   // Output slot #1 writes into a TList container
 181   DefineOutput(1,TList::Class());  //My private output
 182   // Output slot #2 writes cut to private output
 183   //  DefineOutput(2,AliRDHFCutsDplusCorrelationstoKpipi::Class());
 184   DefineOutput(2,TList::Class());
 185   // Output slot #3 writes cut to private output
 186   DefineOutput(3,AliNormalizationCounter::Class());
 187   DefineOutput(4,AliHFAssociatedTrackCuts::Class());
 188
 189
 190 }
 191
 192 //________________________________________________________________________
 193 AliAnalysisTaskSEDplusCorrelations::~AliAnalysisTaskSEDplusCorrelations()
 194 {
 195   //
 196   // Destructor
 197   //
 198   delete fOutput;
 199   delete fListCuts;
 200   delete fRDCutsAnalysis;
 201   delete fCuts;
 202   delete fCounter;
 203   delete fCorrelator;
 204
 205
 206 }
 207 //_________________________________________________________________
 208 void  AliAnalysisTaskSEDplusCorrelations::SetMassLimits(Float_t range){
 209   // set invariant mass limits
 210   Float_t bw=GetBinWidth();
 211   fUpmasslimit = 1.865+range;
 212   fLowmasslimit = 1.865-range;
 213   SetBinWidth(bw);
 214 }
 215 //_________________________________________________________________
 216 void  AliAnalysisTaskSEDplusCorrelations::SetMassLimits(Float_t lowlimit, Float_t uplimit){
 217   // set invariant mass limits
 218   if(uplimit>lowlimit)
 219     {
 220       Float_t bw=GetBinWidth();
 221       fUpmasslimit = lowlimit;
 222       fLowmasslimit = uplimit;
 223       SetBinWidth(bw);
 224     }
 225 }
 226 //________________________________________________________________
 227 void AliAnalysisTaskSEDplusCorrelations::SetBinWidth(Float_t w){
 228   // Define width of mass bins
 229   Float_t width=w;
 230   Int_t nbins=(Int_t)((fUpmasslimit-fLowmasslimit)/width+0.5);
 231   Int_t missingbins=4-nbins%4;
 232   nbins=nbins+missingbins;
 233   width=(fUpmasslimit-fLowmasslimit)/nbins;
 234   if(missingbins!=0){
 235     printf("AliAnalysisTaskSEDplusCorrelations::SetBinWidth: W-bin width of %f will produce histograms not rebinnable by 4. New width set to %f\n",w,width);
 236   }
 237   else{
 238     if(fDebug>1) printf("AliAnalysisTaskSEDplusCorrelations::SetBinWidth: width set to %f\n",width);
 239   }
 240   fBinWidth=width;
 241 }
 242 //_________________________________________________________________
 243 Int_t AliAnalysisTaskSEDplusCorrelations::GetNBinsHistos(){
 244   // Compute number of mass bins
 245   return (Int_t)((fUpmasslimit-fLowmasslimit)/fBinWidth+0.5);
 246 }
 247
 248
 249 //__________________________________________
 250 void AliAnalysisTaskSEDplusCorrelations::Init(){
 251   //
 252   // Initialization
 253   //
 254   if(fDebug > 1) printf("AnalysisTaskSEDplusCorrelations::Init() \n");
 255
 256   //PostData(2,fRDCutsloose);//we should then put those cuts in a tlist if we have more than 1
 257   fListCuts=new TList();
 258   // fListCutsAsso=new TList();
 259
 260   AliRDHFCutsDplustoKpipi *analysis = new AliRDHFCutsDplustoKpipi(*fRDCutsAnalysis);
 261   analysis->SetName("AnalysisCuts");
 262
 263   // AliHFAssociatedTrackCuts *trkcuts = new AliHFAssociatedTrackCuts(*fCuts);
 264   //trkcuts->SetName("Assotrkcuts");
 265
 266   fListCuts->Add(analysis);
 267   //fListCuts->Add(trkcuts);
 268
 269
 270
 271   PostData(2,fListCuts);
 272   PostData(4,fCuts);
 273
 274   return;
 275 }
 276
 277 //________________________________________________________________________
 278 void AliAnalysisTaskSEDplusCorrelations::UserCreateOutputObjects()
 279 {
 280   // Create the output container
 281   //
 282   if(fDebug > 1) printf("AnalysisTaskSEDplusCorrelations::UserCreateOutputObjects() \n");
 283   // correlator creation and definition
 284
 285   Double_t Pi = TMath::Pi();
 286   fCorrelator = new AliHFCorrelator("Correlator",fCuts,fSystem); // fCuts is the hadron cut object, fSystem to switch between pp or PbPb
 287   fCorrelator->SetDeltaPhiInterval((-0.5-1./32)*Pi,(1.5-1./32)*Pi); // set correct phi interval
 288   //fCorrelator->SetDeltaPhiInterval(-1.57,4.71);
 289   fCorrelator->SetEventMixing(fMixing); //set kFALSE/kTRUE for mixing Off/On
 290   fCorrelator->SetAssociatedParticleType(fSelect); // set 1/2/3 for hadron/kaons/kzeros
 291   fCorrelator->SetApplyDisplacementCut(fDisplacement); //set kFALSE/kTRUE for using the displacement cut
 292   fCorrelator->SetUseMC(fReadMC);
 293
 294   Bool_t pooldef = fCorrelator->DefineEventPool();
 295
 296   if(!pooldef) AliInfo("Warning:: Event pool not defined properly");
 297
 298
 299   // Several histograms are more conveniently managed in a TList
 300   fOutput = new TList();
 301   fOutput->SetOwner();
 302   fOutput->SetName("OutputHistos");
 303
 304   TString hisname;
 305   Int_t index=0;
 306   Int_t nbins=GetNBinsHistos();
 307
 308   for(Int_t i=0;i<fNPtBins;i++){
 309
 310     index=GetHistoIndex(i);
 311
 312
 313     hisname.Form("hMassVsdPhiHad%d",i);
 314     fMassVsdPhiHistHad[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,25,(-1)*TMath::Pi()/2,3.*TMath::Pi()/2.);
 315     fMassVsdPhiHistHad[index]->Sumw2();
 316
 317     hisname.Form("hMassVsdEtaHad%d",i);
 318     fMassVsdEtaHistHad[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,25,(-1)*TMath::Pi()/2,3.*TMath::Pi()/2.);
 319     fMassVsdEtaHistHad[index]->Sumw2();
 320
 321     hisname.Form("hMassVsdPhiKaon%d",i);
 322     fMassVsdPhiHistKaon[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,25,(-1)*TMath::Pi()/2,3.*TMath::Pi()/2.);
 323     fMassVsdPhiHistKaon[index]->Sumw2();
 324
 325     hisname.Form("hMassVsdEtaKaon%d",i);
 326     fMassVsdEtaHistKaon[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,25,(-1)*TMath::Pi()/2,3.*TMath::Pi()/2.);
 327     fMassVsdEtaHistKaon[index]->Sumw2();
 328
 329     hisname.Form("hMassVsdPhiK0%d",i);
 330     fMassVsdPhiHistKshort[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,25,(-1)*TMath::Pi()/2,3.*TMath::Pi()/2.);
 331     fMassVsdPhiHistKshort[index]->Sumw2();
 332
 333     hisname.Form("hMassVsdEtaK0%d",i);
 334     fMassVsdEtaHistKshort[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,25,(-1)*TMath::Pi()/2,3.*TMath::Pi()/2.);
 335     fMassVsdEtaHistKshort[index]->Sumw2();
 336
 337     hisname.Form("hMassVsdPhiLeadHad%d",i);
 338     fMassVsdPhiHistLeadHad[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,25,(-1)*TMath::Pi()/2,3.*TMath::Pi()/2.);
 339     fMassVsdPhiHistLeadHad[index]->Sumw2();
 340
 341     hisname.Form("hMassVsdEtaLeadHad%d",i);
 342     fMassVsdEtaHistLeadHad[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,25,(-1)*TMath::Pi()/2,3.*TMath::Pi()/2.);
 343     fMassVsdEtaHistLeadHad[index]->Sumw2();
 344
 345     hisname.Form("hMassPtK0S%d",i);
 346     fMassHistK0S[index]=new TH1F(hisname.Data(),hisname.Data(),500,0.3,0.8);
 347     fMassHistK0S[index]->Sumw2();
 348
 349     hisname.Form("hLeadPt%d",i);
 350     fLeadPt[index]=new TH1F(hisname.Data(),hisname.Data(),500,0.0,50);
 351     fLeadPt[index]->Sumw2();
 352
 353
 354     hisname.Form("hMassPt%d",i);
 355     fMassHist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
 356     fMassHist[index]->Sumw2();
 357
 358
 359     hisname.Form("hMassPt%dTC",i);
 360     fMassHistTC[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
 361     fMassHistTC[index]->Sumw2();
 362
 363     hisname.Form("hMassPt%dTCPlus",i);
 364     fMassHistTCPlus[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
 365     fMassHistTCPlus[index]->Sumw2();
 366
 367     hisname.Form("hMassPt%dTCMinus",i);
 368     fMassHistTCMinus[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
 369     fMassHistTCMinus[index]->Sumw2();
 370
 371
 372
 373     index=GetSignalHistoIndex(i);
 374
 375     hisname.Form("hMassVsdPhiHadSig%d",i);
 376     fMassVsdPhiHistHad[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,25,(-1)*TMath::Pi()/2,3.*TMath::Pi()/2.);
 377     fMassVsdPhiHistHad[index]->Sumw2();
 378
 379     hisname.Form("hMassVsdEtaHadSig%d",i);
 380     fMassVsdEtaHistHad[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,25,(-1)*TMath::Pi()/2,3.*TMath::Pi()/2.);
 381     fMassVsdEtaHistHad[index]->Sumw2();
 382
 383     hisname.Form("hMassVsdPhiKaonSig%d",i);
 384     fMassVsdPhiHistKaon[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,25,(-1)*TMath::Pi()/2,3.*TMath::Pi()/2.);
 385     fMassVsdPhiHistKaon[index]->Sumw2();
 386
 387     hisname.Form("hMassVsdEtaKaonSig%d",i);
 388     fMassVsdEtaHistKaon[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,25,(-1)*TMath::Pi()/2,3.*TMath::Pi()/2.);
 389     fMassVsdEtaHistKaon[index]->Sumw2();
 390
 391     hisname.Form("hMassVsdPhiK0Sig%d",i);
 392     fMassVsdPhiHistKshort[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,25,(-1)*TMath::Pi()/2,3.*TMath::Pi()/2.);
 393     fMassVsdPhiHistKshort[index]->Sumw2();
 394
 395     hisname.Form("hMassVsdEtaK0Sig%d",i);
 396     fMassVsdEtaHistKshort[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,25,(-1)*TMath::Pi()/2,3.*TMath::Pi()/2.);
 397     fMassVsdEtaHistKshort[index]->Sumw2();
 398
 399     hisname.Form("hMassVsdPhiLeadHadSig%d",i);
 400     fMassVsdPhiHistLeadHad[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,25,(-1)*TMath::Pi()/2,3.*TMath::Pi()/2.);
 401     fMassVsdPhiHistLeadHad[index]->Sumw2();
 402
 403     hisname.Form("hMassVsdEtaLeadHadSig%d",i);
 404     fMassVsdEtaHistLeadHad[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,25,(-1)*TMath::Pi()/2,3.*TMath::Pi()/2.);
 405     fMassVsdEtaHistLeadHad[index]->Sumw2();
 406
 407
 408     hisname.Form("hSigPtK0S%d",i);
 409     fMassHistK0S[index]=new TH1F(hisname.Data(),hisname.Data(),500,0.3,0.8);
 410     fMassHistK0S[index]->Sumw2();
 411
 412     hisname.Form("hSigLeadPt%d",i);
 413     fLeadPt[index]=new TH1F(hisname.Data(),hisname.Data(),500,0.0,50);
 414     fLeadPt[index]->Sumw2();
 415
 416     hisname.Form("hSigPt%d",i);
 417     fMassHist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
 418     fMassHist[index]->Sumw2();
 419
 420     hisname.Form("hSigPt%dTC",i);
 421     fMassHistTC[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
 422     fMassHistTC[index]->Sumw2();
 423     hisname.Form("hSigPt%dTCPlus",i);
 424     fMassHistTCPlus[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
 425     fMassHistTCPlus[index]->Sumw2();
 426     hisname.Form("hSigPt%dTCMinus",i);
 427     fMassHistTCMinus[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
 428     fMassHistTCMinus[index]->Sumw2();
 429
 430
 431
 432     index=GetBackgroundHistoIndex(i);
 433
 434     hisname.Form("hMassVsdPhiBkgHad%d",i);
 435     fMassVsdPhiHistHad[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,25,(-1)*TMath::Pi()/2,3.*TMath::Pi()/2.);
 436     fMassVsdPhiHistHad[index]->Sumw2();
 437
 438     hisname.Form("hMassVsdEtaBkgHad%d",i);
 439     fMassVsdEtaHistHad[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,25,(-1)*TMath::Pi()/2,3.*TMath::Pi()/2.);
 440     fMassVsdEtaHistHad[index]->Sumw2();
 441
 442     hisname.Form("hMassVsdPhiBkgKaon%d",i);
 443     fMassVsdPhiHistKaon[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,25,(-1)*TMath::Pi()/2,3.*TMath::Pi()/2.);
 444     fMassVsdPhiHistKaon[index]->Sumw2();
 445
 446     hisname.Form("hMassVsdEtaBkgKaon%d",i);
 447     fMassVsdEtaHistKaon[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,25,(-1)*TMath::Pi()/2,3.*TMath::Pi()/2.);
 448     fMassVsdEtaHistKaon[index]->Sumw2();
 449
 450     hisname.Form("hMassVsdPhiBkgKshort%d",i);
 451     fMassVsdPhiHistKshort[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,25,(-1)*TMath::Pi()/2,3.*TMath::Pi()/2.);
 452     fMassVsdPhiHistKshort[index]->Sumw2();
 453
 454
 455     hisname.Form("hMassVsdPhiBkgKshort%d",i);
 456     fMassVsdEtaHistKshort[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,25,(-1)*TMath::Pi()/2,3.*TMath::Pi()/2.);
 457     fMassVsdEtaHistKshort[index]->Sumw2();
 458
 459     hisname.Form("hMassVsdPhiBkgLeadHad%d",i);
 460     fMassVsdPhiHistLeadHad[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,25,(-1)*TMath::Pi()/2,3.*TMath::Pi()/2.);
 461     fMassVsdPhiHistLeadHad[index]->Sumw2();
 462
 463     hisname.Form("hMassVsdPhiBkgKshort%d",i);
 464     fMassVsdEtaHistLeadHad[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,25,(-1)*TMath::Pi()/2,3.*TMath::Pi()/2.);
 465     fMassVsdEtaHistLeadHad[index]->Sumw2();
 466
 467     hisname.Form("hBkgPtK0S%d",i);
 468     fMassHistK0S[index]=new TH1F(hisname.Data(),hisname.Data(),500,0.3,0.8);
 469     fMassHistK0S[index]->Sumw2();
 470
 471     hisname.Form("hLeadBkgPt%d",i);
 472     fLeadPt[index]=new TH1F(hisname.Data(),hisname.Data(),500,0.0,50);
 473     fLeadPt[index]->Sumw2();
 474
 475     hisname.Form("hBkgPt%dTC",i);
 476     fMassHistTC[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
 477     fMassHistTC[index]->Sumw2();
 478
 479     hisname.Form("hBkgPt%dTCPlus",i);
 480     fMassHistTCPlus[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
 481     fMassHistTCPlus[index]->Sumw2();
 482
 483     hisname.Form("hBkgPt%dTCMinus",i);
 484     fMassHistTCMinus[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
 485     fMassHistTCMinus[index]->Sumw2();
 486   }
 487
 488
 489   for(Int_t i=0; i<3*fNPtBins; i++){
 490     fOutput->Add(fMassVsdPhiHistHad[i]);
 491     fOutput->Add(fMassVsdEtaHistHad[i]);
 492     fOutput->Add(fMassVsdPhiHistKaon[i]);
 493     fOutput->Add(fMassVsdEtaHistKaon[i]);
 494     fOutput->Add(fMassVsdPhiHistKshort[i]);
 495     fOutput->Add(fMassVsdEtaHistKshort[i]);
 496     fOutput->Add(fMassVsdPhiHistLeadHad[i]);
 497     fOutput->Add(fMassVsdEtaHistLeadHad[i]);
 498     // fOutput->Add(fMassVsdPhiHist[i]);
 499     fOutput->Add(fMassHistK0S[i]);
 500     fOutput->Add(fLeadPt[i]);
 501     fOutput->Add(fMassHist[i]);
 502     fOutput->Add(fMassHistTC[i]);
 503     fOutput->Add(fMassHistTCPlus[i]);
 504     fOutput->Add(fMassHistTCMinus[i]);
 505
 506
 507
 508   }
 509   fInvMassK0S = new TH2F("K0S","K0S", 500,0.3,0.8,500,0,50);
 510   fInvMassK0S->GetXaxis()->SetTitle("Invariant Mass (#pi #pi) (GeV/c^{2})");
 511   fInvMassK0S->GetYaxis()->SetTitle("K0S pt (GeV/c)");
 512   fOutput->Add(fInvMassK0S);
 513
 514   fMCSources = new TH1F("MCSources","Origin of associated particles in MC", 10, -0.5, 9.5);
 515   fMCSources->GetXaxis()->SetBinLabel(1,"All ");
 516   fMCSources->GetXaxis()->SetBinLabel(2," from Heavy flavour");
 517   fMCSources->GetXaxis()->SetBinLabel(3," from c->D");
 518   fMCSources->GetXaxis()->SetBinLabel(4," from b->D");
 519   fMCSources->GetXaxis()->SetBinLabel(5," from b->B");
 520   if(fReadMC) fOutput->Add(fMCSources);
 521
 522   fK0Origin = new TH1F("K0Origin","Origin of K0 ", 10, -0.5, 5.5);
 523   fK0Origin->GetXaxis()->SetBinLabel(1,"All K0s");
 524   fK0Origin->GetXaxis()->SetBinLabel(2,"K0s from Heavy flavour");
 525   fK0Origin->GetXaxis()->SetBinLabel(3,"K0s from Charm");
 526   fK0Origin->GetXaxis()->SetBinLabel(4,"K0s from Beauty");
 527   if(fReadMC) fOutput->Add(fK0Origin);
 528
 529   fKaonOrigin = new TH1F("K0Origin","Origin of Kaon ", 10, -0.5, 5.5);
 530   fKaonOrigin->GetXaxis()->SetBinLabel(1,"All Kaons");
 531   fKaonOrigin->GetXaxis()->SetBinLabel(2,"Kaons from Heavy flavour");
 532   fKaonOrigin->GetXaxis()->SetBinLabel(3,"Kaons from Charm");
 533   fKaonOrigin->GetXaxis()->SetBinLabel(4,"Kaons from Beauty");
 534   if(fReadMC) fOutput->Add(fKaonOrigin);
 535
 536
 537
 538   fEventMix = new TH2F("EventMixingCheck","EventMixingCheck",5,-0.5,4.5,7,-0.5,6.5);
 539   if(fMixing)fOutput->Add(fEventMix);
 540
 541
 542   fHistNEvents = new TH1F("fHistNEvents", "number of events ",11,-0.5,10.5);
 543   fHistNEvents->GetXaxis()->SetBinLabel(1,"nEventsAnal");
 544   fHistNEvents->GetXaxis()->SetBinLabel(2,"nEvents accepted");
 545   fHistNEvents->GetXaxis()->SetBinLabel(3,"Rejected due to trigger");
 546   fHistNEvents->GetXaxis()->SetBinLabel(4,"Rejected pileup events");
 547   fHistNEvents->GetXaxis()->SetBinLabel(5,"Rejected due to centrality");
 548   fHistNEvents->GetXaxis()->SetBinLabel(6,"Rejected due to vtxz");
 549   fHistNEvents->GetXaxis()->SetBinLabel(7,"Rejected due to Physics Sel");
 550   fHistNEvents->GetXaxis()->SetBinLabel(8,"Total no. of candidate");
 551   fHistNEvents->GetXaxis()->SetBinLabel(9,"no. of cand wo bitmask");
 552   fHistNEvents->GetXaxis()->SetBinLabel(10,"D+ after loose cuts");
 553   fHistNEvents->GetXaxis()->SetBinLabel(11,"D+ after tight cuts");
 554
 555   fHistNEvents->GetXaxis()->SetNdivisions(1,kFALSE);
 556   fHistNEvents->Sumw2();
 557   fHistNEvents->SetMinimum(0);
 558   fOutput->Add(fHistNEvents);
 559
 560   fPtVsMass=new TH2F("hPtVsMass","PtVsMass (prod. cuts)",nbins,fLowmasslimit,fUpmasslimit,200,0.,20.);
 561   fPtVsMassTC=new TH2F("hPtVsMassTC","PtVsMass (analysis cuts)",nbins,fLowmasslimit,fUpmasslimit,200,0.,20.);
 562   fYVsPt=new TH2F("hYVsPt","YvsPt (prod. cuts)",40,0.,20.,80,-2.,2.);
 563   fYVsPtTC=new TH2F("hYVsPtTC","YvsPt (analysis cuts)",40,0.,20.,80,-2.,2.);
 564   fYVsPtSig=new TH2F("hYVsPtSig","YvsPt (MC, only sig., prod. cuts)",40,0.,20.,80,-2.,2.);
 565   fYVsPtSigTC=new TH2F("hYVsPtSigTC","YvsPt (MC, only Sig, analysis cuts)",40,0.,20.,80,-2.,2.);
 566
 567   fOutput->Add(fPtVsMass);
 568   fOutput->Add(fPtVsMassTC);
 569   fOutput->Add(fYVsPt);
 570   fOutput->Add(fYVsPtTC);
 571   fOutput->Add(fYVsPtSig);
 572   fOutput->Add(fYVsPtSigTC);
 573
 574
 575   // Counter for Normalization
 576   TString normName="NormalizationCounter";
 577   AliAnalysisDataContainer *cont = GetOutputSlot(3)->GetContainer();
 578   if(cont)normName=(TString)cont->GetName();
 579   fCounter = new AliNormalizationCounter(normName.Data());
 580   fCounter->Init();
 581
 582
 583
 584   PostData(1,fOutput);
 585   PostData(3,fCounter);
 586   return;
 587 }
 588
 589 //________________________________________________________________________
 590 void AliAnalysisTaskSEDplusCorrelations::UserExec(Option_t */*option*/)
 591 {
 592   // Do the analysis
 593   AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
 594   TClonesArray *array3Prong = 0;
 595
 596   if(!fMixing){
 597     if(fSelect==1) cout << "TASK::Correlation with hadrons on SE "<< endl;
 598     if(fSelect==2) cout << "TASK::Correlation with kaons on SE "<< endl;
 599     if(fSelect==3) cout << "TASK::Correlation with kzeros on SE "<< endl;
 600   }
 601   if(fMixing){
 602     if(fSelect==1) cout << "TASK::Correlation with hadrons on ME "<< endl;
 603     if(fSelect==2) cout << "TASK::Correlation with kaons on ME "<< endl;
 604     if(fSelect==3) cout << "TASK::Correlation with kzeros on ME "<< endl;
 605   }
 606
 607
 608   if(!aod && AODEvent() && IsStandardAOD()) {
 609
 610     aod = dynamic_cast<AliAODEvent*> (AODEvent());
 611     AliAODHandler* aodHandler = (AliAODHandler*)
 612       ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
 613     if(aodHandler->GetExtensions()) {
 614       AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
 615       AliAODEvent *aodFromExt = ext->GetAOD();
 616       array3Prong=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm3Prong");
 617     }
 618   } else if(aod) {
 619     array3Prong=(TClonesArray*)aod->GetList()->FindObject("Charm3Prong");
 620   }
 621
 622   if(!aod || !array3Prong) {
 623     printf("AliAnalysisTaskSEDplusCorrelations::UserExec: Charm3Prong branch not found!\n");
 624     return;
 625   }
 626
 627
 628   // the AODs with null vertex pointer didn't pass the PhysSel
 629   if(!aod->GetPrimaryVertex()||TMath::Abs(aod->GetMagneticField())<0.001) return;
 630   fCounter->StoreEvent(aod,fRDCutsAnalysis,fReadMC);
 631   fHistNEvents->Fill(0); // count event
 632
 633
 634   Bool_t isEvSel=fRDCutsAnalysis->IsEventSelected(aod);
 635
 636
 637
 638   PostData(1,fOutput);
 639   if(!isEvSel)return;
 640   fHistNEvents->Fill(1);
 641
 642
 643
 644   TClonesArray *arrayMC=0;
 645   AliAODMCHeader *mcHeader=0;
 646   // AOD primary vertex
 647   AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex();
 648   //    vtx1->Print();
 649   TString primTitle = vtx1->GetTitle();
 650   //if(primTitle.Contains("VertexerTracks") && vtx1->GetNContributors()>0)fHistNEvents->Fill(2);
 651
 652   // load MC particles
 653   if(fReadMC){
 654
 655     arrayMC =  (TClonesArray*)aod->GetList()->FindObject(AliAODMCParticle::StdBranchName());
 656     if(!arrayMC) {
 657       printf("AliAnalysisTaskSEDplusCorrelations::UserExec: MC particles branch not found!\n");
 658       return;
 659     }
 660
 661     // load MC header
 662     mcHeader =  (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
 663     if(!mcHeader) {
 664       printf("AliAnalysisTaskSEDplusCorrelations::UserExec: MC header branch not found!\n");
 665       return;
 666     }
 667   }
 668
 669   //HFCorrelators initialization (for this event)
 670
 671   fCorrelator->SetAODEvent(aod); // set the event to be processedfCorrelator->
 672   Bool_t correlatorON = fCorrelator->Initialize(); //define the pool for mixing
 673   if(!correlatorON) {
 674     AliInfo("AliHFCorrelator didn't initialize the pool correctly or processed a bad event");
 675     return;
 676   }
 677   if(fReadMC) fCorrelator->SetMCArray(arrayMC);
 678
 679
 680
 681   Int_t n3Prong = array3Prong->GetEntriesFast();
 682   // printf("Number of D+->Kpipi: %d and of tracks: %d\n",n3Prong,aod->GetNumberOfTracks());
 683   Int_t nOS=0;
 684   Int_t index;
 685   Int_t pdgDgDplustoKpipi[3]={321,211,211};
 686   Int_t nSelectedloose=0,nSelectedtight=0;
 687
 688
 689
 690   for (Int_t i3Prong = 0; i3Prong < n3Prong; i3Prong++) {
 691     AliAODRecoDecayHF3Prong *d = (AliAODRecoDecayHF3Prong*)array3Prong->UncheckedAt(i3Prong);
 692     fHistNEvents->Fill(7);
 693
 694     if(fUseBit && !d->HasSelectionBit(AliRDHFCuts::kDplusCuts)){
 695       fHistNEvents->Fill(8);
 696       continue;
 697     }
 698
 699     Int_t passTightCuts=fRDCutsAnalysis->IsSelected(d,AliRDHFCuts::kAll,aod);
 700
 701     if(!fRDCutsAnalysis->GetIsSelectedCuts()) continue;
 702
 703     Bool_t unsetvtx=kFALSE;
 704     if(!d->GetOwnPrimaryVtx()){
 705       d->SetOwnPrimaryVtx(vtx1);
 706       unsetvtx=kTRUE;
 707     }
 708
 709
 710     Double_t ptCand = d->Pt();
 711     Int_t iPtBin = fRDCutsAnalysis->PtBin(ptCand);
 712     Bool_t recVtx=kFALSE;
 713     AliAODVertex *origownvtx=0x0;
 714     if(fRDCutsAnalysis->GetIsPrimaryWithoutDaughters()){
 715       if(d->GetOwnPrimaryVtx()) origownvtx=new AliAODVertex(*d->GetOwnPrimaryVtx());
 716       if(fRDCutsAnalysis->RecalcOwnPrimaryVtx(d,aod))recVtx=kTRUE;
 717       else fRDCutsAnalysis->CleanOwnPrimaryVtx(d,aod,origownvtx);
 718     }
 719
 720
 721
 722
 723     Int_t labDp=-1;
 724     Bool_t isDplus = kFALSE;
 725
 726     if(fReadMC){
 727       labDp = d->MatchToMC(411,arrayMC,3,pdgDgDplustoKpipi);
 728       if(labDp>=0){
 729         isDplus = kTRUE;
 730       }
 731     }
 732
 733
 734     Double_t invMass=d->InvMassDplus();
 735     Double_t rapid=d->YDplus();
 736     fYVsPt->Fill(ptCand,rapid);
 737     if(passTightCuts>0.) {fYVsPtTC->Fill(ptCand,rapid);nOS++;}
 738     //printf("****************: %d and of tracks: %d\n",n3Prong,nOS);
 739     Bool_t isFidAcc=fRDCutsAnalysis->IsInFiducialAcceptance(ptCand,rapid);
 740     if(isFidAcc){
 741       fPtVsMass->Fill(invMass,ptCand);
 742       if(passTightCuts>0) fPtVsMassTC->Fill(invMass,ptCand);
 743     }
 744
 745
 746     if(iPtBin>=0){
 747       index=GetHistoIndex(iPtBin);
 748       if(isFidAcc){
 749         fHistNEvents->Fill(9);
 750         nSelectedloose++;
 751         fMassHist[index]->Fill(invMass);
 752
 753         // loop for tight cuts
 754         if(passTightCuts>0){
 755           fHistNEvents->Fill(10);
 756           nSelectedtight++;
 757           fMassHistTC[index]->Fill(invMass);
 758
 759           //Dplus info
 760
 761           Double_t phiDplus = fCorrelator->SetCorrectPhiRange(d->Phi());
 762           fCorrelator->SetTriggerParticleProperties(d->Pt(),phiDplus,d->Eta());
 763
 764           Int_t nIDs[3] = {-9999999};
 765           nIDs[0] = ((AliAODTrack*)d->GetDaughter(0))->GetID();
 766           nIDs[1] = ((AliAODTrack*)d->GetDaughter(1))->GetID();
 767           nIDs[2] = ((AliAODTrack*)d->GetDaughter(2))->GetID();
 768
 769           Double_t ptlead = 0;
 770           Double_t philead = 0;
 771           Double_t refpt = 0;
 772
 773
 774
 775           Bool_t execPool = fCorrelator->ProcessEventPool();
 776
 777           //     printf("*************: %d\n",execPool);
 778           if(fMixing && !execPool) {
 779             AliInfo("Mixed event analysis: pool is not ready");
 780             continue;
 781           }
 782           Int_t NofEventsinPool = 1;
 783           if(fMixing) {
 784             NofEventsinPool = fCorrelator->GetNofEventsInPool();
 785           }
 786
 787
 788
 789           for (Int_t jMix =0; jMix < NofEventsinPool; jMix++){// loop on events in the pool; if it is SE analysis, stops at one
 790
 791             Bool_t analyzetracks = fCorrelator->ProcessAssociatedTracks(jMix);
 792             if(!analyzetracks) {
 793               AliInfo("AliHFCorrelator::Cannot process the track array");
 794               continue;
 795             }
 796
 797             //start the track loop
 798
 799             // Int_t NofTracks = fCorrelator->GetNofTracks();
 800
 801             for (Int_t iTrack = 0;iTrack<fCorrelator->GetNofTracks();iTrack++) {
 802               Bool_t runcorrelation = fCorrelator->Correlate(iTrack);
 803
 804               if(!runcorrelation) continue;
 805               Double_t DeltaPhi = fCorrelator->GetDeltaPhi();
 806               Double_t DeltaEta = fCorrelator->GetDeltaEta();
 807
 808               AliReducedParticle* redpart = fCorrelator->GetAssociatedParticle();
 809               Double_t phiHad=redpart->Phi();
 810               Double_t ptHad=redpart->Pt();
 811               Int_t label = redpart->GetLabel();
 812               Int_t trackid = redpart->GetID();
 813               phiHad = fCorrelator->SetCorrectPhiRange(phiHad);
 814
 815
 816               //  discard the dplus daughters
 817               if (!fMixing){
 818                 if( trackid == nIDs[0] || trackid == nIDs[1] || trackid == nIDs[2]) continue;
 819               }
 820               // discard the negative id tracks
 821               if(trackid < 0) continue;
 822
 823
 824               FillCorrelations(d,DeltaPhi,DeltaEta,index,fSelect);
 825
 826               // For leading particle
 827
 828               if (ptHad > refpt) {
 829                 refpt = ptHad; ptlead = ptHad;
 830                 philead = d->Phi() - phiHad;
 831                 if (philead < (-1)*TMath::Pi()/2) philead += 2*TMath::Pi();
 832                 if (philead > 3*TMath::Pi()/2) philead -= 2*TMath::Pi();
 833
 834               }
 835
 836               // montecarlo
 837
 838               if(fReadMC && isDplus) {
 839
 840                 index=GetSignalHistoIndex(iPtBin);
 841
 842                 Int_t PartSource = fCuts->IsMCpartFromHF(label,arrayMC); // check source of associated particle (hadron/kaon/K0)
 843                 FillMCCorrelations(d,DeltaPhi,DeltaEta,index,PartSource,fSelect);
 844
 845
 846               } // readMC
 847
 848             }//count good tracks
 849
 850             // For leading particle
 851             fMassVsdPhiHistLeadHad[index]->Fill(invMass,philead);
 852             fLeadPt[index]->Fill(ptlead);
 853
 854             if(fReadMC && isDplus) {
 855               index=GetSignalHistoIndex(iPtBin);
 856               fMassVsdPhiHistLeadHad[index]->Fill(invMass,philead);
 857               fLeadPt[index]->Fill(ptlead);
 858
 859             }
 860
 861           }//jmix
 862
 863         }// tc
 864       }//fid
 865
 866
 867
 868
 869
 870     }
 871     if(recVtx)fRDCutsAnalysis->CleanOwnPrimaryVtx(d,aod,origownvtx);
 872
 873     if(unsetvtx) d->UnsetOwnPrimaryVtx();
 874
 875   }
 876
 877   if(fMixing){
 878     Bool_t updated = fCorrelator->PoolUpdate();
 879
 880     if(!updated) AliInfo("Pool was not updated");
 881   }
 882   fCounter->StoreCandidates(aod,nSelectedloose,kTRUE);
 883   fCounter->StoreCandidates(aod,nSelectedtight,kFALSE);
 884
 885
 886   PostData(1,fOutput);
 887   PostData(2,fListCuts);
 888   PostData(3,fCounter);
 889   return;
 890 }
 891
 892 //________________________________________________________________________
 893 void AliAnalysisTaskSEDplusCorrelations::Terminate(Option_t */*option*/)
 894 {
 895   // Terminate analysis
 896   //
 897   if(fDebug > 1) printf("AnalysisTaskSEDplusCorrelations: Terminate() \n");
 898
 899   fOutput = dynamic_cast<TList*> (GetOutputData(1));
 900   if (!fOutput) {
 901     printf("ERROR: fOutput not available\n");
 902     return;
 903   }
 904   fHistNEvents = dynamic_cast<TH1F*>(fOutput->FindObject("fHistNEvents"));
 905
 906   TString hisname;
 907   Int_t index=0;
 908
 909   for(Int_t i=0;i<fNPtBins;i++){
 910     index=GetHistoIndex(i);
 911
 912     hisname.Form("hMassPt%dTC",i);
 913     fMassHistTC[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
 914   }
 915
 916   TCanvas *c1=new TCanvas("c1","D+ invariant mass distribution",500,500);
 917   c1->cd();
 918   TH1F *hMassPt=(TH1F*)fOutput->FindObject("hMassPt3TC");
 919   hMassPt->SetLineColor(kBlue);
 920   hMassPt->SetXTitle("M[GeV/c^{2}]");
 921   hMassPt->Draw();
 922
 923   return;
 924 }
 925
 926
 927 //________________________________________________________________________
 928 void AliAnalysisTaskSEDplusCorrelations::FillCorrelations(AliAODRecoDecayHF3Prong* d, Double_t deltaPhi, Double_t deltaEta, Int_t ind, Int_t sel) const{
 929   //Filling histogams
 930
 931   Double_t invMass=d->InvMassDplus();
 932
 933
 934   if(sel==1){
 935     fMassVsdPhiHistHad[ind]->Fill(invMass,deltaPhi);
 936     fMassVsdEtaHistHad[ind]->Fill(invMass,deltaEta);
 937   }
 938   if(sel==2){
 939     fMassVsdPhiHistKaon[ind]->Fill(invMass,deltaPhi);
 940     fMassVsdEtaHistKaon[ind]->Fill(invMass,deltaEta);
 941   }
 942   if(sel==3){
 943     fMassVsdPhiHistKshort[ind]->Fill(invMass,deltaPhi);
 944     fMassVsdEtaHistKshort[ind]->Fill(invMass,deltaEta);
 945   }
 946
 947   return;
 948 }
 949
 950
 951
 952 //________________________________________________________________________
 953 void AliAnalysisTaskSEDplusCorrelations::FillMCCorrelations(AliAODRecoDecayHF3Prong* d, Double_t deltaPhi, Double_t deltaEta, Int_t ind,Int_t mcSource, Int_t sel) const{
 954   // Filling histos with MC information
 955
 956   Double_t invMass=d->InvMassDplus();
 957
 958
 959   if(sel==1){
 960     fMassVsdPhiHistHad[ind]->Fill(invMass,deltaPhi);
 961     fMassVsdEtaHistHad[ind]->Fill(invMass,deltaEta);
 962
 963
 964     fMCSources->Fill(0);
 965
 966     if (mcSource==44){ // is from charm ->D
 967       fMCSources->Fill(1);
 968       fMCSources->Fill(2);
 969     }
 970     if (mcSource==54){ // is from beauty -> D
 971       fMCSources->Fill(1);
 972       fMCSources->Fill(3);
 973     }
 974     if (mcSource==55){ // is from beauty ->B
 975       fMCSources->Fill(1);
 976       fMCSources->Fill(4);
 977     }
 978   }
 979
 980   if(sel==2){
 981     fMassVsdPhiHistKaon[ind]->Fill(invMass,deltaPhi);
 982     fMassVsdEtaHistKaon[ind]->Fill(invMass,deltaEta);
 983     fKaonOrigin->Fill(0);
 984     if (mcSource==44){ // is from charm ->D
 985       fKaonOrigin->Fill(1);
 986       fKaonOrigin->Fill(2);
 987     }
 988     if (mcSource==54){ // is from beauty -> D
 989       fKaonOrigin->Fill(1);
 990       fKaonOrigin->Fill(3);
 991     }
 992     if (mcSource==55){ // is from beauty ->B
 993       fKaonOrigin->Fill(1);
 994       fKaonOrigin->Fill(4);
 995     }
 996   }
 997   if(sel==3){
 998     fMassVsdPhiHistKshort[ind]->Fill(invMass,deltaPhi);
 999     fMassVsdEtaHistKshort[ind]->Fill(invMass,deltaEta);
1000     fK0Origin->Fill(0);
1001     if (mcSource==44){ // is from charm ->D
1002       fK0Origin->Fill(1);
1003       fK0Origin->Fill(2);
1004     }
1005     if (mcSource==54){ // is from beauty -> D
1006       fK0Origin->Fill(1);
1007       fK0Origin->Fill(3);
1008     }
1009     if (mcSource==55){ // is from beauty ->B
1010       fK0Origin->Fill(1);
1011       fK0Origin->Fill(4);
1012     }
1013
1014   }
1015
1016   return;
1017 }
1018
1019
1020
1021
1022
1023
1024
1025
1026
1027