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