PWGGA/GammaConv/macros/AddTask_GammaConvCalo_PbPb.C

   1 void AddTask_GammaConvCalo_PbPb(        Int_t trainConfig = 1,  //change different set of cuts
   2                                                                         Bool_t isMC   = kFALSE, //run MC
   3                                                                         Int_t enableQAMesonTask = 0, //enable QA in AliAnalysisTaskGammaConvV1
   4                                                                         Int_t enableQAPhotonTask = 0, // enable additional QA task
   5                                                                         TString fileNameInputForWeighting = "MCSpectraInput.root", // path to file for weigting input
   6                                                                         Int_t headerSelectionInt = 0,  // 1 pi0 header, 2 eta header, 3 both (only for "named" boxes)
   7                                                                         TString cutnumberAODBranch = "1000000060084000001500000",
   8                                                                         TString periodName = "LHC13d2",  //name of the period for added signals and weighting
   9                                     Bool_t doWeighting = kFALSE,  //enable Weighting
  10                                     Bool_t enableExtendedMatching = kFALSE //enable or disable extended matching histograms for conversion electrons <-> cluster
  11                                                                 ) {
  12
  13         // ================= Load Librariers =================================
  14         gSystem->Load("libCore.so");
  15         gSystem->Load("libTree.so");
  16         gSystem->Load("libGeom.so");
  17         gSystem->Load("libVMC.so");
  18         gSystem->Load("libPhysics.so");
  19         gSystem->Load("libMinuit");
  20         gSystem->Load("libSTEERBase");
  21         gSystem->Load("libESD");
  22         gSystem->Load("libAOD");
  23         gSystem->Load("libANALYSIS");
  24         gSystem->Load("libANALYSISalice");
  25         gSystem->Load("libPWGGAGammaConv.so");
  26         gSystem->Load("libCDB.so");
  27         gSystem->Load("libSTEER.so");
  28         gSystem->Load("libSTEERBase.so");
  29         gSystem->Load("libTENDER.so");
  30         gSystem->Load("libTENDERSupplies.so");
  31
  32         Int_t isHeavyIon = 1;
  33
  34         // ================== GetAnalysisManager ===============================
  35         AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
  36         if (!mgr) {
  37                 Error(Form("AddTask_GammaConvV1_%i",trainConfig), "No analysis manager found.");
  38                 return ;
  39         }
  40
  41         // ================== GetInputEventHandler =============================
  42         AliVEventHandler *inputHandler=mgr->GetInputEventHandler();
  43
  44         //========= Add PID Reponse to ANALYSIS manager ====
  45         if(!(AliPIDResponse*)mgr->GetTask("PIDResponseTask")){
  46                 gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskPIDResponse.C");
  47                 AddTaskPIDResponse(isMC);
  48         }
  49
  50         //=========  Set Cutnumber for V0Reader ================================
  51         TString cutnumberPhoton = "000084001001500000000";
  52         TString cutnumberEvent = "1000000";
  53         AliAnalysisDataContainer *cinput = mgr->GetCommonInputContainer();
  54
  55         //========= Add V0 Reader to  ANALYSIS manager if not yet existent =====
  56         if( !(AliV0ReaderV1*)mgr->GetTask("V0ReaderV1") ){
  57                 AliV0ReaderV1 *fV0ReaderV1 = new AliV0ReaderV1("V0ReaderV1");
  58
  59                 fV0ReaderV1->SetUseOwnXYZCalculation(kTRUE);
  60                 fV0ReaderV1->SetCreateAODs(kFALSE);// AOD Output
  61                 fV0ReaderV1->SetUseAODConversionPhoton(kTRUE);
  62
  63                 if (!mgr) {
  64                         Error("AddTask_V0ReaderV1", "No analysis manager found.");
  65                         return;
  66                 }
  67
  68                 AliConvEventCuts *fEventCuts=NULL;
  69                 if(cutnumberEvent!=""){
  70                         fEventCuts= new AliConvEventCuts(cutnumberEvent.Data(),cutnumberEvent.Data());
  71                         fEventCuts->SetPreSelectionCutFlag(kTRUE);
  72                         if(fEventCuts->InitializeCutsFromCutString(cutnumberEvent.Data())){
  73                                 fV0ReaderV1->SetEventCuts(fEventCuts);
  74                                 fEventCuts->SetFillCutHistograms("",kTRUE);
  75                         }
  76                 }
  77
  78
  79                 // Set AnalysisCut Number
  80                 AliConversionPhotonCuts *fCuts=NULL;
  81                 if(cutnumberPhoton!=""){
  82                         fCuts= new AliConversionPhotonCuts(cutnumberPhoton.Data(),cutnumberPhoton.Data());
  83                         fCuts->SetPreSelectionCutFlag(kTRUE);
  84                         fCuts->SetIsHeavyIon(isHeavyIon);
  85                         if(fCuts->InitializeCutsFromCutString(cutnumberPhoton.Data())){
  86                                 fV0ReaderV1->SetConversionCuts(fCuts);
  87                                 fCuts->SetFillCutHistograms("",kTRUE);
  88                         }
  89                 }
  90
  91                 if(inputHandler->IsA()==AliAODInputHandler::Class()){
  92                 // AOD mode
  93                         fV0ReaderV1->SetDeltaAODBranchName(Form("GammaConv_%s_gamma",cutnumberAODBranch.Data()));
  94                 }
  95                 fV0ReaderV1->Init();
  96
  97                 AliLog::SetGlobalLogLevel(AliLog::kFatal);
  98
  99                 //connect input V0Reader
 100                 mgr->AddTask(fV0ReaderV1);
 101                 mgr->ConnectInput(fV0ReaderV1,0,cinput);
 102
 103         }
 104
 105         //================================================
 106         //========= Add task to the ANALYSIS manager =====
 107         //================================================
 108         AliAnalysisTaskGammaConvCalo *task=NULL;
 109         task= new AliAnalysisTaskGammaConvCalo(Form("GammaConvCalo_%i",trainConfig));
 110         task->SetIsHeavyIon(isHeavyIon);
 111         task->SetIsMC(isMC);
 112         // Cut Numbers to use in Analysis
 113         Int_t numberOfCuts = 5;
 114
 115         TString *eventCutArray = new TString[numberOfCuts];
 116         TString *photonCutArray = new TString[numberOfCuts];
 117         TString *clusterCutArray = new TString[numberOfCuts];
 118         TString *mesonCutArray = new TString[numberOfCuts];
 119
 120         // meson cuts
 121         // meson type (Dalitz or not), BG scheme, pool depth, rotation degrees, rapidity cut, radius cut, alpha, chi2, shared electrons, reject to close v0, MC smearing, dca, dca, dca
 122
 123         if (trainConfig == 1){ // EMCAL clusters
 124                 eventCutArray[ 0] = "6010001"; photonCutArray[ 0] = "002092970028250400000"; clusterCutArray[0] = "10000040022030000"; mesonCutArray[ 0] = "01525065000000"; // 0-5%
 125                 eventCutArray[ 1] = "6120001"; photonCutArray[ 1] = "002092970028250400000"; clusterCutArray[1] = "10000040022030000"; mesonCutArray[ 1] = "01525065000000"; // 5-10%
 126                 eventCutArray[ 2] = "5010001"; photonCutArray[ 2] = "002092970028250400000"; clusterCutArray[2] = "10000040022030000"; mesonCutArray[ 2] = "01525065000000"; // 0-10%
 127                 eventCutArray[ 3] = "5240001"; photonCutArray[ 3] = "002092970028250400000"; clusterCutArray[3] = "10000040022030000"; mesonCutArray[ 3] = "01525065000000"; // 20-40%
 128                 eventCutArray[ 4] = "5250001"; photonCutArray[ 4] = "002092970028250400000"; clusterCutArray[4] = "10000040022030000"; mesonCutArray[ 4] = "01525065000000"; // 20-50%
 129         } else if (trainConfig == 2){ // PHOS clusters
 130                 eventCutArray[ 0] = "6010001"; photonCutArray[ 0] = "002092970028250400000"; clusterCutArray[0] = "20000030022000000"; mesonCutArray[ 0] = "01525065000000"; // 0-5%
 131                 eventCutArray[ 1] = "6120001"; photonCutArray[ 1] = "002092970028250400000"; clusterCutArray[1] = "20000030022000000"; mesonCutArray[ 1] = "01525065000000"; // 5-10%
 132                 eventCutArray[ 2] = "5010001"; photonCutArray[ 2] = "002092970028250400000"; clusterCutArray[2] = "20000030022000000"; mesonCutArray[ 2] = "01525065000000"; // 0-10%
 133                 eventCutArray[ 3] = "5240001"; photonCutArray[ 3] = "002092970028250400000"; clusterCutArray[3] = "20000030022000000"; mesonCutArray[ 3] = "01525065000000"; // 20-40%
 134                 eventCutArray[ 4] = "5250001"; photonCutArray[ 4] = "002092970028250400000"; clusterCutArray[4] = "20000030022000000"; mesonCutArray[ 4] = "01525065000000"; // 20-50%
 135
 136
 137         } else {
 138                 Error(Form("GammaConvCalo_%i",trainConfig), "wrong trainConfig variable no cuts have been specified for the configuration");
 139                 return;
 140         }
 141
 142         TList *EventCutList = new TList();
 143         TList *ConvCutList = new TList();
 144         TList *ClusterCutList = new TList();
 145         TList *MesonCutList = new TList();
 146
 147         TList *HeaderList = new TList();
 148         if (periodName.CompareTo("LHC13d2")==0){
 149                 TObjString *Header1 = new TObjString("pi0_1");
 150                 HeaderList->Add(Header1);
 151         //    TObjString *Header3 = new TObjString("eta_2");
 152         //    HeaderList->Add(Header3);
 153
 154         } else if (periodName.CompareTo("LHC12a17x_fix")==0){
 155                 TObjString *Header1 = new TObjString("PARAM");
 156                 HeaderList->Add(Header1);
 157         } else if (periodName.CompareTo("LHC14a1a")==0){
 158                 if (headerSelectionInt == 1){
 159                         TObjString *Header1 = new TObjString("pi0_1");
 160                         HeaderList->Add(Header1);
 161                 } else if (headerSelectionInt == 2){
 162                         TObjString *Header1 = new TObjString("eta_2");
 163                         HeaderList->Add(Header1);
 164                 } else {
 165                         TObjString *Header1 = new TObjString("pi0_1");
 166                         HeaderList->Add(Header1);
 167                         TObjString *Header2 = new TObjString("eta_2");
 168                         HeaderList->Add(Header2);
 169                 }
 170         } else if (periodName.CompareTo("LHC14a1b")==0 || periodName.CompareTo("LHC14a1c")==0){
 171                 TObjString *Header1 = new TObjString("BOX");
 172                 HeaderList->Add(Header1);
 173         }
 174
 175         EventCutList->SetOwner(kTRUE);
 176         AliConvEventCuts **analysisEventCuts = new AliConvEventCuts*[numberOfCuts];
 177         ConvCutList->SetOwner(kTRUE);
 178         AliConversionPhotonCuts **analysisCuts = new AliConversionPhotonCuts*[numberOfCuts];
 179         ClusterCutList->SetOwner(kTRUE);
 180         AliCaloPhotonCuts **analysisClusterCuts = new AliCaloPhotonCuts*[numberOfCuts];
 181         MesonCutList->SetOwner(kTRUE);
 182         AliConversionMesonCuts **analysisMesonCuts = new AliConversionMesonCuts*[numberOfCuts];
 183
 184         for(Int_t i = 0; i<numberOfCuts; i++){
 185
 186                 analysisEventCuts[i] = new AliConvEventCuts();
 187                 if ( trainConfig == 1){
 188                         if (periodName.CompareTo("LHC14a1a") ==0 || periodName.CompareTo("LHC14a1b") ==0 || periodName.CompareTo("LHC14a1c") ==0 ){
 189                                 if ( i == 0 && doWeighting)  analysisEventCuts[i]->SetUseReweightingWithHistogramFromFile(kTRUE, kTRUE, kFALSE,fileNameInputForWeighting, Form("Pi0_Hijing_%s_PbPb_2760GeV_0005TPC",periodName.Data()), Form("Eta_Hijing_%s_PbPb_2760GeV_0005TPC",periodName.Data()), "","Pi0_Fit_Data_PbPb_2760GeV_0005V0M","Eta_Fit_Data_PbPb_2760GeV_0005V0M");
 190                                 if ( i == 1 && doWeighting)  analysisEventCuts[i]->SetUseReweightingWithHistogramFromFile(kTRUE, kTRUE, kFALSE,fileNameInputForWeighting, Form("Pi0_Hijing_%s_PbPb_2760GeV_0510TPC",periodName.Data()), Form("Eta_Hijing_%s_PbPb_2760GeV_0510TPC",periodName.Data()), "","Pi0_Fit_Data_PbPb_2760GeV_0510V0M","Eta_Fit_Data_PbPb_2760GeV_0510V0M");
 191                                 if ( i == 2 && doWeighting)  analysisEventCuts[i]->SetUseReweightingWithHistogramFromFile(kTRUE, kTRUE, kFALSE,fileNameInputForWeighting, Form("Pi0_Hijing_%s_PbPb_2760GeV_0010TPC",periodName.Data()), Form("Eta_Hijing_%s_PbPb_2760GeV_0010TPC",periodName.Data()), "","Pi0_Fit_Data_PbPb_2760GeV_0010V0M","Eta_Fit_Data_PbPb_2760GeV_0010V0M");
 192                                 if ( i == 3 && doWeighting)  analysisEventCuts[i]->SetUseReweightingWithHistogramFromFile(kTRUE, kTRUE, kFALSE,fileNameInputForWeighting, Form("Pi0_Hijing_%s_PbPb_2760GeV_2040TPC",periodName.Data()), Form("Eta_Hijing_%s_PbPb_2760GeV_2040TPC",periodName.Data()), "","Pi0_Fit_Data_PbPb_2760GeV_2040V0M","Eta_Fit_Data_PbPb_2760GeV_2040V0M");
 193                                 if ( i == 4 && doWeighting)  analysisEventCuts[i]->SetUseReweightingWithHistogramFromFile(kTRUE, kTRUE, kFALSE,fileNameInputForWeighting, Form("Pi0_Hijing_%s_PbPb_2760GeV_2050TPC",periodName.Data()), Form("Eta_Hijing_%s_PbPb_2760GeV_2050TPC",periodName.Data()), "","Pi0_Fit_Data_PbPb_2760GeV_2050V0M","Eta_Fit_Data_PbPb_2760GeV_2050V0M");
 194                         }
 195                 }
 196                 analysisEventCuts[i]->InitializeCutsFromCutString(eventCutArray[i].Data());
 197                 if (periodName.CompareTo("LHC14a1b") ==0 || periodName.CompareTo("LHC14a1c") ==0 ){
 198                         if (headerSelectionInt == 1) analysisEventCuts[i]->SetAddedSignalPDGCode(111);
 199                         if (headerSelectionInt == 2) analysisEventCuts[i]->SetAddedSignalPDGCode(221);
 200                 }
 201                 EventCutList->Add(analysisEventCuts[i]);
 202                 analysisEventCuts[i]->SetFillCutHistograms("",kFALSE);
 203
 204                 analysisCuts[i] = new AliConversionPhotonCuts();
 205                 analysisCuts[i]->InitializeCutsFromCutString(photonCutArray[i].Data());
 206                 ConvCutList->Add(analysisCuts[i]);
 207                 analysisCuts[i]->SetFillCutHistograms("",kFALSE);
 208
 209                 analysisClusterCuts[i] = new AliCaloPhotonCuts();
 210                 analysisClusterCuts[i]->InitializeCutsFromCutString(clusterCutArray[i].Data());
 211                 ClusterCutList->Add(analysisClusterCuts[i]);
 212         analysisClusterCuts[i]->SetExtendedMatching(enableExtendedMatching);
 213                 analysisClusterCuts[i]->SetFillCutHistograms("");
 214
 215                 analysisMesonCuts[i] = new AliConversionMesonCuts();
 216                 analysisMesonCuts[i]->InitializeCutsFromCutString(mesonCutArray[i].Data());
 217                 MesonCutList->Add(analysisMesonCuts[i]);
 218                 analysisMesonCuts[i]->SetFillCutHistograms("");
 219                 analysisEventCuts[i]->SetAcceptedHeader(HeaderList);
 220
 221         }
 222
 223         task->SetEventCutList(numberOfCuts,EventCutList);
 224         task->SetConversionCutList(numberOfCuts,ConvCutList);
 225         task->SetCaloCutList(numberOfCuts,ClusterCutList);
 226         task->SetMesonCutList(numberOfCuts,MesonCutList);
 227         task->SetMoveParticleAccordingToVertex(kTRUE);
 228         task->SetDoMesonAnalysis(kTRUE);
 229         task->SetDoMesonQA(enableQAMesonTask); //Attention new switch for Pi0 QA
 230         task->SetDoPhotonQA(enableQAPhotonTask);  //Attention new switch small for Photon QA
 231         task->SetDoClusterQA(1);  //Attention new switch small for Cluster QA
 232
 233         //connect containers
 234         AliAnalysisDataContainer *coutput =
 235                 mgr->CreateContainer(Form("GammaConvCalo_%i",trainConfig), TList::Class(),
 236                                                         AliAnalysisManager::kOutputContainer,Form("GammaConvCalo_%i.root",trainConfig));
 237
 238         mgr->AddTask(task);
 239         mgr->ConnectInput(task,0,cinput);
 240         mgr->ConnectOutput(task,1,coutput);
 241
 242         return;
 243
 244 }