PWGHF/correlationHF/macros/makeTFile4CutsD0Correlations.C

   1 #include <Riostream.h>
   2 #include <TFile.h>
   3 #include <AliRDHFCutsD0toKpi.h>
   4 #include <TClonesArray.h>
   5 #include <TParameter.h>
   6 #include "AliHFAssociatedTrackCuts.h"
   7
   8 //Use:
   9 //Set hard coded commentet with //set this!!
  10 // root[] .L makeInputD0tasks.C++
  11 // root[] makeInputAliAnalysisTaskSED0Correlations()
  12 //similar macros for the other D mesons
  13
  14 //Author: Fabio Colamaria, fabio.colamaria@ba.infn.it
  15
  16 //macro to make a .root file which contains an AliRDHFCutsD0toKpi for AliAnalysisTaskSED0Mass task
  17
  18 void makeInputAliAnalysisTaskSED0Correlations(){
  19
  20 //____________________________________________________
  21
  22 // Cuts for D0 cuts
  23
  24   AliRDHFCutsD0toKpi* RDHFD0Corr=new AliRDHFCutsD0toKpi();
  25   RDHFD0Corr->SetName("D0toKpiCuts");
  26   RDHFD0Corr->SetTitle("Cuts for D0 analysis");
  27
  28   RDHFD0Corr->SetMinVtxContr(1);
  29
  30   //Quality tracks for daughters
  31   AliESDtrackCuts* esdTrackCuts=new AliESDtrackCuts();
  32   esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
  33   esdTrackCuts->SetRequireTPCRefit(kTRUE);
  34   esdTrackCuts->SetRequireITSRefit(kTRUE);
  35   //esdTrackCuts->SetMinNClustersITS(4); // default is 5
  36   //esdTrackCuts->SetMinNClustersTPC(120);
  37   esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD, AliESDtrackCuts::kAny); // default is kBoth, otherwise kAny
  38   esdTrackCuts->SetMinDCAToVertexXY(0.);
  39   esdTrackCuts->SetEtaRange(-0.8,0.8);
  40   esdTrackCuts->SetPtRange(0.3,1.e10);
  41
  42   RDHFD0Corr->AddTrackCuts(esdTrackCuts);
  43
  44   //D0 selection topological cuts
  45   const Int_t nptbins =14;
  46   const Double_t ptmax = 9999.;
  47   const Int_t nvars=11;
  48   Float_t ptbins[nptbins+1];
  49   ptbins[0]=0.;
  50   ptbins[1]=0.5;
  51   ptbins[2]=1.;
  52   ptbins[3]=2.;
  53   ptbins[4]=3.;
  54   ptbins[5]=4.;
  55   ptbins[6]=5.;
  56   ptbins[7]=6.;
  57   ptbins[8]=7.;
  58   ptbins[9]=8.;
  59   ptbins[10]=12.;
  60   ptbins[11]=16.;
  61   ptbins[12]=20.;
  62   ptbins[13]=24.;
  63   ptbins[14]=ptmax;
  64
  65   RDHFD0Corr->SetGlobalIndex(nvars,nptbins);
  66   RDHFD0Corr->SetPtBins(nptbins+1,ptbins);
  67
  68   Float_t cutsMatrixD0toKpiStand[nptbins][nvars]={{0.400,350.*1E-4,0.8,0.5,0.5,1000.*1E-4,1000.*1E-4,-5000.*1E-8,0.80,0.,0.},/* pt<0.5*/
  69                                                   {0.400,350.*1E-4,0.8,0.5,0.5,1000.*1E-4,1000.*1E-4,-5000.*1E-8,0.80,0.,0.},/* 0.5<pt<1*/
  70                                                   {0.400,300.*1E-4,0.8,0.4,0.4,1000.*1E-4,1000.*1E-4,-25000.*1E-8,0.80,0.,0.},/* 1<pt<2 */
  71                                                   {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-15000.*1E-8,0.85,0.,0.},/* 2<pt<3 */
  72                                                   {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-8000.*1E-8,0.85,0.,0.},/* 3<pt<4 */
  73                                                   {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-8000.*1E-8,0.85,0.,0.},/* 4<pt<5 */
  74                                                   {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-8000.*1E-8,0.85,0.,0.},/* 5<pt<6 */
  75                                                   {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-8000.*1E-8,0.85,0.,0.},/* 6<pt<7 */
  76                                                   {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-7000.*1E-8,0.85,0.,0.},/* 7<pt<8 */
  77                                                   {0.400,300.*1E-4,0.9,0.7,0.7,1000.*1E-4,1000.*1E-4,-5000.*1E-8,0.85,0.,0.},/* 8<pt<12 */
  78                                                   {0.400,300.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,10000.*1E-8,0.85,0.,0.},/* 12<pt<16 */
  79                                                   {0.400,300.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,999999.*1E-8,0.85,0.,0.},/* 16<pt<20 */
  80                                                   {0.400,300.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,999999.*1E-8,0.85,0.,0.},/* 20<pt<24 */
  81                                                   {0.400,300.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,999999.*1E-8,0.85,0.,0.}};/* pt>24 */
  82
  83
  84   //CREATE TRANSPOSE MATRIX...REVERSE INDICES as required by AliRDHFCuts
  85   Float_t **cutsMatrixTransposeStand=new Float_t*[nvars];
  86   for(Int_t iv=0;iv<nvars;iv++)cutsMatrixTransposeStand[iv]=new Float_t[nptbins];
  87
  88   for (Int_t ibin=0;ibin<nptbins;ibin++){
  89     for (Int_t ivar = 0; ivar<nvars; ivar++){
  90       cutsMatrixTransposeStand[ivar][ibin]=cutsMatrixD0toKpiStand[ibin][ivar];
  91     }
  92   }
  93
  94   RDHFD0Corr->SetCuts(nvars,nptbins,cutsMatrixTransposeStand);
  95   RDHFD0Corr->SetUseSpecialCuts(kTRUE);
  96   RDHFD0Corr->SetRemoveDaughtersFromPrim(kTRUE);
  97
  98   for(Int_t iv=0;iv<nvars;iv++) delete [] cutsMatrixTransposeStand[iv];
  99   delete [] cutsMatrixTransposeStand;
 100   cutsMatrixTransposeStand=NULL;
 101
 102   //D0 pid settings
 103   Bool_t pidflag=kTRUE;
 104   RDHFD0Corr->SetUsePID(pidflag);
 105   if(pidflag) cout<<"PID is used"<<endl;
 106   else cout<<"PID is not used"<<endl;
 107
 108   AliAODPidHF* pidObj=new AliAODPidHF();
 109   Int_t mode=1;
 110   const Int_t nlims=2;
 111   Double_t plims[nlims]={0.6,0.8}; //TPC limits in momentum [GeV/c]
 112   Bool_t compat=kTRUE; //effective only for this mode
 113   Bool_t asym=kTRUE;
 114   Double_t sigmas[5]={2.,1.,0.,3.,0.}; //to be checked and to be modified with new implementation of setters by Rossella
 115   pidObj->SetAsym(asym);// if you want to use the asymmetric bands in TPC
 116   pidObj->SetMatch(mode);
 117   pidObj->SetPLimit(plims,nlims);
 118   pidObj->SetSigma(sigmas);
 119   pidObj->SetCompat(compat);
 120   pidObj->SetPCompatTOF(1.5);
 121   pidObj->SetSigmaForTPCCompat(3.);
 122   pidObj->SetSigmaForTOFCompat(3.);
 123   pidObj->SetTPC(kTRUE);
 124   pidObj->SetTOF(kTRUE);
 125   RDHFD0Corr->SetPidHF(pidObj);
 126   RDHFD0Corr->SetUsePID(kTRUE);
 127   RDHFD0Corr->SetUseDefaultPID(kFALSE); //to use the AliAODPidHF
 128   RDHFD0Corr->SetLowPt(kFALSE);
 129
 130   //activate pileup rejection (for pp)
 131   RDHFD0Corr->SetOptPileup(AliRDHFCuts::kRejectPileupEvent);
 132
 133   TString cent="";
 134   //centrality selection (Pb-Pb)
 135   Float_t minc=0,maxc=100;
 136   RDHFD0Corr->SetMinCentrality(minc);
 137   RDHFD0Corr->SetMaxCentrality(maxc);
 138   cent=Form("%.0f%.0f",minc,maxc);
 139   RDHFD0Corr->SetUseCentrality(AliRDHFCuts::kCentOff); //kCentOff,kCentV0M,kCentTRK,kCentTKL,kCentCL1,kCentInvalid
 140
 141   //temporary
 142   //RDHFD0Corr->SetFixRefs();
 143
 144   cout<<"This is the object I'm going to save:"<<endl;
 145   RDHFD0Corr->PrintAll();
 146   TFile* fout=new TFile("D0toKpiCuts.root","recreate");   //set this!!
 147
 148   fout->cd();
 149   RDHFD0Corr->Write();
 150   fout->Close();
 151
 152 //____________________________________________________
 153
 154   // Cuts for correlated tracks
 155
 156   AliHFAssociatedTrackCuts* HFCorrelationCuts=new AliHFAssociatedTrackCuts();
 157   HFCorrelationCuts->SetName("AssociatedTrkCuts");
 158   HFCorrelationCuts->SetTitle("Cuts for associated track");
 159   Float_t eta = 0.9;
 160
 161   // Set quality cuts on tracks
 162   AliESDtrackCuts *esdHadrCuts = new AliESDtrackCuts("AliESDHadrCuts","default");
 163   esdHadrCuts->SetRequireSigmaToVertex(kFALSE);
 164   esdHadrCuts->SetRequireTPCRefit(kTRUE);
 165   esdHadrCuts->SetRequireITSRefit(kTRUE);
 166   esdHadrCuts->SetMinNClustersITS(2); //as for D*
 167   esdHadrCuts->SetMinNClustersTPC(80); //as for D*
 168 //esdHadrCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny);
 169   esdHadrCuts->SetMinDCAToVertexXY(0.);
 170   esdHadrCuts->SetEtaRange(-eta,eta);
 171   esdHadrCuts->SetPtRange(0.,1.e10);
 172   HFCorrelationCuts->AddTrackCuts(esdHadrCuts);
 173
 174   // Set kinematics cuts for AOD track
 175   const int nofcuts = 4;
 176   Float_t* trackcutsarray;
 177   trackcutsarray=new Float_t[nofcuts];
 178   trackcutsarray[0] = 0.3;//track min pt
 179   trackcutsarray[1] = 10000.;//track max pt
 180   trackcutsarray[2] = -99999999.;//track min impact parameter. DON'T put 0 since default value is -999999. and it would skip all tracks if d0 is not calculated!
 181   trackcutsarray[3] = 10000.;//track max impact parameter
 182   HFCorrelationCuts->SetNVarsTrack(nofcuts);
 183   HFCorrelationCuts->SetAODTrackCuts(trackcutsarray);
 184
 185   HFCorrelationCuts->SetCharge(0); // -1/+1 to look for opposite/same charge, 0 no charge selection
 186   HFCorrelationCuts->SetFilterBit(0); // set 0 for analysis with AOD from 2010
 187
 188   // Set kinematics cuts for AOD v0
 189   const int nofcuts2 = 7;
 190   Float_t* vzerocutsarray;
 191   vzerocutsarray=new Float_t[nofcuts2];
 192   vzerocutsarray[0] = 0.2; // max dca between two daugters (cm)
 193   vzerocutsarray[1] = 2; //  max chi square
 194   vzerocutsarray[2] = 2.; // min decay length (cm)
 195   vzerocutsarray[3] = 15; // max decay length (cm)
 196   vzerocutsarray[4] = 100.; // max opening angle between two daugters
 197   vzerocutsarray[5] = 0; // min pt of k0 (GeV/c)
 198   vzerocutsarray[6] = 0.9; // set eta acceptance
 199   HFCorrelationCuts->SetNVarsVzero(nofcuts2);
 200   HFCorrelationCuts->SetAODvZeroCuts(vzerocutsarray);
 201
 202   // Set PID
 203   Int_t mode =1;
 204   Double_t ptlimit[2] = {0.6,0.8};
 205   AliAODPidHF* pidObj=new AliAODPidHF();
 206   pidObj->SetMatch(1);  //A.Rossi mode
 207   pidObj->SetAsym(kTRUE);
 208   pidObj->SetPLimit(ptlimit);
 209   pidObj->SetSigma(0,2.);  //TPC sigma, in three pT ranges
 210   pidObj->SetSigma(1,1.);
 211   pidObj->SetSigma(2,0.);
 212   pidObj->SetSigma(3,3.);  //TOF sigma, whole pT range
 213   pidObj->SetPCompatTOF(1.5);
 214   pidObj->SetSigmaForTPCCompat(3.);
 215   pidObj->SetSigmaForTOFCompat(3.);
 216   pidObj->SetTPC(kTRUE);
 217   pidObj->SetTOF(kTRUE);
 218   pidObj->SetCompat(kTRUE);
 219   HFCorrelationCuts->SetPidHF(pidObj);
 220
 221   //Event Mixing settings
 222   HFCorrelationCuts->SetMaxNEventsInPool(200);
 223   HFCorrelationCuts->SetMinNTracksInPool(1000);
 224   HFCorrelationCuts->SetMinEventsToMix(8);
 225   HFCorrelationCuts->SetNofPoolBins(5,5);
 226
 227   Double_t MBins[]={0,20,40,60,80,500};
 228   Double_t * MultiplicityBins = MBins;
 229   Double_t ZBins[]={-10,-5,-2.5,2.5,5,10};
 230   Double_t *ZVrtxBins = ZBins;
 231
 232   HFCorrelationCuts->SetPoolBins(ZVrtxBins,MultiplicityBins);
 233
 234   // Save to *.root file
 235   HFCorrelationCuts->PrintAll();
 236   TFile* fout=new TFile("AssocPartCuts.root","recreate");   //set this!!
 237   fout->cd();
 238   HFCorrelationCuts->Write();
 239   fout->Close();
 240
 241 }
 242