git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
coverity fixes (Saehanseul Oh <saehanseul.oh@cern.ch>)
[u/mrichter/AliRoot.git] / PWGCF / Correlations / DPhi / FourierDecomposition / AliMuonEffMC.cxx
CommitLineData
7e990b20 1// MUON track QA referring AliMuonEffMC.cxx
2// Author : Saehanseul Oh
3
4#include "AliMuonEffMC.h"
5
6#include <TList.h>
7#include <TH1D.h>
8#include <TH2F.h>
9#include <TH3F.h>
10#include <THn.h>
11#include <TChain.h>
12#include <TFile.h>
13
15b0d1bb 14#include "TParticle.h"
15#include "AliStack.h"
7e990b20 16#include "AliAnalysisManager.h"
17#include "AliAnalysisTask.h"
18#include "AliESDVertex.h"
19#include "AliESDEvent.h"
20#include "AliAODEvent.h"
21#include "AliESDMuonTrack.h"
b673a083 22#include "AliAODTrack.h"
7e990b20 23#include "AliESDVertex.h"
24#include "AliAODVertex.h"
25#include "AliCentrality.h"
26#include "AliVParticle.h"
27#include "AliMCParticle.h"
28#include "AliAnalysisHelperJetTasks.h"
29#include "AliMCEvent.h"
fe76c976 30#include "AliAODMCParticle.h"
7e990b20 31
32using std::cout;
33using std::endl;
34
35ClassImp(AliMuonEffMC)
36
37//________________________________________________________________________
38AliMuonEffMC::AliMuonEffMC() :
15b0d1bb 39 AliAnalysisTaskSE(), fESD(0), fAOD(0), fMC(0), fStack(0), fCentrality(99), fZVertex(99), fOutputList(0x0),
40 fHEventStat(0), fHEvt(0x0), fIsMc(kFALSE), fMDProcess(kFALSE), fCentralityEstimator("V0M"),
41 fNEtaBins(15), fNpTBins(100), fNCentBins(5), fNZvtxBins(10), fNPhiBins(36),
42 fHMuonParGen(0x0), fHMuonDetGen(0x0), fHMuonDetRec(0x0), fHEtcDetRec(0x0), fHFXu(0), fHFXantiu(0), fHFXd(0), fHFXantid(0), fHFXg(0), fHFXetc(0), fHFXmuonP(0), fHFXmuonM(0)
7e990b20 43{
44 // Constructor
a3cd18eb 45 for(Int_t i=0; i<4; i++)
46 {
47 fHMuMotherGenPt[i] = NULL;
48 fHMuMotherRecPt[i] = NULL;
49 fHMuMotherGenPhi[i] = NULL;
50 fHMuMotherRecPhi[i] = NULL;
51 fHMuMotherGenEta[i] = NULL;
52 fHMuMotherRecEta[i] = NULL;
53 fHMuDCA[i] = 0x0;
54 for(Int_t j=0; j<3; j++)
55 {
56 fHMuMohterPhiDifGen[i][j] = NULL;
57 fHMuMohterPhiDifRec[i][j] = NULL;
58 fHMuMohterEtaDifGen[i][j] = NULL;
59 fHMuMohterEtaDifRec[i][j] = NULL;
60 }
61 }
7e990b20 62}
63
64//________________________________________________________________________
65AliMuonEffMC::AliMuonEffMC(const char *name) :
15b0d1bb 66 AliAnalysisTaskSE(name), fESD(0), fAOD(0), fMC(0), fStack(0), fCentrality(99), fZVertex(99), fOutputList(0x0),
67 fHEventStat(0), fHEvt(0x0), fIsMc(kFALSE), fMDProcess(kFALSE), fCentralityEstimator("V0M"),
68 fNEtaBins(15), fNpTBins(100), fNCentBins(5), fNZvtxBins(10), fNPhiBins(36),
69 fHMuonParGen(0x0), fHMuonDetGen(0x0), fHMuonDetRec(0x0), fHEtcDetRec(0x0), fHFXu(0), fHFXantiu(0), fHFXd(0), fHFXantid(0), fHFXg(0), fHFXetc(0), fHFXmuonP(0), fHFXmuonM(0)
7e990b20 70{
71 // Constructor
a3cd18eb 72 for(Int_t i=0; i<4; i++)
73 {
74 fHMuMotherGenPt[i] = NULL;
75 fHMuMotherRecPt[i] = NULL;
76 fHMuMotherGenPhi[i] = NULL;
77 fHMuMotherRecPhi[i] = NULL;
78 fHMuMotherGenEta[i] = NULL;
79 fHMuMotherRecEta[i] = NULL;
80 fHMuDCA[i] = 0x0;
81 for(Int_t j=0; j<3; j++)
82 {
83 fHMuMohterPhiDifGen[i][j] = NULL;
84 fHMuMohterPhiDifRec[i][j] = NULL;
85 fHMuMohterEtaDifGen[i][j] = NULL;
86 fHMuMohterEtaDifRec[i][j] = NULL;
87 }
88 }
7e990b20 89 DefineInput(0, TChain::Class());
90 DefineOutput(1, TList::Class());
fe76c976 91
7e990b20 92}
93
94//________________________________________________________________________
95AliMuonEffMC::~AliMuonEffMC()
96{
97 //Destructor
9370528f 98 if(fOutputList) delete fOutputList;
7e990b20 99}
100
101//________________________________________________________________________
102void AliMuonEffMC::UserCreateOutputObjects()
103{
104 // Create histograms
105 // Called once (per slave on PROOF!)
fe76c976 106
107 fOutputList = new TList();
108 fOutputList->SetOwner(1);
7e990b20 109
7e990b20 110 fHEventStat = new TH1D("fHEventStat","Event statistics for analysis",18,0,18);
111 fHEventStat->GetXaxis()->SetBinLabel(1,"Event");
112 fHEventStat->GetXaxis()->SetBinLabel(2,"SelectedEvent");
113 fHEventStat->GetXaxis()->SetBinLabel(3,"File");
114 fHEventStat->GetXaxis()->SetBinLabel(4,"fSPHighpt"); //!Global Trigger Single plus High p_T
115 fHEventStat->GetXaxis()->SetBinLabel(5,"fSPAllpt"); //!Global Trigger Single plus All p_T
116 fHEventStat->GetXaxis()->SetBinLabel(6,"fSMLowpt"); //!Global Trigger Single minus Low p_T
117 fHEventStat->GetXaxis()->SetBinLabel(7,"fSMHighpt"); //!Global Trigger Single minus High p_T
118 fHEventStat->GetXaxis()->SetBinLabel(8,"fSMAllpt"); //!Global Trigger Single minus All p_T
119 fHEventStat->GetXaxis()->SetBinLabel(9,"fSULowpt"); //!Global Trigger Single undefined Low p_T
120 fHEventStat->GetXaxis()->SetBinLabel(10,"fSUHighpt"); //!Global Trigger Single undefined High p_T
121 fHEventStat->GetXaxis()->SetBinLabel(11,"fSUAllpt"); //!Global Trigger Single undefined All p_T
122 fHEventStat->GetXaxis()->SetBinLabel(12,"fUSLowpt"); //!Global Trigger UnlikeSign pair Low p_T
123 fHEventStat->GetXaxis()->SetBinLabel(13,"fUSHighpt"); //!Global Trigger UnlikeSign pair High p_T
124 fHEventStat->GetXaxis()->SetBinLabel(14,"fUSAllpt"); //!Global Trigger UnlikeSign pair All p_T
125 fHEventStat->GetXaxis()->SetBinLabel(15,"fLSLowpt"); //!Global Trigger LikeSign pair pair Low p_T
126 fHEventStat->GetXaxis()->SetBinLabel(16,"fLSHighpt"); //!Global Trigger LikeSign pair pair High p_T
127 fHEventStat->GetXaxis()->SetBinLabel(17,"fLSAllpt"); //!Global Trigger LikeSign pair pair All p_T
128 fHEventStat->GetXaxis()->SetBinLabel(18,"fSPLowpt"); //!Global Trigger Single plus Low p_T
129 fOutputList->Add(fHEventStat);
130
131 fHEvt = new TH2F("fHEvt", "Event-level variables; Zvtx; Cent", 30, -15, 15, 103, -2, 101);
132 fOutputList->Add(fHEvt);
133
134
135 Int_t iTrackBin[5];
136 Double_t* trackBins[5];
137 const char* trackAxisTitle[5];
138
139 // eta
140 Double_t etaBins[fNEtaBins+1];
141 for (Int_t i=0; i<=fNEtaBins; i++) { etaBins[i] = (Double_t)(-4.0 + 1.5/fNEtaBins*i); }
142 iTrackBin[0] = fNEtaBins;
143 trackBins[0] = etaBins;
144 trackAxisTitle[0] = "#eta";
145
146 // p_T
147 Double_t pTBins[fNpTBins+1];
148 for (Int_t i=0; i<=fNpTBins; i++) { pTBins[i] = (Double_t)(10.0/fNpTBins*i); }
149 iTrackBin[1] = fNpTBins;
150 trackBins[1] = pTBins;
151 trackAxisTitle[1] = "p_{T} (GeV/c)";
152
153 // centrality
154 Double_t CentBins[fNCentBins+1];
155 for (Int_t i=0; i<=fNCentBins; i++) { CentBins[i] = (Double_t)(100.0/fNCentBins*i); }
156 iTrackBin[2] = fNCentBins;
157 trackBins[2] = CentBins;
158 trackAxisTitle[2] = "Cent";
159
160 // Z-vertex
161 Double_t ZvtxBins[fNZvtxBins+1];
162 for (Int_t i=0; i<=fNZvtxBins; i++) { ZvtxBins[i] = (Double_t)(-10.0 + 20.0/fNZvtxBins*i); }
163 iTrackBin[3] = fNZvtxBins;
164 trackBins[3] = ZvtxBins;
165 trackAxisTitle[3] = "Zvtx";
166
167 // phi
168 Double_t phiBins[fNPhiBins+1];
169 for (Int_t i=0; i<=fNPhiBins; i++) { phiBins[i] = (Double_t)(TMath::TwoPi()/fNPhiBins * i); }
170 iTrackBin[4] = fNPhiBins;
171 trackBins[4] = phiBins;
172 trackAxisTitle[4] = "#phi";
173
174 // THn for tracking efficiency
175 fHMuonParGen = new THnF("fHMuonParGen", "", 5, iTrackBin, 0, 0);
176 for (Int_t i=0; i<5; i++)
177 {
178 fHMuonParGen->SetBinEdges(i, trackBins[i]);
179 fHMuonParGen->GetAxis(i)->SetTitle(trackAxisTitle[i]);
180 }
181 fHMuonParGen->Sumw2();
182 fOutputList->Add(fHMuonParGen);
183
184 fHMuonDetGen = (THnF*) fHMuonParGen->Clone("fHMuonDetGen");
185 fHMuonDetGen->Sumw2();
186 fOutputList->Add(fHMuonDetGen);
187
188 fHMuonDetRec = (THnF*) fHMuonParGen->Clone("fHMuonDetRec");
189 fHMuonDetRec->Sumw2();
190 fOutputList->Add(fHMuonDetRec);
191
192 fHEtcDetRec = (THnF*) fHMuonParGen->Clone("fHEtcDetRec");
193 fHEtcDetRec->Sumw2();
194 fOutputList->Add(fHEtcDetRec);
195
15b0d1bb 196 if(fIsMc)
7e990b20 197 {
198 const char* MotherSpecies[4] = {"Pion","Kaon","D", "Etc"};
15b0d1bb 199 const char *MuPt[3] = {"0510","1020","2040"};
200 if(fMDProcess)
7e990b20 201 {
15b0d1bb 202 for(Int_t i=0; i<4; i++)
203 {
204 fHMuMotherGenPt[i] = new TH2F(Form("fHMuMotherGenPt_%s",MotherSpecies[i]),";p_{T,muon}^{gen} (GeV/c);p_{T,mother}^{Truth} (GeV/c);",500, 0, 50, 500, 0, 50);
205 fOutputList->Add(fHMuMotherGenPt[i]);
206 fHMuMotherRecPt[i] = new TH2F(Form("fHMuMotherRecPt_%s",MotherSpecies[i]),";p_{T,muon}^{rec} (GeV/c);p_{T,mother}^{Truth} (GeV/c);",500, 0, 50, 500, 0, 50);
207 fOutputList->Add(fHMuMotherRecPt[i]);
208 fHMuMotherGenPhi[i] = new TH2F(Form("fHMuMotherGenPhi_%s",MotherSpecies[i]),";#phi_{gen};mother #phi;",100, 0, TMath::TwoPi(), 100, 0, TMath::TwoPi());
209 fOutputList->Add(fHMuMotherGenPhi[i]);
210 fHMuMotherRecPhi[i] = new TH2F(Form("fHMuMotherRecPhi_%s",MotherSpecies[i]),";#phi_{rec};mother #phi;",100, 0, TMath::TwoPi(), 100, 0, TMath::TwoPi());
211 fOutputList->Add(fHMuMotherRecPhi[i]);
212 fHMuMotherGenEta[i] = new TH2F(Form("fHMuMotherGenEta_%s",MotherSpecies[i]),";#eta_{gen};mother #eta;",100, -5., -1., 100, -5., -1.);
213 fOutputList->Add(fHMuMotherGenEta[i]);
214 fHMuMotherRecEta[i] = new TH2F(Form("fHMuMotherRecEta_%s",MotherSpecies[i]),";#eta_{rec};mother #eta;",100, -5., -1., 100, -5., -1.);
215 fOutputList->Add(fHMuMotherRecEta[i]);
216 fHMuDCA[i] = new TH1F(Form("fHMuDCA_%s",MotherSpecies[i]), ";DCA", 100, 0, 50);
217 fOutputList->Add(fHMuDCA[i]);
218
219 for(Int_t j=0; j<3; j++)
220 {
221 fHMuMohterPhiDifGen[i][j] = new TH1F(Form("fHMuMohterPhiDifGen_%s_%s",MotherSpecies[i], MuPt[j]),";#Delta#phi",100, -1.0*TMath::Pi(), TMath::Pi());
222 fOutputList->Add(fHMuMohterPhiDifGen[i][j]);
223 fHMuMohterPhiDifRec[i][j] = new TH1F(Form("fHMuMohterPhiDifRec_%s_%s",MotherSpecies[i], MuPt[j]),";#Delta#phi",100, -1.0*TMath::Pi(), TMath::Pi());
224 fOutputList->Add(fHMuMohterPhiDifRec[i][j]);
225 fHMuMohterEtaDifGen[i][j] = new TH1F(Form("fHMuMohterEtaDifGen_%s_%s",MotherSpecies[i], MuPt[j]),";#Delta#eta",100, -5.0, 5.0);
226 fOutputList->Add(fHMuMohterEtaDifGen[i][j]);
227 fHMuMohterEtaDifRec[i][j] = new TH1F(Form("fHMuMohterEtaDifRec_%s_%s",MotherSpecies[i], MuPt[j]),";#Delta#eta",100, -5.0, 5.0);
228 fOutputList->Add(fHMuMohterEtaDifRec[i][j]);
229 }
230 }
7e990b20 231 }
15b0d1bb 232
233 fHFXu = new TH1F("fHFXu",";x_{F}",1000, 0.0, 1.0);
234 fOutputList->Add(fHFXu);
235
236 fHFXantiu = new TH1F("fHFXantiu",";x_{F}",1000, 0.0, 1.0);
237 fOutputList->Add(fHFXantiu);
238
239 fHFXd = new TH1F("fHFXd",";x_{F}",1000, 0.0, 1.0);
240 fOutputList->Add(fHFXd);
241
242 fHFXantid = new TH1F("fHFXantid",";x_{F}",1000, 0.0, 1.0);
243 fOutputList->Add(fHFXantid);
244
245 fHFXg = new TH1F("fHFXg",";x_{F}",1000, 0.0, 1.0);
246 fOutputList->Add(fHFXg);
247
248 fHFXetc = new TH1F("fHFXetc",";x_{F}",1000, 0.0, 1.0);
249 fOutputList->Add(fHFXetc);
250
251 fHFXmuonP = new TH1F("fHFXmuonP",";x_{F}",1000, 0.0, 1.0);
252 fOutputList->Add(fHFXmuonP);
253
254 fHFXmuonM = new TH1F("fHFXmuonM",";x_{F}",1000, 0.0, 1.0);
255 fOutputList->Add(fHFXmuonM);
7e990b20 256 }
257 PostData(1, fOutputList);
258}
259
260//________________________________________________________________________
261void AliMuonEffMC::UserExec(Option_t *)
262{
263 // Main loop, Called for each event
b673a083 264 Int_t ntrks = 0; // number of tracks in an event
265
266 if(((TString)InputEvent()->IsA()->GetName())=="AliAODEvent")
7e990b20 267 {
268 fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
b673a083 269 if (!fAOD) { AliError("AOD event not found. Nothing done!"); return; }
270 ntrks = fAOD->GetNTracks();
271 }
272 else
7e990b20 273 {
b673a083 274 fESD = dynamic_cast<AliESDEvent*>(InputEvent());
275 if (!fESD) { AliError("ESD event not found. Nothing done!"); return; }
276 ntrks = fESD->GetNumberOfMuonTracks();
7e990b20 277 }
278
279 fHEventStat->Fill(0.5);
280 if(fIsMc)
281 {
282 fMC = MCEvent();
7e990b20 283 }
284
b673a083 285 // Centrality, vertex, other event variables...
286 if(fAOD)
287 {
288 const AliAODVertex* vertex = fAOD->GetPrimaryVertex();
289 fZVertex = vertex->GetZ();
fe76c976 290 if(fAOD->GetCentrality()) fCentrality = fAOD->GetCentrality()->GetCentralityPercentile(fCentralityEstimator);
7e990b20 291 }
b673a083 292 else if(fESD)
293 {
294 const AliESDVertex* vertex = fESD->GetPrimaryVertex();
295 fZVertex = vertex->GetZ();
fe76c976 296 if(fESD->GetCentrality()) fCentrality = fESD->GetCentrality()->GetCentralityPercentile(fCentralityEstimator);
b673a083 297 }
7e990b20 298
15b0d1bb 299 if ((fESD && !VertexOk(fESD)) || (fAOD && !VertexOk(fAOD))) { //AliInfo(Form("Event REJECTED. z = %.1f", fZVertex));
300 return; }
301 if (fCentrality > 100. || fCentrality < -1.5) { //AliInfo(Form("Event REJECTED. fCentrality = %.1f", fCentrality));
302 return; }
7e990b20 303
15b0d1bb 304 if(fCentrality < 0) fCentrality = 0.0; //ad hoc centrality for pp
b673a083 305 // Fill Event histogram
306 fHEvt->Fill(fZVertex, fCentrality);
7e990b20 307 fHEventStat->Fill(1.5);
308
b673a083 309 ULong64_t trigword = 0;
310 if(fAOD) trigword=fAOD->GetTriggerMask();
311 else if(fESD) trigword=fESD->GetTriggerMask();
7e990b20 312
313 if (trigword & 0x01) fHEventStat->Fill(17.5);
314 if (trigword & 0x02) fHEventStat->Fill(3.5);
315 if (trigword & 0x04) fHEventStat->Fill(4.5);
316 if (trigword & 0x08) fHEventStat->Fill(5.5);
317 if (trigword & 0x010) fHEventStat->Fill(6.5);
318 if (trigword & 0x020) fHEventStat->Fill(7.5);
319 if (trigword & 0x040) fHEventStat->Fill(8.5);
320 if (trigword & 0x080) fHEventStat->Fill(9.5);
321 if (trigword & 0x100) fHEventStat->Fill(10.5);
322 if (trigword & 0x200) fHEventStat->Fill(11.5);
323 if (trigword & 0x400) fHEventStat->Fill(12.5);
324 if (trigword & 0x800) fHEventStat->Fill(13.5);
325 if (trigword & 0x1000) fHEventStat->Fill(14.5);
326 if (trigword & 0x2000) fHEventStat->Fill(15.5);
327 if (trigword & 0x4000) fHEventStat->Fill(16.5);
fe76c976 328
329 // generated level loop
15b0d1bb 330 if(fIsMc)
7e990b20 331 {
15b0d1bb 332 for (Int_t ipart=0; ipart<fMC->GetNumberOfTracks(); ipart++)
fe76c976 333 {
15b0d1bb 334 if(fAOD)
7e990b20 335 {
15b0d1bb 336 AliAODMCParticle *AodMcParticle = (AliAODMCParticle*)fMC->GetTrack(ipart);
337 if(TMath::Abs(AodMcParticle->PdgCode())==13)
338 {
339 Double_t fillArrayParGen[5] = { AodMcParticle->Eta(), AodMcParticle->Pt(), fCentrality, fZVertex, AodMcParticle->Phi() };
340 fHMuonParGen->Fill(fillArrayParGen);
341 }
342 }
343 else if(fESD)
fe76c976 344 {
15b0d1bb 345 AliMCParticle *McParticle = (AliMCParticle*)fMC->GetTrack(ipart);
346 if(TMath::Abs(McParticle->PdgCode())==13)
347 {
348 Double_t fillArrayParGen[5] = { McParticle->Eta(), McParticle->Pt(), fCentrality, fZVertex, McParticle->Phi() };
349 fHMuonParGen->Fill(fillArrayParGen);
350 }
b673a083 351 }
fe76c976 352 }
353 }
354
355 // reconstructed level loop
356 for (Int_t iTrack = 0; iTrack<ntrks; iTrack++)
357 {
358 Int_t label = 0;
359 Double_t trackpt = 0;
360 Double_t tracketa = 0;
361 Double_t trackphi = 0;
362 Double_t dcavalue = 0;
363
364 Double_t mcpt = 0;
365 Double_t mceta = 0;
366 Double_t mcphi = 0;
367
368 Int_t motherlabel =0;
369 Int_t motherpdg = 0;
370 Double_t motherpt = 0;
371 Double_t mothereta = 0;
372 Double_t motherphi = 0;
373
374 if(fAOD)
375 {
376 AliAODTrack* muonTrack = (AliAODTrack*)fAOD->GetTrack(iTrack);
377 if(muonTrack)
b673a083 378 {
fe76c976 379 if(!(IsGoodMUONtrack(*muonTrack)) || !(muonTrack->IsMuonTrack())) continue;
380 trackpt = muonTrack->Pt();
381 tracketa = muonTrack->Eta();
382 trackphi = muonTrack->Phi();
383 label = TMath::Abs(muonTrack->GetLabel());
b673a083 384 }
fe76c976 385 }
386 else if(fESD)
387 {
388 AliESDMuonTrack* muonTrack = fESD->GetMuonTrack(iTrack);
389 if(muonTrack)
390 {
391 if(!IsGoodMUONtrack(*muonTrack)) continue;
392 trackpt = muonTrack->Pt();
393 tracketa = muonTrack->Eta();
394 trackphi = muonTrack->Phi();
395 label = TMath::Abs(muonTrack->GetLabel());
396 dcavalue = muonTrack->GetDCA();
397 }
398 }
fe76c976 399 Double_t fillArrayDetRec[5] = { tracketa, trackpt, fCentrality, fZVertex, trackphi };
400 fHMuonDetRec->Fill(fillArrayDetRec);
401
15b0d1bb 402 if(fIsMc)
fe76c976 403 {
15b0d1bb 404 if(fAOD)
b673a083 405 {
15b0d1bb 406 AliAODMCParticle *aodMcParticle = (AliAODMCParticle*)fMC->GetTrack(label);
407 if(TMath::Abs(aodMcParticle->PdgCode()) != 13)
408 {
409 fHEtcDetRec->Fill(fillArrayDetRec);
410 continue;
411 }
412 mcpt = aodMcParticle->Pt();
413 mceta = aodMcParticle->Eta();
414 mcphi = aodMcParticle->Phi();
415 motherlabel = aodMcParticle->GetMother();
416 if(motherlabel > 0)
417 {
418 AliAODMCParticle *aodMotherParticle = (AliAODMCParticle*)fMC->GetTrack(motherlabel);
419 motherpdg = TMath::Abs(aodMotherParticle->PdgCode());
420 motherpt = aodMotherParticle->Pt();
421 mothereta = aodMotherParticle->Eta();
422 motherphi = aodMotherParticle->Phi();
423 }
fe76c976 424 }
15b0d1bb 425 else if(fESD)
7e990b20 426 {
15b0d1bb 427 AliMCParticle *McParticle = (AliMCParticle*)fMC->GetTrack(label);
428 if(TMath::Abs(McParticle->PdgCode()) != 13)
fe76c976 429 {
15b0d1bb 430 fHEtcDetRec->Fill(fillArrayDetRec);
431 continue;
fe76c976 432 }
15b0d1bb 433 mcpt = McParticle->Pt();
434 mceta = McParticle->Eta();
435 mcphi = McParticle->Phi();
436 motherlabel = McParticle->GetMother();
437 if(motherlabel > 0)
438 {
439 AliMCParticle *MotherParticle = (AliMCParticle*)fMC->GetTrack(motherlabel);
440 motherpdg = TMath::Abs(MotherParticle->PdgCode());
441 motherpt = MotherParticle->Pt();
442 mothereta = MotherParticle->Eta();
443 motherphi = MotherParticle->Phi();
fe76c976 444 }
15b0d1bb 445 }
446 Double_t fillArrayDetGen[5] = { mceta, mcpt, fCentrality, fZVertex, mcphi };
447 fHMuonDetGen->Fill(fillArrayDetGen);
448
449 // mother-daughter kinematic relation
450 if(fMDProcess && motherlabel>0) MDProcess(motherpdg, mcpt, mcphi, mceta, trackpt, trackphi, tracketa, motherpt, motherphi, mothereta, dcavalue);
451 }
452 }
7e990b20 453 PostData(1, fOutputList);
454 return;
455}
456
457//________________________________________________________________________
458void AliMuonEffMC::Terminate(Option_t *)
459{
460 // Draw result to the screen
461 // Called once at the end of the query
7e990b20 462}
7e990b20 463//________________________________________________________________________
464Bool_t AliMuonEffMC::VertexOk(TObject* obj) const
465{
466 // Modified from AliAnalyseLeadingTrackUE::VertexSelection()
467
468 Int_t nContributors = 0;
469 Double_t zVertex = 999;
470 TString name("");
471
472 if (obj->InheritsFrom("AliESDEvent")) {
473 AliESDEvent* esdevt = (AliESDEvent*) obj;
474 const AliESDVertex* vtx = esdevt->GetPrimaryVertex();
475 if (!vtx)
476 return 0;
477 nContributors = vtx->GetNContributors();
478 zVertex = vtx->GetZ();
479 name = vtx->GetName();
480 }
481 else if (obj->InheritsFrom("AliAODEvent")) {
482 AliAODEvent* aodevt = (AliAODEvent*) obj;
483 if (aodevt->GetNumberOfVertices() < 1)
484 return 0;
485 const AliAODVertex* vtx = aodevt->GetPrimaryVertex();
486 nContributors = vtx->GetNContributors();
487 zVertex = vtx->GetZ();
488 name = vtx->GetName();
489 }
490
491 // Reject if TPC-only vertex
492 if (name.CompareTo("TPCVertex")==0)
493 return kFALSE;
494
495 // Check # contributors and range...
496 if( nContributors < 1 || TMath::Abs(zVertex) > 10 ) {
497 return kFALSE;
498 }
499
500 return kTRUE;
501}
502
503
504//________________________________________________________________________
505Bool_t AliMuonEffMC::IsGoodMUONtrack(AliESDMuonTrack &track)
506{
507 // Applying track cuts for MUON tracks
508 if(!track.ContainTrackerData()) return kFALSE;
7e990b20 509
510 Double_t thetaTrackAbsEnd = TMath::ATan(track.GetRAtAbsorberEnd()/505.) * TMath::RadToDeg();
511 Double_t eta = track.Eta();
512
513 // Theta cut at absorber end
514 if(thetaTrackAbsEnd <= 2. || thetaTrackAbsEnd >= 10.) return kFALSE;
515 // Eta cut
516 if(eta <= -4. || eta >= -2.5) return kFALSE;
517 if(track.GetMatchTrigger() <= 0) return kFALSE;
518
519 return kTRUE;
520}
b673a083 521
522//________________________________________________________________________
523Bool_t AliMuonEffMC::IsGoodMUONtrack(AliAODTrack &track)
524{
525 if (!track.IsMuonTrack()) return kFALSE;
526
527 Double_t dThetaAbs = TMath::ATan(track.GetRAtAbsorberEnd()/505.)
528 * TMath::RadToDeg();
529 if ((dThetaAbs<2.) || (dThetaAbs>10.)) return kFALSE;
530
531 Double_t dEta = track.Eta();
532 if ((dEta<-4.) || (dEta>2.5)) return kFALSE;
533
534 if (track.GetMatchTrigger()<0.5) return kFALSE;
535
536 return kTRUE;
537}
15b0d1bb 538
539//________________________________________________________________________
540void AliMuonEffMC::MDProcess(Int_t motherpdg, Double_t mcpt, Double_t mcphi, Double_t mceta, Double_t trackpt, Double_t trackphi, Double_t tracketa, Double_t motherpt, Double_t motherphi, Double_t mothereta, Double_t dcavalue)
541{
542 if(motherpdg==411 || motherpdg==413 || motherpdg==421 || motherpdg==423 || motherpdg==431 || motherpdg==433 || motherpdg==10413 || motherpdg==10411 || motherpdg==10423 || motherpdg==10421 || motherpdg==10433 || motherpdg==10431 || motherpdg==20413 || motherpdg==415 || motherpdg==20423 || motherpdg==425 || motherpdg==20433 || motherpdg==435)
543 {
544 fHMuMotherGenPt[2]->Fill(mcpt, motherpt);
545 fHMuMotherRecPt[2]->Fill(trackpt, motherpt);
546 fHMuMotherGenPhi[2]->Fill(mcphi, motherphi);
547 fHMuMotherRecPhi[2]->Fill(trackphi, motherphi);
548 fHMuMotherGenEta[2]->Fill(mceta, mothereta);
549 fHMuMotherRecEta[2]->Fill(tracketa, mothereta);
550
551 // generated
552 if((0.5 <= mcpt) && (mcpt < 1.0))
553 {
554 fHMuMohterPhiDifGen[2][0]->Fill(deltaphi(motherphi-mcphi));
555 fHMuMohterEtaDifGen[2][0]->Fill(mothereta-mceta);
556 }
557 else if((1.0 <= mcpt) && (mcpt < 2.0))
558 {
559 fHMuMohterPhiDifGen[2][1]->Fill(deltaphi(motherphi-mcphi));
560 fHMuMohterEtaDifGen[2][1]->Fill(mothereta-mceta);
561 }
562 else if((2.0 <= mcpt) && (mcpt < 4.0))
563 {
564 fHMuMohterPhiDifGen[2][2]->Fill(deltaphi(motherphi-mcphi));
565 fHMuMohterEtaDifGen[2][2]->Fill(mothereta-mceta);
566 }
567 // reconstructed
568 if((0.5 <= trackpt) && (trackpt < 1.0))
569 {
570 fHMuMohterPhiDifRec[2][0]->Fill(deltaphi(motherphi-trackphi));
571 fHMuMohterEtaDifRec[2][0]->Fill(mothereta-tracketa);
572 }
573 else if((1.0 <= trackpt) && (trackpt < 2.0))
574 {
575 fHMuMohterPhiDifRec[2][1]->Fill(deltaphi(motherphi-trackphi));
576 fHMuMohterEtaDifRec[2][1]->Fill(mothereta-tracketa);
577 }
578 else if((2.0 <= trackpt) && (trackpt < 4.0))
579 {
580 fHMuMohterPhiDifRec[2][2]->Fill(deltaphi(motherphi-trackphi));
581 fHMuMohterEtaDifRec[2][2]->Fill(mothereta-tracketa);
582 }
583
584 // DCA
585 if(fESD) fHMuDCA[2]->Fill(dcavalue);
586 }
587
588 else if(motherpdg==211)
589 {
590 fHMuMotherGenPt[0]->Fill(mcpt, motherpt);
591 fHMuMotherRecPt[0]->Fill(trackpt, motherpt);
592 fHMuMotherGenPhi[0]->Fill(mcphi, motherphi);
593 fHMuMotherRecPhi[0]->Fill(trackphi, motherphi);
594 fHMuMotherGenEta[0]->Fill(mceta, mothereta);
595 fHMuMotherRecEta[0]->Fill(tracketa, mothereta);
596 // generated
597 if((0.5 <= mcpt) && (mcpt < 1.0))
598 {
599 fHMuMohterPhiDifGen[0][0]->Fill(deltaphi(motherphi-mcphi));
600 fHMuMohterEtaDifGen[0][0]->Fill(mothereta-mceta);
601 }
602 else if((1.0 <= mcpt) && (mcpt < 2.0))
603 {
604 fHMuMohterPhiDifGen[0][1]->Fill(deltaphi(motherphi-mcphi));
605 fHMuMohterEtaDifGen[0][1]->Fill(mothereta-mceta);
606 }
607 else if((2.0 <= mcpt) && (mcpt < 4.0))
608 {
609 fHMuMohterPhiDifGen[0][2]->Fill(deltaphi(motherphi-mcphi));
610 fHMuMohterEtaDifGen[0][2]->Fill(mothereta-mceta);
611 }
612 // reconstructed
613 if((0.5 <= trackpt) && (trackpt < 1.0))
614 {
615 fHMuMohterPhiDifRec[0][0]->Fill(deltaphi(motherphi-trackphi));
616 fHMuMohterEtaDifRec[0][0]->Fill(mothereta-tracketa);
617 }
618 else if((1.0 <= trackpt) && (trackpt < 2.0))
619 {
620 fHMuMohterPhiDifRec[0][1]->Fill(deltaphi(motherphi-trackphi));
621 fHMuMohterEtaDifRec[0][1]->Fill(mothereta-tracketa);
622 }
623 else if((2.0 <= trackpt) && (trackpt < 4.0))
624 {
625 fHMuMohterPhiDifRec[0][2]->Fill(deltaphi(motherphi-trackphi));
626 fHMuMohterEtaDifRec[0][2]->Fill(mothereta-tracketa);
627 }
628 if(fESD) fHMuDCA[0]->Fill(dcavalue);
629 }
630 else if(motherpdg==321)
631 {
632 fHMuMotherGenPt[1]->Fill(mcpt, motherpt);
633 fHMuMotherRecPt[1]->Fill(trackpt, motherpt);
634 fHMuMotherGenPhi[1]->Fill(mcphi, motherphi);
635 fHMuMotherRecPhi[1]->Fill(trackphi, motherphi);
636 fHMuMotherGenEta[1]->Fill(mceta, mothereta);
637 fHMuMotherRecEta[1]->Fill(tracketa, mothereta);
638 // generated
639 if((0.5 <= mcpt) && (mcpt < 1.0))
640 {
641 fHMuMohterPhiDifGen[1][0]->Fill(deltaphi(motherphi-mcphi));
642 fHMuMohterEtaDifGen[1][0]->Fill(mothereta-mceta);
643 }
644 else if((1.0 <= mcpt) && (mcpt < 2.0))
645 {
646 fHMuMohterPhiDifGen[1][1]->Fill(deltaphi(motherphi-mcphi));
647 fHMuMohterEtaDifGen[1][1]->Fill(mothereta-mceta);
648 }
649 else if((2.0 <= mcpt) && (mcpt < 4.0))
650 {
651 fHMuMohterPhiDifGen[1][2]->Fill(deltaphi(motherphi-mcphi));
652 fHMuMohterEtaDifGen[1][2]->Fill(mothereta-mceta);
653 }
654 // reconstructed
655 if((0.5 <= trackpt) && (trackpt < 1.0))
656 {
657 fHMuMohterPhiDifRec[1][0]->Fill(deltaphi(motherphi-trackphi));
658 fHMuMohterEtaDifRec[1][0]->Fill(mothereta-tracketa);
659 }
660 else if((1.0 <= trackpt) && (trackpt < 2.0))
661 {
662 fHMuMohterPhiDifRec[1][1]->Fill(deltaphi(motherphi-trackphi));
663 fHMuMohterEtaDifRec[1][1]->Fill(mothereta-tracketa);
664 }
665 else if((2.0 <= trackpt) && (trackpt < 4.0))
666 {
667 fHMuMohterPhiDifRec[1][2]->Fill(deltaphi(motherphi-trackphi));
668 fHMuMohterEtaDifRec[1][2]->Fill(mothereta-tracketa);
669 }
670 if(fESD) fHMuDCA[1]->Fill(dcavalue);
671 }
672 else
673 {
674 fHMuMotherGenPt[3]->Fill(mcpt, motherpt);
675 fHMuMotherRecPt[3]->Fill(trackpt, motherpt);
676 fHMuMotherGenPhi[3]->Fill(mcphi, motherphi);
677 fHMuMotherRecPhi[3]->Fill(trackphi, motherphi);
678 fHMuMotherGenEta[3]->Fill(mceta, mothereta);
679 fHMuMotherRecEta[3]->Fill(tracketa, mothereta);
680 // generated
681 if((0.5 <= mcpt) && (mcpt < 1.0))
682 {
683 fHMuMohterPhiDifGen[3][0]->Fill(deltaphi(motherphi-mcphi));
684 fHMuMohterEtaDifGen[3][0]->Fill(mothereta-mceta);
685 }
686 else if((1.0 <= mcpt) && (mcpt < 2.0))
687 {
688 fHMuMohterPhiDifGen[3][1]->Fill(deltaphi(motherphi-mcphi));
689 fHMuMohterEtaDifGen[3][1]->Fill(mothereta-mceta);
690 }
691 else if((2.0 <= mcpt) && (mcpt < 4.0))
692 {
693 fHMuMohterPhiDifGen[3][2]->Fill(deltaphi(motherphi-mcphi));
694 fHMuMohterEtaDifGen[3][2]->Fill(mothereta-mceta);
695 }
696 // reconstructed
697 if((0.5 <= trackpt) && (trackpt < 1.0))
698 {
699 fHMuMohterPhiDifRec[3][0]->Fill(deltaphi(motherphi-trackphi));
700 fHMuMohterEtaDifRec[3][0]->Fill(mothereta-tracketa);
701 }
702 else if((1.0 <= trackpt) && (trackpt < 2.0))
703 {
704 fHMuMohterPhiDifRec[3][1]->Fill(deltaphi(motherphi-trackphi));
705 fHMuMohterEtaDifRec[3][1]->Fill(mothereta-tracketa);
706 }
707 else if((2.0 <= trackpt) && (trackpt < 4.0))
708 {
709 fHMuMohterPhiDifRec[3][2]->Fill(deltaphi(motherphi-trackphi));
710 fHMuMohterEtaDifRec[3][2]->Fill(mothereta-tracketa);
711 }
712 if(fESD) fHMuDCA[3]->Fill(dcavalue);
713 }
714}
715//________________________________________________________________________
716Double_t AliMuonEffMC::deltaphi(Double_t phi)
717{
718 if(phi < -1.0*TMath::Pi()) { return (phi + TMath::TwoPi()); }
719 else if(phi > TMath::Pi()) { return (phi - TMath::TwoPi()); }
720 else { return phi; }
721}