]>
Commit | Line | Data |
---|---|---|
0f67a57c | 1 | /* $Id$ */ |
2 | ||
3 | #include "AlidNdEtaCorrectionTask.h" | |
4 | ||
54b096ef | 5 | #include <TROOT.h> |
0f67a57c | 6 | #include <TCanvas.h> |
7 | #include <TChain.h> | |
8 | #include <TFile.h> | |
9 | #include <TH1F.h> | |
54b096ef | 10 | #include <TH2F.h> |
1c15d51a | 11 | #include <TH3F.h> |
54b096ef | 12 | #include <TProfile.h> |
0f67a57c | 13 | #include <TParticle.h> |
1c15d51a | 14 | #include <TParticlePDG.h> |
0f67a57c | 15 | |
16 | #include <AliLog.h> | |
17 | #include <AliESDVertex.h> | |
18 | #include <AliESDEvent.h> | |
19 | #include <AliStack.h> | |
20 | #include <AliHeader.h> | |
21 | #include <AliGenEventHeader.h> | |
22 | #include <AliMultiplicity.h> | |
23 | #include <AliAnalysisManager.h> | |
24 | #include <AliMCEventHandler.h> | |
25 | #include <AliMCEvent.h> | |
26 | #include <AliESDInputHandler.h> | |
27 | ||
745d6088 | 28 | #include "AliESDtrackCuts.h" |
0f67a57c | 29 | #include "AliPWG0Helper.h" |
30 | //#include "AliCorrection.h" | |
31 | //#include "AliCorrectionMatrix3D.h" | |
32 | #include "dNdEta/dNdEtaAnalysis.h" | |
33 | #include "dNdEta/AlidNdEtaCorrection.h" | |
54b096ef | 34 | #include "AliVertexerTracks.h" |
0f67a57c | 35 | |
36 | ClassImp(AlidNdEtaCorrectionTask) | |
37 | ||
1c15d51a | 38 | AlidNdEtaCorrectionTask::AlidNdEtaCorrectionTask() : |
39 | AliAnalysisTask(), | |
40 | fESD(0), | |
41 | fOutput(0), | |
42 | fOption(), | |
43 | fAnalysisMode(AliPWG0Helper::kTPC), | |
44 | fTrigger(AliPWG0Helper::kMB1), | |
45 | fSignMode(0), | |
46 | fOnlyPrimaries(kFALSE), | |
47 | fEsdTrackCuts(0), | |
48 | fdNdEtaCorrection(0), | |
49 | fdNdEtaAnalysisMC(0), | |
50 | fdNdEtaAnalysisESD(0), | |
51 | fPIDParticles(0), | |
52 | fPIDTracks(0), | |
53 | fVertexCorrelation(0), | |
54 | fVertexProfile(0), | |
55 | fVertexShift(0), | |
56 | fVertexShiftNorm(0), | |
57 | fEtaCorrelation(0), | |
58 | fEtaProfile(0), | |
59 | fEtaResolution(0), | |
60 | fpTResolution(0), | |
61 | fEsdTrackCutsPrim(0), | |
62 | fEsdTrackCutsSec(0), | |
63 | fTemp1(0), | |
64 | fTemp2(0), | |
65 | fMultAll(0), | |
66 | fMultTr(0), | |
67 | fMultVtx(0), | |
68 | fEventStats(0) | |
69 | { | |
70 | // | |
71 | // Constructor. Initialization of pointers | |
72 | // | |
73 | ||
74 | for (Int_t i=0; i<2; i++) | |
75 | fdNdEtaCorrectionProcessType[i] = 0; | |
76 | ||
77 | for (Int_t i=0; i<8; i++) | |
78 | fDeltaPhi[i] = 0; | |
79 | } | |
80 | ||
0f67a57c | 81 | AlidNdEtaCorrectionTask::AlidNdEtaCorrectionTask(const char* opt) : |
82 | AliAnalysisTask("AlidNdEtaCorrectionTask", ""), | |
83 | fESD(0), | |
84 | fOutput(0), | |
85 | fOption(opt), | |
770a1f1d | 86 | fAnalysisMode(AliPWG0Helper::kTPC), |
0fc41645 | 87 | fTrigger(AliPWG0Helper::kMB1), |
0f67a57c | 88 | fSignMode(0), |
3d7758c1 | 89 | fOnlyPrimaries(kFALSE), |
0f67a57c | 90 | fEsdTrackCuts(0), |
91 | fdNdEtaCorrection(0), | |
92 | fdNdEtaAnalysisMC(0), | |
93 | fdNdEtaAnalysisESD(0), | |
94 | fPIDParticles(0), | |
96fcc8a7 | 95 | fPIDTracks(0), |
54b096ef | 96 | fVertexCorrelation(0), |
97 | fVertexProfile(0), | |
1c15d51a | 98 | fVertexShift(0), |
54b096ef | 99 | fVertexShiftNorm(0), |
3d7758c1 | 100 | fEtaCorrelation(0), |
101 | fEtaProfile(0), | |
1c15d51a | 102 | fEtaResolution(0), |
103 | fpTResolution(0), | |
104 | fEsdTrackCutsPrim(0), | |
105 | fEsdTrackCutsSec(0), | |
106 | fTemp1(0), | |
107 | fTemp2(0), | |
4ef701a0 | 108 | fMultAll(0), |
109 | fMultTr(0), | |
50ec344d | 110 | fMultVtx(0), |
111 | fEventStats(0) | |
0f67a57c | 112 | { |
113 | // | |
114 | // Constructor. Initialization of pointers | |
115 | // | |
116 | ||
117 | // Define input and output slots here | |
118 | DefineInput(0, TChain::Class()); | |
119 | DefineOutput(0, TList::Class()); | |
96fcc8a7 | 120 | |
121 | for (Int_t i=0; i<2; i++) | |
122 | fdNdEtaCorrectionProcessType[i] = 0; | |
4ef701a0 | 123 | |
3d7758c1 | 124 | for (Int_t i=0; i<8; i++) |
4ef701a0 | 125 | fDeltaPhi[i] = 0; |
0f67a57c | 126 | } |
127 | ||
128 | AlidNdEtaCorrectionTask::~AlidNdEtaCorrectionTask() | |
129 | { | |
130 | // | |
131 | // Destructor | |
132 | // | |
133 | ||
134 | // histograms are in the output list and deleted when the output | |
135 | // list is deleted by the TSelector dtor | |
136 | ||
137 | if (fOutput) { | |
138 | delete fOutput; | |
139 | fOutput = 0; | |
140 | } | |
141 | } | |
142 | ||
143 | //________________________________________________________________________ | |
144 | void AlidNdEtaCorrectionTask::ConnectInputData(Option_t *) | |
145 | { | |
146 | // Connect ESD | |
147 | // Called once | |
148 | ||
149 | Printf("AlidNdEtaCorrectionTask::ConnectInputData called"); | |
150 | ||
151 | TTree* tree = dynamic_cast<TTree*> (GetInputData(0)); | |
152 | if (!tree) { | |
153 | Printf("ERROR: Could not read tree from input slot 0"); | |
154 | } else { | |
155 | // Disable all branches and enable only the needed ones | |
567160d6 | 156 | tree->SetBranchStatus("*", 0); |
0f67a57c | 157 | |
567160d6 | 158 | tree->SetBranchStatus("AliESDHeader.*", 1); |
159 | tree->SetBranchStatus("*Vertex*", 1); | |
0f67a57c | 160 | |
770a1f1d | 161 | if (fAnalysisMode == AliPWG0Helper::kSPD) |
567160d6 | 162 | tree->SetBranchStatus("AliMultiplicity*", 1); |
0f67a57c | 163 | |
770a1f1d | 164 | if (fAnalysisMode == AliPWG0Helper::kTPC || fAnalysisMode == AliPWG0Helper::kTPCITS) { |
567160d6 | 165 | //AliESDtrackCuts::EnableNeededBranches(tree); |
166 | tree->SetBranchStatus("Tracks*", 1); | |
0f67a57c | 167 | } |
168 | ||
0f67a57c | 169 | AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()); |
170 | ||
171 | if (!esdH) { | |
172 | Printf("ERROR: Could not get ESDInputHandler"); | |
173 | } else | |
174 | fESD = esdH->GetEvent(); | |
175 | } | |
96fcc8a7 | 176 | |
177 | // disable info messages of AliMCEvent (per event) | |
178 | AliLog::SetClassDebugLevel("AliMCEvent", AliLog::kWarning - AliLog::kDebug + 1); | |
0f67a57c | 179 | } |
180 | ||
181 | void AlidNdEtaCorrectionTask::CreateOutputObjects() | |
182 | { | |
183 | // create result objects and add to output list | |
184 | ||
185 | if (fOption.Contains("only-positive")) | |
186 | { | |
187 | Printf("INFO: Processing only positive particles."); | |
188 | fSignMode = 1; | |
189 | } | |
190 | else if (fOption.Contains("only-negative")) | |
191 | { | |
192 | Printf("INFO: Processing only negative particles."); | |
193 | fSignMode = -1; | |
194 | } | |
195 | ||
0f67a57c | 196 | fOutput = new TList; |
197 | fOutput->SetOwner(); | |
198 | ||
770a1f1d | 199 | fdNdEtaCorrection = new AlidNdEtaCorrection("dndeta_correction", "dndeta_correction", fAnalysisMode); |
0f67a57c | 200 | fOutput->Add(fdNdEtaCorrection); |
201 | ||
1c15d51a | 202 | fPIDParticles = new TH1F("fPIDParticles", "PID of generated primary particles", 10001, -5000.5, 5000.5); |
0f67a57c | 203 | fOutput->Add(fPIDParticles); |
204 | ||
1c15d51a | 205 | fPIDTracks = new TH1F("fPIDTracks", "MC PID of reconstructed tracks", 10001, -5000.5, 5000.5); |
0f67a57c | 206 | fOutput->Add(fPIDTracks); |
207 | ||
770a1f1d | 208 | fdNdEtaAnalysisMC = new dNdEtaAnalysis("dndetaMC", "dndetaMC", fAnalysisMode); |
0f67a57c | 209 | fOutput->Add(fdNdEtaAnalysisMC); |
210 | ||
770a1f1d | 211 | fdNdEtaAnalysisESD = new dNdEtaAnalysis("dndetaESD", "dndetaESD", fAnalysisMode); |
0f67a57c | 212 | fOutput->Add(fdNdEtaAnalysisESD); |
96fcc8a7 | 213 | |
1c15d51a | 214 | if (fEsdTrackCuts) |
215 | { | |
216 | fEsdTrackCutsPrim = dynamic_cast<AliESDtrackCuts*> (fEsdTrackCuts->Clone("fEsdTrackCutsPrim")); | |
217 | fEsdTrackCutsSec = dynamic_cast<AliESDtrackCuts*> (fEsdTrackCuts->Clone("fEsdTrackCutsSec")); | |
218 | fOutput->Add(fEsdTrackCutsPrim); | |
219 | fOutput->Add(fEsdTrackCutsSec); | |
220 | } | |
221 | ||
96fcc8a7 | 222 | if (fOption.Contains("process-types")) { |
54b096ef | 223 | fdNdEtaCorrectionProcessType[0] = new AlidNdEtaCorrection("dndeta_correction_ND", "dndeta_correction_ND", fAnalysisMode); |
224 | fdNdEtaCorrectionProcessType[1] = new AlidNdEtaCorrection("dndeta_correction_SD", "dndeta_correction_SD", fAnalysisMode); | |
225 | fdNdEtaCorrectionProcessType[2] = new AlidNdEtaCorrection("dndeta_correction_DD", "dndeta_correction_DD", fAnalysisMode); | |
96fcc8a7 | 226 | |
227 | fOutput->Add(fdNdEtaCorrectionProcessType[0]); | |
228 | fOutput->Add(fdNdEtaCorrectionProcessType[1]); | |
229 | fOutput->Add(fdNdEtaCorrectionProcessType[2]); | |
230 | } | |
231 | ||
1c15d51a | 232 | /* |
233 | fTemp1 = new TH3F("fTemp1", "fTemp1 prim;b0;b1;1/pT", 200, 0, 10, 200, 0, 10, 200, 0, 10); | |
234 | fOutput->Add(fTemp1); | |
235 | fTemp2 = new TH3F("fTemp2", "fTemp2 sec;b0;b1;1/pT", 200, 0, 10, 200, 0, 10, 200, 0, 10); | |
236 | fOutput->Add(fTemp2); | |
237 | */ | |
54b096ef | 238 | |
0fc41645 | 239 | fVertexCorrelation = new TH2F("fVertexCorrelation", "fVertexCorrelation;MC z-vtx;ESD z-vtx", 120, -30, 30, 120, -30, 30); |
240 | fOutput->Add(fVertexCorrelation); | |
54b096ef | 241 | fVertexProfile = new TProfile("fVertexProfile", "fVertexProfile;MC z-vtx;MC z-vtx - ESD z-vtx", 40, -20, 20); |
0fc41645 | 242 | fOutput->Add(fVertexProfile); |
1c15d51a | 243 | fVertexShift = new TH1F("fVertexShift", "fVertexShift;(MC z-vtx - ESD z-vtx);Entries", 201, -2, 2); |
244 | fOutput->Add(fVertexShift); | |
54b096ef | 245 | fVertexShiftNorm = new TH1F("fVertexShiftNorm", "fVertexShiftNorm;(MC z-vtx - ESD z-vtx) / #sigma_{ESD z-vtx};Entries", 200, -100, 100); |
0fc41645 | 246 | fOutput->Add(fVertexShiftNorm); |
247 | ||
3d7758c1 | 248 | fEtaCorrelation = new TH2F("fEtaCorrelation", "fEtaCorrelation;MC #eta;ESD #eta", 120, -3, 3, 120, -3, 3); |
0fc41645 | 249 | fOutput->Add(fEtaCorrelation); |
3d7758c1 | 250 | fEtaProfile = new TProfile("fEtaProfile", "fEtaProfile;MC #eta;MC #eta - ESD #eta", 120, -3, 3); |
0fc41645 | 251 | fOutput->Add(fEtaProfile); |
1c15d51a | 252 | fEtaResolution = new TH1F("fEtaResolution", "fEtaResolution;MC #eta - ESD #eta", 201, -0.2, 0.2); |
253 | fOutput->Add(fEtaResolution); | |
254 | ||
255 | fpTResolution = new TH1F("fpTResolution", "fpTResolution;MC p_{T} - ESD p_{T}", 201, -0.2, 0.2); | |
256 | fOutput->Add(fpTResolution); | |
3d7758c1 | 257 | |
4ef701a0 | 258 | fMultAll = new TH1F("fMultAll", "fMultAll", 500, -0.5, 499.5); |
0fc41645 | 259 | fOutput->Add(fMultAll); |
4ef701a0 | 260 | fMultVtx = new TH1F("fMultVtx", "fMultVtx", 500, -0.5, 499.5); |
0fc41645 | 261 | fOutput->Add(fMultVtx); |
4ef701a0 | 262 | fMultTr = new TH1F("fMultTr", "fMultTr", 500, -0.5, 499.5); |
0fc41645 | 263 | fOutput->Add(fMultTr); |
4ef701a0 | 264 | |
3d7758c1 | 265 | for (Int_t i=0; i<8; i++) |
0fc41645 | 266 | { |
3d7758c1 | 267 | fDeltaPhi[i] = new TH1F(Form("fDeltaPhi_%d", i), ";#Delta phi;Entries", 2000, -0.1, 0.1); |
0fc41645 | 268 | fOutput->Add(fDeltaPhi[i]); |
269 | } | |
50ec344d | 270 | |
ea441adf | 271 | fEventStats = new TH2F("fEventStats", "fEventStats;event type;status;count", 104, -3.5, 100.5, 4, -0.5, 3.5); |
0fc41645 | 272 | fOutput->Add(fEventStats); |
ea441adf | 273 | fEventStats->GetXaxis()->SetBinLabel(1, "INEL"); // x = -3 |
274 | fEventStats->GetXaxis()->SetBinLabel(2, "NSD"); // x = -2 | |
275 | for (Int_t i=0; i<100; i++) | |
276 | fEventStats->GetXaxis()->SetBinLabel(4+i, Form("%d", i)); | |
50ec344d | 277 | |
278 | fEventStats->GetYaxis()->SetBinLabel(1, "nothing"); | |
279 | fEventStats->GetYaxis()->SetBinLabel(2, "trg"); | |
280 | fEventStats->GetYaxis()->SetBinLabel(3, "vtx"); | |
281 | fEventStats->GetYaxis()->SetBinLabel(4, "trgvtx"); | |
1c15d51a | 282 | |
283 | if (fEsdTrackCuts) | |
284 | { | |
285 | fEsdTrackCuts->SetName("fEsdTrackCuts"); | |
286 | fOutput->Add(fEsdTrackCuts); | |
287 | } | |
0f67a57c | 288 | } |
289 | ||
290 | void AlidNdEtaCorrectionTask::Exec(Option_t*) | |
291 | { | |
292 | // process the event | |
293 | ||
294 | // Check prerequisites | |
295 | if (!fESD) | |
296 | { | |
297 | AliDebug(AliLog::kError, "ESD branch not available"); | |
298 | return; | |
299 | } | |
300 | ||
3d7758c1 | 301 | if (fOnlyPrimaries) |
302 | Printf("WARNING: Processing only primaries. For systematical studies only!"); | |
303 | ||
0f67a57c | 304 | // trigger definition |
0fc41645 | 305 | Bool_t eventTriggered = AliPWG0Helper::IsEventTriggered(fESD->GetTriggerMask(), fTrigger); |
54b096ef | 306 | |
3d7758c1 | 307 | if (!eventTriggered) |
308 | Printf("No trigger"); | |
309 | ||
54b096ef | 310 | // post the data already here |
311 | PostData(0, fOutput); | |
312 | ||
313 | // MC info | |
314 | AliMCEventHandler* eventHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler()); | |
315 | if (!eventHandler) { | |
316 | Printf("ERROR: Could not retrieve MC event handler"); | |
317 | return; | |
318 | } | |
319 | ||
320 | AliMCEvent* mcEvent = eventHandler->MCEvent(); | |
321 | if (!mcEvent) { | |
322 | Printf("ERROR: Could not retrieve MC event"); | |
323 | return; | |
324 | } | |
325 | ||
326 | AliStack* stack = mcEvent->Stack(); | |
327 | if (!stack) | |
328 | { | |
329 | AliDebug(AliLog::kError, "Stack not available"); | |
330 | return; | |
331 | } | |
332 | ||
333 | AliHeader* header = mcEvent->Header(); | |
334 | if (!header) | |
335 | { | |
336 | AliDebug(AliLog::kError, "Header not available"); | |
337 | return; | |
338 | } | |
339 | ||
6b62a9c7 | 340 | // get process type; |
341 | Int_t processType = AliPWG0Helper::GetEventProcessType(header); | |
342 | AliDebug(AliLog::kDebug+1, Form("Found process type %d", processType)); | |
54b096ef | 343 | |
6b62a9c7 | 344 | if (processType == AliPWG0Helper::kInvalidProcess) |
ea441adf | 345 | AliDebug(AliLog::kError, "Unknown process type."); |
54b096ef | 346 | |
347 | // get the MC vertex | |
348 | AliGenEventHeader* genHeader = header->GenEventHeader(); | |
349 | TArrayF vtxMC(3); | |
350 | genHeader->PrimaryVertex(vtxMC); | |
351 | ||
352 | // get the ESD vertex | |
353 | const AliESDVertex* vtxESD = AliPWG0Helper::GetVertex(fESD, fAnalysisMode); | |
0f67a57c | 354 | |
770a1f1d | 355 | Bool_t eventVertex = kFALSE; |
54b096ef | 356 | Double_t vtx[3]; |
357 | if (vtxESD) | |
358 | { | |
359 | vtxESD->GetXYZ(vtx); | |
770a1f1d | 360 | eventVertex = kTRUE; |
54b096ef | 361 | |
362 | Double_t diff = vtxMC[2] - vtx[2]; | |
1c15d51a | 363 | fVertexShift->Fill(diff); |
364 | if (vtxESD->GetZRes() > 0) | |
54b096ef | 365 | fVertexShiftNorm->Fill(diff / vtxESD->GetZRes()); |
567160d6 | 366 | } |
54b096ef | 367 | else |
368 | Printf("No vertex found"); | |
0f67a57c | 369 | |
50ec344d | 370 | // fill process type |
371 | Int_t biny = (Int_t) eventTriggered + 2 * (Int_t) eventVertex; | |
372 | // INEL | |
ea441adf | 373 | fEventStats->Fill(-3, biny); |
567160d6 | 374 | // NSD |
ea441adf | 375 | if (processType != AliPWG0Helper::kSD) |
376 | fEventStats->Fill(-2, biny); | |
50ec344d | 377 | |
0f67a57c | 378 | // create list of (label, eta, pt) tuples |
379 | Int_t inputCount = 0; | |
380 | Int_t* labelArr = 0; | |
3d7758c1 | 381 | Int_t* labelArr2 = 0; // only for case of SPD |
0f67a57c | 382 | Float_t* etaArr = 0; |
383 | Float_t* ptArr = 0; | |
4ef701a0 | 384 | Float_t* deltaPhiArr = 0; |
770a1f1d | 385 | if (fAnalysisMode == AliPWG0Helper::kSPD) |
0f67a57c | 386 | { |
387 | // get tracklets | |
388 | const AliMultiplicity* mult = fESD->GetMultiplicity(); | |
389 | if (!mult) | |
390 | { | |
391 | AliDebug(AliLog::kError, "AliMultiplicity not available"); | |
392 | return; | |
393 | } | |
394 | ||
395 | labelArr = new Int_t[mult->GetNumberOfTracklets()]; | |
3d7758c1 | 396 | labelArr2 = new Int_t[mult->GetNumberOfTracklets()]; |
0f67a57c | 397 | etaArr = new Float_t[mult->GetNumberOfTracklets()]; |
398 | ptArr = new Float_t[mult->GetNumberOfTracklets()]; | |
4ef701a0 | 399 | deltaPhiArr = new Float_t[mult->GetNumberOfTracklets()]; |
0f67a57c | 400 | |
1c15d51a | 401 | // get multiplicity from SPD tracklets |
0f67a57c | 402 | for (Int_t i=0; i<mult->GetNumberOfTracklets(); ++i) |
403 | { | |
404 | //printf("%d %f %f %f\n", i, mult->GetTheta(i), mult->GetPhi(i), mult->GetDeltaPhi(i)); | |
405 | ||
4ef701a0 | 406 | Float_t deltaPhi = mult->GetDeltaPhi(i); |
407 | // prevent values to be shifted by 2 Pi() | |
408 | if (deltaPhi < -TMath::Pi()) | |
409 | deltaPhi += TMath::Pi() * 2; | |
410 | if (deltaPhi > TMath::Pi()) | |
411 | deltaPhi -= TMath::Pi() * 2; | |
412 | ||
413 | if (TMath::Abs(deltaPhi) > 1) | |
414 | printf("WARNING: Very high Delta Phi: %d %f %f %f\n", i, mult->GetTheta(i), mult->GetPhi(i), deltaPhi); | |
0f67a57c | 415 | |
3d7758c1 | 416 | if (fOnlyPrimaries) |
417 | if (mult->GetLabel(i, 0) < 0 || mult->GetLabel(i, 0) != mult->GetLabel(i, 1) || !stack->IsPhysicalPrimary(mult->GetLabel(i, 0))) | |
418 | continue; | |
419 | ||
0f67a57c | 420 | etaArr[inputCount] = mult->GetEta(i); |
3d7758c1 | 421 | labelArr[inputCount] = mult->GetLabel(i, 0); |
422 | labelArr2[inputCount] = mult->GetLabel(i, 1); | |
0f67a57c | 423 | ptArr[inputCount] = 0; // no pt for tracklets |
4ef701a0 | 424 | deltaPhiArr[inputCount] = deltaPhi; |
0f67a57c | 425 | ++inputCount; |
426 | } | |
427 | } | |
770a1f1d | 428 | else if (fAnalysisMode == AliPWG0Helper::kTPC || fAnalysisMode == AliPWG0Helper::kTPCITS) |
0f67a57c | 429 | { |
430 | if (!fEsdTrackCuts) | |
431 | { | |
432 | AliDebug(AliLog::kError, "fESDTrackCuts not available"); | |
433 | return; | |
434 | } | |
435 | ||
436 | // get multiplicity from ESD tracks | |
1c15d51a | 437 | TObjArray* list = fEsdTrackCuts->GetAcceptedTracks(fESD, (fAnalysisMode == AliPWG0Helper::kTPC)); |
0f67a57c | 438 | Int_t nGoodTracks = list->GetEntries(); |
439 | ||
54b096ef | 440 | Printf("Accepted %d tracks", nGoodTracks); |
441 | ||
0f67a57c | 442 | labelArr = new Int_t[nGoodTracks]; |
3d7758c1 | 443 | labelArr2 = new Int_t[nGoodTracks]; |
0f67a57c | 444 | etaArr = new Float_t[nGoodTracks]; |
445 | ptArr = new Float_t[nGoodTracks]; | |
4ef701a0 | 446 | deltaPhiArr = new Float_t[nGoodTracks]; |
0f67a57c | 447 | |
448 | // loop over esd tracks | |
449 | for (Int_t i=0; i<nGoodTracks; i++) | |
450 | { | |
451 | AliESDtrack* esdTrack = dynamic_cast<AliESDtrack*> (list->At(i)); | |
452 | if (!esdTrack) | |
453 | { | |
454 | AliDebug(AliLog::kError, Form("ERROR: Could not retrieve track %d.", i)); | |
455 | continue; | |
456 | } | |
457 | ||
458 | etaArr[inputCount] = esdTrack->Eta(); | |
459 | labelArr[inputCount] = TMath::Abs(esdTrack->GetLabel()); | |
3d7758c1 | 460 | labelArr2[inputCount] = labelArr[inputCount]; // no second label for tracks |
0f67a57c | 461 | ptArr[inputCount] = esdTrack->Pt(); |
4ef701a0 | 462 | deltaPhiArr[inputCount] = 0; // no delta phi for tracks |
0f67a57c | 463 | ++inputCount; |
464 | } | |
0f67a57c | 465 | |
54b096ef | 466 | delete list; |
1c15d51a | 467 | |
468 | if (eventTriggered && vtxESD) | |
469 | { | |
470 | // collect values for primaries and secondaries | |
471 | for (Int_t iTrack = 0; iTrack < fESD->GetNumberOfTracks(); iTrack++) | |
472 | { | |
473 | AliESDtrack* track = 0; | |
474 | ||
475 | if (fAnalysisMode == AliPWG0Helper::kTPC) | |
476 | track = AliESDtrackCuts::GetTPCOnlyTrack(fESD, iTrack); | |
477 | else if (fAnalysisMode == AliPWG0Helper::kTPCITS) | |
478 | track = fESD->GetTrack(iTrack); | |
479 | ||
480 | if (!track) | |
481 | continue; | |
482 | ||
483 | Int_t label = TMath::Abs(track->GetLabel()); | |
484 | if (!stack->Particle(label) || !stack->Particle(label)->GetPDG()) | |
485 | { | |
486 | Printf("WARNING: No particle for %d", label); | |
487 | if (stack->Particle(label)) | |
488 | stack->Particle(label)->Print(); | |
489 | continue; | |
490 | } | |
491 | ||
492 | if (stack->Particle(label)->GetPDG()->Charge() == 0) | |
493 | continue; | |
494 | ||
495 | if (stack->IsPhysicalPrimary(label)) | |
496 | { | |
497 | // primary | |
498 | if (fEsdTrackCutsPrim->AcceptTrack(track)) { | |
499 | if (track->Pt() > 0) { | |
500 | Float_t b0, b1; | |
501 | track->GetImpactParameters(b0, b1); | |
502 | if (fTemp1) | |
503 | ((TH3*) fTemp1)->Fill(b0, b1, 1.0 / track->Pt()); | |
504 | //fTemp1->Fill(TMath::Sqrt(b0*b0 + b1*b1), 1.0 / track->Pt()); | |
505 | } | |
506 | ||
507 | if (AliESDtrackCuts::GetSigmaToVertex(track) > 900) | |
508 | { | |
509 | Printf("Track %d has nsigma of %f. Printing track and vertex...", iTrack, AliESDtrackCuts::GetSigmaToVertex(track)); | |
510 | Float_t b[2]; | |
511 | Float_t r[3]; | |
512 | track->GetImpactParameters(b, r); | |
513 | Printf("Impact parameter %f %f and resolution: %f %f %f", b[0], b[1], r[0], r[1], r[2]); | |
514 | track->Print(""); | |
515 | if (vtxESD) | |
516 | vtxESD->Print(); | |
517 | } | |
518 | } | |
519 | } | |
520 | else | |
521 | { | |
522 | // secondary | |
523 | if (fEsdTrackCutsSec->AcceptTrack(track)) { | |
524 | if (track->Pt() > 0) { | |
525 | Float_t b0, b1; | |
526 | track->GetImpactParameters(b0, b1); | |
527 | if (fTemp2) | |
528 | ((TH3*) fTemp2)->Fill(b0, b1, 1.0 / track->Pt()); | |
529 | //fTemp2->Fill(TMath::Sqrt(b0*b0 + b1*b1), 1.0 / track->Pt()); | |
530 | } | |
531 | } | |
532 | } | |
533 | ||
534 | // TODO mem leak in the continue statements above | |
535 | if (fAnalysisMode == AliPWG0Helper::kTPC) | |
536 | delete track; | |
537 | } | |
538 | } | |
0f67a57c | 539 | } |
54b096ef | 540 | else |
0f67a57c | 541 | return; |
0f67a57c | 542 | |
54b096ef | 543 | if (eventTriggered && eventVertex) |
0f67a57c | 544 | { |
54b096ef | 545 | fVertexCorrelation->Fill(vtxMC[2], vtx[2]); |
546 | fVertexProfile->Fill(vtxMC[2], vtxMC[2] - vtx[2]); | |
0f67a57c | 547 | } |
548 | ||
0f67a57c | 549 | |
550 | // loop over mc particles | |
551 | Int_t nPrim = stack->GetNprimary(); | |
4ef701a0 | 552 | Int_t nAccepted = 0; |
0f67a57c | 553 | |
554 | for (Int_t iMc = 0; iMc < nPrim; ++iMc) | |
555 | { | |
556 | //Printf("Getting particle %d", iMc); | |
557 | TParticle* particle = stack->Particle(iMc); | |
558 | ||
559 | if (!particle) | |
560 | continue; | |
561 | ||
562 | if (AliPWG0Helper::IsPrimaryCharged(particle, nPrim) == kFALSE) | |
1c15d51a | 563 | { |
564 | //if (TMath::Abs(particle->GetPdgCode()) > 3000 && TMath::Abs(particle->Eta()) < 1.0) | |
565 | // fPIDParticles->Fill(particle->GetPdgCode()); | |
0f67a57c | 566 | continue; |
1c15d51a | 567 | } |
0f67a57c | 568 | |
569 | if (SignOK(particle->GetPDG()) == kFALSE) | |
570 | continue; | |
571 | ||
1c15d51a | 572 | if (fPIDParticles && TMath::Abs(particle->Eta()) < 1.0) |
573 | fPIDParticles->Fill(particle->GetPdgCode()); | |
574 | ||
0f67a57c | 575 | Float_t eta = particle->Eta(); |
576 | Float_t pt = particle->Pt(); | |
577 | ||
578 | fdNdEtaCorrection->FillMCParticle(vtxMC[2], eta, pt, eventTriggered, eventVertex, processType); | |
579 | ||
96fcc8a7 | 580 | if (fdNdEtaCorrectionProcessType[0]) |
581 | { | |
582 | // non diffractive | |
6b62a9c7 | 583 | if (processType==AliPWG0Helper::kND) |
96fcc8a7 | 584 | fdNdEtaCorrectionProcessType[0]->FillMCParticle(vtxMC[2], eta, pt, eventTriggered, eventVertex, processType); |
585 | ||
586 | // single diffractive | |
6b62a9c7 | 587 | if (processType==AliPWG0Helper::kSD) |
96fcc8a7 | 588 | fdNdEtaCorrectionProcessType[1]->FillMCParticle(vtxMC[2], eta, pt, eventTriggered, eventVertex, processType); |
589 | ||
590 | // double diffractive | |
6b62a9c7 | 591 | if (processType==AliPWG0Helper::kDD) |
96fcc8a7 | 592 | fdNdEtaCorrectionProcessType[2]->FillMCParticle(vtxMC[2], eta, pt, eventTriggered, eventVertex, processType); |
593 | } | |
594 | ||
0f67a57c | 595 | if (eventTriggered) |
596 | if (eventVertex) | |
597 | fdNdEtaAnalysisMC->FillTrack(vtxMC[2], eta, pt); | |
4ef701a0 | 598 | |
0fc41645 | 599 | // TODO this value might be needed lower for the SPD study (only used for control histograms anyway) |
ea441adf | 600 | if (TMath::Abs(eta) < 1.5) // && pt > 0.2) |
4ef701a0 | 601 | nAccepted++; |
602 | } | |
603 | ||
ea441adf | 604 | if (nAccepted == 0) |
605 | fEventStats->Fill(AliPWG0Helper::GetLastProcessType(), biny); | |
606 | ||
4ef701a0 | 607 | fMultAll->Fill(nAccepted); |
608 | if (eventTriggered) { | |
609 | fMultTr->Fill(nAccepted); | |
610 | if (eventVertex) | |
611 | fMultVtx->Fill(nAccepted); | |
0f67a57c | 612 | } |
613 | ||
3d7758c1 | 614 | Int_t processed = 0; |
615 | ||
0f67a57c | 616 | for (Int_t i=0; i<inputCount; ++i) |
617 | { | |
618 | Int_t label = labelArr[i]; | |
3d7758c1 | 619 | Int_t label2 = labelArr2[i]; |
0f67a57c | 620 | |
621 | if (label < 0) | |
622 | { | |
0fc41645 | 623 | Printf("WARNING: cannot find corresponding mc particle for track(let) %d with label %d.", i, label); |
0f67a57c | 624 | continue; |
625 | } | |
626 | ||
627 | TParticle* particle = stack->Particle(label); | |
628 | if (!particle) | |
629 | { | |
630 | AliDebug(AliLog::kError, Form("UNEXPECTED: particle with label %d not found in stack (track loop).", label)); | |
631 | continue; | |
632 | } | |
1c15d51a | 633 | |
634 | fPIDTracks->Fill(particle->GetPdgCode()); | |
3d7758c1 | 635 | |
636 | // find particle that is filled in the correction map | |
637 | // this should be the particle that has been reconstructed | |
638 | // for TPC this is clear and always identified by the track's label | |
639 | // for SPD the labels might not be identical. In this case the particle closest to the beam line that is a primary is taken: | |
640 | // L1 L2 (P = primary, S = secondary) | |
641 | // P P' : --> P | |
642 | // P S : --> P | |
643 | // S P : --> P | |
644 | // S S' : --> S | |
645 | ||
646 | if (label != label2 && label2 >= 0) | |
647 | { | |
648 | if (stack->IsPhysicalPrimary(label) == kFALSE && stack->IsPhysicalPrimary(label2)) | |
649 | { | |
650 | particle = stack->Particle(label2); | |
651 | if (!particle) | |
652 | { | |
653 | AliDebug(AliLog::kError, Form("UNEXPECTED: particle with label %d not found in stack (track loop).", label2)); | |
654 | continue; | |
655 | } | |
656 | } | |
657 | } | |
0f67a57c | 658 | |
659 | if (eventTriggered && eventVertex) | |
660 | { | |
3d7758c1 | 661 | if (SignOK(particle->GetPDG()) == kFALSE) |
662 | continue; | |
663 | ||
664 | processed++; | |
665 | ||
1c15d51a | 666 | // resolutions |
667 | fEtaResolution->Fill(particle->Eta() - etaArr[i]); | |
668 | fpTResolution->Fill(particle->Pt() - ptArr[i]); | |
669 | ||
3d7758c1 | 670 | Bool_t firstIsPrim = stack->IsPhysicalPrimary(label); |
671 | // in case of primary the MC values are filled, otherwise (background) the reconstructed values | |
672 | if (label == label2 && firstIsPrim) | |
673 | { | |
674 | fdNdEtaCorrection->FillTrackedParticle(vtxMC[2], particle->Eta(), particle->Pt()); | |
675 | } | |
676 | else | |
677 | { | |
678 | fdNdEtaCorrection->FillTrackedParticle(vtxMC[2], etaArr[i], ptArr[i]); | |
3d7758c1 | 679 | } |
745d6088 | 680 | |
681 | fEtaProfile->Fill(particle->Eta(), particle->Eta() - etaArr[i]); | |
0f67a57c | 682 | fdNdEtaAnalysisESD->FillTrack(vtxMC[2], particle->Eta(), particle->Pt()); |
3d7758c1 | 683 | |
684 | fEtaCorrelation->Fill(etaArr[i], particle->Eta()); | |
685 | ||
0f67a57c | 686 | if (particle->Pt() > 0.1 && particle->Pt() < 0.2) |
687 | { | |
688 | fPIDTracks->Fill(particle->GetPdgCode()); | |
689 | } | |
96fcc8a7 | 690 | |
691 | if (fdNdEtaCorrectionProcessType[0]) | |
692 | { | |
693 | // non diffractive | |
6b62a9c7 | 694 | if (processType == AliPWG0Helper::kND) |
96fcc8a7 | 695 | fdNdEtaCorrectionProcessType[0]->FillTrackedParticle(vtxMC[2], particle->Eta(), particle->Pt()); |
696 | ||
697 | // single diffractive | |
1c15d51a | 698 | if (processType == AliPWG0Helper::kSD) |
96fcc8a7 | 699 | fdNdEtaCorrectionProcessType[1]->FillTrackedParticle(vtxMC[2], particle->Eta(), particle->Pt()); |
700 | ||
701 | // double diffractive | |
6b62a9c7 | 702 | if (processType == AliPWG0Helper::kDD) |
96fcc8a7 | 703 | fdNdEtaCorrectionProcessType[2]->FillTrackedParticle(vtxMC[2], particle->Eta(), particle->Pt()); |
704 | } | |
4ef701a0 | 705 | |
745d6088 | 706 | // control histograms |
3d7758c1 | 707 | Int_t hist = -1; |
708 | if (label == label2) | |
4ef701a0 | 709 | { |
3d7758c1 | 710 | if (firstIsPrim) |
711 | { | |
712 | hist = 0; | |
713 | } | |
714 | else | |
715 | hist = 1; | |
716 | } | |
717 | else if (label2 >= 0) | |
718 | { | |
719 | Bool_t secondIsPrim = stack->IsPhysicalPrimary(label2); | |
720 | if (firstIsPrim && secondIsPrim) | |
721 | { | |
722 | hist = 2; | |
723 | } | |
724 | else if (firstIsPrim && !secondIsPrim) | |
725 | { | |
726 | hist = 3; | |
727 | ||
728 | // check if secondary is caused by the primary or it is a fake combination | |
729 | //Printf("PS case --> Label 1 is %d, label 2 is %d", label, label2); | |
730 | Int_t mother = label2; | |
731 | while (stack->Particle(mother) && stack->Particle(mother)->GetMother(0) >= 0) | |
732 | { | |
733 | //Printf(" %d created by %d, %d", mother, stack->Particle(mother)->GetMother(0), stack->Particle(mother)->GetMother(1)); | |
734 | if (stack->Particle(mother)->GetMother(0) == label) | |
735 | { | |
736 | /*Printf("The untraceable decay was:"); | |
737 | Printf(" from:"); | |
738 | particle->Print(); | |
739 | Printf(" to (status code %d):", stack->Particle(mother)->GetStatusCode()); | |
740 | stack->Particle(mother)->Print();*/ | |
741 | hist = 4; | |
742 | } | |
743 | mother = stack->Particle(mother)->GetMother(0); | |
744 | } | |
745 | } | |
746 | else if (!firstIsPrim && secondIsPrim) | |
747 | { | |
748 | hist = 5; | |
749 | } | |
750 | else if (!firstIsPrim && !secondIsPrim) | |
751 | { | |
752 | hist = 6; | |
753 | } | |
754 | ||
4ef701a0 | 755 | } |
756 | else | |
3d7758c1 | 757 | hist = 7; |
758 | ||
759 | fDeltaPhi[hist]->Fill(deltaPhiArr[i]); | |
0f67a57c | 760 | } |
761 | } | |
762 | ||
3d7758c1 | 763 | if (processed < inputCount) |
764 | Printf("Only %d out of %d track(let)s were used", processed, inputCount); | |
765 | ||
0f67a57c | 766 | if (eventTriggered && eventVertex) |
96fcc8a7 | 767 | { |
0f67a57c | 768 | fdNdEtaAnalysisMC->FillEvent(vtxMC[2], inputCount); |
96fcc8a7 | 769 | fdNdEtaAnalysisESD->FillEvent(vtxMC[2], inputCount); |
770 | } | |
0f67a57c | 771 | |
96fcc8a7 | 772 | // stuff regarding the vertex reco correction and trigger bias correction |
0f67a57c | 773 | fdNdEtaCorrection->FillEvent(vtxMC[2], inputCount, eventTriggered, eventVertex, processType); |
96fcc8a7 | 774 | |
775 | if (fdNdEtaCorrectionProcessType[0]) | |
776 | { | |
777 | // non diffractive | |
6b62a9c7 | 778 | if (processType == AliPWG0Helper::kND ) |
96fcc8a7 | 779 | fdNdEtaCorrectionProcessType[0]->FillEvent(vtxMC[2], inputCount, eventTriggered, eventVertex, processType); |
780 | ||
781 | // single diffractive | |
6b62a9c7 | 782 | if (processType == AliPWG0Helper::kSD) |
96fcc8a7 | 783 | fdNdEtaCorrectionProcessType[1]->FillEvent(vtxMC[2], inputCount, eventTriggered, eventVertex, processType); |
784 | ||
785 | // double diffractive | |
6b62a9c7 | 786 | if (processType == AliPWG0Helper::kDD) |
96fcc8a7 | 787 | fdNdEtaCorrectionProcessType[2]->FillEvent(vtxMC[2], inputCount, eventTriggered, eventVertex, processType); |
0f67a57c | 788 | } |
789 | ||
790 | delete[] etaArr; | |
791 | delete[] labelArr; | |
3d7758c1 | 792 | delete[] labelArr2; |
0f67a57c | 793 | delete[] ptArr; |
4ef701a0 | 794 | delete[] deltaPhiArr; |
0f67a57c | 795 | } |
796 | ||
797 | void AlidNdEtaCorrectionTask::Terminate(Option_t *) | |
798 | { | |
799 | // The Terminate() function is the last function to be called during | |
800 | // a query. It always runs on the client, it can be used to present | |
801 | // the results graphically or save the results to file. | |
802 | ||
803 | fOutput = dynamic_cast<TList*> (GetOutputData(0)); | |
804 | if (!fOutput) { | |
805 | Printf("ERROR: fOutput not available"); | |
806 | return; | |
807 | } | |
808 | ||
809 | fdNdEtaCorrection = dynamic_cast<AlidNdEtaCorrection*> (fOutput->FindObject("dndeta_correction")); | |
810 | fdNdEtaAnalysisMC = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndetaMC")); | |
811 | fdNdEtaAnalysisESD = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndetaESD")); | |
812 | if (!fdNdEtaCorrection || !fdNdEtaAnalysisMC || !fdNdEtaAnalysisESD) | |
813 | { | |
814 | AliDebug(AliLog::kError, "Could not read object from output list"); | |
815 | return; | |
816 | } | |
817 | ||
818 | fdNdEtaCorrection->Finish(); | |
819 | ||
820 | TString fileName; | |
821 | fileName.Form("correction_map%s.root", fOption.Data()); | |
822 | TFile* fout = new TFile(fileName, "RECREATE"); | |
823 | ||
1c15d51a | 824 | fEsdTrackCuts = dynamic_cast<AliESDtrackCuts*> (fOutput->FindObject("fEsdTrackCuts")); |
0f67a57c | 825 | if (fEsdTrackCuts) |
826 | fEsdTrackCuts->SaveHistograms("esd_track_cuts"); | |
1c15d51a | 827 | |
828 | fEsdTrackCutsPrim = dynamic_cast<AliESDtrackCuts*> (fOutput->FindObject("fEsdTrackCutsPrim")); | |
829 | if (fEsdTrackCutsPrim) | |
830 | fEsdTrackCutsPrim->SaveHistograms("esd_track_cuts_primaries"); | |
831 | ||
832 | fEsdTrackCutsSec = dynamic_cast<AliESDtrackCuts*> (fOutput->FindObject("fEsdTrackCutsSec")); | |
833 | if (fEsdTrackCutsSec) | |
834 | fEsdTrackCutsSec->SaveHistograms("esd_track_cuts_secondaries"); | |
835 | ||
0f67a57c | 836 | fdNdEtaCorrection->SaveHistograms(); |
837 | fdNdEtaAnalysisMC->SaveHistograms(); | |
838 | fdNdEtaAnalysisESD->SaveHistograms(); | |
839 | ||
96fcc8a7 | 840 | fdNdEtaCorrectionProcessType[0] = dynamic_cast<AlidNdEtaCorrection*> (fOutput->FindObject("dndeta_correction_ND")); |
841 | fdNdEtaCorrectionProcessType[1] = dynamic_cast<AlidNdEtaCorrection*> (fOutput->FindObject("dndeta_correction_SD")); | |
842 | fdNdEtaCorrectionProcessType[2] = dynamic_cast<AlidNdEtaCorrection*> (fOutput->FindObject("dndeta_correction_DD")); | |
843 | for (Int_t i=0; i<3; ++i) | |
844 | if (fdNdEtaCorrectionProcessType[i]) | |
845 | fdNdEtaCorrectionProcessType[i]->SaveHistograms(); | |
846 | ||
1c15d51a | 847 | fTemp1 = dynamic_cast<TH1*> (fOutput->FindObject("fTemp1")); |
848 | if (fTemp1) | |
849 | fTemp1->Write(); | |
850 | fTemp2 = dynamic_cast<TH1*> (fOutput->FindObject("fTemp2")); | |
851 | if (fTemp2) | |
852 | fTemp2->Write(); | |
96fcc8a7 | 853 | |
0fc41645 | 854 | fVertexCorrelation = dynamic_cast<TH2F*> (fOutput->FindObject("fVertexCorrelation")); |
54b096ef | 855 | if (fVertexCorrelation) |
856 | fVertexCorrelation->Write(); | |
0fc41645 | 857 | fVertexProfile = dynamic_cast<TProfile*> (fOutput->FindObject("fVertexProfile")); |
54b096ef | 858 | if (fVertexProfile) |
859 | fVertexProfile->Write(); | |
1c15d51a | 860 | fVertexShift = dynamic_cast<TH1F*> (fOutput->FindObject("fVertexShift")); |
861 | if (fVertexShift) | |
862 | fVertexShift->Write(); | |
0fc41645 | 863 | fVertexShiftNorm = dynamic_cast<TH1F*> (fOutput->FindObject("fVertexShiftNorm")); |
54b096ef | 864 | if (fVertexShiftNorm) |
865 | fVertexShiftNorm->Write(); | |
3d7758c1 | 866 | |
0fc41645 | 867 | fEtaCorrelation = dynamic_cast<TH2F*> (fOutput->FindObject("fEtaCorrelation")); |
3d7758c1 | 868 | if (fEtaCorrelation) |
869 | fEtaCorrelation->Write(); | |
0fc41645 | 870 | fEtaProfile = dynamic_cast<TProfile*> (fOutput->FindObject("fEtaProfile")); |
3d7758c1 | 871 | if (fEtaProfile) |
872 | fEtaProfile->Write(); | |
1c15d51a | 873 | fEtaResolution = dynamic_cast<TH1F*> (fOutput->FindObject("fEtaResolution")); |
874 | if (fEtaResolution) | |
875 | fEtaResolution->Write(); | |
876 | fpTResolution = dynamic_cast<TH1F*> (fOutput->FindObject("fpTResolution")); | |
877 | if (fpTResolution) | |
878 | fpTResolution->Write(); | |
3d7758c1 | 879 | |
0fc41645 | 880 | fMultAll = dynamic_cast<TH1F*> (fOutput->FindObject("fMultAll")); |
4ef701a0 | 881 | if (fMultAll) |
882 | fMultAll->Write(); | |
883 | ||
0fc41645 | 884 | fMultTr = dynamic_cast<TH1F*> (fOutput->FindObject("fMultTr")); |
4ef701a0 | 885 | if (fMultTr) |
886 | fMultTr->Write(); | |
0fc41645 | 887 | |
888 | fMultVtx = dynamic_cast<TH1F*> (fOutput->FindObject("fMultVtx")); | |
4ef701a0 | 889 | if (fMultVtx) |
890 | fMultVtx->Write(); | |
891 | ||
3d7758c1 | 892 | for (Int_t i=0; i<8; ++i) |
0fc41645 | 893 | { |
894 | fDeltaPhi[i] = dynamic_cast<TH1F*> (fOutput->FindObject(Form("fDeltaPhi_%d", i))); | |
4ef701a0 | 895 | if (fDeltaPhi[i]) |
896 | fDeltaPhi[i]->Write(); | |
0fc41645 | 897 | } |
54b096ef | 898 | |
0fc41645 | 899 | fEventStats = dynamic_cast<TH2F*> (fOutput->FindObject("fEventStats")); |
50ec344d | 900 | if (fEventStats) |
901 | fEventStats->Write(); | |
902 | ||
1c15d51a | 903 | fPIDParticles = dynamic_cast<TH1F*> (fOutput->FindObject("fPIDParticles")); |
904 | if (fPIDParticles) | |
905 | fPIDParticles->Write(); | |
906 | ||
907 | fPIDTracks = dynamic_cast<TH1F*> (fOutput->FindObject("fPIDTracks")); | |
908 | if (fPIDTracks) | |
909 | fPIDTracks->Write(); | |
0f67a57c | 910 | |
1c15d51a | 911 | //fdNdEtaCorrection->DrawHistograms(); |
0f67a57c | 912 | |
96fcc8a7 | 913 | Printf("Writting result to %s", fileName.Data()); |
914 | ||
0f67a57c | 915 | if (fPIDParticles && fPIDTracks) |
916 | { | |
0f67a57c | 917 | TDatabasePDG* pdgDB = new TDatabasePDG; |
918 | ||
919 | for (Int_t i=0; i <= fPIDParticles->GetNbinsX()+1; ++i) | |
920 | if (fPIDParticles->GetBinContent(i) > 0) | |
1c15d51a | 921 | { |
922 | TObject* pdgParticle = pdgDB->GetParticle((Int_t) fPIDParticles->GetBinCenter(i)); | |
923 | printf("PDG = %d (%s): generated: %d, reconstructed: %d, ratio: %f\n", (Int_t) fPIDParticles->GetBinCenter(i), (pdgParticle) ? pdgParticle->GetName() : "not found", (Int_t) fPIDParticles->GetBinContent(i), (Int_t) fPIDTracks->GetBinContent(i), ((fPIDTracks->GetBinContent(i) > 0) ? fPIDParticles->GetBinContent(i) / fPIDTracks->GetBinContent(i) : -1)); | |
924 | } | |
0f67a57c | 925 | |
926 | delete pdgDB; | |
927 | pdgDB = 0; | |
928 | } | |
1c15d51a | 929 | |
930 | fout->Write(); | |
931 | fout->Close(); | |
0f67a57c | 932 | } |
933 | ||
934 | Bool_t AlidNdEtaCorrectionTask::SignOK(TParticlePDG* particle) | |
935 | { | |
936 | // returns if a particle with this sign should be counted | |
937 | // this is determined by the value of fSignMode, which should have the same sign | |
938 | // as the charge | |
939 | // if fSignMode is 0 all particles are counted | |
940 | ||
941 | if (fSignMode == 0) | |
942 | return kTRUE; | |
943 | ||
944 | if (!particle) | |
945 | { | |
946 | printf("WARNING: not counting a particle that does not have a pdg particle\n"); | |
947 | return kFALSE; | |
948 | } | |
949 | ||
950 | Double_t charge = particle->Charge(); | |
951 | ||
952 | if (fSignMode > 0) | |
953 | if (charge < 0) | |
954 | return kFALSE; | |
955 | ||
956 | if (fSignMode < 0) | |
957 | if (charge > 0) | |
958 | return kFALSE; | |
959 | ||
960 | return kTRUE; | |
961 | } | |
96fcc8a7 | 962 |