]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG0/dNdEta/AlidNdEtaTask.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
- macro for creating OCDB configuration objects of the GlobalHisto component
[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>
68fa248f 12#include <TH1D.h>
0f67a57c 13#include <TH2F.h>
14#include <TH3F.h>
15#include <TParticle.h>
16#include <TRandom.h>
17#include <TNtuple.h>
18#include <TObjString.h>
19#include <TF1.h>
69b09e3b 20#include <TGraph.h>
0f67a57c 21
22#include <AliLog.h>
23#include <AliESDVertex.h>
24#include <AliESDEvent.h>
25#include <AliStack.h>
26#include <AliHeader.h>
27#include <AliGenEventHeader.h>
28#include <AliMultiplicity.h>
29#include <AliAnalysisManager.h>
30#include <AliMCEventHandler.h>
31#include <AliMCEvent.h>
32#include <AliESDInputHandler.h>
1d107532 33#include <AliESDInputHandlerRP.h>
69b09e3b 34#include <AliESDHeader.h>
0f67a57c 35
745d6088 36#include "AliESDtrackCuts.h"
0f67a57c 37#include "AliPWG0Helper.h"
38#include "AliCorrection.h"
39#include "AliCorrectionMatrix3D.h"
40#include "dNdEta/dNdEtaAnalysis.h"
70fdd197 41#include "AliTriggerAnalysis.h"
12bb57f1 42#include "AliPhysicsSelection.h"
43
44//#define FULLALIROOT
45
46#ifdef FULLALIROOT
47 #include "../ITS/AliITSRecPoint.h"
48 #include "AliCDBManager.h"
49 #include "AliCDBEntry.h"
50 #include "AliGeomManager.h"
51 #include "TGeoManager.h"
52#endif
0f67a57c 53
54ClassImp(AlidNdEtaTask)
55
56AlidNdEtaTask::AlidNdEtaTask(const char* opt) :
76532b17 57 AliAnalysisTaskSE("AlidNdEtaTask"),
0f67a57c 58 fESD(0),
59 fOutput(0),
60 fOption(opt),
a7f69e56 61 fAnalysisMode((AliPWG0Helper::AnalysisMode) (AliPWG0Helper::kTPC | AliPWG0Helper::kFieldOn)),
70fdd197 62 fTrigger(AliTriggerAnalysis::kMB1),
69b09e3b 63 fFillPhi(kFALSE),
64 fDeltaPhiCut(-1),
0f67a57c 65 fReadMC(kFALSE),
54b096ef 66 fUseMCVertex(kFALSE),
c17301f3 67 fOnlyPrimaries(kFALSE),
3d7758c1 68 fUseMCKine(kFALSE),
1d107532 69 fSymmetrize(kFALSE),
81be4ee8 70 fMultAxisEta1(kFALSE),
71 fDiffTreatment(AliPWG0Helper::kMCFlags),
0f67a57c 72 fEsdTrackCuts(0),
73 fdNdEtaAnalysisESD(0),
74 fMult(0),
75 fMultVtx(0),
76 fEvents(0),
54b096ef 77 fVertexResolution(0),
0f67a57c 78 fdNdEtaAnalysis(0),
69b09e3b 79 fdNdEtaAnalysisND(0),
0fc41645 80 fdNdEtaAnalysisNSD(0),
81be4ee8 81 fdNdEtaAnalysisOnePart(0),
0f67a57c 82 fdNdEtaAnalysisTr(0),
83 fdNdEtaAnalysisTrVtx(0),
84 fdNdEtaAnalysisTracks(0),
54b096ef 85 fPartPt(0),
51f6de65 86 fVertex(0),
1d107532 87 fVertexVsMult(0),
ea441adf 88 fPhi(0),
69b09e3b 89 fRawPt(0),
ea441adf 90 fEtaPhi(0),
12bb57f1 91 fModuleMap(0),
69b09e3b 92 fDeltaPhi(0),
a7f69e56 93 fDeltaTheta(0),
69b09e3b 94 fFiredChips(0),
1d107532 95 fTrackletsVsClusters(0),
96 fTrackletsVsUnassigned(0),
1d107532 97 fStats(0),
68fa248f 98 fStats2(0),
99 fPtMin(0.15),
100 fEta(0x0),
101 fEtaMC(0x0),
102 fHistEvents(0),
103 fHistEventsMC(0),
104 fTrigEffNum(0),
105 fTrigEffDen(0),
106 fVtxEffNum(0),
107 fVtxEffDen(0),
108 fVtxTrigEffNum(0),
109 fVtxTrigEffDen(0)
0f67a57c 110{
111 //
112 // Constructor. Initialization of pointers
113 //
114
115 // Define input and output slots here
76532b17 116 //DefineInput(0, TChain::Class());
117 DefineOutput(1, TList::Class());
69b09e3b 118
119 fZPhi[0] = 0;
120 fZPhi[1] = 0;
567160d6 121
122 AliLog::SetClassDebugLevel("AlidNdEtaTask", AliLog::kWarning);
0f67a57c 123}
124
125AlidNdEtaTask::~AlidNdEtaTask()
126{
127 //
128 // Destructor
129 //
130
131 // histograms are in the output list and deleted when the output
132 // list is deleted by the TSelector dtor
133
1103419a 134 if (fOutput && !AliAnalysisManager::GetAnalysisManager()->IsProofMode()) {
0f67a57c 135 delete fOutput;
136 fOutput = 0;
137 }
138}
139
76532b17 140Bool_t AlidNdEtaTask::UserNotify()
567160d6 141{
142 static Int_t count = 0;
143 count++;
69b09e3b 144 Printf("Processing %d. file: %s", count, ((TTree*) GetInputData(0))->GetCurrentFile()->GetName());
567160d6 145 return kTRUE;
146}
147
0f67a57c 148//________________________________________________________________________
76532b17 149void AlidNdEtaTask::ConnectInputData(Option_t *opt)
0f67a57c 150{
151 // Connect ESD
152 // Called once
76532b17 153
154 AliAnalysisTaskSE::ConnectInputData(opt);
0f67a57c 155
156 Printf("AlidNdEtaTask::ConnectInputData called");
157
76532b17 158 /*
567160d6 159 AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
160
161 if (!esdH) {
162 Printf("ERROR: Could not get ESDInputHandler");
0f67a57c 163 } else {
567160d6 164 fESD = esdH->GetEvent();
a7f69e56 165
166 TString branches("AliESDHeader Vertex ");
0f67a57c 167
70fdd197 168 if (fAnalysisMode & AliPWG0Helper::kSPD || fTrigger & AliTriggerAnalysis::kOfflineFlag)
a7f69e56 169 branches += "AliMultiplicity ";
170
171 if (fAnalysisMode & AliPWG0Helper::kTPC || fAnalysisMode & AliPWG0Helper::kTPCITS)
172 branches += "Tracks ";
173
c17301f3 174 // Enable only the needed branches
a7f69e56 175 esdH->SetActiveBranches(branches);
0f67a57c 176 }
76532b17 177 */
54b096ef 178
179 // disable info messages of AliMCEvent (per event)
180 AliLog::SetClassDebugLevel("AliMCEvent", AliLog::kWarning - AliLog::kDebug + 1);
12bb57f1 181
182 #ifdef FULLALIROOT
02e1bcdb 183 AliCDBManager::Instance()->SetDefaultStorage("raw://");
184 //AliCDBManager::Instance()->SetDefaultStorage("MC", "Residual");
12bb57f1 185 AliCDBManager::Instance()->SetRun(0);
186
187 AliCDBManager* mgr = AliCDBManager::Instance();
188 AliCDBEntry* obj = mgr->Get(AliCDBPath("GRP", "Geometry", "Data"));
189 AliGeomManager::SetGeometry((TGeoManager*) obj->GetObject());
190
191 AliGeomManager::GetNalignable("ITS");
192 AliGeomManager::ApplyAlignObjsFromCDB("ITS");
193 #endif
0f67a57c 194}
195
76532b17 196void AlidNdEtaTask::UserCreateOutputObjects()
0f67a57c 197{
198 // create result objects and add to output list
199
54b096ef 200 Printf("AlidNdEtaTask::CreateOutputObjects");
1103419a 201 //AliLog::SetClassDebugLevel("AliPhysicsSelection", AliLog::kDebug);
54b096ef 202
5a6310fe 203 if (fOnlyPrimaries)
204 Printf("WARNING: Processing only primaries (MC information used). This is for systematical checks only.");
205
206 if (fUseMCKine)
207 Printf("WARNING: Using MC kine information. This is for systematical checks only.");
208
209 if (fUseMCVertex)
210 Printf("WARNING: Replacing vertex by MC vertex. This is for systematical checks only.");
211
0f67a57c 212 fOutput = new TList;
213 fOutput->SetOwner();
214
770a1f1d 215 fdNdEtaAnalysisESD = new dNdEtaAnalysis("fdNdEtaAnalysisESD", "fdNdEtaAnalysisESD", fAnalysisMode);
0f67a57c 216 fOutput->Add(fdNdEtaAnalysisESD);
217
218 fMult = new TH1F("fMult", "fMult;Ntracks;Count", 201, -0.5, 200.5);
219 fOutput->Add(fMult);
220
221 fMultVtx = new TH1F("fMultVtx", "fMultVtx;Ntracks;Count", 201, -0.5, 200.5);
222 fOutput->Add(fMultVtx);
223
224 for (Int_t i=0; i<3; ++i)
225 {
81be4ee8 226 fPartEta[i] = new TH1F(Form("dndeta_check_%d", i), Form("dndeta_check_%d", i), 60, -3, 3);
0f67a57c 227 fPartEta[i]->Sumw2();
228 fOutput->Add(fPartEta[i]);
229 }
230
51f6de65 231 fEvents = new TH1F("dndeta_check_vertex", "dndeta_check_vertex", 800, -40, 40);
0f67a57c 232 fOutput->Add(fEvents);
233
81be4ee8 234 fVertexResolution = new TH2F("dndeta_vertex_resolution_z", ";z resolution;#Delta #phi", 1000, 0, 2, 100, 0, 0.2);
54b096ef 235 fOutput->Add(fVertexResolution);
236
ea441adf 237 fPhi = new TH1F("fPhi", "fPhi;#phi in rad.;count", 720, 0, 2 * TMath::Pi());
238 fOutput->Add(fPhi);
239
81be4ee8 240 fEtaPhi = new TH2F("fEtaPhi", "fEtaPhi;#eta;#phi in rad.;count", 80, -4, 4, 18*20, 0, 2 * TMath::Pi());
ea441adf 241 fOutput->Add(fEtaPhi);
12bb57f1 242
1d107532 243 fStats = new TH1F("fStats", "fStats", 5, 0.5, 5.5);
69b09e3b 244 fStats->GetXaxis()->SetBinLabel(1, "vertexer 3d");
245 fStats->GetXaxis()->SetBinLabel(2, "vertexer z");
a7f69e56 246 fStats->GetXaxis()->SetBinLabel(3, "trigger");
1d107532 247 fStats->GetXaxis()->SetBinLabel(4, "physics events");
248 fStats->GetXaxis()->SetBinLabel(5, "physics events after veto");
69b09e3b 249 fOutput->Add(fStats);
1d107532 250
76532b17 251 fStats2 = new TH2F("fStats2", "fStats2", 7, -0.5, 6.5, 12, -0.5, 11.5);
12bb57f1 252
253 fStats2->GetXaxis()->SetBinLabel(1, "No trigger");
81be4ee8 254 fStats2->GetXaxis()->SetBinLabel(2, "No Vertex");
255 fStats2->GetXaxis()->SetBinLabel(3, "Vtx res. too bad");
256 fStats2->GetXaxis()->SetBinLabel(4, "|z-vtx| > 15");
257 fStats2->GetXaxis()->SetBinLabel(5, "|z-vtx| > 10");
258 fStats2->GetXaxis()->SetBinLabel(6, "0 tracklets");
12bb57f1 259 fStats2->GetXaxis()->SetBinLabel(7, "Selected");
260
1d107532 261 fStats2->GetYaxis()->SetBinLabel(1, "n/a");
12bb57f1 262 fStats2->GetYaxis()->SetBinLabel(2, "empty");
263 fStats2->GetYaxis()->SetBinLabel(3, "BBA");
264 fStats2->GetYaxis()->SetBinLabel(4, "BBC");
265 fStats2->GetYaxis()->SetBinLabel(5, "BBA BBC");
1d107532 266 fStats2->GetYaxis()->SetBinLabel(6, "BGA");
267 fStats2->GetYaxis()->SetBinLabel(7, "BGC");
12bb57f1 268 fStats2->GetYaxis()->SetBinLabel(8, "BGA BGC");
269 fStats2->GetYaxis()->SetBinLabel(9, "BBA BGC");
270 fStats2->GetYaxis()->SetBinLabel(10, "BGA BBC");
76532b17 271 fStats2->GetYaxis()->SetBinLabel(11, "GFO");
272 fStats2->GetYaxis()->SetBinLabel(12, "NO GFO");
1d107532 273 fOutput->Add(fStats2);
69b09e3b 274
12bb57f1 275 fTrackletsVsClusters = new TH2F("fTrackletsVsClusters", ";tracklets;clusters in ITS", 50, -0.5, 49.5, 1000, -0.5, 999.5);
276 fOutput->Add(fTrackletsVsClusters);
277
a7f69e56 278 if (fAnalysisMode & AliPWG0Helper::kSPD)
69b09e3b 279 {
81be4ee8 280 fDeltaPhi = new TH1F("fDeltaPhi", "fDeltaPhi;#Delta #phi;Entries", 500, -0.16, 0.16);
69b09e3b 281 fOutput->Add(fDeltaPhi);
81be4ee8 282 fDeltaTheta = new TH1F("fDeltaTheta", "fDeltaTheta;#Delta #theta;Entries", 500, -0.05, 0.05);
a7f69e56 283 fOutput->Add(fDeltaTheta);
69b09e3b 284 fFiredChips = new TH2F("fFiredChips", "fFiredChips;Chips L1 + L2;tracklets", 1201, -0.5, 1201, 50, -0.5, 49.5);
285 fOutput->Add(fFiredChips);
1d107532 286 fTrackletsVsUnassigned = new TH2F("fTrackletsVsUnassigned", ";tracklets;unassigned clusters in L0", 50, -0.5, 49.5, 200, -0.5, 199.5);
287 fOutput->Add(fTrackletsVsUnassigned);
69b09e3b 288 for (Int_t i=0; i<2; i++)
289 {
12bb57f1 290 fZPhi[i] = new TH2F(Form("fZPhi_%d", i), Form("fZPhi Layer %d;z (cm);#phi (rad.)", i), 200, -15, 15, 200, 0, TMath::Pi() * 2);
69b09e3b 291 fOutput->Add(fZPhi[i]);
292 }
12bb57f1 293
294 fModuleMap = new TH1F("fModuleMap", "fModuleMap;module number;cluster count", 240, -0.5, 239.5);
295 fOutput->Add(fModuleMap);
69b09e3b 296 }
297
a7f69e56 298 if (fAnalysisMode & AliPWG0Helper::kTPC || fAnalysisMode & AliPWG0Helper::kTPCITS)
69b09e3b 299 {
300 fRawPt = new TH1F("fRawPt", "raw pt;p_{T};Count", 2000, 0, 100);
301 fOutput->Add(fRawPt);
302 }
ea441adf 303
1d107532 304 fVertex = new TH3F("vertex_check", "vertex_check;x;y;z", 100, -1, 1, 100, -1, 1, 100, -30, 30);
51f6de65 305 fOutput->Add(fVertex);
1d107532 306
307 fVertexVsMult = new TH3F("fVertexVsMult", "fVertexVsMult;x;y;multiplicity", 100, -1, 1, 100, -1, 1, 100, -0.5, 99.5);
308 fOutput->Add(fVertexVsMult);
51f6de65 309
0f67a57c 310 if (fReadMC)
311 {
770a1f1d 312 fdNdEtaAnalysis = new dNdEtaAnalysis("dndeta", "dndeta", fAnalysisMode);
0f67a57c 313 fOutput->Add(fdNdEtaAnalysis);
314
69b09e3b 315 fdNdEtaAnalysisND = new dNdEtaAnalysis("dndetaND", "dndetaND", fAnalysisMode);
316 fOutput->Add(fdNdEtaAnalysisND);
317
0fc41645 318 fdNdEtaAnalysisNSD = new dNdEtaAnalysis("dndetaNSD", "dndetaNSD", fAnalysisMode);
319 fOutput->Add(fdNdEtaAnalysisNSD);
320
81be4ee8 321 fdNdEtaAnalysisOnePart = new dNdEtaAnalysis("dndetaOnePart", "dndetaOnePart", fAnalysisMode);
322 fOutput->Add(fdNdEtaAnalysisOnePart);
323
770a1f1d 324 fdNdEtaAnalysisTr = new dNdEtaAnalysis("dndetaTr", "dndetaTr", fAnalysisMode);
0f67a57c 325 fOutput->Add(fdNdEtaAnalysisTr);
326
770a1f1d 327 fdNdEtaAnalysisTrVtx = new dNdEtaAnalysis("dndetaTrVtx", "dndetaTrVtx", fAnalysisMode);
0f67a57c 328 fOutput->Add(fdNdEtaAnalysisTrVtx);
329
770a1f1d 330 fdNdEtaAnalysisTracks = new dNdEtaAnalysis("dndetaTracks", "dndetaTracks", fAnalysisMode);
0f67a57c 331 fOutput->Add(fdNdEtaAnalysisTracks);
332
0f67a57c 333 fPartPt = new TH1F("dndeta_check_pt", "dndeta_check_pt", 1000, 0, 10);
334 fPartPt->Sumw2();
335 fOutput->Add(fPartPt);
336 }
1c15d51a 337
338 if (fEsdTrackCuts)
339 {
340 fEsdTrackCuts->SetName("fEsdTrackCuts");
341 fOutput->Add(fEsdTrackCuts);
342 }
68fa248f 343
344 fEta = new TH1D("fEta", "Eta;#eta;count", 80, -4, 4);
345 fOutput->Add(fEta);
346
347 fEtaMC = new TH1D("fEtaMC", "Eta, MC;#eta;count", 80, -4, 4);
348 fOutput->Add(fEtaMC);
349
350 fHistEvents = new TH1D("fHistEvents", "N. of Events;accepted;count", 2, 0, 2);
351 fOutput->Add(fHistEvents);
352
353 fHistEventsMC = new TH1D("fHistEventsMC", "N. of MC Events;accepted;count", 2, 0, 2);
354 fOutput->Add(fHistEventsMC);
355
356 fTrigEffNum = new TH1D("fTrigEffNum", "N. of triggered events", 100,0,100);
357 fOutput->Add(fTrigEffNum);
358 fTrigEffDen = new TH1D("fTrigEffDen", "N. of true events", 100,0,100);
359 fOutput->Add(fTrigEffDen);
360 fVtxEffNum = new TH1D("fVtxEffNum", "N. of events with vtx", 100,0,100);
361 fOutput->Add(fVtxEffNum);
362 fVtxEffDen = new TH1D("fVtxEffDen", "N. of true events", 100,0,100);
363 fOutput->Add(fVtxEffDen);
364 fVtxTrigEffNum = new TH1D("fVtxTrigEffNum", "N. of triggered events with vtx", 100,0,100);
365 fOutput->Add(fVtxTrigEffNum);
366 fVtxTrigEffDen = new TH1D("fVtxTrigEffDen", "N. of triggered true events", 100,0,100);
367 fOutput->Add(fVtxTrigEffDen);
12bb57f1 368
12bb57f1 369 //AliLog::SetClassDebugLevel("AliPhysicsSelection", AliLog::kDebug);
0f67a57c 370}
371
76532b17 372Bool_t AlidNdEtaTask::IsEventInBinZero()
373{
374 // checks if the event goes to the 0 bin
375
68fa248f 376 Bool_t isZeroBin = kTRUE;
76532b17 377 fESD = (AliESDEvent*) fInputEvent;
378
379 AliInputEventHandler* inputHandler = static_cast<AliInputEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
380 if (!inputHandler)
381 {
382 Printf("ERROR: Could not receive input handler");
383 return kFALSE;
384 }
385
386 static AliTriggerAnalysis* triggerAnalysis = 0;
387 if (!triggerAnalysis)
388 {
389 AliPhysicsSelection* physicsSelection = dynamic_cast<AliPhysicsSelection*> (inputHandler->GetEventSelection());
390 if (physicsSelection)
391 triggerAnalysis = physicsSelection->GetTriggerAnalysis();
392 }
393
394 if (!triggerAnalysis)
395 {
396 Printf("ERROR: Could not receive trigger analysis object");
397 return kFALSE;
398 }
399
400 if (!triggerAnalysis->IsTriggerFired(fESD, fTrigger))
401 return kFALSE;
402
68fa248f 403 // 0 bin check - from Michele
76532b17 404
68fa248f 405 const AliMultiplicity* mult = fESD->GetMultiplicity();
406 if (!mult){
407 Printf("AlidNdEtaTask::IsBinZero: Can't get multiplicity object from ESDs");
408 return kFALSE;
409 }
410 Int_t ntracklet = mult->GetNumberOfTracklets();
411 const AliESDVertex * vtxESD = fESD->GetPrimaryVertexSPD();
412 if(vtxESD) {
413 // If there is a vertex from vertexer z with delta phi > 0.02 we
414 // don't consider it rec (we keep the event in bin0). If quality
415 // is good eneough we check the number of tracklets
416 // if the vertex is more than 15 cm away, this is autamatically bin0
417 if( TMath::Abs(vtxESD->GetZ()) <= 15 ) {
418 if (vtxESD->IsFromVertexerZ()) {
419 if (vtxESD->GetDispersion()<=0.02 ) {
420 if(ntracklet>0) isZeroBin = kFALSE;
421 }
422 } else if(ntracklet>0) isZeroBin = kFALSE; // if the event is not from Vz we chek the n of tracklets
423 }
424 }
425 return isZeroBin;
76532b17 426}
427
428void AlidNdEtaTask::UserExec(Option_t*)
0f67a57c 429{
430 // process the event
431
76532b17 432 fESD = (AliESDEvent*) fInputEvent;
433
1c15d51a 434 // these variables are also used in the MC section below; however, if ESD is off they stay with their default values
435 Bool_t eventTriggered = kFALSE;
436 const AliESDVertex* vtxESD = 0;
437
51f6de65 438 // post the data already here
76532b17 439 PostData(1, fOutput);
51f6de65 440
1c15d51a 441 // ESD analysis
442 if (fESD)
0f67a57c 443 {
69b09e3b 444 AliESDHeader* esdHeader = fESD->GetHeader();
445 if (!esdHeader)
446 {
447 Printf("ERROR: esdHeader could not be retrieved");
448 return;
449 }
ff8c4f30 450
81be4ee8 451 AliInputEventHandler* inputHandler = (AliInputEventHandler*) AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
452 if (!inputHandler)
1d107532 453 {
81be4ee8 454 Printf("ERROR: Could not receive input handler");
455 return;
12bb57f1 456 }
457
1103419a 458 // TODO use flags here!
81be4ee8 459 eventTriggered = inputHandler->IsEventSelected();
460
76532b17 461 static AliTriggerAnalysis* triggerAnalysis = 0;
462 if (!triggerAnalysis)
81be4ee8 463 {
464 AliPhysicsSelection* physicsSelection = dynamic_cast<AliPhysicsSelection*> (inputHandler->GetEventSelection());
465 if (physicsSelection)
76532b17 466 triggerAnalysis = physicsSelection->GetTriggerAnalysis();
81be4ee8 467 }
468
469 if (eventTriggered)
76532b17 470 eventTriggered = triggerAnalysis->IsTriggerFired(fESD, fTrigger);
12bb57f1 471
76532b17 472 AliTriggerAnalysis::V0Decision v0A = triggerAnalysis->V0Trigger(fESD, AliTriggerAnalysis::kASide, kFALSE);
473 AliTriggerAnalysis::V0Decision v0C = triggerAnalysis->V0Trigger(fESD, AliTriggerAnalysis::kCSide, kFALSE);
474 Bool_t fastORHW = (triggerAnalysis->SPDFiredChips(fESD, 1) > 0);
12bb57f1 475
476 Int_t vZero = 0;
477 if (v0A != AliTriggerAnalysis::kV0Invalid && v0C != AliTriggerAnalysis::kV0Invalid)
478 {
479 if (v0A == AliTriggerAnalysis::kV0Empty && v0C == AliTriggerAnalysis::kV0Empty) vZero = 1;
480 if (v0A == AliTriggerAnalysis::kV0BB && v0C == AliTriggerAnalysis::kV0Empty) vZero = 2;
481 if (v0A == AliTriggerAnalysis::kV0Empty && v0C == AliTriggerAnalysis::kV0BB) vZero = 3;
482 if (v0A == AliTriggerAnalysis::kV0BB && v0C == AliTriggerAnalysis::kV0BB) vZero = 4;
483 if (v0A == AliTriggerAnalysis::kV0BG && v0C == AliTriggerAnalysis::kV0Empty) vZero = 5;
484 if (v0A == AliTriggerAnalysis::kV0Empty && v0C == AliTriggerAnalysis::kV0BG) vZero = 6;
485 if (v0A == AliTriggerAnalysis::kV0BG && v0C == AliTriggerAnalysis::kV0BG) vZero = 7;
486 if (v0A == AliTriggerAnalysis::kV0BB && v0C == AliTriggerAnalysis::kV0BG) vZero = 8;
487 if (v0A == AliTriggerAnalysis::kV0BG && v0C == AliTriggerAnalysis::kV0BB) vZero = 9;
488 }
81be4ee8 489
490 if (vZero == 0)
491 Printf("VZERO: %d %d %d %d", vZero, eventTriggered, v0A, v0C);
1d107532 492
12bb57f1 493 Bool_t filled = kFALSE;
494
12bb57f1 495 if (!eventTriggered)
496 {
497 fStats2->Fill(0.0, vZero);
76532b17 498 fStats2->Fill(0.0, (fastORHW) ? 10 : 11);
12bb57f1 499 filled = kTRUE;
500 }
501
12bb57f1 502 if (eventTriggered)
12bb57f1 503 fStats->Fill(3);
12bb57f1 504
81be4ee8 505 /*
506 // ITS cluster tree
507 AliESDInputHandlerRP* handlerRP = dynamic_cast<AliESDInputHandlerRP*> (AliAnalysisManager::GetAnalys
508isManager()->GetInputEventHandler());
509 if (handlerRP)
12bb57f1 510 {
81be4ee8 511 TTree* itsClusterTree = handlerRP->GetTreeR("ITS");
512 if (!itsClusterTree)
513 return;
514
515 TClonesArray* itsClusters = new TClonesArray("AliITSRecPoint");
516 TBranch* itsClusterBranch=itsClusterTree->GetBranch("ITSRecPoints");
517
518 itsClusterBranch->SetAddress(&itsClusters);
519
520 Int_t nItsSubs = (Int_t)itsClusterTree->GetEntries();
521
522 Int_t totalClusters = 0;
523
524 // loop over the its subdetectors
525 for (Int_t iIts=0; iIts < nItsSubs; iIts++) {
526
527 if (!itsClusterTree->GetEvent(iIts))
528 continue;
529
530 Int_t nClusters = itsClusters->GetEntriesFast();
531 totalClusters += nClusters;
532
533 #ifdef FULLALIROOT
534 if (fAnalysisMode & AliPWG0Helper::kSPD)
12bb57f1 535 {
81be4ee8 536 // loop over clusters
537 while (nClusters--) {
538 AliITSRecPoint* cluster = (AliITSRecPoint*) itsClusters->UncheckedAt(nClusters);
539
540 Int_t layer = cluster->GetLayer();
541
542 if (layer > 1)
543 continue;
544
545 Float_t xyz[3] = {0., 0., 0.};
546 cluster->GetGlobalXYZ(xyz);
547
548 Float_t phi = TMath::Pi() + TMath::ATan2(-xyz[1], -xyz[0]);
549 Float_t z = xyz[2];
550
551 fZPhi[layer]->Fill(z, phi);
552 fModuleMap->Fill(layer * 80 + cluster->GetDetectorIndex());
12bb57f1 553 }
1d107532 554 }
81be4ee8 555 #endif
1d107532 556 }
12bb57f1 557 }
81be4ee8 558 */
559
12bb57f1 560 vtxESD = AliPWG0Helper::GetVertex(fESD, fAnalysisMode);
561 if (!vtxESD)
562 {
563 if (!filled)
564 {
81be4ee8 565 fStats2->Fill(1, vZero);
76532b17 566 fStats2->Fill(1, (fastORHW) ? 10 : 11);
12bb57f1 567 filled = kTRUE;
568 }
569 }
570 else
571 {
81be4ee8 572 if (!AliPWG0Helper::TestVertex(vtxESD, fAnalysisMode))
573 {
574 if (!filled)
575 {
576 fStats2->Fill(2, vZero);
76532b17 577 fStats2->Fill(2, (fastORHW) ? 10 : 11);
81be4ee8 578 filled = kTRUE;
579 }
580 }
581
12bb57f1 582 Double_t vtx[3];
583 vtxESD->GetXYZ(vtx);
1d107532 584
81be4ee8 585 // try to compare spd vertex and vertexer from tracks
586 // remove vertices outside +- 15 cm
587 if (TMath::Abs(vtx[2]) > 15)
12bb57f1 588 {
589 if (!filled)
590 {
591 fStats2->Fill(3, vZero);
76532b17 592 fStats2->Fill(3, (fastORHW) ? 10 : 11);
12bb57f1 593 filled = kTRUE;
594 }
595 }
81be4ee8 596
597 if (TMath::Abs(vtx[2]) > 10)
598 {
599 if (!filled)
600 {
601 fStats2->Fill(4, vZero);
76532b17 602 fStats2->Fill(4, (fastORHW) ? 10 : 11);
81be4ee8 603 filled = kTRUE;
604 }
605 }
1d107532 606
12bb57f1 607 const AliMultiplicity* mult = fESD->GetMultiplicity();
4ab23a37 608 if (!mult)
1103419a 609 {
68fa248f 610 Printf("Returning, no Multiplicity found");
12bb57f1 611 return;
1103419a 612 }
1d107532 613
12bb57f1 614 if (mult->GetNumberOfTracklets() == 0)
1d107532 615 {
12bb57f1 616 if (!filled)
1d107532 617 {
81be4ee8 618 fStats2->Fill(5, vZero);
76532b17 619 fStats2->Fill(5, (fastORHW) ? 10 : 11);
12bb57f1 620 filled = kTRUE;
1d107532 621 }
622 }
12bb57f1 623 }
624
12bb57f1 625 if (!filled)
81be4ee8 626 {
12bb57f1 627 fStats2->Fill(6, vZero);
76532b17 628 fStats2->Fill(6, (fastORHW) ? 10 : 11);
81be4ee8 629 //Printf("File: %s, IEV: %d, TRG: ---, Orbit: 0x%x, Period: %d, BC: %d", ((TTree*) GetInputData(0))->GetCurrentFile()->GetName(), fESD->GetEventNumberInFile(), fESD->GetOrbitNumber(),fESD->GetPeriodNumber(),fESD->GetBunchCrossNumber());
630 }
1d107532 631
a7f69e56 632 if (eventTriggered)
633 fStats->Fill(3);
12bb57f1 634
635 fStats->Fill(5);
636
1c15d51a 637 // get the ESD vertex
638 vtxESD = AliPWG0Helper::GetVertex(fESD, fAnalysisMode);
ff8c4f30 639
4ab23a37 640 Double_t vtx[3];
3d7758c1 641
1c15d51a 642 // fill z vertex resolution
643 if (vtxESD)
51f6de65 644 {
81be4ee8 645 if (strcmp(vtxESD->GetTitle(), "vertexer: Z") == 0)
646 fVertexResolution->Fill(vtxESD->GetZRes(), vtxESD->GetDispersion());
647 else
648 fVertexResolution->Fill(vtxESD->GetZRes(), 0);
649
1d107532 650 //if (strcmp(vtxESD->GetTitle(), "vertexer: 3D") == 0)
651 {
652 fVertex->Fill(vtxESD->GetXv(), vtxESD->GetYv(), vtxESD->GetZv());
653 }
654
51f6de65 655 if (AliPWG0Helper::TestVertex(vtxESD, fAnalysisMode))
656 {
657 vtxESD->GetXYZ(vtx);
69b09e3b 658
659 // vertex stats
660 if (strcmp(vtxESD->GetTitle(), "vertexer: 3D") == 0)
69b09e3b 661 fStats->Fill(1);
81be4ee8 662
663 if (strcmp(vtxESD->GetTitle(), "vertexer: Z") == 0)
69b09e3b 664 fStats->Fill(2);
81be4ee8 665
666 // remove vertices outside +- 15 cm
667 if (TMath::Abs(vtx[2]) > 15)
668 vtxESD = 0;
51f6de65 669 }
670 else
671 vtxESD = 0;
672 }
3d7758c1 673
1c15d51a 674 // needed for syst. studies
675 AliStack* stack = 0;
676 TArrayF vtxMC(3);
677
c17301f3 678 if (fUseMCVertex || fUseMCKine || fOnlyPrimaries || fReadMC) {
1c15d51a 679 if (!fReadMC) {
c17301f3 680 Printf("ERROR: fUseMCVertex or fUseMCKine or fOnlyPrimaries set without fReadMC set!");
1c15d51a 681 return;
682 }
683
684 AliMCEventHandler* eventHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
685 if (!eventHandler) {
686 Printf("ERROR: Could not retrieve MC event handler");
687 return;
688 }
689
690 AliMCEvent* mcEvent = eventHandler->MCEvent();
691 if (!mcEvent) {
692 Printf("ERROR: Could not retrieve MC event");
693 return;
694 }
695
696 AliHeader* header = mcEvent->Header();
697 if (!header)
698 {
699 AliDebug(AliLog::kError, "Header not available");
700 return;
701 }
3d7758c1 702
3d7758c1 703 // get the MC vertex
704 AliGenEventHeader* genHeader = header->GenEventHeader();
1c15d51a 705 if (!genHeader)
706 {
707 AliDebug(AliLog::kError, "Could not retrieve genHeader from Header");
708 return;
709 }
3d7758c1 710 genHeader->PrimaryVertex(vtxMC);
711
1c15d51a 712 if (fUseMCVertex)
1c15d51a 713 vtx[2] = vtxMC[2];
3d7758c1 714
c2b08c81 715 stack = mcEvent->Stack();
716 if (!stack)
3d7758c1 717 {
c2b08c81 718 AliDebug(AliLog::kError, "Stack not available");
719 return;
3d7758c1 720 }
721 }
81be4ee8 722
1c15d51a 723 // create list of (label, eta, pt) tuples
724 Int_t inputCount = 0;
725 Int_t* labelArr = 0;
726 Float_t* etaArr = 0;
69b09e3b 727 Float_t* thirdDimArr = 0;
a7f69e56 728 if (fAnalysisMode & AliPWG0Helper::kSPD)
0f67a57c 729 {
81be4ee8 730 if (vtxESD)
1c15d51a 731 {
81be4ee8 732 // get tracklets
733 const AliMultiplicity* mult = fESD->GetMultiplicity();
734 if (!mult)
1d107532 735 {
81be4ee8 736 AliDebug(AliLog::kError, "AliMultiplicity not available");
737 return;
738 }
739
740 Int_t arrayLength = mult->GetNumberOfTracklets();
741 if (fAnalysisMode & AliPWG0Helper::kSPDOnlyL0)
742 arrayLength += mult->GetNumberOfSingleClusters();
743
744 labelArr = new Int_t[arrayLength];
745 etaArr = new Float_t[arrayLength];
746 thirdDimArr = new Float_t[arrayLength];
747
748 // get multiplicity from SPD tracklets
749 for (Int_t i=0; i<mult->GetNumberOfTracklets(); ++i)
750 {
751 //printf("%d %f %f %f\n", i, mult->GetTheta(i), mult->GetPhi(i), mult->GetDeltaPhi(i));
752
753 if (fOnlyPrimaries)
754 if (mult->GetLabel(i, 0) < 0 || mult->GetLabel(i, 0) != mult->GetLabel(i, 1) || !stack->IsPhysicalPrimary(mult->GetLabel(i, 0)))
755 continue;
756
757 Float_t deltaPhi = mult->GetDeltaPhi(i);
758 // prevent values to be shifted by 2 Pi()
759 if (deltaPhi < -TMath::Pi())
760 deltaPhi += TMath::Pi() * 2;
761 if (deltaPhi > TMath::Pi())
762 deltaPhi -= TMath::Pi() * 2;
763
764 if (TMath::Abs(deltaPhi) > 1)
765 printf("WARNING: Very high Delta Phi: %d %f %f %f\n", i, mult->GetTheta(i), mult->GetPhi(i), deltaPhi);
766
767 Int_t label = mult->GetLabel(i, 0);
768 Float_t eta = mult->GetEta(i);
769
770 // control histograms
771 Float_t phi = mult->GetPhi(i);
772 if (phi < 0)
773 phi += TMath::Pi() * 2;
774
775 // TEST exclude probably inefficient phi region
776 //if (phi > 5.70 || phi < 0.06)
777 // continue;
778
779 fPhi->Fill(phi);
780
781 if (vtxESD && TMath::Abs(vtx[2]) < 10)
782 {
783 fEtaPhi->Fill(eta, phi);
784 fDeltaPhi->Fill(deltaPhi);
785 fDeltaTheta->Fill(mult->GetDeltaTheta(i));
786 }
787
788 if (fDeltaPhiCut > 0 && (TMath::Abs(deltaPhi) > fDeltaPhiCut || TMath::Abs(mult->GetDeltaTheta(i)) > fDeltaPhiCut / 0.08 * 0.025))
789 continue;
790
791 if (fUseMCKine)
792 {
793 if (label > 0)
794 {
795 TParticle* particle = stack->Particle(label);
796 eta = particle->Eta();
797 phi = particle->Phi();
798 }
799 else
800 Printf("WARNING: fUseMCKine set without fOnlyPrimaries and no label found");
801 }
802
803 if (fSymmetrize)
804 eta = TMath::Abs(eta);
805
806 etaArr[inputCount] = eta;
807 labelArr[inputCount] = label;
808 thirdDimArr[inputCount] = phi;
809 ++inputCount;
1d107532 810 }
811
81be4ee8 812 if (fAnalysisMode & AliPWG0Helper::kSPDOnlyL0)
51f6de65 813 {
81be4ee8 814 // get additional clusters from L0
815 for (Int_t i=0; i<mult->GetNumberOfSingleClusters(); ++i)
51f6de65 816 {
81be4ee8 817 etaArr[inputCount] = -TMath::Log(TMath::Tan(mult->GetThetaSingle(i)/2.));
818 labelArr[inputCount] = -1;
819 thirdDimArr[inputCount] = mult->GetPhiSingle(i);
820
821 ++inputCount;
51f6de65 822 }
51f6de65 823 }
69b09e3b 824
81be4ee8 825 if (!fFillPhi)
826 {
827 // fill multiplicity in third axis
828 for (Int_t i=0; i<inputCount; ++i)
829 thirdDimArr[i] = inputCount;
830 }
831
832 Int_t firedChips = mult->GetNumberOfFiredChips(0) + mult->GetNumberOfFiredChips(1);
833 fFiredChips->Fill(firedChips, inputCount);
834 Printf("Accepted %d tracklets (%d fired chips)", inputCount, firedChips);
835
836 fTrackletsVsUnassigned->Fill(inputCount, mult->GetNumberOfSingleClusters());
69b09e3b 837 }
0f67a57c 838 }
a7f69e56 839 else if (fAnalysisMode & AliPWG0Helper::kTPC || fAnalysisMode & AliPWG0Helper::kTPCITS)
0f67a57c 840 {
1c15d51a 841 if (!fEsdTrackCuts)
842 {
843 AliDebug(AliLog::kError, "fESDTrackCuts not available");
844 return;
845 }
0f67a57c 846
81be4ee8 847 Bool_t foundInEta10 = kFALSE;
848
69b09e3b 849 if (vtxESD)
0f67a57c 850 {
69b09e3b 851 // get multiplicity from ESD tracks
a7f69e56 852 TObjArray* list = fEsdTrackCuts->GetAcceptedTracks(fESD, fAnalysisMode & AliPWG0Helper::kTPC);
69b09e3b 853 Int_t nGoodTracks = list->GetEntries();
a7f69e56 854 Printf("Accepted %d tracks out of %d total ESD tracks", nGoodTracks, fESD->GetNumberOfTracks());
69b09e3b 855
856 labelArr = new Int_t[nGoodTracks];
857 etaArr = new Float_t[nGoodTracks];
858 thirdDimArr = new Float_t[nGoodTracks];
859
860 // loop over esd tracks
861 for (Int_t i=0; i<nGoodTracks; i++)
1c15d51a 862 {
69b09e3b 863 AliESDtrack* esdTrack = dynamic_cast<AliESDtrack*> (list->At(i));
864 if (!esdTrack)
c17301f3 865 {
69b09e3b 866 AliDebug(AliLog::kError, Form("ERROR: Could not retrieve track %d.", i));
867 continue;
c17301f3 868 }
69b09e3b 869
870 Float_t phi = esdTrack->Phi();
871 if (phi < 0)
872 phi += TMath::Pi() * 2;
873
874 Float_t eta = esdTrack->Eta();
875 Int_t label = TMath::Abs(esdTrack->GetLabel());
876 Float_t pT = esdTrack->Pt();
a7f69e56 877
878 // force pT to fixed value without B field
879 if (!(fAnalysisMode & AliPWG0Helper::kFieldOn))
880 pT = 1;
69b09e3b 881
882 fPhi->Fill(phi);
883 fEtaPhi->Fill(eta, phi);
884 if (eventTriggered && vtxESD)
885 fRawPt->Fill(pT);
886
68fa248f 887 if (esdTrack->Pt() < fPtMin) // even if the pt cut is already in the esdTrackCuts....
81be4ee8 888 continue;
889
a7f69e56 890 if (fOnlyPrimaries)
891 {
68fa248f 892 if (label == 0)
893 continue;
894
895 if (stack->IsPhysicalPrimary(label) == kFALSE)
896 continue;
a7f69e56 897 }
68fa248f 898
899 // 2 types of INEL>0 trigger - choose one
900
901 // 1. HL
902 // INEL>0 trigger
903 // if (TMath::Abs(esdTrack->Eta()) < 1 && esdTrack->Pt() > 0.15)
904 // foundInEta10 = kTRUE;
905
906 // 2. MB working group
907 if (TMath::Abs(esdTrack->Eta()) < 0.8 && esdTrack->Pt() > fPtMin){ // this should be in the trigger selection as well, so one particle should always be found
908 foundInEta10 = kTRUE;
909 }
69b09e3b 910
911 if (fUseMCKine)
912 {
68fa248f 913 if (label > 0)
69b09e3b 914 {
915 TParticle* particle = stack->Particle(label);
916 eta = particle->Eta();
917 pT = particle->Pt();
68fa248f 918 // check when using INEL>0, MB Working Group definition
919 if (TMath::Abs(eta) >=0.8 || pT <= fPtMin){
920 AliDebug(2,Form("WARNING ************* USING KINE: eta = %f, pt = %f",eta,pT));
921 }
69b09e3b 922 }
923 else
924 Printf("WARNING: fUseMCKine set without fOnlyPrimaries and no label found");
68fa248f 925 }
69b09e3b 926
1d107532 927 if (fSymmetrize)
928 eta = TMath::Abs(eta);
69b09e3b 929 etaArr[inputCount] = eta;
930 labelArr[inputCount] = TMath::Abs(esdTrack->GetLabel());
931 thirdDimArr[inputCount] = pT;
932 ++inputCount;
c17301f3 933 }
69b09e3b 934
81be4ee8 935 //if (inputCount > 30)
936 // Printf("Event with %d accepted TPC tracks. Period number: %d Orbit number: %x Bunch crossing number: %d", inputCount, fESD->GetPeriodNumber(), fESD->GetOrbitNumber(), fESD->GetBunchCrossNumber());
12bb57f1 937
69b09e3b 938 // TODO restrict inputCount used as measure for the multiplicity to |eta| < 1
939
940 delete list;
1c15d51a 941 }
81be4ee8 942
68fa248f 943 if (!foundInEta10){
81be4ee8 944 eventTriggered = kFALSE;
68fa248f 945 fHistEvents->Fill(0.1);
946 AliDebug(3,"Event rejected");
947 }
948 else{
949 if (eventTriggered){
950 fHistEvents->Fill(1.1);
951 AliDebug(3,Form("Event Accepted, with inputcount = %d",inputCount));
952 }
953 else{
954 AliDebug(3,"Event has foundInEta10 but was not triggered");
955 }
956 }
1c15d51a 957 }
958 else
959 return;
7307d52c 960
1c15d51a 961 // Processing of ESD information (always)
962 if (eventTriggered)
dd367a14 963 {
964 // control hist
1c15d51a 965 fMult->Fill(inputCount);
966 fdNdEtaAnalysisESD->FillTriggeredEvent(inputCount);
0f67a57c 967
1c15d51a 968 if (vtxESD)
dd367a14 969 {
1c15d51a 970 // control hist
1d107532 971
972 if (strcmp(vtxESD->GetTitle(), "vertexer: 3D") == 0)
973 fVertexVsMult->Fill(vtxESD->GetXv(), vtxESD->GetYv(), inputCount);
974
81be4ee8 975 Int_t part05 = 0;
976 Int_t part10 = 0;
1c15d51a 977 for (Int_t i=0; i<inputCount; ++i)
dd367a14 978 {
1c15d51a 979 Float_t eta = etaArr[i];
69b09e3b 980 Float_t thirdDim = thirdDimArr[i];
df974e92 981 fEta->Fill(eta);
0f67a57c 982
69b09e3b 983 fdNdEtaAnalysisESD->FillTrack(vtx[2], eta, thirdDim);
1c15d51a 984
69b09e3b 985 if (TMath::Abs(vtx[2]) < 10)
1c15d51a 986 {
987 fPartEta[0]->Fill(eta);
988
989 if (vtx[2] < 0)
990 fPartEta[1]->Fill(eta);
991 else
992 fPartEta[2]->Fill(eta);
993 }
81be4ee8 994
995 if (TMath::Abs(eta) < 0.5)
996 part05++;
68fa248f 997 if (TMath::Abs(eta) < 1.0) // in the INEL>0, MB WG definition, this is in any case equivalent to <0.8, due to the EsdTrackCuts
81be4ee8 998 part10++;
dd367a14 999 }
81be4ee8 1000
1001 Int_t multAxis = inputCount;
1002 if (fMultAxisEta1)
1003 multAxis = part10;
1004
1005 fMultVtx->Fill(multAxis);
1006 //if (TMath::Abs(vtx[2]) < 10)
1007 // fMultVtx->Fill(part05);
0f67a57c 1008
1c15d51a 1009 // for event count per vertex
81be4ee8 1010 fdNdEtaAnalysisESD->FillEvent(vtx[2], multAxis);
0f67a57c 1011
1c15d51a 1012 // control hist
81be4ee8 1013 if (multAxis > 0)
12bb57f1 1014 fEvents->Fill(vtx[2]);
1c15d51a 1015
1016 if (fReadMC)
1017 {
1018 // from tracks is only done for triggered and vertex reconstructed events
1019 for (Int_t i=0; i<inputCount; ++i)
1020 {
1021 Int_t label = labelArr[i];
1022
1023 if (label < 0)
1024 continue;
1025
1026 //Printf("Getting particle of track %d", label);
1027 TParticle* particle = stack->Particle(label);
1028
1029 if (!particle)
1030 {
1031 AliDebug(AliLog::kError, Form("ERROR: Could not retrieve particle %d.", label));
1032 continue;
1033 }
1034
69b09e3b 1035 Float_t thirdDim = -1;
a7f69e56 1036 if (fAnalysisMode & AliPWG0Helper::kSPD)
69b09e3b 1037 {
1038 if (fFillPhi)
1039 {
1040 thirdDim = particle->Phi();
1041 }
1042 else
81be4ee8 1043 thirdDim = multAxis;
69b09e3b 1044 }
1045 else
1046 thirdDim = particle->Pt();
1047
1d107532 1048 Float_t eta = particle->Eta();
1049 if (fSymmetrize)
1050 eta = TMath::Abs(eta);
1051 fdNdEtaAnalysisTracks->FillTrack(vtxMC[2], eta, thirdDim);
1c15d51a 1052 } // end of track loop
1053
1054 // for event count per vertex
81be4ee8 1055 fdNdEtaAnalysisTracks->FillEvent(vtxMC[2], multAxis);
1c15d51a 1056 }
1057 }
dd367a14 1058 }
69b09e3b 1059 if (etaArr)
1060 delete[] etaArr;
1061 if (labelArr)
1062 delete[] labelArr;
1063 if (thirdDimArr)
1064 delete[] thirdDimArr;
0f67a57c 1065 }
1066
1067 if (fReadMC) // Processing of MC information (optional)
1068 {
1069 AliMCEventHandler* eventHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
1070 if (!eventHandler) {
1071 Printf("ERROR: Could not retrieve MC event handler");
1072 return;
1073 }
1074
1075 AliMCEvent* mcEvent = eventHandler->MCEvent();
1076 if (!mcEvent) {
1077 Printf("ERROR: Could not retrieve MC event");
1078 return;
1079 }
1080
1c15d51a 1081 AliStack* stack = mcEvent->Stack();
0f67a57c 1082 if (!stack)
1083 {
1084 AliDebug(AliLog::kError, "Stack not available");
1085 return;
1086 }
1087
1088 AliHeader* header = mcEvent->Header();
1089 if (!header)
1090 {
1091 AliDebug(AliLog::kError, "Header not available");
1092 return;
1093 }
1094
1095 // get the MC vertex
1096 AliGenEventHeader* genHeader = header->GenEventHeader();
567160d6 1097 if (!genHeader)
1098 {
1099 AliDebug(AliLog::kError, "Could not retrieve genHeader from Header");
1100 return;
1101 }
1102
0f67a57c 1103 TArrayF vtxMC(3);
1104 genHeader->PrimaryVertex(vtxMC);
81be4ee8 1105
6b62a9c7 1106 // get process type
81be4ee8 1107 Int_t processType = AliPWG0Helper::GetEventProcessType(fESD, header, stack, fDiffTreatment);
6b62a9c7 1108 AliDebug(AliLog::kDebug+1, Form("Found process type %d", processType));
0fc41645 1109
81be4ee8 1110 if (processType == AliPWG0Helper::kInvalidProcess)
6b62a9c7 1111 AliDebug(AliLog::kError, Form("Unknown process type %d.", processType));
0fc41645 1112
0f67a57c 1113 // loop over mc particles
1114 Int_t nPrim = stack->GetNprimary();
1115
1116 Int_t nAcceptedParticles = 0;
81be4ee8 1117 Bool_t oneParticleEvent = kFALSE;
0f67a57c 1118
68fa248f 1119 Int_t nMCparticlesinRange = 0; // number of true particles in my range of interest
1120
51f6de65 1121 // count particles first, then fill
0f67a57c 1122 for (Int_t iMc = 0; iMc < nPrim; ++iMc)
1123 {
1c15d51a 1124 if (!stack->IsPhysicalPrimary(iMc))
1125 continue;
1126
0f67a57c 1127 //Printf("Getting particle %d", iMc);
1128 TParticle* particle = stack->Particle(iMc);
1129
1130 if (!particle)
1131 continue;
1132
1c15d51a 1133 if (particle->GetPDG()->Charge() == 0)
0f67a57c 1134 continue;
1135
1c15d51a 1136 //AliDebug(AliLog::kDebug+1, Form("Accepted primary %d, unique ID: %d", iMc, particle->GetUniqueID()));
0f67a57c 1137 Float_t eta = particle->Eta();
0f67a57c 1138
68fa248f 1139 AliDebug(2,Form("particle %d: eta = %f, pt = %f",iMc,particle->Eta(),particle->Pt()));
1140
1141 // INEL>0: choose one
1142 // 1.HL
1143 // if (TMath::Abs(eta) < 1.0){
1144 // oneParticleEvent = kTRUE;
1145 // nMCparticlesinRange++;
1146 //}
1147 // 2.MB Working Group definition
1148 if (TMath::Abs(eta) < 0.8 && particle->Pt() > fPtMin){
1149 oneParticleEvent = kTRUE;
1150 nMCparticlesinRange++;
1151 }
81be4ee8 1152
1153 // make a rough eta cut (so that nAcceptedParticles is not too far off the true measured value (NB: these histograms are only gathered for comparison))
ea441adf 1154 if (TMath::Abs(eta) < 1.5) // && pt > 0.3)
2b26296f 1155 nAcceptedParticles++;
51f6de65 1156 }
1157
68fa248f 1158 if (TMath::Abs(vtxMC[2]) < 10.){
1159 // if MC vertex is inside 10
1160 if (vtxESD){
1161 fVtxEffNum->Fill(nMCparticlesinRange);
1162 }
1163 fVtxEffDen->Fill(nMCparticlesinRange);
1164 if (eventTriggered){
1165 if (vtxESD){
1166 fVtxTrigEffNum->Fill(nMCparticlesinRange);
1167 }
1168 fVtxTrigEffDen->Fill(nMCparticlesinRange);
1169 fTrigEffNum->Fill(nMCparticlesinRange);
1170 }
1171 fTrigEffDen->Fill(nMCparticlesinRange);
1172 }
1173
1174 if (oneParticleEvent){
1175 fHistEventsMC->Fill(1.1);
1176 }
1177 else{
1178 fHistEventsMC->Fill(0.1);
1179 }
1180
51f6de65 1181 for (Int_t iMc = 0; iMc < nPrim; ++iMc)
1182 {
1183 if (!stack->IsPhysicalPrimary(iMc))
1184 continue;
2b26296f 1185
51f6de65 1186 //Printf("Getting particle %d", iMc);
1187 TParticle* particle = stack->Particle(iMc);
1188
1189 if (!particle)
1190 continue;
1191
1192 if (particle->GetPDG()->Charge() == 0)
1193 continue;
1194
1195 Float_t eta = particle->Eta();
1d107532 1196 if (fSymmetrize)
1197 eta = TMath::Abs(eta);
1198
69b09e3b 1199 Float_t thirdDim = -1;
1200
a7f69e56 1201 if (fAnalysisMode & AliPWG0Helper::kSPD)
69b09e3b 1202 {
1203 if (fFillPhi)
1204 {
1205 thirdDim = particle->Phi();
1206 }
1207 else
1208 thirdDim = nAcceptedParticles;
1209 }
1210 else
1211 thirdDim = particle->Pt();
51f6de65 1212
1213 fdNdEtaAnalysis->FillTrack(vtxMC[2], eta, thirdDim);
0f67a57c 1214
69b09e3b 1215 if (processType != AliPWG0Helper::kSD)
51f6de65 1216 fdNdEtaAnalysisNSD->FillTrack(vtxMC[2], eta, thirdDim);
0fc41645 1217
69b09e3b 1218 if (processType == AliPWG0Helper::kND)
1219 fdNdEtaAnalysisND->FillTrack(vtxMC[2], eta, thirdDim);
81be4ee8 1220
df974e92 1221 if (oneParticleEvent){
1222 AliDebug(3,Form("filling dNdEtaAnalysis object:: vtx = %f, eta = %f, pt = %f",vtxMC[2], eta, thirdDim));
1223 fdNdEtaAnalysisOnePart->FillTrack(vtxMC[2], eta, thirdDim);
1224 }
69b09e3b 1225
0f67a57c 1226 if (eventTriggered)
1227 {
51f6de65 1228 fdNdEtaAnalysisTr->FillTrack(vtxMC[2], eta, thirdDim);
770a1f1d 1229 if (vtxESD)
51f6de65 1230 fdNdEtaAnalysisTrVtx->FillTrack(vtxMC[2], eta, thirdDim);
0f67a57c 1231 }
1232
ff8c4f30 1233 if (TMath::Abs(eta) < 1.0 && particle->Pt() > 0 && particle->P() > 0)
1234 {
1d107532 1235 //Float_t value = 1. / TMath::TwoPi() / particle->Pt() * particle->Energy() / particle->P();
1236 Float_t value = 1;
ff8c4f30 1237 fPartPt->Fill(particle->Pt(), value);
68fa248f 1238 if (TMath::Abs(eta) < 0.8 && particle->Pt() > fPtMin){
1239 fEtaMC->Fill(eta);
1240 }
ff8c4f30 1241 }
0f67a57c 1242 }
1243
1244 fdNdEtaAnalysis->FillEvent(vtxMC[2], nAcceptedParticles);
ea441adf 1245 if (processType != AliPWG0Helper::kSD)
0fc41645 1246 fdNdEtaAnalysisNSD->FillEvent(vtxMC[2], nAcceptedParticles);
69b09e3b 1247 if (processType == AliPWG0Helper::kND)
1248 fdNdEtaAnalysisND->FillEvent(vtxMC[2], nAcceptedParticles);
81be4ee8 1249 if (oneParticleEvent)
1250 fdNdEtaAnalysisOnePart->FillEvent(vtxMC[2], nAcceptedParticles);
0fc41645 1251
0f67a57c 1252 if (eventTriggered)
1253 {
1254 fdNdEtaAnalysisTr->FillEvent(vtxMC[2], nAcceptedParticles);
770a1f1d 1255 if (vtxESD)
0f67a57c 1256 fdNdEtaAnalysisTrVtx->FillEvent(vtxMC[2], nAcceptedParticles);
1257 }
0f67a57c 1258 }
0f67a57c 1259}
1260
1261void AlidNdEtaTask::Terminate(Option_t *)
1262{
1263 // The Terminate() function is the last function to be called during
1264 // a query. It always runs on the client, it can be used to present
1265 // the results graphically or save the results to file.
1266
76532b17 1267 fOutput = dynamic_cast<TList*> (GetOutputData(1));
1c15d51a 1268 if (!fOutput)
0f67a57c 1269 Printf("ERROR: fOutput not available");
0f67a57c 1270
68fa248f 1271 TH1D* dNdEta = new TH1D("dNdEta","#eta;counts",80, -4, 4);
1272 TH1D* dNdEtaMC = new TH1D("dNdEtaMC","#eta,MC;counts",80, -4, 4);
1c15d51a 1273 if (fOutput)
1274 {
1275 fdNdEtaAnalysisESD = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("fdNdEtaAnalysisESD"));
1276 fMult = dynamic_cast<TH1F*> (fOutput->FindObject("fMult"));
1277 fMultVtx = dynamic_cast<TH1F*> (fOutput->FindObject("fMultVtx"));
1278 for (Int_t i=0; i<3; ++i)
1279 fPartEta[i] = dynamic_cast<TH1F*> (fOutput->FindObject(Form("dndeta_check_%d", i)));
1280 fEvents = dynamic_cast<TH1F*> (fOutput->FindObject("dndeta_check_vertex"));
81be4ee8 1281 fVertexResolution = dynamic_cast<TH2F*> (fOutput->FindObject("dndeta_vertex_resolution_z"));
ea441adf 1282
1d107532 1283 fVertex = dynamic_cast<TH3F*> (fOutput->FindObject("vertex_check"));
1284 fVertexVsMult = dynamic_cast<TH3F*> (fOutput->FindObject("fVertexVsMult"));
1c15d51a 1285 fPhi = dynamic_cast<TH1F*> (fOutput->FindObject("fPhi"));
69b09e3b 1286 fRawPt = dynamic_cast<TH1F*> (fOutput->FindObject("fRawPt"));
1c15d51a 1287 fEtaPhi = dynamic_cast<TH2F*> (fOutput->FindObject("fEtaPhi"));
69b09e3b 1288 for (Int_t i=0; i<2; ++i)
1289 fZPhi[i] = dynamic_cast<TH2F*> (fOutput->FindObject(Form("fZPhi_%d", i)));
12bb57f1 1290 fModuleMap = dynamic_cast<TH1F*> (fOutput->FindObject("fModuleMap"));
1c15d51a 1291 fDeltaPhi = dynamic_cast<TH1F*> (fOutput->FindObject("fDeltaPhi"));
a7f69e56 1292 fDeltaTheta = dynamic_cast<TH1F*> (fOutput->FindObject("fDeltaTheta"));
69b09e3b 1293 fFiredChips = dynamic_cast<TH2F*> (fOutput->FindObject("fFiredChips"));
1d107532 1294 fTrackletsVsClusters = dynamic_cast<TH2F*> (fOutput->FindObject("fTrackletsVsClusters"));
1295 fTrackletsVsUnassigned = dynamic_cast<TH2F*> (fOutput->FindObject("fTrackletsVsUnassigned"));
69b09e3b 1296 fStats = dynamic_cast<TH1F*> (fOutput->FindObject("fStats"));
1d107532 1297 fStats2 = dynamic_cast<TH2F*> (fOutput->FindObject("fStats2"));
1c15d51a 1298
1299 fEsdTrackCuts = dynamic_cast<AliESDtrackCuts*> (fOutput->FindObject("fEsdTrackCuts"));
68fa248f 1300
1301 fEta = dynamic_cast<TH1D*>(fOutput->FindObject("fEta"));
1302 fEtaMC = dynamic_cast<TH1D*>(fOutput->FindObject("fEtaMC"));
1303 fHistEvents = dynamic_cast<TH1D*>(fOutput->FindObject("fHistEvents"));
1304 fHistEventsMC = dynamic_cast<TH1D*>(fOutput->FindObject("fHistEventsMC"));
1305 fVtxEffDen = dynamic_cast<TH1D*>(fOutput->FindObject("fVtxEffDen"));
1306 fVtxEffNum = dynamic_cast<TH1D*>(fOutput->FindObject("fVtxEffNum"));
1307 fTrigEffDen = dynamic_cast<TH1D*>(fOutput->FindObject("fTrigEffDen"));
1308 fTrigEffNum = dynamic_cast<TH1D*>(fOutput->FindObject("fTrigEffNum"));
1309 fVtxTrigEffDen = dynamic_cast<TH1D*>(fOutput->FindObject("fVtxTrigEffDen"));
1310 fVtxTrigEffNum = dynamic_cast<TH1D*>(fOutput->FindObject("fVtxTrigEffNum"));
1311 Printf("Events selected = %f", fHistEvents->GetBinContent(fHistEvents->GetXaxis()->FindBin(1.1)));
1312 Printf("Events selected frm MC = %f", fHistEventsMC->GetBinContent(fHistEventsMC->GetXaxis()->FindBin(1.1)));
1313 Float_t nevents = fHistEvents->GetBinContent(fHistEvents->GetXaxis()->FindBin(1.1));
1314 Float_t neventsMC = fHistEventsMC->GetBinContent(fHistEventsMC->GetXaxis()->FindBin(1.1));
1315 for (Int_t ibin = 1; ibin <= fEta->GetNbinsX(); ibin++){
1316 Float_t eta = fEta->GetBinContent(ibin)/nevents/fEta->GetBinWidth(ibin);
1317 Float_t etaerr = fEta->GetBinError(ibin)/nevents/fEta->GetBinWidth(ibin);
1318 Float_t etaMC = fEtaMC->GetBinContent(ibin)/neventsMC/fEtaMC->GetBinWidth(ibin);
1319 Float_t etaerrMC = fEtaMC->GetBinError(ibin)/neventsMC/fEtaMC->GetBinWidth(ibin);
1320 dNdEta->SetBinContent(ibin,eta);
1321 dNdEta->SetBinError(ibin,etaerr);
1322 dNdEtaMC->SetBinContent(ibin,etaMC);
1323 dNdEtaMC->SetBinError(ibin,etaerrMC);
1324 }
1325 new TCanvas("eta", " eta ",50, 50, 550, 550) ;
1326 fEta->Draw();
1327 new TCanvas("etaMC", " etaMC ",50, 50, 550, 550) ;
1328 fEtaMC->Draw();
1329 new TCanvas("dNdEta", "#eta;dNdEta ",50, 50, 550, 550) ;
1330 dNdEta->Draw();
1331 new TCanvas("dNdEtaMC", "#eta,MC;dNdEta ",50, 50, 550, 550) ;
1332 dNdEtaMC->Draw();
1333 new TCanvas("Events", "Events;Events ",50, 50, 550, 550) ;
1334 fHistEvents->Draw();
1335 new TCanvas("Events, MC", "Events, MC;Events ",50, 50, 550, 550) ;
1336 fHistEventsMC->Draw();
1337
1338 TFile* outputFileCheck = new TFile("histogramsCheck.root", "RECREATE");
1339 fEta->Write();
1340 fEtaMC->Write();
1341 dNdEta->Write();
1342 dNdEtaMC->Write();
1343 outputFileCheck->Write();
1344 outputFileCheck->Close();
1345
1c15d51a 1346 }
0f67a57c 1347
1348 if (!fdNdEtaAnalysisESD)
1349 {
1350 AliDebug(AliLog::kError, "ERROR: fdNdEtaAnalysisESD not available");
0f67a57c 1351 }
51f6de65 1352 else
0f67a57c 1353 {
51f6de65 1354 if (fMult && fMultVtx)
1355 {
1356 new TCanvas;
1357 fMult->Draw();
1358 fMultVtx->SetLineColor(2);
1359 fMultVtx->Draw("SAME");
1360 }
0f67a57c 1361
69b09e3b 1362 if (fFiredChips)
1363 {
1364 new TCanvas;
1365 fFiredChips->Draw("COLZ");
1366 }
1367
6714bf63 1368 if (fPartEta[0] && fPartEta[1] && fPartEta[2] && fEvents && fMultVtx && fMult)
51f6de65 1369 {
81be4ee8 1370 Int_t events1 = (Int_t) fEvents->Integral(fEvents->GetXaxis()->FindBin(-9.9), fEvents->GetXaxis()->FindBin(-0.001));
1371 Int_t events2 = (Int_t) fEvents->Integral(fEvents->GetXaxis()->FindBin(0.001), fEvents->GetXaxis()->FindBin(9.9));
0f67a57c 1372
51f6de65 1373 Printf("%d events with vertex used", events1 + events2);
81be4ee8 1374 Printf("%d total events with vertex, %d total triggered events, %d triggered events with 0 multiplicity", (Int_t) fMultVtx->Integral(), (Int_t) fMult->Integral(), (Int_t) fMult->GetBinContent(1));
54b096ef 1375
51f6de65 1376 if (events1 > 0 && events2 > 0)
1377 {
1378 fPartEta[0]->Scale(1.0 / (events1 + events2));
1379 fPartEta[1]->Scale(1.0 / events1);
1380 fPartEta[2]->Scale(1.0 / events2);
0f67a57c 1381
51f6de65 1382 for (Int_t i=0; i<3; ++i)
1383 fPartEta[i]->Scale(1.0 / fPartEta[i]->GetBinWidth(1));
0f67a57c 1384
51f6de65 1385 new TCanvas("control", "control", 500, 500);
1386 for (Int_t i=0; i<3; ++i)
1387 {
1388 fPartEta[i]->SetLineColor(i+1);
1389 fPartEta[i]->Draw((i==0) ? "" : "SAME");
1390 }
1391 }
0f67a57c 1392 }
51f6de65 1393
1394 if (fEvents)
1395 {
1396 new TCanvas("control3", "control3", 500, 500);
1397 fEvents->Draw();
54b096ef 1398 }
12bb57f1 1399
1400 if (fStats2)
1401 {
1402 new TCanvas;
1403 fStats2->Draw("TEXT");
1404 }
0f67a57c 1405
51f6de65 1406 TFile* fout = new TFile("analysis_esd_raw.root", "RECREATE");
0f67a57c 1407
51f6de65 1408 if (fdNdEtaAnalysisESD)
1409 fdNdEtaAnalysisESD->SaveHistograms();
0f67a57c 1410
51f6de65 1411 if (fEsdTrackCuts)
a7f69e56 1412 fEsdTrackCuts->SaveHistograms("esd_track_cuts");
0f67a57c 1413
51f6de65 1414 if (fMult)
1415 fMult->Write();
0f67a57c 1416
51f6de65 1417 if (fMultVtx)
1418 fMultVtx->Write();
0f67a57c 1419
51f6de65 1420 for (Int_t i=0; i<3; ++i)
1421 if (fPartEta[i])
1422 fPartEta[i]->Write();
0f67a57c 1423
51f6de65 1424 if (fEvents)
1425 fEvents->Write();
0f67a57c 1426
51f6de65 1427 if (fVertexResolution)
81be4ee8 1428 {
51f6de65 1429 fVertexResolution->Write();
81be4ee8 1430 fVertexResolution->ProjectionX();
1431 fVertexResolution->ProjectionY();
1432 }
0f67a57c 1433
51f6de65 1434 if (fDeltaPhi)
1435 fDeltaPhi->Write();
54b096ef 1436
a7f69e56 1437 if (fDeltaTheta)
1438 fDeltaTheta->Write();
1439
51f6de65 1440 if (fPhi)
1441 fPhi->Write();
54b096ef 1442
69b09e3b 1443 if (fRawPt)
1444 fRawPt->Write();
1445
51f6de65 1446 if (fEtaPhi)
1447 fEtaPhi->Write();
ea441adf 1448
69b09e3b 1449 for (Int_t i=0; i<2; ++i)
1450 if (fZPhi[i])
1451 fZPhi[i]->Write();
1452
12bb57f1 1453 if (fModuleMap)
1454 fModuleMap->Write();
1455
69b09e3b 1456 if (fFiredChips)
1457 fFiredChips->Write();
1458
1d107532 1459 if (fTrackletsVsClusters)
1460 fTrackletsVsClusters->Write();
1461
1462 if (fTrackletsVsUnassigned)
1463 fTrackletsVsUnassigned->Write();
1464
69b09e3b 1465 if (fStats)
1466 fStats->Write();
1467
1d107532 1468 if (fStats2)
1469 fStats2->Write();
1470
51f6de65 1471 if (fVertex)
1472 fVertex->Write();
ea441adf 1473
1d107532 1474 if (fVertexVsMult)
1475 fVertexVsMult->Write();
1476
68fa248f 1477 if (fVtxEffDen) fVtxEffDen->Write();
1478 else Printf("fVtxEffDen not found");
1479
1480 if (fVtxEffNum) fVtxEffNum->Write();
1481 else Printf("fVtxEffNum not found");
1482
1483 if (fVtxTrigEffNum) fVtxTrigEffNum->Write();
1484 else Printf("fVtxTrigEffNum not found");
1485
1486 if (fVtxTrigEffDen) fVtxTrigEffDen->Write();
1487 else Printf("fVtxTrigEffDen not found");
1488
1489 if (fTrigEffNum) fTrigEffNum->Write();
1490 else Printf("fTrigEffNum not found");
1491
1492 if (fTrigEffDen) fTrigEffDen->Write();
1493 else Printf("fTrigEffDen not found");
1494
51f6de65 1495 fout->Write();
1496 fout->Close();
0f67a57c 1497
51f6de65 1498 Printf("Writting result to analysis_esd_raw.root");
1499 }
0f67a57c 1500
1501 if (fReadMC)
1502 {
1c15d51a 1503 if (fOutput)
1504 {
1505 fdNdEtaAnalysis = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndeta"));
69b09e3b 1506 fdNdEtaAnalysisND = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndetaND"));
1c15d51a 1507 fdNdEtaAnalysisNSD = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndetaNSD"));
81be4ee8 1508 fdNdEtaAnalysisOnePart = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndetaOnePart"));
1c15d51a 1509 fdNdEtaAnalysisTr = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndetaTr"));
1510 fdNdEtaAnalysisTrVtx = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndetaTrVtx"));
1511 fdNdEtaAnalysisTracks = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndetaTracks"));
1512 fPartPt = dynamic_cast<TH1F*> (fOutput->FindObject("dndeta_check_pt"));
1c15d51a 1513 }
0f67a57c 1514
1515 if (!fdNdEtaAnalysis || !fdNdEtaAnalysisTr || !fdNdEtaAnalysisTrVtx || !fPartPt)
1516 {
1517 AliDebug(AliLog::kError, Form("ERROR: Histograms not available %p %p", (void*) fdNdEtaAnalysis, (void*) fPartPt));
1518 return;
1519 }
1520
1521 fdNdEtaAnalysis->Finish(0, -1, AlidNdEtaCorrection::kNone);
69b09e3b 1522 fdNdEtaAnalysisND->Finish(0, -1, AlidNdEtaCorrection::kNone);
0fc41645 1523 fdNdEtaAnalysisNSD->Finish(0, -1, AlidNdEtaCorrection::kNone);
81be4ee8 1524 fdNdEtaAnalysisOnePart->Finish(0, -1, AlidNdEtaCorrection::kNone);
0f67a57c 1525 fdNdEtaAnalysisTr->Finish(0, -1, AlidNdEtaCorrection::kNone);
1526 fdNdEtaAnalysisTrVtx->Finish(0, -1, AlidNdEtaCorrection::kNone);
1527 fdNdEtaAnalysisTracks->Finish(0, -1, AlidNdEtaCorrection::kNone);
1528
1529 Int_t events = (Int_t) fdNdEtaAnalysis->GetData()->GetEventCorrection()->GetMeasuredHistogram()->Integral();
1530 fPartPt->Scale(1.0/events);
1531 fPartPt->Scale(1.0/fPartPt->GetBinWidth(1));
1532
51f6de65 1533 TFile* fout = new TFile("analysis_mc.root","RECREATE");
0f67a57c 1534
1535 fdNdEtaAnalysis->SaveHistograms();
69b09e3b 1536 fdNdEtaAnalysisND->SaveHistograms();
0fc41645 1537 fdNdEtaAnalysisNSD->SaveHistograms();
81be4ee8 1538 fdNdEtaAnalysisOnePart->SaveHistograms();
0f67a57c 1539 fdNdEtaAnalysisTr->SaveHistograms();
1540 fdNdEtaAnalysisTrVtx->SaveHistograms();
1541 fdNdEtaAnalysisTracks->SaveHistograms();
ea441adf 1542
1543 if (fPartPt)
1544 fPartPt->Write();
1545
0f67a57c 1546 fout->Write();
1547 fout->Close();
1548
1549 Printf("Writting result to analysis_mc.root");
1550
1551 if (fPartPt)
1552 {
1553 new TCanvas("control2", "control2", 500, 500);
770a1f1d 1554 fPartPt->Draw();
0f67a57c 1555 }
1556 }
1557}