CbmRoot
Loading...
Searching...
No Matches
CbmFitGlobalTracksQa.cxx
Go to the documentation of this file.
1/* Copyright (C) 2006-2015 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt
2 SPDX-License-Identifier: GPL-3.0-only
3 Authors: Dmytro Kresan, Denis Bertini [committer], Florian Uhlig */
4
5//-------------------------------------------------------------------------
6//----- CbmFitGlobalTracksQa -----
7//----- Created 07/03/06 by D. Kresan -----
8//-------------------------------------------------------------------------
10
11#include "CbmGlobalTrack.h"
12#include "CbmStsHit.h"
13#include "CbmStsPoint.h"
14#include "CbmStsTrack.h"
15#include "CbmTrdHit.h"
16#include "CbmTrdPoint.h"
17#include "CbmTrdTrack.h"
18#include "FairRootManager.h"
19#include "TClonesArray.h"
20#include "TH1.h"
21#include "TMath.h"
22#include "TVector3.h"
23
24#include <iostream>
25
26using std::cout;
27using std::endl;
28
29//________________________________________________________________________
30//
31// CbmFitGlobalTracksQa
32//
33// Task for performance analysis of the global track fit
34//
35
36
37// -----------------------------------------------------------------------
39 : FairTask()
40 , fArrayStsPoint(NULL)
41 , fArrayTrdPoint(NULL)
42 , fArrayStsHit(NULL)
43 , fArrayTrdHit(NULL)
44 , fArrayStsTrack(NULL)
45 , fArrayTrdTrack(NULL)
46 , fArrayGlbTrack(NULL)
47 , fEvents(0)
48 , fh_first_resx(NULL)
49 , fh_first_resy(NULL)
50 , fh_first_restx(NULL)
51 , fh_first_resty(NULL)
52 , fh_first_resz(NULL)
53 , fh_last_resx(NULL)
54 , fh_last_resy(NULL)
55 , fh_last_restx(NULL)
56 , fh_last_resty(NULL)
57 , fh_last_resz(NULL)
58 , fh_first_pullx(NULL)
59 , fh_first_pully(NULL)
60 , fh_first_pulltx(NULL)
61 , fh_first_pullty(NULL)
62 , fh_last_pullx(NULL)
63 , fh_last_pully(NULL)
64 , fh_last_pulltx(NULL)
65 , fh_last_pullty(NULL)
66 , fh_chi2ndf(NULL)
67{
69}
70// -----------------------------------------------------------------------
71
72
73// -----------------------------------------------------------------------
74CbmFitGlobalTracksQa::CbmFitGlobalTracksQa(const char* name, Int_t verbose)
75 : FairTask(name, verbose)
76 , fArrayStsPoint(NULL)
77 , fArrayTrdPoint(NULL)
78 , fArrayStsHit(NULL)
79 , fArrayTrdHit(NULL)
80 , fArrayStsTrack(NULL)
81 , fArrayTrdTrack(NULL)
82 , fArrayGlbTrack(NULL)
83 , fEvents(0)
84 , fh_first_resx(NULL)
85 , fh_first_resy(NULL)
86 , fh_first_restx(NULL)
87 , fh_first_resty(NULL)
88 , fh_first_resz(NULL)
89 , fh_last_resx(NULL)
90 , fh_last_resy(NULL)
91 , fh_last_restx(NULL)
92 , fh_last_resty(NULL)
93 , fh_last_resz(NULL)
94 , fh_first_pullx(NULL)
95 , fh_first_pully(NULL)
96 , fh_first_pulltx(NULL)
97 , fh_first_pullty(NULL)
98 , fh_last_pullx(NULL)
99 , fh_last_pully(NULL)
100 , fh_last_pulltx(NULL)
101 , fh_last_pullty(NULL)
102 , fh_chi2ndf(NULL)
103{
105}
106// -----------------------------------------------------------------------
107
108
109// -----------------------------------------------------------------------
114// -----------------------------------------------------------------------
115
116
117// -----------------------------------------------------------------------
119{
120 // Initialisation of the task
121
122 // Get pointer to the ROOT I/O manager
123 FairRootManager* rootMgr = FairRootManager::Instance();
124 if (NULL == rootMgr) {
125 cout << "-E- CbmFitGlobalTracksQa::Init :"
126 << " ROOT manager is not instantiated" << endl;
127 return kFATAL;
128 }
129 // Get MC points
130 fArrayStsPoint = (TClonesArray*) rootMgr->GetObject("StsPoint");
131 if (NULL == fArrayStsPoint) {
132 cout << "-E- CbmFitGlobalTracksQa::Init : "
133 << "no STS point array!" << endl;
134 }
135 fArrayTrdPoint = (TClonesArray*) rootMgr->GetObject("TrdPoint");
136 if (NULL == fArrayTrdPoint) {
137 cout << "-E- CbmFitGlobalTracksQa::Init : "
138 << "no TRD point array!" << endl;
139 }
140 // Get hit arrays
141 fArrayStsHit = (TClonesArray*) rootMgr->GetObject("StsHit");
142 if (NULL == fArrayStsHit) {
143 cout << "-W- CbmFitGlobalTracksQa::Init :"
144 << " no STS hit array" << endl;
145 }
146 fArrayTrdHit = (TClonesArray*) rootMgr->GetObject("TrdHit");
147 if (NULL == fArrayTrdHit) {
148 cout << "-W- CbmFitGlobalTracksQa::Init :"
149 << " no TRD hit array" << endl;
150 }
151 // Get track arrays
152 fArrayStsTrack = (TClonesArray*) rootMgr->GetObject("StsTrack");
153 if (NULL == fArrayStsTrack) {
154 cout << "-W- CbmFitGlobalTracksQa::Init : "
155 << "no STS track array!" << endl;
156 }
157 fArrayTrdTrack = (TClonesArray*) rootMgr->GetObject("TrdTrack");
158 if (NULL == fArrayTrdTrack) {
159 cout << "-W- CbmFitGlobalTracksQa::Init : "
160 << "no TRD track array!" << endl;
161 }
162 fArrayGlbTrack = (TClonesArray*) rootMgr->GetObject("GlobalTrack");
163 if (NULL == fArrayGlbTrack) {
164 cout << "-W- CbmFitGlobalTracksQa::Init : "
165 << "no global track array!" << endl;
166 }
167 return kSUCCESS;
168}
169// -----------------------------------------------------------------------
170
171
172// -----------------------------------------------------------------------
174{
175 // Task execution
176 if (NULL == fArrayStsPoint || NULL == fArrayTrdPoint || NULL == fArrayStsHit || NULL == fArrayTrdHit
177 || NULL == fArrayStsTrack || NULL == fArrayTrdTrack || NULL == fArrayGlbTrack)
178 return;
179
180 // Loop over the global tracks
181 CbmStsPoint* stsPoint;
182 CbmTrdPoint* trdPoint;
183 CbmStsHit* stsHit;
184 CbmTrdHit* trdHit;
185 CbmStsTrack* stsTrack;
186 CbmTrdTrack* trdTrack;
187 CbmGlobalTrack* glbTrack;
188 Int_t stsTrackIndex;
189 Int_t trdTrackIndex;
190 Int_t stsHitIndex;
191 Int_t trdHitIndex;
192 Int_t stsPointIndex;
193 Int_t trdPointIndex;
194 TVector3 pos;
195 TVector3 mom;
196
197 for (Int_t iGlbTrack = 0; iGlbTrack < fArrayGlbTrack->GetEntriesFast(); iGlbTrack++) {
198
199 // Get pointer to the global track
200 glbTrack = (CbmGlobalTrack*) fArrayGlbTrack->At(iGlbTrack);
201 if (NULL == glbTrack) continue;
202
203
204 // Chi2/NDF
205 fh_chi2ndf->Fill(glbTrack->GetChi2() / (Double_t) glbTrack->GetNDF());
206
207
208 // Performance at the first plane
209 stsTrackIndex = glbTrack->GetStsTrackIndex();
210 if (stsTrackIndex >= 0) {
211 stsTrack = (CbmStsTrack*) fArrayStsTrack->At(stsTrackIndex);
212 if (NULL != stsTrack) {
213 stsPointIndex = -1;
214 if (stsTrack->GetNofStsHits()) {
215 stsHitIndex = stsTrack->GetStsHitIndex(0);
216 stsHit = (CbmStsHit*) fArrayStsHit->At(stsHitIndex);
217 if (NULL != stsHit) {
218 stsPointIndex = stsHit->GetRefId();
219 }
220 }
221 else {
222 cout << "-W- CbmFitGlobalTracksQa::Exec : "
223 << "STS track " << stsTrackIndex << " is empty!" << endl;
224 }
225
226 if (stsPointIndex >= 0) {
227 stsPoint = (CbmStsPoint*) fArrayStsPoint->At(stsPointIndex);
228 if (NULL != stsPoint) {
229 stsPoint->Position(pos);
230 stsPoint->Momentum(mom);
231 fh_first_resx->Fill(glbTrack->GetParamFirst()->GetX() - pos.X());
232 fh_first_resy->Fill(glbTrack->GetParamFirst()->GetY() - pos.Y());
233 fh_first_resz->Fill(glbTrack->GetParamFirst()->GetZ() - pos.Z());
234 fh_first_restx->Fill(glbTrack->GetParamFirst()->GetTx() - mom.X() / mom.Z());
235 fh_first_resty->Fill(glbTrack->GetParamFirst()->GetTy() - mom.Y() / mom.Z());
236 fh_first_pullx->Fill((glbTrack->GetParamFirst()->GetX() - pos.X())
237 / TMath::Sqrt(glbTrack->GetParamFirst()->GetCovariance(0, 0)));
238 fh_first_pully->Fill((glbTrack->GetParamFirst()->GetY() - pos.Y())
239 / TMath::Sqrt(glbTrack->GetParamFirst()->GetCovariance(1, 1)));
240 fh_first_pulltx->Fill((glbTrack->GetParamFirst()->GetTx() - mom.X() / mom.Z())
241 / TMath::Sqrt(glbTrack->GetParamFirst()->GetCovariance(2, 2)));
242 fh_first_pullty->Fill((glbTrack->GetParamFirst()->GetTy() - mom.Y() / mom.Z())
243 / TMath::Sqrt(glbTrack->GetParamFirst()->GetCovariance(3, 3)));
244 }
245 }
246 }
247 }
248 else {
249 trdTrackIndex = glbTrack->GetTrdTrackIndex();
250 trdTrack = (CbmTrdTrack*) fArrayTrdTrack->At(trdTrackIndex);
251 if (NULL != trdTrack) {
252 trdPointIndex = -1;
253 if (trdTrack->GetNofHits()) {
254 trdHitIndex = trdTrack->GetHitIndex(0);
255 trdHit = (CbmTrdHit*) fArrayTrdHit->At(trdHitIndex);
256 if (NULL != trdHit) {
257 trdPointIndex = trdHit->GetRefId();
258 }
259 }
260 else {
261 cout << "-W- CbmFitGlobalTracksQa::Exec : "
262 << "TRD track " << trdTrackIndex << " is empty!" << endl;
263 }
264
265 if (trdPointIndex >= 0) {
266 trdPoint = (CbmTrdPoint*) fArrayTrdPoint->At(trdPointIndex);
267 if (NULL != trdPoint) {
268 trdPoint->Position(pos);
269 trdPoint->Momentum(mom);
270 fh_first_resx->Fill(glbTrack->GetParamFirst()->GetX() - pos.X());
271 fh_first_resy->Fill(glbTrack->GetParamFirst()->GetY() - pos.Y());
272 fh_first_resz->Fill(glbTrack->GetParamFirst()->GetZ() - pos.Z());
273 fh_first_restx->Fill(glbTrack->GetParamFirst()->GetTx() - mom.X() / mom.Z());
274 fh_first_resty->Fill(glbTrack->GetParamFirst()->GetTy() - mom.Y() / mom.Z());
275 fh_first_pullx->Fill((glbTrack->GetParamFirst()->GetX() - pos.X())
276 / TMath::Sqrt(glbTrack->GetParamFirst()->GetCovariance(0, 0)));
277 fh_first_pully->Fill((glbTrack->GetParamFirst()->GetY() - pos.Y())
278 / TMath::Sqrt(glbTrack->GetParamFirst()->GetCovariance(1, 1)));
279 fh_first_pulltx->Fill((glbTrack->GetParamFirst()->GetTx() - mom.X() / mom.Z())
280 / TMath::Sqrt(glbTrack->GetParamFirst()->GetCovariance(2, 2)));
281 fh_first_pullty->Fill((glbTrack->GetParamFirst()->GetTy() - mom.Y() / mom.Z())
282 / TMath::Sqrt(glbTrack->GetParamFirst()->GetCovariance(3, 3)));
283 }
284 }
285 }
286 }
287
288
289 // Performance at the last plane
290 trdTrackIndex = glbTrack->GetTrdTrackIndex();
291 if (trdTrackIndex >= 0) {
292 trdTrack = (CbmTrdTrack*) fArrayTrdTrack->At(trdTrackIndex);
293 if (NULL != trdTrack) {
294 trdPointIndex = -1;
295 if (trdTrack->GetNofHits()) {
296 trdHitIndex = trdTrack->GetHitIndex(trdTrack->GetNofHits() - 1);
297 trdHit = (CbmTrdHit*) fArrayTrdHit->At(trdHitIndex);
298 if (NULL != trdHit) {
299 trdPointIndex = trdHit->GetRefId();
300 }
301 }
302 else {
303 cout << "-W- CbmFitGlobalTracksQa::Exec : "
304 << "TRD track " << trdTrackIndex << " is empty!" << endl;
305 }
306
307 if (trdPointIndex >= 0) {
308 trdPoint = (CbmTrdPoint*) fArrayTrdPoint->At(trdPointIndex);
309 if (NULL != trdPoint) {
310 trdPoint->Position(pos);
311 trdPoint->Momentum(mom);
312 fh_last_resx->Fill(glbTrack->GetParamLast()->GetX() - pos.X());
313 fh_last_resy->Fill(glbTrack->GetParamLast()->GetY() - pos.Y());
314 fh_last_resz->Fill(glbTrack->GetParamLast()->GetZ() - pos.Z());
315 fh_last_restx->Fill(glbTrack->GetParamLast()->GetTx() - mom.X() / mom.Z());
316 fh_last_resty->Fill(glbTrack->GetParamLast()->GetTy() - mom.Y() / mom.Z());
317 fh_last_pullx->Fill((glbTrack->GetParamLast()->GetX() - pos.X())
318 / TMath::Sqrt(glbTrack->GetParamLast()->GetCovariance(0, 0)));
319 fh_last_pully->Fill((glbTrack->GetParamLast()->GetY() - pos.Y())
320 / TMath::Sqrt(glbTrack->GetParamLast()->GetCovariance(1, 1)));
321 fh_last_pulltx->Fill((glbTrack->GetParamLast()->GetTx() - mom.X() / mom.Z())
322 / TMath::Sqrt(glbTrack->GetParamLast()->GetCovariance(2, 2)));
323 fh_last_pullty->Fill((glbTrack->GetParamLast()->GetTy() - mom.Y() / mom.Z())
324 / TMath::Sqrt(glbTrack->GetParamLast()->GetCovariance(3, 3)));
325 }
326 }
327 }
328 }
329 else {
330 stsTrackIndex = glbTrack->GetStsTrackIndex();
331 stsTrack = (CbmStsTrack*) fArrayStsTrack->At(stsTrackIndex);
332 if (NULL != stsTrack) {
333 stsPointIndex = -1;
334 if (stsTrack->GetNofStsHits()) {
335 stsHitIndex = stsTrack->GetStsHitIndex(stsTrack->GetNofStsHits() - 1);
336 stsHit = (CbmStsHit*) fArrayStsHit->At(stsHitIndex);
337 if (NULL != stsHit) {
338 stsPointIndex = stsHit->GetRefId();
339 }
340 }
341 else {
342 cout << "-W- CbmFitGlobalTracksQa::Exec : "
343 << "STS track " << stsTrackIndex << " is empty!" << endl;
344 }
345
346 if (stsPointIndex >= 0) {
347 stsPoint = (CbmStsPoint*) fArrayStsPoint->At(stsPointIndex);
348 if (NULL != stsPoint) {
349 stsPoint->Position(pos);
350 stsPoint->Momentum(mom);
351 fh_last_resx->Fill(glbTrack->GetParamLast()->GetX() - pos.X());
352 fh_last_resy->Fill(glbTrack->GetParamLast()->GetY() - pos.Y());
353 fh_last_resz->Fill(glbTrack->GetParamLast()->GetZ() - pos.Z());
354 fh_last_restx->Fill(glbTrack->GetParamLast()->GetTx() - mom.X() / mom.Z());
355 fh_last_resty->Fill(glbTrack->GetParamLast()->GetTy() - mom.Y() / mom.Z());
356 fh_last_pullx->Fill((glbTrack->GetParamLast()->GetX() - pos.X())
357 / TMath::Sqrt(glbTrack->GetParamLast()->GetCovariance(0, 0)));
358 fh_last_pully->Fill((glbTrack->GetParamLast()->GetY() - pos.Y())
359 / TMath::Sqrt(glbTrack->GetParamLast()->GetCovariance(1, 1)));
360 fh_last_pulltx->Fill((glbTrack->GetParamLast()->GetTx() - mom.X() / mom.Z())
361 / TMath::Sqrt(glbTrack->GetParamLast()->GetCovariance(2, 2)));
362 fh_last_pullty->Fill((glbTrack->GetParamLast()->GetTy() - mom.Y() / mom.Z())
363 / TMath::Sqrt(glbTrack->GetParamLast()->GetCovariance(3, 3)));
364 }
365 }
366 }
367 }
368 }
369
370
371 fEvents += 1;
372}
373
374
375// -----------------------------------------------------------------------
377{
378 // Finish of the task execution
380}
381// -----------------------------------------------------------------------
382
383
384// -----------------------------------------------------------------------
386{
387 // Create control histogramms
388
389 // Residuals at first TRD layer
390 fh_first_resx = new TH1F("h_glb_first_resx", "x residual at first TRD layer", 200, -1., 1.);
391 fh_first_resy = new TH1F("h_glb_first_resy", "y residual at first TRD layer", 200, -1., 1.);
392 fh_first_restx = new TH1F("h_glb_first_restx", "t_{x} residual at first TRD layer", 200, -0.01, 0.01);
393 fh_first_resty = new TH1F("h_glb_first_resty", "t_{y} residual at first TRD layer", 200, -0.01, 0.01);
394 fh_first_resz = new TH1F("h_glb_first_resz", "z residual at first TRD layer", 100, -0.5, 0.5);
395
396 // Residuals at last TRD layer
397 fh_last_resx = new TH1F("h_glb_last_resx", "x residual at last TRD layer", 200, -1., 1.);
398 fh_last_resy = new TH1F("h_glb_last_resy", "y residual at last TRD layer", 200, -1., 1.);
399 fh_last_restx = new TH1F("h_glb_last_restx", "t_{x} residual at last TRD layer", 200, -0.01, 0.01);
400 fh_last_resty = new TH1F("h_glb_last_resty", "t_{y} residual at last TRD layer", 200, -0.01, 0.01);
401 fh_last_resz = new TH1F("h_glb_last_resz", "z residual at last TRD layer", 100, -0.5, 0.5);
402
403 // Pulls at first TRD layer
404 fh_first_pullx = new TH1F("h_glb_first_pullx", "x pull at first TRD layer", 200, -10., 10.);
405 fh_first_pully = new TH1F("h_glb_first_pully", "y pull at first TRD layer", 200, -10., 10.);
406 fh_first_pulltx = new TH1F("h_glb_first_pulltx", "t_{x} pull at first TRD layer", 200, -10., 10.);
407 fh_first_pullty = new TH1F("h_glb_first_pullty", "t_{y} pull at first TRD layer", 200, -10., 10.);
408
409 // Pulls at last TRD layer
410 fh_last_pullx = new TH1F("h_glb_last_pullx", "x pull at last TRD layer", 200, -10., 10.);
411 fh_last_pully = new TH1F("h_glb_last_pully", "y pull at last TRD layer", 200, -10., 10.);
412 fh_last_pulltx = new TH1F("h_glb_last_pulltx", "t_{x} pull at last TRD layer", 200, -10., 10.);
413 fh_last_pullty = new TH1F("h_glb_last_pullty", "t_{y} pull at last TRD layer", 200, -10., 10.);
414
415 // Chi2/NDF of track fit
416 fh_chi2ndf = new TH1F("h_glb_chi2ndf", "#chi^{2}/NDF of track fit", 500, 0., 50.);
417}
418// -----------------------------------------------------------------------
419
420
421// -----------------------------------------------------------------------
423{
424 // Write control histogramms to file
425
426 // Residuals at first and last TRD layer
427 fh_first_resx->Scale(1. / fEvents);
428 fh_first_resy->Scale(1. / fEvents);
429 fh_first_restx->Scale(1. / fEvents);
430 fh_first_resty->Scale(1. / fEvents);
431 fh_last_resx->Scale(1. / fEvents);
432 fh_last_resy->Scale(1. / fEvents);
433 fh_last_restx->Scale(1. / fEvents);
434 fh_last_resty->Scale(1. / fEvents);
435 fh_first_resx->Write();
436 fh_first_resy->Write();
437 fh_first_restx->Write();
438 fh_first_resty->Write();
439 fh_first_resz->Write();
440 fh_last_resx->Write();
441 fh_last_resy->Write();
442 fh_last_restx->Write();
443 fh_last_resty->Write();
444 fh_last_resz->Write();
445
446 // Pulls at first and last TRD layer
447 fh_first_pullx->Scale(1. / fEvents);
448 fh_first_pully->Scale(1. / fEvents);
449 fh_first_pulltx->Scale(1. / fEvents);
450 fh_first_pullty->Scale(1. / fEvents);
451 fh_last_pullx->Scale(1. / fEvents);
452 fh_last_pully->Scale(1. / fEvents);
453 fh_last_pulltx->Scale(1. / fEvents);
454 fh_last_pullty->Scale(1. / fEvents);
455 fh_first_pullx->Write();
456 fh_first_pully->Write();
457 fh_first_pulltx->Write();
458 fh_first_pullty->Write();
459 fh_last_pullx->Write();
460 fh_last_pully->Write();
461 fh_last_pulltx->Write();
462 fh_last_pullty->Write();
463
464 // Chi2/NDF of fit
465 fh_chi2ndf->Scale(1. / fEvents);
466 fh_chi2ndf->Write();
467}
468// -----------------------------------------------------------------------
469
470
static FairRootManager * rootMgr
ClassImp(CbmFitGlobalTracksQa)
Data class for a reconstructed hit in the STS.
Data class for STS tracks.
Class for hits in TRD detector.
virtual InitStatus Init()
virtual void Exec(Option_t *option="")
const FairTrackParam * GetParamLast() const
int32_t GetStsTrackIndex() const
double GetChi2() const
int32_t GetNDF() const
const FairTrackParam * GetParamFirst() const
int32_t GetTrdTrackIndex() const
int32_t GetRefId() const
Definition CbmHit.h:73
data class for a reconstructed 3-d hit in the STS
Definition CbmStsHit.h:35
int32_t GetStsHitIndex(int32_t iHit) const
int32_t GetNofStsHits() const
Definition CbmStsTrack.h:93
virtual int32_t GetNofHits() const
Definition CbmTrack.h:58
int32_t GetHitIndex(int32_t iHit) const
Definition CbmTrack.h:59
data class for a reconstructed Energy-4D measurement in the TRD
Definition CbmTrdHit.h:40