2 #include "AliAODEvent.h"
3 #include "AliAODMCParticle.h"
4 #include "AliAODMCHeader.h"
5 #include "AliAnalysisTaskMFTExample.h"
6 #include "TDatabasePDG.h"
7 #include "AliMUONTrackParam.h"
8 #include "AliMUONTrackExtrap.h"
9 #include "AliAODDimuon.h"
10 #include "AliAODTrack.h"
11 #include "AliAODHeader.h"
12 #include "TClonesArray.h"
13 #include "AliMFTConstants.h"
14 #include "AliMFTAnalysisTools.h"
18 ClassImp(AliAnalysisTaskMFTExample)
20 //====================================================================================================================================================
22 AliAnalysisTaskMFTExample::AliAnalysisTaskMFTExample() :
26 fHistPtSingleMuons(0),
27 fHistPtSingleMuonsFromJpsi(0),
29 fHistMassDimuonsOS(0),
30 fHistPtDimuonsJpsi(0),
31 fHistMassDimuonsJpsi(0),
32 fHistResidualXVtxJpsi(0),
33 fHistResidualYVtxJpsi(0),
34 fHistResidualZVtxJpsi(0) {
36 // Default constructor
38 for (Int_t i=0; i<3; i++) fPrimaryVertex[i] = 0;
39 fVtxResolutionITS[0] = 5.e-4;
40 fVtxResolutionITS[1] = 5.e-4;
41 fVtxResolutionITS[2] = 4.e-4;
45 //====================================================================================================================================================
47 AliAnalysisTaskMFTExample::AliAnalysisTaskMFTExample(const Char_t *name) :
48 AliAnalysisTaskSE(name),
51 fHistPtSingleMuons(0),
52 fHistPtSingleMuonsFromJpsi(0),
54 fHistMassDimuonsOS(0),
55 fHistPtDimuonsJpsi(0),
56 fHistMassDimuonsJpsi(0),
57 fHistResidualXVtxJpsi(0),
58 fHistResidualYVtxJpsi(0),
59 fHistResidualZVtxJpsi(0) {
63 for (Int_t i=0; i<3; i++) fPrimaryVertex[i] = 0;
64 fVtxResolutionITS[0] = 5.e-4;
65 fVtxResolutionITS[1] = 5.e-4;
66 fVtxResolutionITS[2] = 4.e-4;
68 // Define input and output slots here
69 DefineOutput(1, TList::Class());
73 //====================================================================================================================================================
75 void AliAnalysisTaskMFTExample::UserCreateOutputObjects() {
79 fHistogramList = new TList();
80 fHistogramList->SetOwner(kTRUE);
82 fHistPtSingleMuons = new TH1D("fHistPtSingleMuons","p_{T} of single muons (All)", 100, 0, 10);
83 fHistPtSingleMuons -> SetXTitle("p_{T} [GeV/c]");
84 fHistPtSingleMuons -> Sumw2();
86 fHistPtSingleMuonsFromJpsi = new TH1D("fHistPtSingleMuonsFromJpsi", "p_{T} of single muons (from J/#psi)", 100, 0, 10);
87 fHistPtSingleMuonsFromJpsi -> SetXTitle("p_{T} [GeV/c]");
88 fHistPtSingleMuonsFromJpsi -> Sumw2();
90 fHistMassDimuonsOS = new TH1D("fHistMassDimuonsOS", "Mass of OS dimuons", 500, 0, 10);
91 fHistMassDimuonsOS -> SetXTitle("Mass [GeV/c^{2}]");
92 fHistMassDimuonsOS -> Sumw2();
94 fHistMassDimuonsJpsi = new TH1D("fHistMassDimuonsJpsi", "Mass of J/#psi dimuons", 500, 0, 10);
95 fHistMassDimuonsJpsi -> SetXTitle("Mass [GeV/c^{2}]");
96 fHistMassDimuonsJpsi -> Sumw2();
98 fHistPtDimuonsOS = new TH1D("fHistPtDimuonsOS", "p_{T} of OS dimuons", 100, 0, 10);
99 fHistPtDimuonsOS -> SetXTitle("p_{T} [GeV/c]");
100 fHistPtDimuonsOS -> Sumw2();
102 fHistPtDimuonsJpsi = new TH1D("fHistPtDimuonsJpsi", "p_{T} of J/#psi dimuons", 100, 0, 10);
103 fHistPtDimuonsJpsi -> SetXTitle("p_{T} [GeV/c]");
104 fHistPtDimuonsJpsi -> Sumw2();
106 fHistResidualXVtxJpsi = new TH1D("fHistResidualXVtxJpsi", "J/#psi vertex residuals along x", 100, -100, 100);
107 fHistResidualXVtxJpsi -> SetXTitle("Residual [#mum]");
108 fHistResidualXVtxJpsi -> Sumw2();
110 fHistResidualYVtxJpsi = new TH1D("fHistResidualYVtxJpsi", "J/#psi vertex residuals along y", 100, -100, 100);
111 fHistResidualYVtxJpsi -> SetXTitle("Residual [#mum]");
112 fHistResidualYVtxJpsi -> Sumw2();
114 fHistResidualZVtxJpsi = new TH1D("fHistResidualZVtxJpsi", "J/#psi vertex residuals along z", 100, -1000, 1000);
115 fHistResidualZVtxJpsi -> SetXTitle("Residual [#mum]");
116 fHistResidualZVtxJpsi -> Sumw2();
118 fHistogramList->Add(fHistPtSingleMuons);
119 fHistogramList->Add(fHistPtSingleMuonsFromJpsi);
120 fHistogramList->Add(fHistMassDimuonsOS);
121 fHistogramList->Add(fHistMassDimuonsJpsi);
122 fHistogramList->Add(fHistPtDimuonsOS);
123 fHistogramList->Add(fHistPtDimuonsJpsi);
124 fHistogramList->Add(fHistResidualXVtxJpsi);
125 fHistogramList->Add(fHistResidualYVtxJpsi);
126 fHistogramList->Add(fHistResidualZVtxJpsi);
128 PostData(1, fHistogramList);
132 //====================================================================================================================================================
134 void AliAnalysisTaskMFTExample::UserExec(Option_t *) {
137 // Called for each event
139 AliAODEvent *aodEv = dynamic_cast<AliAODEvent*>(InputEvent());
142 aodEv->GetHeader()->InitMagneticField();
143 AliMUONTrackExtrap::SetField();
145 TClonesArray *stackMC = (TClonesArray*) (aodEv->GetList()->FindObject(AliAODMCParticle::StdBranchName()));
147 // Getting primary vertex, either from the generation or from the reconstruction -------------------
149 if (fVertexMode == kGenerated) {
150 AliAODMCHeader *mcHeader = (AliAODMCHeader*) (aodEv->GetList()->FindObject(AliAODMCHeader::StdBranchName()));
151 mcHeader->GetVertex(fPrimaryVertex);
152 for (Int_t i=0; i<3; i++) fPrimaryVertex[i] = gRandom->Gaus(fPrimaryVertex[i], fVtxResolutionITS[i]);
154 else if (fVertexMode == kReconstructed) {
155 aodEv->GetPrimaryVertex()->GetXYZ(fPrimaryVertex);
158 // -------------------------------------------------------------------------------------------------
160 AliAODTrack *recMuon1=0, *recMuon2=0;
161 AliAODMCParticle *mcMuon1=0, *mcMuon2=0;
162 AliAODMCParticle *mcMother1=0;
163 Bool_t isMuon1FromJpsi=0;
165 // Loop over MUON+MFT muons
167 //--------------------------------------------------------------------------------------------------
169 for (Int_t iTrack=0; iTrack<aodEv->GetNumberOfTracks(); iTrack++) {
171 if (!(aodEv->GetTrack(iTrack)->IsMuonGlobalTrack())) continue;
173 recMuon1 = aodEv->GetTrack(iTrack);
174 isMuon1FromJpsi = kFALSE;
177 fHistPtSingleMuons -> Fill(recMuon1->Pt());
179 // pt muons from J/psi
180 if (recMuon1->GetLabel() >= 0) {
181 mcMuon1 = (AliAODMCParticle*) stackMC->At(recMuon1->GetLabel());
183 if (mcMuon1->GetMother()>=0) {
184 mcMother1 = (AliAODMCParticle*) stackMC->At(mcMuon1->GetMother());
185 if (mcMother1->PdgCode()==443) {
186 isMuon1FromJpsi = kTRUE;
187 fHistPtSingleMuonsFromJpsi -> Fill(recMuon1->Pt());
193 for (Int_t jTrack=0; jTrack<iTrack; jTrack++) {
195 if (!(aodEv->GetTrack(jTrack)->IsMuonGlobalTrack())) continue;
197 recMuon2 = aodEv->GetTrack(jTrack);
199 AliAODDimuon *dimuon = new AliAODDimuon(recMuon1, recMuon2);
200 if (dimuon->Charge()) continue;
202 // pt and mass all OS dimuons
203 fHistPtDimuonsOS -> Fill(dimuon->Pt());
204 fHistMassDimuonsOS -> Fill(dimuon->Mass());
206 // pt and mass J/psi dimuons
207 if (!isMuon1FromJpsi) continue;
208 if (recMuon2->GetLabel() >= 0) {
209 mcMuon2 = (AliAODMCParticle*) stackMC->At(recMuon2->GetLabel());
211 if (mcMuon2->GetMother() == mcMuon1->GetMother()) {
213 Double_t pcaQuality=0;
214 TClonesArray *muons = new TClonesArray("AliAODTrack",2);
215 muons -> Add(recMuon1);
216 muons -> Add(recMuon2);
217 TLorentzVector kinem(0,0,0,0);
218 if (!AliMFTAnalysisTools::RefitAODDimuonWithCommonVertex(muons, pca, pcaQuality, kinem)) continue;
219 fHistPtDimuonsJpsi -> Fill(kinem.Pt());
220 fHistMassDimuonsJpsi -> Fill(kinem.M());
221 fHistResidualXVtxJpsi -> Fill(1.e4*(pca[0] - fPrimaryVertex[0]));
222 fHistResidualYVtxJpsi -> Fill(1.e4*(pca[1] - fPrimaryVertex[1]));
223 fHistResidualZVtxJpsi -> Fill(1.e4*(pca[2] - fPrimaryVertex[2]));
228 } // end of loop on 2nd muon
230 } // end of loop on 1st muon
232 //--------------------------------------------------------------------------------------------------
234 PostData(1, fHistogramList);
238 //====================================================================================================================================================
240 void AliAnalysisTaskMFTExample::Terminate(Option_t *) {
242 // Draw result to the screen
243 // Called once at the end of the query
247 //====================================================================================================================================================