]>
Commit | Line | Data |
---|---|---|
7b54a483 | 1 | /************************************************************************** |
2 | * Copyright(c) 1998-2009, 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 | // Class AliAnalysisTaskSEBtoJPSItoEle | |
17 | // AliAnalysisTaskSE for the reconstruction | |
18 | // of heavy-flavour decay candidates | |
19 | // Author: C.Di Giglio, carmelo.digiglio@ba.infn.it | |
20 | //************************************************************************* | |
21 | class AliAnalysisTaskSE; | |
22 | class AliESDEvent; | |
23 | class AliVEvent; | |
24 | class iostream; | |
25 | class TSystem; | |
26 | class TROOT; | |
27 | #include <TFile.h> | |
28 | #include <TClonesArray.h> | |
29 | #include <TDatabasePDG.h> | |
30 | ||
b557eb43 | 31 | #include "AliAnalysisManager.h" |
32 | #include "AliAODHandler.h" | |
7b54a483 | 33 | #include "AliAODEvent.h" |
34 | #include "AliAODMCParticle.h" | |
35 | #include "AliAODMCHeader.h" | |
36 | #include "AliAODRecoDecayHF2Prong.h" | |
37 | ||
38 | #include "AliAnalysisBtoJPSItoEle.h" | |
39 | #include "AliAnalysisTaskSEBtoJPSItoEle.h" | |
40 | ||
41 | ClassImp(AliAnalysisTaskSEBtoJPSItoEle) | |
42 | ||
43 | //______________________________________________________________________________ | |
44 | AliAnalysisTaskSEBtoJPSItoEle::AliAnalysisTaskSEBtoJPSItoEle() : | |
45 | AliAnalysisTaskSE(), | |
46 | fOutput(0), | |
47 | fCdfFit(0), | |
48 | fNtupleJPSI(0), | |
49 | fhDecayTimeMCjpsifromB(0), | |
50 | fhDecayTime(0), | |
51 | fhInvMass(0), | |
52 | fhD0(0), | |
53 | fhD0D0(0), | |
54 | fhCosThetaStar(0), | |
55 | fhCosThetaPointing(0), | |
56 | fhCtsVsD0D0(0), | |
57 | fOkAODMC(kFALSE), | |
58 | fNameMCfile(""), | |
59 | fOkMinimize(kFALSE) | |
60 | { | |
61 | // default constructor | |
7b54a483 | 62 | } |
63 | //_________________________________________________________________________________________________ | |
64 | AliAnalysisTaskSEBtoJPSItoEle::AliAnalysisTaskSEBtoJPSItoEle(const char* name) : | |
65 | AliAnalysisTaskSE(name), | |
66 | fOutput(0), | |
67 | fCdfFit(0), | |
68 | fNtupleJPSI(0), | |
69 | fhDecayTimeMCjpsifromB(0), | |
70 | fhDecayTime(0), | |
71 | fhInvMass(0), | |
72 | fhD0(0), | |
73 | fhD0D0(0), | |
74 | fhCosThetaStar(0), | |
75 | fhCosThetaPointing(0), | |
76 | fhCtsVsD0D0(0), | |
77 | fOkAODMC(kFALSE), | |
78 | fNameMCfile(""), | |
79 | fOkMinimize(kFALSE) | |
80 | { | |
9d6d35b0 | 81 | // standard constructor |
7b54a483 | 82 | |
83 | // Output slot #1 writes into a TList container | |
84 | DefineOutput(1,TList::Class()); //My private output | |
85 | } | |
86 | //_________________________________________________________________________________________________ | |
87 | AliAnalysisTaskSEBtoJPSItoEle::~AliAnalysisTaskSEBtoJPSItoEle() | |
88 | { | |
89 | // destructor | |
90 | ||
91 | if (fOutput) { | |
92 | delete fOutput; | |
93 | fOutput = 0; | |
94 | } | |
95 | } | |
96 | //_________________________________________________________________________________________________ | |
97 | void AliAnalysisTaskSEBtoJPSItoEle::Init() | |
98 | { | |
99 | // Initialization | |
7b54a483 | 100 | |
101 | if(fDebug > 1) printf("AliAnalysisTaskSEBtoJPSItoEle::Init() \n"); | |
102 | ||
103 | return; | |
104 | ||
105 | } | |
106 | //_________________________________________________________________________________________________ | |
107 | void AliAnalysisTaskSEBtoJPSItoEle::UserCreateOutputObjects() | |
108 | { | |
109 | // Create the output container | |
110 | ||
111 | if(fDebug > 1) printf("AliAnalysisTaskSEBtoJPSItoEle::UserCreateOutputObjects() \n"); | |
112 | ||
113 | // Several histograms are more conveniently managed in a TList | |
114 | fOutput = new TList(); | |
115 | fOutput->SetOwner(); | |
116 | ||
117 | if(fOkAODMC){ | |
118 | fhDecayTimeMCjpsifromB = new TH1F("fhDecayTimeMCjpsifromB", "Secondary J/#Psi Monte Carlo pseudo proper decay time; X [#mu m]; Entries",100,-4000.,4000.); | |
119 | fhDecayTimeMCjpsifromB->Sumw2(); | |
120 | fhDecayTimeMCjpsifromB->SetMinimum(0); | |
121 | fOutput->Add(fhDecayTimeMCjpsifromB); | |
122 | } | |
123 | ||
f9fd1412 | 124 | fhDecayTime = new TH1F("fhDecayTime", "J/#Psi pseudo proper decay time; X [#mu m]; Entries",200,-2000.,2000.); |
7b54a483 | 125 | fhDecayTime->Sumw2(); |
126 | fhDecayTime->SetMinimum(0); | |
127 | fOutput->Add(fhDecayTime); | |
128 | ||
129 | fhInvMass = new TH1F("fhInvMass", "J/#Psi invariant mass; M [GeV]; Entries",100,0.,3.2); | |
130 | fhInvMass->Sumw2(); | |
131 | fhInvMass->SetMinimum(0); | |
132 | fOutput->Add(fhInvMass); | |
133 | ||
134 | fhD0 = new TH1F("fhD0", "Impact parameter; D0 [#mu m]; Entries",100,-5000.,5000.); | |
135 | fhD0->Sumw2(); | |
136 | fhD0->SetMinimum(0); | |
137 | fOutput->Add(fhD0); | |
138 | ||
139 | fhD0D0 = new TH1F("fhD0D0", "Product of impact parameters; D0D0 [#mu m^2]; Entries",100,-100000.,100000.); | |
140 | fhD0D0->Sumw2(); | |
141 | fhD0D0->SetMinimum(0); | |
142 | fOutput->Add(fhD0D0); | |
143 | ||
144 | fhCosThetaStar = new TH1F("fhCosThetaStar", "Cosine of decay angle; Cosine Theta Star; Entries",50,-1.2,1.2); | |
145 | fhCosThetaStar->Sumw2(); | |
146 | fhCosThetaStar->SetMinimum(0); | |
147 | fOutput->Add(fhCosThetaStar); | |
148 | ||
149 | fhCosThetaPointing = new TH1F("fhCosThetaPointing", "Cosine of pointing angle; Cosine Pointing Angle; Entries",100,-1,1); | |
150 | fhCosThetaPointing->Sumw2(); | |
151 | fhCosThetaPointing->SetMinimum(0); | |
152 | fOutput->Add(fhCosThetaPointing); | |
153 | ||
154 | fhCtsVsD0D0 = new TH2F("fhCtsVsD0D0", "Cosine theta star Vs. D0; Cosine theta star; D0D0",50,-1,1,100,-100000,100000); | |
155 | fhCtsVsD0D0->Sumw2(); | |
156 | fhCtsVsD0D0->SetMinimum(0); | |
157 | fOutput->Add(fhCtsVsD0D0); | |
158 | ||
159 | fNtupleJPSI = new TNtuple("fNtupleJPSI","J/#Psi pseudo-proper decay time & invariant mass","Xdecaytime:Mass"); | |
160 | fOutput->Add(fNtupleJPSI); | |
161 | ||
162 | return; | |
163 | } | |
164 | //_________________________________________________________________________________________________ | |
165 | void AliAnalysisTaskSEBtoJPSItoEle::UserExec(Option_t */*option*/) | |
166 | { | |
167 | // Execute analysis for current event: | |
168 | ||
169 | AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent()); | |
170 | ||
b557eb43 | 171 | TClonesArray *inputArrayJPSItoEle = 0; |
172 | ||
173 | if(!aod && AODEvent() && IsStandardAOD()) { | |
174 | // In case there is an AOD handler writing a standard AOD, use the AOD | |
175 | // event in memory rather than the input (ESD) event. | |
176 | aod = dynamic_cast<AliAODEvent*> (AODEvent()); | |
177 | // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root) | |
178 | // have to taken from the AOD event hold by the AliAODExtension | |
179 | AliAODHandler* aodHandler = (AliAODHandler*) | |
180 | ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler()); | |
181 | if(aodHandler->GetExtensions()) { | |
182 | AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root"); | |
183 | AliAODEvent *aodFromExt = ext->GetAOD(); | |
184 | // load Jpsi candidates | |
185 | inputArrayJPSItoEle=(TClonesArray*)aodFromExt->GetList()->FindObject("JPSItoEle"); | |
186 | } | |
187 | } else { | |
188 | // load Jpsi candidates | |
189 | inputArrayJPSItoEle=(TClonesArray*)aod->GetList()->FindObject("JPSItoEle"); | |
190 | } | |
8931c313 | 191 | |
7b54a483 | 192 | if(!inputArrayJPSItoEle) { |
193 | printf("AliAnalysisTaskSECompareHF::UserExec: JPSItoEle branch not found!\n"); | |
194 | return; | |
195 | } | |
196 | ||
7c23877d | 197 | // fix for temporary bug in ESDfilter |
198 | // the AODs with null vertex pointer didn't pass the PhysSel | |
5806c290 | 199 | if(!aod->GetPrimaryVertex() || TMath::Abs(aod->GetMagneticField())<0.001) return; |
7c23877d | 200 | |
201 | ||
b557eb43 | 202 | // AOD primary vertex |
203 | AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex(); | |
204 | ||
f9fd1412 | 205 | // load MC particles |
206 | TClonesArray* mcArray = | |
207 | dynamic_cast<TClonesArray*>(aod->FindListObject(AliAODMCParticle::StdBranchName())); | |
208 | if (!mcArray) AliError("Could not find Monte-Carlo in AOD"); | |
209 | AliInfo(Form("MC particles found in mcArray ---> %d ",mcArray->GetEntriesFast())); | |
210 | ||
7b54a483 | 211 | // Read AOD Monte Carlo information |
212 | if(fOkAODMC) ReadAODMCInfo(aod,inputArrayJPSItoEle); | |
213 | ||
214 | // loop over J/Psi->ee candidates | |
215 | Int_t nInJPSItoEle = inputArrayJPSItoEle->GetEntriesFast(); | |
f9fd1412 | 216 | printf("Number of JPSI->e+e- candidates: %d\n",nInJPSItoEle); |
7b54a483 | 217 | |
218 | for (Int_t iJPSItoEle = 0; iJPSItoEle < nInJPSItoEle; iJPSItoEle++) { | |
219 | AliAODRecoDecayHF2Prong *d = (AliAODRecoDecayHF2Prong*)inputArrayJPSItoEle->UncheckedAt(iJPSItoEle); | |
f9fd1412 | 220 | |
221 | Int_t mcLabel = d->MatchToMC(443,mcArray) ; // select only reco JPSI that are true JPSI | |
222 | ||
7b54a483 | 223 | //Secondary vertex |
224 | Double_t vtxSec[3] = {0.,0.,0.}; | |
225 | Double_t vtxPrim[3] = {0.,0.,0.}; | |
226 | d->GetSecondaryVtx(vtxSec); | |
227 | Bool_t unsetvtx=kFALSE; | |
228 | if(!d->GetOwnPrimaryVtx()) { | |
229 | d->SetOwnPrimaryVtx(vtx1); // needed to compute all variables | |
230 | unsetvtx=kTRUE; | |
231 | } | |
232 | //Int_t okJPSI=0; | |
233 | // here analyze only if cuts are passed | |
f9fd1412 | 234 | |
7b54a483 | 235 | if(d->Pt() < fPtCuts[1] && d->Pt() > fPtCuts[0]){ // apply pt cut only |
f9fd1412 | 236 | if (mcLabel == -1) |
237 | { | |
238 | AliDebug(2,"No MC particle found"); | |
239 | } | |
240 | else { | |
241 | ||
7b54a483 | 242 | fhInvMass->Fill(d->InvMassJPSIee()); |
243 | fhD0->Fill(10000*d->ImpParXY()); | |
244 | fhD0D0->Fill(1e8*d->Prodd0d0()); | |
245 | fhCosThetaStar->Fill(d->CosThetaStarJPSI()); | |
246 | fhCtsVsD0D0->Fill(d->CosThetaStarJPSI(),1e8*d->Prodd0d0()); | |
247 | fhCosThetaPointing->Fill(d->CosPointingAngle()); | |
248 | ||
249 | // compute X variable | |
250 | AliAODVertex* primVertex = d->GetOwnPrimaryVtx(); | |
251 | vtxPrim[0] = primVertex->GetX(); | |
252 | vtxPrim[1] = primVertex->GetY(); | |
253 | vtxPrim[2] = primVertex->GetZ(); | |
254 | Double_t lxy = ((vtxSec[0]-vtxPrim[0])*(d->Px()) + (vtxSec[1]-vtxPrim[1])*(d->Py()))/d->Pt(); | |
255 | Double_t pseudoProperDecayTime = lxy*(TDatabasePDG::Instance()->GetParticle(443)->Mass())/d->Pt(); | |
256 | ||
257 | fhDecayTime->Fill(10000*pseudoProperDecayTime); | |
258 | // Post the data already here | |
259 | PostData(1,fOutput); | |
260 | ||
261 | fNtupleJPSI->Fill(pseudoProperDecayTime,d->InvMassJPSIee()); // fill N-tuple with X,M values | |
262 | ||
f9fd1412 | 263 | }// |
264 | ||
7b54a483 | 265 | } // end of JPSItoEle candidates selection according to cuts |
266 | ||
267 | if(unsetvtx) d->UnsetOwnPrimaryVtx(); | |
268 | ||
269 | }// end loop on JPSI to ele candidates | |
270 | ||
271 | } | |
272 | //_________________________________________________________________________________________________ | |
273 | void AliAnalysisTaskSEBtoJPSItoEle::Terminate(Option_t */*option*/) | |
274 | { | |
275 | // | |
276 | // Terminate analysis | |
277 | // | |
278 | ||
279 | if(fDebug > 1) printf("AliAnalysisTaskSEBtoJPSItoEle: Terminate() \n"); | |
280 | ||
281 | fOutput = dynamic_cast<TList*> (GetOutputData(1)); | |
282 | if (!fOutput) { | |
283 | printf("ERROR: fOutput not available\n"); | |
284 | return; | |
285 | } | |
286 | ||
287 | if(fOkAODMC) fhDecayTimeMCjpsifromB = dynamic_cast<TH1F*>(fOutput->FindObject("fhDecayTimeMCjpsifromB")); | |
288 | fhDecayTime = dynamic_cast<TH1F*>(fOutput->FindObject("fhDecayTime")); | |
289 | fhInvMass = dynamic_cast<TH1F*>(fOutput->FindObject("fhInvMass")); | |
290 | fhD0 = dynamic_cast<TH1F*>(fOutput->FindObject("fhD0")); | |
291 | fhD0D0 = dynamic_cast<TH1F*>(fOutput->FindObject("fhD0D0")); | |
292 | fhCosThetaStar = dynamic_cast<TH1F*>(fOutput->FindObject("fhCosThetaStar")); | |
293 | fhCosThetaPointing = dynamic_cast<TH1F*>(fOutput->FindObject("fhCosThetaPointing")); | |
294 | fhCtsVsD0D0 = dynamic_cast<TH2F*>(fOutput->FindObject("fhCtsVsD0D0")); | |
295 | fNtupleJPSI = dynamic_cast<TNtuple*>(fOutput->FindObject("fNtupleJPSI")); | |
296 | ||
297 | if(fOkMinimize){ | |
298 | Double_t* pseudoProperValues=0x0; | |
299 | Double_t* invMassValues=0x0; | |
300 | Int_t ncand=0; | |
301 | fCdfFit = new AliAnalysisBtoJPSItoEle(); | |
302 | ||
303 | printf("+++\n+++ AliAnalysisBtoJPSItoEle object instantiated ---> OK\n+++\n"); | |
304 | fCdfFit->ReadCandidates(fNtupleJPSI,pseudoProperValues,invMassValues,ncand); | |
305 | printf("+++\n+++ X and M vectors filled starting from N-tuple ---> OK\n+++\n"); | |
306 | printf("+++\n+++ Number of candidates ---> %d J/psi \n+++\n",ncand); | |
307 | fCdfFit->SetPtBin(0); | |
308 | printf("+++\n+++ Pt bin setted n. ---> %d \n+++\n",fCdfFit->GetPtBin()); | |
309 | ||
310 | if(fOkAODMC) {fCdfFit->CloneMCtemplate(fhDecayTimeMCjpsifromB);} | |
311 | else { | |
312 | //SetPathMCfile("."); | |
313 | TH1F* hLocal = OpenLocalMCFile(fNameMCfile,fCdfFit->GetPtBin()); | |
314 | fCdfFit->CloneMCtemplate(hLocal); | |
315 | } | |
316 | ||
317 | printf("+++\n+++ MC template histo copied ---> OK\n+++\n"); | |
f9fd1412 | 318 | //fCdfFit->DoMinimization(pseudoProperValues,invMassValues,ncand); |
319 | fCdfFit->DoMinimization(); | |
e11ae259 | 320 | |
321 | if(pseudoProperValues) delete [] pseudoProperValues; | |
322 | pseudoProperValues=NULL; | |
323 | if(invMassValues) delete [] invMassValues; | |
324 | invMassValues=NULL; | |
325 | ||
326 | if(fCdfFit) delete [] fCdfFit; | |
327 | fCdfFit=NULL; | |
328 | ||
329 | } | |
7b54a483 | 330 | |
331 | return; | |
332 | } | |
333 | //_________________________________________________________________________________________________ | |
334 | void AliAnalysisTaskSEBtoJPSItoEle::SetPtCuts(const Double_t ptCuts[2]) | |
335 | { | |
336 | // | |
337 | // apply Pt cuts (lower cut, upper cut) | |
338 | // | |
339 | for(Int_t i = 0; i < 2; i++) fPtCuts[i] = ptCuts[i]; | |
340 | } | |
341 | //_________________________________________________________________________________________________ | |
342 | void AliAnalysisTaskSEBtoJPSItoEle::SetCutsD0(const Double_t cuts[9]) | |
343 | { | |
344 | // | |
345 | // apply D0 and JPSI cuts | |
346 | // | |
347 | for(Int_t i = 0; i < 9; i++) fCuts[i] = cuts[i]; | |
348 | } | |
349 | //_________________________________________________________________________________________________ | |
f9fd1412 | 350 | void AliAnalysisTaskSEBtoJPSItoEle::ReadAODMCInfo(AliAODEvent* aodEvent, const TClonesArray* inputArray) |
7b54a483 | 351 | { |
352 | // | |
353 | // Reads MC information if needed (for example in case of parametrized PID) | |
354 | // | |
355 | ||
356 | // load MC particles | |
f9fd1412 | 357 | TClonesArray* mcArray = |
358 | dynamic_cast<TClonesArray*>(aodEvent->FindListObject(AliAODMCParticle::StdBranchName())); | |
359 | if (!mcArray) AliError("Could not find Monte-Carlo in AOD"); | |
360 | AliInfo(Form("MC particles found in mcArray ---> %d ",mcArray->GetEntriesFast())); | |
7b54a483 | 361 | |
362 | // load MC header | |
363 | AliAODMCHeader* mcHeader = | |
364 | (AliAODMCHeader*)aodEvent->GetList()->FindObject(AliAODMCHeader::StdBranchName()); | |
365 | if(!mcHeader){ | |
366 | printf("AliAnalysisTaskSEBtoJPSItoEle::UserExec: MC AOD header branch not found!\n"); | |
367 | } | |
368 | ||
f9fd1412 | 369 | AliAODVertex *vtx1 = (AliAODVertex*)aodEvent->GetPrimaryVertex(); |
370 | ||
371 | //Double_t rmax = 0.000005; | |
7b54a483 | 372 | Double_t mcPrimVtx[3]; |
373 | ||
374 | // get MC primary Vtx | |
375 | mcHeader->GetVertex(mcPrimVtx); | |
376 | ||
377 | Int_t nInCandidates = inputArray->GetEntriesFast(); | |
378 | printf("Number of Candidates for MC analysis: %d\n",nInCandidates); | |
379 | ||
380 | Int_t lab0, lab1, pdgMother, labMother, pdgJPSI, pdg0, pdg1, labJPSIMother, pdgJPSIMother; | |
381 | Int_t labJPSIdaugh0=0; | |
382 | Int_t labJPSIdaugh1=0; | |
383 | ||
384 | for (Int_t iCandidate = 0; iCandidate < nInCandidates; iCandidate++) { | |
f9fd1412 | 385 | |
7b54a483 | 386 | AliAODRecoDecayHF2Prong *dd = (AliAODRecoDecayHF2Prong*)inputArray->UncheckedAt(iCandidate); |
f9fd1412 | 387 | Int_t mcLabelSec = dd->MatchToMC(443,mcArray) ; // select only reco JPSI that are true JPSI |
388 | ||
389 | Double_t vtxPrim[3] = {0.,0.,0.}; | |
390 | Double_t vtxSec[3] = {0.,0.,0.}; | |
391 | dd->GetSecondaryVtx(vtxSec); | |
392 | Bool_t unsetvtx=kFALSE; | |
393 | if(!dd->GetOwnPrimaryVtx()) { | |
394 | dd->SetOwnPrimaryVtx(vtx1); // needed to compute all variables | |
395 | unsetvtx=kTRUE; | |
396 | } | |
397 | ||
7b54a483 | 398 | if(dd->Pt() < fPtCuts[1] && dd->Pt() > fPtCuts[0]){ // apply pt cut only |
399 | ||
400 | AliAODTrack *trk0 = (AliAODTrack*)dd->GetDaughter(0); | |
401 | AliAODTrack *trk1 = (AliAODTrack*)dd->GetDaughter(1); | |
402 | lab0 = trk0->GetLabel(); | |
403 | lab1 = trk1->GetLabel(); | |
404 | ||
405 | AliAODMCParticle* part0 = (AliAODMCParticle*)mcArray->At(lab0); | |
406 | if(!part0) { | |
407 | printf("no MC particle\n"); | |
408 | continue; | |
409 | } | |
410 | ||
411 | while(part0->GetMother()>=0) { | |
412 | labMother=part0->GetMother(); | |
413 | part0 = (AliAODMCParticle*)mcArray->At(labMother); | |
414 | if(!part0) { | |
415 | printf("no MC mother particle\n"); | |
416 | break; | |
417 | } | |
418 | pdgMother = TMath::Abs(part0->GetPdgCode()); | |
419 | if(pdgMother==443) {//this for JPSI | |
420 | labJPSIdaugh0=labMother; | |
421 | break; | |
422 | } | |
423 | } | |
424 | ||
425 | AliAODMCParticle* part1 = (AliAODMCParticle*)mcArray->At(lab1); | |
426 | if(!part1) { | |
427 | printf("no MC particle\n"); | |
428 | continue; | |
429 | } | |
430 | ||
431 | while(part1->GetMother()>=0) { | |
432 | labMother=part1->GetMother(); | |
433 | part1 = (AliAODMCParticle*)mcArray->At(labMother); | |
434 | if(!part1) { | |
435 | printf("no MC mother particle\n"); | |
436 | break; | |
437 | } | |
438 | pdgMother = TMath::Abs(part1->GetPdgCode()); | |
439 | if(pdgMother==443) {//this for JPSI | |
440 | labJPSIdaugh1=labMother; | |
441 | break; | |
442 | } | |
443 | } | |
444 | ||
f9fd1412 | 445 | if (mcLabelSec == -1) |
446 | { | |
447 | AliDebug(2,"No MC particle found"); | |
448 | } | |
449 | else { | |
450 | ||
7b54a483 | 451 | if(labJPSIdaugh0>=0 && labJPSIdaugh1>=0 && labJPSIdaugh0==labJPSIdaugh1) { // check if JPSI is signal |
f9fd1412 | 452 | AliAODMCParticle* partJPSI = (AliAODMCParticle*)mcArray->At(labJPSIdaugh0); |
453 | pdgJPSI = partJPSI->GetPdgCode(); | |
454 | Int_t pdaugh0 = partJPSI->GetDaughter(0); | |
455 | Int_t pdaugh1 = partJPSI->GetDaughter(1); | |
456 | AliAODMCParticle* partDaugh0 = (AliAODMCParticle*)mcArray->At(pdaugh0); | |
457 | AliAODMCParticle* partDaugh1 = (AliAODMCParticle*)mcArray->At(pdaugh1); | |
458 | pdg0 = partDaugh0->GetPdgCode(); | |
459 | pdg1 = partDaugh1->GetPdgCode(); | |
460 | if(TMath::Abs(pdg0) == 11 && TMath::Abs(pdg1) == 11){ // this is for MC JPSI -> ee | |
461 | labJPSIMother = partJPSI->GetMother(); | |
462 | AliAODMCParticle* partJPSIMother = (AliAODMCParticle*)mcArray->At(labJPSIMother); | |
463 | pdgJPSIMother = partJPSIMother->GetPdgCode(); | |
464 | ||
465 | // chech if JPSI comes from B | |
466 | if( pdgJPSIMother==511 || pdgJPSIMother==521 || | |
467 | pdgJPSIMother==10511 || pdgJPSIMother==10521 || | |
468 | pdgJPSIMother==513 || pdgJPSIMother==523 || | |
469 | pdgJPSIMother==10513 || pdgJPSIMother==10523 || | |
470 | pdgJPSIMother==20513 || pdgJPSIMother==20523 || | |
471 | pdgJPSIMother==515 || pdgJPSIMother==525 || | |
472 | pdgJPSIMother==531 || pdgJPSIMother==10531 || | |
473 | pdgJPSIMother==533 || pdgJPSIMother==10533 || | |
474 | pdgJPSIMother==20533 || pdgJPSIMother==535 || | |
475 | pdgJPSIMother==541 || pdgJPSIMother==10541 || | |
476 | pdgJPSIMother==543 || pdgJPSIMother==10543 || | |
477 | pdgJPSIMother==20543 || pdgJPSIMother==545) | |
478 | { //this is for MC JPSI -> ee from B-hadrons | |
479 | ||
480 | fhInvMass->Fill(dd->InvMassJPSIee()); | |
481 | fhD0->Fill(10000*dd->ImpParXY()); | |
482 | fhD0D0->Fill(1e8*dd->Prodd0d0()); | |
483 | fhCosThetaStar->Fill(dd->CosThetaStarJPSI()); | |
484 | fhCtsVsD0D0->Fill(dd->CosThetaStarJPSI(),1e8*dd->Prodd0d0()); | |
485 | fhCosThetaPointing->Fill(dd->CosPointingAngle()); | |
486 | ||
487 | // compute X variable | |
488 | AliAODVertex* primVertex = dd->GetOwnPrimaryVtx(); | |
489 | vtxPrim[0] = primVertex->GetX(); | |
490 | vtxPrim[1] = primVertex->GetY(); | |
491 | vtxPrim[2] = primVertex->GetZ(); | |
492 | Double_t lxy = ((vtxSec[0]-vtxPrim[0])*(dd->Px()) + (vtxSec[1]-vtxPrim[1])*(dd->Py()))/dd->Pt(); | |
493 | Double_t pseudoProperDecayTime = lxy*(TDatabasePDG::Instance()->GetParticle(443)->Mass())/dd->Pt(); | |
494 | ||
495 | fhDecayTime->Fill(10000*pseudoProperDecayTime); | |
496 | // Post the data already here | |
497 | PostData(1,fOutput); | |
498 | ||
499 | fNtupleJPSI->Fill(pseudoProperDecayTime,dd->InvMassJPSIee()); // fill N-tuple with X,M values | |
7b54a483 | 500 | |
501 | Double_t mcLxy = ((partJPSI->Xv()-mcPrimVtx[0])*(partJPSI->Px()) + (partJPSI->Yv()-mcPrimVtx[1])*(partJPSI->Py()))/partJPSI->Pt(); | |
502 | Double_t mcPseudoProperDecayTime = mcLxy*(TDatabasePDG::Instance()->GetParticle(443)->Mass())/partJPSI->Pt(); | |
503 | fhDecayTimeMCjpsifromB->Fill(10000*mcPseudoProperDecayTime,1); | |
504 | ||
505 | // Post the data already here | |
506 | PostData(1,fOutput); | |
507 | ||
f9fd1412 | 508 | } //this is for MC JPSI -> ee from B-hadrons |
509 | } //this is for MC JPSI -> ee | |
510 | } | |
511 | } // end of check if JPSI is signal | |
7b54a483 | 512 | } |
513 | } | |
514 | ||
515 | } | |
516 | //_________________________________________________________________________________________________ | |
517 | TH1F* AliAnalysisTaskSEBtoJPSItoEle::OpenLocalMCFile(TString datadir, Int_t binNum) | |
518 | { | |
519 | // | |
520 | // Open a local file with MC x distribution for JPSI from B-hadron | |
521 | // | |
522 | ||
523 | TH1F* hCsiMC = new TH1F(); | |
524 | TString useFile = datadir.Data(); | |
525 | useFile.Append("/CsiMCfromKine_10PtBins.root"); | |
526 | TFile *f = new TFile(useFile); | |
527 | Double_t scale = 0.; | |
528 | char processCsiMCXHisto[1024]; | |
529 | if(binNum == 0) sprintf(processCsiMCXHisto,"hPseudoProper"); | |
530 | else sprintf(processCsiMCXHisto,"hPseudoProper%d",binNum); | |
531 | hCsiMC = (TH1F*)f->Get(processCsiMCXHisto); | |
532 | scale=1/hCsiMC->Integral(); | |
533 | hCsiMC->Scale(scale); | |
534 | printf ("Opening Histo with Csi_MC(x) template with n. %f entries", hCsiMC->GetEntries()); | |
535 | ||
536 | return hCsiMC; | |
537 | ||
538 | } | |
539 |