PWGHF/hfe/macros/configs/pPb/ConfigHFEpPbTRD.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])
  10     {
  11 //      new TCanvas;
  12 //        functions[0]->Draw();
  13         printf("Config for centrality class %d found\n", icent);
  14     }
  15     else{
  16       printf("Config for the centrality class %d not found\n", icent);
  17       status = kFALSE;
  18     }
  19   }
  20   delete in;
  21   return status;
  22 }
  23
  24 AliAnalysisTaskHFE* ConfigHFEpPbTRD(Bool_t useMC, Bool_t isAOD, TString appendix,
  25                                  UChar_t TPCcl=70, UChar_t TPCclPID = 80,
  26                                  UChar_t ITScl=3, Double_t DCAxy=1000., Double_t DCAz=1000.,
  27                                  Double_t* tpcdEdxcutlow=NULL,Double_t* tpcdEdxcuthigh=NULL,
  28                                  Double_t TOFs=3., Int_t TOFmis=0,
  29                                  Int_t itshitpixel = 0, Int_t icent=1,
  30                                  Double_t etami=-0.8, Double_t etama=0.8,
  31                                  Int_t TRDtrigger=1, Int_t trdpidmethod=1, Int_t TRDtl=6, Int_t TRDeff=4,
  32                                  TString detector){
  33
  34     Bool_t kAnalyseTaggedTracks = kTRUE;
  35
  36     // TRD settings
  37     Float_t eeff[6] = {0.7, 0.75, 0.8, 0.85, 0.9, 0.95};
  38     Int_t eeffint[6] = {70, 75, 80, 85, 90, 95};
  39     //  if(TRDeff >= 6 || TRDtl < 4 || TRDtl > 6) return NULL;
  40     if(TRDeff >= 6 || TRDtl > 6) return NULL;
  41     printf("TRD settings: %i %f \n",TRDtl, eeff[TRDeff]);
  42
  43   //***************************************//
  44   //        Setting up the HFE cuts        //
  45   //***************************************//
  46
  47   AliHFEcuts *hfecuts = new AliHFEcuts(appendix,"HFE cuts for pPb with TRD");
  48
  49  // hfecuts->SetQAOn();
  50   hfecuts->CreateStandardCuts();
  51   hfecuts->SetMinNClustersTPC(TPCcl);
  52   hfecuts->SetMinNClustersTPCPID(TPCclPID);
  53   hfecuts->SetMinNClustersITS(ITScl);
  54   hfecuts->SetMinRatioTPCclusters(0.6);
  55   hfecuts->SetTPCmodes(AliHFEextraCuts::kFound, AliHFEextraCuts::kFoundOverFindable);
  56   hfecuts->SetCutITSpixel(itshitpixel);
  57   hfecuts->SetCheckITSLayerStatus(kFALSE);
  58 //  hfecuts->SetEtaRange(etami,etama);
  59   if(TRDtl>0) hfecuts->SetMinNTrackletsTRD(TRDtl, kTRUE);   // number of trd tracklets
  60   if(isAOD) hfecuts->SetAODFilterBit(4);
  61
  62   //if((iPixelAny==AliHFEextraCuts::kAny) || (iPixelAny==AliHFEextraCuts::kSecond))
  63   //hfecuts->SetProductionVertex(0,7,0,7);
  64
  65   hfecuts->SetMaxImpactParam(DCAxy,DCAz);
  66   hfecuts->SetUseMixedVertex(kTRUE);
  67   hfecuts->SetVertexRange(10.);
  68   // New pPb cuts (February 2013)
  69   hfecuts->SetUseCorrelationVertex();
  70   hfecuts->SetSPDVtxResolutionCut();
  71   hfecuts->SetpApileupCut();
  72
  73   // TOF settings:
  74   Int_t usetof=0;
  75   Bool_t kTOFmis=kFALSE;
  76   if (TOFs>0.){
  77     usetof = 1;
  78     printf("CONFIGURATION FILE: TOF is used \n");
  79     hfecuts->SetTOFPIDStep(kTRUE);
  80     printf("CONFIGURATION FILE: TOF PID step is requested !!!! \n");
  81     if (TOFmis>0){
  82       kTOFmis = kTRUE;
  83       printf("CONFIGURATION FILE: TOF mismatch rejection is set ON \n");
  84     }
  85   }
  86
  87   //***************************************//
  88   //        Setting up the task            //
  89   //***************************************//
  90
  91   AliAnalysisTaskHFE *task = new AliAnalysisTaskHFE(Form("HFEtask%s",appendix.Data()));
  92   printf("task %p\n", task);
  93   task->SetpPbAnalysis();
  94   if(!isAOD) task->SetRemoveFirstEventInChunk();
  95   task->SetRemovePileUp(kFALSE);
  96   task->SetHFECuts(hfecuts);
  97   task->GetPIDQAManager()->SetHighResolutionHistos();
  98
  99   // Determine the centrality estimator
 100   task->SetCentralityEstimator("V0A");
 101   if (icent == 2) task->SetCentralityEstimator("V0M");
 102   else if (icent == 3) task->SetCentralityEstimator("CL1");
 103   else if (icent == 4) task->SetCentralityEstimator("ZNA");
 104
 105   Bool_t activateTRDTrigger=kFALSE;
 106   if(TRDtrigger>0) activateTRDTrigger=kTRUE;
 107   task->SetTRDTrigger(activateTRDTrigger,TRDtrigger);
 108
 109   //***************************************//
 110   //          Variable manager             //
 111   //***************************************//
 112   // Define Variables
 113   Double_t ptbinning[41] = {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., 22., 24., 26., 28., 30.};
 114   //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.};
 115   //Double_t ptbinning[53] = {0., 0.1, 0.2, 0.22, 0.24, 0.26, 0.28, 0.3, 0.32, 0.34, 0.36, 0.38, 0.4, 0.42, 0.44, 0.46, 0.48, 0.5, 0.55, 0.6, 0.65, 0.7, 0.75, 0.8, 0.85, 0.9, 0.95, 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.};
 116
 117   //Double_t etabinning[9] = {-0.8, -0.6, -0.4, -0.2, 0., 0.2, 0.4, 0.6, 0.8};
 118   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};
 119   //Double_t etabinning[33] = {-0.8,-0.75,-0.7,-0.65,-0.6,-0.55,-0.5,-0.45,-0.4,-0.35,-0.3,-0.25,-0.2,-0.15,-0.1,-0.05,
 120   //0.,0.05,0.1,0.15,0.2,0.25,0.3,0.35,0.4,0.45,0.5,0.55,0.6,0.65,0.7,0.75,0.8};
 121
 122   Int_t sizept=(sizeof(ptbinning)/sizeof(double))-1;
 123   Int_t sizeeta=(sizeof(etabinning)/sizeof(double))-1;
 124
 125   AliHFEvarManager *vm = task->GetVarManager();
 126   vm->AddVariable("pt", sizept, ptbinning);
 127   vm->AddVariable("eta", sizeeta, -0.8,0.8);
 128   vm->AddVariable("phi",21, -0, 2*TMath::Pi());
 129   vm->AddVariable("charge");
 130   vm->AddVariable("source");
 131   vm->AddVariable("centrality");
 132
 133   // For the moment, remove the part dedicated to the background subtraction.
 134   // It will be implemented in a different way, reading it from a root file.
 135
 136   //***************************************//
 137   //          Configure the PID            //
 138   //***************************************//
 139
 140   // Define PID
 141   AliHFEpid *pid = task->GetPID();
 142   if(useMC) pid->SetHasMCData(kTRUE);
 143   pid->SetDetectorsForAnalysis(detector);
 144
 145   // Configure TPC PID
 146   // do the identical thing in data and MC
 147   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};
 148   if(tpcdEdxcutlow) memcpy(paramsTPCdEdxcutlow,tpcdEdxcutlow,sizeof(paramsTPCdEdxcutlow));
 149
 150   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};
 151   if(tpcdEdxcuthigh) memcpy(paramsTPCdEdxcuthigh,tpcdEdxcuthigh,sizeof(paramsTPCdEdxcuthigh));
 152
 153   char *cutmodel;
 154   cutmodel="pol0";
 155
 156   for(Int_t a=0;a<11;a++)
 157     {
 158       //   cout << a << " " << paramsTPCdEdxcut[a] << endl;
 159       Double_t tpcparamlow[1]={paramsTPCdEdxcutlow[a]};
 160       Float_t tpcparamhigh=paramsTPCdEdxcuthigh[a];
 161       pid->ConfigureTPCcentralityCut(a,cutmodel,tpcparamlow,tpcparamhigh);
 162     }
 163   pid->ConfigureTPCdefaultCut(cutmodel,paramsTPCdEdxcutlow,paramsTPCdEdxcuthigh[0]);
 164
 165
 166
 167   // Apply eta correction
 168   AliHFEpidTPC *tpcpid = pid->GetDetPID(AliHFEpid::kTPCpid);
 169   // Setting to cut on absolute dEdx (not sigma)
 170   //tpcpid->UsedEdx();
 171   //TF1 *etacorrection = GetEtaCorrection();
 172   //if(etacorrection) tpcpid->SetEtaCorrection(etacorrection);
 173
 174   // Configure TOF PID
 175   if (usetof){
 176     pid->ConfigureTOF(TOFs);
 177     AliHFEpidTOF *tofpid = pid->GetDetPID(AliHFEpid::kTOFpid);
 178     if (kTOFmis){
 179       tofpid->SetRejectTOFmismatch();
 180     }
 181   }
 182
 183   // To make different upper TOF cut to see contamination effect
 184   // The below two lines should be removed after this check
 185   //AliHFEpidTOF *tofpid = pid->GetDetPID(AliHFEpid::kTOFpid);
 186   //if(TOFs<3.) tofpid->SetTOFnSigmaBand(-3,TOFs); //only to check the assymmetric tof cut
 187
 188   if(TRDtl>0){
 189       AliHFEpidTRD *trdpid = pid->GetDetPID(AliHFEpid::kTRDpid);
 190       if(trdpidmethod==2) trdpid->SetTRD2DPID();
 191       trdpid->SetElectronEfficiency(eeff[TRDeff]);   // efficiency
 192       trdpid->SetNTracklets(TRDtl);      // ntracklets threshold
 193       trdpid->SetCutNTracklets(TRDtl, kTRUE);
 194   }
 195
 196
 197   // Load hadron background
 198   if(!useMC){
 199     Bool_t status = kTRUE;
 200     TF1 *hBackground[12];
 201     //status = ReadContaminationFunctions("hadronContamination_pPbTPCTOF_forwardEta.root", hBackground, tpcdEdxcutlow[0]);
 202     //status = ReadContaminationFunctions("hadronContamination_pPbTPCTOF_eta66.root", hBackground, tpcdEdxcutlow[0]);
 203     if(TRDtl==5){
 204         status = ReadContaminationFunctions("hadroncontamination_TRDTPC_pPb_eta06_5tracklets.root", hBackground, tpcdEdxcutlow[0]);
 205         printf("hadroncont5tracklets \n");
 206     }
 207     if(TRDtl==6)
 208     {
 209         status = ReadContaminationFunctions("hadroncontamination_TRDTPC_pPb_eta06_6tracklets.root", hBackground, tpcdEdxcutlow[0]);
 210         printf("hadroncont6tracklets \n");
 211     }
 212     for(Int_t a=0;a<12;a++) {
 213       //printf("back %f \n",hBackground[a]);
 214       if(status) task->SetBackGroundFactorsFunction(hBackground[a],a);
 215       else printf("not all background functions found\n");
 216     }
 217   }
 218
 219   //***************************************//
 220   //          V0 tagged tracks             //
 221   //***************************************//
 222   if(kAnalyseTaggedTracks){
 223     AliHFEcuts *v0trackCuts = new AliHFEcuts("V0trackCuts", "Track Cuts for tagged track Analysis");
 224     v0trackCuts->CreateStandardCuts();
 225     v0trackCuts->SetMinNClustersTPC(TPCcl);
 226     v0trackCuts->SetMinNClustersTPCPID(TPCclPID);
 227     v0trackCuts->SetMinRatioTPCclusters(0.6);
 228     v0trackCuts->SetTPCmodes(AliHFEextraCuts::kFound, AliHFEextraCuts::kFoundOverFindable);
 229     v0trackCuts->SetMinNClustersITS(1);
 230     v0trackCuts->SetCutITSpixel(AliHFEextraCuts::kFirst);
 231     v0trackCuts->SetCheckITSLayerStatus(kFALSE);
 232     v0trackCuts->UnsetVertexRequirement();
 233     //hfecuts->SetSigmaToVertex(10);
 234     if(usetof) v0trackCuts->SetTOFPIDStep(kTRUE);
 235     if(TRDtl>0)v0trackCuts->SetMinNTrackletsTRD(TRDtl,kTRUE); // condition for TRD tracklets
 236     v0trackCuts->SetQAOn();
 237
 238     task->SwitchOnPlugin(AliAnalysisTaskHFE::kTaggedTrackAnalysis);
 239     task->SetTaggedTrackCuts(v0trackCuts);
 240     task->SetCleanTaggedTrack(kTRUE);
 241   }
 242
 243   // QA
 244   printf("task %p\n", task);
 245   task->SetQAOn(AliAnalysisTaskHFE::kPIDqa);
 246   task->SetQAOn(AliAnalysisTaskHFE::kMCqa);
 247   task->SwitchOnPlugin(AliAnalysisTaskHFE::kDEstep);
 248
 249   printf("*************************************\n");
 250   printf("Configuring standard Task:\n");
 251   task->PrintStatus();
 252   pid->PrintStatus();
 253   printf("*************************************\n");
 254   return task;
 255 }