PWGGA/GammaConv/macros/AddTask_GammaConvFlow_PbPb2.C

   1 /////////////////////////////////////////////////////////////////////////////////////////////
   2 //
   3 // AddTaskGammaConvFlow_PbPb2 macro
   4 // Author: Andrea Dubla, Redmer A. Bertens, Friederike Bock
   5 //         Commented where necessary
   6 /////////////////////////////////////////////////////////////////////////////////////////////
   7
   8 class AliAnalysisDataContainer;
   9 class AliFlowTrackCuts;
  10 class AliFlowTrackSimpleCuts;
  11 class AliFlowEventCuts;
  12 class AliFlowEventSimpleCuts;
  13
  14 void AddTask_GammaConvFlow_PbPb2(
  15                                TString uniqueID = "",
  16                                Int_t harmonic=2,
  17                                Int_t trainConfig = 1,  //change different set of cuts
  18                                //Bool_t isMC   = kFALSE, //run MC
  19                                Int_t enableQAMesonTask = 0, //enable QA in AddTask_GammaConvFlow_PbPb2
  20                                Int_t enableQAPhotonTask = 0, // enable additional QA task
  21                                //TString fileNameInputForWeighting = "MCSpectraInput.root", // path to file for weigting input
  22                                Bool_t doWeighting = kFALSE,  //enable Weighting
  23                                TString cutnumberAODBranch = "1000000060084000001500000",
  24                                Bool_t BasicHistoSP = kFALSE,
  25                                Bool_t debug = kFALSE
  26                                ) {
  27
  28     // ================= Load Librariers =================================
  29     gSystem->Load("libCore.so");
  30     gSystem->Load("libTree.so");
  31     gSystem->Load("libGeom.so");
  32     gSystem->Load("libVMC.so");
  33     gSystem->Load("libPhysics.so");
  34     gSystem->Load("libMinuit");
  35     gSystem->Load("libSTEERBase");
  36     gSystem->Load("libESD");
  37     gSystem->Load("libAOD");
  38     gSystem->Load("libANALYSIS");
  39     gSystem->Load("libANALYSISalice");
  40     gSystem->Load("libCDB.so");
  41     gSystem->Load("libSTEER.so");
  42     gSystem->Load("libSTEERBase.so");
  43     gSystem->Load("libTENDER.so");
  44     gSystem->Load("libTENDERSupplies.so");
  45         gSystem->Load("libPWGflowBase.so");
  46         gSystem->Load("libPWGflowTasks.so");
  47     gSystem->Load("libPWGGAGammaConv.so");
  48
  49     Int_t isHeavyIon = 1;
  50
  51     // ================== GetAnalysisManager ===============================
  52     AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
  53     if (!mgr) {
  54         Error(Form("AddTask_GammaConvV1_%i",trainConfig), "No analysis manager found.");
  55         return ;
  56     }
  57
  58     // ================== GetInputEventHandler =============================
  59     AliVEventHandler *inputHandler=mgr->GetInputEventHandler();
  60
  61     //========= Add PID Reponse to ANALYSIS manager ====
  62     if(!(AliPIDResponse*)mgr->GetTask("PIDResponseTask")){
  63         gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskPIDResponse.C");
  64         AddTaskPIDResponse(isMC);
  65     }
  66
  67     //=========  Set Cutnumber for V0Reader ================================
  68     TString cutnumberPhoton = "000000084001001500000000";
  69     TString cutnumberEvent = "1000000";
  70     AliAnalysisDataContainer *cinput = mgr->GetCommonInputContainer();
  71     //========= Add V0 Reader to  ANALYSIS manager if not yet existent =====
  72     if( !(AliV0ReaderV1*)mgr->GetTask("V0ReaderV1") ){
  73         AliV0ReaderV1 *fV0ReaderV1 = new AliV0ReaderV1("V0ReaderV1");
  74
  75         fV0ReaderV1->SetUseOwnXYZCalculation(kTRUE);
  76         fV0ReaderV1->SetCreateAODs(kFALSE);// AOD Output
  77         fV0ReaderV1->SetUseAODConversionPhoton(kTRUE);
  78
  79         if (!mgr) {
  80             Error("AddTask_V0ReaderV1", "No analysis manager found.");
  81             return;
  82         }
  83
  84         AliConvEventCuts *fEventCuts=NULL;
  85         if(cutnumberEvent!=""){
  86             fEventCuts= new AliConvEventCuts(cutnumberEvent.Data(),cutnumberEvent.Data());
  87             fEventCuts->SetPreSelectionCutFlag(kTRUE);
  88             if(fEventCuts->InitializeCutsFromCutString(cutnumberEvent.Data())){
  89                 fV0ReaderV1->SetEventCuts(fEventCuts);
  90                 fEventCuts->SetFillCutHistograms("",kTRUE);
  91             }
  92         }
  93
  94         // Set AnalysisCut Number
  95         AliConversionPhotonCuts *fCuts=NULL;
  96         if(cutnumberPhoton!=""){
  97             fCuts= new AliConversionPhotonCuts(cutnumberPhoton.Data(),cutnumberPhoton.Data());
  98             fCuts->SetPreSelectionCutFlag(kTRUE);
  99             fCuts->SetIsHeavyIon(isHeavyIon);
 100             if(fCuts->InitializeCutsFromCutString(cutnumberPhoton.Data())){
 101                 fV0ReaderV1->SetConversionCuts(fCuts);
 102                 fCuts->SetFillCutHistograms("",kTRUE);
 103             }
 104         }
 105
 106         if(inputHandler->IsA()==AliAODInputHandler::Class()){
 107             // AOD mode
 108             fV0ReaderV1->SetDeltaAODBranchName(Form("GammaConv_%s_gamma",cutnumberAODBranch.Data()));
 109         }
 110         fV0ReaderV1->Init();
 111
 112         AliLog::SetGlobalLogLevel(AliLog::kInfo);
 113
 114         //connect input V0Reader
 115         mgr->AddTask(fV0ReaderV1);
 116         mgr->ConnectInput(fV0ReaderV1,0,cinput);
 117     }
 118
 119     //================================================
 120     //========= Add task to the ANALYSIS manager =====
 121     //================================================
 122     AliAnalysisTaskGammaConvFlow *task=NULL;
 123     task= new AliAnalysisTaskGammaConvFlow(Form("GammaConvV1_%i_v%d",trainConfig,harmonic));
 124     task->SetIsHeavyIon(isHeavyIon);
 125     //task->SetIsMC(isMC);
 126
 127     cutsRP = new AliFlowTrackCuts(Form("RFPcuts%s",uniqueID));
 128     if(!cutsRP) {
 129         if(debug) cout << " Fatal error: no RP cuts found, could be a library problem! " << endl;
 130         return 0x0;
 131     }
 132     cutsRP = cutsRP->GetStandardVZEROOnlyTrackCuts(); // select vzero tracks
 133     cutsRP->SetVZEROgainEqualizationPerRing(kFALSE);
 134     cutsRP->SetApplyRecentering(kTRUE);
 135
 136     task->SetRPCuts(cutsRP);
 137
 138     AliFlowTrackSimpleCuts *POIfilterVZERO = new AliFlowTrackSimpleCuts();
 139     POIfilterVZERO->SetEtaMin(-0.8);
 140     POIfilterVZERO->SetEtaMax(0.8);
 141     POIfilterVZERO->SetMassMin(263731); POIfilterVZERO->SetMassMax(263733);
 142
 143     if(debug) cout << " === RECEIVED REQUEST FOR FLOW ANALYSIS === " << endl;
 144     AliAnalysisDataContainer *flowEvent = mgr->CreateContainer(Form("FlowContainer_%s",uniqueID.Data()), AliFlowEventSimple::Class(), AliAnalysisManager::kExchangeContainer);
 145     mgr->ConnectOutput(task, 2, flowEvent);
 146     if(debug) cout << "    --> Created IO containers " << flowEvent << endl;
 147     AddSPmethod(Form("SPVZEROQa_in_%s", uniqueID.Data()), "Qa", harmonic, flowEvent, debug,uniqueID, POIfilterVZERO, trainConfig,BasicHistoSP);
 148     if(debug) cout << "    --> Hanging SP Qa task ... succes!" << endl;
 149     AddSPmethod(Form("SPVZEROQb_in_%s", uniqueID.Data()), "Qb", harmonic, flowEvent, debug,uniqueID, POIfilterVZERO, trainConfig,BasicHistoSP);
 150     if(debug) cout << "    --> Hanging SP Qb task ... succes!"<< endl;
 151
 152
 153     // Cut Numbers to use in Analysis
 154     Int_t numberOfCuts = 1;
 155
 156     TString *eventCutArray = new TString[numberOfCuts];
 157     TString *photonCutArray = new TString[numberOfCuts];
 158     TString *mesonCutArray = new TString[numberOfCuts];
 159
 160     if (trainConfig == 1){
 161         eventCutArray[ 0] = "6010001"; photonCutArray[ 0] = "042000092970023220000000"; mesonCutArray[ 0] = "01522045009000";
 162     } else if (trainConfig == 2) {
 163         eventCutArray[ 0] = "6120001"; photonCutArray[ 0] = "042000092970023220000000"; mesonCutArray[ 0] = "01522045009000";
 164     } else if (trainConfig == 3) {
 165         eventCutArray[ 0] = "5010001"; photonCutArray[ 0] = "042000092970023220000000"; mesonCutArray[ 0] = "01522045009000";
 166     } else if (trainConfig == 4) {
 167         eventCutArray[ 0] = "5020001"; photonCutArray[ 0] = "042000092970023220000000"; mesonCutArray[ 0] = "01522045009000";
 168     } else if (trainConfig == 5) {
 169         eventCutArray[ 0] = "5120001"; photonCutArray[ 0] = "042000092970023220000000"; mesonCutArray[ 0] = "01522045009000";
 170     } else if (trainConfig == 6) {
 171         eventCutArray[ 0] = "5240001"; photonCutArray[ 0] = "042000092970023220000000"; mesonCutArray[ 0] = "01522045009000";
 172     } else if (trainConfig == 7) {
 173         eventCutArray[ 0] = "5460001"; photonCutArray[ 0] = "042000092970023220000000"; mesonCutArray[ 0] = "01522065009000";
 174     } else if (trainConfig == 8) {
 175         eventCutArray[ 0] = "5480001"; photonCutArray[ 0] = "042000092970023220000000"; mesonCutArray[ 0] = "01522065009000";
 176     } else if (trainConfig == 9) {
 177         eventCutArray[ 0] = "5450001"; photonCutArray[ 0] = "042000092970023220000000"; mesonCutArray[ 0] = "01522065009000";
 178     } else if (trainConfig == 10) {
 179         eventCutArray[ 0] = "5560001"; photonCutArray[ 0] = "042000092970023220000000"; mesonCutArray[ 0] = "01522065009000";
 180     } else if (trainConfig == 11) {
 181         eventCutArray[ 0] = "5680001"; photonCutArray[ 0] = "042000092970023220000000"; mesonCutArray[ 0] = "01522065009000";
 182     } else if (trainConfig == 12) {
 183         eventCutArray[ 0] = "5670001"; photonCutArray[ 0] = "042000092970023220000000"; mesonCutArray[ 0] = "01522065009000";
 184     } else if (trainConfig == 13) {
 185         eventCutArray[ 0] = "5780001"; photonCutArray[ 0] = "042000092970023220000000"; mesonCutArray[ 0] = "01522065009000";
 186     } else if (trainConfig == 14) {
 187         eventCutArray[ 0] = "4690001"; photonCutArray[ 0] = "042000092970023220000000"; mesonCutArray[ 0] = "01522065009000";
 188     } else if (trainConfig == 15) {
 189         eventCutArray[ 0] = "5890001"; photonCutArray[ 0] = "042000092970023220000000"; mesonCutArray[ 0] = "01522065009000";
 190     } else  if (trainConfig == 16){
 191         eventCutArray[ 0] = "6010001"; photonCutArray[ 0] = "002000092970028250400000"; mesonCutArray[ 0] = "01525065000000";
 192     } else if (trainConfig == 17) {
 193         eventCutArray[ 0] = "6120001"; photonCutArray[ 0] = "002000092970028250400000"; mesonCutArray[ 0] = "01525065000000";
 194     } else if (trainConfig == 18) {
 195         eventCutArray[ 0] = "5010001"; photonCutArray[ 0] = "002000092970028250400000"; mesonCutArray[ 0] = "01525065000000";
 196     } else if (trainConfig == 19) {
 197         eventCutArray[ 0] = "5020001"; photonCutArray[ 0] = "002000092970028250400000"; mesonCutArray[ 0] = "01525065000000";
 198     } else if (trainConfig == 20) {
 199         eventCutArray[ 0] = "5120001"; photonCutArray[ 0] = "002000092970028250400000"; mesonCutArray[ 0] = "01525065000000";
 200     } else if (trainConfig == 21) {
 201         eventCutArray[ 0] = "5240001"; photonCutArray[ 0] = "002000092970028250400000"; mesonCutArray[ 0] = "01525065000000";
 202     } else if (trainConfig == 22) {
 203         eventCutArray[ 0] = "5460001"; photonCutArray[ 0] = "002000092970028250400000"; mesonCutArray[ 0] = "01525065000000";
 204     } else if (trainConfig == 23) {
 205         eventCutArray[ 0] = "5480001"; photonCutArray[ 0] = "002000092970028250400000"; mesonCutArray[ 0] = "01525065000000";
 206     } else if (trainConfig == 24) {
 207         eventCutArray[ 0] = "5450001"; photonCutArray[ 0] = "002000092970028250400000"; mesonCutArray[ 0] = "01525065000000";
 208     } else if (trainConfig == 25) {
 209         eventCutArray[ 0] = "5560001"; photonCutArray[ 0] = "002000092970028250400000"; mesonCutArray[ 0] = "01525065000000";
 210     } else if (trainConfig == 26) {
 211         eventCutArray[ 0] = "5680001"; photonCutArray[ 0] = "002000092970028250400000"; mesonCutArray[ 0] = "01525065000000";
 212     } else if (trainConfig == 27) {
 213         eventCutArray[ 0] = "5670001"; photonCutArray[ 0] = "002000092970028250400000"; mesonCutArray[ 0] = "01525065000000";
 214     } else if (trainConfig == 28) {
 215         eventCutArray[ 0] = "5780001"; photonCutArray[ 0] = "002000092970028250400000"; mesonCutArray[ 0] = "01525065000000";
 216     } else if (trainConfig == 29) {
 217         eventCutArray[ 0] = "4690001"; photonCutArray[ 0] = "002000092970028250400000"; mesonCutArray[ 0] = "01525065000000";
 218     } else if (trainConfig == 30) {
 219         eventCutArray[ 0] = "5890001"; photonCutArray[ 0] = "002000092970028250400000"; mesonCutArray[ 0] = "01525065000000";
 220     } else if (trainConfig == 31) {
 221         eventCutArray[ 0] = "5080001"; photonCutArray[ 0] = "002000092970028250400000"; mesonCutArray[ 0] = "01525065000000";
 222     } else if (trainConfig == 32) {
 223         eventCutArray[ 0] = "5250001"; photonCutArray[ 0] = "002000092970028250400000"; mesonCutArray[ 0] = "01525065000000";
 224     } else if (trainConfig == 33) {
 225         eventCutArray[ 0] = "5350001"; photonCutArray[ 0] = "002000092970028250400000"; mesonCutArray[ 0] = "01525065000000";
 226     } else if (trainConfig == 34) {
 227         eventCutArray[ 0] = "5450001"; photonCutArray[ 0] = "002000092970028250400000"; mesonCutArray[ 0] = "01525065000000";
 228     } else if (trainConfig == 35) {
 229         eventCutArray[ 0] = "5340001"; photonCutArray[ 0] = "002000092970028250400000"; mesonCutArray[ 0] = "01525065000000";
 230     } else if (trainConfig == 36) {
 231         eventCutArray[ 0] = "5230001"; photonCutArray[ 0] = "002000092970028250400000"; mesonCutArray[ 0] = "01525065000000";
 232     } else {
 233         Error(Form("GammaConvV1_%i",trainConfig), "wrong trainConfig variable no cuts have been specified for the configuration");
 234         return;
 235     }
 236
 237     TList *EventCutList = new TList();
 238     TList *ConvCutList = new TList();
 239     TList *MesonCutList = new TList();
 240
 241     TList *HeaderList = new TList();
 242     TObjString *Header1 = new TObjString("BOX");
 243     HeaderList->Add(Header1);
 244     //    TObjString *Header3 = new TObjString("eta_2");
 245     //    HeaderList->Add(Header3);
 246
 247     EventCutList->SetOwner(kTRUE);
 248     AliConvEventCuts **analysisEventCuts = new AliConvEventCuts*[numberOfCuts];
 249     ConvCutList->SetOwner(kTRUE);
 250     AliConversionPhotonCuts **analysisCuts = new AliConversionPhotonCuts*[numberOfCuts];
 251     MesonCutList->SetOwner(kTRUE);
 252     AliConversionMesonCuts **analysisMesonCuts = new AliConversionMesonCuts*[numberOfCuts];
 253
 254     for(Int_t i = 0; i<numberOfCuts; i++){
 255         analysisEventCuts[i] = new AliConvEventCuts();
 256         analysisEventCuts[i]->InitializeCutsFromCutString(eventCutArray[i].Data());
 257         EventCutList->Add(analysisEventCuts[i]);
 258         analysisEventCuts[i]->SetFillCutHistograms("",kFALSE);
 259
 260         analysisCuts[i] = new AliConversionPhotonCuts();
 261         analysisCuts[i]->InitializeCutsFromCutString(photonCutArray[i].Data());
 262         ConvCutList->Add(analysisCuts[i]);
 263         analysisCuts[i]->SetFillCutHistograms("",kFALSE);
 264
 265         analysisMesonCuts[i] = new AliConversionMesonCuts();
 266         analysisMesonCuts[i]->InitializeCutsFromCutString(mesonCutArray[i].Data());
 267         MesonCutList->Add(analysisMesonCuts[i]);
 268         analysisMesonCuts[i]->SetFillCutHistograms("");
 269
 270         analysisEventCuts[i]->SetAcceptedHeader(HeaderList);
 271     }
 272
 273     task->SetEventCutList(numberOfCuts,EventCutList);
 274     task->SetConversionCutList(numberOfCuts,ConvCutList);
 275 //    task->SetMesonCutList(numberOfCuts,MesonCutList);
 276 //    task->SetMoveParticleAccordingToVertex(kTRUE);
 277     task->SetDoMesonAnalysis(kFALSE);
 278     task->SetDoMesonQA(enableQAMesonTask); //Attention new switch for Pi0 QA
 279     task->SetDoPhotonQA(enableQAPhotonTask);  //Attention new switch small for Photon QA
 280
 281     //connect containers
 282     AliAnalysisDataContainer *coutput =
 283         mgr->CreateContainer(Form("GammaConvV1_%i_v%d",trainConfig,harmonic), TList::Class(),
 284                          AliAnalysisManager::kOutputContainer,Form("GammaConvFlow_%i.root",trainConfig));
 285
 286     mgr->AddTask(task);
 287     mgr->ConnectInput(task,0,cinput);
 288     mgr->ConnectOutput(task,1,coutput);
 289
 290     return;
 291
 292 }
 293
 294 //_____________________________________________________________________________
 295 void AddSPmethod(char *name, char *Qvector, int harmonic, AliAnalysisDataContainer *flowEvent, bool debug, TString uniqueID,  AliFlowTrackSimpleCuts* POIfilter,Int_t trainConfig, bool BasicHistoSP = kTRUE)
 296 {
 297     // add sp task and invm filter tasks
 298     if(debug)  cout << " ******* Switching to SP task ******* " << endl;
 299     TString fileName = Form("GammaConvFlow_%i.root:SP_V0",trainConfig);
 300     if(debug) cout << "    --> fileName " << fileName << endl;
 301     TString myFolder = fileName;
 302     if(debug) cout << "    --> myFolder " << myFolder << endl;
 303     TString myNameSP;
 304     myNameSP = Form("%sSPv%d%s", name, harmonic, Qvector);
 305     if(debug) cout << " myNameSP " << myNameSP << endl;
 306     AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
 307     AliAnalysisDataContainer *flowEventOut = mgr->CreateContainer(Form("Filter_%s",myNameSP.Data()),AliFlowEventSimple::Class(),AliAnalysisManager::kExchangeContainer);
 308     AliAnalysisTaskFilterFE *tskFilter = new AliAnalysisTaskFilterFE(Form("TaskFilter_%s", myNameSP.Data()), NULL, POIfilter);
 309    // tskFilter->SetSubeventEtaRange(minEtaA, maxEtaA, minEtaB, maxEtaB);
 310     tskFilter->SetSubeventEtaRange(-10, -1, 1, 10);
 311     mgr->AddTask(tskFilter);
 312     mgr->ConnectInput(tskFilter, 0, flowEvent);
 313     mgr->ConnectOutput(tskFilter, 1, flowEventOut);
 314     AliAnalysisDataContainer *outSP = mgr->CreateContainer(myNameSP.Data(), TList::Class(), AliAnalysisManager::kOutputContainer, fileName);
 315     AliAnalysisTaskScalarProduct *tskSP = new AliAnalysisTaskScalarProduct(Form("TaskScalarProduct_%s", myNameSP.Data()), kFALSE);
 316     tskSP->SetApplyCorrectionForNUA(kTRUE);
 317     tskSP->SetHarmonic(harmonic);
 318     tskSP->SetTotalQvector(Qvector);
 319     //if (bEP) tskSP->SetBehaveAsEP();
 320     tskSP->SetBookOnlyBasicCCH(BasicHistoSP);
 321     mgr->AddTask(tskSP);
 322     mgr->ConnectInput(tskSP, 0, flowEventOut);
 323     mgr->ConnectOutput(tskSP, 1, outSP);
 324 }
 325 //_____________________________________________________________________________
 326
 327
 328
 329
 330