Moving PbPb multiplicity in the new directory structure
[u/mrichter/AliRoot.git] / PWGLF / SPECTRA / ChargedHadrons / multPbPb / AliAnalysisTaskMultPbTracks.cxx
CommitLineData
a23f7c97 1// AliAnalysisTaskMultPbTracks
2
3// Author: Michele Floris, CERN
4// TODO:
5// - Add chi2/cluster plot for primary, secondaries and fakes
a23f7c97 6
7#include "AliAnalysisTaskMultPbTracks.h"
8#include "AliESDInputHandler.h"
9#include "AliAnalysisMultPbTrackHistoManager.h"
10#include "AliAnalysisManager.h"
11#include "AliESDtrackCuts.h"
a23f7c97 12#include "AliMCEvent.h"
13#include "AliStack.h"
14#include "TH1I.h"
15#include "TH3D.h"
bcc49144 16#include "TH2D.h"
a23f7c97 17#include "AliMCParticle.h"
18#include "AliGenEventHeader.h"
9ef42f6c 19#include "AliCentrality.h"
7f5f2e0c 20#include "AliMultiplicity.h"
abd808b9 21#include "TFile.h"
a23f7c97 22#include <iostream>
2bbfd8c6 23#include "AliAnalysisMultPbCentralitySelector.h"
bcf2601a 24#include "AliTriggerAnalysis.h"
2b42cee2 25#include "AliPIDResponse.h"
a23f7c97 26
27using namespace std;
28
29ClassImp(AliAnalysisTaskMultPbTracks)
30
31AliAnalysisTaskMultPbTracks::AliAnalysisTaskMultPbTracks()
54b31d6a 32: AliAnalysisTaskSE("TaskMultPbTracks"),
2b42cee2 33 fESD(0),fHistoManager(0),fCentrSelector(0),fTrackCuts(0),fTrackCutsNoDCA(0),fOfflineTrigger(0), fIsMC(0),fIsTPCOnly(0), fTriggerAnalysis(0),fPIDResponse(0)
a23f7c97 34{
35 // constructor
36
37 DefineOutput(1, AliAnalysisMultPbTrackHistoManager::Class());
54b31d6a 38 DefineOutput(2, AliESDtrackCuts::Class());
2bbfd8c6 39 DefineOutput(3, AliAnalysisMultPbCentralitySelector::Class());
a23f7c97 40
41 fHistoManager = new AliAnalysisMultPbTrackHistoManager("histoManager","Hitograms, Multiplicity, Track analysis");
42 if(fIsMC) fHistoManager->SetSuffix("MC");
43
44}
45AliAnalysisTaskMultPbTracks::AliAnalysisTaskMultPbTracks(const char * name)
54b31d6a 46 : AliAnalysisTaskSE(name),
2b42cee2 47 fESD(0),fHistoManager(0),fCentrSelector(0),fTrackCuts(0),fTrackCutsNoDCA(0),fOfflineTrigger(0),fIsMC(0),fIsTPCOnly(0), fTriggerAnalysis(0),fPIDResponse(0)
a23f7c97 48{
49 //
50 // Standard constructur which should be used
51 //
52
53 DefineOutput(1, AliAnalysisMultPbTrackHistoManager::Class());
54b31d6a 54 DefineOutput(2, AliESDtrackCuts::Class());
2bbfd8c6 55 DefineOutput(3, AliAnalysisMultPbCentralitySelector::Class());
a23f7c97 56
57 fHistoManager = new AliAnalysisMultPbTrackHistoManager("histoManager","Hitograms, Multiplicity, Track analysis");
58 if(fIsMC) fHistoManager->SetSuffix("MC");
59
60}
61
62AliAnalysisTaskMultPbTracks::AliAnalysisTaskMultPbTracks(const AliAnalysisTaskMultPbTracks& obj) :
2b42cee2 63 AliAnalysisTaskSE(obj) ,fESD (0), fHistoManager(0), fCentrSelector(0), fTrackCuts(0),fTrackCutsNoDCA(0),fOfflineTrigger(0),fIsMC(0),fIsTPCOnly(0), fTriggerAnalysis(0),fPIDResponse(0)
a23f7c97 64{
65 //copy ctor
66 fESD = obj.fESD ;
67 fHistoManager= obj.fHistoManager;
a23f7c97 68 fTrackCuts = obj.fTrackCuts;
69 fTrackCutsNoDCA = obj.fTrackCutsNoDCA;
2bbfd8c6 70 fCentrSelector = obj.fCentrSelector;
e9c445f5 71 fOfflineTrigger = obj.fOfflineTrigger;
9441cdaa 72 fIsMC = obj.fIsMC;
73 fIsTPCOnly = obj.fIsTPCOnly;
bcf2601a 74 fTriggerAnalysis = obj.fTriggerAnalysis;
a23f7c97 75
76}
77
78AliAnalysisTaskMultPbTracks::~AliAnalysisTaskMultPbTracks(){
79 // destructor
80
81 if(!AliAnalysisManager::GetAnalysisManager()->IsProofMode()) {
82 if(fHistoManager) {
83 delete fHistoManager;
84 fHistoManager = 0;
bcf2601a 85 delete fTriggerAnalysis;
86 fTriggerAnalysis=0;
a23f7c97 87 }
88 }
89 // Histo list should not be destroyed: fListWrapper is owner!
90
91}
92void AliAnalysisTaskMultPbTracks::UserCreateOutputObjects()
93{
94 // Called once
95
96 // For the DCA distributions, we want to use exactly the same cuts
97 // as for the other distributions, with the exception of the DCA cut
98 fTrackCutsNoDCA = new AliESDtrackCuts(*fTrackCuts); // clone cuts
99 // Reset all DCA cuts; FIXME: is this all?
100 fTrackCutsNoDCA->SetMaxDCAToVertexXY();
101 fTrackCutsNoDCA->SetMaxDCAToVertexZ ();
102 fTrackCutsNoDCA->SetMaxDCAToVertexXYPtDep();
103 fTrackCutsNoDCA->SetMaxDCAToVertexZPtDep();
104
bcf2601a 105 fTriggerAnalysis = new AliTriggerAnalysis();
bfae2924 106 fTriggerAnalysis->SetAnalyzeMC(fIsMC);
2b42cee2 107
108 //The common PID object can then be retrieved from the input handler. This can naturally be done in the UserCreateOutputObjects:
109 AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
110 AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
111 fPIDResponse = inputHandler->GetPIDResponse();
112
a23f7c97 113}
114
115
116void AliAnalysisTaskMultPbTracks::UserExec(Option_t *)
117{
118 // User code
119
e9c445f5 120
a23f7c97 121 /* PostData(0) is taken care of by AliAnalysisTaskSE */
122 PostData(1,fHistoManager);
54b31d6a 123 PostData(2,fTrackCuts);
2bbfd8c6 124 PostData(3,fCentrSelector);
a23f7c97 125
abd808b9 126 // Cache (some) histogram pointers to avoid loop in the histo manager
e0376287 127 static TH3D * hTracks [AliAnalysisMultPbTrackHistoManager::kNHistos];
bcc49144 128 static TH2D * hDCA [AliAnalysisMultPbTrackHistoManager::kNHistos];
e0376287 129 static TH1D * hNTracks [AliAnalysisMultPbTrackHistoManager::kNHistos];
a23f7c97 130 hTracks[AliAnalysisMultPbTrackHistoManager::kHistoGen] = fHistoManager->GetHistoPtEtaVz(AliAnalysisMultPbTrackHistoManager::kHistoGen );
131 hTracks[AliAnalysisMultPbTrackHistoManager::kHistoRec] = fHistoManager->GetHistoPtEtaVz(AliAnalysisMultPbTrackHistoManager::kHistoRec );
132 hTracks[AliAnalysisMultPbTrackHistoManager::kHistoRecPrim] = fHistoManager->GetHistoPtEtaVz(AliAnalysisMultPbTrackHistoManager::kHistoRecPrim );
133 hTracks[AliAnalysisMultPbTrackHistoManager::kHistoRecFake] = fHistoManager->GetHistoPtEtaVz(AliAnalysisMultPbTrackHistoManager::kHistoRecFake );
134 hTracks[AliAnalysisMultPbTrackHistoManager::kHistoRecSecMat] = fHistoManager->GetHistoPtEtaVz(AliAnalysisMultPbTrackHistoManager::kHistoRecSecMat );
135 hTracks[AliAnalysisMultPbTrackHistoManager::kHistoRecSecWeak] = fHistoManager->GetHistoPtEtaVz(AliAnalysisMultPbTrackHistoManager::kHistoRecSecWeak);
136
137 hDCA[AliAnalysisMultPbTrackHistoManager::kHistoGen] = fHistoManager->GetHistoDCA(AliAnalysisMultPbTrackHistoManager::kHistoGen );
138 hDCA[AliAnalysisMultPbTrackHistoManager::kHistoRec] = fHistoManager->GetHistoDCA(AliAnalysisMultPbTrackHistoManager::kHistoRec );
139 hDCA[AliAnalysisMultPbTrackHistoManager::kHistoRecPrim] = fHistoManager->GetHistoDCA(AliAnalysisMultPbTrackHistoManager::kHistoRecPrim );
140 hDCA[AliAnalysisMultPbTrackHistoManager::kHistoRecFake] = fHistoManager->GetHistoDCA(AliAnalysisMultPbTrackHistoManager::kHistoRecFake );
141 hDCA[AliAnalysisMultPbTrackHistoManager::kHistoRecSecMat] = fHistoManager->GetHistoDCA(AliAnalysisMultPbTrackHistoManager::kHistoRecSecMat );
142 hDCA[AliAnalysisMultPbTrackHistoManager::kHistoRecSecWeak] = fHistoManager->GetHistoDCA(AliAnalysisMultPbTrackHistoManager::kHistoRecSecWeak);
143
2bbfd8c6 144 hNTracks[AliAnalysisMultPbTrackHistoManager::kHistoGen] = fHistoManager->GetHistoMult(AliAnalysisMultPbTrackHistoManager::kHistoGen );
145 hNTracks[AliAnalysisMultPbTrackHistoManager::kHistoRec] = fHistoManager->GetHistoMult(AliAnalysisMultPbTrackHistoManager::kHistoRec );
a23f7c97 146
abd808b9 147
148 fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRec)->Reset();
149
a23f7c97 150 fESD = dynamic_cast<AliESDEvent*>(fInputEvent);
151 if (strcmp(fESD->ClassName(),"AliESDEvent")) {
152 AliFatal("Not processing ESDs");
153 }
154
a23f7c97 155 fHistoManager->GetHistoStats()->Fill(AliAnalysisMultPbTrackHistoManager::kStatAll);
156
8b628a20 157
54b31d6a 158 // Centrality selection
2bbfd8c6 159 Bool_t isCentralitySelected = fCentrSelector->IsCentralityBinSelected(fESD,fTrackCuts);
160 if(!isCentralitySelected) return;
161
162 // AliESDCentrality *centrality = fESD->GetCentrality();
163 // if(!centrality) {
164 // AliFatal("Centrality object not available");
165 // }
166 // else {
167 // Int_t centrBin = centrality->GetCentralityClass5(fCentralityEstimator.Data()) ;
54b31d6a 168
2bbfd8c6 169 // if (centrBin != fCentrBin && fCentrBin != -1) return;
170 // }
a23f7c97 171
172 fHistoManager->GetHistoStats()->Fill(AliAnalysisMultPbTrackHistoManager::kStatCentr);
173
72491d7c 174 Bool_t isSelected = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & fOfflineTrigger);
175
176 if(!isSelected) return;
177 fHistoManager->GetHistoStats()->Fill(AliAnalysisMultPbTrackHistoManager::kStatPhysSel);
4d0aa70f 178
179
a23f7c97 180 if (fIsMC) {
4d0aa70f 181 // Get MC vertex
182 //FIXME: which vertex do I take for MC?
183 TArrayF vertex;
184 fMCEvent->GenEventHeader()->PrimaryVertex(vertex);
185 Float_t zvGen = vertex[2];
a23f7c97 186
187 if (!fMCEvent) {
188 AliError("No MC info found");
189 } else {
190
e0376287 191 //loop on the MC event, only over primaries, which are always
192 // the first in stack.
bcc49144 193 // WRONG: D&B decays are produced in the transport... Need To Loop over all particles
194 // Int_t nMCTracks = fMCEvent->GetNumberOfPrimaries();
195 Int_t nMCTracks = fMCEvent->GetNumberOfTracks();
e0376287 196 Int_t nPhysicalPrimaries = 0;
197 for (Int_t ipart=0; ipart<nMCTracks; ipart++) {
a23f7c97 198
199 AliMCParticle *mcPart = (AliMCParticle*)fMCEvent->GetTrack(ipart);
200
201 // We don't care about neutrals and non-physical primaries
202 if(mcPart->Charge() == 0) continue;
203
54b31d6a 204 //check if current particle is a physical primary
e0376287 205 if(!IsPhysicalPrimaryAndTransportBit(ipart)) continue;
bcc49144 206 if(TMath::Abs(mcPart->Zv()-zvGen)>1e-6) {
207 // This cures a bug in Hijing
208 // A little hack here: I put those in the underflow bin of the process type to keep track of them
209 fHistoManager->GetHistoProcess(AliAnalysisMultPbTrackHistoManager::kHistoGen)->Fill(-1);
210 continue;
211 }
212
e0376287 213 nPhysicalPrimaries++;
3b8cbf2d 214 // Fill species histo and particle species
215 fHistoManager->GetHistoProcess(AliAnalysisMultPbTrackHistoManager::kHistoGen)->Fill(mcPart->Particle()->GetUniqueID());
216 fHistoManager->FillParticleID(AliAnalysisMultPbTrackHistoManager::kHistoGen, mcPart);
e0376287 217
a23f7c97 218 // Float_t zv = vtxESD->GetZ();
219 // Fill generated histo
4d0aa70f 220 hTracks[AliAnalysisMultPbTrackHistoManager::kHistoGen]->Fill(mcPart->Pt(),mcPart->Eta(),zvGen);
221 Int_t partCode = fHistoManager->GetLocalParticleID(mcPart);
abd808b9 222 // cout << "Part " << partCode << endl;
4d0aa70f 223 if (partCode == AliAnalysisMultPbTrackHistoManager::kPartPiPlus ||
bcc49144 224 partCode == AliAnalysisMultPbTrackHistoManager::kPartPiMinus ||
225 partCode == AliAnalysisMultPbTrackHistoManager::kPartKPlus ||
226 partCode == AliAnalysisMultPbTrackHistoManager::kPartKMinus ||
227 partCode == AliAnalysisMultPbTrackHistoManager::kPartP ||
228 partCode == AliAnalysisMultPbTrackHistoManager::kPartPBar
229 ){
abd808b9 230 // cout << "Found " << partCode << endl;
bcc49144 231
4d0aa70f 232 fHistoManager->GetHistoPtEtaVz(AliAnalysisMultPbTrackHistoManager::kHistoGen, partCode);
bcc49144 233 }
a23f7c97 234 }
e0376287 235 hNTracks[AliAnalysisMultPbTrackHistoManager::kHistoGen]->Fill(nPhysicalPrimaries);
4d0aa70f 236 fHistoManager->GetHistoVzEvent(AliAnalysisMultPbTrackHistoManager::kHistoGen)->Fill(zvGen);
237
238
a23f7c97 239 }
240 }
241
e0376287 242
a23f7c97 243 const AliESDVertex* vtxESD = fESD->GetPrimaryVertex();
244 if(!vtxESD) return;
bcc49144 245 // TODO: leave the cuts here or find a better place?)
246 // Quality cut on vertexer Z, as suggested by Francesco Prino
4d0aa70f 247 if(vtxESD->IsFromVertexerZ()) {
248 if (vtxESD->GetNContributors() <= 0) return;
bcf2601a 249 if (vtxESD->GetDispersion() > 0.04) return;
250 if (vtxESD->GetZRes() > 0.25) return;
4d0aa70f 251 }
55e8e56a 252 // "Beam gas" vertex cut
253 const AliESDVertex * vtxESDTPC= fESD->GetPrimaryVertexTPC();
bcc49144 254 if(vtxESDTPC->GetNContributors()<1) return;
bcf2601a 255 if (vtxESDTPC->GetNContributors()<(-10.+0.25*fESD->GetMultiplicity()->GetNumberOfITSClusters(0))) return;
bcc49144 256
bcf2601a 257 // Fill statistics
a23f7c97 258 fHistoManager->GetHistoStats()->Fill(AliAnalysisMultPbTrackHistoManager::kStatVtx);
bcf2601a 259
260 // ZDC cut, only ZNs
bfae2924 261 Bool_t zdcA = fTriggerAnalysis->ZDCTDCTrigger(fESD, AliTriggerAnalysis::kASide, kTRUE, kFALSE) ;
262 Bool_t zdcC = fTriggerAnalysis->ZDCTDCTrigger(fESD, AliTriggerAnalysis::kCSide, kTRUE, kFALSE) ;
bcf2601a 263
bfae2924 264 if (!(zdcA && zdcC) && (!fIsMC)) return;
bcf2601a 265 fHistoManager->GetHistoStats()->Fill(AliAnalysisMultPbTrackHistoManager::kStatZDCCut);
266
abd808b9 267
268 if(TMath::Abs(vtxESD->GetZ()) > 10) return;
269 fHistoManager->GetHistoStats()->Fill(AliAnalysisMultPbTrackHistoManager::kStatVtxRangeCut);
270
271 // FIXME
272 Float_t ntracksOk = fTrackCuts->CountAcceptedTracks(fESD);
273 const AliMultiplicity * mult = fESD->GetMultiplicity();
274 if (ntracksOk < 1000 && mult->GetNumberOfITSClusters(1) > 4500) {
275 Float_t dummy;
276 Printf("IEV: %d, Orbit: %x, Period: %d, BC: %d\n",fESD->GetEventNumberInFile(), fESD->GetOrbitNumber(),fESD->GetPeriodNumber(),fESD->GetBunchCrossNumber());
277 printf("%s ----> Processing event # %lld\n", CurrentFileName(), Entry());
278 cout << "File " << AliAnalysisManager::GetAnalysisManager()->GetTree()->GetCurrentFile()->GetName() << endl;
279 cout << "Nt " << ntracksOk << ", "<< fESD->GetNumberOfTracks() << " V0 " << fCentrSelector->GetCorrV0(fESD,dummy) << " SPD1 " << mult->GetNumberOfITSClusters(1) << endl;
280 }
281
282
bcf2601a 283
284 // Fill Vertex
285
4d0aa70f 286 fHistoManager->GetHistoVzEvent(AliAnalysisMultPbTrackHistoManager::kHistoRec)->Fill(vtxESD->GetZ());
a23f7c97 287
288 // loop on tracks
289 Int_t ntracks = fESD->GetNumberOfTracks();
2bbfd8c6 290 Int_t acceptedTracks = 0;
a23f7c97 291
292 for (Int_t itrack = 0; itrack<ntracks; itrack++) {
9441cdaa 293 AliESDtrack * esdTrack = fIsTPCOnly ? fTrackCuts->GetTPCOnlyTrack(fESD,itrack) : fESD->GetTrack(itrack);// FIXMEL TrackCuts or TrackCutsNoDCA for the TPC?
294 if (!esdTrack) continue;
295
a23f7c97 296 // Fill DCA distibutions, without the DCA cut, Fill the other stuff, with all cuts!
297 Bool_t accepted = fTrackCuts->AcceptTrack(esdTrack);
298 Bool_t acceptedNoDCA = fTrackCutsNoDCA->AcceptTrack(esdTrack);
299
2b42cee2 300 // accepted = accepted && ((fPIDResponse->NumberOfSigmasTPC(esdTrack,AliPID::kElectron) > 2) ||
301 // (fPIDResponse->NumberOfSigmasTPC(esdTrack,AliPID::kPion ) < 1) ||
302 // (fPIDResponse->NumberOfSigmasTPC(esdTrack,AliPID::kProton ) < 1) ||
303 // (fPIDResponse->NumberOfSigmasTPC(esdTrack,AliPID::kKaon ) < 1)
304 // // (esdTrack->P() > 1.2)
305 // );//FIXME SKIP ELECTRONS below p = 1.2 gev, make configurable. Keep particles if they are in the crossing
306
307 // acceptedNoDCA = acceptedNoDCA && ((fPIDResponse->NumberOfSigmasTPC(esdTrack,AliPID::kElectron) > 2) ||
308 // (fPIDResponse->NumberOfSigmasTPC(esdTrack,AliPID::kPion ) < 1) ||
309 // (fPIDResponse->NumberOfSigmasTPC(esdTrack,AliPID::kProton ) < 1) ||
310 // (fPIDResponse->NumberOfSigmasTPC(esdTrack,AliPID::kKaon ) < 1)
311 // // (esdTrack->P() > 1.2)
312 // );//FIXME SKIP ELECTRONS below p = 1.2 gev, make configurable. Keep particles if they are in the crossing
313
2bbfd8c6 314 if(accepted) acceptedTracks++;
315
a23f7c97 316 // Compute weighted offset
317 // TODO: other choiches for the DCA?
a23f7c97 318 Double_t weightedDCA = 10;
319
320
bcc49144 321#if defined WEIGHTED_DCA
322 Double_t b = fESD->GetMagneticField();
323 Double_t dca[2];
324 Double_t cov[3];
a23f7c97 325 if (esdTrack->PropagateToDCA(vtxESD, b,3., dca, cov)) {
326 Double_t det = cov[0]*cov[2]-cov[1]*cov[1];
327 if (det<=0) {
328 AliError("DCA Covariance matrix is not positive definite!");
329 }
330 else {
331 weightedDCA = (dca[0]*dca[0]*cov[2]+dca[1]*dca[1]*cov[0]-2*dca[0]*dca[1]*cov[1])/det;
332 weightedDCA = weightedDCA>0 ? TMath::Sqrt(weightedDCA) : 10;
333 }
334 // printf("dR:%e dZ%e sigR2:%e sigRZ:%e sigZ2:%e\n",dca[0],dca[1],cov[0],cov[1],cov[2]);
335 }
336
bcc49144 337#elif defined TRANSVERSE_DCA
338 Float_t xz[2];
339 esdTrack->GetDZ(vtxESD->GetX(),vtxESD->GetY(),vtxESD->GetZ(), fESD->GetMagneticField(), xz);
340 weightedDCA = xz[0];
341#endif
a23f7c97 342
343
344 // Alternative: p*DCA
345 // Float_t xz[2];
346 // esdTrack->GetDZ(vtxESD->GetX(),vtxESD->GetY(),vtxESD->GetZ(), fESD->GetMagneticField(), xz);
347 // Float_t dca = xz[0]*esdTrack->P();
348
bcc49144 349 // This cures a bug in Hijing (displaced primaries)
350 if (fIsMC) {
351 Int_t label = TMath::Abs(esdTrack->GetLabel()); // no fakes!!!
352 if (IsPhysicalPrimaryAndTransportBit(label)) {
353 AliMCParticle *mcPart = (AliMCParticle*)fMCEvent->GetTrack(label);
354 if(!mcPart) AliFatal("No particle");
355 TArrayF vertex;
356 fMCEvent->GenEventHeader()->PrimaryVertex(vertex);
357 Float_t zvGen = vertex[2];
358 if(TMath::Abs(mcPart->Zv()-zvGen)>1e-6) continue;
359 }
360 }
361 // for each track
abd808b9 362 if(accepted) {
bcc49144 363 hTracks[AliAnalysisMultPbTrackHistoManager::kHistoRec]->Fill(esdTrack->Pt(),esdTrack->Eta(),vtxESD->GetZ());
abd808b9 364 if (TMath::Abs(esdTrack->Eta())<0.5 && TMath::Abs(vtxESD->GetZ())<7)
365 fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRec)->Fill(esdTrack->Pt());
366 }
bcc49144 367 if(acceptedNoDCA)
368 hDCA[AliAnalysisMultPbTrackHistoManager::kHistoRec]->Fill(weightedDCA,esdTrack->Pt());
a23f7c97 369
370 // Fill separately primaries and secondaries
371 // FIXME: fakes? Is this the correct way to account for them?
372 // Get label and corresponding mcPart;
373 if (fIsMC) {
72491d7c 374 Int_t label = TMath::Abs(esdTrack->GetLabel()); // no fakes!!!
e0376287 375 AliMCParticle *mcPart = label < 0 ? 0 : (AliMCParticle*)fMCEvent->GetTrack(label);
a23f7c97 376 if (!mcPart) {
377 if(accepted)
378 hTracks[AliAnalysisMultPbTrackHistoManager::kHistoRecFake]->Fill(esdTrack->Pt(),esdTrack->Eta(),vtxESD->GetZ());
379 if(acceptedNoDCA)
bcc49144 380 hDCA[AliAnalysisMultPbTrackHistoManager::kHistoRecFake]->Fill(weightedDCA,esdTrack->Pt());
a23f7c97 381 }
382 else {
e0376287 383 if(IsPhysicalPrimaryAndTransportBit(label)) {
384 // Fill species histo
3b8cbf2d 385 fHistoManager->GetHistoProcess(AliAnalysisMultPbTrackHistoManager::kHistoRecPrim)->Fill(mcPart->Particle()->GetUniqueID());
4d0aa70f 386 if(accepted) {
a23f7c97 387 hTracks[AliAnalysisMultPbTrackHistoManager::kHistoRecPrim]->Fill(esdTrack->Pt(),esdTrack->Eta(),vtxESD->GetZ());
4d0aa70f 388 Int_t partCode = fHistoManager->GetLocalParticleID(mcPart);
389 if (partCode == AliAnalysisMultPbTrackHistoManager::kPartPiPlus ||
bcc49144 390 partCode == AliAnalysisMultPbTrackHistoManager::kPartPiMinus ||
391 partCode == AliAnalysisMultPbTrackHistoManager::kPartKPlus ||
392 partCode == AliAnalysisMultPbTrackHistoManager::kPartKMinus ||
393 partCode == AliAnalysisMultPbTrackHistoManager::kPartP ||
394 partCode == AliAnalysisMultPbTrackHistoManager::kPartPBar
395 )
4d0aa70f 396 fHistoManager->GetHistoPtEtaVz(AliAnalysisMultPbTrackHistoManager::kHistoRecPrim, partCode);
397 }
a23f7c97 398 if(acceptedNoDCA)
bcc49144 399 hDCA[AliAnalysisMultPbTrackHistoManager::kHistoRecPrim]->Fill(weightedDCA,esdTrack->Pt());
a23f7c97 400 }
401 else {
402 Int_t mfl=0;
403 Int_t indexMoth=mcPart->Particle()->GetFirstMother();
404 if(indexMoth>=0){
405 TParticle* moth = fMCEvent->Stack()->Particle(indexMoth);
406 Float_t codemoth = TMath::Abs(moth->GetPdgCode());
407 mfl = Int_t (codemoth/ TMath::Power(10, Int_t(TMath::Log10(codemoth))));
408 }
409 if(mfl==3){ // strangeness
3b8cbf2d 410 fHistoManager->GetHistoProcess(AliAnalysisMultPbTrackHistoManager::kHistoRecSecWeak)->Fill(mcPart->Particle()->GetUniqueID());
a23f7c97 411 if(accepted)
412 hTracks[AliAnalysisMultPbTrackHistoManager::kHistoRecSecWeak]->Fill(esdTrack->Pt(),esdTrack->Eta(),vtxESD->GetZ());
413 if(acceptedNoDCA)
bcc49144 414 hDCA[AliAnalysisMultPbTrackHistoManager::kHistoRecSecWeak]->Fill(weightedDCA,esdTrack->Pt());
a23f7c97 415 }else{ // material
3b8cbf2d 416 fHistoManager->GetHistoProcess(AliAnalysisMultPbTrackHistoManager::kHistoRecSecMat)->Fill(mcPart->Particle()->GetUniqueID());
a23f7c97 417 if(accepted)
418 hTracks[AliAnalysisMultPbTrackHistoManager::kHistoRecSecMat]->Fill(esdTrack->Pt(),esdTrack->Eta(),vtxESD->GetZ());
419 if(acceptedNoDCA)
bcc49144 420 hDCA[AliAnalysisMultPbTrackHistoManager::kHistoRecSecMat]->Fill(weightedDCA,esdTrack->Pt());
a23f7c97 421
422 }
423
424
425 }
426 }
427 }
428
429
430 }
2bbfd8c6 431 // cout << acceptedTracks << endl;
432
433 hNTracks[AliAnalysisMultPbTrackHistoManager::kHistoRec] ->Fill(acceptedTracks);
84bad6b8 434
435 // FIXME: this used to be filled with rescaled V0. I'm using raw V0 to test some newer production here.
436 AliESDVZERO* esdV0 = fESD->GetVZEROData();
437 Float_t multV0A=esdV0->GetMTotV0A();
438 Float_t multV0C=esdV0->GetMTotV0C();
439 Float_t multV0 = multV0A+multV0C;
440 fHistoManager->GetHistoV0vsNtracks(AliAnalysisMultPbTrackHistoManager::kHistoRec)->Fill(acceptedTracks, multV0);
441 // FIXME: the old version:
442 // Float_t v0;
443 // fHistoManager->GetHistoV0vsNtracks(AliAnalysisMultPbTrackHistoManager::kHistoRec)->Fill(acceptedTracks, fCentrSelector->GetCorrV0(fESD,v0));
a23f7c97 444
abd808b9 445 // Fill <pt> analysis histos
446 fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRec)->Scale(1,"width");// correct bin width
447 if (fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRec)->GetEntries()>0)
448 fHistoManager->GetHistoMeanPt(AliAnalysisMultPbTrackHistoManager::kHistoRec)->Fill(fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRec)->GetMean());
449 if (fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRec)->GetMean() >
450 fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRecHighestMeanPt)->GetMean()) {
451 // Found a new highest pt: first reset than add it to the container for the highest
452 fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRecHighestMeanPt)->Reset();
453 fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRecHighestMeanPt)
454 ->Add(fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRec));
455 }
456 if ((fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRec)->GetMean() <
457 fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRecLowestMeanPt)->GetMean() &&
458 fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRec)->GetEntries()>0) ||
459 !(fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRecLowestMeanPt)->GetEntries()>0)) {
460 // Found a new lowest pt: first reset than add it to the container for the lowest
461 fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRecLowestMeanPt)->Reset();
462 fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRecLowestMeanPt)
463 ->Add(fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRec));
464 }
a23f7c97 465}
466
467void AliAnalysisTaskMultPbTracks::Terminate(Option_t *){
468 // terminate
469
470}
471
472
e0376287 473Bool_t AliAnalysisTaskMultPbTracks::IsPhysicalPrimaryAndTransportBit(Int_t ipart) {
bcc49144 474 // Check if it is a primary and if it was transported
e0376287 475 Bool_t physprim=fMCEvent->IsPhysicalPrimary(ipart);
476 if (!physprim) return kFALSE;
477 Bool_t transported = ((AliMCParticle*)fMCEvent->GetTrack(ipart))->Particle()->TestBit(kTransportBit);
478 if(!transported) return kFALSE;
479
480 return kTRUE;
481
482}
eef42d18 483
484// void AliAnalysisTaskEvil::PrintProcInfo()
485// {
486// ProcInfo_t info;
487// gSystem->GetProcInfo(&info);
488// AliInfo(Form("fMemResident=%10ld kB fMemVirtual=%10ld kB",info.fMemResident,info.fMemVirtual));
489// }