]>
Commit | Line | Data |
---|---|---|
6419bebb | 1 | // AliTRDonlineTrackletFilter implements an analysis task which reads |
2 | // the tracklets from the TRD.Tracklets.root file and makes them | |
3 | // available to other analysis tasks. This allows to access the full | |
4 | // tracklet information without repeating the techncicalities in each | |
5 | // analysis task. The tracklets are made available in an exchange | |
6 | // container on the output slot 1. | |
7 | ||
b220e595 | 8 | #include "TFile.h" |
9 | #include "TTree.h" | |
10 | #include "TChain.h" | |
11 | #include "TH1F.h" | |
12 | #include "TH2F.h" | |
13 | #include "TProfile.h" | |
14 | #include "TCanvas.h" | |
15 | ||
16 | #include "AliAnalysisManager.h" | |
17 | #include "AliESDEvent.h" | |
18 | #include "AliAODEvent.h" | |
19 | #include "AliMCEvent.h" | |
20 | #include "AliESDInputHandler.h" | |
21 | #include "AliAODHandler.h" | |
22 | #include "AliMCEventHandler.h" | |
23 | #include "AliLog.h" | |
24 | #include "AliESDTrdTrack.h" | |
25 | ||
26 | #include "AliTRDtrackletMCM.h" | |
27 | #include "AliTRDtrackletWord.h" | |
28 | #include "AliVParticle.h" | |
29 | #include "AliMCParticle.h" | |
30 | ||
31 | #include "AliTRDonlineTrackletFilter.h" | |
32 | ||
33 | ClassImp(AliTRDonlineTrackletFilter) | |
34 | ||
35 | AliTRDonlineTrackletFilter::AliTRDonlineTrackletFilter(const char *name) : | |
36 | AliAnalysisTask(name, ""), | |
37 | fESD(0x0), | |
38 | fInputHandler(0x0), | |
39 | fInputEvent(0x0), | |
40 | fOutputAOD(0x0), | |
41 | fMCEvent(0x0), | |
42 | fTrackletsRaw(new TClonesArray("AliTRDtrackletWord")), | |
43 | fTrackletsSim(new TClonesArray("AliTRDtrackletMCM")), | |
54d34aac | 44 | fTrackletTree(0x0), |
b220e595 | 45 | fGeo(new AliTRDgeometry), |
46 | fNevent(0), | |
54d34aac | 47 | fPath(""), |
b220e595 | 48 | fTrackletFile(0x0), |
54d34aac | 49 | fNEventsPerFile(0), |
50 | fEvent(0), | |
51 | fFileNumber(0), | |
b220e595 | 52 | fTrackletTreeSim(0x0), |
54d34aac | 53 | fTrackletTreeRaw(0x0) |
b220e595 | 54 | { |
55 | // ctor | |
56 | ||
57 | DefineInput(0, TChain::Class()); | |
58 | ||
59 | DefineOutput(0, TTree::Class()); | |
60 | DefineOutput(1, TTree::Class()); | |
61 | } | |
62 | ||
63 | AliTRDonlineTrackletFilter::~AliTRDonlineTrackletFilter() | |
64 | { | |
65 | // dtor | |
66 | ||
67 | delete fTrackletsRaw; | |
68 | delete fTrackletsSim; | |
69 | delete fGeo; | |
70 | } | |
71 | ||
54d34aac | 72 | void AliTRDonlineTrackletFilter::ConnectInputData(const Option_t */* option */) |
b220e595 | 73 | { |
6419bebb | 74 | // connect all the input handlers to access the data |
75 | ||
b220e595 | 76 | fInputHandler = (AliInputEventHandler*) AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler(); |
77 | if (fInputHandler) | |
78 | fInputEvent = fInputHandler->GetEvent(); | |
79 | ||
80 | AliMCEventHandler *mcH = (AliMCEventHandler*) AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler(); | |
81 | if (mcH) | |
82 | fMCEvent = mcH->MCEvent(); | |
83 | } | |
84 | ||
85 | void AliTRDonlineTrackletFilter::CreateOutputObjects() | |
86 | { | |
6419bebb | 87 | // create the output objects used to make the tracklets available to |
88 | // other tasks | |
89 | ||
b220e595 | 90 | OpenFile(1); |
91 | ||
92 | fTrackletTree = new TTree("tracklets", "on-line tracklets"); | |
93 | fTrackletTree->Branch("tracklets_sim", fTrackletsSim); | |
94 | fTrackletTree->Branch("tracklets_raw", fTrackletsRaw); | |
95 | } | |
96 | ||
97 | Bool_t AliTRDonlineTrackletFilter::Notify() | |
98 | { | |
6419bebb | 99 | // we use notify to figure out which TRD.Tracklets.root file we have |
100 | // to use for accessing the tracklets | |
b220e595 | 101 | |
102 | TString filename(AliAnalysisManager::GetAnalysisManager()->GetTree()->GetCurrentFile()->GetName()); | |
103 | ||
104 | AliInfo(Form("Now reading from %s", filename.Data())); | |
105 | ||
106 | if (filename.Contains("AliAOD.root")) | |
107 | filename.ReplaceAll("AliAOD.root", ""); | |
108 | else if (filename.Contains("AliESDs.root")) | |
109 | filename.ReplaceAll("AliESDs.root", ""); | |
110 | else if (filename.Contains("galice.root")) | |
111 | filename.ReplaceAll("galice.root", ""); | |
112 | else if (filename.BeginsWith("root:")) | |
113 | filename.Append("?ZIP="); | |
114 | ||
115 | fPath = filename; | |
116 | ||
117 | fTrackletFile = TFile::Open(Form("%sTRD.Tracklets.root", fPath.Data())); | |
118 | ||
119 | if (!fTrackletFile) { | |
120 | AliError("No tracklet file"); | |
121 | return kFALSE; | |
122 | } | |
123 | ||
124 | fNEventsPerFile = fTrackletFile->GetNkeys() - fTrackletFile->GetNProcessIDs(); //??? | |
125 | ||
126 | fEvent = -1; | |
127 | fFileNumber = 0; | |
54d34aac | 128 | |
129 | return kTRUE; | |
b220e595 | 130 | } |
131 | ||
132 | ||
54d34aac | 133 | void AliTRDonlineTrackletFilter::Exec(const Option_t * /* option */) |
b220e595 | 134 | { |
135 | // execute this for each event | |
136 | ||
137 | if (!LoadEvent()) | |
138 | return; | |
139 | ||
140 | // ----- simulated tracklets ----- | |
141 | if (fTrackletTreeSim) { | |
142 | AliTRDtrackletMCM *trkl = 0x0; | |
143 | ||
144 | TBranch *br = fTrackletTreeSim->GetBranch("mcmtrklbranch"); | |
145 | br->SetAddress(&trkl); | |
146 | ||
147 | for (Int_t iTracklet = 0; iTracklet < br->GetEntries(); iTracklet++) { | |
148 | br->GetEntry(iTracklet); | |
149 | new ((*fTrackletsSim)[fTrackletsSim->GetEntriesFast()]) AliTRDtrackletMCM(*trkl); | |
150 | } | |
151 | } | |
152 | ||
153 | // ----- raw tracklets ----- | |
154 | if (fTrackletTreeRaw) { | |
155 | Int_t hc; | |
156 | TClonesArray *trklArray = 0x0; | |
157 | fTrackletTreeRaw->SetBranchAddress("hc", &hc); | |
158 | fTrackletTreeRaw->SetBranchAddress("trkl", &trklArray); | |
159 | for (Int_t iHCidx = 0; iHCidx < fTrackletTreeRaw->GetEntries(); iHCidx++) { | |
160 | fTrackletTreeRaw->GetEntry(iHCidx); | |
161 | for (Int_t iTracklet = 0; iTracklet < trklArray->GetEntries(); iTracklet++) { | |
162 | AliTRDtrackletWord *trklWord = (AliTRDtrackletWord*) ((*trklArray)[iTracklet]); | |
163 | trklWord->SetDetector(hc/2); | |
164 | new ((*fTrackletsRaw)[fTrackletsRaw->GetEntriesFast()]) AliTRDtrackletWord(*trklWord); | |
165 | } | |
166 | } | |
167 | } | |
168 | ||
169 | AliInfo(Form("%i tracklets", fTrackletsSim->GetEntriesFast())); | |
170 | fTrackletTree->SetBranchAddress("tracklets_sim", &fTrackletsSim); | |
171 | fTrackletTree->SetBranchAddress("tracklets_raw", &fTrackletsRaw); | |
172 | fTrackletTree->Fill(); | |
173 | PostData(1, fTrackletTree); | |
174 | } | |
175 | ||
176 | void AliTRDonlineTrackletFilter::LocalInit() | |
177 | { | |
178 | ||
179 | } | |
180 | ||
54d34aac | 181 | void AliTRDonlineTrackletFilter::Terminate(const Option_t * /* option */) |
b220e595 | 182 | { |
183 | ||
184 | } | |
185 | ||
186 | Bool_t AliTRDonlineTrackletFilter::LoadEvent() | |
187 | { | |
188 | // load tracklets for the current event | |
189 | ||
190 | // ----- cleaning ----- | |
191 | fTrackletsSim->Delete(); | |
192 | fTrackletsRaw->Delete(); | |
193 | ||
194 | // ----- initialization ----- | |
195 | if (!fInputEvent) { | |
196 | AliError("No event found!"); | |
197 | return kFALSE; | |
198 | } | |
199 | fESD = dynamic_cast<AliESDEvent*> (fInputEvent); | |
200 | ||
201 | fEvent++; | |
202 | Int_t inew = fEvent / fNEventsPerFile; | |
203 | if ( inew != fFileNumber) { | |
204 | fFileNumber++; | |
205 | ||
206 | delete fTrackletFile; | |
207 | fTrackletFile = TFile::Open(Form("%sTRD.Tracklets%d.root", fPath.Data(), fFileNumber)); | |
208 | ||
209 | if (!fTrackletFile) { | |
210 | AliError("No tracklet file"); | |
211 | return kFALSE; | |
212 | } | |
213 | } | |
214 | ||
215 | if (!fTrackletFile) { | |
216 | AliError("no tracklet file"); | |
217 | return kFALSE; | |
218 | } | |
219 | ||
220 | // tracklets from simulation | |
221 | char treename[30]; | |
69503e16 | 222 | snprintf(treename, 30, "Event%d/tracklets", fEvent); |
b220e595 | 223 | |
224 | fTrackletTreeSim = (TTree*) fTrackletFile->Get(treename); | |
225 | ||
226 | // tracklets from raw | |
227 | char treenameRaw[30]; | |
69503e16 | 228 | snprintf(treenameRaw, 30, "Event%d/tracklets-raw", fEvent); |
b220e595 | 229 | fTrackletTreeRaw = (TTree*) fTrackletFile->Get(treenameRaw); |
230 | ||
231 | return kTRUE; | |
232 | } | |
233 |