]>
Commit | Line | Data |
---|---|---|
3a219f60 | 1 | /************************************************************************** |
2 | * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. * | |
3 | * * | |
4 | * Author: The ALICE Off-line Project. * | |
5 | * Contributors are mentioned in the code where appropriate. * | |
6 | * * | |
7 | * Permission to use, copy, modify and distribute this software and its * | |
8 | * documentation strictly for non-commercial purposes is hereby granted * | |
9 | * without fee, provided that the above copyright notice appears in all * | |
10 | * copies and that both the copyright notice and this permission notice * | |
11 | * appear in the supporting documentation. The authors make no claims * | |
12 | * about the suitability of this software for any purpose. It is * | |
13 | * provided "as is" without express or implied warranty. * | |
14 | **************************************************************************/ | |
15 | ||
16 | ///////////////////////////////////////////////////////////// | |
17 | // | |
18 | // AliAnalysisTaskSE for the comparison of heavy flavor | |
19 | // decay candidates with the MC truth. | |
20 | // | |
21 | // Author: A.Dainese, andrea.dainese@lnl.infn.it | |
22 | ///////////////////////////////////////////////////////////// | |
23 | ||
24 | #include <TClonesArray.h> | |
25 | #include <TNtuple.h> | |
26 | #include <TList.h> | |
27 | #include <TH1F.h> | |
28 | ||
b557eb43 | 29 | #include "AliAnalysisManager.h" |
30 | #include "AliAODHandler.h" | |
3a219f60 | 31 | #include "AliAODEvent.h" |
32 | #include "AliAODVertex.h" | |
607719af | 33 | #include "AliESDtrack.h" |
3a219f60 | 34 | #include "AliAODTrack.h" |
35 | #include "AliAODMCHeader.h" | |
36 | #include "AliAODMCParticle.h" | |
37 | #include "AliAODRecoDecayHF2Prong.h" | |
c4cdf105 | 38 | #include "AliAODRecoDecayHF3Prong.h" |
39 | #include "AliAODRecoDecayHF4Prong.h" | |
b3999999 | 40 | #include "AliAODRecoCascadeHF.h" |
2ae9c7cb | 41 | #include "AliAnalysisVertexingHF.h" |
3a219f60 | 42 | #include "AliAnalysisTaskSE.h" |
43 | #include "AliAnalysisTaskSECompareHF.h" | |
44 | ||
45 | ClassImp(AliAnalysisTaskSECompareHF) | |
46 | ||
47 | ||
48 | //________________________________________________________________________ | |
49 | AliAnalysisTaskSECompareHF::AliAnalysisTaskSECompareHF(): | |
50 | AliAnalysisTaskSE(), | |
51 | fOutput(0), | |
2a041947 | 52 | fNtupleCmp(0), |
12bfc069 | 53 | fHistMass(0), |
a273622f | 54 | fHistNEvents(0), |
12bfc069 | 55 | fVHF(0) |
3a219f60 | 56 | { |
57 | // Default constructor | |
3a219f60 | 58 | |
9d6d35b0 | 59 | // NO DefineOutput() HERE (ONLY IN STANDARD CONSTRUCTOR) |
3a219f60 | 60 | } |
61 | ||
62 | //________________________________________________________________________ | |
63 | AliAnalysisTaskSECompareHF::AliAnalysisTaskSECompareHF(const char *name): | |
64 | AliAnalysisTaskSE(name), | |
65 | fOutput(0), | |
2a041947 | 66 | fNtupleCmp(0), |
12bfc069 | 67 | fHistMass(0), |
a273622f | 68 | fHistNEvents(0), |
12bfc069 | 69 | fVHF(0) |
3a219f60 | 70 | { |
9d6d35b0 | 71 | // Standard constructor |
3a219f60 | 72 | |
73 | // Output slot #1 writes into a TList container | |
74 | DefineOutput(1,TList::Class()); //My private output | |
10bdd1ae | 75 | // Output slot #2 writes into a TNtuple container |
76 | DefineOutput(2,TNtuple::Class()); //My private output | |
3a219f60 | 77 | } |
78 | ||
79 | //________________________________________________________________________ | |
80 | AliAnalysisTaskSECompareHF::~AliAnalysisTaskSECompareHF() | |
81 | { | |
82 | // Destructor | |
83 | if (fOutput) { | |
84 | delete fOutput; | |
85 | fOutput = 0; | |
86 | } | |
12bfc069 | 87 | if (fVHF) { |
88 | delete fVHF; | |
89 | fVHF = 0; | |
90 | } | |
3a219f60 | 91 | } |
92 | ||
93 | //________________________________________________________________________ | |
94 | void AliAnalysisTaskSECompareHF::Init() | |
95 | { | |
96 | // Initialization | |
97 | ||
98 | if(fDebug > 1) printf("AnalysisTaskSECompareHF::Init() \n"); | |
9d6d35b0 | 99 | |
12bfc069 | 100 | gROOT->LoadMacro("ConfigVertexingHF.C"); |
101 | ||
102 | fVHF = (AliAnalysisVertexingHF*)gROOT->ProcessLine("ConfigVertexingHF()"); | |
103 | fVHF->PrintStatus(); | |
9d6d35b0 | 104 | |
3a219f60 | 105 | return; |
106 | } | |
107 | ||
108 | //________________________________________________________________________ | |
109 | void AliAnalysisTaskSECompareHF::UserCreateOutputObjects() | |
110 | { | |
111 | // Create the output container | |
112 | // | |
113 | if(fDebug > 1) printf("AnalysisTaskSECompareHF::UserCreateOutputObjects() \n"); | |
114 | ||
115 | // Several histograms are more conveniently managed in a TList | |
116 | fOutput = new TList(); | |
117 | fOutput->SetOwner(); | |
118 | ||
119 | fHistMass = new TH1F("fHistMass", "D^{0} invariant mass; M [GeV]; Entries",200,1.765,1.965); | |
120 | fHistMass->Sumw2(); | |
121 | fHistMass->SetMinimum(0); | |
122 | fOutput->Add(fHistMass); | |
123 | ||
a273622f | 124 | fHistNEvents = new TH1F("fHistNEvents", "Number of processed events; ; Events",3,-1.5,1.5); |
125 | fHistNEvents->Sumw2(); | |
126 | fHistNEvents->SetMinimum(0); | |
127 | fOutput->Add(fHistNEvents); | |
128 | ||
2a5c77c6 | 129 | OpenFile(2); // 2 is the slot number of the ntuple |
607719af | 130 | fNtupleCmp = new TNtuple("fNtupleCmp","Charm comparison","pdg:nprongs:VxRec:VxTrue:ErrVx:VyRec:VyTrue:ErrVy:VzRec:VzTrue:ErrVz:Chi2toNDF:PtRec:Mrec:CPta:Prodd0"); |
3a219f60 | 131 | |
132 | return; | |
133 | } | |
134 | ||
135 | //________________________________________________________________________ | |
136 | void AliAnalysisTaskSECompareHF::UserExec(Option_t */*option*/) | |
137 | { | |
138 | // Execute analysis for current event: | |
139 | // heavy flavor candidates association to MC truth | |
9d6d35b0 | 140 | |
3a219f60 | 141 | |
142 | AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent()); | |
8931c313 | 143 | |
b557eb43 | 144 | TClonesArray *inputArrayVertices = 0; |
145 | TClonesArray *inputArrayD0toKpi = 0; | |
146 | TClonesArray *inputArrayDstar = 0; | |
147 | ||
148 | if(!aod && AODEvent() && IsStandardAOD()) { | |
149 | // In case there is an AOD handler writing a standard AOD, use the AOD | |
150 | // event in memory rather than the input (ESD) event. | |
151 | aod = dynamic_cast<AliAODEvent*> (AODEvent()); | |
152 | // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root) | |
153 | // have to taken from the AOD event hold by the AliAODExtension | |
154 | AliAODHandler* aodHandler = (AliAODHandler*) | |
155 | ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler()); | |
156 | if(aodHandler->GetExtensions()) { | |
157 | AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root"); | |
158 | AliAODEvent *aodFromExt = ext->GetAOD(); | |
159 | // load HF vertices | |
160 | inputArrayVertices = (TClonesArray*)aodFromExt->GetList()->FindObject("VerticesHF"); | |
161 | // load D0->Kpi candidates | |
162 | inputArrayD0toKpi = (TClonesArray*)aodFromExt->GetList()->FindObject("D0toKpi"); | |
163 | // load D*+ candidates | |
164 | inputArrayDstar = (TClonesArray*)aodFromExt->GetList()->FindObject("Dstar"); | |
165 | } | |
166 | } else { | |
167 | // load HF vertices | |
168 | inputArrayVertices = (TClonesArray*)aod->GetList()->FindObject("VerticesHF"); | |
169 | // load D0->Kpi candidates | |
170 | inputArrayD0toKpi = (TClonesArray*)aod->GetList()->FindObject("D0toKpi"); | |
171 | // load D*+ candidates | |
172 | inputArrayDstar = (TClonesArray*)aod->GetList()->FindObject("Dstar"); | |
173 | } | |
174 | ||
2a041947 | 175 | |
2a041947 | 176 | if(!inputArrayVertices) { |
177 | printf("AliAnalysisTaskSECompareHF::UserExec: Vertices branch not found!\n"); | |
178 | return; | |
179 | } | |
3a219f60 | 180 | if(!inputArrayD0toKpi) { |
181 | printf("AliAnalysisTaskSECompareHF::UserExec: D0toKpi branch not found!\n"); | |
182 | return; | |
183 | } | |
b3999999 | 184 | if(!inputArrayDstar) { |
185 | printf("AliAnalysisTaskSECompareHF::UserExec: Dstar branch not found!\n"); | |
186 | return; | |
187 | } | |
9d6d35b0 | 188 | |
b3999999 | 189 | |
b557eb43 | 190 | fHistNEvents->Fill(0); // count event |
191 | // Post the data already here | |
192 | PostData(1,fOutput); | |
193 | ||
3a219f60 | 194 | // AOD primary vertex |
195 | AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex(); | |
196 | //vtx1->Print(); | |
197 | ||
198 | // load MC particles | |
199 | TClonesArray *mcArray = | |
200 | (TClonesArray*)aod->GetList()->FindObject(AliAODMCParticle::StdBranchName()); | |
201 | if(!mcArray) { | |
202 | printf("AliAnalysisTaskSECompareHF::UserExec: MC particles branch not found!\n"); | |
203 | return; | |
204 | } | |
205 | ||
206 | // load MC header | |
207 | AliAODMCHeader *mcHeader = | |
208 | (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName()); | |
209 | if(!mcHeader) { | |
210 | printf("AliAnalysisTaskSECompareHF::UserExec: MC header branch not found!\n"); | |
211 | return; | |
212 | } | |
213 | ||
2a041947 | 214 | Int_t nprongs,lab,okD0,okD0bar,pdg; |
215 | Bool_t unsetvtx; | |
90426506 | 216 | Double_t invmass,posRec[3],posTrue[3],covRec[6],errx,erry,errz; |
952f1e33 | 217 | AliAODRecoDecayHF2Prong *d2=0; |
218 | AliAODRecoDecayHF3Prong *d3=0; | |
219 | AliAODRecoDecayHF4Prong *d4=0; | |
2a041947 | 220 | |
c8112d2f | 221 | Int_t pdgDgD0toKpi[2]={321,211}; |
222 | Int_t pdgDgDplustoKpipi[3]={321,211,211}; | |
223 | Int_t pdgDgD0toKpipipi[4]={321,211,211,211}; | |
224 | ||
2a041947 | 225 | // loop over vertices |
226 | Int_t nVertices = inputArrayVertices->GetEntriesFast(); | |
3b1598da | 227 | if(fDebug>1) printf("Number of vertices: %d\n",nVertices); |
2a041947 | 228 | |
229 | for (Int_t iVtx = 0; iVtx < nVertices; iVtx++) { | |
230 | AliAODVertex *vtx = (AliAODVertex*)inputArrayVertices->UncheckedAt(iVtx); | |
231 | ||
90426506 | 232 | vtx->GetXYZ(posRec); |
233 | vtx->GetCovarianceMatrix(covRec); | |
234 | errx=1.; erry=1.; errz=1.; | |
235 | if(covRec[0]>0) errx = TMath::Sqrt(covRec[0]); | |
236 | if(covRec[2]>0) erry = TMath::Sqrt(covRec[2]); | |
237 | if(covRec[5]>0) errz = TMath::Sqrt(covRec[5]); | |
238 | ||
239 | ||
2a041947 | 240 | // get parent AliAODRecoDecayHF |
241 | nprongs= vtx->GetNDaughters(); | |
607719af | 242 | // check that the daughters have kITSrefit and kTPCrefit |
243 | Bool_t allDgOK=kTRUE; | |
244 | for(Int_t idg=0; idg<nprongs; idg++) { | |
245 | AliAODTrack *track = (AliAODTrack*)vtx->GetDaughter(idg); | |
246 | if(!(track->GetStatus()&AliESDtrack::kITSrefit)) allDgOK=kFALSE; | |
247 | if(!(track->GetStatus()&AliESDtrack::kTPCrefit)) allDgOK=kFALSE; | |
248 | } | |
249 | if(!allDgOK) continue; | |
250 | ||
2a041947 | 251 | |
252 | switch(nprongs) { | |
253 | case 2: // look for D0->Kpi | |
952f1e33 | 254 | d2 = (AliAODRecoDecayHF2Prong*)vtx->GetParent(); |
27917274 | 255 | if(d2->IsLikeSign()) continue; |
a273622f | 256 | if(d2->Charge() != 0) continue; // these are D* |
c8112d2f | 257 | lab = d2->MatchToMC(421,mcArray,2,pdgDgD0toKpi); |
2a041947 | 258 | if(lab>=0) { |
259 | unsetvtx=kFALSE; | |
260 | if(!d2->GetOwnPrimaryVtx()) { | |
261 | d2->SetOwnPrimaryVtx(vtx1); // needed to compute all variables | |
262 | unsetvtx=kTRUE; | |
263 | } | |
264 | okD0=0; okD0bar=0; | |
607719af | 265 | if(d2->SelectD0(fVHF->GetD0toKpiCuts(),okD0,okD0bar)) { // beware: cuts may bias the resolution! |
266 | AliAODMCParticle *dMC = (AliAODMCParticle*)mcArray->At(lab); | |
267 | pdg = dMC->GetPdgCode(); | |
268 | invmass = (pdg==421 ? d2->InvMassD0() : d2->InvMassD0bar()); | |
607719af | 269 | // get a daughter for true pos of decay vertex |
270 | AliAODMCParticle *dg0MC = (AliAODMCParticle*)mcArray->At(dMC->GetDaughter(0)); | |
271 | dg0MC->XvYvZv(posTrue); | |
caf4ca8b | 272 | fHistMass->Fill(invmass); |
273 | // Post the data already here | |
274 | PostData(1,fOutput); | |
607719af | 275 | Float_t tmp[16]={(Float_t)pdg,(Float_t)nprongs, |
276 | (Float_t)posRec[0],(Float_t)posTrue[0],(Float_t)errx, | |
277 | (Float_t)posRec[1],(Float_t)posTrue[1],(Float_t)erry, | |
278 | (Float_t)posRec[2],(Float_t)posTrue[2],(Float_t)errz, | |
279 | (Float_t)d2->GetReducedChi2(),(Float_t)d2->Pt(),(Float_t)invmass, | |
280 | (Float_t)d2->CosPointingAngle(),(Float_t)d2->Prodd0d0()}; | |
281 | fNtupleCmp->Fill(tmp); | |
282 | PostData(2,fNtupleCmp); | |
283 | } | |
2a041947 | 284 | if(unsetvtx) d2->UnsetOwnPrimaryVtx(); |
3a219f60 | 285 | } |
2a041947 | 286 | break; |
287 | case 3: // look for D+ | |
952f1e33 | 288 | d3 = (AliAODRecoDecayHF3Prong*)vtx->GetParent(); |
27917274 | 289 | if(d3->IsLikeSign()) continue; |
c8112d2f | 290 | lab = d3->MatchToMC(411,mcArray,3,pdgDgDplustoKpipi); |
2a041947 | 291 | if(lab>=0) { |
292 | unsetvtx=kFALSE; | |
293 | if(!d3->GetOwnPrimaryVtx()) { | |
294 | d3->SetOwnPrimaryVtx(vtx1); // needed to compute all variables | |
295 | unsetvtx=kTRUE; | |
296 | } | |
90426506 | 297 | //if(d3->SelectDplus(fVHF->GetDplusCuts())) { |
298 | AliAODMCParticle *dMC = (AliAODMCParticle*)mcArray->At(lab); | |
299 | pdg = dMC->GetPdgCode(); | |
300 | invmass = d3->InvMassDplus(); | |
301 | // get a daughter for true pos of decay vertex | |
302 | AliAODMCParticle *dg0MC = (AliAODMCParticle*)mcArray->At(dMC->GetDaughter(0)); | |
303 | dg0MC->XvYvZv(posTrue); | |
607719af | 304 | Float_t tmp[16]={(Float_t)pdg,(Float_t)nprongs, |
305 | (Float_t)posRec[0],(Float_t)posTrue[0],(Float_t)errx, | |
306 | (Float_t)posRec[1],(Float_t)posTrue[1],(Float_t)erry, | |
307 | (Float_t)posRec[2],(Float_t)posTrue[2],(Float_t)errz, | |
308 | (Float_t)d3->GetReducedChi2(),(Float_t)d3->Pt(),(Float_t)invmass, | |
309 | (Float_t)d3->CosPointingAngle(),(Float_t)(d3->Getd0Prong(0)*d3->Getd0Prong(1)*d3->Getd0Prong(2))}; | |
310 | fNtupleCmp->Fill(tmp); | |
a7c1934a | 311 | PostData(2,fNtupleCmp); |
90426506 | 312 | //} |
2a041947 | 313 | if(unsetvtx) d3->UnsetOwnPrimaryVtx(); |
314 | } | |
315 | break; | |
316 | case 4: // look for D0->Kpipipi | |
952f1e33 | 317 | d4 = (AliAODRecoDecayHF4Prong*)vtx->GetParent(); |
27917274 | 318 | if(d4->IsLikeSign()) continue; |
c8112d2f | 319 | lab = d4->MatchToMC(421,mcArray,4,pdgDgD0toKpipipi); |
2a041947 | 320 | if(lab>=0) { |
321 | unsetvtx=kFALSE; | |
322 | if(!d4->GetOwnPrimaryVtx()) { | |
323 | d4->SetOwnPrimaryVtx(vtx1); // needed to compute all variables | |
324 | unsetvtx=kTRUE; | |
325 | } | |
326 | okD0=0; okD0bar=0; | |
90426506 | 327 | //if(d4->SelectD0(fVHF->GetD0to4ProngsCuts(),okD0,okD0bar)) { |
328 | AliAODMCParticle *dMC = (AliAODMCParticle*)mcArray->At(lab); | |
329 | pdg = dMC->GetPdgCode(); | |
330 | //invmass = (pdg==421 ? d->InvMassD0() : d->InvMassD0bar()); | |
331 | invmass = 10.; // dummy | |
332 | // get a daughter for true pos of decay vertex | |
333 | AliAODMCParticle *dg0MC = (AliAODMCParticle*)mcArray->At(dMC->GetDaughter(0)); | |
334 | dg0MC->XvYvZv(posTrue); | |
607719af | 335 | Float_t tmp[16]={(Float_t)pdg,(Float_t)nprongs, |
336 | (Float_t)posRec[0],(Float_t)posTrue[0],(Float_t)errx, | |
337 | (Float_t)posRec[1],(Float_t)posTrue[1],(Float_t)erry, | |
338 | (Float_t)posRec[2],(Float_t)posTrue[2],(Float_t)errz, | |
339 | (Float_t)d4->GetReducedChi2(),(Float_t)d4->Pt(),(Float_t)invmass, | |
340 | (Float_t)d4->CosPointingAngle(),(Float_t)(d4->Getd0Prong(0)*d4->Getd0Prong(1)*d4->Getd0Prong(2)*d4->Getd0Prong(3))}; | |
341 | fNtupleCmp->Fill(tmp); | |
a7c1934a | 342 | PostData(2,fNtupleCmp); |
90426506 | 343 | //} |
2a041947 | 344 | if(unsetvtx) d4->UnsetOwnPrimaryVtx(); |
345 | } | |
346 | break; | |
3a219f60 | 347 | } |
3a219f60 | 348 | |
2a041947 | 349 | } // end loop on vertices |
3a219f60 | 350 | |
9d6d35b0 | 351 | |
b3999999 | 352 | // loop over D*+ candidates |
4dd17b11 | 353 | /* |
b3999999 | 354 | for (Int_t iDstar = 0; iDstar < inputArrayDstar->GetEntries(); iDstar++) { |
355 | AliAODRecoCascadeHF *c = (AliAODRecoCascadeHF*)inputArrayDstar->UncheckedAt(iDstar); | |
356 | Int_t labDstar = c->MatchToMC(413,421,mcArray); | |
357 | if(labDstar>=0) printf("GOOD MATCH FOR D*+\n"); | |
358 | } | |
4dd17b11 | 359 | */ |
2a041947 | 360 | |
4dd17b11 | 361 | // Post the data already here |
362 | PostData(1,fOutput); | |
2a041947 | 363 | PostData(2,fNtupleCmp); |
b3999999 | 364 | |
3a219f60 | 365 | return; |
366 | } | |
3a219f60 | 367 | //________________________________________________________________________ |
368 | void AliAnalysisTaskSECompareHF::Terminate(Option_t */*option*/) | |
369 | { | |
370 | // Terminate analysis | |
371 | // | |
372 | if(fDebug > 1) printf("AnalysisTaskSECompareHF: Terminate() \n"); | |
373 | ||
374 | fOutput = dynamic_cast<TList*> (GetOutputData(1)); | |
375 | if (!fOutput) { | |
376 | printf("ERROR: fOutput not available\n"); | |
377 | return; | |
378 | } | |
379 | ||
3a219f60 | 380 | fHistMass = dynamic_cast<TH1F*>(fOutput->FindObject("fHistMass")); |
a273622f | 381 | fHistNEvents = dynamic_cast<TH1F*>(fOutput->FindObject("fHistNEvents")); |
3a219f60 | 382 | |
27917274 | 383 | //fNtupleCmp = dynamic_cast<TNtuple*> (GetOutputData(2)); |
10bdd1ae | 384 | |
3a219f60 | 385 | return; |
386 | } | |
387 |