PWGHF/hfe/macros/configs/pPb/ConfigHFEpPb.C

   1 Bool_t ReadContaminationFunctions(TString filename, TF1 **functions, double sigma){
   2   TFile *in = TFile::Open(Form("$ALICE_ROOT/PWGHF/hfe/macros/configs/pPb/%s", filename.Data()));
   3   gROOT->cd();
   4   int isig = static_cast<int>(sigma * 100.);
   5   printf("Getting hadron background for the sigma cut: %d\n", isig);
   6   bool status = kTRUE;
   7   for(int icent = 0; icent < 12; icent++){
   8     functions[icent] = dynamic_cast<TF1 *>(in->Get(Form("hback_%d_%d", isig, icent)));
   9     if(functions[icent]) printf("Config for centrality class %d found\n", icent);
  10     else{
  11       printf("Config for the centrality class %d not found\n", icent);
  12       status = kFALSE;
  13     }
  14   }
  15   delete in;
  16   return status;
  17 }
  18
  19 AliAnalysisTaskHFE* ConfigHFEpPb(Bool_t useMC, Bool_t isAOD, TString appendix,
  20                                  UChar_t TPCcl=70, UChar_t TPCclPID = 80,
  21                                  UChar_t ITScl=3, Double_t DCAxy=1000., Double_t DCAz=1000.,
  22                                  Double_t* tpcdEdxcutlow=NULL,Double_t* tpcdEdxcuthigh=NULL,
  23                                  Double_t TOFs=3., Int_t TOFmis=0,
  24                                  Int_t itshitpixel = 0, Int_t icent=1,
  25                                  Double_t etami=-0.8, Double_t etama=0.8, Bool_t esdPreselection = kFALSE){
  26
  27   Bool_t kAnalyseTaggedTracks = kTRUE;
  28   Bool_t kApplyPreselection = isAOD ? kFALSE : esdPreselection;
  29
  30   //***************************************//
  31   //        Setting up the HFE cuts        //
  32   //***************************************//
  33
  34   AliHFEcuts *hfecuts = new AliHFEcuts(appendix,"HFE cuts for pPb");
  35   //hfecuts->SetQAOn();
  36   hfecuts->CreateStandardCuts();
  37   hfecuts->SetMinNClustersTPC(TPCcl);
  38   hfecuts->SetMinNClustersTPCPID(TPCclPID);
  39   hfecuts->SetMinNClustersITS(ITScl);
  40   hfecuts->SetMinRatioTPCclusters(0.6);
  41   hfecuts->SetTPCmodes(AliHFEextraCuts::kFound, AliHFEextraCuts::kFoundOverFindable);
  42   hfecuts->SetCutITSpixel(itshitpixel);
  43   hfecuts->SetCheckITSLayerStatus(kFALSE);
  44   hfecuts->SetEtaRange(etami,etama);
  45   hfecuts->SetRejectKinkDaughters();
  46   hfecuts->SetAcceptKinkMothers();
  47   if(isAOD) hfecuts->SetAODFilterBit(4);
  48
  49   //if((iPixelAny==AliHFEextraCuts::kAny) || (iPixelAny==AliHFEextraCuts::kSecond))
  50   //hfecuts->SetProductionVertex(0,7,0,7);
  51
  52   hfecuts->SetMaxImpactParam(DCAxy,DCAz);
  53   hfecuts->SetUseMixedVertex(kTRUE);
  54   hfecuts->SetVertexRange(10.);
  55   // New pPb cuts (February 2013)
  56   hfecuts->SetUseCorrelationVertex();
  57   hfecuts->SetSPDVtxResolutionCut();
  58   hfecuts->SetpApileupCut();
  59
  60   // TOF settings:
  61   Int_t usetof=0;
  62   Bool_t kTOFmis=kFALSE;
  63   if (TOFs>0.){
  64     usetof = 1;
  65     printf("CONFIGURATION FILE: TOF is used \n");
  66     hfecuts->SetTOFPIDStep(kTRUE);
  67     printf("CONFIGURATION FILE: TOF PID step is requested !!!! \n");
  68     if (TOFmis>0){
  69       kTOFmis = kTRUE;
  70       printf("CONFIGURATION FILE: TOF mismatch rejection is set ON \n");
  71     }
  72   }
  73
  74   //***************************************//
  75   //        Setting up the task            //
  76   //***************************************//
  77
  78   AliAnalysisTaskHFE *task = new AliAnalysisTaskHFE(Form("HFEtask%s",appendix.Data()));
  79   printf("task %p\n", task);
  80   task->SetpPbAnalysis();
  81   if(!isAOD) task->SetRemoveFirstEventInChunk();
  82   task->SetRemovePileUp(kFALSE);
  83   task->SetHFECuts(hfecuts);
  84   task->GetPIDQAManager()->SetHighResolutionHistos();
  85   task->SetRejectKinkMother(kFALSE);
  86
  87   // Determine the centrality estimator
  88   task->SetCentralityEstimator("V0A");
  89   if (icent == 2) task->SetCentralityEstimator("V0M");
  90   else if (icent == 3) task->SetCentralityEstimator("CL1");
  91   else if (icent == 4) task->SetCentralityEstimator("ZNA");
  92
  93   //***************************************//
  94   //        Prepare preselection           //
  95   // This mimics the ESD->AOD filter in    //
  96   // case of the ESD analysis and selects  //
  97   // only tracks which will be selected in //
  98   // the AOD analysis with the given filter//
  99   // bit. Not to be applied for AODS.      //
 100   // For pPb the cuts used are (bit 4)     //
 101   // esdTrackCutsHG0 from file $ALICE_ROOT///
 102   // ANALYSIS/macros/AddTaskESDFilter.C    //
 103   //***************************************//
 104
 105   if(kApplyPreselection){
 106     AliESDtrackCuts* esdfilter = AliESDtrackCuts::GetStandardITSTPCTrackCuts2011(kFALSE);
 107     esdfilter->SetMaxDCAToVertexXY(2.4);
 108     esdfilter->SetMaxDCAToVertexZ(3.2);
 109     esdfilter->SetDCAToVertex2D(kTRUE);
 110
 111     task->SetHFECutsPreselect(esdfilter);
 112     printf("Put a preselection cut\n");
 113     task->SetFillNoCuts(kTRUE);
 114   }
 115
 116   //***************************************//
 117   //          Variable manager             //
 118   //***************************************//
 119   // Define Variables
 120   Double_t ptbinning[36] = {0., 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1., 1.1, 1.2, 1.3, 1.4, 1.5, 1.75, 2., 2.25, 2.5, 2.75, 3., 3.5, 4., 4.5, 5., 5.5, 6., 7., 8., 10., 12., 14., 16., 18., 20.};
 121   Double_t etabinning[17] = {-0.8, -0.7, -0.6, -0.5, -0.4, -0.3, -0.2, -0.1, 0., 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8};
 122
 123   Int_t sizept=(sizeof(ptbinning)/sizeof(double))-1;
 124   Int_t sizeeta=(sizeof(etabinning)/sizeof(double))-1;
 125
 126   AliHFEvarManager *vm = task->GetVarManager();
 127   vm->AddVariable("pt", sizept, ptbinning);
 128   vm->AddVariable("eta", sizeeta, -0.8,0.8);
 129   vm->AddVariable("phi",21, -0, 2*TMath::Pi());
 130   vm->AddVariable("charge");
 131   vm->AddVariable("source");
 132   vm->AddVariable("centrality");
 133
 134   // For the moment, remove the part dedicated to the background subtraction.
 135   // It will be implemented in a different way, reading it from a root file.
 136
 137   //***************************************//
 138   //          Configure the PID            //
 139   //***************************************//
 140
 141   // Define PID
 142   AliHFEpid *pid = task->GetPID();
 143   if(useMC) pid->SetHasMCData(kTRUE);
 144
 145   if (usetof){
 146     pid->AddDetector("TOF", 0);
 147     pid->AddDetector("TPC", 1);
 148   } else {
 149     pid->AddDetector("TPC", 0);
 150   }
 151
 152   // Configure TPC PID
 153   // do the identical thing in data and MC
 154   Double_t paramsTPCdEdxcutlow[12] ={0.0, 0.0, 0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0};
 155   if(tpcdEdxcutlow) memcpy(paramsTPCdEdxcutlow,tpcdEdxcutlow,sizeof(paramsTPCdEdxcutlow));
 156
 157   Double_t paramsTPCdEdxcuthigh[12] ={3.0, 3.0, 3.0,3.0,3.0,3.0,3.0,3.0,3.0,3.0,3.0,3.0};
 158   if(tpcdEdxcuthigh) memcpy(paramsTPCdEdxcuthigh,tpcdEdxcuthigh,sizeof(paramsTPCdEdxcuthigh));
 159
 160   char *cutmodel;
 161   cutmodel="pol0";
 162
 163   for(Int_t a=0;a<11;a++){
 164     // Not necessary anymore, since the pPb case is handled similarly to the pp case
 165     //   cout << a << " " << paramsTPCdEdxcut[a] << endl;
 166     Double_t tpcparamlow[1]={paramsTPCdEdxcutlow[a]};
 167     Float_t tpcparamhigh=paramsTPCdEdxcuthigh[a];
 168     pid->ConfigureTPCcentralityCut(a,cutmodel,tpcparamlow,tpcparamhigh);
 169   }
 170   pid->ConfigureTPCdefaultCut(cutmodel,paramsTPCdEdxcutlow,paramsTPCdEdxcuthigh[0]); // After introducing the pPb flag, pPb is merged with pp and this line defines the cut
 171
 172   // Configure TOF PID
 173   if (usetof){
 174     pid->ConfigureTOF(TOFs);
 175     AliHFEpidTOF *tofpid = pid->GetDetPID(AliHFEpid::kTOFpid);
 176     if (kTOFmis){
 177       tofpid->SetRejectTOFmismatch();
 178     }
 179   }
 180
 181   // To make different upper TOF cut to see contamination effect
 182   // The below two lines should be removed after this check
 183   //AliHFEpidTOF *tofpid = pid->GetDetPID(AliHFEpid::kTOFpid);
 184   //if(TOFs<3.) tofpid->SetTOFnSigmaBand(-3,TOFs); //only to check the assymmetric tof cut
 185
 186   // Load hadron background
 187   if(!useMC){
 188     Bool_t status = kTRUE;
 189     TF1 *hBackground[12];
 190     status = ReadContaminationFunctions("hadroncontamination_TOFTPC_pPb_eta06_newsplines_try3.root", hBackground, tpcdEdxcutlow[0]);
 191     for(Int_t a=0;a<12;a++) {
 192       //printf("back %f \n",hBackground[a]);
 193       if(status) task->SetBackGroundFactorsFunction(hBackground[a],a);
 194       else printf("not all background functions found\n");
 195     }
 196   }
 197
 198   //***************************************//
 199   //          V0 tagged tracks             //
 200   //***************************************//
 201
 202   if(kAnalyseTaggedTracks){
 203     AliHFEcuts *v0trackCuts = new AliHFEcuts("V0trackCuts", "Track Cuts for tagged track Analysis");
 204     v0trackCuts->CreateStandardCuts();
 205     v0trackCuts->SetMinNClustersTPC(TPCcl);
 206     v0trackCuts->SetMinNClustersTPCPID(TPCclPID);
 207     v0trackCuts->SetMinRatioTPCclusters(0.6);
 208     v0trackCuts->SetTPCmodes(AliHFEextraCuts::kFound, AliHFEextraCuts::kFoundOverFindable);
 209     v0trackCuts->SetMinNClustersITS(1);
 210     v0trackCuts->SetCutITSpixel(AliHFEextraCuts::kFirst);
 211     v0trackCuts->SetCheckITSLayerStatus(kFALSE);
 212     v0trackCuts->UnsetVertexRequirement();
 213     //hfecuts->SetSigmaToVertex(10);
 214     if(usetof) v0trackCuts->SetTOFPIDStep(kTRUE);
 215     v0trackCuts->SetQAOn();
 216
 217     task->SwitchOnPlugin(AliAnalysisTaskHFE::kTaggedTrackAnalysis);
 218     task->SetTaggedTrackCuts(v0trackCuts);
 219     task->SetCleanTaggedTrack(kTRUE);
 220   }
 221
 222   // QA
 223   printf("task %p\n", task);
 224   task->SetQAOn(AliAnalysisTaskHFE::kPIDqa);
 225   task->SetQAOn(AliAnalysisTaskHFE::kMCqa);
 226   task->SwitchOnPlugin(AliAnalysisTaskHFE::kDEstep);
 227
 228   printf("*************************************\n");
 229   printf("Configuring standard Task:\n");
 230   task->PrintStatus();
 231   pid->PrintStatus();
 232   printf("*************************************\n");
 233   return task;
 234 }