AliPhysics  59e0e03 (59e0e03)
AliAnalysisTaskDmesonJetsDetectorResponse.cxx
Go to the documentation of this file.
1 /*************************************************************************
2 * Copyright(c) 1998-2016, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15 
16 #include "AliAnalysisManager.h"
17 #include "AliVEventHandler.h"
18 #include "AliMCEventHandler.h"
19 #include "AliHFTrackContainer.h"
21 
23 
24 // Definitions of class AliAnalysisTaskDmesonJetsDetectorResponse::AliDmesonMatchInfoSummary
25 
29 
30 // Definitions of class AliAnalysisTaskDmesonJetsDetectorResponse::AliD0MatchInfoSummary
31 
35 
38 {
39  fGenerated.Reset();
41 }
42 
47 {
48  fReconstructed.Set(reco);
49 }
50 
55 {
56  fGenerated.Set(truth);
57 }
58 
59 // Definitions of class AliAnalysisTaskDmesonJetsDetectorResponse::AliDStarMatchInfoSummary
60 
64 
67 {
68  fGenerated.Reset();
70 }
71 
76 {
77  fReconstructed.Set(reco);
78 }
79 
84 {
85  fGenerated.Set(truth);
86 }
87 
88 // Definitions of class AliAnalysisTaskDmesonJetsDetectorResponse::ResponseEngine
89 
93 
96  TObject(),
97  fCandidateType(kD0toKpi),
98  fInhibit(kFALSE),
99  fMaxJetDmesonDistance(0),
100  fName(),
101  fTree(0),
102  fCurrentDmeson(0),
103  fCurrentJetInfoReco(0),
104  fCurrentJetInfoTruth(0),
105  fDataSlotNumber(-1),
106  fRecontructed(0),
107  fGenerated(0),
109 {
110 
111 }
112 
117  TObject(),
118  fCandidateType(type),
119  fInhibit(kFALSE),
121  fName(),
122  fTree(0),
123  fCurrentDmeson(0),
126  fDataSlotNumber(-1),
127  fRecontructed(0),
128  fGenerated(0),
130 {
131 
132 }
133 
138  TObject(source),
140  fInhibit(source.fInhibit),
142  fName(source.fName),
143  fTree(0),
144  fCurrentDmeson(0),
149  fGenerated(source.fGenerated),
151 {
152 }
153 
158 {
159  new (this) ResponseEngine(source);
160  return *this;
161 }
162 
164 {
165  if (lhs.fRecontructed && rhs.fRecontructed) return (*(lhs.fRecontructed) < *(rhs.fRecontructed));
166  if (!lhs.fRecontructed && rhs.fRecontructed) return true;
167  return false;
168 }
169 
170 
172 {
173  if (lhs.fRecontructed && rhs.fRecontructed) return (*(lhs.fRecontructed) == *(rhs.fRecontructed));
174  if (!lhs.fRecontructed && !rhs.fRecontructed) return true;
175  return false;
176 }
177 
182 {
183  fInhibit = (fRecontructed == 0 || fGenerated == 0);
185  return !fInhibit;
186 }
187 
190 {
191 
192 }
193 
199 {
200  TString classname;
201  switch (fCandidateType) {
202  case kD0toKpi:
203  case kD0toKpiLikeSign:
204  classname = "AliAnalysisTaskDmesonJetsDetectorResponse::AliD0MatchInfoSummary";
206  break;
207  case kDstartoKpipi:
208  classname = "AliAnalysisTaskDmesonJetsDetectorResponse::AliDStarMatchInfoSummary";
210  break;
211  }
212  TString treeName = TString::Format("%s_%s", taskName, GetName());
213  fTree = new TTree(treeName, treeName);
214  fTree->Branch("DmesonJet", classname, &fCurrentDmeson);
215 
217  for (Int_t i = 0; i < fGenerated->GetJetDefinitions().size(); i++) {
218  if (fGenerated->GetJetDefinitions()[i].GetRhoName().IsNull()) {
220  TString bname = TString::Format("%s_truth", fGenerated->GetJetDefinitions()[i].GetName());
221  fTree->Branch(bname, "AliAnalysisTaskDmesonJets::AliJetInfoSummary", &fCurrentJetInfoTruth[i]);
222  }
223  else {
225  TString bname = TString::Format("%s_truth", fGenerated->GetJetDefinitions()[i].GetName());
226  fTree->Branch(bname, "AliAnalysisTaskDmesonJets::AliJetInfoPbPbSummary", &fCurrentJetInfoTruth[i]);
227  }
228  }
229 
231  for (Int_t i = 0; i < fRecontructed->GetJetDefinitions().size(); i++) {
232  if (fGenerated->GetJetDefinitions()[i].GetRhoName().IsNull()) {
234  TString bname = TString::Format("%s_reco", fRecontructed->GetJetDefinitions()[i].GetName());
235  fTree->Branch(bname, "AliAnalysisTaskDmesonJets::AliJetInfoSummary", &fCurrentJetInfoReco[i]);
236  }
237  else {
239  TString bname = TString::Format("%s_reco", fRecontructed->GetJetDefinitions()[i].GetName());
240  fTree->Branch(bname, "AliAnalysisTaskDmesonJets::AliJetInfoPbPbSummary", &fCurrentJetInfoReco[i]);
241  }
242  }
243 
244  return fTree;
245 }
246 
252 {
253  std::map<int, AliDmesonJetInfo>& recoDmesons = fRecontructed->GetDmesons();
254  std::map<int, AliDmesonJetInfo>& truthDmesons = fGenerated->GetDmesons();
255 
256  TString hname;
257  // Loop over reconstructed D meson and look for their generated counterparts
258  for (auto& dmeson_reco : recoDmesons) {
259  if (dmeson_reco.second.fMCLabel < 0) {
260  hname = TString::Format("%s/fHistGeneratedDMesonNotFoundRecoPt", GetName());
261  fHistManagerResponse->FillTH1(hname, dmeson_reco.second.fD.Pt());
262 
263  hname = TString::Format("%s/fHistGeneratedDMesonNotFoundRecoEta", GetName());
264  fHistManagerResponse->FillTH1(hname, dmeson_reco.second.fD.Eta());
265 
266  hname = TString::Format("%s/fHistGeneratedDMesonNotFoundRecoPhi", GetName());
267  fHistManagerResponse->FillTH1(hname, dmeson_reco.second.fD.Pt());
268  continue; // it should never happen
269  }
270 
271  // Reset D meson and jet tree branches
273  for (UInt_t ij = 0; ij < fRecontructed->GetJetDefinitions().size(); ij++) {
274  fCurrentJetInfoReco[ij]->Reset();
275  }
276  for (UInt_t ij = 0; ij < fGenerated->GetJetDefinitions().size(); ij++) {
278  }
279 
280  // Copy the reconstructed D meson information in the tree branch
281  fCurrentDmeson->SetReconstructed(dmeson_reco.second);
282 
283  Int_t accRecJets = 0;
284  Int_t accGenJets = 0;
285 
286  // Copy the reconstructed jet information in the tree branch
287  for (UInt_t ij = 0; ij < fRecontructed->GetJetDefinitions().size(); ij++) {
288  AliJetInfo* jet = dmeson_reco.second.GetJet(fRecontructed->GetJetDefinitions()[ij].GetName());
289  if (!jet) continue;
290  if (applyKinCuts && !fRecontructed->GetJetDefinitions()[ij].IsJetInAcceptance(*jet)) continue;
291  fCurrentJetInfoReco[ij]->Set(dmeson_reco.second, fRecontructed->GetJetDefinitions()[ij].GetName());
292  accRecJets++;
293  }
294 
295  // Look for the generated D meson counterpart
296  std::map<int, AliDmesonJetInfo>::iterator it = truthDmesons.find(dmeson_reco.second.fMCLabel);
297  if (it != truthDmesons.end()) {
298  std::pair<const int, AliDmesonJetInfo>& dmeson_truth = (*it);
299  // Flag the generated D meson as reconstructed
300  dmeson_truth.second.fReconstructed = kTRUE;
301  // Copy the generated D meson information in the tree branch
302  fCurrentDmeson->SetGenerated((*it).second);
303 
304  // Copy the generated jet information in the tree branch
305  for (UInt_t ij = 0; ij < fGenerated->GetJetDefinitions().size(); ij++) {
306  AliJetInfo* jet = dmeson_truth.second.GetJet(fGenerated->GetJetDefinitions()[ij].GetName());
307  if (!jet) continue;
308  if (!applyKinCuts || fGenerated->GetJetDefinitions()[ij].IsJetInAcceptance(*jet)) accGenJets++;
309  fCurrentJetInfoTruth[ij]->Set(dmeson_truth.second, fGenerated->GetJetDefinitions()[ij].GetName());
310  }
311  }
312  else {
313  hname = TString::Format("%s/fHistGeneratedDMesonOutsideAccRecoPt", GetName());
314  fHistManagerResponse->FillTH1(hname, dmeson_reco.second.fD.Pt());
315 
316  hname = TString::Format("%s/fHistGeneratedDMesonOutsideAccRecoEta", GetName());
317  fHistManagerResponse->FillTH1(hname, dmeson_reco.second.fD.Eta());
318 
319  hname = TString::Format("%s/fHistGeneratedDMesonOutsideAccRecoPhi", GetName());
320  fHistManagerResponse->FillTH1(hname, dmeson_reco.second.fD.Pt());
321  continue;
322  }
323 
324  // Fill the tree with the current D meson
325  if (accRecJets > 0 || accGenJets > 0) fTree->Fill();
326  }
327 
328  // Reset jet tree branches
329  for (UInt_t ij = 0; ij < fRecontructed->GetJetDefinitions().size(); ij++) fCurrentJetInfoReco[ij]->Reset();
330  for (UInt_t ij = 0; ij < fGenerated->GetJetDefinitions().size(); ij++) fCurrentJetInfoTruth[ij]->Reset();
331 
332  // Loop over generated D meson that were not reconstructed
333  for (auto& dmeson_truth : truthDmesons) {
334  // Skip if the generated D meson was reconstructed
335  if (dmeson_truth.second.fReconstructed) continue;
336 
337  // Reset D meson tree branch
339 
340  // Copy the generated D meson information in the tree branch
341  fCurrentDmeson->SetGenerated(dmeson_truth.second);
342 
343  // Copy the reconstructed jet information in the tree branch
344  Int_t accGenJets = 0;
345  for (UInt_t ij = 0; ij < fGenerated->GetJetDefinitions().size(); ij++) {
347  AliJetInfo* jet = dmeson_truth.second.GetJet(fGenerated->GetJetDefinitions()[ij].GetName());
348  if (!jet) continue;
349  if (!applyKinCuts || fGenerated->GetJetDefinitions()[ij].IsJetInAcceptance(*jet)) accGenJets++;
350  fCurrentJetInfoTruth[ij]->Set(dmeson_truth.second, fGenerated->GetJetDefinitions()[ij].GetName());
351  }
352 
353  // Fill the tree with the current D meson
354  if (accGenJets > 0) fTree->Fill();
355  }
356 
357  for (auto& dmeson_truth : truthDmesons) {
358  dmeson_truth.second.fReconstructed = kFALSE;
359  }
360 
361  return kTRUE;
362 }
363 
364 // Definitions of class AliAnalysisTaskDmesonJetsDetectorResponse
365 
369 
375 {
377 }
378 
383  AliAnalysisTaskDmesonJets(name, nOutputTrees),
384  fHistManagerResponse(TString::Format("%s_QA", name)),
386 {
388 }
389 
392 {
393  ::Info("UserCreateOutputObjects", "CreateOutputObjects of task %s", GetName());
394 
396 
397  for (auto &param : fAnalysisEngines) {
398  if (param.IsInhibit()) continue;
399  if (param.GetMCMode() != kSignalOnly) continue;
400 
401  ResponseEngine resp(param.GetCandidateType());
402  resp.SetReconstructedAnalysisEngine(&param);
403  fResponseEngines.push_back(resp);
404  }
405 
406  for (auto &resp : fResponseEngines) {
407  for (auto &param : fAnalysisEngines) {
408  if (param.IsInhibit()) continue;
409  if (param.GetMCMode() != kMCTruth) continue;
410  if (resp.fRecontructed->GetCandidateType() != param.GetCandidateType()) continue;
411 
412  resp.SetGeneratedAnalysisEngine(&param);
413  }
414  }
415 
416  Int_t treeSlot = 0;
417 
418  TString hname;
419  TString htitle;
420  TH1* h = nullptr;
421 
422  for (auto &resp : fResponseEngines) {
423  if (!resp.CheckInit()) continue;
424 
425  resp.BuildTree(GetName());
426  if (treeSlot < fNOutputTrees) {
427  resp.AssignDataSlot(treeSlot+2);
428  treeSlot++;
430  }
431  else {
432  AliError(Form("Number of data output slots %d not sufficient. Tree of response engine %s will not be posted!", fNOutputTrees, resp.GetName()));
433  }
434 
435  resp.fHistManagerResponse = &fHistManagerResponse;
436 
437  hname = TString::Format("%s/fHistGeneratedDMesonNotFoundRecoPt", resp.GetName());
438  htitle = hname + ";#it{p}_{T,D} (GeV/#it{c});counts";
439  h = fHistManagerResponse.CreateTH1(hname, htitle, 300, 0, 150);
440 
441  hname = TString::Format("%s/fHistGeneratedDMesonNotFoundRecoEta", resp.GetName());
442  htitle = hname + ";#eta_{D};counts";
443  h = fHistManagerResponse.CreateTH1(hname, htitle, 400, -5, 5);
444 
445  hname = TString::Format("%s/fHistGeneratedDMesonNotFoundRecoPhi", resp.GetName());
446  htitle = hname + ";#phi_{D};counts";
447  h = fHistManagerResponse.CreateTH1(hname, htitle, 400, 0, TMath::TwoPi());
448 
449  hname = TString::Format("%s/fHistGeneratedDMesonOutsideAccRecoPt", resp.GetName());
450  htitle = hname + ";#it{p}_{T,D} (GeV/#it{c});counts";
451  h = fHistManagerResponse.CreateTH1(hname, htitle, 300, 0, 150);
452 
453  hname = TString::Format("%s/fHistGeneratedDMesonOutsideAccRecoEta", resp.GetName());
454  htitle = hname + ";#eta_{D};counts";
455  h = fHistManagerResponse.CreateTH1(hname, htitle, 400, -5, 5);
456 
457  hname = TString::Format("%s/fHistGeneratedDMesonOutsideAccRecoPhi", resp.GetName());
458  htitle = hname + ";#phi_{D};counts";
459  h = fHistManagerResponse.CreateTH1(hname, htitle, 400, 0, TMath::TwoPi());
460  }
461 
463 
464  PostData(1, fOutput);
465 }
466 
470 {
472 }
473 
478 {
480 }
481 
486 {
487  for (auto &resp : fResponseEngines) {
488  if (resp.IsInhibit()) continue;
489 
490  resp.FillTree(fApplyKinematicCuts);
491 
493  }
494 
495  for (auto &ana : fAnalysisEngines) {
496  if (ana.IsInhibit()) continue;
497 
498  ana.GetOutputHandler()->FillOutput(fApplyKinematicCuts);
499  }
500 
502  return kTRUE;
503 
504 }
505 
510 {
512  AliWarning("This class only provides a tree output.");
513  }
514 
515  // Always set it to kOnlyQAOutput: base class does not generate any output (other than QA histograms), all the output comes from the derived class
517 }
518 
525 {
526  if (eng.GetDataSlotNumber() >= 0 && eng.GetTree()) {
527  PostData(eng.GetDataSlotNumber(), eng.GetTree());
528  return eng.GetDataSlotNumber();
529  }
530  else {
531  return -1;
532  }
533 }
534 
544 {
545  // Get the pointer to the existing analysis manager via the static access method.
546  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
547  if (!mgr) {
548  ::Error("AddTaskDmesonJetsDetectorResponse", "No analysis manager to connect to.");
549  return 0;
550  }
551 
552  // Check the analysis type using the event handlers connected to the analysis manager.
553  AliVEventHandler* handler = mgr->GetInputEventHandler();
554  if (!handler) {
555  ::Error("AddTaskDmesonJetsDetectorResponse", "This task requires an input event handler");
556  return 0;
557  }
558 
560 
561  if (handler->InheritsFrom("AliESDInputHandler")) {
562  dataType = kESD;
563  }
564  else if (handler->InheritsFrom("AliAODInputHandler")) {
565  dataType = kAOD;
566  }
567 
568  // Init the task and do settings
569  if (trackName == "usedefault") {
570  if (dataType == kESD) {
571  trackName = "Tracks";
572  }
573  else if (dataType == kAOD) {
574  trackName = "tracks";
575  }
576  else {
577  trackName = "";
578  }
579  }
580 
581  if (clusName == "usedefault") {
582  if (dataType == kESD) {
583  clusName = "CaloClusters";
584  }
585  else if (dataType == kAOD) {
586  clusName = "caloClusters";
587  }
588  else {
589  clusName = "";
590  }
591  }
592 
593  if (mcPartName == "usedefault") {
594  mcPartName = "mcparticles"; // Always needs AliAODMCParticle objects
595  }
596 
597  TString name("AliAnalysisTaskDmesonJetsDetectorResponse");
598  if (!suffix.IsNull()) {
599  name += TString::Format("_%s", suffix.Data());
600  }
601 
603 
604  if (!trackName.IsNull()) {
605  AliHFTrackContainer* trackCont = new AliHFTrackContainer(trackName);
606  jetTask->AdoptParticleContainer(trackCont);
607  }
608 
609  if (!mcPartName.IsNull()) {
610  AliMCParticleContainer* partCont = new AliHFAODMCParticleContainer(mcPartName);
611  partCont->SetEtaLimits(-1.5, 1.5);
612  partCont->SetPtLimits(0, 1000);
613  jetTask->AdoptParticleContainer(partCont);
614  }
615 
616  jetTask->AddClusterContainer(clusName.Data());
617 
618  // Final settings, pass to manager and set the containers
619  mgr->AddTask(jetTask);
620 
621  // Create containers for input/output
622  AliAnalysisDataContainer* cinput1 = mgr->GetCommonInputContainer();
623  TString contname1(name);
624  contname1 += "_histos";
625  AliAnalysisDataContainer* coutput1 = mgr->CreateContainer(contname1.Data(),
626  TList::Class(), AliAnalysisManager::kOutputContainer,
627  Form("%s", AliAnalysisManager::GetCommonFileName()));
628 
629  mgr->ConnectInput(jetTask, 0, cinput1);
630  mgr->ConnectOutput(jetTask, 1, coutput1);
631 
632  for (Int_t i = 0; i < nMaxTrees; i++) {
633  TString contname = TString::Format("%s_tree_%d", name.Data(), i);
634  AliAnalysisDataContainer *coutput = mgr->CreateContainer(contname.Data(),
635  TTree::Class(),AliAnalysisManager::kOutputContainer,
636  Form("%s", AliAnalysisManager::GetCommonFileName()));
637  mgr->ConnectOutput(jetTask, 2+i, coutput);
638  }
639  return jetTask;
640 }
641 
std::list< AnalysisEngine > fAnalysisEngines
Array of analysis parameters.
static AliAnalysisTaskDmesonJetsDetectorResponse * AddTaskDmesonJetsDetectorResponse(TString trackName, TString clusName, TString mcPartName, Int_t nMaxTrees, TString suffix)
Analysis task for D meson jets.
Lightweight class that encapsulates matching between reconstructed and generated D0 mesons...
Bool_t fInhibit
!inhibit execution of the task
virtual void SetGenerated(const AliDmesonJetInfo &truth)=0
EDataType_t
Switch for the data type.
Analysis task used to build a detector response for D meson jets.
virtual void Set(const AliDmesonJetInfo &source)
Lightweight class that encapsulates D meson jets.
AliHFAODMCParticleContainer * fMCContainer
! MC particle container
Analysis engine to produce detector response matrix in the D meson jet analysis.
virtual void UserCreateOutputObjects()
Creates the output containers.
void FillPartonLevelHistograms()
Fill histograms with parton-level information.
void ana(anaModes mode=mLocal)
Definition: ana.C:986
std::vector< AliAnalysisTaskDmesonJets::AliHFJetDefinition > & GetJetDefinitions()
virtual void UserCreateOutputObjects()
Creates the output containers.
AliAnalysisTaskDmesonJetsDetectorResponse()
This is the default constructor, used for ROOT I/O purposes.
AnalysisEngine * fRecontructed
! Reconstructed level analysis engine
Lightweight class that encapsulates matching between reconstructed and generated D mesons...
Select tracks based on specific prescriptions of HF analysis.
friend bool operator==(const ResponseEngine &lhs, const ResponseEngine &rhs)
virtual void Set(const AliDmesonJetInfo &source)
int Int_t
Definition: External.C:63
unsigned int UInt_t
Definition: External.C:33
THashList * GetListOfHistograms() const
Get the list of histograms.
Definition: THistManager.h:671
virtual void Reset()
Reset the current object.
Lightweight class that encapsulates D meson jets for PbPb analysis.
TH1 * CreateTH1(const char *name, const char *title, int nbins, double xmin, double xmax, Option_t *opt="")
Create a new TH1 within the container.
Double_t fMaxJetDmesonDistance
Maximum distance between a generated D meson and a reconstructed jet, used for geometrical matching (...
THistManager fHistManagerResponse
Histogram manager for response.
void FillTH1(const char *hname, double x, double weight=1., Option_t *opt="")
Fill a 1D histogram within the container.
virtual void SetReconstructed(const AliDmesonJetInfo &reco)=0
Select MC particles based on specific prescriptions of HF analysis.
Lightweight class that encapsulates matching between reconstructed and generated D* mesons...
std::map< int, AliDmesonJetInfo > & GetDmesons()
Bool_t fApplyKinematicCuts
Apply jet kinematic cuts.
std::vector< ResponseEngine > fResponseEngines
! Response engines
EOutputType_t fOutputType
Output type: none, TTree or THnSparse.
friend bool operator<(const ResponseEngine &lhs, const ResponseEngine &rhs)
bool Bool_t
Definition: External.C:53
Int_t fNOutputTrees
Maximum number of output trees.
Container for MC-true particles within the EMCAL framework.
Definition: External.C:196
void RunAnalysis()
Run the requested analysis for the current event.
virtual void Set(const AliDmesonJetInfo &source, std::string n)