]>
Commit | Line | Data |
---|---|---|
aad6618e | 1 | /************************************************************************** |
2 | * Copyright(c) 1998-2007, 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 | // Implementation of the single muon analysis class | |
18 | // An example of usage can be found in the macro runSingleMuAnalysis.C. | |
19 | //---------------------------------------------------------------------------- | |
20 | ||
21 | #define AliAnalysisTaskSingleMu_cxx | |
22 | ||
23 | // ROOT includes | |
24 | #include "Riostream.h" | |
25 | #include "TChain.h" | |
26 | #include "TH1.h" | |
27 | #include "TCanvas.h" | |
28 | #include "TSystem.h" | |
29 | #include "TROOT.h" | |
30 | #include "TParticle.h" | |
31 | #include "TLorentzVector.h" | |
32 | ||
33 | // STEER includes | |
34 | #include "AliLog.h" | |
35 | ||
36 | #include "AliAODEvent.h" | |
37 | #include "AliAODTrack.h" | |
38 | #include "AliAODVertex.h" | |
39 | #include "AliAODInputHandler.h" | |
40 | ||
41 | #include "AliAnalysisTask.h" | |
42 | #include "AliAnalysisDataSlot.h" | |
43 | #include "AliAnalysisManager.h" | |
44 | #include "AliAnalysisTaskSingleMu.h" | |
45 | ||
46 | ClassImp(AliAnalysisTaskSingleMu) | |
47 | ||
48 | //________________________________________________________________________ | |
49 | AliAnalysisTaskSingleMu::AliAnalysisTaskSingleMu(const char *name) : | |
50 | AliAnalysisTask(name,""), | |
51 | fAOD(0), | |
52 | fOutputContainer(0) | |
53 | { | |
54 | // | |
55 | /// Constructor. | |
56 | // | |
57 | ResetHistos(); | |
58 | // Input slot #0 works with an Ntuple | |
59 | DefineInput(0, TChain::Class()); | |
60 | // Output slot #0 writes into a TObjArray container | |
61 | DefineOutput(0, TObjArray::Class()); | |
62 | } | |
63 | ||
64 | //___________________________________________________________________________ | |
65 | void AliAnalysisTaskSingleMu::ConnectInputData(Option_t *) | |
66 | { | |
67 | // | |
68 | /// Connect AOD here | |
69 | /// Called once | |
70 | // | |
71 | ||
72 | TTree* tree = dynamic_cast<TTree*> (GetInputData(0)); | |
73 | if (!tree) { | |
74 | Printf("ERROR: Could not read chain from input slot 0"); | |
75 | } else { | |
76 | /* | |
77 | // Disable all branches and enable only the needed ones | |
78 | // The next two lines are different when data produced as AliAODEvent is read | |
79 | tree->SetBranchStatus("*", kFALSE); | |
80 | tree->SetBranchStatus("fTracks.*", kTRUE); | |
81 | */ | |
82 | ||
83 | AliAODInputHandler *aodH = dynamic_cast<AliAODInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()); | |
84 | ||
85 | if (!aodH) { | |
86 | Printf("ERROR: Could not get AODInputHandler"); | |
87 | } else | |
88 | printf(" ConnectInputData of task %s\n", GetName()); | |
89 | fAOD = aodH->GetEvent(); | |
90 | } | |
91 | } | |
92 | ||
93 | //___________________________________________________________________________ | |
94 | void AliAnalysisTaskSingleMu::CreateOutputObjects() | |
95 | { | |
96 | // | |
97 | /// Create output histograms | |
98 | // | |
99 | printf(" CreateOutputObjects of task %s\n", GetName()); | |
100 | ||
101 | Int_t ptBins = 60; | |
102 | Float_t ptLow = 0., ptHigh = 30.; | |
103 | Char_t *ptName = "P_{t} (GeV/c)"; | |
104 | ||
105 | Int_t vzBins = 40; | |
106 | Float_t vzLow = -20., vzHigh = 20.; | |
107 | Char_t *vzName = "Vz (cm)"; | |
108 | ||
109 | TString baseName, histoName; | |
110 | fOutputContainer = new TObjArray(fgkNhistos*fgkNTrigCuts); | |
111 | fOutputContainer->SetName("SingleMuAnalysisContainer"); | |
112 | Int_t iHisto = 0; | |
113 | ||
114 | for(Int_t iTrig=0; iTrig<fgkNTrigCuts; iTrig++){ | |
115 | ||
116 | // 2D histos | |
117 | if(!fVzVsPt[iTrig]){ | |
118 | baseName = "fVzVsPt"; | |
119 | histoName = baseName + trigName[iTrig]; | |
120 | fVzVsPt[iTrig] = new TH2F(histoName, histoName, | |
121 | ptBins, ptLow, ptHigh, | |
122 | vzBins, vzLow, vzHigh); | |
123 | fVzVsPt[iTrig]->GetXaxis()->SetTitle(ptName); | |
124 | fVzVsPt[iTrig]->GetYaxis()->SetTitle(vzName); | |
125 | ||
126 | fOutputContainer->AddAt(fVzVsPt[iTrig], iHisto); | |
127 | iHisto++; | |
128 | } | |
129 | } | |
130 | } | |
131 | ||
132 | //________________________________________________________________________ | |
133 | void AliAnalysisTaskSingleMu::Exec(Option_t *) | |
134 | { | |
135 | // | |
136 | /// Main loop | |
137 | /// Called for each event | |
138 | // | |
139 | ||
140 | TTree *tinput = (TTree*)GetInputData(0); | |
141 | tinput->GetReadEntry(); | |
142 | ||
143 | if (!fAOD) { | |
144 | Printf("ERROR: fAOD not available"); | |
145 | return; | |
146 | } | |
147 | ||
148 | ||
149 | // Object declaration | |
150 | AliAODTrack *muonTrack = 0x0; | |
151 | TLorentzVector lorVec; | |
152 | Int_t trigMatch = -1; | |
153 | ||
154 | Int_t nTracks = fAOD->GetNumberOfTracks(); | |
155 | for (Int_t itrack = 0; itrack < nTracks; itrack++) { | |
156 | muonTrack = fAOD->GetTrack(itrack); | |
157 | ||
158 | // Apply cuts | |
159 | if(!MuonPassesCuts(*muonTrack, lorVec, trigMatch)) continue; | |
160 | ||
161 | for(Int_t iTrig=0; iTrig<=trigMatch; iTrig++){ | |
162 | fVzVsPt[iTrig]->Fill(lorVec.Pt(), fAOD->GetPrimaryVertex()->GetZ()); | |
163 | } | |
164 | } | |
165 | ||
166 | // Post final data. It will be written to a file with option "RECREATE" | |
167 | PostData(0, fOutputContainer); | |
168 | } | |
169 | ||
170 | //________________________________________________________________________ | |
171 | void AliAnalysisTaskSingleMu::Terminate(Option_t *) { | |
172 | // | |
173 | /// Draw some histogram at the end. | |
174 | // | |
175 | if (!gROOT->IsBatch()) { | |
176 | TCanvas *c1 = new TCanvas("c1","Vz vs Pt",10,10,310,310); | |
177 | c1->SetFillColor(10); c1->SetHighLightColor(10); | |
178 | c1->SetLeftMargin(0.15); c1->SetBottomMargin(0.15); | |
179 | c1->Divide(2,2); | |
180 | for(Int_t iTrig=0; iTrig<fgkNTrigCuts; iTrig++){ | |
181 | c1->cd(iTrig+1); | |
182 | fVzVsPt[iTrig]->DrawCopy("COLZ"); | |
183 | } | |
184 | } | |
185 | } | |
186 | ||
187 | //________________________________________________________________________ | |
188 | void AliAnalysisTaskSingleMu::ResetHistos() | |
189 | { | |
190 | // | |
191 | /// Reset histograms | |
192 | // | |
193 | for(Int_t iTrig=0; iTrig<fgkNTrigCuts; iTrig++){ | |
194 | fVzVsPt[iTrig] = 0x0; | |
195 | } | |
196 | trigName[kNoMatchTrig] = "NoMatchTrig"; | |
197 | trigName[kAllPtTrig] = "AllPtTrig"; | |
198 | trigName[kLowPtTrig] = "LowPtTrig"; | |
199 | trigName[kHighPtTrig] = "HighPtTrig"; | |
200 | } | |
201 | ||
202 | ||
203 | //________________________________________________________________________ | |
204 | Bool_t AliAnalysisTaskSingleMu::MuonPassesCuts(AliAODTrack &muonTrack, | |
205 | TLorentzVector &lorVec, | |
206 | Int_t &trigMatch) | |
207 | { | |
208 | // | |
209 | /// Fill lorentz vector and check cuts | |
210 | // | |
211 | ||
212 | // Check if track is a muon | |
213 | if(muonTrack.GetMostProbablePID()!=AliAODTrack::kMuon) return kFALSE; | |
214 | ||
215 | // Check if track is triggered | |
216 | trigMatch = kNoMatchTrig; | |
217 | if (muonTrack.MatchTriggerHighPt()) { | |
218 | trigMatch = kHighPtTrig; | |
219 | } else if (muonTrack.MatchTriggerLowPt()) { | |
220 | trigMatch = kLowPtTrig; | |
221 | } else if (muonTrack.MatchTriggerAnyPt()){ | |
222 | trigMatch = kAllPtTrig; | |
223 | } | |
224 | ||
225 | // Fill track parameters | |
226 | Double_t px = muonTrack.Px(); | |
227 | Double_t py = muonTrack.Py(); | |
228 | Double_t pz = muonTrack.Pz(); | |
229 | Double_t p = muonTrack.P(); | |
230 | ||
231 | const Double_t kMuonMass = 0.105658369; | |
232 | ||
233 | Double_t energy = TMath::Sqrt(p*p + kMuonMass*kMuonMass); | |
234 | lorVec.SetPxPyPzE(px,py,pz,energy); | |
235 | ||
236 | return kTRUE; | |
237 | } |