]>
Commit | Line | Data |
---|---|---|
c7b7f445 | 1 | // $Id$ |
2 | //************************************************************************** | |
3 | //* This file is property of and copyright by the ALICE HLT Project * | |
4 | //* ALICE Experiment at CERN, All rights reserved. * | |
5 | //* * | |
6 | //* Primary Authors: Matthias Richter <Matthias.Richter@ift.uib.no> * | |
7 | //* Jochen Thaeder <jochen@thaeder.de> * | |
8 | //* for The ALICE HLT Project. * | |
9 | //* * | |
10 | //* Permission to use, copy, modify and distribute this software and its * | |
11 | //* documentation strictly for non-commercial purposes is hereby granted * | |
12 | //* without fee, provided that the above copyright notice appears in all * | |
13 | //* copies and that both the copyright notice and this permission notice * | |
14 | //* appear in the supporting documentation. The authors make no claims * | |
15 | //* about the suitability of this software for any purpose. It is * | |
16 | //* provided "as is" without express or implied warranty. * | |
17 | //************************************************************************** | |
18 | ||
19 | /// @file AliHLTTRDTriggerComponent.cxx | |
20 | /// @author Felix Rettig, Stefan Kirsch | |
21 | /// @date 2012-08-16 | |
22 | /// @brief | |
23 | ||
24 | #include <cstdlib> | |
25 | #include "TSystem.h" | |
26 | #include "TObjArray.h" | |
27 | #include "TObjString.h" | |
28 | #include "TFile.h" | |
29 | #include "TH1I.h" | |
30 | #include "TH2I.h" | |
31 | #include "AliESDtrack.h" | |
32 | #include "AliESDEvent.h" | |
33 | #include "AliHLTTRDDefinitions.h" | |
34 | #include "AliHLTTriggerDecision.h" | |
35 | #include "AliHLTDomainEntry.h" | |
36 | #include "AliHLTLogging.h" | |
37 | #include "AliHLTErrorGuard.h" | |
38 | #include "AliHLTGlobalBarrelTrack.h" | |
39 | #include "AliESDtrackCuts.h" | |
40 | //#include "AliHLTESDTrackCuts.h" | |
41 | #include "AliGeomManager.h" | |
42 | #include "AliHLTComponentBenchmark.h" | |
43 | #include "AliHLTTRDTriggerComponent.h" | |
44 | #include "AliTRDonlineTrackingDataContainer.h" | |
45 | #include "AliTRDpadPlane.h" | |
46 | ||
47 | ClassImp(AliHLTTRDTriggerComponent) | |
48 | ||
49 | #define LogError( ... ) { HLTError(__VA_ARGS__); if (fDebugLevel >= 1) { DbgLog("ERROR", __VA_ARGS__); } } | |
50 | #define LogInfo( ... ) { HLTInfo(__VA_ARGS__); if (fDebugLevel >= 1) { DbgLog("INFO", __VA_ARGS__); } } | |
51 | #define LogInspect( ... ) { HLTDebug(__VA_ARGS__); if (fDebugLevel >= 1) { DbgLog("INSPECT", __VA_ARGS__); } } | |
52 | #define LogDebug( ... ) { if (fDebugLevel >= 1) { HLTInfo(__VA_ARGS__); DbgLog("DEBUG", __VA_ARGS__); } } | |
53 | ||
54 | AliHLTTRDTriggerComponent::AliHLTTRDTriggerComponent() : | |
55 | AliHLTTrigger(), | |
56 | fName(), | |
57 | fRefTrackSelectionEtaLimit(0.9), | |
58 | fRefTrackSelectionVertexXYLimit(20.), | |
59 | fRefTrackSelectionVertexZLimit(30.), | |
60 | fRefTrackSelectionPtThreshold(0.7), | |
61 | fMatchRatingThreshold(0.25), | |
62 | fElectronTriggerPtThresholdHSE(3.), | |
63 | fElectronTriggerPIDThresholdHSE(144), | |
64 | fElectronTriggerPtThresholdHQU(2.), | |
65 | fElectronTriggerPIDThresholdHQU(164), | |
66 | fApplyRefTrackCuts(kFALSE), | |
67 | fElectronTriggerOnL1TrgOnly(kFALSE), | |
68 | fHistoMode(1), | |
69 | fDebugLevel(0), | |
70 | fExtendedHistos(kFALSE), | |
71 | fEventRendering(kFALSE), | |
72 | fPushHistos(kFALSE), | |
73 | fWriteHistos(kFALSE), | |
74 | fEventId(fgkInvalidEventId), | |
75 | fRunNumber(-1), | |
76 | fChunkId(NULL), | |
77 | fSectorsWithData(0), | |
78 | fIsMinBiasEvent(kFALSE), | |
79 | fIsTRDElectronEvent(kFALSE), | |
80 | fESDtracksPresent(kFALSE), | |
81 | fHLTtracksPresent(kFALSE), | |
82 | fTRDGeometry(NULL), | |
83 | fEsdEvent(NULL), | |
84 | fTrackingData(NULL), | |
85 | fHLTTracks(NULL), | |
86 | fRefTrackCuts(NULL), | |
87 | #ifdef __TRDHLTDEBUG | |
88 | fEventDisplay(NULL), | |
89 | fBenchmark(NULL), | |
90 | #endif | |
91 | fHistArray(NULL), | |
92 | fHistMatchRating(NULL), | |
93 | fHistMatchRatingByPt(NULL), | |
94 | fHistMatchRatingByPid(NULL), | |
95 | fHistTrackPt(NULL), | |
96 | fHistTrackPtMatched(NULL), | |
97 | fHistTrackPtCorr(NULL), | |
98 | fHistTrackPid(NULL), | |
99 | fHistTrackPidMatched(NULL), | |
100 | fHistElectronCandidatePt(NULL), | |
101 | fHistElectronCandidateMatchedPt(NULL), | |
102 | fHistElectronCandidatePid(NULL), | |
103 | fHistElectronCandidateMatchedPid(NULL), | |
104 | fHistRefTrackPid(NULL), | |
105 | fHistMatchedRefTrackPid(NULL), | |
106 | fHistPIDvsTruncPID(NULL), | |
107 | fHistElectronFalsePIDvsTruncPID(NULL), | |
108 | fHistElectronConfirmedPIDvsTruncPID(NULL), | |
109 | fHistTrackMatchingCombinations(NULL), | |
110 | fHistElectronTriggerBaseMinBias(NULL), | |
111 | fHistElectronTriggerBaseTrdL1(NULL) | |
112 | { | |
113 | } | |
114 | ||
115 | const char* AliHLTTRDTriggerComponent::fgkDefaultOCDBEntry = "HLT/ConfigHLT/TRDTrigger"; | |
116 | const char* AliHLTTRDTriggerComponent::fgkTriggerDecisionElectronHSE = "TRD-ELECTRON-HSE"; | |
117 | const char* AliHLTTRDTriggerComponent::fgkTriggerDecisionElectronHQU = "TRD-ELECTRON-HQU"; | |
118 | ||
119 | AliHLTTRDTriggerComponent::~AliHLTTRDTriggerComponent() | |
120 | { | |
121 | } | |
122 | ||
123 | const char* AliHLTTRDTriggerComponent::GetTriggerName() const | |
124 | { | |
125 | if (!fName.IsNull()) | |
126 | return fName.Data(); | |
127 | else | |
128 | return "TRDTriggerComponent"; | |
129 | } | |
130 | ||
131 | AliHLTComponent* AliHLTTRDTriggerComponent::Spawn() | |
132 | { | |
133 | return new AliHLTTRDTriggerComponent; | |
134 | } | |
135 | ||
136 | int AliHLTTRDTriggerComponent::DoInit(int argc, const char** argv) | |
137 | { | |
138 | int iResult = 0; | |
139 | ||
140 | do { | |
141 | ||
142 | fChunkId = new TString("XXXXX"); | |
143 | #ifdef __TRDHLTDEBUG | |
144 | if (fChunkId){ | |
145 | // chunk identification for debug output | |
146 | *fChunkId = gSystem->WorkingDirectory(); | |
147 | fChunkId->Remove(0, fChunkId->Last('/') + 1); | |
148 | if (fChunkId->Contains("hlt_trd_comp")) | |
149 | *fChunkId = "L"; | |
150 | } else { | |
151 | iResult = -ENOMEM; | |
152 | break; | |
153 | } | |
154 | #endif | |
155 | ||
156 | AliGeomManager::LoadGeometry(); | |
157 | ||
158 | fTRDGeometry = new AliTRDgeometry(); | |
159 | if (!fTRDGeometry) { | |
160 | iResult = -ENOMEM; | |
161 | break; | |
162 | } | |
163 | ||
164 | fTrackingData = new AliTRDonlineTrackingDataContainer(); | |
165 | if (!fTrackingData) { | |
166 | iResult = -ENOMEM; | |
167 | break; | |
168 | } | |
169 | ||
170 | fHLTTracks = new vector<AliHLTGlobalBarrelTrack>; | |
171 | if (!fHLTTracks) { | |
172 | iResult = -ENOMEM; | |
173 | break; | |
174 | } | |
175 | ||
176 | fRefTrackCuts = new AliESDtrackCuts("AliESDtrackCuts", "No track cuts"); | |
177 | if (fRefTrackCuts){ | |
178 | // fRefTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2011(kTRUE, 0); | |
179 | fRefTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts(); | |
180 | fRefTrackCuts->SetEtaRange(-0.8, 0.8); | |
181 | } else { | |
182 | iResult = -ENOMEM; | |
183 | break; | |
184 | } | |
185 | ||
186 | fHistArray = new TObjArray(25); | |
187 | if(!fHistArray) | |
188 | return -ENOMEM; | |
189 | fHistArray->SetOwner(kTRUE); | |
190 | ||
191 | fHistMatchRating = new TH1I("trdtrg_match_rating", "Match rating distribution HLT tracks vs. GTU track;match rating;abundance;", | |
192 | 1./0.02, 0., 1.); | |
193 | fHistArray->AddLast(fHistMatchRating); | |
194 | ||
195 | fHistMatchRatingByPt = new TH2I("trdtrg_match_rating_pt", "Match rating distribution HLT tracks vs. GTU track by p_{T};match rating m;p_{T}^{TRD,online} (GeV/c)", | |
196 | 1./0.02, 0., 1., 20./0.02, 0., 20.); | |
197 | fHistArray->AddLast(fHistMatchRatingByPt); | |
198 | ||
199 | fHistMatchRatingByPid = new TH2I("trdtrg_match_rating_pid", "Match rating distribution HLT tracks vs. GTU track by PID;match rating m;PID (a.u.)", | |
200 | 1./0.02, 0., 1., 256, 0, 256); | |
201 | fHistArray->AddLast(fHistMatchRatingByPid); | |
202 | ||
203 | fHistTrackPt = new TH1I("trdtrg_track_pt", "GTU track p_{T};GTU track p_{T} (GeV/c);abundance", 200, -20, 20.); | |
204 | fHistArray->AddLast(fHistTrackPt); | |
205 | ||
206 | fHistTrackPtMatched = new TH1I("trdtrg_matched_track_pt", "matched GTU track p_{T};GTU track p_{T} (GeV/c);abundance", 200, -20, 20.); | |
207 | fHistArray->AddLast(fHistTrackPtMatched); | |
208 | ||
209 | fHistTrackPtCorr = new TH2I("trdtrg_track_pt_corr", "HLT vs. GTU track p_{T};p_{T}^{HLT} (GeV/c);p_{T}^{GTU} (GeV/c)", | |
210 | 400, -40., 40., 400, -40., 40.); | |
211 | fHistArray->AddLast(fHistTrackPtCorr); | |
212 | ||
213 | fHistTrackPid = new TH1I("trdtrg_track_pid", "GTU track PID;GTU track PID (a.u.);abundance", 256, 0, 256); | |
214 | fHistArray->AddLast(fHistTrackPid); | |
215 | ||
216 | fHistTrackPidMatched = new TH1I("trdtrg_matched_track_pid", "matched GTU track PID;GTU track PID (a.u.);abundance", 256, 0, 256); | |
217 | fHistArray->AddLast(fHistTrackPidMatched); | |
218 | ||
219 | fHistTrackMatchingCombinations = new TH2I("trdtrg_matching_combinations", "HLT-GTU track matching combinations;set;combinations to check;", | |
220 | 10, 0, 10, 2000, 0, 10000); | |
221 | fHistArray->AddLast(fHistTrackMatchingCombinations); | |
222 | ||
223 | fHistElectronTriggerBaseMinBias = new TH1I("trdtrg_electron_trigger_base_mb", "min. bias base numbers for electron trigger analysis;set;abundance;", | |
224 | 10, 0, 10); | |
225 | fHistElectronTriggerBaseMinBias->GetXaxis()->SetBinLabel(1, "min. bias events"); | |
226 | fHistElectronTriggerBaseMinBias->GetXaxis()->SetBinLabel(2, "TRD L1 electron triggered"); | |
227 | fHistElectronTriggerBaseMinBias->GetXaxis()->SetBinLabel(3, "TRD HLT electron triggered"); | |
228 | fHistArray->AddLast(fHistElectronTriggerBaseMinBias); | |
229 | ||
230 | fHistElectronTriggerBaseTrdL1 = new TH1I("trdtrg_electron_trigger_base_l1", "TRD L1 base numbers for electron trigger analysis;set;abundance;", | |
231 | 10, 0, 10); | |
232 | fHistElectronTriggerBaseTrdL1->GetXaxis()->SetBinLabel(1, "TRD L1 electron triggered"); | |
233 | fHistElectronTriggerBaseTrdL1->GetXaxis()->SetBinLabel(2, "TRD HLT electron triggered"); | |
234 | fHistArray->AddLast(fHistElectronTriggerBaseTrdL1); | |
235 | ||
236 | fHistElectronCandidatePt = new TH1I("trdtrg_electron_candidate_pt", "GTU electron candidate p_{T};GTU track p_{T} (GeV/c);abundance", 200, -20, 20.); | |
237 | fHistArray->AddLast(fHistElectronCandidatePt); | |
238 | ||
239 | fHistElectronCandidateMatchedPt = new TH1I("trdtrg_electron_candidate_matched_pt", "matching GTU electron candidate p_{T};GTU track p_{T} (GeV/c);abundance", 200, -20, 20.); | |
240 | fHistArray->AddLast(fHistElectronCandidateMatchedPt); | |
241 | ||
242 | fHistElectronCandidatePid = new TH1I("trdtrg_electron_candidate_pid", "GTU electron candidate PID;GTU track PID (a.u.);abundance", 256, 0, 256); | |
243 | fHistArray->AddLast(fHistElectronCandidatePid); | |
244 | ||
245 | fHistElectronCandidateMatchedPid = new TH1I("trdtrg_electron_candidate_matched_pid", "matching GTU electron candidate PID;GTU track PID (a.u.);abundance", 256, 0, 256); | |
246 | fHistArray->AddLast(fHistElectronCandidateMatchedPid); | |
247 | ||
248 | } while (0); | |
249 | ||
250 | if (iResult < 0){ | |
251 | ||
252 | if (fHistArray) delete fHistArray; | |
253 | fHistArray = NULL; | |
254 | ||
255 | if (fTRDGeometry) delete fTRDGeometry; | |
256 | fTRDGeometry = NULL; | |
257 | ||
258 | if (fRefTrackCuts) delete fRefTrackCuts; | |
259 | fRefTrackCuts = NULL; | |
260 | ||
261 | if (fHLTTracks) delete fHLTTracks; | |
262 | fHLTTracks = NULL; | |
263 | ||
264 | if (fTrackingData) delete fTrackingData; | |
265 | fTrackingData = NULL; | |
266 | ||
267 | if (fChunkId) delete fChunkId; | |
268 | fChunkId = NULL; | |
269 | ||
270 | } | |
271 | ||
272 | // check if the -triggername argument is used, the trigger name determines the following initialization | |
273 | vector<const char*> remainingArgs; | |
274 | for (int i = 0; i < argc; ++i) { | |
275 | if (strcmp(argv[i], "-triggername") == 0) { | |
276 | if (++i < argc){ | |
277 | fName = argv[i]; | |
278 | } else { | |
279 | LogError("invalid parameter for argument '-triggername', string expected"); | |
280 | return -EINVAL; | |
281 | } | |
282 | continue; | |
283 | } | |
284 | remainingArgs.push_back(argv[i]); | |
285 | } | |
286 | ||
287 | TString cdbPath; | |
288 | if (!fName.IsNull()) { | |
289 | cdbPath = "HLT/ConfigHLT/"; | |
290 | cdbPath += fName; | |
291 | } else { | |
292 | cdbPath = fgkDefaultOCDBEntry; | |
293 | } | |
294 | ||
295 | LogInfo("cdbPath = <%s>", cdbPath.Data()); | |
296 | iResult = ConfigureFromCDBObject(cdbPath); | |
297 | ||
298 | if (iResult >= 0 && argc > 0) | |
299 | iResult = ConfigureFromArgumentString(remainingArgs.size(), &(remainingArgs[0])); | |
300 | ||
301 | return iResult; | |
302 | } | |
303 | ||
304 | int AliHLTTRDTriggerComponent::DoDeinit() | |
305 | { | |
306 | ||
307 | #ifdef __TRDHLTDEBUG | |
308 | if (fEventDisplay) delete fEventDisplay; | |
309 | fEventDisplay = NULL; | |
310 | #endif | |
311 | ||
312 | if ((fHistoMode == 1) && (fWriteHistos)){ | |
313 | TFile out("trg_out/trg_hists.root", "RECREATE"); | |
314 | if (!out.IsZombie()) { | |
315 | out.cd(); | |
316 | UInt_t numHists = fHistArray->GetEntries(); | |
317 | for (UInt_t iHist = 0; iHist < numHists; ++iHist) | |
318 | if (fHistArray->At(iHist)) | |
319 | fHistArray->At(iHist)->Write(); | |
320 | out.Close(); | |
321 | } | |
322 | } | |
323 | ||
324 | if (fHistArray) delete fHistArray; | |
325 | fHistArray = NULL; | |
326 | ||
327 | if (fRefTrackCuts) delete fRefTrackCuts; | |
328 | fRefTrackCuts = NULL; | |
329 | ||
330 | if (!fHLTTracks) delete fHLTTracks; | |
331 | fHLTTracks = NULL; | |
332 | ||
333 | if (fTrackingData) delete fTrackingData; | |
334 | fTrackingData = NULL; | |
335 | ||
336 | if (fTRDGeometry) delete fTRDGeometry; | |
337 | fTRDGeometry = NULL; | |
338 | ||
339 | if (fChunkId) delete fChunkId; | |
340 | fChunkId = NULL; | |
341 | ||
342 | return 0; | |
343 | } | |
344 | ||
345 | int AliHLTTRDTriggerComponent::Reconfigure(const char* cdbEntry, const char* /*chainId*/) | |
346 | { | |
347 | TString cdbPath; | |
348 | if (!cdbEntry || cdbEntry[0] == 0) { | |
349 | if (!fName.IsNull()) { | |
350 | cdbPath = "HLT/ConfigHLT/"; | |
351 | cdbPath += fName; | |
352 | } else { | |
353 | cdbPath = fgkDefaultOCDBEntry; | |
354 | } | |
355 | } else | |
356 | cdbPath = cdbEntry; | |
357 | ||
358 | return ConfigureFromCDBObject(cdbPath); | |
359 | } | |
360 | ||
361 | int AliHLTTRDTriggerComponent::ReadPreprocessorValues(const char* /*modules*/) | |
362 | { | |
363 | return 0; | |
364 | } | |
365 | ||
366 | Int_t AliHLTTRDTriggerComponent::ConfigureFromCDBObject(TString /*cdbPath*/) | |
367 | { | |
368 | return 0; | |
369 | } | |
370 | ||
371 | int AliHLTTRDTriggerComponent::ScanConfigurationArgument(int argc, const char** argv) | |
372 | { | |
373 | ||
374 | if (argc <= 0) | |
375 | return 0; | |
376 | ||
377 | UShort_t iArg = 0; | |
378 | TString argument(argv[iArg]); | |
379 | ||
380 | if (!argument.CompareTo("-match-sel-eta")){ | |
381 | if (++iArg >= argc) return -EPROTO; | |
382 | argument = argv[iArg]; | |
383 | fRefTrackSelectionEtaLimit = argument.Atof(); | |
384 | LogInfo("ref track selection eta limit set to %.1f.", fRefTrackSelectionEtaLimit); | |
385 | return 2; | |
386 | } | |
387 | ||
388 | if (!argument.CompareTo("-match-sel-vxy")){ | |
389 | if (++iArg >= argc) return -EPROTO; | |
390 | argument = argv[iArg]; | |
391 | fRefTrackSelectionVertexXYLimit = argument.Atof(); | |
392 | LogInfo("ref track selection vertex xy limit set to %.1f.", fRefTrackSelectionVertexXYLimit); | |
393 | return 2; | |
394 | } | |
395 | ||
396 | if (!argument.CompareTo("-match-sel-vz")){ | |
397 | if (++iArg >= argc) return -EPROTO; | |
398 | argument = argv[iArg]; | |
399 | fRefTrackSelectionVertexZLimit = argument.Atof(); | |
400 | LogInfo("ref track selection vertex z limit set to %.1f.", fRefTrackSelectionVertexZLimit); | |
401 | return 2; | |
402 | } | |
403 | ||
404 | if (!argument.CompareTo("-match-sel-pt")){ | |
405 | if (++iArg >= argc) return -EPROTO; | |
406 | argument = argv[iArg]; | |
407 | fRefTrackSelectionPtThreshold = argument.Atof(); | |
408 | LogInfo("ref track selection pt threshold set to %.1f GeV/c.", fRefTrackSelectionPtThreshold); | |
409 | return 2; | |
410 | } | |
411 | ||
412 | if (!argument.CompareTo("-match-rating")){ | |
413 | if (++iArg >= argc) return -EPROTO; | |
414 | argument = argv[iArg]; | |
415 | fMatchRatingThreshold = argument.Atof(); | |
416 | LogInfo("match rating threshold set to %.2f GeV/c.", fMatchRatingThreshold); | |
417 | return 2; | |
418 | } | |
419 | ||
420 | if (!argument.CompareTo("-histo-mode")){ | |
421 | if (++iArg >= argc) return -EPROTO; | |
422 | argument = argv[iArg]; | |
423 | fHistoMode = argument.Atoi(); | |
424 | LogInfo("histo mode set to %d", fHistoMode); | |
425 | return 2; | |
426 | } | |
427 | ||
428 | if (!argument.CompareTo("-trghse")){ | |
429 | LogInfo("TRD HLT electron trigger HSE enabled."); | |
430 | return 1; | |
431 | } | |
432 | ||
433 | if (!argument.CompareTo("-hse-pt")){ | |
434 | if (++iArg >= argc) return -EPROTO; | |
435 | argument = argv[iArg]; | |
436 | fElectronTriggerPtThresholdHSE = argument.Atof(); | |
437 | LogInfo("pt threshold for HSE trigger set to %.1f GeV/c.", fElectronTriggerPtThresholdHSE); | |
438 | return 2; | |
439 | } | |
440 | ||
441 | if (!argument.CompareTo("-hse-pid")){ | |
442 | if (++iArg >= argc) return -EPROTO; | |
443 | argument = argv[iArg]; | |
444 | fElectronTriggerPIDThresholdHSE = argument.Atof(); | |
445 | LogInfo("PID threshold for HSE trigger set to %d.", fElectronTriggerPIDThresholdHSE); | |
446 | return 2; | |
447 | } | |
448 | ||
449 | if (!argument.CompareTo("-trghqu")){ | |
450 | LogInfo("TRD HLT electron trigger HQU enabled."); | |
451 | return 1; | |
452 | } | |
453 | ||
454 | if (!argument.CompareTo("-hqu-pt")){ | |
455 | if (++iArg >= argc) return -EPROTO; | |
456 | argument = argv[iArg]; | |
457 | fElectronTriggerPtThresholdHQU = argument.Atof(); | |
458 | LogInfo("pt threshold for HQU trigger set to %.1f GeV/c.", fElectronTriggerPtThresholdHQU); | |
459 | return 2; | |
460 | } | |
461 | ||
462 | if (!argument.CompareTo("-hqu-pid")){ | |
463 | if (++iArg >= argc) return -EPROTO; | |
464 | argument = argv[iArg]; | |
465 | fElectronTriggerPIDThresholdHQU = argument.Atof(); | |
466 | LogInfo("PID threshold for HQU trigger set to %.1f GeV/c.", fElectronTriggerPIDThresholdHQU); | |
467 | return 2; | |
468 | } | |
469 | ||
470 | if (!argument.CompareTo("-l1-only")){ | |
471 | LogInfo("evaluation of electron trigger only for events with TRD L1 electron trigger enabled."); | |
472 | fElectronTriggerOnL1TrgOnly = kTRUE; | |
473 | return 1; | |
474 | } | |
475 | ||
476 | if (!argument.CompareTo("-histo-ext")){ | |
477 | LogInfo("extended histogramming enabled."); | |
478 | fExtendedHistos = kTRUE; // enable extended histogramming, for debugging/tuning only! | |
479 | ||
480 | if (!fHistRefTrackPid){ | |
481 | fHistRefTrackPid = new TH2I("trdtrg_reftrack_pid", "TPC dE/dx by p_{T} for matching GTU tracks;p_{T} (GeV/c);dE/dx (a. u.)", | |
482 | 500, 0., 10., 500, 0., 500); | |
483 | fHistArray->AddLast(fHistRefTrackPid); | |
484 | } | |
485 | ||
486 | if (!fHistMatchedRefTrackPid){ | |
487 | fHistMatchedRefTrackPid = new TH2I("trdtrg_matched_reftrack_pid", "TPC dE/dx by p_{T} for matching GTU tracks;p_{T} (GeV/c);dE/dx (a. u.)", | |
488 | 500, 0., 10., 500, 0., 500); | |
489 | fHistArray->AddLast(fHistMatchedRefTrackPid); | |
490 | } | |
491 | ||
492 | if (!fHistPIDvsTruncPID){ | |
493 | fHistPIDvsTruncPID = new TH2I("trdtrg_pid_trunc_pid", "GTU track PID vs. truncated PID;PID (a.u.);truncated PID (a. u.)", | |
494 | 256, 0., 256., 256, 0., 256.); | |
495 | fHistArray->AddLast(fHistPIDvsTruncPID); | |
496 | } | |
497 | ||
498 | if (!fHistElectronFalsePIDvsTruncPID){ | |
499 | fHistElectronFalsePIDvsTruncPID = new TH2I("trdtrg_electron_false_pid_trunc_pid", "false electron PID vs. truncated PID;PID (a.u.);truncated PID (a. u.)", | |
500 | 256, 0., 256., 256, 0., 256.); | |
501 | fHistArray->AddLast(fHistElectronFalsePIDvsTruncPID); | |
502 | } | |
503 | ||
504 | if (!fHistElectronConfirmedPIDvsTruncPID){ | |
505 | fHistElectronConfirmedPIDvsTruncPID = new TH2I("trdtrg_electron_confirmed_pid_trunc_pid", "confirmed electron PID vs. truncated PID;PID (a.u.);truncated PID (a. u.)", | |
506 | 256, 0., 256., 256, 0., 256.); | |
507 | fHistArray->AddLast(fHistElectronConfirmedPIDvsTruncPID); | |
508 | } | |
509 | ||
510 | return 1; | |
511 | } | |
512 | ||
513 | if (!argument.CompareTo("-ref-cuts")){ | |
514 | LogInfo("ref track cuts for matching enabled."); | |
515 | fApplyRefTrackCuts = kTRUE; | |
516 | return 1; | |
517 | } | |
518 | ||
519 | if (!argument.CompareTo("-push-histograms")){ | |
520 | LogInfo("pushing of histograms enabled."); | |
521 | fPushHistos = kTRUE; // enable histogram pushing, for debugging/tuning only! | |
522 | return 1; | |
523 | } | |
524 | ||
525 | if (!argument.CompareTo("-write-histograms")){ | |
526 | LogInfo("writing of histograms enabled."); | |
527 | fWriteHistos = kTRUE; // enable histogram writing, for debugging/tuning only! | |
528 | return 1; | |
529 | } | |
530 | ||
531 | if (!argument.CompareTo("-render")){ | |
532 | #ifdef __TRDHLTDEBUG | |
533 | LogInfo("rendering of interesting events enabled."); | |
534 | if (!fEventDisplay) | |
535 | fEventDisplay = new AliTRDtrackingEventDisplay(); | |
536 | LogInfo("event rendering activated. this is for debugging only!"); | |
537 | fEventRendering = kTRUE; // enable event rendering, for debugging/tuning only! | |
538 | #endif | |
539 | return 1; | |
540 | } | |
541 | ||
542 | if (!argument.CompareTo("-debug")){ | |
543 | if (++iArg >= argc) return -EPROTO; | |
544 | argument = argv[iArg]; | |
545 | fDebugLevel = argument.Atoi(); | |
546 | LogInfo("debug level set to %d.", fDebugLevel); | |
547 | return 2; | |
548 | } | |
549 | ||
550 | if (!argument.CompareTo("-ref-pt")){ | |
551 | if (++iArg >= argc) return -EPROTO; | |
552 | argument = argv[iArg]; | |
553 | Float_t minPt, maxPt; | |
554 | fRefTrackCuts->GetPtRange(minPt, maxPt); | |
555 | maxPt = argument.Atof(); | |
556 | fRefTrackCuts->SetPtRange(minPt, maxPt); | |
557 | LogInfo("ref track pt range set to %.1f .. %.1f GeV/c.", minPt, maxPt); | |
558 | return 2; | |
559 | } | |
560 | ||
561 | if (!argument.CompareTo("-ref-tdca")){ | |
562 | if (++iArg >= argc) return -EPROTO; | |
563 | argument = argv[iArg]; | |
564 | fRefTrackCuts->SetMaxDCAToVertexXY(argument.Atof()); | |
565 | LogInfo("ref track DCA transverse threshold set to %.3f", argument.Atof()); | |
566 | return 2; | |
567 | } | |
568 | ||
569 | if (!argument.CompareTo("-ref-ldca")){ | |
570 | if (++iArg >= argc) return -EPROTO; | |
571 | argument = argv[iArg]; | |
572 | fRefTrackCuts->SetMaxDCAToVertexZ(argument.Atof()); | |
573 | LogInfo("ref track longitudinal DCA threshold set to %.3f", argument.Atof()); | |
574 | return 2; | |
575 | } | |
576 | ||
577 | return 1; | |
578 | } | |
579 | ||
580 | void AliHLTTRDTriggerComponent::ScanTriggerClasses(const char* firedTriggerClasses) { | |
581 | ||
582 | fIsMinBiasEvent = kFALSE; | |
583 | ||
584 | TString trg(firedTriggerClasses); | |
585 | if (trg.Index("CINT7WU-S-NOPF-ALL") >= 0) | |
586 | fIsMinBiasEvent = kTRUE; | |
587 | ||
588 | if (trg.Index("CINT8WU-S-NOPF-ALL") >= 0) | |
589 | fIsMinBiasEvent = kTRUE; | |
590 | ||
591 | } | |
592 | ||
593 | Bool_t AliHLTTRDTriggerComponent::CheckRefTrackCuts(AliESDtrack* track){ | |
594 | ||
595 | // cuts by cut class | |
596 | if (fApplyRefTrackCuts) | |
597 | return fRefTrackCuts->AcceptTrack(track) ? kTRUE : kFALSE; | |
598 | ||
599 | // simple custom cuts | |
600 | Float_t dcaToVertexXY, dcaToVertexZ; | |
601 | track->GetImpactParametersTPC(dcaToVertexXY, dcaToVertexZ); | |
602 | // LogDebug("IMPACT TPC %.4f %.4f", dcaToVertexXY, dcaToVertexZ); | |
603 | ||
604 | if ((dcaToVertexXY > fRefTrackSelectionVertexXYLimit) || (dcaToVertexZ >= fRefTrackSelectionVertexZLimit)) | |
605 | return kFALSE; | |
606 | ||
607 | if (TMath::Abs(track->Eta()) > fRefTrackSelectionEtaLimit) | |
608 | return kFALSE; | |
609 | ||
610 | return kTRUE; | |
611 | } | |
612 | ||
613 | ||
614 | void AliHLTTRDTriggerComponent::DbgLog(const char* prefix, ...){ | |
615 | #ifdef __TRDHLTDEBUG | |
616 | AliHLTEventID_t eventNumber = fEventId; | |
617 | printf("TRDHLTGM %s-%s-%s: [TRG] %s", | |
618 | (fRunNumber >= 0) ? Form("%06d", fRunNumber) : "XXXXXX", | |
619 | fChunkId->Data(), | |
620 | (eventNumber != fgkInvalidEventId) ? Form("%05llu", eventNumber) : "XXXXX", | |
621 | (strlen(prefix) > 0) ? Form("<%s> ", prefix) : ""); | |
622 | #endif | |
623 | va_list args; | |
624 | va_start(args, prefix); | |
625 | char* fmt = va_arg(args, char*); | |
626 | vprintf(fmt, args); | |
627 | printf("\n"); | |
628 | va_end(args); | |
629 | } | |
630 | ||
631 | int AliHLTTRDTriggerComponent::PrepareESDData(){ | |
632 | ||
633 | int result = 0; | |
634 | fESDtracksPresent = kFALSE; | |
635 | ||
636 | // check access to ESD event data | |
637 | const TObject* obj = GetFirstInputObject(kAliHLTAllDataTypes, "AliESDEvent"); | |
638 | fEsdEvent = dynamic_cast<AliESDEvent*>(const_cast<TObject*>(obj)); | |
639 | ||
640 | if (fEsdEvent) { | |
641 | fEsdEvent->GetStdContent(); | |
642 | fRunNumber = fEsdEvent->GetRunNumber(); | |
643 | fESDtracksPresent = kTRUE; | |
644 | unsigned int numTracks = fEsdEvent->GetNumberOfTracks(); | |
645 | ||
646 | // process trigger classes | |
647 | ScanTriggerClasses(fEsdEvent->GetFiredTriggerClasses().Data()); | |
648 | ||
649 | LogDebug("ESD event data: %d ESD tracks [event num: %d] [minbias: %d (fired:%s)]", | |
650 | numTracks, fEsdEvent->GetEventNumberInFile(), fIsMinBiasEvent, fEsdEvent->GetFiredTriggerClasses().Data()); | |
651 | ||
652 | // process ESD tracks | |
653 | if ((fDebugLevel >= 2) || fExtendedHistos){ | |
654 | AliESDtrack* esdTrack; | |
655 | TString paramStr(""); | |
656 | for (unsigned int iTrack = 0; iTrack < numTracks; ++iTrack) { | |
657 | esdTrack = fEsdEvent->GetTrack(iTrack); | |
658 | ||
659 | if (fExtendedHistos){ | |
660 | fHistRefTrackPid->Fill(esdTrack->Pt(), esdTrack->GetTPCsignal()); | |
661 | } | |
662 | ||
663 | if (fDebugLevel >= 2){ | |
664 | paramStr = ""; | |
665 | if (esdTrack){ | |
666 | ||
667 | if (esdTrack->GetInnerParam()) | |
668 | paramStr += "I"; | |
669 | ||
670 | if (esdTrack->GetOuterParam()) | |
671 | paramStr += "O"; | |
672 | ||
673 | LogDebug("ESD track %4d - pt: %+.2fGeV/c [params: S%s]", iTrack, esdTrack->GetSignedPt(), paramStr.Data()); | |
674 | } else | |
675 | LogError("ESD data for track %d invalid", iTrack); | |
676 | } | |
677 | ||
678 | } // loop over ESD tracks | |
679 | } | |
680 | ||
681 | fTrackingData->SetGtuPtMultiplierFromMagField(fEsdEvent->GetMagneticField()); // used for sign correction | |
682 | ||
683 | result = 1; | |
684 | ||
685 | } | |
686 | ||
687 | return result; | |
688 | ||
689 | } | |
690 | ||
691 | int AliHLTTRDTriggerComponent::PrepareHLTData(){ | |
692 | ||
693 | int iResult = 0; | |
694 | UInt_t numHLTTracks = 0; | |
695 | fHLTtracksPresent = kFALSE; | |
696 | ||
697 | if (!fHLTTracks){ | |
698 | LogError("HLT track vector instance not available."); | |
699 | return 0; | |
700 | } | |
701 | ||
702 | fHLTTracks->clear(); | |
703 | ||
704 | for (const AliHLTComponentBlockData *pBlock = GetFirstInputBlock(kAliHLTDataTypeTrack | kAliHLTDataOriginTPC); | |
705 | pBlock != NULL; pBlock = GetNextInputBlock()) { | |
706 | LogDebug("#hlttrk - data block with HLT raw tracks received"); | |
707 | fHLTtracksPresent = kTRUE; | |
708 | ||
709 | // vector<AliHLTGlobalBarrelTrack> hltTracks; | |
710 | if ((iResult = AliHLTGlobalBarrelTrack::ConvertTrackDataArray(reinterpret_cast<const AliHLTTracksData*>(pBlock->fPtr), pBlock->fSize, *fHLTTracks)) >= 0) { | |
711 | for (vector<AliHLTGlobalBarrelTrack>::iterator element = fHLTTracks->begin(); | |
712 | element != fHLTTracks->end(); element++) { | |
713 | ||
714 | numHLTTracks++; | |
715 | //## AliHLTGlobalBarrelTrack -> AliKalmanTrack -> AliExternalTrackParam | |
716 | ||
717 | } // loop over HLT tracks | |
718 | } | |
719 | } // loop over data blocks | |
720 | ||
721 | LogDebug("#hlttrk - %d HLT raw tracks found", numHLTTracks); | |
722 | ||
723 | return iResult; | |
724 | ||
725 | } | |
726 | ||
727 | int AliHLTTRDTriggerComponent::PrepareTRDData() { | |
728 | ||
729 | int result = 1; | |
730 | ||
731 | fSectorsWithData = 0; | |
732 | fTrackingData->Clear(); | |
733 | for (const AliHLTComponentBlockData* datablock = GetFirstInputBlock(AliHLTTRDDefinitions::fgkOnlineDataType); | |
734 | datablock != NULL; | |
735 | datablock = GetNextInputBlock()) | |
736 | { | |
737 | fSectorsWithData |= datablock->fSpecification; | |
738 | fTrackingData->Decompress(datablock->fPtr, datablock->fSize, kTRUE); | |
739 | } | |
740 | ||
741 | fTrackingData->SetGtuPtMultiplierFromMagField(GetBz()); // used for sign correction | |
742 | ||
743 | fTrackingData->PrintSummary("trigger component"); | |
744 | ||
745 | return result; | |
746 | ||
747 | } | |
748 | ||
749 | int AliHLTTRDTriggerComponent::MatchTRDTracksESD(){ | |
750 | ||
751 | if (!fEsdEvent) { | |
752 | LogError("ESD event data not available in MatchTRDTracks()."); | |
753 | return 0; | |
754 | } | |
755 | ||
756 | if (fHistoMode == 0){ | |
757 | UInt_t numHists = fHistArray->GetEntries(); | |
758 | for (UInt_t iHist = 0; iHist < numHists; ++iHist) | |
759 | if (fHistArray->At(iHist)) | |
ca892663 | 760 | ((TH1*) (fHistArray->At(iHist)))->Reset(); |
c7b7f445 | 761 | } |
762 | ||
763 | int result = 1; | |
764 | unsigned int numRefTracks = fEsdEvent->GetNumberOfTracks(); | |
765 | unsigned int numGtuTracks; | |
766 | Int_t refTrackIndices[fkMaxRefTracksPerStack]; | |
767 | UInt_t refTrackCount = 0; | |
768 | Bool_t isRelevant; | |
769 | Double_t distY, distZ; | |
770 | Double_t matchRating; | |
771 | Double_t bestMatchRating; | |
772 | Int_t bestMatchRefIndex; | |
773 | AliESDtrack* refTrack = NULL; | |
774 | UInt_t numComparisonsDone = 0; | |
775 | UInt_t numUnmatchedTracks = 0; | |
776 | UInt_t numMatchedTracks = 0; | |
777 | Double_t magField = fEsdEvent->GetMagneticField(); | |
778 | ||
779 | for (UShort_t iStack = 0; iStack < fkTRDStacks; ++iStack) { | |
780 | numGtuTracks = fTrackingData->GetNumTracks(iStack); | |
781 | refTrackCount = 0; | |
782 | ||
783 | // preselect ESD relevant ESD tracks | |
784 | for (UInt_t iRefTrack = 0; iRefTrack < numRefTracks; ++iRefTrack) { | |
785 | refTrack = fEsdEvent->GetTrack(iRefTrack); | |
786 | isRelevant = (refTrack->Pt() >= fRefTrackSelectionPtThreshold); //## implement properly | |
787 | if (isRelevant){ | |
788 | if (refTrackCount < fkMaxRefTracksPerStack) | |
789 | refTrackIndices[refTrackCount++] = iRefTrack; | |
790 | else { | |
791 | LogError("number of HLT tracks exceeding limit of %d. Skipping some tracks.", fkMaxRefTracksPerStack); | |
792 | break; | |
793 | } | |
794 | } | |
795 | } | |
796 | ||
797 | // try to match GTU track with ref tracks | |
798 | for (UInt_t iGtuTrack = 0; iGtuTrack < numGtuTracks; ++iGtuTrack) { | |
799 | bestMatchRating = 0.; | |
800 | bestMatchRefIndex = -1; | |
801 | Double_t gpt = fTrackingData->GetTrackPt(iStack, iGtuTrack); | |
802 | UShort_t gpid = fTrackingData->GetTrackPID(iStack, iGtuTrack); | |
803 | UShort_t layerMask = fTrackingData->GetTrackLayerMask(iStack, iGtuTrack); | |
804 | ||
805 | Float_t trklLocalY[fkTRDLayers]; | |
806 | Int_t trklBinZ[fkTRDLayers]; | |
807 | for (UShort_t iLayer = 0; iLayer < fkTRDLayers; ++iLayer){ | |
808 | if ((layerMask >> iLayer) & 1){ | |
809 | trklLocalY[iLayer] = fTrackingData->GetTrackTrackletLocalY(iStack, iGtuTrack, iLayer); | |
810 | trklBinZ[iLayer] = fTrackingData->GetTrackTrackletBinZ(iStack, iGtuTrack, iLayer); | |
811 | } | |
812 | } | |
813 | ||
814 | for (UInt_t iRefTrack = 0; iRefTrack < refTrackCount; ++iRefTrack) { | |
815 | refTrack = fEsdEvent->GetTrack(refTrackIndices[iRefTrack]); | |
816 | ||
817 | // if (!CheckRefTrackCuts(refTrack)) | |
818 | // continue; | |
819 | ||
820 | numComparisonsDone++; | |
821 | ||
822 | // if ((!refTrack->GetOuterParam()) && (!refTrack->GetInnerParam())) // use tracks with TPC outer param only (all HLT tracks have it anyways) | |
823 | // continue; | |
824 | ||
825 | Int_t distRes = EstimateTrackDistance(refTrack, iStack, layerMask, trklLocalY, trklBinZ, | |
826 | magField, &distY, &distZ); | |
827 | ||
828 | if (fDebugLevel >= 3){ | |
829 | printf("CHECKMATCH = %i distY %.2f distZ %.2f pt: %.2f %.2f\n", | |
830 | distRes, distY, distZ, gpt, refTrack->GetSignedPt()); | |
831 | } | |
832 | ||
833 | if (distRes == 0){ | |
834 | matchRating = RateTrackMatch(distY, distZ, gpt, refTrack->GetSignedPt()); | |
835 | } else { | |
836 | matchRating = 0.; | |
837 | } | |
838 | ||
839 | if ((matchRating >= fMatchRatingThreshold) && (matchRating > bestMatchRating)) { | |
840 | bestMatchRefIndex = refTrackIndices[iRefTrack]; | |
841 | bestMatchRating = matchRating; | |
842 | } else if (matchRating > bestMatchRating) | |
843 | bestMatchRating = matchRating; | |
844 | ||
845 | // DbgLog("", Form("#match - comparing GTU track %d in S%02d-%d with ref track %d: [gpt: %+5.2f rpt: %+5.2f] dy: %.1f dz: %.1f rating: %.1f", | |
846 | // iGtuTrack, iStack/5, iStack%5, refTrackIndices[iRefTrack], | |
847 | // fTrackingData->GetTrackPt(iStack, iGtuTrack), refTrack->GetSignedPt(), | |
848 | // distY, distZ, matchRating)); | |
849 | ||
850 | } // loop over ref tracks in stack | |
851 | ||
852 | fHistMatchRating->Fill(bestMatchRating); | |
853 | fHistMatchRatingByPt->Fill(bestMatchRating, TMath::Abs(gpt)); | |
854 | fHistMatchRatingByPid->Fill(bestMatchRating, gpid); | |
855 | fHistTrackPt->Fill(gpt); | |
856 | fHistTrackPid->Fill(gpid); | |
857 | ||
858 | if (fExtendedHistos){ | |
859 | fHistMatchedRefTrackPid->Fill(refTrack->Pt(), refTrack->GetTPCsignal()); | |
860 | } | |
861 | ||
862 | if (bestMatchRefIndex >= 0){ | |
863 | // GTU track has matching reference track | |
864 | refTrack = fEsdEvent->GetTrack(bestMatchRefIndex); | |
865 | Double_t rpt = refTrack->GetSignedPt(); | |
866 | LogDebug("#match - match found: rating %.2f, gpt: %+5.2f, rpt: %+5.2f (%i) -> diff: %.2f%%", | |
867 | bestMatchRating,gpt, rpt, bestMatchRefIndex, (gpt - rpt)/rpt*100.); | |
868 | fTrackingData->SetTrackAddInfo(iStack, iGtuTrack, bestMatchRefIndex); | |
869 | if (fDebugLevel >= 3){ | |
870 | LogDebug("#match-info rating: %.2f gpt: %+.2f matchref: %d rpt: %+.2f gpid: %d S%02d-%d %d rpid: %.3f", | |
871 | bestMatchRating, gpt, bestMatchRefIndex, rpt, gpid, | |
872 | iStack/5, iStack%5, iGtuTrack, refTrack->GetTPCsignal()); | |
873 | } | |
874 | numMatchedTracks++; | |
875 | ||
876 | fHistTrackPtMatched->Fill(gpt); | |
877 | fHistTrackPtCorr->Fill(rpt, gpt); | |
878 | fHistTrackPidMatched->Fill(gpid); | |
879 | ||
880 | } else { | |
881 | if (fDebugLevel >= 3){ | |
882 | LogDebug("#match-info rating: %.2f gpt: %+.2f no ref matching gpid: %d S%02d-%d %d", | |
883 | bestMatchRating, gpt, gpid, | |
884 | iStack/5, iStack%5, iGtuTrack); | |
885 | } | |
886 | numUnmatchedTracks++; | |
887 | } | |
888 | ||
889 | } // loop over gtu tracks in stack | |
890 | ||
891 | } // loop over stacks | |
892 | ||
893 | LogInfo("#match - %d matched GTU tracks, %d unmatched", | |
894 | numMatchedTracks, numUnmatchedTracks); | |
895 | ||
896 | fHistTrackMatchingCombinations->Fill(0., (Double_t)(numRefTracks * fTrackingData->GetNumTracks())); // index 0: full combinatorics | |
897 | fHistTrackMatchingCombinations->Fill(1., numComparisonsDone); // index 1: track matching comparisons actually done | |
898 | ||
899 | return result; | |
900 | ||
901 | } | |
902 | ||
903 | int AliHLTTRDTriggerComponent::MatchTRDTracksHLT(){ | |
904 | ||
905 | if (!fHLTTracks){ | |
906 | LogError("HLT track data available."); | |
907 | return 0; | |
908 | } | |
909 | ||
910 | Double_t magField = GetBz(); | |
911 | unsigned int numGtuTracks; | |
912 | unsigned int numHLTTracks = 0; | |
913 | Bool_t isRelevant; | |
914 | Double_t distY, distZ; | |
915 | Double_t matchRating; | |
916 | Double_t bestMatchRating; | |
917 | Int_t bestMatchRefIndex; | |
918 | UInt_t numComparisonsDone = 0; | |
919 | UInt_t numUnmatchedTracks = 0; | |
920 | UInt_t numMatchedTracks = 0; | |
921 | ||
922 | Bool_t hltTrackPreSel[fkMaxRefTracks]; | |
923 | UInt_t iHltTrack; | |
924 | ||
925 | for (UShort_t iStack = 0; iStack < fkTRDStacks; ++iStack) { | |
926 | numGtuTracks = fTrackingData->GetNumTracks(iStack); | |
927 | ||
928 | // preselect relevant HLT tracks | |
929 | iHltTrack = 0; | |
930 | for (vector<AliHLTGlobalBarrelTrack>::iterator element = fHLTTracks->begin(); | |
931 | element != fHLTTracks->end(); element++) { | |
932 | numHLTTracks++; | |
933 | isRelevant = (element->Pt() >= fRefTrackSelectionPtThreshold); //## implement properly | |
934 | ||
935 | //## use cuts here | |
936 | hltTrackPreSel[iHltTrack++] = isRelevant; | |
1d94ad26 | 937 | if (iHltTrack >= fkMaxRefTracks) { |
c7b7f445 | 938 | LogError("maximum number of HLT tracks exceeded."); |
1d94ad26 | 939 | break; |
940 | } | |
c7b7f445 | 941 | } // loop over HLT tracks; |
942 | ||
943 | // search for matching HLT track for each GTU track | |
944 | for (UInt_t iGtuTrack = 0; iGtuTrack < numGtuTracks; ++iGtuTrack) { | |
945 | bestMatchRating = 0.; | |
946 | bestMatchRefIndex = -1; | |
947 | Double_t gpt = fTrackingData->GetTrackPt(iStack, iGtuTrack); | |
948 | UShort_t gpid = fTrackingData->GetTrackPID(iStack, iGtuTrack); | |
949 | UShort_t layerMask = fTrackingData->GetTrackLayerMask(iStack, iGtuTrack); | |
950 | ||
951 | Float_t trklLocalY[fkTRDLayers]; | |
952 | Int_t trklBinZ[fkTRDLayers]; | |
953 | for (UShort_t iLayer = 0; iLayer < fkTRDLayers; ++iLayer){ | |
954 | if ((layerMask >> iLayer) & 1){ | |
955 | trklLocalY[iLayer] = fTrackingData->GetTrackTrackletLocalY(iStack, iGtuTrack, iLayer); | |
956 | trklBinZ[iLayer] = fTrackingData->GetTrackTrackletBinZ(iStack, iGtuTrack, iLayer); | |
957 | } | |
958 | } | |
959 | ||
960 | iHltTrack = 0; | |
961 | for (vector<AliHLTGlobalBarrelTrack>::iterator element = fHLTTracks->begin(); | |
962 | element != fHLTTracks->end(); element++) { | |
963 | if (!hltTrackPreSel[iHltTrack]){ | |
964 | iHltTrack++; | |
965 | continue; | |
966 | } | |
967 | ||
968 | // compare GTU track and relevant HLT track | |
969 | numComparisonsDone++; | |
970 | ||
971 | AliHLTGlobalBarrelTrack hltPar(*element); | |
972 | Int_t distRes = EstimateTrackDistance(&hltPar, iStack, layerMask, trklLocalY, trklBinZ, | |
973 | magField, &distY, &distZ); | |
974 | ||
975 | if (fDebugLevel >= 3){ | |
976 | printf("CHECKMATCH = %i distY %.2f distZ %.2f pt: %.2f %.2f\n", | |
977 | distRes, distY, distZ, gpt, element->GetSignedPt()); | |
978 | } | |
979 | ||
980 | if (distRes == 0){ | |
981 | matchRating = RateTrackMatch(distY, distZ, gpt, element->GetSignedPt()); | |
982 | } else { | |
983 | matchRating = 0.; | |
984 | } | |
985 | ||
986 | if ((matchRating >= fMatchRatingThreshold) && (matchRating > bestMatchRating)) { | |
987 | bestMatchRefIndex = iHltTrack; | |
988 | bestMatchRating = matchRating; | |
989 | } else if (matchRating > bestMatchRating) | |
990 | bestMatchRating = matchRating; | |
991 | ||
992 | ||
993 | iHltTrack++; | |
994 | } // loop over HLT tracks; | |
995 | ||
996 | fHistMatchRating->Fill(bestMatchRating); | |
997 | fHistMatchRatingByPt->Fill(bestMatchRating, TMath::Abs(gpt)); | |
998 | fHistMatchRatingByPid->Fill(bestMatchRating, gpid); | |
999 | fHistTrackPt->Fill(gpt); | |
1000 | fHistTrackPid->Fill(gpid); | |
1001 | ||
1002 | if (bestMatchRefIndex >= 0){ | |
1003 | // GTU track has matching reference track | |
1004 | Double_t rpt = fHLTTracks->at(bestMatchRefIndex).GetSignedPt(); | |
1005 | LogDebug("#match - match found: rating %.2f, gpt: %+5.2f, rpt: %+5.2f (%i) -> diff: %.2f%%", | |
1006 | bestMatchRating,gpt, rpt, bestMatchRefIndex, (gpt - rpt)/rpt*100.); | |
1007 | fTrackingData->SetTrackAddInfo(iStack, iGtuTrack, bestMatchRefIndex); | |
1008 | ||
1009 | // if (fExtendedHistos){ | |
1010 | // fHistMatchedRefTrackPid->Fill(element->Pt(), element->GetTPCsignal()); | |
1011 | // } | |
1012 | // | |
1013 | // if (fDebugLevel >= 3){ | |
1014 | // LogDebug("#match-info rating: %.2f gpt: %+.2f matchref: %d rpt: %+.2f gpid: %d S%02d-%d %d rpid: %.3f", | |
1015 | // bestMatchRating, gpt, bestMatchRefIndex, rpt, gpid, | |
1016 | // iStack/5, iStack%5, iGtuTrack, refTrack->GetTPCsignal()); | |
1017 | // } | |
1018 | numMatchedTracks++; | |
1019 | ||
1020 | fHistTrackPtMatched->Fill(gpt); | |
1021 | fHistTrackPtCorr->Fill(rpt, gpt); | |
1022 | fHistTrackPidMatched->Fill(gpid); | |
1023 | ||
1024 | } else { | |
1025 | if (fDebugLevel >= 3){ | |
1026 | LogDebug("#match-info rating: %.2f gpt: %+.2f no ref matching gpid: %d S%02d-%d %d", | |
1027 | bestMatchRating, gpt, gpid, | |
1028 | iStack/5, iStack%5, iGtuTrack); | |
1029 | } | |
1030 | numUnmatchedTracks++; | |
1031 | } | |
1032 | ||
1033 | } // loop over gtu tracks | |
1034 | } // loop over stacks | |
1035 | ||
1036 | LogInfo("#match - %d matched GTU tracks, %d unmatched (%d comparisons)", | |
1037 | numMatchedTracks, numUnmatchedTracks, numComparisonsDone); | |
1038 | ||
1039 | fHistTrackMatchingCombinations->Fill(0., (Double_t)(numHLTTracks * fTrackingData->GetNumTracks())); // index 0: full combinatorics | |
1040 | fHistTrackMatchingCombinations->Fill(1., numComparisonsDone); // index 1: track matching comparisons actually done | |
1041 | ||
1042 | return kTRUE; | |
1043 | } | |
1044 | ||
1045 | int AliHLTTRDTriggerComponent::MatchTRDTracks(){ | |
1046 | ||
1047 | if (!fEsdEvent) { | |
1048 | LogError("ESD event data not available in MatchTRDTracks()."); | |
1049 | return 0; | |
1050 | } | |
1051 | ||
1052 | if (fHistoMode == 0){ | |
1053 | UInt_t numHists = fHistArray->GetEntries(); | |
1054 | for (UInt_t iHist = 0; iHist < numHists; ++iHist) | |
1055 | if (fHistArray->At(iHist)) | |
ca892663 | 1056 | ((TH1*) (fHistArray->At(iHist)))->Reset(); |
c7b7f445 | 1057 | } |
1058 | ||
1059 | int result = 1; | |
1060 | unsigned int numRefTracks = fEsdEvent->GetNumberOfTracks(); | |
1061 | unsigned int numGtuTracks; | |
1062 | Int_t refTrackIndices[fkMaxRefTracksPerStack]; | |
1063 | UInt_t refTrackCount = 0; | |
1064 | Bool_t isRelevant; | |
1065 | Double_t distY, distZ; | |
1066 | Double_t matchRating; | |
1067 | Double_t bestMatchRating; | |
1068 | Int_t bestMatchRefIndex; | |
1069 | AliESDtrack* refTrack; | |
1070 | UInt_t numComparisonsDone = 0; | |
1071 | UInt_t numUnmatchedTracks = 0; | |
1072 | UInt_t numMatchedTracks = 0; | |
1073 | ||
1074 | for (UShort_t iStack = 0; iStack < fkTRDStacks; ++iStack) { | |
1075 | numGtuTracks = fTrackingData->GetNumTracks(iStack); | |
1076 | refTrackCount = 0; | |
1077 | ||
1078 | // preselect ESD relevant ESD tracks | |
1079 | for (UInt_t iRefTrack = 0; iRefTrack < numRefTracks; ++iRefTrack) { | |
1080 | refTrack = fEsdEvent->GetTrack(iRefTrack); | |
1081 | isRelevant = (refTrack->Pt() >= fRefTrackSelectionPtThreshold); | |
1082 | if (isRelevant){ | |
1083 | if (refTrackCount < fkMaxRefTracksPerStack) | |
1084 | refTrackIndices[refTrackCount++] = iRefTrack; | |
1085 | else { | |
1086 | LogError("number of HLT tracks exceeding limit of %d. Skipping some tracks.", fkMaxRefTracksPerStack); | |
1087 | break; | |
1088 | } | |
1089 | } | |
1090 | } | |
1091 | ||
1092 | // try to match GTU track with ref tracks | |
1093 | for (UInt_t iGtuTrack = 0; iGtuTrack < numGtuTracks; ++iGtuTrack) { | |
1094 | bestMatchRating = 0.; | |
1095 | bestMatchRefIndex = -1; | |
1096 | Double_t gpt = fTrackingData->GetTrackPt(iStack, iGtuTrack); | |
1097 | UShort_t gpid = fTrackingData->GetTrackPID(iStack, iGtuTrack); | |
1098 | UShort_t layerMask = fTrackingData->GetTrackLayerMask(iStack, iGtuTrack); | |
1099 | ||
1100 | for (UInt_t iRefTrack = 0; iRefTrack < refTrackCount; ++iRefTrack) { | |
1101 | refTrack = fEsdEvent->GetTrack(refTrackIndices[iRefTrack]); | |
1102 | ||
1103 | if (!CheckRefTrackCuts(refTrack)) | |
1104 | continue; | |
1105 | ||
1106 | numComparisonsDone++; | |
1107 | ||
1108 | // if ((!refTrack->GetOuterParam()) && (!refTrack->GetInnerParam())) // use tracks with TPC outer param only (all HLT tracks have it anyways) | |
1109 | // continue; | |
1110 | ||
1111 | Float_t trklLocalY[fkTRDLayers]; | |
1112 | Int_t trklBinZ[fkTRDLayers]; | |
1113 | for (UShort_t iLayer = 0; iLayer < fkTRDLayers; ++iLayer){ | |
1114 | if ((layerMask >> iLayer) & 1){ | |
1115 | trklLocalY[iLayer] = fTrackingData->GetTrackTrackletLocalY(iStack, iGtuTrack, iLayer); | |
1116 | trklBinZ[iLayer] = fTrackingData->GetTrackTrackletBinZ(iStack, iGtuTrack, iLayer); | |
1117 | } | |
1118 | } | |
1119 | ||
1120 | Int_t distRes = EstimateTrackDistance(refTrack, iStack, layerMask, trklLocalY, trklBinZ, | |
1121 | fEsdEvent->GetMagneticField(), &distY, &distZ); | |
1122 | ||
1123 | if (distRes == 0){ | |
1124 | matchRating = RateTrackMatch(distY, distZ, gpt, refTrack->GetSignedPt()); | |
1125 | } else { | |
1126 | matchRating = 0.; | |
1127 | } | |
1128 | ||
1129 | if ((matchRating >= fMatchRatingThreshold) && (matchRating > bestMatchRating)) { | |
1130 | bestMatchRefIndex = refTrackIndices[iRefTrack]; | |
1131 | bestMatchRating = matchRating; | |
1132 | } else if (matchRating > bestMatchRating) | |
1133 | bestMatchRating = matchRating; | |
1134 | ||
1135 | // DbgLog("", Form("#match - comparing GTU track %d in S%02d-%d with ref track %d: [gpt: %+5.2f rpt: %+5.2f] dy: %.1f dz: %.1f rating: %.1f", | |
1136 | // iGtuTrack, iStack/5, iStack%5, refTrackIndices[iRefTrack], | |
1137 | // fTrackingData->GetTrackPt(iStack, iGtuTrack), refTrack->GetSignedPt(), | |
1138 | // distY, distZ, matchRating)); | |
1139 | ||
1140 | } // loop over ref tracks in stack | |
1141 | ||
1142 | fHistMatchRating->Fill(bestMatchRating); | |
1143 | fHistMatchRatingByPt->Fill(bestMatchRating, TMath::Abs(gpt)); | |
1144 | fHistMatchRatingByPid->Fill(bestMatchRating, gpid); | |
1145 | fHistTrackPt->Fill(gpt); | |
1146 | fHistTrackPid->Fill(gpid); | |
1147 | ||
1148 | if (bestMatchRefIndex >= 0){ | |
1149 | // GTU track has matching reference track | |
1150 | refTrack = fEsdEvent->GetTrack(bestMatchRefIndex); | |
1151 | Double_t rpt = refTrack->GetSignedPt(); | |
1152 | LogDebug("#match - match found: rating %.2f, gpt: %+5.2f, rpt: %+5.2f (%i) -> diff: %.2f%%", | |
1153 | bestMatchRating,gpt, rpt, bestMatchRefIndex, (gpt - rpt)/rpt*100.); | |
1154 | fTrackingData->SetTrackAddInfo(iStack, iGtuTrack, bestMatchRefIndex); | |
1155 | ||
1156 | if (fExtendedHistos){ | |
1157 | fHistMatchedRefTrackPid->Fill(refTrack->Pt(), refTrack->GetTPCsignal()); | |
1158 | } | |
1159 | ||
1160 | if (fDebugLevel >= 3){ | |
1161 | LogDebug("#match-info rating: %.2f gpt: %+.2f matchref: %d rpt: %+.2f gpid: %d S%02d-%d %d rpid: %.3f", | |
1162 | bestMatchRating, gpt, bestMatchRefIndex, rpt, gpid, | |
1163 | iStack/5, iStack%5, iGtuTrack, refTrack->GetTPCsignal()); | |
1164 | } | |
1165 | numMatchedTracks++; | |
1166 | ||
1167 | fHistTrackPtMatched->Fill(gpt); | |
1168 | fHistTrackPtCorr->Fill(rpt, gpt); | |
1169 | fHistTrackPidMatched->Fill(gpid); | |
1170 | ||
1171 | } else { | |
1172 | if (fDebugLevel >= 3){ | |
1173 | LogDebug("#match-info rating: %.2f gpt: %+.2f no ref matching gpid: %d S%02d-%d %d", | |
1174 | bestMatchRating, gpt, gpid, | |
1175 | iStack/5, iStack%5, iGtuTrack); | |
1176 | } | |
1177 | numUnmatchedTracks++; | |
1178 | } | |
1179 | ||
1180 | } // loop over gtu tracks in stack | |
1181 | ||
1182 | } // loop over stacks | |
1183 | ||
1184 | LogInfo("#match - %d matched GTU tracks, %d unmatched", | |
1185 | numMatchedTracks, numUnmatchedTracks); | |
1186 | ||
1187 | fHistTrackMatchingCombinations->Fill(0., (Double_t)(numRefTracks * fTrackingData->GetNumTracks())); // index 0: full combinatorics | |
1188 | fHistTrackMatchingCombinations->Fill(1., numComparisonsDone); // index 1: track matching comparisons actually done | |
1189 | ||
1190 | return result; | |
1191 | } | |
1192 | ||
1193 | void AliHLTTRDTriggerComponent::DumpTrackingData(){ | |
1194 | ||
1195 | if (fTrackingData->GetNumTracklets() + fTrackingData->GetNumTracks() == 0) | |
1196 | return; | |
1197 | ||
1198 | TString trklStr(""); | |
1199 | TString matchStr(""); | |
1200 | UShort_t layerMask; | |
1201 | ||
1202 | // for (UShort_t iSector = 0; iSector < 18; ++iSector){ | |
1203 | // if (fTrackingData->GetSectorTrgWord(iSector) != ((UInt_t)0x2345352 ^ ((UInt_t)iSector + 34))) | |
1204 | // LogError("invalid sector trigger word in sector %02d: trg word is 0x%08x, should be 0x%08x", | |
1205 | // iSector, fTrackingData->GetSectorTrgWord(iSector), ((UInt_t)0x2345352 ^ ((UInt_t)iSector + 34))); | |
1206 | // for (UShort_t iStack = 0; iStack < 5; ++iStack){ | |
1207 | // ULong64_t dwl = ((ULong64_t)0xaffe << 16) | (ULong64_t)iSector << 8 | iStack; | |
1208 | // ULong64_t dw = (((ULong64_t)0xbead << 16) | (ULong64_t)iSector << 8 | iStack) | ((ULong64_t)dwl << 32); | |
1209 | // if (fTrackingData->GetStackTrgWord(iSector, iStack) != dw) | |
1210 | // LogError("stack %02d-%d trg word is 0x%016llx, should be 0x%016llx", | |
1211 | // iSector, iStack, fTrackingData->GetStackTrgWord(iSector, iStack), dw); | |
1212 | // } | |
1213 | // } | |
1214 | ||
1215 | for (UShort_t iStack = 0; iStack < fkTRDStacks; ++iStack){ | |
1216 | for (Int_t iTrk = 0; iTrk < fTrackingData->GetNumTracks(iStack); ++iTrk){ | |
1217 | ||
1218 | layerMask = fTrackingData->GetTrackLayerMask(iStack, iTrk); | |
1219 | ||
1220 | trklStr = Form("trkl: "); | |
1221 | for (Short_t iLayer = 5; iLayer >= 0; --iLayer){ | |
1222 | if ((layerMask >> iLayer) & 1) | |
1223 | trklStr += Form("0x%08x (%+8.3f) ", | |
1224 | fTrackingData->GetTrackTrackletWord(iStack, iTrk, iLayer), | |
1225 | fTrackingData->GetTrackTrackletLocalY(iStack, iTrk, iLayer)); | |
1226 | else | |
1227 | trklStr += "--------------------- "; | |
1228 | } // loop over layers | |
1229 | trklStr.Remove(trklStr.Length() - 2, 2); | |
1230 | ||
1231 | if (fTrackingData->GetTrackAddInfo(iStack, iTrk) >= 0){ | |
1232 | Double_t rpt = (fESDtracksPresent) ? fEsdEvent->GetTrack(fTrackingData->GetTrackAddInfo(iStack, iTrk))->GetSignedPt() : | |
1233 | fHLTTracks->at(fTrackingData->GetTrackAddInfo(iStack, iTrk)).GetSignedPt(); | |
1234 | matchStr = Form("mpt: %+7.2f", rpt); | |
1235 | } else | |
1236 | matchStr = "unmatched"; | |
1237 | ||
1238 | if (fDebugLevel >= 3){ | |
1239 | ||
1240 | printf("###DOTDA EV%04llu GTU TRACK - S%02d-%d pt: %+7.2f pid: %3d lm: 0x%02x %s %s\n", | |
1241 | fEventId, | |
1242 | iStack/5, iStack%5, fTrackingData->GetTrackPt(iStack, iTrk), fTrackingData->GetTrackPID(iStack, iTrk), | |
1243 | layerMask, trklStr.Data(), | |
1244 | matchStr.Data()); | |
1245 | ||
1246 | printf("###DOTDB EV%04llu GTU TRACK - S%02d-%d pt: %+7.2f pid: %3d lm: 0x%02x %s\n", | |
1247 | fEventId, | |
1248 | iStack/5, iStack%5, fTrackingData->GetTrackPt(iStack, iTrk), fTrackingData->GetTrackPID(iStack, iTrk), | |
1249 | layerMask, trklStr.Data()); | |
1250 | } | |
1251 | ||
1252 | ||
1253 | // paranoia checks | |
1254 | for (Short_t iLayer = 5; iLayer >= 0; --iLayer){ | |
1255 | if (((layerMask >> iLayer) & 1) && (fTrackingData->GetTrackTrackletWord(iStack, iTrk, iLayer) == 0x10001000)) | |
1256 | LogError("invalid layer mask / tracklet value combination A"); | |
1257 | ||
1258 | if ((((layerMask >> iLayer) & 1) == 0) && (fTrackingData->GetTrackTrackletWord(iStack, iTrk, iLayer) != 0x10001000)) | |
1259 | LogError("invalid layer mask / tracklet value combination B"); | |
1260 | } | |
1261 | ||
1262 | } // loop over tracks in stack | |
1263 | } // loop over stacks | |
1264 | ||
1265 | } | |
1266 | ||
1267 | void AliHLTTRDTriggerComponent::AssignTrackInfo(TString* infoStr, const UInt_t stack, const UInt_t trackIndex, const char* flagStr) { | |
1268 | ||
1269 | TString flags(flagStr); | |
1270 | Bool_t appendRefData = kFALSE; | |
1271 | ||
1272 | if (fTrackingData->GetTrackAddInfo(stack, trackIndex) >= 0){ | |
1273 | appendRefData = kTRUE; | |
1274 | if (flags.First('M') < 0) | |
1275 | flags += "M"; | |
1276 | } | |
1277 | ||
1278 | *infoStr = Form("EXCH-TRK-INFO DS<%s> CH<%s> EV%05llu SEC%02d-%d-%d gpt: %+6.1f gpid: %3d flags: %s", | |
1279 | "TRDHLTTRG", | |
1280 | fChunkId->Data(), fEventId, | |
1281 | stack/5, stack%5, trackIndex, | |
1282 | fTrackingData->GetTrackPt(stack, trackIndex), | |
1283 | fTrackingData->GetTrackPID(stack, trackIndex), flags.Data()); | |
1284 | ||
1285 | if (appendRefData){ | |
1286 | Double_t rpt = (fESDtracksPresent) ? fEsdEvent->GetTrack(fTrackingData->GetTrackAddInfo(stack, trackIndex))->GetSignedPt() : | |
1287 | fHLTTracks->at(fTrackingData->GetTrackAddInfo(stack, trackIndex)).GetSignedPt(); | |
1288 | *infoStr += Form(" rpt: %+6.1f", rpt); | |
1289 | } | |
1290 | ||
1291 | } | |
1292 | ||
1293 | #ifdef __TRDHLTDEBUG | |
1294 | void AliHLTTRDTriggerComponent::RenderEvent(const Bool_t showGtuTracks, const Bool_t showTracklets, const Bool_t showRefTracks) { | |
1295 | ||
1296 | if ((!fEventDisplay) || (!fEsdEvent)) | |
1297 | return; | |
1298 | ||
1299 | LogDebug("rendering event"); | |
1300 | ||
1301 | const Float_t refTrackPtDisplayThreshold = 0.7; | |
1302 | const Float_t trackPtEmphasizeThreshold = 1.8; | |
1303 | ||
1304 | fEventDisplay->Reset(); | |
1305 | fEventDisplay->SetMagField(fEsdEvent->GetMagneticField()); | |
1306 | for (UInt_t iDet = 0; iDet < fkTRDStacks*fkTRDLayers; ++iDet) | |
1307 | fEventDisplay->SetChamberState(iDet/fkTRDLayers, iDet%6, ((fSectorsWithData >> (iDet/30)) & 1) ); | |
1308 | ||
1309 | if (showTracklets){ | |
1310 | for (UShort_t iDet = 0; iDet < fkTRDStacks*fkTRDLayers; ++iDet){ | |
1311 | UInt_t numTrkl = fTrackingData->GetNumTracklets(iDet); | |
1312 | for (UInt_t iTrkl = 0; iTrkl < numTrkl; ++iTrkl){ | |
1313 | fEventDisplay->AddTracklet(fTrackingData->GetTracklet(iDet, iTrkl)); | |
1314 | } | |
1315 | } | |
1316 | } | |
1317 | ||
1318 | if (showGtuTracks){ | |
1319 | for (UShort_t iStack = 0; iStack < fkTRDStacks; ++iStack){ | |
1320 | for (Int_t iTrk = 0; iTrk < fTrackingData->GetNumTracks(iStack); ++iTrk){ | |
1321 | AliESDTrdTrack* track = fTrackingData->GetTrack(iStack, iTrk); | |
1322 | Int_t trkIndex = fEventDisplay->AddTrack(track, (TMath::Abs(track->Pt()) >= trackPtEmphasizeThreshold) ? (kMagenta + 2) : kMagenta, 1, 1); | |
1323 | if (fTrackingData->GetTrackAddInfo(iStack, iTrk) >= 0) | |
1324 | fEventDisplay->SetTrackTrackletStyle(trkIndex, kRed, 12); | |
1325 | else | |
1326 | fEventDisplay->SetTrackTrackletStyle(trkIndex, kViolet - 5, 12); | |
1327 | } | |
1328 | } | |
1329 | } | |
1330 | ||
1331 | unsigned int numRefTracks = fEsdEvent->GetNumberOfTracks(); | |
1332 | unsigned int numRefTracksRendered = 0; | |
1333 | if (showRefTracks){ | |
1334 | // check for some marks for rendering | |
1335 | UShort_t marks[40000]; | |
1336 | memset(marks, 0, sizeof(UShort_t)*40000); | |
1337 | for (UShort_t iStack = 0; iStack < fkTRDStacks; ++iStack){ | |
1338 | for (Int_t iTrk = 0; iTrk < fTrackingData->GetNumTracks(iStack); ++iTrk){ | |
1339 | if (fTrackingData->GetTrackAddInfo(iStack, iTrk) >= 0){ | |
1340 | marks[fTrackingData->GetTrackAddInfo(iStack, iTrk)] |= 1; | |
1341 | if ( | |
1342 | (fTrackingData->GetTrackPID(iStack, iTrk) >= fElectronTriggerPIDThresholdHQU) && | |
1343 | (TMath::Abs(fTrackingData->GetTrackPt(iStack, iTrk)) >= fElectronTriggerPtThresholdHQU)) | |
1344 | marks[fTrackingData->GetTrackAddInfo(iStack, iTrk)] |= 2; | |
1345 | } | |
1346 | } | |
1347 | } | |
1348 | ||
1349 | // add to rendering | |
1350 | for (unsigned int iTrack = 0; iTrack < numRefTracks; ++iTrack){ | |
1351 | AliESDtrack* track = fEsdEvent->GetTrack(iTrack); | |
1352 | if ((track->Pt() >= refTrackPtDisplayThreshold) || (marks[iTrack] != 0)){ | |
1353 | Color_t color = (track->Pt() >= trackPtEmphasizeThreshold) ? (kGray + 2) : kGray; | |
1354 | UShort_t width = 1; | |
1355 | UShort_t style = 1; | |
1356 | if (marks[iTrack] & 0x2){ | |
1357 | color = kRed + 1; | |
1358 | width = 6; | |
1359 | } | |
1360 | if (!CheckRefTrackCuts(track)) | |
1361 | style = 2; | |
1362 | ||
1363 | fEventDisplay->AddTrack(track, color, width, style); | |
1364 | numRefTracksRendered++; | |
1365 | } | |
1366 | } | |
1367 | } | |
1368 | ||
1369 | if (!fEventDisplay->IsEmpty()){ | |
1370 | fEventDisplay->SetTitle(""); | |
1371 | fEventDisplay->SetSetupText("HLT", Form("%.1f#scale[0.5]{ }/#scale[0.5]{ }%.1f", | |
1372 | refTrackPtDisplayThreshold, trackPtEmphasizeThreshold)); | |
1373 | fEventDisplay->SetBottomText(Form("%05i-%s-%05llu", | |
1374 | fRunNumber, fChunkId->Data(), fEventId)); | |
1375 | fEventDisplay->SetBottomTextRight(Form("%d (%d) HLT tracks, %d tracklets, %d GTU tracks", | |
1376 | numRefTracks, | |
1377 | numRefTracksRendered, | |
1378 | fTrackingData->GetNumTracklets(), | |
1379 | fTrackingData->GetNumTracks())); | |
1380 | fEventDisplay->SetLook(AliTRDtrackingEventDisplay::dmMediumLight); | |
1381 | fEventDisplay->SetDisplayMode(AliTRDtrackingEventDisplay::dmFullXY); | |
1382 | fEventDisplay->SaveToFile(Form("display/event-%s-%05llu.eps", | |
1383 | fChunkId->Data(), fEventId)); | |
1384 | } | |
1385 | ||
1386 | } | |
1387 | #endif | |
1388 | ||
1389 | Bool_t AliHLTTRDTriggerComponent::TRDElectronTrigger(const char *ident, const Double_t minPt, const UShort_t minPID){ | |
1390 | ||
1391 | LogDebug("Electron trigger processing (%s: pt>=%.1f, pid>=%d)...", ident, minPt, minPID); | |
1392 | ||
1393 | UInt_t numTracks; | |
1394 | Bool_t highPtElectronSeenGTU = kFALSE; | |
1395 | Bool_t highPtElectronSeen = kFALSE; | |
1396 | UInt_t truncMeanPID = 0; | |
1397 | TString trackExchangeInfo(""); | |
1398 | TString flags(""); | |
1399 | ||
1400 | UInt_t trdSectorTrgContribs = 0; | |
1401 | for (UShort_t iSector = 0; iSector < fkTRDSectors; ++iSector) | |
1402 | trdSectorTrgContribs |= fTrackingData->GetSectorTrgContribs(iSector); | |
1403 | if ((trdSectorTrgContribs >> 5) & 0x3) | |
1404 | fIsTRDElectronEvent = kTRUE; | |
1405 | else | |
1406 | fIsTRDElectronEvent = kFALSE; | |
1407 | ||
1408 | if (fIsMinBiasEvent) | |
1409 | fHistElectronTriggerBaseMinBias->Fill(0.); | |
1410 | ||
1411 | if (fIsTRDElectronEvent) | |
1412 | fHistElectronTriggerBaseTrdL1->Fill(0.); | |
1413 | ||
1414 | if ((fElectronTriggerOnL1TrgOnly) && (!fIsTRDElectronEvent)){ | |
1415 | // evaluate trigger for events with TRD L1 electron trigger fired | |
1416 | DbgLog("skipping %s electron trigger evaluation for event, because no TRD trigger flag set (ctbs: 0x%02x)", | |
1417 | ident, trdSectorTrgContribs); | |
1418 | return 0; | |
1419 | } | |
1420 | ||
1421 | for (UShort_t iStack = 0; iStack < fkTRDStacks; ++iStack){ | |
1422 | numTracks = fTrackingData->GetNumTracks(iStack); | |
1423 | for (UInt_t iTrk = 0; iTrk < numTracks; ++iTrk){ | |
1424 | Double_t gpt = fTrackingData->GetTrackPt(iStack, iTrk); | |
1425 | Bool_t trdElectronCandidate = kFALSE; | |
1426 | Bool_t hltElectronCandidate = kFALSE; | |
1427 | ||
1428 | // re-evaluate GTU only decision for comparison | |
1429 | if ( | |
1430 | (TMath::Abs(gpt) >= minPt) && | |
1431 | (fTrackingData->GetTrackPID(iStack, iTrk) >= minPID) | |
1432 | ) { | |
1433 | trdElectronCandidate = kTRUE; | |
1434 | highPtElectronSeenGTU = kTRUE; | |
1435 | fHistElectronCandidatePt->Fill(gpt); | |
1436 | fHistElectronCandidatePid->Fill(fTrackingData->GetTrackPID(iStack, iTrk)); | |
1437 | ||
1438 | if (fExtendedHistos){ | |
1439 | ||
1440 | // idea to be checked | |
1441 | truncMeanPID = 0; | |
1442 | Short_t minLyr = -1; | |
1443 | UShort_t minValue = 255; | |
1444 | Short_t maxLyr = -1; | |
1445 | UShort_t maxValue = 0; | |
1446 | UShort_t lyrCount = 0; | |
1447 | UInt_t layerMask = fTrackingData->GetTrackLayerMask(iStack, iTrk); | |
1448 | ||
1449 | // scan for min & max values | |
1450 | for (UShort_t iLayer = 0; iLayer < fkTRDLayers; ++iLayer){ | |
1451 | if ((layerMask >> iLayer) & 1){ | |
1452 | if (fTrackingData->GetTrackTrackletPID(iStack, iTrk, iLayer) < minValue){ | |
1453 | minValue = fTrackingData->GetTrackTrackletPID(iStack, iTrk, iLayer); | |
1454 | minLyr = iLayer; | |
1455 | } | |
1456 | if (fTrackingData->GetTrackTrackletPID(iStack, iTrk, iLayer) > maxValue){ | |
1457 | maxValue = fTrackingData->GetTrackTrackletPID(iStack, iTrk, iLayer); | |
1458 | maxLyr = iLayer; | |
1459 | } | |
1460 | } | |
1461 | } // loop over layers | |
1462 | ||
1463 | // calculate trunc mean | |
1464 | for (UShort_t iLayer = 0; iLayer < fkTRDLayers; ++iLayer){ | |
1465 | if (((layerMask >> iLayer) & 1) && (iLayer != minLyr) && (iLayer != maxLyr)){ | |
1466 | truncMeanPID += fTrackingData->GetTrackTrackletPID(iStack, iTrk, iLayer); | |
1467 | lyrCount++; | |
1468 | } | |
1469 | } // loop over layers | |
1470 | truncMeanPID = TMath::Nint((Double_t)(truncMeanPID)/(Double_t)(lyrCount)); | |
1471 | ||
1472 | fHistPIDvsTruncPID->Fill(fTrackingData->GetTrackPID(iStack, iTrk), truncMeanPID); | |
1473 | if (fTrackingData->GetTrackAddInfo(iStack, iTrk) < 0) | |
1474 | fHistElectronFalsePIDvsTruncPID->Fill(fTrackingData->GetTrackPID(iStack, iTrk), truncMeanPID); | |
1475 | ||
1476 | } | |
1477 | ||
1478 | LogInspect("#hlt-trd-trg - GTU flagged %s %s high-pt electron seen in S%02d-%d: pt: %+6.1f pid: %d [id: S%02d-%d-%d] [trunc-pid: %d]", | |
1479 | (fTrackingData->GetTrackAddInfo(iStack, iTrk) < 0) ? "unmatched" : "matched", ident, | |
1480 | iStack/5, iStack%5, | |
1481 | gpt, fTrackingData->GetTrackPID(iStack, iTrk), | |
1482 | iStack/5, iStack%5, iTrk, | |
1483 | truncMeanPID); | |
1484 | } | |
1485 | ||
1486 | // evaluate HLT-level trigger decision using match information | |
1487 | if ( | |
1488 | (fTrackingData->GetTrackAddInfo(iStack, iTrk) >= 0) && | |
1489 | (TMath::Abs(gpt) >= minPt) && | |
1490 | (fTrackingData->GetTrackPID(iStack, iTrk) >= minPID) | |
1491 | ) { | |
1492 | hltElectronCandidate = kTRUE; | |
1493 | highPtElectronSeen = kTRUE; | |
1494 | ||
1495 | fHistElectronCandidateMatchedPt->Fill(gpt); | |
1496 | fHistElectronCandidateMatchedPid->Fill(fTrackingData->GetTrackPID(iStack, iTrk)); | |
1497 | ||
1498 | if (fExtendedHistos){ | |
1499 | fHistElectronConfirmedPIDvsTruncPID->Fill(fTrackingData->GetTrackPID(iStack, iTrk), truncMeanPID); | |
1500 | } | |
1501 | ||
1502 | Double_t rpt = (fESDtracksPresent) ? fEsdEvent->GetTrack(fTrackingData->GetTrackAddInfo(iStack, iTrk))->GetSignedPt() : | |
1503 | fHLTTracks->at(fTrackingData->GetTrackAddInfo(iStack, iTrk)).GetSignedPt(); | |
1504 | LogInspect("#hlt-trd-trg - HLT matched %s high-pt electron seen in S%02d-%d: gpt: %+6.1f gpid: %d rpt: %+6.1f [id: S%02d-%d-%d]", | |
1505 | ident, iStack/5, iStack%5, fTrackingData->GetTrackPt(iStack, iTrk), fTrackingData->GetTrackPID(iStack, iTrk), | |
1506 | rpt, iStack/5, iStack%5, iTrk); | |
1507 | } | |
1508 | ||
1509 | // log output for subsequent offline analysis | |
1510 | if ((fDebugLevel >= 3) && (trdElectronCandidate || hltElectronCandidate)){ | |
1511 | trackExchangeInfo = ""; | |
1512 | flags = ""; | |
1513 | if (trdElectronCandidate) | |
1514 | flags += "G"; | |
1515 | if (hltElectronCandidate) | |
1516 | flags += "H"; | |
1517 | AssignTrackInfo(&trackExchangeInfo, iStack, iTrk, flags.Data()); | |
1518 | LogDebug("%s\n", trackExchangeInfo.Data()); | |
1519 | } | |
1520 | ||
1521 | } // loop over tracks in stack | |
1522 | } // loop over stacks | |
1523 | ||
1524 | if (highPtElectronSeenGTU || highPtElectronSeen){ | |
1525 | LogInspect("#hlt-trd-trg - event triggered by %s electron trigger (TRD L1: %d, TRD HLT: %d)", | |
1526 | ident, highPtElectronSeenGTU, highPtElectronSeen); | |
1527 | } | |
1528 | ||
1529 | if (highPtElectronSeenGTU){ | |
1530 | if (fIsMinBiasEvent) | |
1531 | fHistElectronTriggerBaseMinBias->Fill(1.); | |
1532 | } | |
1533 | ||
1534 | if (highPtElectronSeen){ | |
1535 | if (fIsMinBiasEvent) | |
1536 | fHistElectronTriggerBaseMinBias->Fill(2.); | |
1537 | if (fIsTRDElectronEvent) | |
1538 | fHistElectronTriggerBaseTrdL1->Fill(1.); | |
1539 | } | |
1540 | ||
1541 | return (highPtElectronSeen) ? kTRUE : kFALSE; | |
1542 | } | |
1543 | ||
1544 | ||
1545 | int AliHLTTRDTriggerComponent::DoTrigger() | |
1546 | { | |
1547 | ||
1548 | fEsdEvent = NULL; | |
1549 | fRunNumber = -1; | |
1550 | int iResult = 0; | |
1551 | UShort_t firedTriggers = 0; | |
1552 | ||
1553 | const AliHLTComponentEventData* hltEventData = GetEventData(); | |
1554 | fEventId = hltEventData->fEventID; | |
1555 | LogDebug("### START DoTrigger [event id: %llu, %d blocks, size: %d]", | |
1556 | fEventId, hltEventData->fBlockCnt, hltEventData->fStructSize); | |
1557 | ||
1558 | // LogDebug("### START DoTrigger [event id: %llu, %d blocks, size: %d]", | |
1559 | // fEventId, hltEventData->fBlockCnt, hltEventData->fStructSize); | |
1560 | ||
1561 | if (!IsDataEvent()) { // process data events only | |
1562 | IgnoreEvent(); | |
1563 | LogDebug("### END DoTrigger [event id: %llu, %d blocks, size: %d] (skipped: no data event)", | |
1564 | fEventId, hltEventData->fBlockCnt, hltEventData->fStructSize); | |
1565 | return iResult; | |
1566 | } | |
1567 | ||
1568 | fTrackingData->SetLogPrefix(Form("TRDHLTGM XXXXXX-%05llu: [TRG] {TrkDat} ", fEventId)); | |
1569 | ||
1570 | do { | |
1571 | ||
1572 | // access to TRD specific data from AliHLTTRDPreprocessorComponent | |
1573 | if (!PrepareTRDData()){ | |
1574 | LogError("access to TRD data failed. Skipping event..."); | |
1575 | break; | |
1576 | } | |
1577 | ||
1578 | if (fTrackingData->GetNumTracks() + fTrackingData->GetNumTracklets() == 0) { | |
1579 | LogDebug("no trigger-relevant TRD information, skipping further event processing"); | |
1580 | break; | |
1581 | } | |
1582 | ||
1583 | // access to ESD data | |
1584 | if (!PrepareESDData()){ | |
1585 | LogInfo("access to ESD event data failed."); | |
1586 | } | |
1587 | ||
1588 | // access to alternative HLT data | |
1589 | if (!fESDtracksPresent){ | |
1590 | if (!PrepareHLTData()){ | |
1591 | LogError("access to HLT event data failed."); | |
1592 | } | |
1593 | } | |
1594 | ||
1595 | // match TRD and HLT tracks | |
1596 | if (fESDtracksPresent){ | |
1597 | if (!MatchTRDTracksESD()){ | |
1598 | LogError("matching TRD tracks to ESD tracks failed. Skipping event..."); | |
1599 | break; | |
1600 | } | |
1601 | } else if (fHLTtracksPresent){ | |
1602 | if (!MatchTRDTracksHLT()){ | |
1603 | LogError("matching TRD tracks to HLT tracks failed. Skipping event..."); | |
1604 | break; | |
1605 | } | |
1606 | } else { | |
1607 | LogError("No HLT track information available. Skipping event..."); | |
1608 | break; | |
1609 | } | |
1610 | ||
1611 | // if (!MatchTRDTracks()){ | |
1612 | // LogError("matching TRD tracks to TPC tracks failed. Skipping event..."); | |
1613 | // break; | |
1614 | // } | |
1615 | ||
1616 | if (fDebugLevel >= 1) | |
1617 | DumpTrackingData(); | |
1618 | ||
1619 | // evaluate electron trigger conditions | |
1620 | if (TRDElectronTrigger("HSE", fElectronTriggerPtThresholdHSE, fElectronTriggerPIDThresholdHSE)) | |
1621 | firedTriggers |= fkElectronTriggerHSE; | |
1622 | ||
1623 | if (TRDElectronTrigger("HQU", fElectronTriggerPtThresholdHQU, fElectronTriggerPIDThresholdHQU)) | |
1624 | firedTriggers |= fkElectronTriggerHQU; | |
1625 | ||
1626 | break; | |
1627 | ||
1628 | } while (1); | |
1629 | ||
1630 | ||
1631 | // trigger decision | |
1632 | TString description(""); | |
1633 | if (firedTriggers & fkElectronTriggerHSE){ | |
1634 | if (description.Length() > 0) | |
1635 | description += " "; | |
1636 | description += fgkTriggerDecisionElectronHSE; | |
1637 | } | |
1638 | ||
1639 | if (firedTriggers & fkElectronTriggerHQU){ | |
1640 | if (description.Length() > 0) | |
1641 | description += " "; | |
1642 | description += fgkTriggerDecisionElectronHQU; | |
1643 | } | |
1644 | ||
1645 | SetDescription(description.Data()); | |
1646 | AliHLTTriggerDecision decision((firedTriggers) ? kTRUE : kFALSE, | |
1647 | GetTriggerName(), | |
1648 | GetReadoutList(), | |
1649 | GetDescription() | |
1650 | ); | |
1651 | TriggerEvent(&decision, kAliHLTDataTypeTObject | kAliHLTDataOriginOut); | |
1652 | ||
1653 | if (firedTriggers){ | |
1654 | LogInspect("TRD HLT trigger fired for event: description: >%s<, flags: 0x%04x", | |
1655 | description.Data(), firedTriggers); | |
1656 | #ifdef __TRDHLTDEBUG | |
1657 | if (fEventRendering) | |
1658 | RenderEvent(); | |
1659 | #endif | |
1660 | } else { | |
1661 | LogDebug("TRD HLT trigger did not fire for event"); | |
1662 | } | |
1663 | ||
1664 | if (fPushHistos){ | |
1665 | PushBack(fHistArray, (kAliHLTDataTypeTObjArray | kAliHLTDataOriginTRD), 0x3fffff); | |
1666 | } | |
1667 | ||
1668 | LogDebug("### END DoTrigger [event id: %llu, %d blocks, size: %d]", | |
1669 | fEventId, hltEventData->fBlockCnt, hltEventData->fStructSize); | |
1670 | ||
1671 | return iResult; | |
1672 | } | |
1673 | ||
1674 | Bool_t AliHLTTRDTriggerComponent::TrackPlaneIntersect(AliExternalTrackParam *trk, Double_t pnt[3], Double_t norm[3], Double_t mag){ | |
1675 | ||
1676 | UInt_t its = 0; | |
1677 | Double_t r = 290.; | |
1678 | Double_t dist = 99999, dist_prev = 99999; | |
1679 | Double_t x[3] = {0., 0., 0.}; | |
1680 | Bool_t ret = kTRUE; | |
1681 | ||
1682 | dist = (x[0] - pnt[0]) * norm[0] + (x[1] - pnt[1]) *norm[1] + (x[2] - pnt[2]) * norm[2]; | |
1683 | ||
1684 | while(TMath::Abs(dist) > 0.1) { | |
1685 | ||
1686 | trk->GetXYZAt(r, mag, x); | |
1687 | ||
1688 | if ((x[0] * x[0] + x[1] * x[1]) < 100.) { // extrapolation to radius failed | |
1689 | ret = kFALSE; | |
1690 | break; | |
1691 | } | |
1692 | ||
1693 | //distance between current track position and plane | |
1694 | dist_prev = TMath::Abs(dist); | |
1695 | dist = (x[0] - pnt[0]) * norm[0] + (x[1] - pnt[1]) * norm[1]; | |
1696 | r -= dist; | |
1697 | its++; | |
1698 | if(TMath::Abs(dist) >= dist_prev || | |
1699 | (r > 380.) || (r < 100.)){ | |
1700 | ret = kFALSE; | |
1701 | break; | |
1702 | } | |
1703 | } | |
1704 | ||
1705 | for (Int_t i=0; i<3; i++){ | |
1706 | if(ret) | |
1707 | pnt[i] = x[i]; | |
1708 | else | |
1709 | pnt[i] = 0.; | |
1710 | } | |
1711 | ||
1712 | return kTRUE; | |
1713 | } | |
1714 | ||
1715 | ||
1716 | Double_t AliHLTTRDTriggerComponent::RateTrackMatch(Double_t distY, Double_t distZ, Double_t rpt, Double_t gpt){ | |
1717 | ||
1718 | // maximum limits for spatial distance | |
1719 | if ((distY > 5.) || (distZ > 20.)) | |
1720 | return 0.; | |
1721 | ||
1722 | // same pt sign required | |
1723 | if ((rpt * gpt) < 0.) | |
1724 | return 0.; | |
1725 | ||
1726 | Double_t rating_distY = -0.1 * distY + 1.; | |
1727 | Double_t rating_distZ = -0.025 * distZ + 1.; | |
1728 | Double_t rating_ptDiff = 1. - TMath::Abs((TMath::Abs(rpt) > 0.000001) ? ((gpt-rpt)/rpt) : 0.); | |
1729 | ||
1730 | if (rating_ptDiff < 0.) | |
1731 | rating_ptDiff = 0.2; | |
1732 | ||
1733 | Double_t total = rating_distY * rating_distZ * rating_ptDiff; | |
1734 | ||
1735 | // DbgLog("", Form("#matching: rating: dy: %.3f dz: %.3f dpt: %.3f -> total: %.3f", | |
1736 | // rating_distY, rating_distZ, rating_ptDiff, total)); | |
1737 | ||
1738 | if (total > 1.) | |
1739 | LogError("track match rating exceeds limit of 1.0: %.3f", total); | |
1740 | ||
1741 | return total; | |
1742 | } | |
1743 | ||
1744 | Int_t AliHLTTRDTriggerComponent::EstimateTrackDistance(AliESDtrack *esd_track, | |
1745 | const UShort_t stack, | |
1746 | const UShort_t layerMask, | |
1747 | const Float_t trklLocalY[6], const Int_t trklBinZ[6], | |
1748 | Double_t mag, Double_t *ydist, Double_t *zdist){ | |
1749 | if (!esd_track) | |
1750 | return -3; | |
1751 | ||
1752 | AliExternalTrackParam* refParam = NULL; | |
1753 | if (esd_track->GetOuterParam()) | |
1754 | refParam = new AliExternalTrackParam(*(esd_track->GetOuterParam())); | |
1755 | if (!refParam) | |
1756 | refParam = new AliExternalTrackParam(*(esd_track)); | |
1757 | ||
1758 | Int_t res = EstimateTrackDistance(refParam, stack, layerMask, trklLocalY, trklBinZ, mag, ydist, zdist); | |
1759 | ||
1760 | if (refParam) | |
1761 | delete refParam; | |
1762 | ||
1763 | return res; | |
1764 | ||
1765 | } | |
1766 | ||
1767 | Int_t AliHLTTRDTriggerComponent::EstimateTrackDistance(AliExternalTrackParam *refParam, | |
1768 | const UShort_t stack, | |
1769 | const UShort_t layerMask, | |
1770 | const Float_t trklLocalY[6], const Int_t trklBinZ[6], | |
1771 | Double_t mag, Double_t *ydist, Double_t *zdist){ | |
1772 | ||
1773 | Float_t diff_y = 0; | |
1774 | Float_t diff_z = 0; | |
1775 | Int_t nLayers = 0; | |
1776 | Double_t xtrkl[3]; | |
1777 | Double_t ptrkl[3]; | |
1778 | Double_t ptrkl2[3]; | |
1779 | UInt_t trklDet; | |
1780 | UShort_t trklLayer; | |
1d94ad26 | 1781 | // UInt_t stack_gtu; |
c7b7f445 | 1782 | UShort_t stackInSector; |
1783 | ||
1784 | AliTRDpadPlane* padPlane; | |
1785 | ||
1786 | for (UShort_t iLayer = 0; iLayer < 6; iLayer++){ | |
1787 | if ((layerMask >> iLayer) & 1){ | |
1788 | trklDet = stack*6 + iLayer; | |
1789 | trklLayer = iLayer; | |
1d94ad26 | 1790 | // stack_gtu = stack; |
c7b7f445 | 1791 | stackInSector = stack % 5; |
1792 | ||
1793 | // local coordinates of the outer end point of the tracklet | |
1794 | xtrkl[0] = AliTRDgeometry::AnodePos(); | |
1795 | xtrkl[1] = trklLocalY[iLayer]; | |
1796 | ||
1797 | padPlane = fTRDGeometry->GetPadPlane(trklLayer, stackInSector); | |
1798 | if(stackInSector == 2){ // corrected version by Felix Muecke | |
1799 | xtrkl[2] = padPlane->GetRowPos(trklBinZ[iLayer]) - (padPlane->GetRowSize(trklBinZ[iLayer]))/2. - padPlane->GetRowPos(6); | |
1800 | } else { | |
1801 | xtrkl[2] = padPlane->GetRowPos(trklBinZ[iLayer]) - (padPlane->GetRowSize(trklBinZ[iLayer]))/2. - padPlane->GetRowPos(8); | |
1802 | } | |
1803 | ||
1804 | // transform to global coordinates | |
1805 | TGeoHMatrix *matrix = fTRDGeometry->GetClusterMatrix(trklDet); | |
1806 | if (!matrix){ | |
1807 | LogError("invalid TRD cluster matrix in EstimateTrackDistance for detector %i", trklDet); | |
1808 | return -5; | |
1809 | } | |
1810 | matrix->LocalToMaster(xtrkl, ptrkl); | |
1811 | fTRDGeometry->RotateBack((stack/5) * 30, ptrkl, ptrkl2); // ptrkl2 now contains the global position of the outer end point of the tracklet | |
1812 | ||
1813 | // calculate parameterization of plane representing the tracklets layer | |
1814 | Double_t layer_zero_local[3] = {0., 0., 0.}; | |
1815 | Double_t layer_zero_global[3], layer_zero_global2[3]; | |
1816 | ||
1817 | matrix->LocalToMaster(layer_zero_local, layer_zero_global); | |
1818 | fTRDGeometry->RotateBack(trklDet, layer_zero_global, layer_zero_global2); // layer_zero_global2 points to chamber origin in global coords | |
1819 | ||
1820 | Double_t layer_ref_local[3] = {AliTRDgeometry::AnodePos(), 0., 0.}; | |
1821 | Double_t layer_ref_global[3], layer_ref_global2[3]; | |
1822 | ||
1823 | matrix->LocalToMaster(layer_ref_local, layer_ref_global); | |
1824 | fTRDGeometry->RotateBack(trklDet, layer_ref_global, layer_ref_global2); // layer_ref_global2 points to center anode pos within plane in global coords | |
1825 | ||
1826 | Double_t n0[3] = {layer_ref_global2[0]-layer_zero_global2[0], | |
1827 | layer_ref_global2[1]-layer_zero_global2[1], | |
1828 | layer_ref_global2[2]-layer_zero_global2[2]}; | |
1829 | ||
1830 | Double_t n_len = TMath::Sqrt(n0[0]*n0[0] + n0[1]*n0[1] + n0[2]*n0[2]); | |
1831 | if (n_len == 0.){ // This should never happen | |
1832 | //printf("<ERROR> divison by zero in estimate_track_distance!"); | |
1833 | n_len = 1.; | |
1834 | } | |
1835 | Double_t n[3] = {n0[0]/n_len, n0[1]/n_len, n0[2]/n_len}; // normal vector of plane | |
1836 | ||
1837 | Bool_t isects = TrackPlaneIntersect(refParam, layer_ref_global2, n, mag); // find intersection point between track and TRD layer | |
1838 | ||
1839 | if (isects == kFALSE){ // extrapolation fails, because track never reaches the TRD radius | |
1840 | return -1; | |
1841 | } | |
1842 | ||
1843 | Double_t m[2] = {ptrkl2[0] - layer_ref_global2[0], ptrkl2[1] - layer_ref_global2[1]}; | |
1844 | Double_t len_m = TMath::Sqrt(m[0]*m[0] + m[1]*m[1]); | |
1845 | diff_y += len_m; | |
1846 | diff_z += TMath::Abs(ptrkl2[2] - layer_ref_global2[2]); | |
1847 | nLayers++; | |
1848 | } | |
1849 | } | |
1850 | ||
1851 | if (nLayers > 0){ | |
1852 | *ydist = diff_y / nLayers; | |
1853 | *zdist = diff_z / nLayers; | |
1854 | return 0; | |
1855 | } else { | |
1856 | LogError("invalid number of contributing layers (%d) in EstimateTrackDistance()", nLayers); | |
1857 | return -4; | |
1858 | } | |
1859 | } |