2 // Jet mass analysis task.
6 #include <TClonesArray.h>
10 #include <THnSparse.h>
12 #include <TLorentzVector.h>
19 #include "AliVCluster.h"
20 #include "AliVTrack.h"
21 #include "AliEmcalJet.h"
22 #include "AliRhoParameter.h"
24 #include "AliEmcalParticle.h"
25 #include "AliMCEvent.h"
26 #include "AliGenPythiaEventHeader.h"
27 #include "AliAODMCHeader.h"
28 #include "AliMCEvent.h"
29 #include "AliAnalysisManager.h"
30 #include "AliJetContainer.h"
32 #include "AliAODEvent.h"
34 #include "AliAnalysisTaskEmcalJetMass.h"
36 ClassImp(AliAnalysisTaskEmcalJetMass)
38 //________________________________________________________________________
39 AliAnalysisTaskEmcalJetMass::AliAnalysisTaskEmcalJetMass() :
40 AliAnalysisTaskEmcalJet("AliAnalysisTaskEmcalJetMass", kTRUE),
42 fMinFractionShared(0),
43 fh2PtJet1VsLeadPtAllSel(0),
44 fh2PtJet1VsLeadPtTagged(0),
45 fh2PtVsMassJet1All(0),
46 fh2PtVsMassJet1Tagged(0),
48 fpPtVsMassJet1Tagged(0),
49 fh2MassVsAreaJet1All(0),
50 fh2MassVsAreaJet1Tagged(0),
51 fh2MassVsNConstJet1All(0),
52 fh2MassVsNConstJet1Tagged(0),
55 // Default constructor.
57 fh2PtJet1VsLeadPtAllSel = new TH2F*[fNcentBins];
58 fh2PtJet1VsLeadPtTagged = new TH2F*[fNcentBins];
59 fh2PtVsMassJet1All = new TH2F*[fNcentBins];
60 fh2PtVsMassJet1Tagged = new TH2F*[fNcentBins];
61 fpPtVsMassJet1All = new TProfile*[fNcentBins];
62 fpPtVsMassJet1Tagged = new TProfile*[fNcentBins];
63 fh2MassVsAreaJet1All = new TH2F*[fNcentBins];
64 fh2MassVsAreaJet1Tagged = new TH2F*[fNcentBins];
65 fh2MassVsNConstJet1All = new TH2F*[fNcentBins];
66 fh2MassVsNConstJet1Tagged = new TH2F*[fNcentBins];
67 fh2EtMassOverEtRSq = new TH2F*[fNcentBins];
69 for (Int_t i = 0; i < fNcentBins; i++) {
70 fh2PtJet1VsLeadPtAllSel[i] = 0;
71 fh2PtJet1VsLeadPtTagged[i] = 0;
72 fh2PtVsMassJet1All[i] = 0;
73 fh2PtVsMassJet1Tagged[i] = 0;
74 fpPtVsMassJet1All[i] = 0;
75 fpPtVsMassJet1Tagged[i] = 0;
76 fh2MassVsAreaJet1All[i] = 0;
77 fh2MassVsAreaJet1Tagged[i] = 0;
78 fh2MassVsNConstJet1All[i] = 0;
79 fh2MassVsNConstJet1Tagged[i] = 0;
80 fh2EtMassOverEtRSq[i] = 0;
83 SetMakeGeneralHistograms(kTRUE);
87 //________________________________________________________________________
88 AliAnalysisTaskEmcalJetMass::AliAnalysisTaskEmcalJetMass(const char *name) :
89 AliAnalysisTaskEmcalJet(name, kTRUE),
91 fMinFractionShared(0),
92 fh2PtJet1VsLeadPtAllSel(0),
93 fh2PtJet1VsLeadPtTagged(0),
94 fh2PtVsMassJet1All(0),
95 fh2PtVsMassJet1Tagged(0),
97 fpPtVsMassJet1Tagged(0),
98 fh2MassVsAreaJet1All(0),
99 fh2MassVsAreaJet1Tagged(0),
100 fh2MassVsNConstJet1All(0),
101 fh2MassVsNConstJet1Tagged(0),
102 fh2EtMassOverEtRSq(0)
104 // Standard constructor.
106 fh2PtJet1VsLeadPtAllSel = new TH2F*[fNcentBins];
107 fh2PtJet1VsLeadPtTagged = new TH2F*[fNcentBins];
108 fh2PtVsMassJet1All = new TH2F*[fNcentBins];
109 fh2PtVsMassJet1Tagged = new TH2F*[fNcentBins];
110 fpPtVsMassJet1All = new TProfile*[fNcentBins];
111 fpPtVsMassJet1Tagged = new TProfile*[fNcentBins];
112 fh2MassVsAreaJet1All = new TH2F*[fNcentBins];
113 fh2MassVsAreaJet1Tagged = new TH2F*[fNcentBins];
114 fh2MassVsNConstJet1All = new TH2F*[fNcentBins];
115 fh2MassVsNConstJet1Tagged = new TH2F*[fNcentBins];
116 fh2EtMassOverEtRSq = new TH2F*[fNcentBins];
118 for (Int_t i = 0; i < fNcentBins; i++) {
119 fh2PtJet1VsLeadPtAllSel[i] = 0;
120 fh2PtJet1VsLeadPtTagged[i] = 0;
121 fh2PtVsMassJet1All[i] = 0;
122 fh2PtVsMassJet1Tagged[i] = 0;
123 fpPtVsMassJet1All[i] = 0;
124 fpPtVsMassJet1Tagged[i] = 0;
125 fh2MassVsAreaJet1All[i] = 0;
126 fh2MassVsAreaJet1Tagged[i] = 0;
127 fh2MassVsNConstJet1All[i] = 0;
128 fh2MassVsNConstJet1Tagged[i] = 0;
129 fh2EtMassOverEtRSq[i] = 0;
132 SetMakeGeneralHistograms(kTRUE);
135 //________________________________________________________________________
136 AliAnalysisTaskEmcalJetMass::~AliAnalysisTaskEmcalJetMass()
141 //________________________________________________________________________
142 void AliAnalysisTaskEmcalJetMass::UserCreateOutputObjects()
144 // Create user output.
146 AliAnalysisTaskEmcalJet::UserCreateOutputObjects();
148 Bool_t oldStatus = TH1::AddDirectoryStatus();
149 TH1::AddDirectory(kFALSE);
151 const Int_t nBinsPt = 250;
152 const Double_t minPt = -50.;
153 const Double_t maxPt = 200.;
155 const Int_t nBinsArea = 100;
156 const Double_t minArea = 0.;
157 const Double_t maxArea = 1.;
159 const Int_t nBinsNConst = 100;
160 const Double_t minNConst = 0.;
161 const Double_t maxNConst = 500.;
163 TString histName = "";
164 TString histTitle = "";
165 for (Int_t i = 0; i < fNcentBins; i++) {
166 histName = TString::Format("fh2PtJet1VsLeadPtAllSel_%d",i);
167 histTitle = TString::Format("%s;#it{p}_{T,jet1};#it{p}_{T,lead trk}",histName.Data());
168 fh2PtJet1VsLeadPtAllSel[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,20,0.,20.);
169 fOutput->Add(fh2PtJet1VsLeadPtAllSel[i]);
171 histName = TString::Format("fh2PtJet1VsLeadPtTagged_%d",i);
172 histTitle = TString::Format("%s;#it{p}_{T,jet1};#it{p}_{T,lead trk}",histName.Data());
173 fh2PtJet1VsLeadPtTagged[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,20,0.,20.);
174 fOutput->Add(fh2PtJet1VsLeadPtTagged[i]);
176 histName = TString::Format("fh2PtVsMassJet1All_%d",i);
177 histTitle = TString::Format("%s;#it{p}_{T,jet1};#it{M}_{jet1}",histName.Data());
178 fh2PtVsMassJet1All[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,nBinsPt,minPt,maxPt);
179 fOutput->Add(fh2PtVsMassJet1All[i]);
181 histName = TString::Format("fh2PtVsMassJet1Tagged_%d",i);
182 histTitle = TString::Format("%s;#it{p}_{T,jet1};#it{M}_{jet1}",histName.Data());
183 fh2PtVsMassJet1Tagged[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,nBinsPt,minPt,maxPt);
184 fOutput->Add(fh2PtVsMassJet1Tagged[i]);
186 histName = TString::Format("fpPtVsMassJet1All_%d",i);
187 histTitle = TString::Format("%s;#it{p}_{T,jet1};Avg #it{M}_{jet1}",histName.Data());
188 fpPtVsMassJet1All[i] = new TProfile(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt);
189 fOutput->Add(fpPtVsMassJet1All[i]);
191 histName = TString::Format("fpPtVsMassJet1Tagged_%d",i);
192 histTitle = TString::Format("%s;#it{p}_{T,jet1};Avg #it{M}_{jet1}",histName.Data());
193 fpPtVsMassJet1Tagged[i] = new TProfile(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt);
194 fOutput->Add(fpPtVsMassJet1Tagged[i]);
196 histName = TString::Format("fh2MassVsAreaJet1All_%d",i);
197 histTitle = TString::Format("%s;#it{M}_{jet1};#it{A}",histName.Data());
198 fh2MassVsAreaJet1All[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,nBinsArea,minArea,maxArea);
199 fOutput->Add(fh2MassVsAreaJet1All[i]);
201 histName = TString::Format("fh2MassVsAreaJet1Tagged_%d",i);
202 histTitle = TString::Format("%s;#it{M}_{jet1};#it{A}",histName.Data());
203 fh2MassVsAreaJet1Tagged[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,nBinsArea,minArea,maxArea);
204 fOutput->Add(fh2MassVsAreaJet1Tagged[i]);
206 histName = TString::Format("fh2MassVsNConstJet1All_%d",i);
207 histTitle = TString::Format("%s;#it{M}_{jet1};#it{N}_{constituents}",histName.Data());
208 fh2MassVsNConstJet1All[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,nBinsNConst,minNConst,maxNConst);
209 fOutput->Add(fh2MassVsNConstJet1All[i]);
211 histName = TString::Format("fh2MassVsNConstJet1Tagged_%d",i);
212 histTitle = TString::Format("%s;#it{M}_{jet1};#it{N}_{constituents}",histName.Data());
213 fh2MassVsNConstJet1Tagged[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,nBinsNConst,minNConst,maxNConst);
214 fOutput->Add(fh2MassVsNConstJet1Tagged[i]);
216 histName = TString::Format("fh2EtMassOverEtRSq_%d",i);
217 histTitle = TString::Format("%s;#it{E}_{T};(#it{M}/(#it{E}_{T}#it{R}))^{2}",histName.Data());
218 fh2EtMassOverEtRSq[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,100,0.,1.);
219 fOutput->Add(fh2EtMassOverEtRSq[i]);
222 // =========== Switch on Sumw2 for all histos ===========
223 for (Int_t i=0; i<fOutput->GetEntries(); ++i) {
224 TH1 *h1 = dynamic_cast<TH1*>(fOutput->At(i));
229 THnSparse *hn = dynamic_cast<THnSparse*>(fOutput->At(i));
233 TH1::AddDirectory(oldStatus);
235 PostData(1, fOutput); // Post data for ALL output slots > 0 here.
238 //________________________________________________________________________
239 Bool_t AliAnalysisTaskEmcalJetMass::Run()
241 // Run analysis code here, if needed. It will be executed before FillHistograms().
246 //________________________________________________________________________
247 Bool_t AliAnalysisTaskEmcalJetMass::FillHistograms()
251 AliEmcalJet* jet1 = NULL;
253 AliJetContainer *jetCont = GetJetContainer(fContainerBase);
255 jetCont->ResetCurrentID();
256 while((jet1 = jetCont->GetNextAcceptJet())) {
257 Double_t fraction = jetCont->GetFractionSharedPt(jet1);
258 if(fMinFractionShared>0. && fraction<fMinFractionShared) continue;
260 Double_t ptJet1 = jet1->Pt() - GetRhoVal(fContainerBase)*jet1->Area();//jetCont->GetJetPtCorr(jetCont->GetCurrentID());//
261 fh2PtJet1VsLeadPtAllSel[fCentBin]->Fill(ptJet1,jet1->MaxTrackPt());
262 fh2PtVsMassJet1All[fCentBin]->Fill(ptJet1,jet1->M());
263 fpPtVsMassJet1All[fCentBin]->Fill(ptJet1,jet1->M());
264 fh2MassVsAreaJet1All[fCentBin]->Fill(jet1->M(),jet1->Area());
265 fh2MassVsNConstJet1All[fCentBin]->Fill(jet1->M(),jet1->GetNumberOfConstituents());
267 if(jet1->GetTagStatus()<1 || !jet1->GetTaggedJet())
270 fh2PtJet1VsLeadPtTagged[fCentBin]->Fill(ptJet1,jet1->MaxTrackPt());
271 fh2PtVsMassJet1Tagged[fCentBin]->Fill(ptJet1,jet1->M());
272 fpPtVsMassJet1Tagged[fCentBin]->Fill(ptJet1,jet1->M());
273 fh2MassVsAreaJet1Tagged[fCentBin]->Fill(jet1->M(),jet1->Area());
274 fh2MassVsNConstJet1Tagged[fCentBin]->Fill(jet1->M(),jet1->GetNumberOfConstituents());
276 Double_t Et2 = jet1->M()*jet1->M() + jet1->Pt()*jet1->Pt();
277 Double_t Et = 0.; Double_t massOverEtR = 0.;
278 if(Et2>0.) Et = TMath::Sqrt(Et2);
279 if((Et*jetCont->GetJetRadius())>0.)
280 massOverEtR = jet1->M()/(Et*jetCont->GetJetRadius());
281 fh2EtMassOverEtRSq[fCentBin]->Fill(Et,massOverEtR*massOverEtR);
288 //________________________________________________________________________
289 Bool_t AliAnalysisTaskEmcalJetMass::RetrieveEventObjects() {
291 // retrieve event objects
294 if (!AliAnalysisTaskEmcalJet::RetrieveEventObjects())
301 //_______________________________________________________________________
302 void AliAnalysisTaskEmcalJetMass::Terminate(Option_t *)
304 // Called once at the end of the analysis.