fixing warnings
[u/mrichter/AliRoot.git] / PWG0 / dNdEta / AlidNdEtaTask.cxx
CommitLineData
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
745d6088 32#include "AliESDtrackCuts.h"
0f67a57c 33#include "AliPWG0Helper.h"
34#include "AliCorrection.h"
35#include "AliCorrectionMatrix3D.h"
36#include "dNdEta/dNdEtaAnalysis.h"
37
38ClassImp(AlidNdEtaTask)
39
40AlidNdEtaTask::AlidNdEtaTask(const char* opt) :
41 AliAnalysisTask("AlidNdEtaTask", ""),
42 fESD(0),
43 fOutput(0),
44 fOption(opt),
770a1f1d 45 fAnalysisMode(AliPWG0Helper::kTPC),
0fc41645 46 fTrigger(AliPWG0Helper::kMB1),
0f67a57c 47 fReadMC(kFALSE),
54b096ef 48 fUseMCVertex(kFALSE),
3d7758c1 49 fUseMCKine(kFALSE),
0f67a57c 50 fEsdTrackCuts(0),
51 fdNdEtaAnalysisESD(0),
52 fMult(0),
53 fMultVtx(0),
54 fEvents(0),
54b096ef 55 fVertexResolution(0),
0f67a57c 56 fdNdEtaAnalysis(0),
0fc41645 57 fdNdEtaAnalysisNSD(0),
0f67a57c 58 fdNdEtaAnalysisTr(0),
59 fdNdEtaAnalysisTrVtx(0),
60 fdNdEtaAnalysisTracks(0),
61 fVertex(0),
54b096ef 62 fPartPt(0),
63 fDeltaPhi(0)
0f67a57c 64{
65 //
66 // Constructor. Initialization of pointers
67 //
68
69 // Define input and output slots here
70 DefineInput(0, TChain::Class());
71 DefineOutput(0, TList::Class());
72}
73
74AlidNdEtaTask::~AlidNdEtaTask()
75{
76 //
77 // Destructor
78 //
79
80 // histograms are in the output list and deleted when the output
81 // list is deleted by the TSelector dtor
82
83 if (fOutput) {
84 delete fOutput;
85 fOutput = 0;
86 }
87}
88
89//________________________________________________________________________
90void AlidNdEtaTask::ConnectInputData(Option_t *)
91{
92 // Connect ESD
93 // Called once
94
95 Printf("AlidNdEtaTask::ConnectInputData called");
96
97 TTree* tree = dynamic_cast<TTree*> (GetInputData(0));
98 if (!tree) {
99 Printf("ERROR: Could not read tree from input slot 0");
100 } else {
101 // Disable all branches and enable only the needed ones
102 //tree->SetBranchStatus("*", 0);
103
0fc41645 104 tree->SetBranchStatus("TriggerMask", 1);
105 tree->SetBranchStatus("SPDVertex*", 1);
106 tree->SetBranchStatus("PrimaryVertex*", 1);
107 tree->SetBranchStatus("TPCVertex*", 1);
0f67a57c 108
54b096ef 109 if (fAnalysisMode == AliPWG0Helper::kSPD) {
0fc41645 110 tree->SetBranchStatus("AliMultiplicity*", 1);
54b096ef 111 }
0fc41645 112
770a1f1d 113 if (fAnalysisMode == AliPWG0Helper::kTPC || fAnalysisMode == AliPWG0Helper::kTPCITS) {
0fc41645 114 //AliESDtrackCuts::EnableNeededBranches(tree);
115 tree->SetBranchStatus("Tracks*", 1);
0f67a57c 116 }
117
0f67a57c 118 AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
119
120 if (!esdH) {
121 Printf("ERROR: Could not get ESDInputHandler");
122 } else
123 fESD = esdH->GetEvent();
124 }
54b096ef 125
126 // disable info messages of AliMCEvent (per event)
127 AliLog::SetClassDebugLevel("AliMCEvent", AliLog::kWarning - AliLog::kDebug + 1);
0f67a57c 128}
129
130void AlidNdEtaTask::CreateOutputObjects()
131{
132 // create result objects and add to output list
133
54b096ef 134 Printf("AlidNdEtaTask::CreateOutputObjects");
135
0f67a57c 136 fOutput = new TList;
137 fOutput->SetOwner();
138
770a1f1d 139 fdNdEtaAnalysisESD = new dNdEtaAnalysis("fdNdEtaAnalysisESD", "fdNdEtaAnalysisESD", fAnalysisMode);
0f67a57c 140 fOutput->Add(fdNdEtaAnalysisESD);
141
142 fMult = new TH1F("fMult", "fMult;Ntracks;Count", 201, -0.5, 200.5);
143 fOutput->Add(fMult);
144
145 fMultVtx = new TH1F("fMultVtx", "fMultVtx;Ntracks;Count", 201, -0.5, 200.5);
146 fOutput->Add(fMultVtx);
147
148 for (Int_t i=0; i<3; ++i)
149 {
150 fPartEta[i] = new TH1F(Form("dndeta_check_%d", i), Form("dndeta_check_%d", i), 60, -6, 6);
151 fPartEta[i]->Sumw2();
152 fOutput->Add(fPartEta[i]);
153 }
154
0fc41645 155 fEvents = new TH1F("dndeta_check_vertex", "dndeta_check_vertex", 160, -40, 40);
0f67a57c 156 fOutput->Add(fEvents);
157
54b096ef 158 fVertexResolution = new TH1F("dndeta_vertex_resolution_z", "dndeta_vertex_resolution_z", 1000, 0, 10);
159 fOutput->Add(fVertexResolution);
160
0f67a57c 161 if (fReadMC)
162 {
770a1f1d 163 fdNdEtaAnalysis = new dNdEtaAnalysis("dndeta", "dndeta", fAnalysisMode);
0f67a57c 164 fOutput->Add(fdNdEtaAnalysis);
165
0fc41645 166 fdNdEtaAnalysisNSD = new dNdEtaAnalysis("dndetaNSD", "dndetaNSD", fAnalysisMode);
167 fOutput->Add(fdNdEtaAnalysisNSD);
168
770a1f1d 169 fdNdEtaAnalysisTr = new dNdEtaAnalysis("dndetaTr", "dndetaTr", fAnalysisMode);
0f67a57c 170 fOutput->Add(fdNdEtaAnalysisTr);
171
770a1f1d 172 fdNdEtaAnalysisTrVtx = new dNdEtaAnalysis("dndetaTrVtx", "dndetaTrVtx", fAnalysisMode);
0f67a57c 173 fOutput->Add(fdNdEtaAnalysisTrVtx);
174
770a1f1d 175 fdNdEtaAnalysisTracks = new dNdEtaAnalysis("dndetaTracks", "dndetaTracks", fAnalysisMode);
0f67a57c 176 fOutput->Add(fdNdEtaAnalysisTracks);
177
178 fVertex = new TH3F("vertex_check", "vertex_check", 50, -50, 50, 50, -50, 50, 50, -50, 50);
179 fOutput->Add(fVertex);
180
181 fPartPt = new TH1F("dndeta_check_pt", "dndeta_check_pt", 1000, 0, 10);
182 fPartPt->Sumw2();
183 fOutput->Add(fPartPt);
184 }
54b096ef 185
186 if (fAnalysisMode == AliPWG0Helper::kSPD)
187 fDeltaPhi = new TH1F("fDeltaPhi", "fDeltaPhi;#Delta #phi;Entries", 1000, -3.14, 3.14);
0f67a57c 188}
189
190void AlidNdEtaTask::Exec(Option_t*)
191{
192 // process the event
193
194 // Check prerequisites
195 if (!fESD)
196 {
197 AliDebug(AliLog::kError, "ESD branch not available");
198 return;
199 }
200
201 // trigger definition
0fc41645 202 Bool_t eventTriggered = AliPWG0Helper::IsEventTriggered(fESD->GetTriggerMask(), fTrigger);
0f67a57c 203
0f67a57c 204 // get the ESD vertex
770a1f1d 205 const AliESDVertex* vtxESD = AliPWG0Helper::GetVertex(fESD, fAnalysisMode);
206
0f67a57c 207 Double_t vtx[3];
54b096ef 208 if (vtxESD) {
770a1f1d 209 vtxESD->GetXYZ(vtx);
54b096ef 210 }
211 else
212 Printf("No vertex");
213
214 // fill z vertex resolution
215 if (vtxESD)
216 fVertexResolution->Fill(vtxESD->GetZRes());
0f67a57c 217
218 // post the data already here
219 PostData(0, fOutput);
220
3d7758c1 221 // needed for syst. studies
222 AliStack* stack = 0;
223
224 if (fUseMCVertex || fUseMCKine) {
225 if (!fReadMC) {
226 Printf("ERROR: fUseMCVertex or fUseMCKine set without fReadMC set!");
227 return;
228 }
229
230 AliMCEventHandler* eventHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
231 if (!eventHandler) {
232 Printf("ERROR: Could not retrieve MC event handler");
233 return;
234 }
235
236 AliMCEvent* mcEvent = eventHandler->MCEvent();
237 if (!mcEvent) {
238 Printf("ERROR: Could not retrieve MC event");
239 return;
240 }
241
242 AliHeader* header = mcEvent->Header();
243 if (!header)
244 {
245 AliDebug(AliLog::kError, "Header not available");
246 return;
247 }
248
249 if (fUseMCVertex)
250 {
251 Printf("WARNING: Replacing vertex by MC vertex. This is for systematical checks only.");
252 // get the MC vertex
253 AliGenEventHeader* genHeader = header->GenEventHeader();
254 TArrayF vtxMC(3);
255 genHeader->PrimaryVertex(vtxMC);
256
257 vtx[2] = vtxMC[2];
258 }
259
260 if (fUseMCKine)
261 {
3d7758c1 262 stack = mcEvent->Stack();
263 if (!stack)
264 {
265 AliDebug(AliLog::kError, "Stack not available");
266 return;
267 }
268 }
269 }
270
0f67a57c 271 // create list of (label, eta, pt) tuples
272 Int_t inputCount = 0;
273 Int_t* labelArr = 0;
274 Float_t* etaArr = 0;
275 Float_t* ptArr = 0;
770a1f1d 276 if (fAnalysisMode == AliPWG0Helper::kSPD)
0f67a57c 277 {
278 // get tracklets
279 const AliMultiplicity* mult = fESD->GetMultiplicity();
280 if (!mult)
281 {
282 AliDebug(AliLog::kError, "AliMultiplicity not available");
283 return;
284 }
285
286 labelArr = new Int_t[mult->GetNumberOfTracklets()];
287 etaArr = new Float_t[mult->GetNumberOfTracklets()];
288 ptArr = new Float_t[mult->GetNumberOfTracklets()];
289
3d7758c1 290 if (fUseMCKine && stack)
291 Printf("Processing only primaries (MC information used). This is for systematical checks only.");
292
0f67a57c 293 // get multiplicity from ITS tracklets
294 for (Int_t i=0; i<mult->GetNumberOfTracklets(); ++i)
295 {
296 //printf("%d %f %f %f\n", i, mult->GetTheta(i), mult->GetPhi(i), mult->GetDeltaPhi(i));
297
3d7758c1 298 if (fUseMCKine && stack)
299 if (mult->GetLabel(i, 0) < 0 || mult->GetLabel(i, 0) != mult->GetLabel(i, 1) || !stack->IsPhysicalPrimary(mult->GetLabel(i, 0)))
300 continue;
301
302 Float_t deltaPhi = mult->GetDeltaPhi(i);
303 // prevent values to be shifted by 2 Pi()
304 if (deltaPhi < -TMath::Pi())
305 deltaPhi += TMath::Pi() * 2;
306 if (deltaPhi > TMath::Pi())
307 deltaPhi -= TMath::Pi() * 2;
308
309 if (TMath::Abs(deltaPhi) > 1)
310 printf("WARNING: Very high Delta Phi: %d %f %f %f\n", i, mult->GetTheta(i), mult->GetPhi(i), deltaPhi);
0f67a57c 311
3d7758c1 312 fDeltaPhi->Fill(deltaPhi);
54b096ef 313
0f67a57c 314 etaArr[inputCount] = mult->GetEta(i);
3d7758c1 315 labelArr[inputCount] = mult->GetLabel(i, 0);
0f67a57c 316 ptArr[inputCount] = 0; // no pt for tracklets
317 ++inputCount;
318 }
54b096ef 319
320 Printf("Accepted %d tracklets", inputCount);
0f67a57c 321 }
770a1f1d 322 else if (fAnalysisMode == AliPWG0Helper::kTPC || fAnalysisMode == AliPWG0Helper::kTPCITS)
0f67a57c 323 {
324 if (!fEsdTrackCuts)
325 {
326 AliDebug(AliLog::kError, "fESDTrackCuts not available");
327 return;
328 }
329
330 // get multiplicity from ESD tracks
331 TObjArray* list = fEsdTrackCuts->GetAcceptedTracks(fESD);
332 Int_t nGoodTracks = list->GetEntries();
333
334 labelArr = new Int_t[nGoodTracks];
335 etaArr = new Float_t[nGoodTracks];
336 ptArr = new Float_t[nGoodTracks];
337
338 // loop over esd tracks
339 for (Int_t i=0; i<nGoodTracks; i++)
340 {
341 AliESDtrack* esdTrack = dynamic_cast<AliESDtrack*> (list->At(i));
342 if (!esdTrack)
343 {
344 AliDebug(AliLog::kError, Form("ERROR: Could not retrieve track %d.", i));
345 continue;
346 }
347
348 etaArr[inputCount] = esdTrack->Eta();
349 labelArr[inputCount] = TMath::Abs(esdTrack->GetLabel());
350 ptArr[inputCount] = esdTrack->Pt();
351 ++inputCount;
352 }
2b26296f 353
54b096ef 354 Printf("Accepted %d tracks", nGoodTracks);
355
2b26296f 356 delete list;
0f67a57c 357 }
358 else
359 return;
360
361 // Processing of ESD information (always)
362 if (eventTriggered)
363 {
dd367a14 364 // control hist
0f67a57c 365 fMult->Fill(inputCount);
7307d52c 366 fdNdEtaAnalysisESD->FillTriggeredEvent(inputCount);
367
770a1f1d 368 if (vtxESD)
dd367a14 369 {
370 // control hist
0f67a57c 371 fMultVtx->Fill(inputCount);
372
dd367a14 373 for (Int_t i=0; i<inputCount; ++i)
374 {
375 Float_t eta = etaArr[i];
376 Float_t pt = ptArr[i];
0f67a57c 377
dd367a14 378 fdNdEtaAnalysisESD->FillTrack(vtx[2], eta, pt);
0f67a57c 379
dd367a14 380 if (TMath::Abs(vtx[2]) < 20)
381 {
382 fPartEta[0]->Fill(eta);
0f67a57c 383
dd367a14 384 if (vtx[2] < 0)
385 fPartEta[1]->Fill(eta);
386 else
387 fPartEta[2]->Fill(eta);
388 }
0f67a57c 389 }
0f67a57c 390
dd367a14 391 // for event count per vertex
392 fdNdEtaAnalysisESD->FillEvent(vtx[2], inputCount);
0f67a57c 393
54b096ef 394 // control hist
dd367a14 395 fEvents->Fill(vtx[2]);
396 }
0f67a57c 397 }
398
399 if (fReadMC) // Processing of MC information (optional)
400 {
401 AliMCEventHandler* eventHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
402 if (!eventHandler) {
403 Printf("ERROR: Could not retrieve MC event handler");
404 return;
405 }
406
407 AliMCEvent* mcEvent = eventHandler->MCEvent();
408 if (!mcEvent) {
409 Printf("ERROR: Could not retrieve MC event");
410 return;
411 }
412
745d6088 413 stack = mcEvent->Stack();
0f67a57c 414 if (!stack)
415 {
416 AliDebug(AliLog::kError, "Stack not available");
417 return;
418 }
419
420 AliHeader* header = mcEvent->Header();
421 if (!header)
422 {
423 AliDebug(AliLog::kError, "Header not available");
424 return;
425 }
426
427 // get the MC vertex
428 AliGenEventHeader* genHeader = header->GenEventHeader();
429 TArrayF vtxMC(3);
430 genHeader->PrimaryVertex(vtxMC);
431
0fc41645 432 // get process type; NB: this only works for Pythia
433 Int_t processType = AliPWG0Helper::GetPythiaEventProcessType(header);
434 AliDebug(AliLog::kDebug+1, Form("Found pythia process type %d", processType));
435
436 if (processType<0)
437 AliDebug(AliLog::kError, Form("Unknown Pythia process type %d.", processType));
438
0f67a57c 439 // loop over mc particles
440 Int_t nPrim = stack->GetNprimary();
441
442 Int_t nAcceptedParticles = 0;
443
444 for (Int_t iMc = 0; iMc < nPrim; ++iMc)
445 {
446 //Printf("Getting particle %d", iMc);
447 TParticle* particle = stack->Particle(iMc);
448
449 if (!particle)
450 continue;
451
452 if (AliPWG0Helper::IsPrimaryCharged(particle, nPrim) == kFALSE)
453 continue;
454
455 AliDebug(AliLog::kDebug+1, Form("Accepted primary %d, unique ID: %d", iMc, particle->GetUniqueID()));
0f67a57c 456 Float_t eta = particle->Eta();
457 Float_t pt = particle->Pt();
458
54b096ef 459 // make a rough eta cut (so that nAcceptedParticles is not too far off the true measured value (NB: this histograms are only gathered for comparison))
460 if (TMath::Abs(eta) < 1.5 && pt > 0.3)
2b26296f 461 nAcceptedParticles++;
462
0f67a57c 463 fdNdEtaAnalysis->FillTrack(vtxMC[2], eta, pt);
464 fVertex->Fill(particle->Vx(), particle->Vy(), particle->Vz());
465
0fc41645 466 if (processType != 92 && processType != 93)
467 fdNdEtaAnalysisNSD->FillTrack(vtxMC[2], eta, pt);
468
0f67a57c 469 if (eventTriggered)
470 {
471 fdNdEtaAnalysisTr->FillTrack(vtxMC[2], eta, pt);
770a1f1d 472 if (vtxESD)
0f67a57c 473 fdNdEtaAnalysisTrVtx->FillTrack(vtxMC[2], eta, pt);
474 }
475
476 if (TMath::Abs(eta) < 0.8)
477 fPartPt->Fill(pt);
478 }
479
480 fdNdEtaAnalysis->FillEvent(vtxMC[2], nAcceptedParticles);
0fc41645 481 if (processType != 92 && processType != 93)
482 fdNdEtaAnalysisNSD->FillEvent(vtxMC[2], nAcceptedParticles);
483
0f67a57c 484 if (eventTriggered)
485 {
486 fdNdEtaAnalysisTr->FillEvent(vtxMC[2], nAcceptedParticles);
770a1f1d 487 if (vtxESD)
0f67a57c 488 fdNdEtaAnalysisTrVtx->FillEvent(vtxMC[2], nAcceptedParticles);
489 }
490
770a1f1d 491 if (eventTriggered && vtxESD)
0f67a57c 492 {
493 // from tracks is only done for triggered and vertex reconstructed events
494
495 for (Int_t i=0; i<inputCount; ++i)
496 {
497 Int_t label = labelArr[i];
498
499 if (label < 0)
500 continue;
501
502 //Printf("Getting particle of track %d", label);
503 TParticle* particle = stack->Particle(label);
504
505 if (!particle)
506 {
507 AliDebug(AliLog::kError, Form("ERROR: Could not retrieve particle %d.", label));
508 continue;
509 }
510
511 fdNdEtaAnalysisTracks->FillTrack(vtxMC[2], particle->Eta(), particle->Pt());
512 } // end of track loop
513
514 // for event count per vertex
515 fdNdEtaAnalysisTracks->FillEvent(vtxMC[2], inputCount);
516 }
517 }
518
519 delete[] etaArr;
520 delete[] labelArr;
521 delete[] ptArr;
522}
523
524void AlidNdEtaTask::Terminate(Option_t *)
525{
526 // The Terminate() function is the last function to be called during
527 // a query. It always runs on the client, it can be used to present
528 // the results graphically or save the results to file.
529
530 fOutput = dynamic_cast<TList*> (GetOutputData(0));
531 if (!fOutput) {
532 Printf("ERROR: fOutput not available");
533 return;
534 }
535
536 fdNdEtaAnalysisESD = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("fdNdEtaAnalysisESD"));
537 fMult = dynamic_cast<TH1F*> (fOutput->FindObject("fMult"));
538 fMultVtx = dynamic_cast<TH1F*> (fOutput->FindObject("fMultVtx"));
539 fEvents = dynamic_cast<TH1F*> (fOutput->FindObject("dndeta_check_vertex"));
54b096ef 540 fVertexResolution = dynamic_cast<TH1F*> (fOutput->FindObject("dndeta_vertex_resolution_z"));
0f67a57c 541 for (Int_t i=0; i<3; ++i)
542 fPartEta[i] = dynamic_cast<TH1F*> (fOutput->FindObject(Form("dndeta_check_%d", i)));
543
544 if (!fdNdEtaAnalysisESD)
545 {
546 AliDebug(AliLog::kError, "ERROR: fdNdEtaAnalysisESD not available");
547 return;
548 }
549
550 if (fMult && fMultVtx)
551 {
552 new TCanvas;
553 fMult->Draw();
554 fMultVtx->SetLineColor(2);
770a1f1d 555 fMultVtx->Draw("SAME");
0f67a57c 556 }
557
558 if (fPartEta[0])
559 {
560 Int_t events1 = (Int_t) fEvents->Integral(fEvents->GetXaxis()->FindBin(-19.9), fEvents->GetXaxis()->FindBin(-0.001));
561 Int_t events2 = (Int_t) fEvents->Integral(fEvents->GetXaxis()->FindBin(0.001), fEvents->GetXaxis()->FindBin(19.9));
562
54b096ef 563 Printf("%d events with vertex used", events1 + events2);
564
565 if (events1 > 0 && events2 > 0)
566 {
0f67a57c 567 fPartEta[0]->Scale(1.0 / (events1 + events2));
568 fPartEta[1]->Scale(1.0 / events1);
569 fPartEta[2]->Scale(1.0 / events2);
570
571 for (Int_t i=0; i<3; ++i)
572 fPartEta[i]->Scale(1.0 / fPartEta[i]->GetBinWidth(1));
573
574 new TCanvas("control", "control", 500, 500);
575 for (Int_t i=0; i<3; ++i)
576 {
577 fPartEta[i]->SetLineColor(i+1);
578 fPartEta[i]->Draw((i==0) ? "" : "SAME");
579 }
54b096ef 580 }
0f67a57c 581 }
582
583 if (fEvents)
584 {
585 new TCanvas("control3", "control3", 500, 500);
770a1f1d 586 fEvents->Draw();
0f67a57c 587 }
588
589 TFile* fout = new TFile("analysis_esd_raw.root", "RECREATE");
590
591 if (fdNdEtaAnalysisESD)
592 fdNdEtaAnalysisESD->SaveHistograms();
593
594 if (fEsdTrackCuts)
595 fEsdTrackCuts->SaveHistograms("esd_tracks_cuts");
596
597 if (fMult)
598 fMult->Write();
599
600 if (fMultVtx)
601 fMultVtx->Write();
602
603 for (Int_t i=0; i<3; ++i)
604 if (fPartEta[i])
605 fPartEta[i]->Write();
606
607 if (fEvents)
608 fEvents->Write();
609
54b096ef 610 if (fVertexResolution)
611 fVertexResolution->Write();
612
613 if (fDeltaPhi)
614 fDeltaPhi->Write();
615
0f67a57c 616 fout->Write();
617 fout->Close();
618
619 Printf("Writting result to analysis_esd_raw.root");
620
621 if (fReadMC)
622 {
623 fdNdEtaAnalysis = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndeta"));
0fc41645 624 fdNdEtaAnalysisNSD = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndetaNSD"));
0f67a57c 625 fdNdEtaAnalysisTr = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndetaTr"));
626 fdNdEtaAnalysisTrVtx = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndetaTrVtx"));
627 fdNdEtaAnalysisTracks = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndetaTracks"));
628 fPartPt = dynamic_cast<TH1F*> (fOutput->FindObject("dndeta_check_pt"));
629
630 if (!fdNdEtaAnalysis || !fdNdEtaAnalysisTr || !fdNdEtaAnalysisTrVtx || !fPartPt)
631 {
632 AliDebug(AliLog::kError, Form("ERROR: Histograms not available %p %p", (void*) fdNdEtaAnalysis, (void*) fPartPt));
633 return;
634 }
635
636 fdNdEtaAnalysis->Finish(0, -1, AlidNdEtaCorrection::kNone);
0fc41645 637 fdNdEtaAnalysisNSD->Finish(0, -1, AlidNdEtaCorrection::kNone);
0f67a57c 638 fdNdEtaAnalysisTr->Finish(0, -1, AlidNdEtaCorrection::kNone);
639 fdNdEtaAnalysisTrVtx->Finish(0, -1, AlidNdEtaCorrection::kNone);
640 fdNdEtaAnalysisTracks->Finish(0, -1, AlidNdEtaCorrection::kNone);
641
642 Int_t events = (Int_t) fdNdEtaAnalysis->GetData()->GetEventCorrection()->GetMeasuredHistogram()->Integral();
643 fPartPt->Scale(1.0/events);
644 fPartPt->Scale(1.0/fPartPt->GetBinWidth(1));
645
745d6088 646 fout = new TFile("analysis_mc.root","RECREATE");
0f67a57c 647
648 fdNdEtaAnalysis->SaveHistograms();
0fc41645 649 fdNdEtaAnalysisNSD->SaveHistograms();
0f67a57c 650 fdNdEtaAnalysisTr->SaveHistograms();
651 fdNdEtaAnalysisTrVtx->SaveHistograms();
652 fdNdEtaAnalysisTracks->SaveHistograms();
653 fPartPt->Write();
654
655 fout->Write();
656 fout->Close();
657
658 Printf("Writting result to analysis_mc.root");
659
660 if (fPartPt)
661 {
662 new TCanvas("control2", "control2", 500, 500);
770a1f1d 663 fPartPt->Draw();
0f67a57c 664 }
665 }
666}