AliPhysics  master (3d17d9d)
AliAnalysisTaskEmcalJetHCorrelations.cxx
Go to the documentation of this file.
1 //Measure Jet-hadron correlations
3 //Does event Mixing using AliEventPoolManager
5 
7 
8 #include <bitset>
9 
10 #include <TH1F.h>
11 #include <TH2F.h>
12 #include <TH3F.h>
13 #include <THnSparse.h>
14 #include <TVector3.h>
15 #include <TFile.h>
16 #include <TGrid.h>
17 #include <TList.h>
18 
19 #include <AliAnalysisManager.h>
20 #include <AliInputEventHandler.h>
21 #include <AliEventPoolManager.h>
22 #include <AliLog.h>
23 #include <AliVAODHeader.h>
24 #include <AliVTrack.h>
25 
26 #include "AliEmcalJet.h"
27 #include "AliTLorentzVector.h"
28 #include "AliBasicParticle.h"
29 #include "AliEmcalContainerUtils.h"
30 #include "AliClusterContainer.h"
31 #include "AliTrackContainer.h"
32 #include "AliJetContainer.h"
35 
37 
38 namespace PWGJE {
39 namespace EMCALJetTasks {
40 
45  AliAnalysisTaskEmcalJet("AliAnalysisTaskEmcalJetHCorrelations", kFALSE),
46  fYAMLConfig(),
47  fConfigurationInitialized(false),
48  fTrackBias(5),
49  fClusterBias(5),
50  fDoEventMixing(kFALSE),
51  fNMixingTracks(50000), fMinNTracksMixedEvents(5000), fMinNEventsMixedEvents(5), fNCentBinsMixedEvent(10),
52  fPoolMgr(nullptr),
53  fTriggerType(AliVEvent::kEMCEJE), fMixingEventType(AliVEvent::kMB | AliVEvent::kCentral | AliVEvent::kSemiCentral),
54  fDisableFastPartition(kFALSE),
55  fRandom(0),
56  fEfficiencyPeriodIdentifier(AliAnalysisTaskEmcalJetHUtils::kDisableEff),
57  fArtificialTrackInefficiency(1.0),
58  fNoMixedEventJESCorrection(kFALSE),
59  fJESCorrectionHist(nullptr),
60  fDoLessSparseAxes(kFALSE), fDoWiderTrackBin(kFALSE),
61  fRequireMatchedJetWhenEmbedding(kTRUE),
62  fMinSharedMomentumFraction(0.),
63  fRequireMatchedPartLevelJet(false),
64  fMaxMatchedJetDistance(-1),
65  fHistManager(),
66  fHistJetHTrackPt(nullptr),
67  fHistJetEtaPhi(nullptr),
68  fHistJetHEtaPhi(nullptr),
69  fhnMixedEvents(nullptr),
70  fhnJH(nullptr),
71  fhnTrigger(nullptr)
72 {
73  // Default Constructor
75 }
76 
81  AliAnalysisTaskEmcalJet(name, kTRUE),
82  fYAMLConfig(),
84  fTrackBias(5),
85  fClusterBias(5),
86  fDoEventMixing(kFALSE),
89  fTriggerType(AliVEvent::kEMCEJE), fMixingEventType(AliVEvent::kMB | AliVEvent::kCentral | AliVEvent::kSemiCentral),
90  fDisableFastPartition(kFALSE),
91  fRandom(0),
96  fDoLessSparseAxes(kFALSE), fDoWiderTrackBin(kFALSE),
101  fHistManager(name),
106  fhnJH(nullptr),
108 {
109  // Constructor
111  // Ensure that additional general histograms are created
113 }
114 
119 {
120  for(Int_t trackPtBin = 0; trackPtBin < kMaxTrackPtBins; trackPtBin++){
121  fHistTrackEtaPhi[trackPtBin] = nullptr;
122  }
123  for(Int_t centralityBin = 0; centralityBin < kMaxCentralityBins; ++centralityBin){
124  fHistJetPt[centralityBin] = nullptr;
125  fHistJetPtBias[centralityBin] = nullptr;
126  }
127 }
128 
133 {
135 
136  // Ensure that we have at least one configuration in the YAML config.
137  if (fYAMLConfig.DoesConfigurationExist(0) == false) {
138  // No configurations exist. Return immediately.
140  }
141 
142  // Always initialize for streaming purposes
144 
145  // Setup task based on the properties defined in the YAML config
146  AliDebugStream(2) << "Configuring task from the YAML configuration.\n";
148  AliDebugStream(2) << "Finished configuring via the YAML configuration.\n";
149 
150  // Print the results of the initialization
151  // Print outside of the ALICE Log system to ensure that it is always available!
152  std::cout << *this;
153 
156 }
157 
162 {
163  // Base class options
164  // Task physics (trigger) selection.
165  std::string baseName = "eventCuts";
166  std::vector<std::string> physicsSelection;
167  bool res = fYAMLConfig.GetProperty(std::vector<std::string>({"eventCuts", "physicsSelection"}), physicsSelection, false);
168  if (res) {
169  fOfflineTriggerMask = AliEmcalContainerUtils::DeterminePhysicsSelectionFromYAML(physicsSelection);
170  }
171 
172  // Event cuts
173  // If event cuts are enabled (which they exceptionally are by default), then we want to configure them here.
174  // If the event cuts are explicitly disabled, then we invert that value to enable the AliAnylsisTaskEmcal
175  // builtin event selection.
176  bool tempBool;
177  fYAMLConfig.GetProperty({baseName, "enabled"}, tempBool, false);
178  fUseBuiltinEventSelection = !tempBool;
179  if (fUseBuiltinEventSelection == false) {
180  // Need to include the namespace so that AliDebug will work properly...
181  std::string taskName = "PWGJE::EMCALJetTasks::";
182  taskName += GetName();
183  AliAnalysisTaskEmcalJetHUtils::ConfigureEventCuts(fAliEventCuts, fYAMLConfig, fOfflineTriggerMask, baseName, taskName);
184  }
185 
186  // General task options
187  baseName = "general";
188  fYAMLConfig.GetProperty({baseName, "nCentBins"}, fNcentBins, false);
189 
190  // Efficiency
191  std::string tempStr = "";
192  baseName = "efficiency";
193  res = fYAMLConfig.GetProperty({baseName, "periodIdentifier"}, tempStr, false);
194  if (res) {
196  }
197 }
198 
203 {
204  // Called once
206 
207  // Check that the task was initialized
209  AliFatal("Task was not initialized. Please ensure that Initialize() was called!");
210  }
211  // Reinitialize the YAML configuration
213 
214  // Create histograms
215  fHistJetHTrackPt = new TH1F("fHistJetHTrackPt", "P_{T} distribution", 1000, 0.0, 100.0);
216  fHistJetEtaPhi = new TH2F("fHistJetEtaPhi","Jet eta-phi",900,-1.8,1.8,720,-3.2,3.2);
217  fHistJetHEtaPhi = new TH2F("fHistJetHEtaPhi","Jet-Hadron deta-dphi",900,-1.8,1.8,720,-1.6,4.8);
218 
220  fOutput->Add(fHistJetEtaPhi);
221  fOutput->Add(fHistJetHEtaPhi);
222 
223  TString name;
224  for(Int_t trackPtBin = 0; trackPtBin < kMaxTrackPtBins; ++trackPtBin){
225  name = Form("fHistTrackEtaPhi_%i", trackPtBin);
226  fHistTrackEtaPhi[trackPtBin] = new TH2F(name,name,400,-1,1,720,0.0,2.0*TMath::Pi());
227  fOutput->Add(fHistTrackEtaPhi[trackPtBin]);
228  }
229 
230  for(Int_t centralityBin = 0; centralityBin < kMaxCentralityBins; ++centralityBin){
231  name = Form("fHistJetPt_%i",centralityBin);
232  fHistJetPt[centralityBin] = new TH1F(name,name,240,-40,200);
233  fOutput->Add(fHistJetPt[centralityBin]);
234 
235  name = Form("fHistJetPtBias_%i",centralityBin);
236  fHistJetPtBias[centralityBin] = new TH1F(name,name,240,-40,200);
237  fOutput->Add(fHistJetPtBias[centralityBin]);
238  }
239 
240  // Jet matching cuts
241  // Only need if we actually jet matching
243  std::vector<std::string> binLabels = {"noMatch", "matchedJet", "sharedMomentumFraction", "partLevelMatchedJet", "jetDistance", "passedAllCuts"};
244  for (auto histName : std::vector<std::string>({"SameEvent", "MixedEvent"})) {
245  name = std::string("fHistJetMatching") + histName.c_str() + "Cuts";
246  std::string title = std::string("Jets which passed matching jet cuts for ") + histName;
247  auto histMatchedJetCuts = fHistManager.CreateTH1(name, title.c_str(), binLabels.size(), 0, binLabels.size());
248  // Set label names
249  for (unsigned int i = 1; i <= binLabels.size(); i++) {
250  histMatchedJetCuts->GetXaxis()->SetBinLabel(i, binLabels.at(i-1).c_str());
251  }
252  histMatchedJetCuts->GetYaxis()->SetTitle("Number of jets");
253  }
254  }
255 
256  // NOTE: The bit encoding doesn't preserve the order defined here. It's
257  // just using the bit values.
258  UInt_t cifras = 0; // bit coded, see GetDimParams() below
259  if(fDoLessSparseAxes) {
260  cifras = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<7;
261  } else {
262  cifras = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<7 | 1<<8 | 1<<9;
263  }
264  fhnJH = NewTHnSparseF("fhnJH", cifras);
265  fhnJH->Sumw2();
266  fOutput->Add(fhnJH);
267 
268  if(fDoEventMixing){
269  // The event plane angle does not need to be included because the semi-central determined that the EP angle didn't change
270  // significantly for any of the EP orientations. However, it will be included so this can be demonstrated for the central
271  // analysis if so desired.
272  if(fDoLessSparseAxes) {
273  cifras = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<7;
274  } else {
275  cifras = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<7 | 1<<8 | 1<<9;
276  }
277  fhnMixedEvents = NewTHnSparseF("fhnMixedEvents", cifras);
278  fhnMixedEvents->Sumw2();
279  fOutput->Add(fhnMixedEvents);
280  }
281 
282  // Trigger THnSparse
283  cifras = 1<<0 | 1<<1 | 1<<7;
284  fhnTrigger = NewTHnSparseF("fhnTrigger", cifras);
285  fhnTrigger->Sumw2();
286  fOutput->Add(fhnTrigger);
287 
288  // Store hist manager output in the output list
289  TIter next(fHistManager.GetListOfHistograms());
290  TObject* obj = 0;
291  while ((obj = next())) {
292  fOutput->Add(obj);
293  }
294 
295  PostData(1, fOutput);
296 
297  // Event Mixing
298  Int_t poolSize = -1; // Maximum number of events. Set to -1 to avoid limits on number of events
299  // ZVertex
300  Int_t nZVertexBins = 10;
301  Double_t* zVertexBins = GenerateFixedBinArray(nZVertexBins, -10, 10);
302  // Event activity (centrality of multiplicity)
303  Int_t nEventActivityBins = 8;
304  Double_t* eventActivityBins = 0;
305  // +1 to accomodate the fact that we define bins rather than array entries.
306  Double_t multiplicityBins[kMixedEventMultiplicityBins+1] = {0., 4., 9., 15., 25., 35., 55., 100., 500.};
307 
308  // Cannot use GetBeamType() since it is not available until UserExec()
309  if (fForceBeamType != AliAnalysisTaskEmcal::kpp ) { //all besides pp
310  // Event Activity is centrality in AA, pA
311  nEventActivityBins = fNCentBinsMixedEvent;
312  eventActivityBins = GenerateFixedBinArray(nEventActivityBins, 0, 100);
313  }
314  else if (fForceBeamType == AliAnalysisTaskEmcal::kpp) { //for pp only
315  // Event Activity is multiplicity in pp
316  eventActivityBins = multiplicityBins;
317  }
318 
319  fPoolMgr = new AliEventPoolManager(poolSize, fNMixingTracks, nEventActivityBins, eventActivityBins, nZVertexBins, zVertexBins);
320 
321  // Print pool properties
322  fPoolMgr->Validate();
323 }
324 
329 {
330  // Base class.
332 
333  // Ensure that the random number generator is seeded in each job.
334  fRandom.SetSeed(0);
335 }
336 
344 {
345  Int_t ptBin = -1;
346  if (pt < 0.5) ptBin = 0;
347  else if (pt < 1 ) ptBin = 1;
348  else if (pt < 2 ) ptBin = 2;
349  else if (pt < 3 ) ptBin = 3;
350  else if (pt < 5 ) ptBin = 4;
351  else if (pt < 8 ) ptBin = 5;
352  else if (pt < 20 ) ptBin = 6;
353 
354  return ptBin;
355 }
356 
363 {
364  UInt_t eventTrigger = 0;
365  if (fIsEmbedded) {
366  auto embeddingHelper = AliAnalysisTaskEmcalEmbeddingHelper::GetInstance();
367  if (embeddingHelper) {
368  auto aodHeader = dynamic_cast<AliVAODHeader *>(embeddingHelper->GetEventHeader());
369  if (aodHeader) {
370  AliDebugStream(5) << "Retrieving trigger mask from embedded event\n";
371  eventTrigger = aodHeader->GetOfflineTrigger();
372  }
373  else {
374  AliErrorStream() << "Failed to retrieve requested AOD header from embedding helper\n";
375  }
376  }
377  else {
378  AliErrorStream() << "Failed to retrieve requested embedding helper\n";
379  }
380  }
381  else {
382  AliDebugStream(5) << "Retrieving trigger mask from internal event\n";
383  eventTrigger = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
384  }
385 
386  return eventTrigger;
387 }
388 
393 {
394  // NOTE: Clusters are never used directly in the task, so the container is neither created not retrieved
395  // Retrieve tracks
396  AliTrackContainer * tracks = static_cast<AliTrackContainer * >(GetParticleContainer("tracksForCorrelations"));
397  if (!tracks) {
398  AliError(Form("%s: Unable to retrieve tracks!", GetName()));
399  return kFALSE;
400  }
401 
402  // Retrieve jets
403  AliJetContainer * jets = GetJetContainer(0);
404  if (!jets) {
405  AliError(Form("%s: Unable to retrieve jets!", GetName()));
406  return kFALSE;
407  }
408 
409  // Keep track of the tracks which are rejected with an aritificial track inefficiency
410  std::vector<unsigned int> rejectedTrackIndices;
411  bool useListOfRejectedIndices = false;
412 
413  // Get z vertex
414  Double_t zVertex=fVertex[2];
415  // Flags
416  Bool_t biasedJet = kFALSE;
417  Bool_t leadJet = kFALSE;
418  // Jet pt
419  Double_t jetPt = 0;
420  // Rho
421  // NOTE: Defaults to 0 if rho was not specified.
422  Double_t rhoVal = jets->GetRhoVal();
423  // Relative angles and distances
424  Double_t deltaPhi = 0;
425  Double_t deltaEta = 0;
426  Double_t deltaR = 0;
427  Double_t epAngle = 0;
428  // Event activity (centrality or multipilicity)
429  Double_t eventActivity = 0;
430  // Efficiency correction
431  Double_t efficiency = -999;
432  // For comparison to the current jet
433  AliEmcalJet * leadingJet = jets->GetLeadingJet();
434  // For getting the proper properties of tracks
435  AliTLorentzVector track;
436 
437  // Determine the trigger for the current event
438  UInt_t eventTrigger = RetrieveTriggerMask();
439 
440  AliDebugStream(5) << "Beginning main processing. Number of jets: " << jets->GetNJets() << ", accepted jets: " << jets->GetNAcceptedJets() << "\n";
441 
442  // Handle fast partition if selected
443  if ((eventTrigger & AliVEvent::kFastOnly) && fDisableFastPartition) {
444  AliDebugStream(4) << GetName() << ": Fast partition disabled\n";
445  if (fGeneralHistograms) {
446  fHistEventRejection->Fill("Fast Partition", 1);
447  }
448  return kFALSE;
449  }
450 
451  for (auto jet : jets->accepted()) {
452  // Selects only events that we are interested in (ie triggered)
453  if (!(eventTrigger & fTriggerType)) {
454  // The two bitwise and to zero yet are still equal when both are 0, so we allow for that possibility
455  if (eventTrigger == fTriggerType && eventTrigger == 0) {
456  AliDebugStream(5) << "Event accepted because the physics selection is \"0\".\n";
457  }
458  else {
459  AliDebugStream(5) << "Rejected jets due to physics selection. Phys sel: " << std::bitset<32>(eventTrigger) << ", requested triggers: " << std::bitset<32>(fTriggerType) << " \n";
460  // We can break here - the physics selection is not going to change within an event.
461  break;
462  }
463  }
464 
465  AliDebugStream(5) << "Jet passed event selection!\nJet: " << jet->toString().Data() << "\n";
466 
467  // Require the found jet to be matched
468  // This match should be between detector and particle level MC
470  bool foundMatchedJet = CheckForMatchedJet(jets, jet, "fHistJetMatchingSameEventCuts");
471  if (foundMatchedJet == false) {
472  continue;
473  }
474  }
475 
476  // Determine event activity
477  if (fBeamType == kAA || fBeamType == kpA) {
478  eventActivity = fCent;
479  }
480  else if (fBeamType == kpp) {
481  eventActivity = static_cast<Double_t>(tracks->GetNTracks());
482  }
483 
484  // Jet properties
485  jetPt = AliAnalysisTaskEmcalJetHUtils::GetJetPt(jet, rhoVal);
486  // Determine if we have the lead jet
487  leadJet = kFALSE;
488  if (jet == leadingJet) leadJet = kTRUE;
489  biasedJet = BiasedJet(jet);
491 
492  // Fill jet properties
493  fHistJetEtaPhi->Fill(jet->Eta(), jet->Phi());
494  FillHist(fHistJetPt[fCentBin], jetPt);
495  if (biasedJet == kTRUE) {
496  FillHist(fHistJetPtBias[fCentBin], jetPt);
497 
498  const double triggerInfo[] = {eventActivity, jetPt, epAngle};
499  fhnTrigger->Fill(triggerInfo);
500  }
501 
502  // Cut on jet pt of 15 to reduce the size of the sparses
503  if (jetPt > 15) {
504 
505  AliDebugStream(4) << "Passed min jet pt cut of 15. Jet: " << jet->toString().Data() << "\n";
506  AliDebugStream(4) << "N accepted tracks: " << tracks->GetNAcceptedTracks() << "\n";
507  auto tracksIter = tracks->accepted_momentum();
508  for (auto trackIter = tracksIter.begin(); trackIter != tracksIter.end(); trackIter++ ) {
509  // Get proper track properties
510  track.Clear();
511  track = trackIter->first;
512 
513  // Artificial inefficiency
514  // Note that we already randomly rejected tracks so that the same tracks will be rejected for the mixed events
515  bool rejectParticle = CheckArtificialTrackEfficiency(trackIter.current_index(), rejectedTrackIndices, useListOfRejectedIndices);
516  if (rejectParticle) {
517  AliDebugStream(4) << "Track rejected in signal correlation loop.\n";
518  continue;
519  }
520 
521  // Determine relative angles and distances and set the respective variables
522  GetDeltaEtaDeltaPhiDeltaR(track, jet, deltaEta, deltaPhi, deltaR);
523 
524  // Fill track properties
525  fHistJetHTrackPt->Fill(track.Pt());
526  fHistJetHEtaPhi->Fill(deltaEta, deltaPhi);
527 
528  // Calculate single particle tracking efficiency for correlations
529  efficiency = EffCorrection(track.Eta(), track.Pt());
530  AliDebugStream(6) << "track eta: " << track.Eta() << ", track pt: " << track.Pt() << ", efficiency: " << efficiency << "\n";
531 
532  if (biasedJet == kTRUE) {
533  if(fDoLessSparseAxes) { // check if we want all dimensions
534  double triggerEntries[] = {eventActivity, jetPt, track.Pt(), deltaEta, deltaPhi, static_cast<Double_t>(leadJet), epAngle};
535  FillHist(fhnJH, triggerEntries, 1.0/efficiency);
536  } else {
537  double triggerEntries[] = {eventActivity, jetPt, track.Pt(), deltaEta, deltaPhi, static_cast<Double_t>(leadJet), epAngle, zVertex, deltaR};
538  FillHist(fhnJH, triggerEntries, 1.0/efficiency);
539  }
540  }
541 
542  } //track loop
543 
544  // After one jet (and looping over whatever tracks are available in this event), we want to use the list of rejected indices,
545  // both for the next possible signal jet in the event and for the mixed events
546  AliDebugStream(4) << "Switching to list of rejected track indices. Number of indices: " << rejectedTrackIndices.size() << "\n";
547  useListOfRejectedIndices = true;
548 
549  }//jet pt cut
550  }//jet loop
551 
552  //Prepare to do event mixing
553 
554  // create a list of reduced objects. This speeds up processing and reduces memory consumption for the event pool
555  TObjArray* tracksClone = 0;
556 
557  if(fDoEventMixing == kTRUE){
558 
559  // event mixing
560 
561  // 1. First get an event pool corresponding in mult (cent) and
562  // zvertex to the current event. Once initialized, the pool
563  // should contain nMix (reduced) events. This routine does not
564  // pre-scan the chain. The first several events of every chain
565  // will be skipped until the needed pools are filled to the
566  // specified depth. If the pool categories are not too rare, this
567  // should not be a problem. If they are rare, you could lose
568  // statistics.
569 
570  // 2. Collect the whole pool's content of tracks into one TObjArray
571  // (bgTracks), which is effectively a single background super-event.
572 
573  // 3. The reduced and bgTracks arrays must both be passed into
574  // FillCorrelations(). Also nMix should be passed in, so a weight
575  // of 1./nMix can be applied.
576 
577  AliEventPool *pool = 0;
578  if (fBeamType == kAA || fBeamType == kpA) {//everything but pp
579  pool = fPoolMgr->GetEventPool(fCent, zVertex);
580  }
581  else if (fBeamType == kpp) {//pp only
582  pool = fPoolMgr->GetEventPool(static_cast<Double_t>(tracks->GetNTracks()), zVertex);
583  }
584 
585  if (!pool){
586  if (fBeamType == kAA || fBeamType == kpA) AliFatal(Form("No pool found for centrality = %f, zVertex = %f", fCent, zVertex));
587  else if (fBeamType == kpp) AliFatal(Form("No pool found for ntracks_pp = %d, zVertex = %f", tracks->GetNTracks(), zVertex));
588  return kTRUE;
589  }
590 
591  // The number of events in the pool
592  Int_t nMix = pool->GetCurrentNEvents();
593 
594  // The two bitwise and to zero yet are still equal when both are 0, so we allow for that possibility
595  if((eventTrigger & fTriggerType) || eventTrigger == fTriggerType) {
596  // check for a trigger jet
597  if (pool->IsReady() || pool->NTracksInPool() >= fMinNTracksMixedEvents || nMix >= fMinNEventsMixedEvents) {
598 
599  for (auto jet : jets->accepted()) {
600  // Require the found jet to be matched
601  // This match should be between detector and particle level MC
603  bool foundMatchedJet = CheckForMatchedJet(jets, jet, "fHistJetMatchingMixedEventCuts");
604  if (foundMatchedJet == false) {
605  continue;
606  }
607  }
608 
609  if (fBeamType == kAA || fBeamType == kpA) { //pA and AA
610  eventActivity = fCent;
611  }
612  else if (fBeamType == kpp) {
613  eventActivity = static_cast<Double_t>(tracks->GetNTracks());
614  }
615 
616  // Jet properties
617  jetPt = AliAnalysisTaskEmcalJetHUtils::GetJetPt(jet, rhoVal);
618  // Determine if we have the lead jet
619  leadJet = kFALSE;
620  if (jet == leadingJet) { leadJet = kTRUE; }
621  biasedJet = BiasedJet(jet);
623 
624  // Make sure event contains a biased jet above our threshold (reduce stats of sparse)
625  if (jetPt < 15 || biasedJet == kFALSE) continue;
626 
627  // Fill mixed-event histos here
628  for (Int_t jMix=0; jMix < nMix; jMix++) {
629  TObjArray* bgTracks = pool->GetEvent(jMix);
630 
631  for (Int_t ibg=0; ibg < bgTracks->GetEntries(); ibg++){
632  AliBasicParticle *bgTrack = static_cast<AliBasicParticle*>(bgTracks->At(ibg));
633  if(!bgTrack) {
634  AliError(Form("%s:Failed to retrieve tracks from mixed events", GetName()));
635  }
636 
637  // NOTE: We don't need to apply the artificial track inefficiency here because we already applied
638  // it when will filling into the event pool (in CloneAndReduceTrackList()).
639 
640  // Fill into TLorentzVector for use with functions below
641  track.Clear();
642  track.SetPtEtaPhiE(bgTrack->Pt(), bgTrack->Eta(), bgTrack->Phi(), 0);
643 
644  // Calculate single particle tracking efficiency of mixed events for correlations
645  efficiency = EffCorrection(track.Eta(), track.Pt());
646 
647  // Phi is [-0.5*TMath::Pi(), 3*TMath::Pi()/2.]
648  GetDeltaEtaDeltaPhiDeltaR(track, jet, deltaEta, deltaPhi, deltaR);
649 
650  if (fDoLessSparseAxes) { // check if we want all the axis filled
651  double triggerEntries[] = {eventActivity, jetPt, track.Pt(), deltaEta, deltaPhi, static_cast<Double_t>(leadJet), epAngle};
652  FillHist(fhnMixedEvents, triggerEntries, 1./(nMix*efficiency), fNoMixedEventJESCorrection);
653  } else {
654  double triggerEntries[] = {eventActivity, jetPt, track.Pt(), deltaEta, deltaPhi, static_cast<Double_t>(leadJet), epAngle, zVertex, deltaR};
655  FillHist(fhnMixedEvents, triggerEntries, 1./(nMix*efficiency), fNoMixedEventJESCorrection);
656  }
657  }
658  }
659  }
660  }
661  }
662 
663  // The two bitwise and to zero yet are still equal when both are 0, so we allow for that possibility
664  if ((eventTrigger & fMixingEventType) || eventTrigger == fMixingEventType) {
665  tracksClone = CloneAndReduceTrackList(rejectedTrackIndices, useListOfRejectedIndices);
666 
667  //update pool if jet in event or not
668  pool->UpdatePool(tracksClone);
669  }
670 
671  } // end of event mixing
672 
673  return kTRUE;
674 }
675 
683 {
684  if ((jet->MaxTrackPt() > fTrackBias) || (jet->MaxClusterPt() > fClusterBias))
685  {
686  return kTRUE;
687  }
688  return kFALSE;
689 }
690 
700 void AliAnalysisTaskEmcalJetHCorrelations::GetDeltaEtaDeltaPhiDeltaR(AliTLorentzVector & particleOne, AliVParticle * particleTwo, Double_t & deltaEta, Double_t & deltaPhi, Double_t & deltaR)
701 {
702  // Define dPhi = jet.Phi() - particle.Phi() and similarly for dEta
703  // Returns deltaPhi in symmetric range so that we can calculate DeltaR.
704  deltaPhi = DeltaPhi(particleOne.Phi(), particleTwo->Phi(), -1.0*TMath::Pi(), TMath::Pi());
705  deltaEta = particleTwo->Eta() - particleOne.Eta();
706  deltaR = TMath::Sqrt(deltaPhi*deltaPhi + deltaEta*deltaEta);
707 
708  // Adjust to the normal range after the DeltaR caluclation
709  deltaPhi = DeltaPhi(particleTwo->Phi(), particleOne.Phi(), -0.5*TMath::Pi(), 3*TMath::Pi()/2.);
710 }
711 
720 bool AliAnalysisTaskEmcalJetHCorrelations::CheckArtificialTrackEfficiency(unsigned int trackIndex, std::vector<unsigned int> & rejectedTrackIndices, bool useRejectedList)
721 {
722  bool returnValue = false;
723  if (fArtificialTrackInefficiency < 1.) {
724  if (useRejectedList) {
725  if (std::find(rejectedTrackIndices.begin(), rejectedTrackIndices.end(), trackIndex) != rejectedTrackIndices.end()) {
726  AliDebugStream(4) << "Track " << trackIndex << " rejected due to artificial tracking inefficiency (from list)\n";
727  returnValue = true;
728  }
729  }
730  else {
731  // Rejet randomly
732  Double_t rnd = fRandom.Rndm();
733  if (fArtificialTrackInefficiency < rnd) {
734  // Store index so we can reject it again if it is also filled for mixed events
735  rejectedTrackIndices.push_back(trackIndex);
736  AliDebugStream(4) << "Track " << trackIndex << " rejected due to artificial tracking inefficiency (from random)\n";
737  returnValue = true;
738  }
739  }
740  }
741 
742  return returnValue;
743 }
744 
764 {
765  bool returnValue = false;
766  if (jet->ClosestJet()) {
767  fHistManager.FillTH1(histName.c_str(), "matchedJet");
768  returnValue = true;
769  // TODO: Can it be merged with the function in JetHPerformance?
770  AliDebugStream(4) << "Jet is matched!\nJet: " << jet->toString() << "\n";
771  // Check shared momentum fraction
772  // We explicitly want to use indices instead of geometric matching
773  double sharedFraction = jets->GetFractionSharedPt(jet, nullptr);
774  if (sharedFraction < fMinSharedMomentumFraction) {
775  AliDebugStream(4) << "Jet rejected due to shared momentum fraction of " << sharedFraction << ", which is smaller than the min momentum fraction of " << fMinSharedMomentumFraction << "\n";
776  returnValue = false;
777  }
778  else {
779  AliDebugStream(4) << "Passed shared momentum fraction with value of " << sharedFraction << "\n";
780  fHistManager.FillTH1(histName.c_str(), "sharedMomentumFraction");
781  }
782 
784  AliEmcalJet * detLevelJet = jet->ClosestJet();
785  AliEmcalJet * partLevelJet = detLevelJet->ClosestJet();
786  if (!partLevelJet) {
787  AliDebugStream(4) << "Jet rejected due to no matching part level jet.\n";
788  returnValue = false;
789  }
790  else {
791  AliDebugStream(4) << "Det level jet has a required match to a part level jet.\n" << "Part level jet: " << partLevelJet->toString() << "\n";
792  fHistManager.FillTH1(histName.c_str(), "partLevelMatchedJet");
793  }
794  }
795 
796  // Only check matched jet distance if a value has been set
797  if (fMaxMatchedJetDistance > 0) {
798  double matchedJetDistance = jet->ClosestJetDistance();
799  if (matchedJetDistance > fMaxMatchedJetDistance) {
800  AliDebugStream(4) << "Jet rejected due to matching distance of " << matchedJetDistance << ", which is larger than the max distance of " << fMaxMatchedJetDistance << "\n";
801  returnValue = false;
802  }
803  else {
804  AliDebugStream(4) << "Jet passed distance cut with distance of " << matchedJetDistance << "\n";
805  fHistManager.FillTH1(histName.c_str(), "jetDistance");
806  }
807  }
808 
809  // Record all cuts passed
810  if (returnValue == true) {
811  fHistManager.FillTH1(histName.c_str(), "passedAllCuts");
812  }
813  }
814  else {
815  AliDebugStream(5) << "Rejected jet because it was not matched to a external event jet.\n";
816  fHistManager.FillTH1(histName.c_str(), "noMatch");
817  returnValue = false;
818  }
819 
820  return returnValue;
821 }
822 
830 THnSparse* AliAnalysisTaskEmcalJetHCorrelations::NewTHnSparseF(const char* name, UInt_t entries)
831 {
832  Int_t count = 0;
833  UInt_t tmp = entries;
834  while(tmp!=0){
835  count++;
836  tmp = tmp &~ -tmp; // clear lowest bit
837  }
838 
839  TString hnTitle(name);
840  const Int_t dim = count;
841  Int_t nbins[dim];
842  Double_t xmin[dim];
843  Double_t xmax[dim];
844 
845  Int_t i=0;
846  Int_t c=0;
847  while(c<dim && i<32){
848  if(entries&(1<<i)){
849 
850  TString label("");
851  GetDimParams(i, label, nbins[c], xmin[c], xmax[c]);
852  hnTitle += Form(";%s",label.Data());
853  c++;
854  }
855 
856  i++;
857  }
858  hnTitle += ";";
859 
860  return new THnSparseF(name, hnTitle.Data(), dim, nbins, xmin, xmax);
861 }
862 
873 {
874  const Double_t pi = TMath::Pi();
875 
876  switch(iEntry){
877 
878  case 0:
879  label = "V0 centrality (%)";
880  nbins = 10;
881  xmin = 0.;
882  xmax = 100.;
883  // Adjust for pp, since we are retrieving multiplicity instead
885  label = "Multiplicity";
886  xmax = 200.;
887  }
888  break;
889 
890  case 1:
891  label = "Jet p_{T}";
892  nbins = 24;
893  xmin = -40.;
894  xmax = 200.;
895  break;
896 
897  case 2:
898  if(fDoWiderTrackBin) {
899  label = "Track p_{T}";
900  nbins = 40;
901  xmin = 0.;
902  xmax = 10.;
903  } else {
904  label = "Track p_{T}";
905  nbins = 100;
906  xmin = 0.;
907  xmax = 10;
908  }
909  break;
910 
911  case 3:
912  label = "#Delta#eta";
913  nbins = 28;
914  xmin = -1.4;
915  xmax = 1.4;
916  break;
917 
918  case 4:
919  label = "#Delta#phi";
920  nbins = 72;
921  xmin = -0.5*pi;
922  xmax = 1.5*pi;
923  break;
924 
925  case 5:
926  label = "Leading Jet";
927  nbins = 3;
928  xmin = -0.5;
929  xmax = 2.5;
930  break;
931 
932  case 6:
933  label = "Trigger track";
934  nbins = 10;
935  xmin = 0;
936  xmax = 50;
937  break;
938 
939  case 7:
940  label = "Event plane angle";
941  nbins = 3;
942  xmin = 0;
943  xmax = TMath::Pi()/2.;
944  break;
945 
946  case 8:
947  label = "Z vertex (cm)";
948  nbins = 10;
949  xmin = -10;
950  xmax = 10;
951  break;
952 
953  case 9:
954  label = "deltaR";
955  nbins = 10;
956  xmin = 0.;
957  xmax = 5.0;
958  break;
959 
960  case 10:
961  label = "Leading track";
962  nbins = 20;
963  xmin = 0;
964  xmax = 50;
965  break;
966  }
967 }
968 
976 TObjArray* AliAnalysisTaskEmcalJetHCorrelations::CloneAndReduceTrackList(std::vector<unsigned int> & rejectedTrackIndices, const bool useRejectedList)
977 {
978  // clones a track list by using AliBasicTrack which uses much less memory (used for event mixing)
979  TObjArray* tracksClone = new TObjArray;
980  tracksClone->SetOwner(kTRUE);
981 
982  // Loop over all tracks
983  AliVParticle * particle = 0;
984  AliBasicParticle * clone = 0;
985  AliTrackContainer * tracks = GetTrackContainer("tracksForCorrelations");
986 
987  auto particlesIter = tracks->accepted_momentum();
988  for (auto particleIter = particlesIter.begin(); particleIter != particlesIter.end(); particleIter++)
989  {
990  // Retrieve the particle
991  particle = particleIter->second;
992 
993  // Artificial inefficiency
994  bool rejectParticle = CheckArtificialTrackEfficiency(particleIter.current_index(), rejectedTrackIndices, useRejectedList);
995  if (rejectParticle) {
996  AliDebugStream(4) << "Track rejected in CloneAndReduceTrackList()\n";
997  continue;
998  }
999 
1000  // Fill some QA information about the tracks
1001  Int_t trackPtBin = GetTrackPtBin(particle->Pt());
1002  if(trackPtBin > -1) fHistTrackEtaPhi[trackPtBin]->Fill(particle->Eta(),particle->Phi());
1003 
1004  // Create new particle
1005  clone = new AliBasicParticle(particle->Eta(), particle->Phi(), particle->Pt(), particle->Charge());
1006  // Set so that we can do comparisons using the IsEqual() function.
1007  clone->SetUniqueID(particle->GetUniqueID());
1008 
1009  tracksClone->Add(clone);
1010  }
1011 
1012  return tracksClone;
1013 }
1014 
1025  trackPt, trackEta, fCentBin, fEfficiencyPeriodIdentifier,
1026  "PWGJE::EMCALJetTasks::AliAnalysisTaskEmcalJetHCorrelations");
1027 }
1028 
1038 void AliAnalysisTaskEmcalJetHCorrelations::FillHist(TH1 * hist, Double_t fillValue, Double_t weight, Bool_t noCorrection)
1039 {
1040  if (fJESCorrectionHist == 0 || noCorrection == kTRUE)
1041  {
1042  AliDebugStream(3) << GetName() << ":" << hist->GetName() << ": " << std::boolalpha << "Using normal weights: JESHist: " << (fJESCorrectionHist ? fJESCorrectionHist->GetName() : "Null") << ", noCorrection: " << noCorrection << std::endl;
1043  hist->Fill(fillValue, weight);
1044  }
1045  else
1046  {
1047  // Determine where to get the values in the correction hist
1048  Int_t xBin = fJESCorrectionHist->GetXaxis()->FindBin(fillValue);
1049 
1050  std::vector <Double_t> yBinsContent;
1051  AliDebug(3, TString::Format("%s: Attempt to access weights from JES correction hist %s with jet pt %f!", GetName(), hist->GetName(), fillValue));
1052  AccessSetOfYBinValues(fJESCorrectionHist, xBin, yBinsContent);
1053  AliDebug(3, TString::Format("weights size: %zd", yBinsContent.size()));
1054 
1055  // Loop over all possible bins to contribute.
1056  // If content is 0 then calling Fill won't make a difference
1057  for (Int_t index = 1; index <= fJESCorrectionHist->GetYaxis()->GetNbins(); index++)
1058  {
1059  // Don't bother trying to fill in the weight is 0
1060  if (yBinsContent.at(index-1) > 0) {
1061  // Determine the value to fill based on the center of the bins.
1062  // This in principle allows the binning between the correction and hist to be different
1063  Double_t fillLocation = fJESCorrectionHist->GetYaxis()->GetBinCenter(index);
1064  AliDebug(4, TString::Format("fillLocation: %f, weight: %f", fillLocation, yBinsContent.at(index-1)));
1065  // minus 1 since loop starts at 1
1066  hist->Fill(fillLocation, weight*yBinsContent.at(index-1));
1067  }
1068  }
1069 
1070  //TEMP
1071  //hist->Draw();
1072  //END TEMP
1073  }
1074 }
1075 
1087 void AliAnalysisTaskEmcalJetHCorrelations::FillHist(THnSparse * hist, Double_t *fillValue, Double_t weight, Bool_t noCorrection)
1088 {
1089  if (fJESCorrectionHist == 0 || noCorrection == kTRUE)
1090  {
1091  AliDebugStream(3) << GetName() << ":" << hist->GetName() << ": " << std::boolalpha << "Using normal weights: JESHist: " << (fJESCorrectionHist ? fJESCorrectionHist->GetName() : "Null") << ", noCorrection: " << noCorrection << std::endl;
1092  hist->Fill(fillValue, weight);
1093  }
1094  else
1095  {
1096  // Jet pt is always located in the second position
1097  Double_t jetPt = fillValue[1];
1098 
1099  // Determine where to get the values in the correction hist
1100  Int_t xBin = fJESCorrectionHist->GetXaxis()->FindBin(jetPt);
1101 
1102  std::vector <Double_t> yBinsContent;
1103  AliDebug(3, TString::Format("%s: Attempt to access weights from JES correction hist %s with jet pt %f!", GetName(), hist->GetName(), jetPt));
1104  AccessSetOfYBinValues(fJESCorrectionHist, xBin, yBinsContent);
1105  AliDebug(3, TString::Format("weights size: %zd", yBinsContent.size()));
1106 
1107  // Loop over all possible bins to contribute.
1108  // If content is 0 then calling Fill won't make a difference
1109  for (Int_t index = 1; index <= fJESCorrectionHist->GetYaxis()->GetNbins(); index++)
1110  {
1111  // Don't bother trying to fill in the weight is 0
1112  if (yBinsContent.at(index-1) > 0) {
1113  // Determine the value to fill based on the center of the bins.
1114  // This in principle allows the binning between the correction and hist to be different
1115  fillValue[1] = fJESCorrectionHist->GetYaxis()->GetBinCenter(index);
1116  AliDebug(4,TString::Format("fillValue[1]: %f, weight: %f", fillValue[1], yBinsContent.at(index-1)));
1117  // minus 1 since loop starts at 1
1118  hist->Fill(fillValue, weight*yBinsContent.at(index-1));
1119  }
1120  }
1121  }
1122 }
1123 
1133 void AliAnalysisTaskEmcalJetHCorrelations::AccessSetOfYBinValues(TH2D * hist, Int_t xBin, std::vector <Double_t> & yBinsContent, Double_t scaleFactor)
1134 {
1135  for (Int_t index = 1; index <= hist->GetYaxis()->GetNbins(); index++)
1136  {
1137  //yBinsContent[index-1] = hist->GetBinContent(hist->GetBin(xBin,index));
1138  yBinsContent.push_back(hist->GetBinContent(hist->GetBin(xBin,index)));
1139 
1140  if (scaleFactor >= 0)
1141  {
1142  // -1 since index starts at 1
1143  hist->SetBinContent(hist->GetBin(xBin,index), yBinsContent.at(index-1)/scaleFactor);
1144  }
1145  }
1146 }
1147 
1165 {
1166  // Initialize grid connection if necessary
1167  if (filename.Contains("alien://") && !gGrid) {
1168  TGrid::Connect("alien://");
1169  }
1170 
1171  // Setup hist name if a track or cluster bias was defined.
1172  // NOTE: This can always be disabled by setting kDisableBias.
1173  // We arbitrarily add 0.1 to test since the values are doubles and cannot be
1174  // tested directly for equality. If we are still less than disable bins, then
1175  // it has been set and we should format it.
1176  // NOTE: To ensure we can disable, we don't just take the member values!
1177  // NOTE: The histBaseName will be attempted if the formatted name cannot be found.
1178  TString histBaseName = histName;
1180  histName = TString::Format("%s_Track%.2f", histName.Data(), trackBias);
1181  }
1182  if (clusterBias + 0.1 < AliAnalysisTaskEmcalJetHCorrelations::kDisableBias) {
1183  histName = TString::Format("%s_Clus%.2f", histName.Data(), clusterBias);
1184  }
1185 
1186  // Open file containing the correction
1187  TFile * jesCorrectionFile = TFile::Open(filename);
1188  if (!jesCorrectionFile || jesCorrectionFile->IsZombie()) {
1189  AliError(TString::Format("%s: Could not open JES correction file %s", GetName(), filename.Data()));
1190  return kFALSE;
1191  }
1192 
1193  // Retrieve the histogram containing the correction and safely add it to the task.
1194  TH2D * JESCorrectionHist = dynamic_cast<TH2D*>(jesCorrectionFile->Get(histName.Data()));
1195  if (JESCorrectionHist) {
1196  AliInfo(TString::Format("%s: JES correction hist name \"%s\" loaded from file %s.", GetName(), histName.Data(), filename.Data()));
1197  }
1198  else {
1199  AliError(TString::Format("%s: JES correction hist name \"%s\" not found in file %s.", GetName(), histName.Data(), filename.Data()));
1200 
1201  // Attempt the base name instead of the formatted hist name
1202  JESCorrectionHist = dynamic_cast<TH2D*>(jesCorrectionFile->Get(histBaseName.Data()));
1203  if (JESCorrectionHist) {
1204  AliInfo(TString::Format("%s: JES correction hist name \"%s\" loaded from file %s.", GetName(), histBaseName.Data(), filename.Data()));
1205  histName = histBaseName;
1206  }
1207  else
1208  {
1209  AliError(TString::Format("%s: JES correction with base hist name %s not found in file %s.", GetName(), histBaseName.Data(), filename.Data()));
1210  return kFALSE;
1211  }
1212  }
1213 
1214  // Clone to ensure that the hist is available
1215  TH2D * tempHist = static_cast<TH2D *>(JESCorrectionHist->Clone());
1216  tempHist->SetDirectory(0);
1217  SetJESCorrectionHist(tempHist);
1218 
1219  // Close file
1220  jesCorrectionFile->Close();
1221 
1222  // Append to task name for clarity
1223  // Unfortunately, this doesn't change the name of the output list (it would need to be
1224  // changed in the AnalysisManager output container), so the suffix is still important
1225  // if this correction is manually configured!
1226  TString tempName = GetName();
1227  TString tag = "_JESCorr";
1228  // Append the tag if it isn't already included
1229  if (tempName.Index(tag) == -1) {
1230  // Insert before the suffix
1231  Ssiz_t suffixLocation = tempName.Last('_');
1232  tempName.Insert(suffixLocation, tag.Data());
1233 
1234  // Set the new name
1235  AliDebug(3, TString::Format("%s: Setting task name to %s", GetName(), tempName.Data()));
1236  SetName(tempName.Data());
1237  }
1238 
1239  // Successful
1240  return kTRUE;
1241 }
1242 
1248  const char *nTracks,
1249  const char *nCaloClusters,
1250  // Jet options
1251  const Double_t trackBias,
1252  const Double_t clusterBias,
1253  // Mixed event options
1254  const Int_t nTracksMixedEvent, // Additionally acts as a switch for enabling mixed events
1255  const Int_t minNTracksMixedEvent,
1256  const Int_t minNEventsMixedEvent,
1257  const UInt_t nCentBinsMixedEvent,
1258  // Triggers
1259  UInt_t trigEvent,
1260  UInt_t mixEvent,
1261  // Options
1262  const Bool_t lessSparseAxes,
1263  const Bool_t widerTrackBin,
1264  // Corrections
1265  const Bool_t JESCorrection,
1266  const char * JESCorrectionFilename,
1267  const char * JESCorrectionHistName,
1268  const char *suffix
1269  )
1270 {
1271  // Get the pointer to the existing analysis manager via the static access method.
1272  //==============================================================================
1273  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
1274  if (!mgr)
1275  {
1276  AliErrorClass("No analysis manager to connect to.");
1277  return nullptr;
1278  }
1279 
1280  //-------------------------------------------------------
1281  // Init the task and do settings
1282  //-------------------------------------------------------
1283 
1284  // Determine cluster and track names
1285  TString trackName(nTracks);
1286  TString clusName(nCaloClusters);
1287 
1288  if (trackName == "usedefault") {
1290  }
1291 
1292  if (clusName == "usedefault") {
1294  }
1295 
1296  TString name("AliAnalysisTaskJetH");
1297  if (!trackName.IsNull()) {
1298  name += TString::Format("_%s", trackName.Data());
1299  }
1300  if (!clusName.IsNull()) {
1301  name += TString::Format("_%s", clusName.Data());
1302  }
1303  if (strcmp(suffix, "") != 0) {
1304  name += TString::Format("_%s", suffix);
1305  }
1306 
1308  // Set jet bias
1309  correlationTask->SetTrackBias(trackBias);
1310  correlationTask->SetClusterBias(clusterBias);
1311  // Mixed events
1312  correlationTask->SetEventMixing(static_cast<Bool_t>(nTracksMixedEvent));
1313  correlationTask->SetNumberOfMixingTracks(nTracksMixedEvent);
1314  correlationTask->SetMinNTracksForMixedEvents(minNTracksMixedEvent);
1315  correlationTask->SetMinNEventsForMixedEvents(minNEventsMixedEvent);
1316  correlationTask->SetNCentBinsMixedEvent(nCentBinsMixedEvent);
1317  // Triggers
1318  correlationTask->SetTriggerType(trigEvent);
1319  correlationTask->SetMixedEventTriggerType(mixEvent);
1320  // Options
1321  correlationTask->SetNCentBins(5);
1322  correlationTask->SetDoLessSparseAxes(lessSparseAxes);
1323  correlationTask->SetDoWiderTrackBin(widerTrackBin);
1324  // Corrections
1325  if (JESCorrection == kTRUE)
1326  {
1327  Bool_t result = correlationTask->RetrieveAndInitializeJESCorrectionHist(JESCorrectionFilename, JESCorrectionHistName, correlationTask->GetTrackBias(), correlationTask->GetClusterBias());
1328  if (!result) {
1329  AliErrorClass("Failed to successfully retrieve and initialize the JES correction! Task initialization continuing without JES correction (can be set manually later).");
1330  }
1331  }
1332 
1333  //-------------------------------------------------------
1334  // Final settings, pass to manager and set the containers
1335  //-------------------------------------------------------
1336 
1337  mgr->AddTask(correlationTask);
1338 
1339  // Create containers for input/output
1340  mgr->ConnectInput (correlationTask, 0, mgr->GetCommonInputContainer() );
1341  AliAnalysisDataContainer* cojeth =
1342  mgr->CreateContainer(correlationTask->GetName(), TList::Class(), AliAnalysisManager::kOutputContainer,
1343  Form("%s", AliAnalysisManager::GetCommonFileName()));
1344  mgr->ConnectOutput(correlationTask, 1, cojeth);
1345 
1346  return correlationTask;
1347 }
1348 
1355  std::string clusName,
1356  const double jetConstituentPtCut,
1357  const double trackEta,
1358  const double jetRadius)
1359 {
1360  bool returnValue = false;
1361  AliInfoStream() << "Configuring Jet-H Correlations task for a standard analysis.\n";
1362 
1363  // Add Containers
1364  // Clusters
1365  if (clusName == "usedefault") {
1367  }
1368  // For jet finding
1369  AliClusterContainer * clustersForJets = new AliClusterContainer(clusName.c_str());
1370  clustersForJets->SetName("clustersForJets");
1371  clustersForJets->SetMinE(jetConstituentPtCut);
1372 
1373  // Tracks
1374  // For jet finding
1375  if (trackName == "usedefault") {
1377  }
1379  particlesForJets->SetName("particlesForJets");
1380  particlesForJets->SetMinPt(jetConstituentPtCut);
1381  particlesForJets->SetEtaLimits(-1.0*trackEta, trackEta);
1382  // Don't need to adopt the container - we'll just use it to find the right jet collection
1383  // For correlations
1384  AliParticleContainer * particlesForCorrelations = AliAnalysisTaskEmcalJetHUtils::CreateParticleOrTrackContainer(trackName.c_str());
1385  if (particlesForCorrelations)
1386  {
1387  particlesForCorrelations->SetName("tracksForCorrelations");
1388  particlesForCorrelations->SetMinPt(0.15);
1389  particlesForCorrelations->SetEtaLimits(-1.0*trackEta, trackEta);
1390  // Adopt the container
1391  this->AdoptParticleContainer(particlesForCorrelations);
1392  }
1393  else {
1394  AliWarningStream() << "No particle container was successfully created!\n";
1395  }
1396 
1397  // Jets
1401  jetRadius,
1403  particlesForJets,
1404  clustersForJets);
1405  // 0.6 * jet area
1406  jetContainer->SetJetAreaCut(jetRadius * jetRadius * TMath::Pi() * 0.6);
1407  jetContainer->SetMaxTrackPt(100);
1408  jetContainer->SetJetPtCut(0.1);
1409 
1410  // Successfully configured
1411  returnValue = true;
1412 
1413  return returnValue;
1414 }
1415 
1422  std::string clusName,
1423  const double jetConstituentPtCut,
1424  const double trackEta,
1425  const double jetRadius,
1426  const std::string & jetTag,
1427  const std::string & correlationsTracksCutsPeriod)
1428 {
1429  bool returnValue = false;
1430  AliInfoStream() << "Configuring Jet-H Correlations task for an embedding analysis.\n";
1431 
1432  // Set the task to know it that is embedded
1433  this->SetIsEmbedded(true);
1434 
1435  // Add Containers
1436  // Clusters
1437  if (clusName == "usedefault") {
1439  }
1440  // For jet finding
1441  AliClusterContainer * clustersForJets = new AliClusterContainer(clusName.c_str());
1442  clustersForJets->SetName("clustersForJets");
1443  clustersForJets->SetMinE(jetConstituentPtCut);
1444  // We need the combined clusters, which should be available in the internal event.
1445  // However, we don't need to adopt the container - we'll just use it to find the right jet collection
1446  // For correlations
1447  /*AliClusterContainer * clustersforCorrelations = new AliClusterContainer("usedefault");
1448  clustersForCorrelations->SetName("clustersForCorrelations");
1449  clustersForCorrelations->SetMinE(0.30);
1450  clustersForCorrelations->SetIsEmbedding(true);
1451  this->AdoptClusterContainer(clustersForCorrelations);*/
1452 
1453  // Tracks
1454  // For jet finding
1455  if (trackName == "usedefault") {
1457  }
1459  particlesForJets->SetName("particlesForJets");
1460  particlesForJets->SetMinPt(jetConstituentPtCut);
1461  particlesForJets->SetEtaLimits(-1.0*trackEta, trackEta);
1462  // Don't need to adopt the container - we'll just use it to find the right jet collection
1463  // For correlations
1464  AliParticleContainer * particlesForCorrelations = AliAnalysisTaskEmcalJetHUtils::CreateParticleOrTrackContainer(trackName.c_str());
1465  // Ensure that we don't operate on a null pointer
1466  if (particlesForCorrelations)
1467  {
1468  particlesForCorrelations->SetName("tracksForCorrelations");
1469  particlesForCorrelations->SetMinPt(0.15);
1470  particlesForCorrelations->SetEtaLimits(-1.0*trackEta, trackEta);
1471  particlesForCorrelations->SetIsEmbedding(true);
1472  AliTrackContainer * trackCont = dynamic_cast<AliTrackContainer *>(particlesForCorrelations);
1473  if (trackCont) {
1474  // This option only exists for track containers
1475  trackCont->SetTrackCutsPeriod(correlationsTracksCutsPeriod.c_str());
1476  }
1477  // Adopt the container
1478  this->AdoptParticleContainer(particlesForCorrelations);
1479  }
1480  else {
1481  AliWarningStream() << "No particle container was successfully created!\n";
1482  }
1483 
1484  // Jets
1485  // The tag "hybridLevelJets" is defined in the jet finder
1489  jetRadius,
1491  particlesForJets,
1492  clustersForJets,
1493  jetTag);
1494  // 0.6 * jet area
1495  jetContainer->SetJetAreaCut(jetRadius * jetRadius * TMath::Pi() * 0.6);
1496  jetContainer->SetMaxTrackPt(100);
1497  jetContainer->SetJetPtCut(0.1);
1498 
1499  // Successfully configured
1500  returnValue = true;
1501 
1502  return returnValue;
1503 }
1504 
1511 {
1512  std::stringstream tempSS;
1513  tempSS << std::boolalpha;
1514  tempSS << "Jet collections:\n";
1515  TIter next(&fJetCollArray);
1516  AliJetContainer * jetCont;
1517  while ((jetCont = static_cast<AliJetContainer *>(next()))) {
1518  tempSS << "\t" << jetCont->GetName() << ": " << jetCont->GetArrayName() << "\n";
1519  }
1520  tempSS << "Event selection\n";
1521  tempSS << "\tUse AliEventCuts: " << !fUseBuiltinEventSelection << "\n";
1522  tempSS << "\tTrigger event selection: " << std::bitset<32>(fTriggerType) << "\n";
1523  tempSS << "\tMixed event selection: " << std::bitset<32>(fMixingEventType) << "\n";
1524  tempSS << "\tEnabled only for non-fast partition: " << fDisableFastPartition << "\n";
1525  tempSS << "Jet settings:\n";
1526  tempSS << "\tTrack bias: " << fTrackBias << "\n";
1527  tempSS << "\tCluster bias: " << fClusterBias << "\n";
1528  tempSS << "Event mixing:\n";
1529  tempSS << "\tEnabled: " << fDoEventMixing << "\n";
1530  tempSS << "\tN mixed tracks: " << fNMixingTracks << "\n";
1531  tempSS << "\tMin number of tracks for mixing: " << fMinNTracksMixedEvents << "\n";
1532  tempSS << "\tMin number of events for mixing: " << fMinNEventsMixedEvents << "\n";
1533  tempSS << "\tNumber of centrality bins for mixing: " << fNCentBinsMixedEvent << "\n";
1534  tempSS << "Histogramming options:\n";
1535  tempSS << "\tLess sparse axes: " << fDoLessSparseAxes << "\n";
1536  tempSS << "\tWider associated track pt bins: " << fDoWiderTrackBin << "\n";
1537  tempSS << "Jet matching options for embedding:\n";
1538  tempSS << "\tRequire the jet to match to a embedded jets: " << fRequireMatchedJetWhenEmbedding << "\n";
1539  tempSS << "\tRequire an additional match to a part level jet: " << fRequireMatchedPartLevelJet << "\n";
1540  tempSS << "\tMinimum shared momentum fraction: " << fMinSharedMomentumFraction << "\n";
1541  tempSS << "\tMax matched jet distance: " << fMaxMatchedJetDistance << "\n";
1542  tempSS << "Efficiency\n";
1543  tempSS << "\tSingle track efficiency identifier: " << fEfficiencyPeriodIdentifier << "\n";
1544  tempSS << "\tArtifical track inefficiency: " << fArtificialTrackInefficiency << "\n";
1545 
1546  return tempSS.str();
1547 }
1548 
1555 std::ostream & AliAnalysisTaskEmcalJetHCorrelations::Print(std::ostream & in) const {
1556  in << toString();
1557  return in;
1558 }
1559 
1567 {
1568  Printf("%s", toString().c_str());
1569 }
1570 
1571 } /* namespace EMCALJetTasks */
1572 } /* namespace PWGJE */
1573 
1581 std::ostream & operator<<(std::ostream & in, const PWGJE::EMCALJetTasks::AliAnalysisTaskEmcalJetHCorrelations & myTask)
1582 {
1583  std::ostream & result = myTask.Print(in);
1584  return result;
1585 }
1586 
bool CheckArtificialTrackEfficiency(unsigned int trackIndex, std::vector< unsigned int > &rejectedTrackIndices, bool useRejectedList)
Bool_t fNoMixedEventJESCorrection
True if the jet energy scale correction should be applied to mixed event histograms.
void AccessSetOfYBinValues(TH2D *hist, Int_t xBin, std::vector< Double_t > &yBinsContent, Double_t scaleFactor=-1.0)
bool ConfigureForEmbeddingAnalysis(std::string trackName="usedefault", std::string clusName="caloClustersCombined", const double jetConstituentPtCut=3, const double trackEta=0.8, const double jetRadius=0.2, const std::string &jetTag="hybridLevelJets", const std::string &correlationsTracksCutsPeriod="lhc11a")
const char * filename
Definition: TestFCM.C:1
void SetTrackCutsPeriod(const char *period)
Double_t GetRhoVal() const
double Double_t
Definition: External.C:58
static double DetermineTrackingEfficiency(const double trackPt, const double trackEta, const int centralityBin, const EEfficiencyPeriodIdentifier_t efficiencyPeriodIdentifier, const std::string &taskName)
static Double_t DeltaPhi(Double_t phia, Double_t phib, Double_t rMin=-TMath::Pi()/2, Double_t rMax=3 *TMath::Pi()/2)
Definition: External.C:236
Int_t GetNTracks() const
const char * title
Definition: MakeQAPdf.C:27
bool CheckForMatchedJet(AliJetContainer *jets, AliEmcalJet *jet, const std::string &histName)
AliEmcalJet * ClosestJet() const
Definition: AliEmcalJet.h:341
AliJetContainer * GetJetContainer(Int_t i=0) const
void SetName(const char *n)
Set the name of the class of the objets inside the underlying array.
Double_t ClosestJetDistance() const
Definition: AliEmcalJet.h:342
void AdoptParticleContainer(AliParticleContainer *cont)
void SetMinE(Double_t min)
TH2D * fJESCorrectionHist
Histogram containing the jet energy scale correction.
AliAnalysisTaskEmcalJetHUtils::EEfficiencyPeriodIdentifier_t fEfficiencyPeriodIdentifier
Identifies the period for determining the efficiency correction to apply.
Container with name, TClonesArray and cuts for particles.
virtual void SetTrackBias(Double_t b)
Require a track with pt > b in jet.
Declaration of class AliTLorentzVector.
TH1 * fHistJetPtBias[6]
! Jet pt spectrum of jets which meet the constituent bias criteria (the array corresponds to centrali...
Double_t fEPV0
!event plane V0
void SetMinPt(Double_t min)
Bool_t fGeneralHistograms
whether or not it should fill some general histograms
Declaration of class AliAnalysisTaskEmcalEmbeddingHelper.
TCanvas * c
Definition: TestFitELoss.C:172
Bool_t fDoLessSparseAxes
True if there should be fewer THnSparse axes.
AliJetContainer * AddJetContainer(const char *n, TString defaultCutType, Float_t jetRadius=0.4)
virtual void UserExecOnce()
Task initializations handled in user tasks.
Int_t fCentBin
!event centrality bin
bool GetProperty(std::vector< std::string > propertyPath, const std::string &propertyName, T &property, const bool requiredProperty) const
Double_t fMaxMatchedJetDistance
Maximum distance between two matched jets.
Bool_t fDoWiderTrackBin
True if the track pt bins in the THnSparse should be wider.
virtual Double_t Pt() const
UInt_t fTriggerType
Event selection for jets (ie triggered events).
void SetEtaLimits(Double_t min, Double_t max)
Bool_t fUseBuiltinEventSelection
Use builtin event selection of the AliAnalysisTaskEmcal instead of AliEventCuts.
bool ConfigureForStandardAnalysis(std::string trackName="usedefault", std::string clusName="usedefault", const double jetConstituentPtCut=3, const double trackEta=0.8, const double jetRadius=0.2)
virtual Double_t Eta() const
bool DoesConfigurationExist(const int i) const
Container for particles within the EMCAL framework.
TObjArray * CloneAndReduceTrackList(std::vector< unsigned int > &rejectedTrackIndices, const bool useRejectedList)
Bool_t fIsEmbedded
trigger, embedded signal
void SetIsEmbedding(Bool_t b)
Set embedding status.
AliParticleContainer * GetParticleContainer(Int_t i=0) const
Get particle container attached to this task.
static void ConfigureEventCuts(AliEventCuts &eventCuts, PWG::Tools::AliYAMLConfiguration &yamlConfig, const UInt_t offlineTriggerMask, const std::string &baseName, const std::string &taskName)
AliEmcalJet * GetLeadingJet(const char *opt="")
static UInt_t DeterminePhysicsSelectionFromYAML(const std::vector< std::string > &selections)
virtual THnSparse * NewTHnSparseF(const char *name, UInt_t entries)
int Int_t
Definition: External.C:63
virtual void SetClusterBias(Double_t b)
Require a cluster with pt > b in jet.
void SetJetPtCut(Float_t cut)
unsigned int UInt_t
Definition: External.C:33
THashList * GetListOfHistograms() const
Get the list of histograms.
Definition: THistManager.h:671
Int_t GetNJets() const
TString toString() const
TH1 * CreateTH1(const char *name, const char *title, int nbins, double xmin, double xmax, Option_t *opt="")
Create a new TH1 within the container.
BeamType fForceBeamType
forced beam type
Definition: External.C:228
Int_t fNcentBins
how many centrality bins
Double_t MaxTrackPt() const
Definition: AliEmcalJet.h:155
void FillHist(TH1 *hist, Double_t fillValue, Double_t weight=1.0, Bool_t noCorrection=kFALSE)
BeamType fBeamType
!event beam type
static double GetJetPt(const AliEmcalJet *jet, const double rho)
Double_t fCent
!event centrality
virtual void SetMixedEventTriggerType(UInt_t me)
Set the mixed event trigger selection.
Bool_t RetrieveAndInitializeJESCorrectionHist(TString filename, TString histName, Double_t trackBias=AliAnalysisTaskEmcalJetHCorrelations::kDisableBias, Double_t clusterBias=AliAnalysisTaskEmcalJetHCorrelations::kDisableBias)
static AliParticleContainer * CreateParticleOrTrackContainer(const std::string &collectionName)
static const std::map< std::string, EEfficiencyPeriodIdentifier_t > fgkEfficiencyPeriodIdentifier
! Map from name to efficiency period identifier for use with the YAML config
void FillTH1(const char *hname, double x, double weight=1., Option_t *opt="")
Fill a 1D histogram within the container.
static double RelativeEPAngle(double jetAngle, double epAngle)
TObjArray fJetCollArray
jet collection array
Double_t fMinSharedMomentumFraction
Minimum shared momentum with matched jet.
virtual void SetNCentBins(Int_t n)
static AliAnalysisTaskEmcalJetHCorrelations * AddTaskEmcalJetHCorrelations(const char *nTracks="usedefault", const char *nCaloClusters="usedefault", const Double_t trackBias=5, const Double_t clusterBias=5, const Int_t nTracksMixedEvent=0, const Int_t minNTracksMixedEvent=5000, const Int_t minNEventsMixedEvent=5, const UInt_t nCentBinsMixedEvent=10, UInt_t trigEvent=AliVEvent::kAny, UInt_t mixEvent=AliVEvent::kAny, const Bool_t lessSparseAxes=kFALSE, const Bool_t widerTrackBin=kFALSE, const Bool_t JESCorrection=kFALSE, const char *JESCorrectionFilename="alien:///alice/cern.ch/user/r/rehlersi/JESCorrection.root", const char *JESCorrectionHistName="JESCorrection", const char *suffix="biased")
const TString & GetArrayName() const
const char * GetName() const
Bool_t fRequireMatchedJetWhenEmbedding
True if jets are required to be matched (ie. jet->MatchedJet() != nullptr)
static Double_t * GenerateFixedBinArray(Int_t n, Double_t min, Double_t max)
AliEventCuts fAliEventCuts
Event cuts (run2 defaults)
AliEmcalList * fOutput
!output list
PWG::Tools::AliYAMLConfiguration fYAMLConfig
YAML configuration file.
TH2 * fHistTrackEtaPhi[7]
! Track eta-phi distribution (the array corresponds to track pt)
bool fConfigurationInitialized
True if the task configuration has been successfully initialized.
TH1 * fHistJetPt[6]
! Jet pt spectrum (the array corresponds to centrality bins)
Double_t fVertex[3]
!event vertex
Bool_t fDisableFastPartition
True if task should be disabled for the fast partition, where the EMCal is not included.
AliTrackContainer * GetTrackContainer(Int_t i=0) const
TH1 * fHistEventRejection
!book keep reasons for rejecting event
void SetMakeGeneralHistograms(Bool_t g)
Enable general histograms.
Base task in the EMCAL jet framework.
Double_t MaxClusterPt() const
Definition: AliEmcalJet.h:154
Represent a jet reconstructed using the EMCal jet framework.
Definition: AliEmcalJet.h:51
static std::string DetermineUseDefaultName(InputObject_t objType)
Double_t GetFractionSharedPt(const AliEmcalJet *jet, AliParticleContainer *cont2=0x0) const
const AliTrackIterableMomentumContainer accepted_momentum() const
const char Option_t
Definition: External.C:48
TH2 * fHistJetHEtaPhi
! Eta-phi distribution of jets which are in jet-hadron correlations
void GetDeltaEtaDeltaPhiDeltaR(AliTLorentzVector &particleOne, AliVParticle *particleTwo, Double_t &deltaEta, Double_t &deltaPhi, Double_t &deltaR)
void UserCreateOutputObjects()
Main initialization function on the worker.
const Double_t pi
const Int_t nbins
const AliJetIterableContainer accepted() const
bool Bool_t
Definition: External.C:53
friend std::ostream & operator<<(std::ostream &in, const AliAnalysisTaskEmcalJetHCorrelations &myTask)
virtual void SetTriggerType(UInt_t te)
Set the trigger event trigger selection.
void SetMaxTrackPt(Float_t b)
virtual Double_t Phi() const
Container structure for EMCAL clusters.
EMCal fiducial acceptance (each eta, phi edge narrowed by jet R)
Definition: AliEmcalJet.h:70
! Arbitrarily large value which can be used to disable the constituent bias. Can be used for either t...
virtual void GetDimParams(Int_t iEntry, TString &label, Int_t &nbins, Double_t &xmin, Double_t &xmax)
Container for jet within the EMCAL jet framework.
bool fRequireMatchedPartLevelJet
True if matched jets are required to be matched to a particle level jet.
Definition: External.C:196
Double_t fArtificialTrackInefficiency
Artificial track inefficiency. Enabled if < 1.0.
void SetJetAreaCut(Float_t cut)
static const AliAnalysisTaskEmcalEmbeddingHelper * GetInstance()
TRandom3 fRandom
! Random number generator for artificial track inefficiency.