From 4873fa4ad43d89d43b31da6961660aa17ba669a6 Mon Sep 17 00:00:00 2001 From: pulvir Date: Tue, 8 Mar 2011 13:22:36 +0000 Subject: [PATCH] First version of tasks for central train --- .../macros/train/AddRsnAnalysisTask.C | 66 +++++ .../macros/train/AddRsnAnalysisTaskEffPhi.C | 210 ++++++++++++++++ .../macros/train/AddRsnEventComputations.C | 99 ++++++++ .../macros/train/AddTenderSupplies.C | 54 ++++ PWG2/RESONANCES/macros/train/AnalysisSetup.C | 135 ++++++++++ .../RESONANCES/macros/train/CPhiCutsAndAxes.C | 209 ++++++++++++++++ PWG2/RESONANCES/macros/train/RsnConfigPhi.C | 231 ++++++++++++++++++ 7 files changed, 1004 insertions(+) create mode 100644 PWG2/RESONANCES/macros/train/AddRsnAnalysisTask.C create mode 100644 PWG2/RESONANCES/macros/train/AddRsnAnalysisTaskEffPhi.C create mode 100644 PWG2/RESONANCES/macros/train/AddRsnEventComputations.C create mode 100644 PWG2/RESONANCES/macros/train/AddTenderSupplies.C create mode 100644 PWG2/RESONANCES/macros/train/AnalysisSetup.C create mode 100644 PWG2/RESONANCES/macros/train/CPhiCutsAndAxes.C create mode 100644 PWG2/RESONANCES/macros/train/RsnConfigPhi.C diff --git a/PWG2/RESONANCES/macros/train/AddRsnAnalysisTask.C b/PWG2/RESONANCES/macros/train/AddRsnAnalysisTask.C new file mode 100644 index 00000000000..2af5a9c9905 --- /dev/null +++ b/PWG2/RESONANCES/macros/train/AddRsnAnalysisTask.C @@ -0,0 +1,66 @@ +// +// This macro serves to add the RSN analysis task to the steering macro. +// +// Inputs: +// - dataLabel = a string with informations about the type of data +// which could be needed to be ported to the config macro +// to set up some cuts +// - configMacro = macro which configures the analysis; it has *ALWAYS* +// defined inside a function named 'RsnConfigTask()', +// whatever the name of the macro itself, whose first two +// arguments must have to be the task and the 'dataLabel' argument. +// +Bool_t AddRsnAnalysisTask +( + Bool_t isMC, + Bool_t isMix, + const char *options, + const char *path = "$(ALICE_ROOT)/PWG2/RESONANCES/macros/train", + const char *taskName = "RSNtask" +) +{ + // retrieve analysis manager + AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager(); + + // create the task and connect with physics selection + AliRsnAnalysisTask *task = new AliRsnAnalysisTask(taskName); + task->SetZeroEventPercentWarning(100.0); + task->SelectCollisionCandidates(); + task->SetMixing(isMix); + ::Info("AddRsnAnalysisTask.C", "Mixing: %s", (task->IsMixing() ? "YES" : "NO")); + ::Info("AddRsnAnalysisTask.C", "MC: %s", (isMC ? "YES" : "NO")); + + // add the task to manager + mgr->AddTask(task); + + // add the event computations with the options to eventually select centrality cut + gROOT->LoadMacro(Form("%s/AddRsnEventComputations.C", path)); + AddRsnEventComputations(isMC, options); + + // add all configs for phi + gROOT->LoadMacro(Form("%s/RsnConfigPhi.C", path)); + RsnConfigPhi(isMC, "tpcpid_tofpid", path, taskName); + //RsnConfigPhi(isMC, "itspid_tpcpid_tofpid", path, taskName); + + // in case of MC, add efficiency tasks + if (isMC) { + gROOT->LoadMacro(Form("%s/AddRsnAnalysisTaskEffPhi.C", path)); + AddRsnAnalysisTaskEffPhi("tpcpid_tofpid"); + //AddRsnAnalysisTaskEffPhi("itspid_tpcpid_tofpid"); + } + + // connect input container according to source choice + mgr->ConnectInput(task, 0, mgr->GetCommonInputContainer()); + + // create paths for the output in the common file + Char_t commonPath[500]; + sprintf(commonPath, "%s", AliAnalysisManager::GetCommonFileName()); + + // create containers for output + AliAnalysisDataContainer *outputInfo = mgr->CreateContainer("RsnInfo", TList::Class(), AliAnalysisManager::kOutputContainer, commonPath); + AliAnalysisDataContainer *outputHist = mgr->CreateContainer("RsnHist", TList::Class(), AliAnalysisManager::kOutputContainer, commonPath); + mgr->ConnectOutput(task, 1, outputInfo); + mgr->ConnectOutput(task, 2, outputHist); + + return kTRUE; +} diff --git a/PWG2/RESONANCES/macros/train/AddRsnAnalysisTaskEffPhi.C b/PWG2/RESONANCES/macros/train/AddRsnAnalysisTaskEffPhi.C new file mode 100644 index 00000000000..5598a7e5b36 --- /dev/null +++ b/PWG2/RESONANCES/macros/train/AddRsnAnalysisTaskEffPhi.C @@ -0,0 +1,210 @@ +// +// This macro add an analysis task for computing efficiency. +// It will have as output an AliCFContainer with several steps: +// +// 0) all resonances in MC which decay in the pair specified +// 1) subset of (0) whose daughters are in acceptance +// 2) subset of (1) whose daughters satisfy quality track cuts (covariance, chi square && nTPCclusters) +// 3) subset of (2) whose daughters satisfy primary track cuts (nsigma to vertex, no kink daughters) +// 4) subset of (3) whose daughters satisty the BB TPC compatibility cut at 3 sigma +// +Bool_t AddRsnAnalysisTaskEffPhi +( + const char *options, + const char *path = "$(ALICE_ROOT)/PWG2/RESONANCES/macros/train" +) +{ + // ================================================================================================================== + // == OPTIONS ======================================================================================================= + // ================================================================================================================== + + TString opt(options); + opt.ToUpper(); + opt.ReplaceAll(" ", ""); + + Bool_t addITS = opt.Contains("ITS"); + Bool_t addTPC = opt.Contains("TPC"); + Bool_t useITS = opt.Contains("ITSPID"); + Bool_t useTPC = opt.Contains("TPCPID"); + Bool_t useTOF = opt.Contains("TOFPID"); + + // correct options when needed + if (!addITS) useITS = kFALSE; + if (!addTPC) useTPC = useTOF = kFALSE; + + // ================================================================================================================== + // == DEFINITIONS =================================================================================================== + // ================================================================================================================== + + // We put here the definitions of all objects which are needed in the following, in order to have then + // a more readable code in the points where these objects are added to the analysis manager. + + // pair definition: + // phi --> K+ K- + AliRsnPairDef *pairPhi = new AliRsnPairDef(AliPID::kKaon, '+', AliPID::kKaon, '-', 333, 1.019455); + + // for cuts and axes, load the support macro + gROOT->LoadMacro("$(ALICE_ROOT)/PWG2/RESONANCES/macros/train/PhiCutsAndAxes.C"); + + AliRsnCutTrackQuality *cutQualityITS = TrackQualityITS(); + AliRsnCutTrackQuality *cutQualityTPC = TrackQualityTPC(); + AliRsnCutPIDITS *cutPIDITSkaon = PIDITS(isMC, AliRsnDaughter::kKaon, 0.0, 1E20, 3.0, "cutPIDITSkaon"); + AliRsnCutPIDTPC *cutPIDTPCkaonLow = PIDTPC(isMC, AliRsnDaughter::kKaon, 0.0, 0.350, 5.0, "cutPIDTPCkaonLow"); + AliRsnCutPIDTPC *cutPIDTPCkaonHigh = PIDTPC(isMC, AliRsnDaughter::kKaon, 0.350, 1E20, 3.0, "cutPIDTPCkaonHigh"); + AliRsnCutPIDTOF *cutPIDTOFkaon = PIDTOF(AliRsnDaughter::kKaon, 3.0, !useTPC, "cutPIDTOFkaon"); + + AliRsnValue *axisIM = AxisIM(); + AliRsnValue *axisPt = AxisPt(); + AliRsnValue *axisRes = AxisRes(); + AliRsnValue *axisMultSPD = AxisMultSPD(); + AliRsnValue *axisMultMC = AxisMultMC(); + + // ================================================================================================================== + // == PRELIMINARY OPERATIONS ======================================================================================== + // ================================================================================================================== + + // retrieve analysis manager + AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager(); + + // create task + AliRsnAnalysisEffSE *task = new AliRsnAnalysisEffSE(Form("RsnEff_%s", options)); + task->SelectCollisionCandidates(); + + // add pair definition, to choose the checked resonance + task->AddPairDef(pairPhi); + + // add the output histogram axis + task->AddAxis(axisIM); + task->AddAxis(axisRes); + task->AddAxis(axisPt); + task->AddAxis(axisY); + task->AddAxis(axisMultSPD); + task->AddAxis(axisMultMC); + + // ================================================================================================================== + // == EVENT CUTS ==================================================================================================== + // ================================================================================================================== + + gROOT->LoadMacro(Form("%s/AddEventStuff.C", path)); + AddEventStuff(task->GetName(), options); + + // ================================================================================================================== + // == STEPS ========================================================================================================= + // ================================================================================================================== + + Char_t qualityITS[255], qualityTPC[255]; + Char_t pidITS[255], pidTPC[255], pidTOF[255]; + Char_t schemeITS[255], schemeTPC[255], scheme[255]; + + sprintf(qualityITS, "%s" , conf.cutQualityITS->GetName()); + sprintf(qualityTPC, "%s" , conf.cutQualityTPC->GetName()); + sprintf(pidITS , "%s" , conf.cutPIDITSkaon->GetName()); + sprintf(pidTPC , "(%s|%s)", conf.cutPIDTPCkaonHigh->GetName(), conf.cutPIDTPCkaonLow->GetName()); + sprintf(pidTOF , "%s" , conf.cutPIDTOFkaon->GetName()); + sprintf(schemeITS , ""); + sprintf(schemeTPC , ""); + sprintf(scheme , ""); + + // + // *** STEP 0 - All resonances which decay in the specified pair + // + // This step does not need any kind of definition, since + // its requirement is automatically checked during execution, + // but to avoid segfaults, it is needed to initialize a cut manager. + // + AliRsnCutManager *mgr_step0 = new AliRsnCutManager("mc_step0", ""); + + // + // *** STEP 1 - Track quality + // + // All resonances whose daughters were reconstructed + // and pass quality track cuts will enter this step + // + AliRsnCutManager *mgr_step1 = new AliRsnCutManager("rec_step1", ""); + AliRsnCutSet *set_step1 = mgr_step1->GetCommonDaughterCuts(); + + if (addTPC && addITS) { + set_step1->AddCut(conf.cutQualityTPC); + set_step1->AddCut(conf.cutQualityITS); + set_step1->SetCutScheme(Form("%s|%s", qualityTPC, qualityITS)); + } else if (addTPC) { + set_step1->AddCut(conf.cutQualityTPC); + set_step1->SetCutScheme(qualityTPC); + } else if (addITS) { + set_step1->AddCut(conf.cutQualityITS); + set_step1->SetCutScheme(qualityITS); + } else { + ::Error("Need to ad at least one between ITS and TPC tracks"); + return kFALSE; + } + ::Info("AddRsnAnalysisTaskEffPhi", "Cut scheme for step #1: %s", set_step1->GetCutScheme().Data()); + + // + // *** STEP 2 - PID + // + // Add all TPC cuts, according to options + // + AliRsnCutManager *mgr_step2 = new AliRsnCutManager("esd_step2", ""); + AliRsnCutSet *set_step2 = mgr_step2->GetCommonDaughterCuts(); + + if (addITS && useITS) { + sprintf(schemeITS, "%s & %s", qualityITS, pidITS); + set_step2->AddCut(conf.cutPIDITSkaon); + } + if (addTPC) { + if (useTPC && useTOF) { + set_step2->AddCut(conf.cutPIDTPCkaonLow); + set_step2->AddCut(conf.cutPIDTPCkaonHigh); + set_step2->AddCut(conf.cutPIDTOFkaon); + sprintf(schemeTPC, "%s & %s", pidTPC, pidTOF); + } else if (useTPC) { + set_step2->AddCut(conf.cutPIDTPCkaonLow); + set_step2->AddCut(conf.cutPIDTPCkaonHigh); + sprintf(schemeTPC, "%s & %s", pidTPC); + } else if (useTOF) { + set_step2->AddCut(conf.cutPIDTOFkaon); + sprintf(schemeTPC, "%s & %s", pidTOF); + } + } + + // final scheme depends on what of the above were added + // in case both ITS-SA and TPC tracks are added, the scheme + // is the OR of the cuts for the first and those for the second + // category of tracks + if (strlen(schemeITS) > 0 && strlen(schemeTPC) > 0) { + sprintf(scheme, "(%s) | (%s)", schemeITS, schemeTPC); + } else if (strlen(schemeITS) > 0) { + sprintf(scheme, "%s", schemeITS); + } else if (strlen(schemeTPC) > 0) { + sprintf(scheme, "%s", schemeTPC); + } else { + ::Error("Scheme is empty!"); + return kFALSE; + } + // check scheme + set_step2->SetCutScheme(scheme); + ::Info("AddRsnAnalysisTaskEffPhi", "Cut scheme for step #1: %s", set_step2->GetCutScheme().Data()); + + // add all steps to the task: + // - first step computed on MC + // - all other steps computed on reconstruction + task->AddStepMC(mgr_step0); + task->AddStepESD(mgr_step1); + task->AddStepESD(mgr_step2); + + // add the task to the manager and connect to input + mgr->AddTask(task); + mgr->ConnectInput(task, 0, mgr->GetCommonInputContainer()); + + // create paths for the output in the common file + TString infoname(task->GetName()); + TString histname(task->GetName()); + infoname.Append("Info"); + histname.Append("Hist"); + AliAnalysisDataContainer *outputInfo = mgr->CreateContainer(infoname.Data(), TList::Class(), AliAnalysisManager::kOutputContainer, mgr->GetCommonFileName()); + AliAnalysisDataContainer *outputHist = mgr->CreateContainer(histname.Data(), TList::Class(), AliAnalysisManager::kOutputContainer, mgr->GetCommonFileName()); + mgr->ConnectOutput(task, 1, outputInfo); + mgr->ConnectOutput(task, 2, outputHist); + + return kTRUE; +} diff --git a/PWG2/RESONANCES/macros/train/AddRsnEventComputations.C b/PWG2/RESONANCES/macros/train/AddRsnEventComputations.C new file mode 100644 index 00000000000..d69645516a6 --- /dev/null +++ b/PWG2/RESONANCES/macros/train/AddRsnEventComputations.C @@ -0,0 +1,99 @@ +// +// Configuration script for monitor task with 2010 runs +// +// It contains a class definition where the cuts for each object +// are defined separately, functions are initialized and so on. +// This is used in the main function (named after the file name), +// which is called by the 'AddTask' function. +// + +#if !defined(__CINT__) || defined(__MAKECINT__) + + #include "TString.h" + #include "AliAnalysisManager.h" + #include "AliRsnValue.h" + #include "AliRsnFunction.h" + #include "AliRsnCutValue.h" + #include "AliRsnCutPrimaryVertex.h" + #include "AliRsnAnalysisTask.h" + +#endif + +Bool_t AddRsnEventComputations(Bool_t isMC, const char *options = "", const char *taskName = "RSNtask") +{ + // ================================================================================================================== + // == PRELIMINARY OPERATIONS ======================================================================================== + // ================================================================================================================== + + // retrieve task from manager, using its name + AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager(); + AliRsnAnalysisTask *task = (AliRsnAnalysisTask*)mgr->GetTask(taskName); + if (!task) { + Error("RsnConfigMonitor", "Task not found"); + return kFALSE; + } + + TString opt(options); + opt.ToUpper(); + opt.ReplaceAll(" ", ""); + + Bool_t central = opt.Contains("CENT"); + Bool_t peripheral = opt.Contains("PERI"); + + // ================================================================================================================== + // == EVENT CUTS ==================================================================================================== + // ================================================================================================================== + + // event cuts are added directly to a cutSet in the task + // we create all and then add thos which are needed + + // primary vertex: + // - 2nd argument --> |Vz| range + // - 3rd argument --> minimum required number of contributors + // - 4th argument --> tells if TPC stand-alone vertexes must be accepted + // we switch on the check for pileup + AliRsnCutPrimaryVertex *cutVertex = new AliRsnCutPrimaryVertex("cutVertex", 10.0, 0, kFALSE); + cutVertex->SetCheckPileUp(kTRUE); + + // centrality: + // - 2nd argument --> one of the centrality evaluation methods + // - 3rd, 4th argument --> centrality ranges in percentile (0-10 for central, 60-70 for peripheral) + AliRsnCutValue *cutCentrality = 0x0; + if (central) + cutCentrality = new AliRsnCutValue("cutCentral", AliRsnValue::kEventCentralityV0, 0.0, 10.0); + else if (peripheral) + cutCentrality = new AliRsnCutValue("cutPeripheral", AliRsnValue::kEventCentralityV0, 60.0, 70.0); + + // primary vertex is always used + task->GetEventCuts()->AddCut(cutVertex); + + // set cut scheme as AND of primary vertex and centrality, if initialized + if (cutCentrality) { + task->GetEventCuts()->AddCut(cutCentrality); + task->GetEventCuts()->SetCutScheme(Form("%s & %s", cutVertex->GetName(), cutCentrality->GetName())); + } else { + task->GetEventCuts()->SetCutScheme(cutVertex->GetName()); + } + ::Info("AddEventStuff", "Scheme for event cuts: %s", task->GetEventCuts()->GetCutScheme().Data()); + + // ================================================================================================================== + // == EVENT FUNCTIONS =============================================================================================== + // ================================================================================================================== + + // we want to add an AliRsnFunction to compute multiplicity distribution + // it is needed in order to know how many events we have in each multiplicity bin + + // axes + AliRsnValue *axisEvMultSPD = new AliRsnValue("MultSPD", AliRsnValue::kEventMultSPD, 0.0, 150.0, 1.0); + AliRsnValue *axisEvMultMC = new AliRsnValue("MultMC" , AliRsnValue::kEventMultMC , 0.0, 150.0, 1.0); + + // create function and add axis + AliRsnFunction *fcnEv = new AliRsnFunction; + if (!fcnEv->AddAxis(axisEvMultSPD)) return kFALSE; + if (isMC && !fcnEv->AddAxis(axisEvMultMC)) return kFALSE; + + // add functions to pairs + task->GetInfo()->AddEventFunction(fcnEv); + + return kTRUE; +} diff --git a/PWG2/RESONANCES/macros/train/AddTenderSupplies.C b/PWG2/RESONANCES/macros/train/AddTenderSupplies.C new file mode 100644 index 00000000000..733c805d2fb --- /dev/null +++ b/PWG2/RESONANCES/macros/train/AddTenderSupplies.C @@ -0,0 +1,54 @@ +AliAnalysisTask *AddTenderSupplies +( + Float_t tofres = 80, + Bool_t corrExpTimes = kFALSE, + Bool_t applyT0 = kFALSE +) +{ + // get the current analysis manager + AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager(); + if (!mgr) { + Error("AddTask_tender_Tender", "No analysis manager found."); + exit(0); + return 0; + } + + // + // === Add tender to the ANALYSIS manager and set default storage ===== + // + AliTender *tender = new AliTender("AnalysisTender"); + tender->SetCheckEventSelection(kFALSE); + //tender->SetDefaultCDBStorage("raw://"); + tender->SetDefaultCDBStorage("alien://folder=/alice/data/2010/OCDB"); + mgr->AddTask(tender); + + // + // === Attach VZERO supply ============================================ + // + AliVZEROTenderSupply *VZEROtender = new AliVZEROTenderSupply("VZEROtender"); + tender->AddSupply(VZEROtender); + + // + // === Attach TPC supply ============================================== + // + AliTPCTenderSupply *TPCtender = new AliTPCTenderSupply("TPCtender"); + tender->AddSupply(TPCtender); + + // + // === Attach TOF supply ============================================== + // + AliTOFTenderSupply *TOFtender = new AliTOFTenderSupply("TOFtender"); + TOFtender->SetTOFres(tofres); + TOFtender->SetApplyT0(applyT0); + TOFtender->SetCorrectExpTimes(corrExpTimes); + tender->AddSupply(TOFtender); + + // + // === Define output containers, please use 'username'_'somename' ===== + // + AliAnalysisDataContainer *coutput1 = mgr->CreateContainer("tender_event", AliESDEvent::Class(), AliAnalysisManager::kExchangeContainer, "default_tender"); + mgr->ConnectInput(tender, 0, mgr->GetCommonInputContainer()); + mgr->ConnectOutput(tender, 1, coutput1); + + return tender; +} diff --git a/PWG2/RESONANCES/macros/train/AnalysisSetup.C b/PWG2/RESONANCES/macros/train/AnalysisSetup.C new file mode 100644 index 00000000000..e6543e1075f --- /dev/null +++ b/PWG2/RESONANCES/macros/train/AnalysisSetup.C @@ -0,0 +1,135 @@ +// +// This macro sets all the aspects of configuration of an Analysis Train run +// which are always the same for all kinds of analysis (local, PROOF, AliEn) +// +// Inputs: +// +// - taskList = a string containin all the 'add-task' macros to be used +// - options = a set of keywords which drive some configurations +// - outputFileName = name of file produced by train +// - configPath = a path where all required config macros are stored +// +// Notes: +// +// - in case the source is an ESD, and if inputs are a MC production +// the MC input handler is created by default +// +Bool_t AnalysisSetup +( + Bool_t isMC, + const char *options, + const char *outputFileName, + const char *configPath +) +{ + // + // === EXAMINE OPTIONS ========================================================================== + // + + // this is done using the utility 'RsnOptions.C' + // which provides a unique way to interpret them + + TString opt(options); + opt.ToUpper(); + opt.ReplaceAll(" ", ""); + + Bool_t isMix = opt.Contains("MIX"); + Bool_t isESD = opt.Contains("ESD"); + Bool_t isAOD = opt.Contains("AOD"); + Bool_t central = opt.Contains("CEN"); + Bool_t peripheral = opt.Contains("PER"); + Bool_t useTender = opt.Contains("TENDER"); + Bool_t usePhysSel = opt.Contains("PHYS"); + Bool_t noV0 = opt.Contains("NOV0"); + + // + // === LOAD LIBRARIES =========================================================================== + // + + gSystem->Load("libVMC.so"); + gSystem->Load("libTree.so"); + gSystem->Load("libPhysics.so"); + gSystem->Load("libMatrix.so"); + gSystem->Load("libMinuit.so"); + gSystem->Load("libXMLParser.so"); + gSystem->Load("libGui.so"); + gSystem->Load("libSTEERBase.so"); + gSystem->Load("libESD.so"); + gSystem->Load("libAOD.so"); + gSystem->Load("libANALYSIS.so"); + gSystem->Load("libANALYSISalice.so"); + gSystem->Load("libEventMixing.so"); + gSystem->Load("libCORRFW.so"); + + if (useTender) { + ::Info("AnalysisSetup", "Loading tender libraries"); + gSystem->Load("libTENDER.so"); + gSystem->Load("libTENDERSupplies.so"); + } + + if (!AliAnalysisAlien::SetupPar("PWG2resonances.par")) return kFALSE; + + // + // === ANALYSIS MANAGER CONFIGURATION =========================================================== + // + + // create analysis manager + AliAnalysisManager *mgr = new AliAnalysisManager("RsnAnalysisManager"); + mgr->SetCommonFileName(outputFileName); + ::Info("AnalysisSetup", "Common file name: %s", outputFileName); + + // + // === INPUT / OUTPUT HANDLER CONFIGURATION ===================================================== + // + + // create input handler + // since there is an exit point above if the job + // isn't either ESD or AOD, here we don't recheck that + if (isESD) { + ::Info("AnalysisSetup", "Configuring for ESD"); + AliESDInputHandler *esdHandler = new AliESDInputHandler(); + mgr->SetInputEventHandler(esdHandler); + // if possible, create also MC handler + if (isMC) { + ::Info("AnalysisSetup", "Creating MC handler"); + AliMCEventHandler *mcHandler = new AliMCEventHandler(); + mgr->SetMCtruthEventHandler(mcHandler); + } + } else if (isAOD) { + ::Info("AnalysisSetup", "Configuring for AOD"); + AliAODInputHandler *aodHandler = new AliAODInputHandler(); + mgr->SetInputEventHandler(aodHandler); + } else { + ::Error("AnalysisSetup", "Require ESD or AOD"); + return kFALSE; + } + + // + // === CONFIGURE AND INSERT PHYSICS SELECTION & TENDER SUPPLY =================================== + // + + // add event selection for data if running ESD + if (isESD) { + // add tender supply for TOF + if (useTender) { + ::Info("AnalysisSetup", "options '%s' require to add tender", options); + gROOT->LoadMacro(Form("%s/AddTenderSupplies.C", configPath)); + AddTenderSupplies(100.0, kTRUE, kFALSE); + } + + // add event selection for data + // and swtich off VZERO if tender is not used + if (usePhysSel) { + ::Info("AnalysisSetup", "options '%s' require to add physics selection", options); + gROOT->LoadMacro("$(ALICE_ROOT)/ANALYSIS/macros/AddTaskPhysicsSelection.C"); + AliPhysicsSelectionTask* physSelTask = AddTaskPhysicsSelection(isMC); + if (noV0) { + ::Info("AnalysisSetup", "options '%s' require to skip V0 info", options); + physSelTask->GetPhysicsSelection()->SetSkipV0(kTRUE); + } + } + } + + ::Info("AnalysisSetup", "Setup successful"); + return kTRUE; +} diff --git a/PWG2/RESONANCES/macros/train/CPhiCutsAndAxes.C b/PWG2/RESONANCES/macros/train/CPhiCutsAndAxes.C new file mode 100644 index 00000000000..1b469b0a474 --- /dev/null +++ b/PWG2/RESONANCES/macros/train/CPhiCutsAndAxes.C @@ -0,0 +1,209 @@ +// +// Initialization of some example cuts of common use +// + +#if !defined(__CINT__) || defined(__MAKECINT__) + + #include "AliRsnCutPrimaryVertex.h" + #include "AliRsnCutTrackQuality.h" + #include "AliRsnCutPIDTPC.h" + #include "AliRsnCutPIDTOF.h" + #include "AliRsnCutPIDITS.h" + #include "AliRsnCutValue.h" + + #include "AliRsnValue.h" + #include "AliRsnPairDef.h" + +#endif + +//__________________________________________________________________________________________________ +// +// Track quality for ITS standalone: +// this cut is used to select tracks of good quality, irrespective of the PID. +// When adding status flags, the second argument tells if each considered flag +// must be active or not in the track status, since the ITS-SA tracks need that +// some of them are OFF (e.g.: kTPCin) +// +AliRsnCutTrackQuality* TrackQualityITS() +{ + AliRsnCutTrackQuality *cutQualityITS = new AliRsnCutTrackQuality("cutQualityITS"); + cutQualityITS->AddStatusFlag(AliESDtrack::kITSin , kTRUE); + cutQualityITS->AddStatusFlag(AliESDtrack::kTPCin , kFALSE); + cutQualityITS->AddStatusFlag(AliESDtrack::kITSrefit , kTRUE); + cutQualityITS->AddStatusFlag(AliESDtrack::kTPCrefit , kFALSE); + cutQualityITS->AddStatusFlag(AliESDtrack::kITSpureSA, kFALSE); + cutQualityITS->AddStatusFlag(AliESDtrack::kITSpid , kTRUE); + cutQualityITS->SetPtRange(0.15, 1E+20); + cutQualityITS->SetEtaRange(-0.8, 0.8); + cutQualityITS->SetDCARPtFormula("0.0595+0.0182/pt^1.55"); + cutQualityITS->SetDCAZmax(2.0); + cutQualityITS->SetSPDminNClusters(1); + cutQualityITS->SetITSminNClusters(4); + cutQualityITS->SetITSmaxChi2(2.0); + cutQualityITS->SetTPCminNClusters(0); + cutQualityITS->SetTPCmaxChi2(1E+10); + cutQualityITS->SetRejectKinkDaughters(); + + return cutQualityITS; +} + +//__________________________________________________________________________________________________ +// +// Track quality for TPC+ITS: +// works exactly like the one above, but has settings for selecting TPC+ITS tracks +// in this case, the flags required are all necessary, so here the procedure is simpler +// + +AliRsnCutTrackQuality *TrackQualityTPC() +{ + AliRsnCutTrackQuality *cutQualityTPC = new AliRsnCutTrackQuality("cutQualityTPC"); + cutQualityTPC->AddStatusFlag(AliESDtrack::kTPCin , kTRUE); + cutQualityTPC->AddStatusFlag(AliESDtrack::kTPCrefit, kTRUE); + cutQualityTPC->AddStatusFlag(AliESDtrack::kITSrefit, kTRUE); + cutQualityTPC->SetPtRange(0.15, 1E+20); + cutQualityTPC->SetEtaRange(-0.8, 0.8); + cutQualityTPC->SetDCARPtFormula("0.0182+0.0350/pt^1.01"); + cutQualityTPC->SetDCAZmax(2.0); + cutQualityTPC->SetSPDminNClusters(1); + cutQualityTPC->SetITSminNClusters(0); + cutQualityTPC->SetITSmaxChi2(1E+20); + cutQualityTPC->SetTPCminNClusters(70); + cutQualityTPC->SetTPCmaxChi2(4.0); + cutQualityTPC->SetRejectKinkDaughters(); + + return cutQualityTPC; +} + +//__________________________________________________________________________________________________ +// +// The ITS PID is done with a 3sigma band in the dE/dx vs. Bethe-Bloch comparison, +// in the whole momentum range, with a set of parameters which are different in +// data and MonteCarlo, so we need to know on what we are running. +// +AliRsnCutPIDITS *PIDITS +(Bool_t isMC, AliPID::EParticleType pid, Double_t pmin, Double_t pmax, Double_t nsigma, const char *name) +{ + AliRsnCutPIDITS *cutPIDITS = new AliRsnCutPIDITS(name, pid, -nsigma, nsigma); + cutPIDITS->SetMC(isMC); + cutPIDITS->SetMomentumRange(pmin, pmax); + + return cutPIDITS; +} + +//__________________________________________________________________________________________________ +// +// For the TPC PID we define two zones: one below 350 MeV/c in momentum (at the TPC barrel) +// where we apply a 5sigma cut, and another above 350 MeV/c where the cut is tightened +// to the range of 3sigma. +// Even in this case we have different parameters for data and MC, but here we have to set +// them manually (they are passed to the AliESDpid object which is in the cut) +// +AliRsnCutPIDTPC *PIDTPC +(Bool_t isMC, AliPID::EParticleType pid, Double_t pmin, Double_t pmax, Double_t nsigma, const char *name) +{ + AliRsnCutPIDTPC *cutPIDTPC = new AliRsnCutPIDTPC(name, pid, -nsigma, nsigma); + + // BB parameterization depends on data sample (MC, data) + // the momentum range is passed and tracks outside it are rejected + Double_t bbPar[5]; + if (isMC) { + bbPar[0] = 2.15898 / 50.0; + bbPar[1] = 1.75295E1; + bbPar[2] = 3.40030E-9; + bbPar[3] = 1.96178; + bbPar[4] = 3.91720; + } else { + bbPar[0] = 1.41543 / 50.0; + bbPar[1] = 2.63394E1; + bbPar[2] = 5.0411E-11; + bbPar[3] = 2.12543; + bbPar[4] = 4.88663; + } + cutPIDTPC->SetBBParam(bbPar); + cutPIDTPC->SetMomentumRange(pmin, pmax); + cutPIDTPC->SetRejectOutside(kTRUE); + + return cutPIDTPC; +} + +//__________________________________________________________________________________________________ +// +// The TOF PID is simply a 3sigma cut. +// Since it is not doing a comparison with a reference function, +// we don't need to know if we are on data or MC. +// It is important to choose if this cut must reject tracks not matched in TOF. +// Usually, if TPC pid is used, we can accept them, otherwise we must reject. +// +AliRsnCutPIDTOF *PIDTOF +(AliPID::EParticleType pid, Double_t nsigma, Bool_t rejectUnmatched, const char *name) +{ + AliRsnCutPIDTOF *cutPIDTOF = new AliRsnCutPIDTOF(name, pid, -nsigma, nsigma); + cutPIDTOF->SetRejectUnmatched(rejectUnmatched); + + return cutPIDTOF; +} + +//__________________________________________________________________________________________________ +// +// A cut is applied on the rapidity of the pair. +// This is a vaolue which can be computed for filling histograms, +// and all such values can be cut using the generic AliRsnCutValue class. +// Since rapidity needs a mass hypothesis, we need to add a support object +// which is a reference pair definition (we can use any of those defined above) +// +AliRsnCutValue *RapidityRange +(AliRsnPairDef *def, Double_t min, Double_t max, const char *name = "cutY") +{ + AliRsnCutValue *cutRapidity = new AliRsnCutValue(name, AliRsnValue::kPairY, min, max); + cutRapidity->GetValueObj()->SetSupportObject(def); + + return cutRapidity; +} + +//__________________________________________________________________________________________________ +// +// Define all axes +// +AliRsnValue *AxisIM(Double_t min = 0.9, Double_t max = 1.4, Double_t step = 0.001) +{ + AliRsnValue *axis = new AliRsnValue("IM", AliRsnValue::kPairInvMass, min, max, step); + return axis; +} + +AliRsnValue *AxisRes(Double_t min = -0.5, Double_t max = 0.5, Double_t step = 0.001) +{ + AliRsnValue *axis = new AliRsnValue("IM", AliRsnValue::kPairInvMassRes, min, max, step); + return axis; +} + +AliRsnValue *AxisPt(Double_t min = 0.0, Double_t max = 5.0, Double_t step = 0.1) +{ + AliRsnValue *axis = new AliRsnValue("PT", AliRsnValue::kPairPt, min, max, step); + return axis; +} + +AliRsnValue *AxisY(Double_t min = -1.5, Double_t max = 1.5, Double_t step = 0.1) +{ + AliRsnValue *axis = new AliRsnValue("Y", AliRsnValue::kPairPt, min, max, step); + return axis; +} + +AliRsnValue *AxisMultSPD() +{ + Double_t mult[] = { 0., 1., 2., 3., 4., 5., 6., 7., 8., 9., 10., 11., 12., 13., 14., 15., 16., 17., 18., 19., + 20., 21., 22., 23., 24., 25., 30., 35., 40., 50., 60., 70., 80., 90., 100., 120., 140., 160., 180., 200., 500.}; + Int_t nmult = sizeof(mult) / sizeof(mult[0]); + + AliRsnValue *axis = new AliRsnValue("MSPD", AliRsnValue::kEventMultSPD, nmult, mult); + return axis; +} + +AliRsnValue *AxisMultMC() +{ + Double_t mult[] = { 0., 1., 2., 3., 4., 5., 6., 7., 8., 9., 10., 11., 12., 13., 14., 15., 16., 17., 18., 19., + 20., 21., 22., 23., 24., 25., 30., 35., 40., 50., 60., 70., 80., 90., 100., 120., 140., 160., 180., 200., 500.}; + Int_t nmult = sizeof(mult) / sizeof(mult[0]); + + AliRsnValue *axis = new AliRsnValue("MSPD", AliRsnValue::kEventMultMC, nmult, mult); + return axis; +} diff --git a/PWG2/RESONANCES/macros/train/RsnConfigPhi.C b/PWG2/RESONANCES/macros/train/RsnConfigPhi.C new file mode 100644 index 00000000000..2f7d121f0ca --- /dev/null +++ b/PWG2/RESONANCES/macros/train/RsnConfigPhi.C @@ -0,0 +1,231 @@ +// +// *** Configuration script for phi->KK analysis with 2010 runs *** +// +// A configuration script for RSN package needs to define the followings: +// +// (1) decay tree of each resonance to be studied, which is needed to select +// true pairs and to assign the right mass to all candidate daughters +// (2) cuts at all levels: single daughters, tracks, events +// (3) output objects: histograms or trees +// +Bool_t RsnConfigPhi +( + Bool_t isMC, + const char *options, + const char *path = "$(ALICE_ROOT)/PWG2/RESONANCES/macros/train", + const char *taskName = "RSNtask" +) +{ + // ================================================================================================================== + // == PRELIMINARY OPERATIONS ======================================================================================== + // ================================================================================================================== + + // retrieve task from manager, using its name + AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager(); + AliRsnAnalysisTask *task = (AliRsnAnalysisTask*)mgr->GetTask(taskName); + if (!task) { + Error("RsnConfigMonitor", "Task not found"); + return kFALSE; + } + + // ================================================================================================================== + // == OPTIONS ======================================================================================================= + // ================================================================================================================== + + // Instead or getting confused with plenty of arguments in the macro (with default values), + // we use a unique string of options with a set of conventional strings to set up the job: + // -- "MC"/"DATA" --> what kind of sample + // -- "ITS"/"TPC" --> what tracks to use (ITS standalone and/or TPC+ITS) + // -- "xxxPID" --> add the PID cuts for the detector xxx. + // + // In this point, these options are converted into boolean variables. + + TString opt(options); + opt.ToUpper(); + opt.ReplaceAll(" ", ""); + + Bool_t addITS = opt.Contains("ITS"); + Bool_t addTPC = opt.Contains("TPC"); + Bool_t useITS = opt.Contains("ITSPID"); + Bool_t useTPC = opt.Contains("TPCPID"); + Bool_t useTOF = opt.Contains("TOFPID"); + + // correct options when needed + if (!addITS) useITS = kFALSE; + if (!addTPC) useTPC = useTOF = kFALSE; + + // ================================================================================================================== + // == DEFINITIONS =================================================================================================== + // ================================================================================================================== + + // We put here the definitions of all objects which are needed in the following, in order to have then + // a more readable code in the points where these objects are added to the analysis manager. + + // pair definitions --> decay channels: + // in our case, unlike-charged KK pairs for the signal, and like-charged ones for background + AliRsnPairDef *pairDef[3]; + pairDef[0] = new AliRsnPairDef(AliRsnDaughter::kKaon, '+', AliRsnDaughter::kKaon, '-', 333, 1.019455); // unlike + pairDef[1] = new AliRsnPairDef(AliRsnDaughter::kKaon, '+', AliRsnDaughter::kKaon, '+', 333, 1.019455); // like ++ + pairDef[2] = new AliRsnPairDef(AliRsnDaughter::kKaon, '-', AliRsnDaughter::kKaon, '-', 333, 1.019455); // like -- + + // computation objects: + // these are the objects which drive the computations, whatever it is (histos or tree filling) + // and all tracks passed to them will be given the mass according to the reference pair definition + // we create two unlike-sign pair computators, one for all pairs and another for true pairs (useful in MC) + AliRsnPairFunctions *pair[4]; + pair[0] = new AliRsnPairFunctions(Form("%s_kaonP_kaonM_phi", opt.Data()), pairDef[0]); // unlike - true + pair[1] = new AliRsnPairFunctions(Form("%s_kaonP_kaonM_all", opt.Data()), pairDef[0]); // unlike - all + pair[2] = new AliRsnPairFunctions(Form("%s_kaonP_kaonM_all", opt.Data()), pairDef[1]); // like ++ + pair[3] = new AliRsnPairFunctions(Form("%s_kaonP_kaonM_all", opt.Data()), pairDef[2]); // like -- + + // set additional option for true pairs (slot [0]) + pair[0]->SetOnlyTrue(kTRUE); + pair[0]->SetCheckDecay(kTRUE); + + // for cuts and axes, load the support macro + gSystem->AddIncludePath("-I$ALICE_ROOT/include -I$ALICE_ROOT/PWG2/RESONANCES"); + gROOT->LoadMacro(Form("%s/CPhiCutsAndAxes.C++", path)); + + // ================================================================================================================== + // == SINGLE DAUGHTER CUTS ========================================================================================== + // ================================================================================================================== + + // we need to define cuts for track quality and for PID + // which one will be used, depend on the analysis type + + // use functions in the loaded macro to create cuts + AliRsnCutTrackQuality *cutQualityITS = TrackQualityITS(); + AliRsnCutTrackQuality *cutQualityTPC = TrackQualityTPC(); + AliRsnCutPIDITS *cutPIDITSkaon = PIDITS(isMC, AliRsnDaughter::kKaon, 0.0, 1E20, 3.0, "cutPIDITSkaon"); + AliRsnCutPIDTPC *cutPIDTPCkaonLow = PIDTPC(isMC, AliRsnDaughter::kKaon, 0.0, 0.350, 5.0, "cutPIDTPCkaonLow"); + AliRsnCutPIDTPC *cutPIDTPCkaonHigh = PIDTPC(isMC, AliRsnDaughter::kKaon, 0.350, 1E20, 3.0, "cutPIDTPCkaonHigh"); + AliRsnCutPIDTOF *cutPIDTOFkaon = PIDTOF(AliRsnDaughter::kKaon, 3.0, !useTPC, "cutPIDTOFkaon"); + + Char_t qualityITS[255], qualityTPC[255]; + Char_t pidITS[255], pidTPC[255], pidTOF[255]; + Char_t schemeITS[255], schemeTPC[255], scheme[255]; + + sprintf(qualityITS, "%s" , cutQualityITS->GetName()); + sprintf(qualityTPC, "%s" , cutQualityTPC->GetName()); + sprintf(pidITS , "%s" , cutPIDITSkaon->GetName()); + sprintf(pidTPC , "(%s|%s)", cutPIDTPCkaonHigh->GetName(), cutPIDTPCkaonLow->GetName()); + sprintf(pidTOF , "%s" , cutPIDTOFkaon->GetName()); + sprintf(schemeITS , ""); + sprintf(schemeTPC , ""); + + // choose cuts to add depending on used tracks + for (Int_t i = 0; i < 4; i++) { + // assign cuts to common manager for daughters + AliRsnCutSet *cutSet = pair[i]->GetCutManager()->GetCommonDaughterCuts(); + // add cuts and define schemes depending on options + if (addITS) { + // if adding ITS standalone, + // by default we use the quality check + // and if required we add the PID (option "ITSPID") + cutSet->AddCut(cutQualityITS); + sprintf(schemeITS, "%s", qualityITS); + if (useITS) { + sprintf(schemeITS, "%s & %s", qualityITS, pidITS); + cutSet->AddCut(cutPIDITSkaon); + } + } + if (addTPC) { + // if adding TPC+ITS tracks, + // by default we use the quality check + // and then we can use no PID, or both, + // or only one among TPC and TOF + // the TPC PID cut is always an 'OR' of + // the two ones, in order to cover the full momentum range + cutSet->AddCut(cutQualityTPC); + sprintf(schemeTPC, "%s", qualityTPC); + if (useTPC && useTOF) { + cutSet->AddCut(cutPIDTPCkaonLow); + cutSet->AddCut(cutPIDTPCkaonHigh); + cutSet->AddCut(cutPIDTOFkaon); + sprintf(schemeTPC, "%s & %s & %s", qualityTPC, pidTPC, pidTOF); + } else if (useTPC) { + cutSet->AddCut(cutPIDTPCkaonLow); + cutSet->AddCut(cutPIDTPCkaonHigh); + sprintf(schemeTPC, "%s & %s", qualityTPC, pidTPC); + } else if (useTOF) { + cutSet->AddCut(cutPIDTOFkaon); + sprintf(schemeTPC, "%s & %s", qualityTPC, pidTOF); + } + } + + // final scheme depends on what of the above were added + // in case both ITS-SA and TPC tracks are added, the scheme + // is the OR of the cuts for the first and those for the second + // category of tracks + if (strlen(schemeITS) > 0 && strlen(schemeTPC) > 0) { + sprintf(scheme, "(%s) | (%s)", schemeITS, schemeTPC); + } else if (strlen(schemeITS) > 0) { + sprintf(scheme, "%s", schemeITS); + } else if (strlen(schemeTPC) > 0) { + sprintf(scheme, "%s", schemeTPC); + } else { + ::Error("Scheme is empty!"); + return kFALSE; + } + cutSet->SetCutScheme(scheme); + ::Info("RsnConfigPhi", "Scheme for daughter cuts: %s", cutSet->GetCutScheme().Data()); + } + + // ================================================================================================================== + // == PAIR CUTS ===================================================================================================== + // ================================================================================================================== + + // we define here our rapidity window + AliRsnCutValue *cutRapidity = RapidityRange(pairDef[0], -0.5, 0.5); + + // in this case, we add the cut to the specific cut sets of all pairs + // and we must then loop over all pairs, to add cuts to the related sets + for (Int_t ipair = 0; ipair < 4; ipair++) { + pair[ipair]->GetCutManager()->GetMotherCuts()->AddCut(cutRapidity); + pair[ipair]->GetCutManager()->GetMotherCuts()->SetCutScheme(cutRapidity->GetName()); + ::Info("RsnConfigPhi", "Scheme for pair cuts: %s", pair[ipair]->GetCutManager()->GetMotherCuts()->GetCutScheme().Data()); + } + + // ================================================================================================================== + // == FUNCTIONS FOR HISTOGRAMS ====================================================================================== + // ================================================================================================================== + + // we define a histogram with inv mass and other support binning values + // and, when possible, a resolution on inv. mass + + AliRsnValue *axisIM = AxisIM(); + AliRsnValue *axisPt = AxisPt(); + AliRsnValue *axisRes = AxisRes(); + AliRsnValue *axisMultSPD = AxisMultSPD(); + AliRsnValue *axisMultMC = AxisMultMC(); + + // create function for inv. mass and add axes + AliRsnFunction *fcnIM = new AliRsnFunction; + if (!fcnIM->AddAxis(axisIM) ) return kFALSE; + if (!fcnIM->AddAxis(axisPt) ) return kFALSE; + if (!fcnIM->AddAxis(axisMultSPD)) return kFALSE; + fcnIM->UseSparse(); + + // create function for inv. mass and add axes + AliRsnFunction *fcnIMRes = new AliRsnFunction; + if (!fcnIMRes->AddAxis(axisIM) ) return kFALSE; + if (!fcnIMRes->AddAxis(axisPt) ) return kFALSE; + if (!fcnIMRes->AddAxis(axisMultSPD)) return kFALSE; + if (!fcnIMRes->AddAxis(axisRes) ) return kFALSE; + + // add functions to pairs + pair[0]->AddFunction(fcnIMRes); + for (Int_t ipair = 1; ipair < 4; ipair++) pair[ipair]->AddFunction(fcnIM); + + // ================================================================================================================== + // == CONCLUSION ==================================================================================================== + // ================================================================================================================== + + // add all created AliRsnPair objects to the AliRsnAnalysisManager in the task + task->GetAnalysisManager()->Add(pair[1]); + task->GetAnalysisManager()->Add(pair[2]); + task->GetAnalysisManager()->Add(pair[3]); + if (isMC) task->GetAnalysisManager()->Add(pair[0]); + + return kTRUE; +} -- 2.43.0