PWG2/RESONANCES/macros/AddAnalysisTaskRsn.C

   1 //
   2 // Macro to create the full analysis manager for Resonances
   3 //
   4
   5 static Double_t  cov11 = 2;
   6 static Double_t  cov22 = 2;
   7 static Double_t  cov33 = 0.5;
   8 static Double_t  cov44 = 0.5;
   9 static Double_t  cov55 = 2;
  10 static Double_t  nSigmaToVertex = 4;
  11 static Double_t  dcaToVertex = 3.0;
  12 static Double_t  maxChi2PerClusterTPC = 3.5;
  13 static Bool_t    requireTPCRefit = kTRUE;
  14 static Bool_t    requireSigmaToVertex = kTRUE;
  15 static Bool_t    acceptKinkDaughters = kFALSE;
  16 static Int_t     minNClustersTPC = 50;
  17
  18 Bool_t AddAnalysisTaskRsn(const char *outFile = "rsn.root", Bool_t useKine = kTRUE)
  19 {
  20   // Retrieve analysis manager.
  21   // Since it is usually created in the steering macro,
  22   // then we don't ever need to initialize a new one
  23   AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
  24
  25   // no need to  create handlers for input and MonteCarlo
  26   // Done in the general steering macro
  27   //  AliMCEventHandler  *mcHandler  = new AliMCEventHandler();
  28   //  AliESDInputHandler *esdHandler = new AliESDInputHandler();
  29   //AliAODHandler      *aodHandler = new AliAODHandler();
  30   //aodHandler->SetOutputFileName("AliAOD.root");
  31   //  if (useKine) mgr->SetMCtruthEventHandler(mcHandler);
  32   //  mgr->SetInputEventHandler(esdHandler);
  33   //mgr->SetOutputEventHandler(aodHandler);
  34
  35   // initialize task
  36   AliRsnAnalysisSE *task = new AliRsnAnalysisSE("AliRsnAnalysisSE");
  37
  38   // set prior probabilities for PID
  39   task->SetPriorProbability(AliPID::kElectron, 0.02);
  40   task->SetPriorProbability(AliPID::kMuon,     0.02);
  41   task->SetPriorProbability(AliPID::kPion,     0.83);
  42   task->SetPriorProbability(AliPID::kKaon,     0.07);
  43   task->SetPriorProbability(AliPID::kProton,   0.06);
  44   task->DumpPriors();
  45
  46   // initialize analysis manager with pairs from config
  47   AliRsnAnalysisManager *anaMgr = task->GetAnalysisManager("MyAnalysisSE");
  48   anaMgr->Add(CreatePairsNoPID("PHI_NoPID_0"   , 333, AliPID::kKaon, AliPID::kKaon, 10000.0));
  49   anaMgr->Add(CreatePairsNoPID("PHI_NoPID_10"  , 333, AliPID::kKaon, AliPID::kKaon,     0.1));
  50   anaMgr->Add(CreatePairsNoPID("PHI_NoPID_20"  , 333, AliPID::kKaon, AliPID::kKaon,     0.2));
  51   anaMgr->Add(CreatePairsPID  ("PHI_PID"       , 333, AliPID::kKaon, AliPID::kKaon));
  52   anaMgr->Add(CreatePairsNoPID("KSTAR_NoPID_0" , 313, AliPID::kKaon, AliPID::kPion, 10000.0));
  53   anaMgr->Add(CreatePairsNoPID("KSTAR_NoPID_10", 313, AliPID::kKaon, AliPID::kPion,     0.1));
  54   anaMgr->Add(CreatePairsNoPID("KSTAR_NoPID_20", 313, AliPID::kKaon, AliPID::kPion,     0.2));
  55   anaMgr->Add(CreatePairsPID  ("KSTAR_PID"     , 313, AliPID::kKaon, AliPID::kPion));
  56
  57   // setup cuts
  58   AliESDtrackCuts *esdCuts = new AliESDtrackCuts;
  59   esdCuts->SetMaxCovDiagonalElements(cov11, cov22, cov33, cov44, cov55);
  60   esdCuts->SetRequireSigmaToVertex(requireSigmaToVertex);
  61   if (requireSigmaToVertex) esdCuts->SetMaxNsigmaToVertex(nSigmaToVertex);
  62   else
  63   {
  64     esdCuts->SetDCAToVertexZ(dcaToVertex);
  65     esdCuts->SetDCAToVertexXY(dcaToVertex);
  66   }
  67   esdCuts->SetRequireTPCRefit(requireTPCRefit);
  68   esdCuts->SetAcceptKinkDaughters(acceptKinkDaughters);
  69   esdCuts->SetMinNClustersTPC(minNClustersTPC);
  70   esdCuts->SetMaxChi2PerClusterTPC(maxChi2PerClusterTPC);
  71   task->SetESDtrackCuts(esdCuts);
  72
  73   // add the task to manager
  74   mgr->AddTask(task);
  75
  76   // initialize container for the output
  77   // (input container is common)
  78   AliAnalysisDataContainer *out = mgr->CreateContainer("RSN", TList::Class(), AliAnalysisManager::kOutputContainer, outFile);
  79
  80   // connect input and output slots for this task
  81   mgr->ConnectInput(task, 0, mgr->GetCommonInputContainer());
  82   mgr->ConnectOutput(task, 1, out);
  83 }
  84
  85 void AddRsnFunctionsToPair(AliRsnPair *pair)
  86 {
  87 //
  88 // In general, for all processed pairs in one analysis the same functions are computed.
  89 // Then, they are defined separately here and added in the same way to all pairs.
  90 //
  91
  92   // define histogram templates
  93   AliRsnHistoDef *hdIM   = new AliRsnHistoDef(1000,  0.0,   2.0);
  94   AliRsnHistoDef *hdPt   = new AliRsnHistoDef(  10,  0.0,  10.0);
  95   AliRsnHistoDef *hdMult = new AliRsnHistoDef(   8,  0.0, 200.0);
  96
  97   // define functions:
  98   // -- one computed w.r. to Pt and Mult
  99   // -- one computed w.r. to Eta and Mult
 100   AliRsnFunction *fcnIM1 = new AliRsnFunction(AliRsnFunction::kInvMass, AliRsnFunction::kBinPairPt, AliRsnFunction::kBinEventMult, hdIM, hdPt, hdMult);
 101
 102   // add to pair
 103   pair->AddFunction(fcnIM1);
 104 }
 105
 106 AliRsnPairManager* CreatePairsNoPID
 107 (
 108   const char            *pairMgrName,    // name for the pair manager
 109   Int_t                  resonancePDG,   // PDG code of resonance (for true pairs)
 110   AliPID::EParticleType  type1,          // PID of one member of decay (+)
 111   AliPID::EParticleType  type2,          // PID of other member of decay (-)
 112   Double_t               bbCut = 10000.0 // cut on Bethe-Bloch
 113 )
 114 {
 115 //
 116 // Creates an AliRsnPairMgr for a specified resonance, which contains:
 117 // - signal (inv. mass)
 118 // - event mixing (inv. mass)
 119 // - like-sign (inv. mass)
 120 // - true pairs (inv. mass, resolution)
 121 //
 122 // For all pairs, a binning in Pt and Eta is provided, and a cut in multiplicity
 123 // which defines a multiplicity bin where the analysis is computed.
 124 //
 125 // Arguments define how the pair manager must be built, and are explained above
 126 //
 127
 128   AliRsnPairManager  *pairMgr  = new AliRsnPairManager(pairMgrName);
 129
 130   // === PAIR DEFINITIONS =========================================================================
 131
 132   // if particle #1 and #2 are different, two histograms must be built
 133   // for each scheme (signal, true, mixed, like-sign) exchanging both particles and signs
 134   Int_t i, j, nArray = 1;
 135   if (type1 != type2) nArray = 2;
 136
 137   AliRsnPairDef *defUnlike[2] = {0, 0};
 138   AliRsnPairDef *defLikePP[2] = {0, 0};
 139   AliRsnPairDef *defLikeMM[2] = {0, 0};
 140
 141   defUnlike[0] = new AliRsnPairDef(type1, '+', type2, '-', resonancePDG);
 142   defLikePP[0] = new AliRsnPairDef(type1, '+', type2, '+', resonancePDG);
 143   defLikeMM[0] = new AliRsnPairDef(type1, '-', type2, '-', resonancePDG);
 144
 145   defUnlike[1] = new AliRsnPairDef(type2, '+', type1, '-', resonancePDG);
 146   defLikePP[1] = new AliRsnPairDef(type2, '+', type1, '+', resonancePDG);
 147   defLikeMM[1] = new AliRsnPairDef(type2, '-', type1, '-', resonancePDG);
 148
 149   // === PAIR ANALYSIS ENGINES ====================================================================
 150
 151   // define null (dummy) objects and initialize only the ones which are needed,
 152   // depending again on particle types;
 153   // array is organized as follows:
 154   // [0] - true pairs
 155   // [1] - signal
 156   // [2] - like PP
 157   // [3] - like MM
 158   AliRsnPair *noPID[2][4] = {0,0,0,0,0,0,0,0};
 159
 160   for (i = 0; i < nArray; i++) {
 161     noPID[i][0] = new AliRsnPair(AliRsnPair::kNoPID, defUnlike[i]);
 162     noPID[i][1] = new AliRsnPair(AliRsnPair::kNoPID, defUnlike[i]);
 163     noPID[i][2] = new AliRsnPair(AliRsnPair::kNoPID, defLikePP[i]);
 164     noPID[i][3] = new AliRsnPair(AliRsnPair::kNoPID, defLikeMM[i]);
 165   }
 166
 167   // === CUTS =====================================================================================
 168
 169   // cuts for tracks:
 170   // -- Bethe-Bloch & required kaon PID
 171   AliRsnCutStd *cutReqKaon = new AliRsnCutStd("cutReqKaon", AliRsnCutStd::kRequiredPID, (Int_t)AliPID::kKaon);
 172   AliRsnCutBetheBloch *cutKaonBB = new AliRsnCutBetheBloch("cutKaon", bbCut, AliPID::kKaon);
 173   cutKaonBB->SetCalibConstant(0, 0.76176e-1);
 174   cutKaonBB->SetCalibConstant(1, 10.632);
 175   cutKaonBB->SetCalibConstant(2, 0.13279e-4);
 176   cutKaonBB->SetCalibConstant(3, 1.8631);
 177   cutKaonBB->SetCalibConstant(4, 1.9479);
 178
 179   // cuts on pairs:
 180   // -- true daughters of a phi resonance (only for true pairs histogram)cutSetPairTrue->AddCut(cutTrue);
 181   AliRsnCutStd *cutTruePair = new AliRsnCutStd("cutTrue", AliRsnCutStd::kTruePair, resonancePDG);
 182
 183   // cut set definition for all pairs
 184   AliRsnCutSet *cutSetParticle = new AliRsnCutSet("trackCuts");
 185   cutSetParticle->AddCut(cutKaonBB);
 186   cutSetParticle->AddCut(cutReqKaon);
 187   cutSetParticle->SetCutScheme("(cutReqKaon&cutKaonBB)|(!cutReqKaon)");
 188
 189   // cut set definition for true pairs
 190   AliRsnCutSet *cutSetPairTrue = new AliRsnCutSet("truePairs");
 191   cutSetPairTrue->AddCut(cutTruePair);
 192   cutSetPairTrue->SetCutScheme("cutTrue");
 193
 194   // cut manager for all pairs
 195   AliRsnCutMgr *cutMgrAll = new AliRsnCutMgr("std", "All");
 196   cutMgrAll->SetCutSet(AliRsnCut::kParticle, cutSetParticle);
 197
 198   // cut manager for all pairs
 199   AliRsnCutMgr *cutMgrTrue = new AliRsnCutMgr("true", "True");
 200   cutMgrTrue->SetCutSet(AliRsnCut::kParticle, cutSetParticle);
 201   cutMgrTrue->SetCutSet(AliRsnCut::kPair, cutSetPairTrue);
 202
 203   for (i = 0; i < nArray; i++) {
 204     noPID[i][0]->SetCutMgr(cutMgrTrue);
 205     for (j = 1; j < 4; j++) {
 206       noPID[i][j]->SetCutMgr(cutMgrAll);
 207     }
 208   }
 209
 210   // === FUNCTIONS ================================================================================
 211
 212   for (i = 0; i < nArray; i++) {
 213     for (j = 0; j < 4; j++) {
 214       AddRsnFunctionsToPair(noPID[i][j]);
 215     }
 216   }
 217
 218   // === ADD TO PAIR MANAGER ======================================================================
 219
 220   for (i = 0; i < nArray; i++) {
 221     for (j = 0; j < 4; j++) {
 222       pairMgr->AddPair(noPID[i][j]);
 223     }
 224   }
 225
 226   return pairMgr;
 227 }
 228
 229 AliRsnPairManager* CreatePairsPID
 230 (
 231   const char            *pairMgrName,    // name for the pair manager
 232   Int_t                  resonancePDG,   // PDG code of resonance (for true pairs)
 233   AliPID::EParticleType  type1,          // PID of one member of decay (+)
 234   AliPID::EParticleType  type2           // PID of other member of decay (-)
 235 )
 236 {
 237 //
 238 // Creates an AliRsnPairMgr for a specified resonance, which contains:
 239 // - signal (inv. mass)
 240 // - event mixing (inv. mass)
 241 // - like-sign (inv. mass)
 242 // - true pairs (inv. mass, resolution)
 243 //
 244 // For all pairs, a binning in Pt and Eta is provided, and a cut in multiplicity
 245 // which defines a multiplicity bin where the analysis is computed.
 246 //
 247 // Arguments define how the pair manager must be built, and are explained above
 248 //
 249
 250   AliRsnPairManager  *pairMgr  = new AliRsnPairManager(pairMgrName);
 251
 252   // === PAIR DEFINITIONS =========================================================================
 253
 254   // if particle #1 and #2 are different, two histograms must be built
 255   // for each scheme (signal, true, mixed, like-sign) exchanging both particles and signs
 256   Int_t i, j, nArray = 1;
 257   if (type1 != type2) nArray = 2;
 258
 259   AliRsnPairDef *defUnlike[2] = {0, 0};
 260   AliRsnPairDef *defLikePP[2] = {0, 0};
 261   AliRsnPairDef *defLikeMM[2] = {0, 0};
 262
 263   defUnlike[0] = new AliRsnPairDef(type1, '+', type2, '-', resonancePDG);
 264   defLikePP[0] = new AliRsnPairDef(type1, '+', type2, '+', resonancePDG);
 265   defLikeMM[0] = new AliRsnPairDef(type1, '-', type2, '-', resonancePDG);
 266
 267   defUnlike[1] = new AliRsnPairDef(type2, '+', type1, '-', resonancePDG);
 268   defLikePP[1] = new AliRsnPairDef(type2, '+', type1, '+', resonancePDG);
 269   defLikeMM[1] = new AliRsnPairDef(type2, '-', type1, '-', resonancePDG);
 270
 271   // === PAIR ANALYSIS ENGINES ====================================================================
 272
 273   // define null (dummy) objects and initialize only the ones which are needed,
 274   // depending again on particle types;
 275   // array is organized as follows:
 276   // [0] - true pairs
 277   // [1] - signal
 278   // [2] - like PP
 279   // [3] - like MM
 280   AliRsnPair *perfectPID[2][4]   = {0,0,0,0,0,0,0,0};
 281   AliRsnPair *realisticPID[2][4] = {0,0,0,0,0,0,0,0};
 282
 283   for (i = 0; i < nArray; i++) {
 284     perfectPID[i][0] = new AliRsnPair(AliRsnPair::kPerfectPID, defUnlike[i]);
 285     perfectPID[i][1] = new AliRsnPair(AliRsnPair::kPerfectPID, defUnlike[i]);
 286     perfectPID[i][2] = new AliRsnPair(AliRsnPair::kPerfectPID, defLikePP[i]);
 287     perfectPID[i][3] = new AliRsnPair(AliRsnPair::kPerfectPID, defLikeMM[i]);
 288
 289     realisticPID[i][0] = new AliRsnPair(AliRsnPair::kRealisticPID, defUnlike[i]);
 290     realisticPID[i][1] = new AliRsnPair(AliRsnPair::kRealisticPID, defUnlike[i]);
 291     realisticPID[i][2] = new AliRsnPair(AliRsnPair::kRealisticPID, defLikePP[i]);
 292     realisticPID[i][3] = new AliRsnPair(AliRsnPair::kRealisticPID, defLikeMM[i]);
 293   }
 294
 295   // === CUTS =====================================================================================
 296
 297   // cuts for tracks:
 298   // -- nothing
 299
 300   // cuts on pairs:
 301   // -- true daughters of a phi resonance (only for true pairs histogram)cutSetPairTrue->AddCut(cutTrue);
 302   AliRsnCutStd *cutTruePair = new AliRsnCutStd("cutTrue", AliRsnCutStd::kTruePair, resonancePDG);
 303
 304   // cut set definition for true pairs
 305   AliRsnCutSet *cutSetPairTrue = new AliRsnCutSet("truePairs");
 306   cutSetPairTrue->AddCut(cutTruePair);
 307   cutSetPairTrue->SetCutScheme("cutTrue");
 308
 309   // cut manager for all pairs
 310   AliRsnCutMgr *cutMgrAll = new AliRsnCutMgr("std", "All");
 311
 312   // cut manager for all pairs
 313   AliRsnCutMgr *cutMgrTrue = new AliRsnCutMgr("true", "True");
 314   cutMgrTrue->SetCutSet(AliRsnCut::kPair, cutSetPairTrue);
 315
 316   for (i = 0; i < nArray; i++) {
 317     perfectPID[i][0]->SetCutMgr(cutMgrTrue);
 318     realisticPID[i][0]->SetCutMgr(cutMgrTrue);
 319     for (j = 1; j < 4; j++) {
 320       perfectPID[i][j]->SetCutMgr(cutMgrAll);
 321       realisticPID[i][j]->SetCutMgr(cutMgrAll);
 322     }
 323   }
 324
 325   // === FUNCTIONS ================================================================================
 326
 327   for (i = 0; i < nArray; i++) {
 328     for (j = 0; j < 4; j++) {
 329       AddRsnFunctionsToPair(perfectPID[i][j]);
 330       AddRsnFunctionsToPair(realisticPID[i][j]);
 331     }
 332   }
 333
 334   // === ADD TO PAIR MANAGER ======================================================================
 335
 336   for (i = 0; i < nArray; i++) {
 337     for (j = 0; j < 4; j++) {
 338       pairMgr->AddPair(perfectPID[i][j]);
 339       pairMgr->AddPair(realisticPID[i][j]);
 340     }
 341   }
 342
 343   return pairMgr;
 344 }
 345