]>
Commit | Line | Data |
---|---|---|
0f67a57c | 1 | /* $Id$ */ |
2 | ||
3 | #include "AlidNdEtaTask.h" | |
4 | ||
5 | #include <TStyle.h> | |
6 | #include <TSystem.h> | |
7 | #include <TCanvas.h> | |
8 | #include <TVector3.h> | |
9 | #include <TChain.h> | |
10 | #include <TFile.h> | |
11 | #include <TH1F.h> | |
12 | #include <TH2F.h> | |
13 | #include <TH3F.h> | |
14 | #include <TParticle.h> | |
15 | #include <TRandom.h> | |
16 | #include <TNtuple.h> | |
17 | #include <TObjString.h> | |
18 | #include <TF1.h> | |
19 | ||
20 | #include <AliLog.h> | |
21 | #include <AliESDVertex.h> | |
22 | #include <AliESDEvent.h> | |
23 | #include <AliStack.h> | |
24 | #include <AliHeader.h> | |
25 | #include <AliGenEventHeader.h> | |
26 | #include <AliMultiplicity.h> | |
27 | #include <AliAnalysisManager.h> | |
28 | #include <AliMCEventHandler.h> | |
29 | #include <AliMCEvent.h> | |
30 | #include <AliESDInputHandler.h> | |
31 | ||
32 | #include "esdTrackCuts/AliESDtrackCuts.h" | |
33 | #include "AliPWG0Helper.h" | |
34 | #include "AliCorrection.h" | |
35 | #include "AliCorrectionMatrix3D.h" | |
36 | #include "dNdEta/dNdEtaAnalysis.h" | |
37 | ||
38 | ClassImp(AlidNdEtaTask) | |
39 | ||
40 | AlidNdEtaTask::AlidNdEtaTask(const char* opt) : | |
41 | AliAnalysisTask("AlidNdEtaTask", ""), | |
42 | fESD(0), | |
43 | fOutput(0), | |
44 | fOption(opt), | |
770a1f1d | 45 | fAnalysisMode(AliPWG0Helper::kTPC), |
0f67a57c | 46 | fReadMC(kFALSE), |
47 | fEsdTrackCuts(0), | |
48 | fdNdEtaAnalysisESD(0), | |
49 | fMult(0), | |
50 | fMultVtx(0), | |
51 | fEvents(0), | |
52 | fdNdEtaAnalysis(0), | |
53 | fdNdEtaAnalysisTr(0), | |
54 | fdNdEtaAnalysisTrVtx(0), | |
55 | fdNdEtaAnalysisTracks(0), | |
56 | fVertex(0), | |
57 | fPartPt(0) | |
58 | { | |
59 | // | |
60 | // Constructor. Initialization of pointers | |
61 | // | |
62 | ||
63 | // Define input and output slots here | |
64 | DefineInput(0, TChain::Class()); | |
65 | DefineOutput(0, TList::Class()); | |
66 | } | |
67 | ||
68 | AlidNdEtaTask::~AlidNdEtaTask() | |
69 | { | |
70 | // | |
71 | // Destructor | |
72 | // | |
73 | ||
74 | // histograms are in the output list and deleted when the output | |
75 | // list is deleted by the TSelector dtor | |
76 | ||
77 | if (fOutput) { | |
78 | delete fOutput; | |
79 | fOutput = 0; | |
80 | } | |
81 | } | |
82 | ||
83 | //________________________________________________________________________ | |
84 | void AlidNdEtaTask::ConnectInputData(Option_t *) | |
85 | { | |
86 | // Connect ESD | |
87 | // Called once | |
88 | ||
89 | Printf("AlidNdEtaTask::ConnectInputData called"); | |
90 | ||
91 | TTree* tree = dynamic_cast<TTree*> (GetInputData(0)); | |
92 | if (!tree) { | |
93 | Printf("ERROR: Could not read tree from input slot 0"); | |
94 | } else { | |
95 | // Disable all branches and enable only the needed ones | |
96 | //tree->SetBranchStatus("*", 0); | |
97 | ||
98 | tree->SetBranchStatus("fTriggerMask", 1); | |
99 | tree->SetBranchStatus("fSPDVertex*", 1); | |
770a1f1d | 100 | // PrimaryVertex also needed |
0f67a57c | 101 | |
770a1f1d | 102 | if (fAnalysisMode == AliPWG0Helper::kSPD) |
0f67a57c | 103 | tree->SetBranchStatus("fSPDMult*", 1); |
104 | ||
770a1f1d | 105 | if (fAnalysisMode == AliPWG0Helper::kTPC || fAnalysisMode == AliPWG0Helper::kTPCITS) { |
0f67a57c | 106 | AliESDtrackCuts::EnableNeededBranches(tree); |
107 | tree->SetBranchStatus("fTracks.fLabel", 1); | |
108 | } | |
109 | ||
0f67a57c | 110 | AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()); |
111 | ||
112 | if (!esdH) { | |
113 | Printf("ERROR: Could not get ESDInputHandler"); | |
114 | } else | |
115 | fESD = esdH->GetEvent(); | |
116 | } | |
117 | } | |
118 | ||
119 | void AlidNdEtaTask::CreateOutputObjects() | |
120 | { | |
121 | // create result objects and add to output list | |
122 | ||
123 | fOutput = new TList; | |
124 | fOutput->SetOwner(); | |
125 | ||
770a1f1d | 126 | fdNdEtaAnalysisESD = new dNdEtaAnalysis("fdNdEtaAnalysisESD", "fdNdEtaAnalysisESD", fAnalysisMode); |
0f67a57c | 127 | fOutput->Add(fdNdEtaAnalysisESD); |
128 | ||
129 | fMult = new TH1F("fMult", "fMult;Ntracks;Count", 201, -0.5, 200.5); | |
130 | fOutput->Add(fMult); | |
131 | ||
132 | fMultVtx = new TH1F("fMultVtx", "fMultVtx;Ntracks;Count", 201, -0.5, 200.5); | |
133 | fOutput->Add(fMultVtx); | |
134 | ||
135 | for (Int_t i=0; i<3; ++i) | |
136 | { | |
137 | fPartEta[i] = new TH1F(Form("dndeta_check_%d", i), Form("dndeta_check_%d", i), 60, -6, 6); | |
138 | fPartEta[i]->Sumw2(); | |
139 | fOutput->Add(fPartEta[i]); | |
140 | } | |
141 | ||
142 | fEvents = new TH1F("dndeta_check_vertex", "dndeta_check_vertex", 40, -20, 20); | |
143 | fOutput->Add(fEvents); | |
144 | ||
145 | if (fReadMC) | |
146 | { | |
770a1f1d | 147 | fdNdEtaAnalysis = new dNdEtaAnalysis("dndeta", "dndeta", fAnalysisMode); |
0f67a57c | 148 | fOutput->Add(fdNdEtaAnalysis); |
149 | ||
770a1f1d | 150 | fdNdEtaAnalysisTr = new dNdEtaAnalysis("dndetaTr", "dndetaTr", fAnalysisMode); |
0f67a57c | 151 | fOutput->Add(fdNdEtaAnalysisTr); |
152 | ||
770a1f1d | 153 | fdNdEtaAnalysisTrVtx = new dNdEtaAnalysis("dndetaTrVtx", "dndetaTrVtx", fAnalysisMode); |
0f67a57c | 154 | fOutput->Add(fdNdEtaAnalysisTrVtx); |
155 | ||
770a1f1d | 156 | fdNdEtaAnalysisTracks = new dNdEtaAnalysis("dndetaTracks", "dndetaTracks", fAnalysisMode); |
0f67a57c | 157 | fOutput->Add(fdNdEtaAnalysisTracks); |
158 | ||
159 | fVertex = new TH3F("vertex_check", "vertex_check", 50, -50, 50, 50, -50, 50, 50, -50, 50); | |
160 | fOutput->Add(fVertex); | |
161 | ||
162 | fPartPt = new TH1F("dndeta_check_pt", "dndeta_check_pt", 1000, 0, 10); | |
163 | fPartPt->Sumw2(); | |
164 | fOutput->Add(fPartPt); | |
165 | } | |
166 | } | |
167 | ||
168 | void AlidNdEtaTask::Exec(Option_t*) | |
169 | { | |
170 | // process the event | |
171 | ||
172 | // Check prerequisites | |
173 | if (!fESD) | |
174 | { | |
175 | AliDebug(AliLog::kError, "ESD branch not available"); | |
176 | return; | |
177 | } | |
178 | ||
179 | // trigger definition | |
180 | Bool_t eventTriggered = AliPWG0Helper::IsEventTriggered(fESD->GetTriggerMask(), AliPWG0Helper::kMB2); | |
181 | ||
0f67a57c | 182 | // get the ESD vertex |
770a1f1d | 183 | const AliESDVertex* vtxESD = AliPWG0Helper::GetVertex(fESD, fAnalysisMode); |
184 | ||
0f67a57c | 185 | Double_t vtx[3]; |
770a1f1d | 186 | if (vtxESD) |
187 | vtxESD->GetXYZ(vtx); | |
0f67a57c | 188 | |
189 | // post the data already here | |
190 | PostData(0, fOutput); | |
191 | ||
192 | // create list of (label, eta, pt) tuples | |
193 | Int_t inputCount = 0; | |
194 | Int_t* labelArr = 0; | |
195 | Float_t* etaArr = 0; | |
196 | Float_t* ptArr = 0; | |
770a1f1d | 197 | if (fAnalysisMode == AliPWG0Helper::kSPD) |
0f67a57c | 198 | { |
199 | // get tracklets | |
200 | const AliMultiplicity* mult = fESD->GetMultiplicity(); | |
201 | if (!mult) | |
202 | { | |
203 | AliDebug(AliLog::kError, "AliMultiplicity not available"); | |
204 | return; | |
205 | } | |
206 | ||
207 | labelArr = new Int_t[mult->GetNumberOfTracklets()]; | |
208 | etaArr = new Float_t[mult->GetNumberOfTracklets()]; | |
209 | ptArr = new Float_t[mult->GetNumberOfTracklets()]; | |
210 | ||
211 | // get multiplicity from ITS tracklets | |
212 | for (Int_t i=0; i<mult->GetNumberOfTracklets(); ++i) | |
213 | { | |
214 | //printf("%d %f %f %f\n", i, mult->GetTheta(i), mult->GetPhi(i), mult->GetDeltaPhi(i)); | |
215 | ||
216 | // This removes non-tracklets in PDC06 data. Very bad solution. New solution is implemented for newer data. Keeping this for compatibility. | |
217 | if (mult->GetDeltaPhi(i) < -1000) | |
218 | continue; | |
219 | ||
220 | etaArr[inputCount] = mult->GetEta(i); | |
221 | labelArr[inputCount] = mult->GetLabel(i); | |
222 | ptArr[inputCount] = 0; // no pt for tracklets | |
223 | ++inputCount; | |
224 | } | |
225 | } | |
770a1f1d | 226 | else if (fAnalysisMode == AliPWG0Helper::kTPC || fAnalysisMode == AliPWG0Helper::kTPCITS) |
0f67a57c | 227 | { |
228 | if (!fEsdTrackCuts) | |
229 | { | |
230 | AliDebug(AliLog::kError, "fESDTrackCuts not available"); | |
231 | return; | |
232 | } | |
233 | ||
234 | // get multiplicity from ESD tracks | |
235 | TObjArray* list = fEsdTrackCuts->GetAcceptedTracks(fESD); | |
236 | Int_t nGoodTracks = list->GetEntries(); | |
237 | ||
238 | labelArr = new Int_t[nGoodTracks]; | |
239 | etaArr = new Float_t[nGoodTracks]; | |
240 | ptArr = new Float_t[nGoodTracks]; | |
241 | ||
242 | // loop over esd tracks | |
243 | for (Int_t i=0; i<nGoodTracks; i++) | |
244 | { | |
245 | AliESDtrack* esdTrack = dynamic_cast<AliESDtrack*> (list->At(i)); | |
246 | if (!esdTrack) | |
247 | { | |
248 | AliDebug(AliLog::kError, Form("ERROR: Could not retrieve track %d.", i)); | |
249 | continue; | |
250 | } | |
251 | ||
252 | etaArr[inputCount] = esdTrack->Eta(); | |
253 | labelArr[inputCount] = TMath::Abs(esdTrack->GetLabel()); | |
254 | ptArr[inputCount] = esdTrack->Pt(); | |
255 | ++inputCount; | |
256 | } | |
2b26296f | 257 | |
258 | delete list; | |
0f67a57c | 259 | } |
260 | else | |
261 | return; | |
262 | ||
263 | // Processing of ESD information (always) | |
264 | if (eventTriggered) | |
265 | { | |
dd367a14 | 266 | // control hist |
0f67a57c | 267 | fMult->Fill(inputCount); |
770a1f1d | 268 | if (vtxESD) |
dd367a14 | 269 | { |
270 | // control hist | |
0f67a57c | 271 | fMultVtx->Fill(inputCount); |
272 | ||
dd367a14 | 273 | for (Int_t i=0; i<inputCount; ++i) |
274 | { | |
275 | Float_t eta = etaArr[i]; | |
276 | Float_t pt = ptArr[i]; | |
0f67a57c | 277 | |
dd367a14 | 278 | fdNdEtaAnalysisESD->FillTrack(vtx[2], eta, pt); |
0f67a57c | 279 | |
dd367a14 | 280 | if (TMath::Abs(vtx[2]) < 20) |
281 | { | |
282 | fPartEta[0]->Fill(eta); | |
0f67a57c | 283 | |
dd367a14 | 284 | if (vtx[2] < 0) |
285 | fPartEta[1]->Fill(eta); | |
286 | else | |
287 | fPartEta[2]->Fill(eta); | |
288 | } | |
0f67a57c | 289 | } |
0f67a57c | 290 | |
dd367a14 | 291 | // for event count per vertex |
292 | fdNdEtaAnalysisESD->FillEvent(vtx[2], inputCount); | |
0f67a57c | 293 | |
dd367a14 | 294 | // control hists |
295 | fEvents->Fill(vtx[2]); | |
296 | } | |
0f67a57c | 297 | } |
298 | ||
299 | if (fReadMC) // Processing of MC information (optional) | |
300 | { | |
301 | AliMCEventHandler* eventHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler()); | |
302 | if (!eventHandler) { | |
303 | Printf("ERROR: Could not retrieve MC event handler"); | |
304 | return; | |
305 | } | |
306 | ||
307 | AliMCEvent* mcEvent = eventHandler->MCEvent(); | |
308 | if (!mcEvent) { | |
309 | Printf("ERROR: Could not retrieve MC event"); | |
310 | return; | |
311 | } | |
312 | ||
313 | AliStack* stack = mcEvent->Stack(); | |
314 | if (!stack) | |
315 | { | |
316 | AliDebug(AliLog::kError, "Stack not available"); | |
317 | return; | |
318 | } | |
319 | ||
320 | AliHeader* header = mcEvent->Header(); | |
321 | if (!header) | |
322 | { | |
323 | AliDebug(AliLog::kError, "Header not available"); | |
324 | return; | |
325 | } | |
326 | ||
327 | // get the MC vertex | |
328 | AliGenEventHeader* genHeader = header->GenEventHeader(); | |
329 | TArrayF vtxMC(3); | |
330 | genHeader->PrimaryVertex(vtxMC); | |
331 | ||
332 | // loop over mc particles | |
333 | Int_t nPrim = stack->GetNprimary(); | |
334 | ||
335 | Int_t nAcceptedParticles = 0; | |
336 | ||
337 | for (Int_t iMc = 0; iMc < nPrim; ++iMc) | |
338 | { | |
339 | //Printf("Getting particle %d", iMc); | |
340 | TParticle* particle = stack->Particle(iMc); | |
341 | ||
342 | if (!particle) | |
343 | continue; | |
344 | ||
345 | if (AliPWG0Helper::IsPrimaryCharged(particle, nPrim) == kFALSE) | |
346 | continue; | |
347 | ||
348 | AliDebug(AliLog::kDebug+1, Form("Accepted primary %d, unique ID: %d", iMc, particle->GetUniqueID())); | |
0f67a57c | 349 | Float_t eta = particle->Eta(); |
350 | Float_t pt = particle->Pt(); | |
351 | ||
2b26296f | 352 | // make a rough eta cut (so that nAcceptedParticles is not too far off) |
353 | if (TMath::Abs(eta) < 1.5) | |
354 | nAcceptedParticles++; | |
355 | ||
0f67a57c | 356 | fdNdEtaAnalysis->FillTrack(vtxMC[2], eta, pt); |
357 | fVertex->Fill(particle->Vx(), particle->Vy(), particle->Vz()); | |
358 | ||
359 | if (eventTriggered) | |
360 | { | |
361 | fdNdEtaAnalysisTr->FillTrack(vtxMC[2], eta, pt); | |
770a1f1d | 362 | if (vtxESD) |
0f67a57c | 363 | fdNdEtaAnalysisTrVtx->FillTrack(vtxMC[2], eta, pt); |
364 | } | |
365 | ||
366 | if (TMath::Abs(eta) < 0.8) | |
367 | fPartPt->Fill(pt); | |
368 | } | |
369 | ||
370 | fdNdEtaAnalysis->FillEvent(vtxMC[2], nAcceptedParticles); | |
371 | if (eventTriggered) | |
372 | { | |
373 | fdNdEtaAnalysisTr->FillEvent(vtxMC[2], nAcceptedParticles); | |
770a1f1d | 374 | if (vtxESD) |
0f67a57c | 375 | fdNdEtaAnalysisTrVtx->FillEvent(vtxMC[2], nAcceptedParticles); |
376 | } | |
377 | ||
770a1f1d | 378 | if (eventTriggered && vtxESD) |
0f67a57c | 379 | { |
380 | // from tracks is only done for triggered and vertex reconstructed events | |
381 | ||
382 | for (Int_t i=0; i<inputCount; ++i) | |
383 | { | |
384 | Int_t label = labelArr[i]; | |
385 | ||
386 | if (label < 0) | |
387 | continue; | |
388 | ||
389 | //Printf("Getting particle of track %d", label); | |
390 | TParticle* particle = stack->Particle(label); | |
391 | ||
392 | if (!particle) | |
393 | { | |
394 | AliDebug(AliLog::kError, Form("ERROR: Could not retrieve particle %d.", label)); | |
395 | continue; | |
396 | } | |
397 | ||
398 | fdNdEtaAnalysisTracks->FillTrack(vtxMC[2], particle->Eta(), particle->Pt()); | |
399 | } // end of track loop | |
400 | ||
401 | // for event count per vertex | |
402 | fdNdEtaAnalysisTracks->FillEvent(vtxMC[2], inputCount); | |
403 | } | |
404 | } | |
405 | ||
406 | delete[] etaArr; | |
407 | delete[] labelArr; | |
408 | delete[] ptArr; | |
409 | } | |
410 | ||
411 | void AlidNdEtaTask::Terminate(Option_t *) | |
412 | { | |
413 | // The Terminate() function is the last function to be called during | |
414 | // a query. It always runs on the client, it can be used to present | |
415 | // the results graphically or save the results to file. | |
416 | ||
417 | fOutput = dynamic_cast<TList*> (GetOutputData(0)); | |
418 | if (!fOutput) { | |
419 | Printf("ERROR: fOutput not available"); | |
420 | return; | |
421 | } | |
422 | ||
423 | fdNdEtaAnalysisESD = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("fdNdEtaAnalysisESD")); | |
424 | fMult = dynamic_cast<TH1F*> (fOutput->FindObject("fMult")); | |
425 | fMultVtx = dynamic_cast<TH1F*> (fOutput->FindObject("fMultVtx")); | |
426 | fEvents = dynamic_cast<TH1F*> (fOutput->FindObject("dndeta_check_vertex")); | |
427 | for (Int_t i=0; i<3; ++i) | |
428 | fPartEta[i] = dynamic_cast<TH1F*> (fOutput->FindObject(Form("dndeta_check_%d", i))); | |
429 | ||
430 | if (!fdNdEtaAnalysisESD) | |
431 | { | |
432 | AliDebug(AliLog::kError, "ERROR: fdNdEtaAnalysisESD not available"); | |
433 | return; | |
434 | } | |
435 | ||
436 | if (fMult && fMultVtx) | |
437 | { | |
438 | new TCanvas; | |
439 | fMult->Draw(); | |
440 | fMultVtx->SetLineColor(2); | |
770a1f1d | 441 | fMultVtx->Draw("SAME"); |
0f67a57c | 442 | } |
443 | ||
444 | if (fPartEta[0]) | |
445 | { | |
446 | Int_t events1 = (Int_t) fEvents->Integral(fEvents->GetXaxis()->FindBin(-19.9), fEvents->GetXaxis()->FindBin(-0.001)); | |
447 | Int_t events2 = (Int_t) fEvents->Integral(fEvents->GetXaxis()->FindBin(0.001), fEvents->GetXaxis()->FindBin(19.9)); | |
448 | ||
449 | fPartEta[0]->Scale(1.0 / (events1 + events2)); | |
450 | fPartEta[1]->Scale(1.0 / events1); | |
451 | fPartEta[2]->Scale(1.0 / events2); | |
452 | ||
453 | for (Int_t i=0; i<3; ++i) | |
454 | fPartEta[i]->Scale(1.0 / fPartEta[i]->GetBinWidth(1)); | |
455 | ||
456 | new TCanvas("control", "control", 500, 500); | |
457 | for (Int_t i=0; i<3; ++i) | |
458 | { | |
459 | fPartEta[i]->SetLineColor(i+1); | |
460 | fPartEta[i]->Draw((i==0) ? "" : "SAME"); | |
461 | } | |
462 | } | |
463 | ||
464 | if (fEvents) | |
465 | { | |
466 | new TCanvas("control3", "control3", 500, 500); | |
770a1f1d | 467 | fEvents->Draw(); |
0f67a57c | 468 | } |
469 | ||
470 | TFile* fout = new TFile("analysis_esd_raw.root", "RECREATE"); | |
471 | ||
472 | if (fdNdEtaAnalysisESD) | |
473 | fdNdEtaAnalysisESD->SaveHistograms(); | |
474 | ||
475 | if (fEsdTrackCuts) | |
476 | fEsdTrackCuts->SaveHistograms("esd_tracks_cuts"); | |
477 | ||
478 | if (fMult) | |
479 | fMult->Write(); | |
480 | ||
481 | if (fMultVtx) | |
482 | fMultVtx->Write(); | |
483 | ||
484 | for (Int_t i=0; i<3; ++i) | |
485 | if (fPartEta[i]) | |
486 | fPartEta[i]->Write(); | |
487 | ||
488 | if (fEvents) | |
489 | fEvents->Write(); | |
490 | ||
491 | fout->Write(); | |
492 | fout->Close(); | |
493 | ||
494 | Printf("Writting result to analysis_esd_raw.root"); | |
495 | ||
496 | if (fReadMC) | |
497 | { | |
498 | fdNdEtaAnalysis = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndeta")); | |
499 | fdNdEtaAnalysisTr = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndetaTr")); | |
500 | fdNdEtaAnalysisTrVtx = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndetaTrVtx")); | |
501 | fdNdEtaAnalysisTracks = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndetaTracks")); | |
502 | fPartPt = dynamic_cast<TH1F*> (fOutput->FindObject("dndeta_check_pt")); | |
503 | ||
504 | if (!fdNdEtaAnalysis || !fdNdEtaAnalysisTr || !fdNdEtaAnalysisTrVtx || !fPartPt) | |
505 | { | |
506 | AliDebug(AliLog::kError, Form("ERROR: Histograms not available %p %p", (void*) fdNdEtaAnalysis, (void*) fPartPt)); | |
507 | return; | |
508 | } | |
509 | ||
510 | fdNdEtaAnalysis->Finish(0, -1, AlidNdEtaCorrection::kNone); | |
511 | fdNdEtaAnalysisTr->Finish(0, -1, AlidNdEtaCorrection::kNone); | |
512 | fdNdEtaAnalysisTrVtx->Finish(0, -1, AlidNdEtaCorrection::kNone); | |
513 | fdNdEtaAnalysisTracks->Finish(0, -1, AlidNdEtaCorrection::kNone); | |
514 | ||
515 | Int_t events = (Int_t) fdNdEtaAnalysis->GetData()->GetEventCorrection()->GetMeasuredHistogram()->Integral(); | |
516 | fPartPt->Scale(1.0/events); | |
517 | fPartPt->Scale(1.0/fPartPt->GetBinWidth(1)); | |
518 | ||
519 | TFile* fout = new TFile("analysis_mc.root","RECREATE"); | |
520 | ||
521 | fdNdEtaAnalysis->SaveHistograms(); | |
522 | fdNdEtaAnalysisTr->SaveHistograms(); | |
523 | fdNdEtaAnalysisTrVtx->SaveHistograms(); | |
524 | fdNdEtaAnalysisTracks->SaveHistograms(); | |
525 | fPartPt->Write(); | |
526 | ||
527 | fout->Write(); | |
528 | fout->Close(); | |
529 | ||
530 | Printf("Writting result to analysis_mc.root"); | |
531 | ||
532 | if (fPartPt) | |
533 | { | |
534 | new TCanvas("control2", "control2", 500, 500); | |
770a1f1d | 535 | fPartPt->Draw(); |
0f67a57c | 536 | } |
537 | } | |
538 | } |