Working on the electron cut
[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"),
5ad99723 33 fESD(0),fHistoManager(0),fCentrSelector(0),fTrackCuts(0),fTrackCutsNoDCA(0),fOfflineTrigger(0), fIsMC(0),fIsTPCOnly(0), fTriggerAnalysis(0),fPIDResponse(0),fRejectElectrons(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),
5ad99723 47 fESD(0),fHistoManager(0),fCentrSelector(0),fTrackCuts(0),fTrackCutsNoDCA(0),fOfflineTrigger(0),fIsMC(0),fIsTPCOnly(0), fTriggerAnalysis(0),fPIDResponse(0),fRejectElectrons(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) :
5ad99723 63 AliAnalysisTaskSE(obj) ,fESD (0), fHistoManager(0), fCentrSelector(0), fTrackCuts(0),fTrackCutsNoDCA(0),fOfflineTrigger(0),fIsMC(0),fIsTPCOnly(0), fTriggerAnalysis(0),fPIDResponse(0),fRejectElectrons(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
5ad99723 113 PostData(1,fHistoManager);
114 PostData(2,fTrackCuts);
115 PostData(3,fCentrSelector);
116
117
a23f7c97 118}
119
120
121void AliAnalysisTaskMultPbTracks::UserExec(Option_t *)
122{
123 // User code
124
e9c445f5 125
a23f7c97 126 /* PostData(0) is taken care of by AliAnalysisTaskSE */
127 PostData(1,fHistoManager);
54b31d6a 128 PostData(2,fTrackCuts);
2bbfd8c6 129 PostData(3,fCentrSelector);
a23f7c97 130
abd808b9 131 // Cache (some) histogram pointers to avoid loop in the histo manager
e0376287 132 static TH3D * hTracks [AliAnalysisMultPbTrackHistoManager::kNHistos];
bcc49144 133 static TH2D * hDCA [AliAnalysisMultPbTrackHistoManager::kNHistos];
e0376287 134 static TH1D * hNTracks [AliAnalysisMultPbTrackHistoManager::kNHistos];
a23f7c97 135 hTracks[AliAnalysisMultPbTrackHistoManager::kHistoGen] = fHistoManager->GetHistoPtEtaVz(AliAnalysisMultPbTrackHistoManager::kHistoGen );
136 hTracks[AliAnalysisMultPbTrackHistoManager::kHistoRec] = fHistoManager->GetHistoPtEtaVz(AliAnalysisMultPbTrackHistoManager::kHistoRec );
137 hTracks[AliAnalysisMultPbTrackHistoManager::kHistoRecPrim] = fHistoManager->GetHistoPtEtaVz(AliAnalysisMultPbTrackHistoManager::kHistoRecPrim );
138 hTracks[AliAnalysisMultPbTrackHistoManager::kHistoRecFake] = fHistoManager->GetHistoPtEtaVz(AliAnalysisMultPbTrackHistoManager::kHistoRecFake );
139 hTracks[AliAnalysisMultPbTrackHistoManager::kHistoRecSecMat] = fHistoManager->GetHistoPtEtaVz(AliAnalysisMultPbTrackHistoManager::kHistoRecSecMat );
140 hTracks[AliAnalysisMultPbTrackHistoManager::kHistoRecSecWeak] = fHistoManager->GetHistoPtEtaVz(AliAnalysisMultPbTrackHistoManager::kHistoRecSecWeak);
141
142 hDCA[AliAnalysisMultPbTrackHistoManager::kHistoGen] = fHistoManager->GetHistoDCA(AliAnalysisMultPbTrackHistoManager::kHistoGen );
143 hDCA[AliAnalysisMultPbTrackHistoManager::kHistoRec] = fHistoManager->GetHistoDCA(AliAnalysisMultPbTrackHistoManager::kHistoRec );
144 hDCA[AliAnalysisMultPbTrackHistoManager::kHistoRecPrim] = fHistoManager->GetHistoDCA(AliAnalysisMultPbTrackHistoManager::kHistoRecPrim );
145 hDCA[AliAnalysisMultPbTrackHistoManager::kHistoRecFake] = fHistoManager->GetHistoDCA(AliAnalysisMultPbTrackHistoManager::kHistoRecFake );
146 hDCA[AliAnalysisMultPbTrackHistoManager::kHistoRecSecMat] = fHistoManager->GetHistoDCA(AliAnalysisMultPbTrackHistoManager::kHistoRecSecMat );
147 hDCA[AliAnalysisMultPbTrackHistoManager::kHistoRecSecWeak] = fHistoManager->GetHistoDCA(AliAnalysisMultPbTrackHistoManager::kHistoRecSecWeak);
148
2bbfd8c6 149 hNTracks[AliAnalysisMultPbTrackHistoManager::kHistoGen] = fHistoManager->GetHistoMult(AliAnalysisMultPbTrackHistoManager::kHistoGen );
150 hNTracks[AliAnalysisMultPbTrackHistoManager::kHistoRec] = fHistoManager->GetHistoMult(AliAnalysisMultPbTrackHistoManager::kHistoRec );
a23f7c97 151
abd808b9 152
153 fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRec)->Reset();
154
a23f7c97 155 fESD = dynamic_cast<AliESDEvent*>(fInputEvent);
156 if (strcmp(fESD->ClassName(),"AliESDEvent")) {
157 AliFatal("Not processing ESDs");
158 }
159
a23f7c97 160 fHistoManager->GetHistoStats()->Fill(AliAnalysisMultPbTrackHistoManager::kStatAll);
161
8b628a20 162
54b31d6a 163 // Centrality selection
2bbfd8c6 164 Bool_t isCentralitySelected = fCentrSelector->IsCentralityBinSelected(fESD,fTrackCuts);
165 if(!isCentralitySelected) return;
166
167 // AliESDCentrality *centrality = fESD->GetCentrality();
168 // if(!centrality) {
169 // AliFatal("Centrality object not available");
170 // }
171 // else {
172 // Int_t centrBin = centrality->GetCentralityClass5(fCentralityEstimator.Data()) ;
54b31d6a 173
2bbfd8c6 174 // if (centrBin != fCentrBin && fCentrBin != -1) return;
175 // }
a23f7c97 176
177 fHistoManager->GetHistoStats()->Fill(AliAnalysisMultPbTrackHistoManager::kStatCentr);
178
72491d7c 179 Bool_t isSelected = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & fOfflineTrigger);
180
181 if(!isSelected) return;
182 fHistoManager->GetHistoStats()->Fill(AliAnalysisMultPbTrackHistoManager::kStatPhysSel);
4d0aa70f 183
184
a23f7c97 185 if (fIsMC) {
4d0aa70f 186 // Get MC vertex
187 //FIXME: which vertex do I take for MC?
188 TArrayF vertex;
189 fMCEvent->GenEventHeader()->PrimaryVertex(vertex);
190 Float_t zvGen = vertex[2];
a23f7c97 191
192 if (!fMCEvent) {
193 AliError("No MC info found");
194 } else {
195
e0376287 196 //loop on the MC event, only over primaries, which are always
197 // the first in stack.
bcc49144 198 // WRONG: D&B decays are produced in the transport... Need To Loop over all particles
199 // Int_t nMCTracks = fMCEvent->GetNumberOfPrimaries();
200 Int_t nMCTracks = fMCEvent->GetNumberOfTracks();
e0376287 201 Int_t nPhysicalPrimaries = 0;
202 for (Int_t ipart=0; ipart<nMCTracks; ipart++) {
a23f7c97 203
204 AliMCParticle *mcPart = (AliMCParticle*)fMCEvent->GetTrack(ipart);
205
206 // We don't care about neutrals and non-physical primaries
207 if(mcPart->Charge() == 0) continue;
208
5ad99723 209 // If we are cutting away electrons, also exclude them here:
210 if(fRejectElectrons)
211 if (TMath::Abs(mcPart->PdgCode())==11) continue;
212
213 //check if current particle is a physical primary and fill map of physical primaries
e0376287 214 if(!IsPhysicalPrimaryAndTransportBit(ipart)) continue;
5ad99723 215 fHistoManager->FillSpeciesMap(mcPart->PdgCode());
216
bcc49144 217 if(TMath::Abs(mcPart->Zv()-zvGen)>1e-6) {
218 // This cures a bug in Hijing
219 // A little hack here: I put those in the underflow bin of the process type to keep track of them
220 fHistoManager->GetHistoProcess(AliAnalysisMultPbTrackHistoManager::kHistoGen)->Fill(-1);
221 continue;
222 }
223
e0376287 224 nPhysicalPrimaries++;
3b8cbf2d 225 // Fill species histo and particle species
226 fHistoManager->GetHistoProcess(AliAnalysisMultPbTrackHistoManager::kHistoGen)->Fill(mcPart->Particle()->GetUniqueID());
227 fHistoManager->FillParticleID(AliAnalysisMultPbTrackHistoManager::kHistoGen, mcPart);
e0376287 228
a23f7c97 229 // Float_t zv = vtxESD->GetZ();
230 // Fill generated histo
4d0aa70f 231 hTracks[AliAnalysisMultPbTrackHistoManager::kHistoGen]->Fill(mcPart->Pt(),mcPart->Eta(),zvGen);
232 Int_t partCode = fHistoManager->GetLocalParticleID(mcPart);
abd808b9 233 // cout << "Part " << partCode << endl;
4d0aa70f 234 if (partCode == AliAnalysisMultPbTrackHistoManager::kPartPiPlus ||
bcc49144 235 partCode == AliAnalysisMultPbTrackHistoManager::kPartPiMinus ||
236 partCode == AliAnalysisMultPbTrackHistoManager::kPartKPlus ||
237 partCode == AliAnalysisMultPbTrackHistoManager::kPartKMinus ||
238 partCode == AliAnalysisMultPbTrackHistoManager::kPartP ||
239 partCode == AliAnalysisMultPbTrackHistoManager::kPartPBar
240 ){
abd808b9 241 // cout << "Found " << partCode << endl;
bcc49144 242
4d0aa70f 243 fHistoManager->GetHistoPtEtaVz(AliAnalysisMultPbTrackHistoManager::kHistoGen, partCode);
bcc49144 244 }
a23f7c97 245 }
e0376287 246 hNTracks[AliAnalysisMultPbTrackHistoManager::kHistoGen]->Fill(nPhysicalPrimaries);
4d0aa70f 247 fHistoManager->GetHistoVzEvent(AliAnalysisMultPbTrackHistoManager::kHistoGen)->Fill(zvGen);
248
249
a23f7c97 250 }
251 }
252
e0376287 253
a23f7c97 254 const AliESDVertex* vtxESD = fESD->GetPrimaryVertex();
255 if(!vtxESD) return;
bcc49144 256 // TODO: leave the cuts here or find a better place?)
257 // Quality cut on vertexer Z, as suggested by Francesco Prino
4d0aa70f 258 if(vtxESD->IsFromVertexerZ()) {
259 if (vtxESD->GetNContributors() <= 0) return;
bcf2601a 260 if (vtxESD->GetDispersion() > 0.04) return;
261 if (vtxESD->GetZRes() > 0.25) return;
4d0aa70f 262 }
55e8e56a 263 // "Beam gas" vertex cut
264 const AliESDVertex * vtxESDTPC= fESD->GetPrimaryVertexTPC();
bcc49144 265 if(vtxESDTPC->GetNContributors()<1) return;
bcf2601a 266 if (vtxESDTPC->GetNContributors()<(-10.+0.25*fESD->GetMultiplicity()->GetNumberOfITSClusters(0))) return;
bcc49144 267
bcf2601a 268 // Fill statistics
a23f7c97 269 fHistoManager->GetHistoStats()->Fill(AliAnalysisMultPbTrackHistoManager::kStatVtx);
bcf2601a 270
271 // ZDC cut, only ZNs
bfae2924 272 Bool_t zdcA = fTriggerAnalysis->ZDCTDCTrigger(fESD, AliTriggerAnalysis::kASide, kTRUE, kFALSE) ;
273 Bool_t zdcC = fTriggerAnalysis->ZDCTDCTrigger(fESD, AliTriggerAnalysis::kCSide, kTRUE, kFALSE) ;
bcf2601a 274
bfae2924 275 if (!(zdcA && zdcC) && (!fIsMC)) return;
bcf2601a 276 fHistoManager->GetHistoStats()->Fill(AliAnalysisMultPbTrackHistoManager::kStatZDCCut);
277
abd808b9 278
279 if(TMath::Abs(vtxESD->GetZ()) > 10) return;
280 fHistoManager->GetHistoStats()->Fill(AliAnalysisMultPbTrackHistoManager::kStatVtxRangeCut);
281
282 // FIXME
283 Float_t ntracksOk = fTrackCuts->CountAcceptedTracks(fESD);
284 const AliMultiplicity * mult = fESD->GetMultiplicity();
285 if (ntracksOk < 1000 && mult->GetNumberOfITSClusters(1) > 4500) {
286 Float_t dummy;
287 Printf("IEV: %d, Orbit: %x, Period: %d, BC: %d\n",fESD->GetEventNumberInFile(), fESD->GetOrbitNumber(),fESD->GetPeriodNumber(),fESD->GetBunchCrossNumber());
288 printf("%s ----> Processing event # %lld\n", CurrentFileName(), Entry());
289 cout << "File " << AliAnalysisManager::GetAnalysisManager()->GetTree()->GetCurrentFile()->GetName() << endl;
290 cout << "Nt " << ntracksOk << ", "<< fESD->GetNumberOfTracks() << " V0 " << fCentrSelector->GetCorrV0(fESD,dummy) << " SPD1 " << mult->GetNumberOfITSClusters(1) << endl;
291 }
292
293
bcf2601a 294
295 // Fill Vertex
296
4d0aa70f 297 fHistoManager->GetHistoVzEvent(AliAnalysisMultPbTrackHistoManager::kHistoRec)->Fill(vtxESD->GetZ());
a23f7c97 298
299 // loop on tracks
300 Int_t ntracks = fESD->GetNumberOfTracks();
2bbfd8c6 301 Int_t acceptedTracks = 0;
a23f7c97 302
303 for (Int_t itrack = 0; itrack<ntracks; itrack++) {
9441cdaa 304 AliESDtrack * esdTrack = fIsTPCOnly ? fTrackCuts->GetTPCOnlyTrack(fESD,itrack) : fESD->GetTrack(itrack);// FIXMEL TrackCuts or TrackCutsNoDCA for the TPC?
305 if (!esdTrack) continue;
306
a23f7c97 307 // Fill DCA distibutions, without the DCA cut, Fill the other stuff, with all cuts!
308 Bool_t accepted = fTrackCuts->AcceptTrack(esdTrack);
309 Bool_t acceptedNoDCA = fTrackCutsNoDCA->AcceptTrack(esdTrack);
310
5ad99723 311 //fHistoManager->GetHistoElectronCutQA()->Fill(esdTrack->Pt(), esdTrack->GetTPCsignal()); FIXME: temporary hack
2b42cee2 312
5ad99723 313 if(fRejectElectrons && acceptedNoDCA) { // don't check this for tracks to be anyways rejected
314 Bool_t isEl = IsElectron(esdTrack);
315 if(isEl) {
316 // cout << " --> FOUND!" << endl;
317 fHistoManager->GetHistoElectronCutQA()->Fill(esdTrack->Pt(), esdTrack->GetTPCsignal());
318 }
319 accepted = accepted && !isEl;
320 acceptedNoDCA = acceptedNoDCA && !isEl;
321 }
2b42cee2 322
2bbfd8c6 323 if(accepted) acceptedTracks++;
324
a23f7c97 325 // Compute weighted offset
326 // TODO: other choiches for the DCA?
a23f7c97 327 Double_t weightedDCA = 10;
328
329
bcc49144 330#if defined WEIGHTED_DCA
331 Double_t b = fESD->GetMagneticField();
332 Double_t dca[2];
333 Double_t cov[3];
a23f7c97 334 if (esdTrack->PropagateToDCA(vtxESD, b,3., dca, cov)) {
335 Double_t det = cov[0]*cov[2]-cov[1]*cov[1];
336 if (det<=0) {
337 AliError("DCA Covariance matrix is not positive definite!");
338 }
339 else {
340 weightedDCA = (dca[0]*dca[0]*cov[2]+dca[1]*dca[1]*cov[0]-2*dca[0]*dca[1]*cov[1])/det;
341 weightedDCA = weightedDCA>0 ? TMath::Sqrt(weightedDCA) : 10;
342 }
343 // printf("dR:%e dZ%e sigR2:%e sigRZ:%e sigZ2:%e\n",dca[0],dca[1],cov[0],cov[1],cov[2]);
344 }
345
bcc49144 346#elif defined TRANSVERSE_DCA
347 Float_t xz[2];
348 esdTrack->GetDZ(vtxESD->GetX(),vtxESD->GetY(),vtxESD->GetZ(), fESD->GetMagneticField(), xz);
349 weightedDCA = xz[0];
350#endif
a23f7c97 351
352
353 // Alternative: p*DCA
354 // Float_t xz[2];
355 // esdTrack->GetDZ(vtxESD->GetX(),vtxESD->GetY(),vtxESD->GetZ(), fESD->GetMagneticField(), xz);
356 // Float_t dca = xz[0]*esdTrack->P();
357
bcc49144 358 // This cures a bug in Hijing (displaced primaries)
359 if (fIsMC) {
360 Int_t label = TMath::Abs(esdTrack->GetLabel()); // no fakes!!!
361 if (IsPhysicalPrimaryAndTransportBit(label)) {
362 AliMCParticle *mcPart = (AliMCParticle*)fMCEvent->GetTrack(label);
363 if(!mcPart) AliFatal("No particle");
364 TArrayF vertex;
365 fMCEvent->GenEventHeader()->PrimaryVertex(vertex);
366 Float_t zvGen = vertex[2];
367 if(TMath::Abs(mcPart->Zv()-zvGen)>1e-6) continue;
368 }
369 }
370 // for each track
abd808b9 371 if(accepted) {
bcc49144 372 hTracks[AliAnalysisMultPbTrackHistoManager::kHistoRec]->Fill(esdTrack->Pt(),esdTrack->Eta(),vtxESD->GetZ());
abd808b9 373 if (TMath::Abs(esdTrack->Eta())<0.5 && TMath::Abs(vtxESD->GetZ())<7)
374 fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRec)->Fill(esdTrack->Pt());
375 }
bcc49144 376 if(acceptedNoDCA)
377 hDCA[AliAnalysisMultPbTrackHistoManager::kHistoRec]->Fill(weightedDCA,esdTrack->Pt());
a23f7c97 378
379 // Fill separately primaries and secondaries
380 // FIXME: fakes? Is this the correct way to account for them?
381 // Get label and corresponding mcPart;
382 if (fIsMC) {
72491d7c 383 Int_t label = TMath::Abs(esdTrack->GetLabel()); // no fakes!!!
e0376287 384 AliMCParticle *mcPart = label < 0 ? 0 : (AliMCParticle*)fMCEvent->GetTrack(label);
a23f7c97 385 if (!mcPart) {
386 if(accepted)
387 hTracks[AliAnalysisMultPbTrackHistoManager::kHistoRecFake]->Fill(esdTrack->Pt(),esdTrack->Eta(),vtxESD->GetZ());
388 if(acceptedNoDCA)
bcc49144 389 hDCA[AliAnalysisMultPbTrackHistoManager::kHistoRecFake]->Fill(weightedDCA,esdTrack->Pt());
a23f7c97 390 }
391 else {
e0376287 392 if(IsPhysicalPrimaryAndTransportBit(label)) {
393 // Fill species histo
3b8cbf2d 394 fHistoManager->GetHistoProcess(AliAnalysisMultPbTrackHistoManager::kHistoRecPrim)->Fill(mcPart->Particle()->GetUniqueID());
4d0aa70f 395 if(accepted) {
a23f7c97 396 hTracks[AliAnalysisMultPbTrackHistoManager::kHistoRecPrim]->Fill(esdTrack->Pt(),esdTrack->Eta(),vtxESD->GetZ());
4d0aa70f 397 Int_t partCode = fHistoManager->GetLocalParticleID(mcPart);
398 if (partCode == AliAnalysisMultPbTrackHistoManager::kPartPiPlus ||
bcc49144 399 partCode == AliAnalysisMultPbTrackHistoManager::kPartPiMinus ||
400 partCode == AliAnalysisMultPbTrackHistoManager::kPartKPlus ||
401 partCode == AliAnalysisMultPbTrackHistoManager::kPartKMinus ||
402 partCode == AliAnalysisMultPbTrackHistoManager::kPartP ||
5ad99723 403 partCode == AliAnalysisMultPbTrackHistoManager::kPartPBar ||
404 partCode == AliAnalysisMultPbTrackHistoManager::kPartLPlus ||
405 partCode == AliAnalysisMultPbTrackHistoManager::kPartKMinus
bcc49144 406 )
4d0aa70f 407 fHistoManager->GetHistoPtEtaVz(AliAnalysisMultPbTrackHistoManager::kHistoRecPrim, partCode);
408 }
a23f7c97 409 if(acceptedNoDCA)
bcc49144 410 hDCA[AliAnalysisMultPbTrackHistoManager::kHistoRecPrim]->Fill(weightedDCA,esdTrack->Pt());
a23f7c97 411 }
412 else {
413 Int_t mfl=0;
414 Int_t indexMoth=mcPart->Particle()->GetFirstMother();
415 if(indexMoth>=0){
416 TParticle* moth = fMCEvent->Stack()->Particle(indexMoth);
417 Float_t codemoth = TMath::Abs(moth->GetPdgCode());
418 mfl = Int_t (codemoth/ TMath::Power(10, Int_t(TMath::Log10(codemoth))));
419 }
420 if(mfl==3){ // strangeness
3b8cbf2d 421 fHistoManager->GetHistoProcess(AliAnalysisMultPbTrackHistoManager::kHistoRecSecWeak)->Fill(mcPart->Particle()->GetUniqueID());
a23f7c97 422 if(accepted)
423 hTracks[AliAnalysisMultPbTrackHistoManager::kHistoRecSecWeak]->Fill(esdTrack->Pt(),esdTrack->Eta(),vtxESD->GetZ());
424 if(acceptedNoDCA)
bcc49144 425 hDCA[AliAnalysisMultPbTrackHistoManager::kHistoRecSecWeak]->Fill(weightedDCA,esdTrack->Pt());
a23f7c97 426 }else{ // material
3b8cbf2d 427 fHistoManager->GetHistoProcess(AliAnalysisMultPbTrackHistoManager::kHistoRecSecMat)->Fill(mcPart->Particle()->GetUniqueID());
a23f7c97 428 if(accepted)
429 hTracks[AliAnalysisMultPbTrackHistoManager::kHistoRecSecMat]->Fill(esdTrack->Pt(),esdTrack->Eta(),vtxESD->GetZ());
430 if(acceptedNoDCA)
bcc49144 431 hDCA[AliAnalysisMultPbTrackHistoManager::kHistoRecSecMat]->Fill(weightedDCA,esdTrack->Pt());
a23f7c97 432
433 }
434
435
436 }
437 }
438 }
439
440
441 }
2bbfd8c6 442 // cout << acceptedTracks << endl;
443
444 hNTracks[AliAnalysisMultPbTrackHistoManager::kHistoRec] ->Fill(acceptedTracks);
84bad6b8 445
446 // FIXME: this used to be filled with rescaled V0. I'm using raw V0 to test some newer production here.
447 AliESDVZERO* esdV0 = fESD->GetVZEROData();
448 Float_t multV0A=esdV0->GetMTotV0A();
449 Float_t multV0C=esdV0->GetMTotV0C();
450 Float_t multV0 = multV0A+multV0C;
451 fHistoManager->GetHistoV0vsNtracks(AliAnalysisMultPbTrackHistoManager::kHistoRec)->Fill(acceptedTracks, multV0);
452 // FIXME: the old version:
453 // Float_t v0;
454 // fHistoManager->GetHistoV0vsNtracks(AliAnalysisMultPbTrackHistoManager::kHistoRec)->Fill(acceptedTracks, fCentrSelector->GetCorrV0(fESD,v0));
a23f7c97 455
abd808b9 456 // Fill <pt> analysis histos
457 fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRec)->Scale(1,"width");// correct bin width
458 if (fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRec)->GetEntries()>0)
459 fHistoManager->GetHistoMeanPt(AliAnalysisMultPbTrackHistoManager::kHistoRec)->Fill(fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRec)->GetMean());
460 if (fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRec)->GetMean() >
461 fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRecHighestMeanPt)->GetMean()) {
462 // Found a new highest pt: first reset than add it to the container for the highest
463 fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRecHighestMeanPt)->Reset();
464 fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRecHighestMeanPt)
465 ->Add(fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRec));
466 }
467 if ((fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRec)->GetMean() <
468 fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRecLowestMeanPt)->GetMean() &&
469 fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRec)->GetEntries()>0) ||
470 !(fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRecLowestMeanPt)->GetEntries()>0)) {
471 // Found a new lowest pt: first reset than add it to the container for the lowest
472 fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRecLowestMeanPt)->Reset();
473 fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRecLowestMeanPt)
474 ->Add(fHistoManager->GetHistoPtEvent(AliAnalysisMultPbTrackHistoManager::kHistoRec));
475 }
a23f7c97 476}
477
478void AliAnalysisTaskMultPbTracks::Terminate(Option_t *){
479 // terminate
480
481}
482
483
e0376287 484Bool_t AliAnalysisTaskMultPbTracks::IsPhysicalPrimaryAndTransportBit(Int_t ipart) {
bcc49144 485 // Check if it is a primary and if it was transported
e0376287 486 Bool_t physprim=fMCEvent->IsPhysicalPrimary(ipart);
487 if (!physprim) return kFALSE;
488 Bool_t transported = ((AliMCParticle*)fMCEvent->GetTrack(ipart))->Particle()->TestBit(kTransportBit);
489 if(!transported) return kFALSE;
490
491 return kTRUE;
492
493}
eef42d18 494
5ad99723 495Bool_t AliAnalysisTaskMultPbTracks::IsElectron(AliESDtrack * esdTrack) {
496 Bool_t isElectron = kFALSE;
497 // ((fPIDResponse->NumberOfSigmasTPC(esdTrack,AliPID::kElectron) < 2 ) &&
498 // (fPIDResponse->NumberOfSigmasTPC(esdTrack,AliPID::kPion ) > 1 ) &&
499 // (fPIDResponse->NumberOfSigmasTPC(esdTrack,AliPID::kProton ) > 1 ) &&
500 // (fPIDResponse->NumberOfSigmasTPC(esdTrack,AliPID::kKaon ) > 1 )
501 // );//FIXME SKIP ELECTRONS below p = 1.2 gev, make configurable. Keep particles if they are in the crossing
502
503 // cout << " - TPC: " << fPIDResponse->NumberOfSigmasTPC(esdTrack,AliPID::kElectron) << ", TOF: " << fPIDResponse->NumberOfSigmasTOF(esdTrack,AliPID::kElectron);
504 // printf("[%f,%f,%d]\n", esdTrack->Pt(), esdTrack->Eta(), esdTrack->GetTPCclusters(0));
505 // if (fPIDResponse->NumberOfSigmasTPC(esdTrack,AliPID::kElectron) <0){
506
507
508 // }
509 isElectron = (TMath::Abs(fPIDResponse->NumberOfSigmasTPC(esdTrack,AliPID::kElectron)) < 2) && (TMath::Abs(fPIDResponse->NumberOfSigmasTOF(esdTrack,AliPID::kElectron) < 3));
510
511 return isElectron;
512}
513
eef42d18 514// void AliAnalysisTaskEvil::PrintProcInfo()
515// {
516// ProcInfo_t info;
517// gSystem->GetProcInfo(&info);
518// AliInfo(Form("fMemResident=%10ld kB fMemVirtual=%10ld kB",info.fMemResident,info.fMemVirtual));
519// }