AliPhysics  e6c8d43 (e6c8d43)
AliAnalysisTaskEmcalTriggerJetsIDcorr.cxx
Go to the documentation of this file.
1 /************************************************************************************
2  * Copyright (C) 2017, Copyright Holders of the ALICE Collaboration *
3  * All rights reserved. *
4  * *
5  * Redistribution and use in source and binary forms, with or without *
6  * modification, are permitted provided that the following conditions are met: *
7  * * Redistributions of source code must retain the above copyright *
8  * notice, this list of conditions and the following disclaimer. *
9  * * Redistributions in binary form must reproduce the above copyright *
10  * notice, this list of conditions and the following disclaimer in the *
11  * documentation and/or other materials provided with the distribution. *
12  * * Neither the name of the <organization> nor the *
13  * names of its contributors may be used to endorse or promote products *
14  * derived from this software without specific prior written permission. *
15  * *
16  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND *
17  * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED *
18  * WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE *
19  * DISCLAIMED. IN NO EVENT SHALL ALICE COLLABORATION BE LIABLE FOR ANY *
20  * DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES *
21  * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; *
22  * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND *
23  * ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT *
24  * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS *
25  * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. *
26  ************************************************************************************/
27 #include <THistManager.h>
28 #include <TLinearBinning.h>
29 #include <TMath.h>
30 #include <TString.h>
31 #include <TVector3.h>
32 
33 #include "AliAnalysisManager.h"
34 #include "AliAODTrack.h"
35 #include "AliClusterContainer.h"
36 #include "AliEmcalJet.h"
37 #include "AliInputEventHandler.h"
38 #include "AliJetContainer.h"
39 #include "AliLog.h"
40 #include "AliPIDResponse.h"
41 #include "AliTrackContainer.h"
43 #include "AliVCluster.h"
44 #include "AliVEvent.h"
45 #include "AliVTrack.h"
46 
47 #include <array>
48 #include <iostream>
49 
53 
54 namespace EmcalTriggerJets {
55 
56 AliAnalysisTaskEmcalTriggerJetsIDcorr::AliAnalysisTaskEmcalTriggerJetsIDcorr():
58  fJetCont(nullptr),
59  fPIDResponse(nullptr),
60  fHistos(nullptr)
61 {
62 
63 }
64 
66  AliAnalysisTaskEmcalJet(name, true),
70 {
71 
72 }
73 
75  delete fHistos;
76 }
77 
80 
81  const std::array<TString, 3> kEmcalTriggers = {"INT7", "EJ1", "EJ2"};
82  const int kNJetPtBins = 9;
83  const int kNJetRadiusBins = 7;
84 
85  TLinearBinning jetptbinning(9, 20, 200), pbinning(300, 0., 30.), massbinning(600., 0., 6.), radiusBinning(10, 0., 1.);
86  const TBinning *binningPID[4] = {&jetptbinning, &pbinning, &radiusBinning, &massbinning};
87  const std::array<TString, 2> kSpecies = {"Proton", "Deuteron"};
88  fHistos = new THistManager("EmcalJetHistos");
89  fHistos->CreateTH1("hTPCdEdxErrors", "Error counter for TPC dE/dx", AliPID::kSPECIESC, -0.5, AliPID::kSPECIESC - 0.5);
90  for(auto t : kEmcalTriggers){
91  fHistos->CreateTH1("hEventCount" + t, "Number of events for trigger " + t, 1, 0.5, 1.5);
92  fHistos->CreateTH1("hPtRawAllJet" + t, " Raw jet pt spectrum for trigger " + t, 1000, 0., 1000.);
93  fHistos->CreateTH1("hPtRawSelJet" + t, " Raw jet pt spectrum for trigger " + t, 1000, 0., 1000.);
94  fHistos->CreateTH1("hNJetsAll" + t, "All reconstructed jets for trigger " + t, 101, -0.5, 100.5);
95  fHistos->CreateTH1("hNJetsSelected" + t, "Selected reconstructed jets for trigger " + t, 101, -0.5, 101.5);
96  for(auto s: kSpecies){
97  fHistos->CreateTH1("hNCandidatesPerEvent" + s + t, "Number of " + s + " candidates per events ", 101, -0.5, 100.5);
98  fHistos->CreateTHnSparse("hPIDAssociateFullJet" + s + t, TString::Format("PID (TOF masss) for particles associated to full jets for R=0.4f in EMCAL for %s candidates for trigger %s", s.Data(), t.Data()), 4, binningPID);
99  }
100  }
101  for(auto h : *(fHistos->GetListOfHistograms())){
102  fOutput->Add(h);
103  }
104  PostData(1, fOutput);
105 }
106 
108  fPIDResponse = fInputHandler->GetPIDResponse();
109  if(!fPIDResponse){
110  AliErrorStream() << "PID Response not available - PID plots will not be filled" << std::endl;
111  }
112  fJetCont = this->GetJetContainer("FullJetsR04EMCAL");
113 }
114 
116  std::vector<TString> triggers, kEmcalTriggers = {"EJ1", "EJ2"};
117  if(fInputHandler->IsEventSelected() & AliVEvent::kINT7) triggers.push_back("INT7");
118  if(fInputHandler->IsEventSelected() & AliVEvent::kEMCEJE){
119  TString fired = fInputEvent->GetFiredTriggerClasses();
120  for(auto e : kEmcalTriggers){
121  if(fired.Contains(e))
122  triggers.push_back(e);
123  }
124  }
125  if(!triggers.size()) return false;
126 
127  for(const auto &t : triggers) fHistos->FillTH1("hEventCount" + t, 1);
128 
129  // Get (rough) TPC candidates for protons and deuterons
130  std::vector<AliVTrack *> protonCandidates , deuteronCandidates;
131  try{
132  protonCandidates = this->GetTPCPIDCandidates(AliPID::kProton);
133  } catch(TPCdEdxException &e) {
134  fHistos->FillTH1("hTPCdEdxErrors", e.GetParticleType());
135  }
136  try{
137  deuteronCandidates = this->GetTPCPIDCandidates(AliPID::kDeuteron);
138  } catch(TPCdEdxException &e) {
139  fHistos->FillTH1("hTPCdEdxErrors", e.GetParticleType());
140  }
141  for(const auto &t : triggers){
142  fHistos->FillTH1("hNCandidatesPerEventProton" + t, protonCandidates.size());
143  fHistos->FillTH1("hNCandidatesPerEventDeuteron" + t, deuteronCandidates.size());
144  }
145 
146  int njetsAll(0), njetsSel(0);
147  for(auto j : fJetCont->all()) {
148  njetsAll++;
149  for(const auto &t : triggers) fHistos->FillTH1("hPtRawAllJet" + t, j->Pt());
150  }
151  for(auto j : fJetCont->accepted()){
152  njetsSel++;
153  for(const auto &t : triggers) fHistos->FillTH1("hPtRawSelJet" + t, j->Pt());
154  TVector3 jetmomentum(j->Px(), j->Py(), j->Pz());
155  // Proton-jet correlations
156  std::vector<CorrParticleInfo> protonsCorrelated = CorrelateCandidatesToJet(jetmomentum, protonCandidates);
157  for(auto c : protonsCorrelated) {
158  Double_t point[4] = {TMath::Abs(j->Pt()), c.fPt, c.fDR, c.fMass*c.fMass};
159  for(auto t : triggers) fHistos->FillTHnSparse("hPIDAssociateFullJetProton" +t, point);
160  }
161 
162  // deuteron-jet correlations
163  std::vector<CorrParticleInfo> deuteronsCorrelated = CorrelateCandidatesToJet(jetmomentum, deuteronCandidates);
164  for(auto c : deuteronsCorrelated) {
165  Double_t point[4] = {TMath::Abs(j->Pt()), c.fPt, c.fDR, c.fMass*c.fMass};
166  for(auto t : triggers) fHistos->FillTHnSparse("hPIDAssociateFullJetDeuteron" +t, point);
167  }
168  }
169  for(const auto &t : triggers){
170  fHistos->FillTH1("hNJetsAll" + t, njetsAll);
171  fHistos->FillTH1("hNJetsSelected" + t, njetsSel);
172  }
173 
174  return true;
175 }
176 
177 std::vector<CorrParticleInfo> AliAnalysisTaskEmcalTriggerJetsIDcorr::CorrelateCandidatesToJet(const TVector3 &jet, std::vector<AliVTrack *> candidates) const {
178  std::vector<CorrParticleInfo> correlatedParticles;
179  for(auto c : candidates) {
180  TVector3 particleMom(c->Px(), c->Py(), c->Pz());
181  double dR = particleMom.DeltaR(jet);
182  if(dR > 1.) continue;
183  double mTOF = -1;
184  try {
185  mTOF = GetTOFMass(c);
186  } catch(TOFMassException &e) {
187  continue;
188  }
189  correlatedParticles.push_back({c->Pt(), dR, mTOF});
190  }
191  return correlatedParticles;
192 }
193 
194 double AliAnalysisTaskEmcalTriggerJetsIDcorr::GetTOFMass(const AliVTrack *const track) const {
195  if(!((track->GetStatus() & AliVTrack::kTOFout) && (track->GetStatus() & AliVTrack::kTIME))) throw TOFMassException();
196  Double_t trtime = (track->GetTOFsignal() - fPIDResponse->GetTOFResponse().GetTimeZero()) * 1e-12;
197  Double_t v = track->GetIntegratedLength()/(100. * trtime);
198  Double_t beta = v / TMath::C(), gamma = 1 / TMath::Sqrt(1-beta*beta);
199  return track->P() / (beta*gamma);
200 }
201 
202 std::vector<AliVTrack *> AliAnalysisTaskEmcalTriggerJetsIDcorr::GetTPCPIDCandidates(AliPID::EParticleType type) const {
203  std::vector<AliVTrack *> result;
204  for(int itrk = 0; itrk < fInputEvent->GetNumberOfTracks(); itrk++){
205  AliVTrack *trk = static_cast<AliVTrack *>(fInputEvent->GetTrack(itrk));
206  if(!trk) continue;
207  if(trk->IsA() == AliAODTrack::Class()){
208  AliAODTrack *aodtrk = static_cast<AliAODTrack *>(trk);
209  if(!(aodtrk->IsHybridGlobalConstrainedGlobal() || aodtrk->IsHybridTPCConstrainedGlobal())) continue;
210  if(aodtrk->GetTPCsignalN() < 30) continue;
211  if(aodtrk->GetTPCSharedMapPtr()->CountBits(0) > 70) continue; // exclude tracks with more than 70 shared clusters ->
212  Double_t nSigmaTPC = fPIDResponse->NumberOfSigmasTPC(trk, type);
213  if(nSigmaTPC <= -999.) continue; //throw TPCdEdxException(type);
214  if(TMath::Abs(nSigmaTPC) > 3) continue;
215  result.push_back(trk);
216  }
217  }
218  return result;
219 }
220 
222  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
223 
225  mgr->AddTask(task);
226 
227  // Adding containers for clusters and tracks
228  AliClusterContainer *clustercont = task->AddClusterContainer("caloClusters");
229  clustercont->SetMinE(0.3);
230  AliTrackContainer *trackcont = task->AddTrackContainer("tracks");
231  trackcont->SetMinPt(0.15);
232 
233  // Adding Jet containers
234  // - Using jet radius 0.4
235  // - Only EMCAL side
236  // - processing only full jets
237 
238 
239  // Full jets, R=0.4, EMCAL
240  AliJetContainer *cont = task->AddJetContainer(
244  0.4,
246  trackcont, clustercont
247  );
248  cont->SetName("FullJetsR04EMCAL");
249  cont->SetJetPtCut(20.);
250  cont->SetJetPtCutMax(200.);
251 
252 
253  // Connect Input / Output containers
254  TString outfilename = mgr->GetCommonFileName();
255  outfilename += ":EmcalTriggerJetsIDcorr";
256  mgr->ConnectInput(task, 0, mgr->GetCommonInputContainer());
257  mgr->ConnectOutput(task, 1, mgr->CreateContainer("HistsEmcalTriggerJetsIDcorr", TList::Class(), AliAnalysisManager::kOutputContainer, outfilename));
258 
259  return task;
260 }
261 
262 } /* namespace EmcalTriggerJets */
double Double_t
Definition: External.C:58
Class creating a linear binning, used in the histogram manager.
AliJetContainer * GetJetContainer(Int_t i=0) const
Container with name, TClonesArray and cuts for particles.
virtual bool Run()
Run function. This is the core function of the analysis and contains the user code. Therefore users have to implement this function.
std::vector< AliVTrack * > GetTPCPIDCandidates(AliPID::EParticleType type) const
virtual void UserExecOnce()
Task initializations handled in user tasks.
TCanvas * c
Definition: TestFitELoss.C:172
AliJetContainer * AddJetContainer(const char *n, TString defaultCutType, Float_t jetRadius=0.4)
Interface for binnings used by the histogram handler.
Definition: TBinning.h:21
AliClusterContainer * AddClusterContainer(const char *n)
Create new cluster container and attach it to the task.
void FillTHnSparse(const char *name, const double *x, double weight=1., Option_t *opt="")
static AliAnalysisTaskEmcalTriggerJetsIDcorr * AddTaskEmcalTriggerJetsIDcorr(const char *name)
void SetJetPtCut(Float_t cut)
THashList * GetListOfHistograms() const
Get the list of histograms.
Definition: THistManager.h:671
TH1 * CreateTH1(const char *name, const char *title, int nbins, double xmin, double xmax, Option_t *opt="")
Create a new TH1 within the container.
void FillTH1(const char *hname, double x, double weight=1., Option_t *opt="")
Fill a 1D histogram within the container.
AliEmcalList * fOutput
!output list
std::vector< CorrParticleInfo > CorrelateCandidatesToJet(const TVector3 &jet, std::vector< AliVTrack * > candidates) const
AliTrackContainer * AddTrackContainer(const char *n)
Create new track container and attach it to the task.
Base task in the EMCAL jet framework.
Container class for histograms.
Definition: THistManager.h:99
void UserCreateOutputObjects()
Main initialization function on the worker.
const AliJetIterableContainer accepted() const
THnSparse * CreateTHnSparse(const char *name, const char *title, int ndim, const int *nbins, const double *min, const double *max, Option_t *opt="")
Create a new THnSparse within the container.
Container structure for EMCAL clusters.
EMCal fiducial acceptance (each eta, phi edge narrowed by jet R)
Definition: AliEmcalJet.h:70
Container for jet within the EMCAL jet framework.
void SetJetPtCutMax(Float_t cut)
const AliJetIterableContainer all() const