PWG2/RESONANCES/macros/train/ConfigTaskRsn.C

   1 //
   2 // CONFIGURATION SCRIPT FOR RESONANCES ANALYSIS TRAIN
   3 //
   4 // The configuration script made of two functions:
   5 //   1) RsnConfigPairManager
   6 //   2) RsnConfigTask
   7 //
   8 // Function 'RsnConfigPairManager' defines a common scheme used
   9 // for axis definition, cuts, and whatever is in common among
  10 // all the output histograms.
  11 // Function 'RsnConfigTask' calls the above function for all
  12 // specific resonance pairs one wants to create.
  13 //
  14
  15 //_________________________________________________________________________________________________
  16 AliRsnPairManager* RsnConfigPairManager
  17 (
  18   const char            *pairMgrName,
  19   Int_t                  resonancePDG,
  20   Double_t               resonanceMass,
  21   AliPID::EParticleType  type1,
  22   AliPID::EParticleType  type2,
  23   Bool_t                 addTruePairs
  24 )
  25 //
  26 // -- RsnConfigPairManager --
  27 //
  28 // This function defines:
  29 //  - what pairs one wants
  30 //  - what cuts have to be applied
  31 //  - what axes the output histogram will have
  32 //
  33 // Output:
  34 //  - an AliRsnPairManager containing all pairs related to one analysis
  35 //   (usually, all pairs in the same manager refer to the same resonance
  36 //    and this is assumed in the philosophy of the function itself)
  37 //
  38 // Arguments:
  39 //  - pairMgrName   = name assigned to output object
  40 //  - resonancePDG  = PDG code of resonance in study
  41 //  - resonanceMass = nominal (PDG) mass of resonance in study (required for mass-based axes, like Y and Mt)
  42 //  - type1         = first daughter type
  43 //  - type2         = second daughter type
  44 //  - addTruePairs  = flag to decide if true pairs histogram must be computed
  45 //                    (requires MC info, so it is not feasible with AOD only)
  46 //
  47 // == NOTE 1 ==
  48 // in this version of the configuration macro, the name given to PairManager must contain
  49 // some keywords which define a choice of cuts which will be added to the analysis:
  50 //  - "NOPID": completely no PID analysis (only primary track cuts are applied)
  51 //  - "BB"   : all tracks are used, but the TPC Bethe-Bloch cut is applied (cut value = 0.2)
  52 //  - "PID"  : realistic PID is used
  53 //
  54 // == NOTE 2 ==
  55 // in this version of the configuration macro, the output histograms have the following axes:
  56 //  0 = inv. mass
  57 //  1 = transverse mass
  58 //  2 = rapidity
  59 //
  60 {
  61   // === NAME DEFINITIONS =========================================================================
  62
  63   AliRsnPairManager  *pairMgr  = new AliRsnPairManager(pairMgrName);
  64
  65   // examines the given name to define details about track selection and cuts
  66   TString               str(pairMgrName);
  67   AliRsnPair::EPairType pidType;
  68   Bool_t                useBBCut;
  69   if (str.Contains("NOPID"))
  70   {
  71     pidType  = AliRsnPair::kNoPID;
  72     useBBCut = kFALSE;
  73     Info("RsnConfig", "PID TYPE = No PID        -- BB CUT: not used");
  74   }
  75   else if (str.Contains("BB"))
  76   {
  77     pidType  = AliRsnPair::kNoPID;
  78     useBBCut = kTRUE;
  79     Info("RsnConfig", "PID TYPE = No PID        -- BB CUT: used");
  80   }
  81   else if (str.Contains("PID"))
  82   {
  83     pidType  = AliRsnPair::kRealisticPID;
  84     useBBCut = kFALSE;
  85     Info("RsnConfig", "PID TYPE = Realistic PID -- BB CUT: not used");
  86   }
  87   else
  88   {
  89     Error("RsnConfig", "Unrecognized keywords in the name. Can't continue");
  90     return 0x0;
  91   }
  92
  93   // === PAIR DEFINITIONS =========================================================================
  94
  95   // if particle #1 and #2 are different, two histograms must be built
  96   // for each scheme (signal, true, mixed, like-sign) exchanging both particles and signs
  97   Int_t i, j, nArray = 1;
  98   if (type1 != type2) nArray = 2;
  99
 100   AliRsnPairDef *defUnlike[2] = {0, 0};
 101   AliRsnPairDef *defLikePP[2] = {0, 0};
 102   AliRsnPairDef *defLikeMM[2] = {0, 0};
 103
 104   defUnlike[0] = new AliRsnPairDef(type1, '+', type2, '-', resonancePDG, resonanceMass);
 105   defLikePP[0] = new AliRsnPairDef(type1, '+', type2, '+', resonancePDG, resonanceMass);
 106   defLikeMM[0] = new AliRsnPairDef(type1, '-', type2, '-', resonancePDG, resonanceMass);
 107
 108   defUnlike[1] = new AliRsnPairDef(type2, '+', type1, '-', resonancePDG, resonanceMass);
 109   defLikePP[1] = new AliRsnPairDef(type2, '+', type1, '+', resonancePDG, resonanceMass);
 110   defLikeMM[1] = new AliRsnPairDef(type2, '-', type1, '-', resonancePDG, resonanceMass);
 111
 112   // === PAIR ANALYSIS ENGINES ====================================================================
 113
 114   // define null (dummy) objects and initialize only the ones which are needed,
 115   // depending again on particle types;
 116   // array is organized as follows:
 117   // [0] - true pairs
 118   // [1] - signal
 119   // [2] - like PP
 120   // [3] - like MM
 121   AliRsnPair *pair[2][4];
 122   for (i = 0; i < 2; i++) for (j = 0; j < 4; j++) pair[i][j] = 0x0;
 123
 124   for (i = 0; i < nArray; i++) {
 125     pair[i][0] = new AliRsnPair(pidType, defUnlike[i]);
 126     pair[i][1] = new AliRsnPair(pidType, defUnlike[i]);
 127     pair[i][2] = new AliRsnPair(pidType, defLikePP[i]);
 128     pair[i][3] = new AliRsnPair(pidType, defLikeMM[i]);
 129   }
 130
 131   // === CUTS =====================================================================================
 132
 133   // cuts for tracks:
 134   // -- primary track quality
 135   AliRsnCutESDPrimary *cutESDPrimary = new AliRsnCutESDPrimary("cutESDPrimary");
 136   cutESDPrimary->GetCuts()->SetMaxCovDiagonalElements(2.0, 2.0, 0.5, 0.5, 2.0);
 137   cutESDPrimary->GetCuts()->SetRequireSigmaToVertex(kTRUE);
 138   cutESDPrimary->GetCuts()->SetMaxNsigmaToVertex(3.0);
 139   cutESDPrimary->GetCuts()->SetRequireTPCRefit(kTRUE);
 140   cutESDPrimary->GetCuts()->SetAcceptKinkDaughters(kFALSE);
 141   cutESDPrimary->GetCuts()->SetMinNClustersTPC(50);
 142   cutESDPrimary->GetCuts()->SetMaxChi2PerClusterTPC(3.5);
 143   // -- Bethe-Bloch with kaon mass hypothesis
 144   Double_t sigmaTPC = 0.065;
 145   AliRsnCutBetheBloch *cutKaonBB = new AliRsnCutBetheBloch("cutKaonBB", 3.0 * sigmaTPC, AliPID::kKaon);
 146   cutKaonBB->SetCalibConstant(0, 0.76176e-1);
 147   cutKaonBB->SetCalibConstant(1, 10.632);
 148   cutKaonBB->SetCalibConstant(2, 0.13279e-4);
 149   cutKaonBB->SetCalibConstant(3, 1.8631);
 150   cutKaonBB->SetCalibConstant(4, 1.9479);
 151
 152   // cuts on pairs:
 153   // -- true daughters of a phi resonance (only for true pairs histogram)
 154   AliRsnCutStd *cutTruePair = new AliRsnCutStd("cutTrue", AliRsnCutStd::kTruePair, resonancePDG);
 155   // -- whole interval in Pt and Eta
 156   //AliRsnCutStd *cutEtaPair = new AliRsnCutStd("cutEtaPair", AliRsnCutStd::kEta, -0.9,  0.9);
 157   //AliRsnCutStd *cutPtPair  = new AliRsnCutStd("cutPtPair" , AliRsnCutStd::kPt ,  0.0, 10.0);
 158   //AliRsnCutStd *cutYPair   = new AliRsnCutStd("cutYPair"  , AliRsnCutStd::kEta, -0.9,  0.9);
 159   //cutYPair->SetMass(resonanceMass);
 160
 161   // cuts on event (specific for this analysis):
 162   // -- whole interval in multiplicity
 163   AliRsnCutStd *cutMultiplicity = new AliRsnCutStd("cutMult", AliRsnCutStd::kMult, 0, 200);
 164
 165   // cut set definition for all pairs
 166   AliRsnCutSet *cutSetParticle = new AliRsnCutSet("trackCuts");
 167   cutSetParticle->AddCut(cutESDPrimary);
 168   if (useBBCut)
 169   {
 170     cutSetParticle->AddCut(cutKaonBB);
 171     cutSetParticle->SetCutScheme("cutKaonBB&cutESDPrimary");
 172   }
 173   else
 174   {
 175     cutSetParticle->SetCutScheme("cutESDPrimary");
 176   }
 177
 178   // cut set definition for true pairs
 179   AliRsnCutSet *cutSetPairTrue = new AliRsnCutSet("truePairs");
 180   //cutSetPairTrue->AddCut(cutPtPair);
 181   //cutSetPairTrue->AddCut(cutEtaPair);
 182   //cutSetPairTrue->AddCut(cutTruePair);
 183   //cutSetPairTrue->SetCutScheme("cutPtPair&cutEtaPair&cutTrue");
 184   cutSetPairTrue->SetCutScheme("cutTrue");
 185
 186   // cut set definition for all pairs
 187   //AliRsnCutSet *cutSetPairAll = new AliRsnCutSet("allPairs");
 188   //cutSetPairAll->AddCut(cutPtPair);
 189   //if (useY)
 190   //{
 191   //  cutSetPairAll->AddCut(cutYPair);
 192   //  cutSetPairAll->SetCutScheme("cutPtPair&cutYPair");
 193   //}
 194   //else
 195   //{
 196   //  cutSetPairAll->AddCut(cutEtaPair);
 197   //  cutSetPairAll->SetCutScheme("cutPtPair&cutEtaPair");
 198   //}
 199
 200   // cut set definition for events
 201   //AliRsnCutSet *cutSetEvent = new AliRsnCutSet("cutSetMult");
 202   //cutSetEvent->AddCut(cutMultiplicity);
 203   //cutSetEvent->SetCutScheme("cutMult");
 204
 205   // cut manager for all pairs
 206   // define a proper name for each mult bin, to avoid omonyme output histos
 207   AliRsnCutMgr *cutMgrAll = new AliRsnCutMgr("std", "All");
 208   cutMgrAll->SetCutSet(AliRsnCut::kParticle, cutSetParticle);
 209   //cutMgrAll->SetCutSet(AliRsnCut::kPair, cutSetPairAll);
 210   //cutMgrAll->SetCutSet(AliRsnCut::kEvent, cutSetEvent);
 211
 212   // cut manager for true pairs
 213   AliRsnCutMgr *cutMgrTrue = new AliRsnCutMgr("true", "True");
 214   cutMgrTrue->SetCutSet(AliRsnCut::kParticle, cutSetParticle);
 215   cutMgrTrue->SetCutSet(AliRsnCut::kPair, cutSetPairTrue);
 216   //cutMgrAll->SetCutSet(AliRsnCut::kEvent, cutSetEvent);
 217
 218   for (i = 0; i < nArray; i++) {
 219     pair[i][0]->SetCutMgr(cutMgrTrue);
 220     pair[i][0]->SetOnlyTrue();
 221     for (j = 1; j < 4; j++) {
 222       pair[i][j]->SetCutMgr(cutMgrAll);
 223     }
 224   }
 225
 226   // === FUNCTIONS ================================================================================
 227
 228   // resonance type defines IM binning
 229   Double_t imMin, imMax, imBin;
 230   Int_t    nbins;
 231   if (str.Contains("PHI")) {
 232     imMin = 0.9;
 233     imMax = 1.6;
 234     imBin = 1.0;
 235   }
 236   else
 237   {
 238     imMin = 0.6;
 239     imMax = 2.1;
 240     imBin = 10.0;
 241   }
 242   nbins = (Int_t)((imMax - imMin) / (imBin / 1000.0));
 243   cout << "MIN  = " << imMin << endl;
 244   cout << "MAX  = " << imMax << endl;
 245   cout << "BIN  = " << imBin << endl;
 246   cout << "NBIN = " << nbins << endl;
 247
 248   // define all binnings
 249   AliRsnFunctionAxis *axisIM   = new AliRsnFunctionAxis(AliRsnFunctionAxis::kPairInvMass,   nbins,  imMin, imMax);
 250   AliRsnFunctionAxis *axisPt   = new AliRsnFunctionAxis(AliRsnFunctionAxis::kPairPt,           50,  0.0,  10.0);
 251   AliRsnFunctionAxis *axisY    = new AliRsnFunctionAxis(AliRsnFunctionAxis::kPairY,            20, -1.0,   1.0);
 252   AliRsnFunctionAxis *axisEta  = new AliRsnFunctionAxis(AliRsnFunctionAxis::kPairEta,          20, -1.0,   1.0);
 253   AliRsnFunctionAxis *axisMult = new AliRsnFunctionAxis(AliRsnFunctionAxis::kEventMult,         8,  0.0, 200.0);
 254
 255   AliRsnFunction *fcnPtY   = new AliRsnFunction();
 256   AliRsnFunction *fcnPtEta = new AliRsnFunction();
 257   AliRsnFunction *fcnMult  = new AliRsnFunction();
 258
 259   fcnPtY  ->AddAxis(axisIM);
 260   fcnPtEta->AddAxis(axisIM);
 261   fcnMult ->AddAxis(axisIM);
 262
 263   fcnPtY  ->AddAxis(axisPt);
 264   fcnPtEta->AddAxis(axisPt);
 265   fcnMult ->AddAxis(axisMult);
 266
 267   fcnPtY  ->AddAxis(axisY);
 268   fcnPtEta->AddAxis(axisEta);
 269
 270   // add functions to pairs
 271   for (i = 0; i < nArray; i++) {
 272     for (j = 0; j < 4; j++) {
 273       pair[i][j]->AddFunction(fcnPtY);
 274       pair[i][j]->AddFunction(fcnPtEta);
 275       pair[i][j]->AddFunction(fcnMult);
 276     }
 277   }
 278
 279   // === ADD TO PAIR MANAGER ======================================================================
 280
 281   //
 282   // true pairs are array [0] element:
 283   // if they must be excluded, this element is skipped everywhere
 284   //
 285
 286   for (i = 0; i < nArray; i++) {
 287     Int_t start = (addTruePairs ? 0 : 1);
 288     for (j = start; j < 4; j++) {
 289       pairMgr->AddPair(pair[i][j]);
 290     }
 291   }
 292
 293   return pairMgr;
 294 }
 295
 296 //_________________________________________________________________________________________________
 297 void RsnConfigTask(AliRsnAnalysisSE* &task)
 298 //
 299 // -- RsnConfigTask --
 300 //
 301 // This function configures the entire task for all resonances the user is interested in.
 302 // It recalls the above function many times to configure all required pair managers.
 303 //
 304 // By default, this function has the used task as argument.
 305 //
 306 {
 307   if (!task)
 308   {
 309     Error("ConfigTaskRsn", "Task not found");
 310     return;
 311   }
 312
 313   // set prior probabilities for PID
 314   task->SetPriorProbability(AliPID::kElectron, 0.02);
 315   task->SetPriorProbability(AliPID::kMuon,     0.02);
 316   task->SetPriorProbability(AliPID::kPion,     0.83);
 317   task->SetPriorProbability(AliPID::kKaon,     0.07);
 318   task->SetPriorProbability(AliPID::kProton,   0.06);
 319   task->DumpPriors();
 320
 321   // initialize analysis manager with pairs from config
 322   AliRsnAnalysisManager *anaMgr = task->GetAnalysisManager(0);
 323
 324   // create pair managers for phi
 325   anaMgr->Add(RsnConfigPairManager("PHI_NOPID", 333, 1.0193, AliPID::kKaon, AliPID::kKaon, kTRUE));
 326   anaMgr->Add(RsnConfigPairManager("PHI_BB"   , 333, 1.0193, AliPID::kKaon, AliPID::kKaon, kTRUE));
 327   anaMgr->Add(RsnConfigPairManager("PHI_PID"  , 333, 1.0193, AliPID::kKaon, AliPID::kKaon, kTRUE));
 328   // create pair managers for kstar
 329   anaMgr->Add(RsnConfigPairManager("KSTAR_NOPID", 313, 0.896, AliPID::kPion, AliPID::kKaon, kTRUE));
 330   anaMgr->Add(RsnConfigPairManager("KSTAR_BB"   , 313, 0.896, AliPID::kPion, AliPID::kKaon, kTRUE));
 331   anaMgr->Add(RsnConfigPairManager("KSTAR_PID"  , 313, 0.896, AliPID::kPion, AliPID::kKaon, kTRUE));
 332
 333   // setup cuts for events (good primary vertex)
 334   AliRsnCutPrimaryVertex *cutVertex   = new AliRsnCutPrimaryVertex("cutVertex", 3);
 335   AliRsnCutSet           *cutSetEvent = new AliRsnCutSet("eventCuts");
 336   cutSetEvent->AddCut(cutVertex);
 337   cutSetEvent->SetCutScheme("cutVertex");
 338   task->SetEventCuts(cutSetEvent);
 339 }