]>
Commit | Line | Data |
---|---|---|
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 | |
27 | using namespace std; | |
28 | ||
29 | ClassImp(AliAnalysisTaskMultPbTracks) | |
30 | ||
31 | AliAnalysisTaskMultPbTracks::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 | } | |
45 | AliAnalysisTaskMultPbTracks::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 | ||
62 | AliAnalysisTaskMultPbTracks::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 | ||
78 | AliAnalysisTaskMultPbTracks::~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 | } | |
92 | void 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 | ||
121 | void 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 | ||
478 | void AliAnalysisTaskMultPbTracks::Terminate(Option_t *){ | |
479 | // terminate | |
480 | ||
481 | } | |
482 | ||
483 | ||
e0376287 | 484 | Bool_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 | 495 | Bool_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 | // } |