1 /**************************************************************************
\r
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
\r
4 * Author: The ALICE Off-line Project. *
\r
5 * Contributors are mentioned in the code where appropriate. *
\r
7 * Permission to use, copy, modify and distribute this software and its *
\r
8 * documentation strictly for non-commercial purposes is hereby granted *
\r
9 * without fee, provided that the above copyright notice appears in all *
\r
10 * copies and that both the copyright notice and this permission notice *
\r
11 * appear in the supporting documentation. The authors make no claims *
\r
12 * about the suitability of this software for any purpose. It is *
\r
13 * provided "as is" without express or implied warranty. *
\r
14 **************************************************************************/
\r
16 //-----------------------------------------------------------------------
\r
17 // Example of task (running locally, on AliEn and CAF),
\r
18 // which provides standard way of calculating acceptance and efficiency
\r
19 // between different steps of the procedure.
\r
20 // The ouptut of the task is a AliCFContainer from which the efficiencies
\r
21 // can be calculated
\r
22 //-----------------------------------------------------------------------
\r
23 // Author : R. Vernet, Consorzio Cometa - Catania (it)
\r
24 //-----------------------------------------------------------------------
\r
25 // Modification done by X. Lopez - LPC Clermont (fr)
\r
26 //-----------------------------------------------------------------------
\r
29 #ifndef ALICFMUONRESTASK1_CXX
\r
30 #define ALICFMUONRESTASK1_CXX
\r
32 #include "AliCFMuonResTask1.h"
\r
33 #include "AliHeader.h"
\r
34 #include "AliESDHeader.h"
\r
35 #include "AliStack.h"
\r
36 #include "TParticle.h"
\r
38 #include "AliMCEvent.h"
\r
39 #include "AliAnalysisManager.h"
\r
40 #include "AliESDEvent.h"
\r
41 #include "AliAODEvent.h"
\r
42 #include "AliCFManager.h"
\r
43 #include "AliCFCutBase.h"
\r
44 #include "AliCFContainer.h"
\r
46 #include "AliESDtrack.h"
\r
48 #include "AliESDMuonTrack.h"
\r
49 #include "AliESDtrack.h"
\r
50 #include "AliESDInputHandler.h"
\r
51 #include "TCanvas.h"
\r
53 ClassImp(AliCFMuonResTask1)
\r
55 //__________________________________________________________________________
\r
56 AliCFMuonResTask1::AliCFMuonResTask1() :
\r
60 fHistEventsProcessed(0x0),
\r
67 //___________________________________________________________________________
\r
68 AliCFMuonResTask1::AliCFMuonResTask1(const Char_t* name) :
\r
69 AliAnalysisTaskSE(name),
\r
73 fHistEventsProcessed(0x0),
\r
77 // Constructor. Initialization of Inputs and Outputs
\r
79 Info("AliCFMuonResTask1","Calling Constructor");
\r
81 fHistEventsProcessed = new TH1I("fHistEventsProcessed","",1,0,1) ;
\r
83 DefineOutput(1,TH1I::Class());
\r
84 DefineOutput(2,AliCFContainer::Class());
\r
88 //___________________________________________________________________________
\r
89 AliCFMuonResTask1& AliCFMuonResTask1::operator=(const AliCFMuonResTask1& c)
\r
92 // Assignment operator
\r
95 AliAnalysisTaskSE::operator=(c) ;
\r
96 fReadAODData = c.fReadAODData ;
\r
97 fCFManager = c.fCFManager;
\r
98 fQAHistList = c.fQAHistList ;
\r
99 fHistEventsProcessed = c.fHistEventsProcessed;
\r
105 //___________________________________________________________________________
\r
106 AliCFMuonResTask1::AliCFMuonResTask1(const AliCFMuonResTask1& c) :
\r
107 AliAnalysisTaskSE(c),
\r
108 fReadAODData(c.fReadAODData),
\r
109 fCFManager(c.fCFManager),
\r
110 fQAHistList(c.fQAHistList),
\r
111 fHistEventsProcessed(c.fHistEventsProcessed),
\r
115 // Copy Constructor
\r
119 //___________________________________________________________________________
\r
120 AliCFMuonResTask1::~AliCFMuonResTask1() {
\r
124 Info("~AliCFMuonResTask1","Calling Destructor");
\r
125 if (fCFManager) delete fCFManager ;
\r
126 if (fHistEventsProcessed) delete fHistEventsProcessed ;
\r
127 if (fQAHistList) {fQAHistList->Clear(); delete fQAHistList;}
\r
130 //_________________________________________________
\r
131 void AliCFMuonResTask1::UserExec(Option_t *)
\r
134 // Main loop function
\r
137 Info("UserExec","") ;
\r
139 Error("UserExec","NO MC EVENT FOUND!");
\r
144 Double_t containerInput[9] ;
\r
150 fCFManager->SetEventInfo(fMCEvent);
\r
151 AliStack* stack = fMCEvent->Stack();
\r
153 // loop on the MC event
\r
154 for (Int_t ipart=0; ipart<fMCEvent->GetNumberOfTracks(); ipart++) {
\r
155 AliMCParticle *mcPart = fMCEvent->GetTrack(ipart);
\r
157 TParticle *part = mcPart->Particle();
\r
158 TParticle *part0 = mcPart->Particle();
\r
159 TParticle *part1 = mcPart->Particle();
\r
161 // Selection of the resonance
\r
162 if (!fCFManager->CheckParticleCuts(AliCFManager::kPartGenCuts,mcPart)) continue;
\r
164 // Mother kinematics
\r
165 Float_t e = part->Energy();
\r
166 Float_t pz = part->Pz();
\r
167 Float_t py = part->Py();
\r
168 Float_t px = part->Px();
\r
169 Float_t rapmc = Rap(e,pz);
\r
170 Float_t mass = part->GetCalcMass();
\r
172 // Decays kinematics
\r
174 Int_t p0 = part->GetDaughter(0);
\r
175 part0 = stack->Particle(p0);
\r
176 // selection of the rapidity for first muon
\r
177 AliMCParticle *mcpart0 = new AliMCParticle(part0);
\r
178 if (!fCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,mcpart0)) continue;
\r
179 Int_t pdg0 = part0->GetPdgCode();
\r
181 Int_t p1 = part->GetDaughter(1);
\r
182 part1 = stack->Particle(p1);
\r
183 // selection of the rapidity for second muon
\r
184 AliMCParticle *mcpart1 = new AliMCParticle(part1);
\r
185 if (!fCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,mcpart1)) continue;
\r
186 Int_t pdg1 = part1->GetPdgCode();
\r
189 Float_t e0, pz0, py0, px0, phi0, rapmc0;
\r
190 Float_t e1, pz1, py1, px1, phi1, rapmc1;
\r
192 // ordering sign: first = mu-
\r
194 e0 = part0->Energy();
\r
198 phi0 = part0->Phi(); // Warning in TParticle Phi = pi + ATan2(Py,Px) = [0,2pi]
\r
199 phi0 = Phideg(phi0);
\r
200 rapmc0=Rap(e0,pz0);
\r
202 e1 = part1->Energy();
\r
206 phi1 = part1->Phi();
\r
207 phi1 = Phideg(phi1);
\r
208 rapmc1=Rap(e1,pz1);
\r
212 e1 = part0->Energy();
\r
216 phi1 = part0->Phi();
\r
217 phi1 = Phideg(phi1);
\r
218 rapmc1=Rap(e1,pz1);
\r
220 e0 = part1->Energy();
\r
224 phi0 = part1->Phi();
\r
225 phi0 = Phideg(phi0);
\r
226 rapmc0=Rap(e0,pz0);
\r
230 if(pdg0==13 || pdg1==13) {
\r
232 Float_t pmc = TMath::Sqrt((px0+px1)*(px0+px1)+(py0+py1)*(py0+py1)+
\r
233 (pz0+pz1)*(pz0+pz1));
\r
234 Float_t ptmc = TMath::Sqrt((px0+px1)*(px0+px1)+(py0+py1)*(py0+py1));
\r
235 Float_t imassmc = Imass(e0,px0,py0,pz0,e1,px1,py1,pz1);
\r
236 Float_t rapmc_check=Rap((e0+e1),(pz0+pz1));
\r
238 containerInput[0] = fNevt ;
\r
239 containerInput[1] = rapmc0 ;
\r
240 containerInput[2] = phi0 ;
\r
241 containerInput[3] = rapmc1 ;
\r
242 containerInput[4] = phi1 ;
\r
243 containerInput[5] = imassmc ;
\r
244 containerInput[6] = rapmc ;
\r
245 containerInput[7] = ptmc;
\r
246 containerInput[8] = pmc;
\r
248 // fill the container at the first step
\r
249 fCFManager->GetParticleContainer()->Fill(containerInput,0);
\r
257 AliESDEvent *fESD;
\r
258 AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*>
\r
259 (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
\r
260 fESD = esdH->GetEvent();
\r
262 Int_t mult1 = fESD->GetNumberOfMuonTracks() ;
\r
264 for (Int_t j = 0; j<mult1; j++) {
\r
265 AliESDMuonTrack* mu1 = new AliESDMuonTrack(*(fESD->GetMuonTrack(j)));
\r
266 Float_t zr1 = mu1->Charge();
\r
269 Float_t pxr1 = mu1->Px();
\r
270 Float_t pyr1 = mu1->Py();
\r
271 Float_t pzr1 = mu1->Pz();
\r
272 Float_t phir1 = mu1->Phi(); // Warning in AliESDMuonTrack Phi = ATan2(Py,Px) = [-pi,pi]
\r
273 phir1 = phir1 * 57.296;
\r
274 Float_t er1 = mu1->E();
\r
275 Float_t rap1=Rap(er1,pzr1);
\r
276 // selection of the rapidity for first muon
\r
277 if (!fCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,mu1)) continue;
\r
279 for (Int_t jj = 0; jj<mult1; jj++) {
\r
280 AliESDMuonTrack* mu2 = new AliESDMuonTrack(*(fESD->GetMuonTrack(jj)));
\r
281 Float_t zr2 = mu2->Charge();
\r
283 if(zr2>0 && jj !=j){
\r
284 Float_t pxr2 = mu2->Px();
\r
285 Float_t pyr2 = mu2->Py();
\r
286 Float_t pzr2 = mu2->Pz();
\r
287 Float_t phir2 = mu2->Phi();
\r
288 phir2 = phir2 * 57.296;
\r
289 Float_t er2 = mu2->E();
\r
290 Float_t rap2=Rap(er2,pzr2);
\r
291 // selection of the rapidity for second muon
\r
292 if (!fCFManager->CheckParticleCuts(AliCFManager::kPartRecCuts,mu2)) continue;
\r
294 Float_t prec = TMath::Sqrt((pxr1+pxr2)*(pxr1+pxr2)+(pyr1+pyr2)*(pyr1+pyr2)+
\r
295 (pzr1+pzr2)*(pzr1+pzr2));
\r
296 Float_t ptrec = TMath::Sqrt((pxr1+pxr2)*(pxr1+pxr2)+(pyr1+pyr2)*(pyr1+pyr2));
\r
297 Float_t raprec=Rap((er1+er2),(pzr1+pzr2));
\r
298 Float_t imassrec = Imass(er1,pxr1,pyr1,pzr1,er2,pxr2,pyr2,pzr2);
\r
300 containerInput[0] = fNevt ;
\r
301 containerInput[1] = rap1 ;
\r
302 containerInput[2] = phir1 ;
\r
303 containerInput[3] = rap2 ;
\r
304 containerInput[4] = phir2 ;
\r
305 containerInput[5] = imassrec ;
\r
306 containerInput[6] = raprec ;
\r
307 containerInput[7] = ptrec;
\r
308 containerInput[8] = prec;
\r
310 // fill the container at the second step
\r
311 fCFManager->GetParticleContainer()->Fill(containerInput,1);
\r
314 } // second mu Loop
\r
319 fHistEventsProcessed->Fill(0);
\r
320 PostData(1,fHistEventsProcessed) ;
\r
321 PostData(2,fCFManager->GetParticleContainer()) ;
\r
323 //________________________________________________________________________
\r
324 const Float_t AliCFMuonResTask1::Imass(Float_t e1, Float_t px1, Float_t py1, Float_t pz1,
\r
325 Float_t e2, Float_t px2, Float_t py2, Float_t pz2)
\r
327 // invariant mass calculation
\r
328 Float_t imassrec = TMath::Sqrt((e1+e2)*(e1+e2)-((px1+px2)*(px1+px2)+
\r
329 (py1+py2)*(py1+py2)+(pz1+pz2)*(pz1+pz2)));
\r
332 //________________________________________________________________________
\r
333 const Float_t AliCFMuonResTask1::Rap(Float_t e, Float_t pz)
\r
335 // calculate rapidity
\r
338 rap = 0.5*TMath::Log((e+pz)/(e-pz));
\r
346 //________________________________________________________________________
\r
347 const Float_t AliCFMuonResTask1::Phideg(Float_t phi)
\r
349 // calculate Phi from TParticle in range [-180,180]
\r
352 phideg = phi-TMath::Pi();
\r
353 phideg = phideg*57.296;
\r
356 //________________________________________________________________________
\r
357 void AliCFMuonResTask1::Terminate(Option_t *)
\r
359 // draw result of the Invariant mass MC and ESD
\r
361 AliCFContainer *cont = dynamic_cast<AliCFContainer*> (GetOutputData(2));
\r
363 TH1D *kmass = cont->ShowProjection(5,0);
\r
364 TH1D *rmass = cont->ShowProjection(5,1);
\r
366 TCanvas *c1 = new TCanvas("AliCFMuonResTask1","JPSI MC & ESD",10,10,510,510);
\r
369 kmass->Draw("HIST");
\r
371 rmass->Draw("HIST");
\r
373 //________________________________________________________________________
\r