3 // Base class for rho calculation.
4 // Calculates parameterized rho for given centrality independent of input.
11 #include <TClonesArray.h>
14 #include "AliRhoParameter.h"
15 #include "AliEmcalJet.h"
17 #include "AliAnalysisTaskRhoBase.h"
19 ClassImp(AliAnalysisTaskRhoBase)
21 //________________________________________________________________________
22 AliAnalysisTaskRhoBase::AliAnalysisTaskRhoBase() :
23 AliAnalysisTaskEmcalJet("AliAnalysisTaskRhoBase", kFALSE),
26 fCompareRhoScaledName(),
29 fInEventSigmaRho(35.83),
30 fAttachToEvent(kTRUE),
35 fHistJetAreavsCent(0),
37 fHistJetPtvsNtrack(0),
38 fHistJetAreavsNtrack(0),
41 fHistRhoScaledvsCent(0),
42 fHistDeltaRhovsCent(0),
43 fHistDeltaRhoScalevsCent(0),
45 fHistRhoScaledvsNtrack(0),
46 fHistDeltaRhovsNtrack(0),
47 fHistDeltaRhoScalevsNtrack(0),
48 fHistRhovsNcluster(0),
49 fHistRhoScaledvsNcluster(0)
53 for (Int_t i = 0; i < 4; i++) {
54 fHistJetNconstVsPt[i] = 0;
56 for (Int_t i = 0; i < 12; i++) {
57 fHistNjUEoverNjVsNj[i] = 0;
61 //________________________________________________________________________
62 AliAnalysisTaskRhoBase::AliAnalysisTaskRhoBase(const char *name, Bool_t histo) :
63 AliAnalysisTaskEmcalJet(name, histo),
66 fCompareRhoScaledName(),
69 fInEventSigmaRho(35.83),
70 fAttachToEvent(kTRUE),
75 fHistJetAreavsCent(0),
77 fHistJetPtvsNtrack(0),
78 fHistJetAreavsNtrack(0),
81 fHistRhoScaledvsCent(0),
82 fHistDeltaRhovsCent(0),
83 fHistDeltaRhoScalevsCent(0),
85 fHistRhoScaledvsNtrack(0),
86 fHistDeltaRhovsNtrack(0),
87 fHistDeltaRhoScalevsNtrack(0),
88 fHistRhovsNcluster(0),
89 fHistRhoScaledvsNcluster(0)
93 for (Int_t i = 0; i < 4; i++) {
94 fHistJetNconstVsPt[i] = 0;
96 for (Int_t i = 0; i < 12; i++) {
97 fHistNjUEoverNjVsNj[i] = 0;
100 SetMakeGeneralHistograms(histo);
103 //________________________________________________________________________
104 void AliAnalysisTaskRhoBase::UserCreateOutputObjects()
106 // User create output objects, called at the beginning of the analysis.
111 AliAnalysisTaskEmcalJet::UserCreateOutputObjects();
113 fHistRhovsCent = new TH2F("RhovsCent", "RhovsCent", 101, -1, 100, fNbins, fMinBinPt, fMaxBinPt*2);
114 fOutput->Add(fHistRhovsCent);
116 if (!fTracksName.IsNull()) {
117 fHistRhovsNtrack = new TH2F("RhovsNtrack", "RhovsNtrack", 125, 0, 4000, fNbins, fMinBinPt, fMaxBinPt*2);
118 fOutput->Add(fHistRhovsNtrack);
121 if (!fCaloName.IsNull()) {
122 fHistRhovsNcluster = new TH2F("RhovsNcluster", "RhovsNcluster", 50, 0, 1500, fNbins, fMinBinPt, fMaxBinPt*2);
123 fOutput->Add(fHistRhovsNcluster);
126 if (!fJetsName.IsNull()) {
127 fHistJetPtvsCent = new TH2F("JetPtvsCent", "JetPtvsCent", 101, -1, 100, fNbins, fMinBinPt, fMaxBinPt);
128 fHistJetAreavsCent = new TH2F("JetAreavsCent", "JetAreavsCent", 101, -1, 100, 30, 0, fJetRadius * fJetRadius * TMath::Pi() * 3);
129 fHistNjetvsCent = new TH2F("NjetvsCent", "NjetvsCent", 101, -1, 100, 150, -0.5, 149.5);
131 fOutput->Add(fHistJetPtvsCent);
132 fOutput->Add(fHistJetAreavsCent);
133 fOutput->Add(fHistNjetvsCent);
135 if (!fTracksName.IsNull()) {
136 fHistJetPtvsNtrack = new TH2F("JetPtvsNtrack", "JetPtvsNtrack", 125, 0, 4000, fNbins, fMinBinPt, fMaxBinPt);
137 fHistJetAreavsNtrack = new TH2F("JetAreavsNtrack", "JetAreavsNtrack", 125, 0, 4000, 30, 0, fJetRadius * fJetRadius * TMath::Pi() * 3);
138 fHistNjetvsNtrack = new TH2F("NjetvsNtrack", "rNjetvsNtrack", 125, 0, 4000, 150, -0.5, 149.5);
140 fOutput->Add(fHistJetPtvsNtrack);
141 fOutput->Add(fHistJetAreavsNtrack);
142 fOutput->Add(fHistNjetvsNtrack);
147 for (Int_t i = 0; i < 4; i++) {
148 name = Form("fHistJetNconstVsPt_%d",i);
149 fHistJetNconstVsPt[i] = new TH2F(name, name, 150, -0.5, 149.5, fNbins, fMinBinPt, fMaxBinPt);
150 fHistJetNconstVsPt[i]->GetXaxis()->SetTitle("# of constituents");
151 fHistJetNconstVsPt[i]->GetYaxis()->SetTitle("p_{T,jet} (GeV/c)");
152 fOutput->Add(fHistJetNconstVsPt[i]);
154 for (Int_t j = 0; j < 3; j++) {
155 name = Form("NjUEoverNjVsNj_%d_%d",i,j+1);
156 fHistNjUEoverNjVsNj[i*3+j] = new TH2F(name, name, 150, -0.5, 149.5, 120, 0.01, 1.21);
157 fHistNjUEoverNjVsNj[i*3+j]->GetXaxis()->SetTitle("N_{jet}");
158 fHistNjUEoverNjVsNj[i*3+j]->GetYaxis()->SetTitle("N_{jet_{UE}} / N_{jet}");
159 fOutput->Add(fHistNjUEoverNjVsNj[i*3+j]);
164 if (!fCompareRhoName.IsNull()) {
165 fHistDeltaRhovsCent = new TH2F("DeltaRhovsCent", "DetlaRhovsCent", 101, -1, 100, fNbins, -fMaxBinPt, fMaxBinPt);
166 fOutput->Add(fHistDeltaRhovsCent);
167 if (!fTracksName.IsNull()) {
168 fHistDeltaRhovsNtrack = new TH2F("DeltaRhovsNtrack", "DeltaRhovsNtrack", 125, 0, 4000, fNbins, -fMaxBinPt, fMaxBinPt);
169 fOutput->Add(fHistDeltaRhovsNtrack);
173 if (fScaleFunction) {
174 fHistRhoScaledvsCent = new TH2F("RhoScaledvsCent", "RhoScalevsCent", 101, -1, 100, fNbins, fMinBinPt , fMaxBinPt*2);
175 fOutput->Add(fHistRhoScaledvsCent);
177 if (!fTracksName.IsNull()) {
178 fHistRhoScaledvsNtrack = new TH2F("RhoScaledvsNtrack", "RhoScaledvsNtrack", 125, 0, 4000, fNbins, fMinBinPt, fMaxBinPt*2);
179 fOutput->Add(fHistRhoScaledvsNtrack);
182 if (!fCaloName.IsNull()) {
183 fHistRhoScaledvsNcluster = new TH2F("RhoScaledvsNcluster", "RhoScaledvsNcluster", 50, 0, 1500, fNbins, fMinBinPt, fMaxBinPt*2);
184 fOutput->Add(fHistRhoScaledvsNcluster);
187 if (!fCompareRhoScaledName.IsNull()) {
188 fHistDeltaRhoScalevsCent = new TH2F("DeltaRhoScalevsCent", "DeltaRhoScalevsCent", 101, -1, 100, fNbins, -fMaxBinPt, fMaxBinPt);
189 fOutput->Add(fHistDeltaRhoScalevsCent);
191 if (!fTracksName.IsNull()) {
192 fHistDeltaRhoScalevsNtrack = new TH2F("DeltaRhoScalevsNtrack", "DeltaRhoScalevsNtrack", 125, 0, 4000, fNbins, -fMaxBinPt, fMaxBinPt);
193 fOutput->Add(fHistDeltaRhoScalevsNtrack);
199 //________________________________________________________________________
200 Bool_t AliAnalysisTaskRhoBase::Run()
204 Double_t rho = GetRhoFactor(fCent);
207 if (fScaleFunction) {
208 Double_t rhoScaled = rho * GetScaleFactor(fCent);
209 fRhoScaled->SetVal(rhoScaled);
215 //________________________________________________________________________
216 Bool_t AliAnalysisTaskRhoBase::FillHistograms()
224 Ntracks = fTracks->GetEntriesFast();
226 Nclusters = fCaloClusters->GetEntriesFast();
229 Int_t Njets = fJets->GetEntries();
231 Int_t NjetUE1Sigma = 0;
232 Int_t NjetUE2Sigma = 0;
233 Int_t NjetUE3Sigma = 0;
234 Double_t rhoPlus1Sigma = fRho->GetVal() + fInEventSigmaRho;
235 Double_t rhoPlus2Sigma = fRho->GetVal() + 2*fInEventSigmaRho;
236 Double_t rhoPlus3Sigma = fRho->GetVal() + 3*fInEventSigmaRho;
238 for (Int_t i = 0; i < Njets; ++i) {
240 AliEmcalJet *jet = static_cast<AliEmcalJet*>(fJets->At(i));
242 AliError(Form("%s: Could not receive jet %d", GetName(), i));
249 fHistJetPtvsCent->Fill(fCent, jet->Pt());
250 fHistJetAreavsCent->Fill(fCent, jet->Area());
253 fHistJetPtvsNtrack->Fill(Ntracks, jet->Pt());
254 fHistJetAreavsNtrack->Fill(Ntracks, jet->Area());
257 fHistJetNconstVsPt[fCentBin]->Fill(jet->GetNumberOfConstituents(), jet->Pt());
259 if (jet->Pt() < rhoPlus1Sigma * jet->Area())
262 if (jet->Pt() < rhoPlus2Sigma * jet->Area())
265 if (jet->Pt() < rhoPlus3Sigma * jet->Area())
272 fHistNjUEoverNjVsNj[fCentBin*3 ]->Fill(NjetAcc,1.*NjetUE1Sigma/NjetAcc);
273 fHistNjUEoverNjVsNj[fCentBin*3+1]->Fill(NjetAcc,1.*NjetUE2Sigma/NjetAcc);
274 fHistNjUEoverNjVsNj[fCentBin*3+2]->Fill(NjetAcc,1.*NjetUE3Sigma/NjetAcc);
277 fHistNjetvsCent->Fill(fCent, NjetAcc);
279 fHistNjetvsNtrack->Fill(Ntracks, NjetAcc);
282 fHistRhovsCent->Fill(fCent, fRho->GetVal());
285 fHistRhovsNtrack->Fill(Ntracks, fRho->GetVal());
287 fHistRhovsNcluster->Fill(Nclusters, fRho->GetVal());
289 fHistDeltaRhovsCent->Fill(fCent, fRho->GetVal() - fCompareRho->GetVal());
291 fHistDeltaRhovsNtrack->Fill(Ntracks, fRho->GetVal() - fCompareRho->GetVal());
295 fHistRhoScaledvsCent->Fill(fCent, fRhoScaled->GetVal());
297 fHistRhoScaledvsNtrack->Fill(Ntracks, fRhoScaled->GetVal());
299 fHistRhoScaledvsNcluster->Fill(Nclusters, fRhoScaled->GetVal());
300 if (fCompareRhoScaled) {
301 fHistDeltaRhoScalevsCent->Fill(fCent, fRhoScaled->GetVal() - fCompareRhoScaled->GetVal());
303 fHistDeltaRhoScalevsNtrack->Fill(Ntracks, fRhoScaled->GetVal() - fCompareRhoScaled->GetVal());
311 //________________________________________________________________________
312 void AliAnalysisTaskRhoBase::ExecOnce()
314 // Init the analysis.
317 fRho = new AliRhoParameter(fRhoName, 0);
319 if (fAttachToEvent) {
320 if (!(InputEvent()->FindListObject(fRhoName))) {
321 InputEvent()->AddObject(fRho);
323 AliFatal(Form("%s: Container with same name %s already present. Aborting", GetName(), fRhoName.Data()));
328 if (fScaleFunction && !fRhoScaled) {
329 fRhoScaled = new AliRhoParameter(fRhoScaledName, 0);
330 if (fAttachToEvent) {
331 if (!(InputEvent()->FindListObject(fRhoScaledName))) {
332 InputEvent()->AddObject(fRhoScaled);
334 AliFatal(Form("%s: Container with same name %s already present. Aborting", GetName(), fRhoScaledName.Data()));
341 if (!fCompareRhoName.IsNull() && !fCompareRho) {
342 fCompareRho = dynamic_cast<AliRhoParameter*>(InputEvent()->FindListObject(fCompareRhoName));
344 AliWarning(Form("%s: Could not retrieve rho %s!", GetName(), fCompareRhoName.Data()));
348 if (!fCompareRhoScaledName.IsNull() && !fCompareRhoScaled) {
349 fCompareRhoScaled = dynamic_cast<AliRhoParameter*>(InputEvent()->FindListObject(fCompareRhoScaledName));
350 if (!fCompareRhoScaled) {
351 AliWarning(Form("%s: Could not retrieve rho %s!", GetName(), fCompareRhoScaledName.Data()));
355 AliAnalysisTaskEmcalJet::ExecOnce();
358 //________________________________________________________________________
359 Double_t AliAnalysisTaskRhoBase::GetRhoFactor(Double_t cent)
361 // Return rho per centrality.
365 rho = fRhoFunction->Eval(cent);
369 //________________________________________________________________________
370 Double_t AliAnalysisTaskRhoBase::GetScaleFactor(Double_t cent)
376 scale = fScaleFunction->Eval(cent);