BmnRoot
Loading...
Searching...
No Matches
CbmL1.cxx
Go to the documentation of this file.
1/*
2 *====================================================================
3 *
4 * CBM Level 1 Reconstruction
5 *
6 * Authors: I.Kisel, S.Gorbunov
7 *
8 * e-mail : ikisel@kip.uni-heidelberg.de
9 *
10 *====================================================================
11 *
12 * L1 main program
13 *
14 *====================================================================
15 */
16#include "CbmL1.h"
17
18#include "CbmGeoStsPar.h"
19#include "CbmKF.h"
20#include "CbmMvdHit.h"
21#include "CbmStsDigiPar.h"
22#include "CbmStsFindTracks.h"
23#include "CbmStsHit.h"
24#include "CbmStsSector.h"
25#include "CbmStsSensor.h" // for field approx.
26#include "CbmStsStation.h"
27#include "FairRootManager.h"
28#include "FairRunAna.h"
29#include "FairRuntimeDb.h"
30#include "L1Algo/L1Algo.h"
31#include "L1Algo/L1Branch.h"
32#include "L1Algo/L1Field.h"
33#include "L1Algo/L1StsHit.h"
34#include "TFile.h"
35#include "TMatrixD.h"
36#include "TProfile2D.h"
37#include "TROOT.h"
38#include "TRandom3.h"
39#include "TVectorD.h"
40
41#include <fstream>
42#include <iomanip>
43#include <iostream>
44
45using namespace std;
46
47static L1Algo algo_static _fvecalignment;
48
49CbmL1* CbmL1::fInstance = 0;
50
51CbmL1::CbmL1()
52 : algo(0)
53 , // for access to L1 Algorithm from L1::Instance
54 PF(0)
55 , vRTracks()
56 , // reconstructed tracks
57 vHitStore()
58 , NStation(0)
59 , NMvdStations(0)
60 , NStsStations(0)
61 , // number of detector stations (all\sts\mvd)
62 fPerformance(0)
63 , fSTAPDataMode(0)
64 , fSTAPDataDir("")
65 ,
66
67 fTrackingLevel(2)
68 , // really doesn't used
69 fMomentumCutOff(0.1)
70 , // really doesn't used
71 fGhostSuppression(1)
72 , // really doesn't used
73 fUseMVD(0)
74 , // really doesn't used
75
76 StsDigi()
77 , PrimVtx()
78 ,
79
80 listMCTracks(0)
81 , listStsPts(0)
82 , listStsDigi(0)
83 , listStsClusters(0)
84 , listStsHits(0)
85 ,
86
87 listMvdPts(0)
88 , listMvdHits(0)
89 , listMvdHitMatches(0)
90 , vStsHits()
91 , vMCPoints()
92 , vMCTracks()
93 , vHitMCRef()
94 , vRParticles()
95 , vMCParticles()
96 , MCtoRParticleId()
97 , RtoMCParticleId()
98 , histodir(0)
99 , fFindParticlesMode()
100 , fMatBudgetFileName("")
101 , fExtrapolateToTheEndOfSTS(false)
102{
103 if (!fInstance)
104 fInstance = this;
105 PF = new CbmL1ParticlesFinder();
106 StsDigi = CbmStsDigiScheme::Instance(); // AZ
107}
108
109CbmL1::CbmL1(const char* name,
110 Int_t iVerbose,
111 Int_t _fPerformance,
112 int fSTAPDataMode_,
113 TString fSTAPDataDir_,
114 int findParticleMode_)
115 : BmnTask(name, iVerbose)
116 , algo(0)
117 , // for access to L1 Algorithm from L1::Instance
118 PF(0)
119 ,
120
121 vRTracks()
122 , // reconstructed tracks
123 vHitStore()
124 , NStation(0)
125 , NMvdStations(0)
126 , NStsStations(0)
127 , // number of detector stations (all\sts\mvd)
128 fPerformance(_fPerformance)
129 , fSTAPDataMode(fSTAPDataMode_)
130 , fSTAPDataDir(fSTAPDataDir_)
131 ,
132
133 fTrackingLevel(2)
134 , // really doesn't used
135 fMomentumCutOff(0.1)
136 , // really doesn't used
137 fGhostSuppression(1)
138 , // really doesn't used
139 fUseMVD(0)
140 , // really doesn't used
141
142 StsDigi()
143 , PrimVtx()
144 ,
145
146 listMCTracks(0)
147 , listStsPts(0)
148 , listStsDigi(0)
149 , listStsClusters(0)
150 , listStsHits(0)
151 ,
152
153 listMvdPts(0)
154 , listMvdHits(0)
155 , listMvdHitMatches(0)
156 , vStsHits()
157 , vMCPoints()
158 , vMCTracks()
159 , vHitMCRef()
160 , vRParticles()
161 , vMCParticles()
162 , MCtoRParticleId()
163 , RtoMCParticleId()
164 , histodir(0)
165 , fFindParticlesMode(findParticleMode_)
166 , fMatBudgetFileName("")
167 , fExtrapolateToTheEndOfSTS(false)
168{
169 if (!fInstance)
170 fInstance = this;
171 PF = new CbmL1ParticlesFinder();
172 StsDigi = CbmStsDigiScheme::Instance(); // AZ
173}
174
175CbmL1::~CbmL1()
176{
177 if (fInstance == this)
178 fInstance = 0;
179 delete PF;
180}
181
182void CbmL1::SetParContainers()
183{
185 // FairRunAna* ana = FairRunAna::Instance();
186 // FairRuntimeDb* rtdb=ana->GetRuntimeDb();
187 FairRuntimeDb* rtdb = FairRuntimeDb::instance();
189 rtdb->getContainer("FairGeoPassivePar");
190 rtdb->getContainer("CbmGeoStsPar");
191 rtdb->getContainer("CbmStsDigiPar");
192}
193
194InitStatus CbmL1::ReInit()
195{
196 // AZ StsDigi.Clear();
197 StsDigi->Clear();
198 SetParContainers();
199 return Init();
200}
201
202InitStatus CbmL1::Init()
203{
204 if (fVerbose > 1) {
205 char y[20] = " [0;33;44m"; // yellow
206 char Y[20] = " [1;33;44m"; // yellow bold
207 char W[20] = " [1;37;44m"; // white bold
208 char o[20] = " [0m\n"; // color off
209 Y[0] = y[0] = W[0] = o[0] = 0x1B; // escape character
210
211 cout << endl << endl;
212 cout << " " << W << " " << o;
213 cout << " " << W << " ===////====================================================== " << o;
214 cout << " " << W << " = = " << o;
215 cout << " " << W << " = " << Y << "L1 on-line finder" << W << " = "
216 << o;
217 cout << " " << W << " = = " << o;
218 cout << " " << W << " = " << W << "Cellular Automaton 3.1 Vector" << y << " with " << W << "KF Quadro" << y
219 << " technology" << W << " = " << o;
220 cout << " " << W << " = = " << o;
221 cout << " " << W << " = " << y << "Designed for CBM collaboration" << W << " = "
222 << o;
223 cout << " " << W << " = " << y << "All rights reserved" << W << " = "
224 << o;
225 cout << " " << W << " = = " << o;
226 cout << " " << W << " ========================================================////= " << o;
227 cout << " " << W << " " << o;
228 cout << endl << endl;
229 }
230
231 FairRootManager* fManger = FairRootManager::Instance();
232
234 // FairRunAna *Run = FairRunAna::Instance();
235 // FairRuntimeDb *RunDB = Run->GetRuntimeDb();
236 FairRuntimeDb* RunDB = FairRuntimeDb::instance();
238 {
239 CbmGeoStsPar* StsPar = L1_DYNAMIC_CAST<CbmGeoStsPar*>(RunDB->findContainer("CbmGeoStsPar"));
240 CbmStsDigiPar* digiPar = L1_DYNAMIC_CAST<CbmStsDigiPar*>(RunDB->findContainer("CbmStsDigiPar"));
241 // AZ StsDigi.Init(StsPar, digiPar);
242 StsDigi->Init(StsPar, digiPar);
243 }
244 {
245 // AZ-040623 fUseMVD = 1;
246 fUseMVD = 0; // AZ-040623
248 // CbmStsFindTracks * FindTask = L1_DYNAMIC_CAST<CbmStsFindTracks*>( Run->GetTask("STSFindTracks") );
249 // if( FindTask ) fUseMVD = FindTask->MvdUsage();
251 }
252
253 histodir = gROOT->mkdir("L1");
254
255 if (fPerformance) {
256 listMCTracks = L1_DYNAMIC_CAST<TClonesArray*>(fManger->GetObject("MCTrack"));
257 listStsPts = L1_DYNAMIC_CAST<TClonesArray*>(fManger->GetObject("StsPoint"));
258 listStsClusters = L1_DYNAMIC_CAST<TClonesArray*>(fManger->GetObject("StsCluster"));
259 listStsDigi = L1_DYNAMIC_CAST<TClonesArray*>(fManger->GetObject("StsDigi"));
260 } else {
261 listMCTracks = 0;
262 listStsPts = 0;
263 listStsClusters = 0;
264 listStsDigi = 0;
265 }
266 listStsHits = L1_DYNAMIC_CAST<TClonesArray*>(fManger->GetObject("StsHit"));
267
268 if (fPerformance) {
269 if (!fUseMVD) {
270 listMvdPts = 0;
271 listMvdHitMatches = 0;
272 } else {
273 listMvdPts = L1_DYNAMIC_CAST<TClonesArray*>(fManger->GetObject("MvdPoint"));
274 listMvdHitMatches = L1_DYNAMIC_CAST<TClonesArray*>(fManger->GetObject("MvdHitMatch"));
275 }
276 } else {
277 listMvdPts = 0;
278 listMvdHitMatches = 0;
279 }
280 if (!fUseMVD) {
281 listMvdHits = 0;
282 } else {
283 listMvdHits = L1_DYNAMIC_CAST<TClonesArray*>(fManger->GetObject("MvdHit"));
284 }
285
286 NMvdStations = 0;
287 NStsStations = 0;
288 NStation = 0;
289 algo = &algo_static;
290
291 fscal geo[10000];
292
293 int ind = 0;
294 for (int i = 0; i < 3; i++) {
295 Double_t point[3] = {0, 0, 2.5 * i};
296 Double_t B[3] = {0, 0, 0};
297 if (CbmKF::Instance()->GetMagneticField())
298 CbmKF::Instance()->GetMagneticField()->GetFieldValue(point, B);
299 geo[ind++] = 2.5 * i;
300 geo[ind++] = B[0];
301 geo[ind++] = B[1];
302 geo[ind++] = B[2];
303 }
304
305 NMvdStations = (fUseMVD) ? CbmKF::Instance()->vMvdMaterial.size() : 0;
306 // AZ NStsStations = StsDigi.GetNStations();
307 NStsStations = StsDigi->GetNStations();
308 NStation = NMvdStations + NStsStations;
309 geo[ind++] = NStation;
310 geo[ind++] = NMvdStations;
311
312 // field
313 const int M = 5; // polinom order
314 const int N = (M + 1) * (M + 2) / 2;
315
316 // { // field at the z=0 plane
317 // const double Xmax = 10, Ymax = 10;
318 // const double z = 0.;
319 // double dx = 1.; // step for the field approximation
320 // double dy = 1.;
321 // if( dx > Xmax/N/2 ) dx = Xmax/N/4.;
322 // if( dy > Ymax/N/2 ) dy = Ymax/N/4.;
323 //
324 // TMatrixD A(N,N);
325 // TVectorD b0(N), b1(N), b2(N);
326 // for( int i=0; i<N; i++){
327 // for( int j=0; j<N; j++) A(i,j)==0.;
328 // b0(i)=b1(i)=b2(i) = 0.;
329 // }
330 // for( double x=-Xmax; x<=Xmax; x+=dx )
331 // for( double y=-Ymax; y<=Ymax; y+=dy )
332 // {
333 // double r = sqrt(fabs(x*x/Xmax/Xmax+y/Ymax*y/Ymax));
334 // if( r>1. ) continue;
335 // Double_t w = 1./(r*r+1);
336 // Double_t p[3] = { x, y, z};
337 // Double_t B[3] = {0.,0.,0.};
338 // CbmKF::Instance()->GetMagneticField()->GetFieldValue(p, B);
339 // TVectorD m(N);
340 // m(0)=1;
341 // for( int i=1; i<=M; i++){
342 // int k = (i-1)*(i)/2;
343 // int l = i*(i+1)/2;
344 // for( int j=0; j<i; j++ ) m(l+j) = x*m(k+j);
345 // m(l+i) = y*m(k+i-1);
346 // }
347 //
348 // TVectorD mt = m;
349 // for( int i=0; i<N; i++){
350 // for( int j=0; j<N;j++) A(i,j)+=w*m(i)*m(j);
351 // b0(i)+=w*B[0]*m(i);
352 // b1(i)+=w*B[1]*m(i);
353 // b2(i)+=w*B[2]*m(i);
354 // }
355 // }
356 // double det;
357 // A.Invert(&det);
358 // TVectorD c0 = A*b0, c1 = A*b1, c2 = A*b2;
359 //
360 // targetFieldSlice = new L1FieldSlice;
361 // for(int i=0; i<N; i++){
362 // targetFieldSlice->cx[i] = c0(i);
363 // targetFieldSlice->cy[i] = c1(i);
364 // targetFieldSlice->cz[i] = c2(i);
365 // }
366 //
367 // } // target field
368
369 for (Int_t ist = 0; ist < NStation; ist++) {
370 double C[3][N];
371 double z = 0;
372 double Xmax, Ymax;
373 if (ist < NMvdStations) {
374 CbmKFTube& t = CbmKF::Instance()->vMvdMaterial[ist];
375 geo[ind++] = t.z;
376 geo[ind++] = t.dz;
377 geo[ind++] = t.r;
378 geo[ind++] = t.R;
379 geo[ind++] = t.RadLength;
380 fscal f_phi = 0, f_sigma = 5.e-4, b_phi = 3.14159265358 / 2., b_sigma = 5.e-4;
381 f_phi = 0;
382 f_sigma = 5.e-4;
383 b_phi = 3.14159265358 / 2.;
384 b_sigma = 5.e-4;
385 geo[ind++] = f_phi;
386 geo[ind++] = f_sigma;
387 geo[ind++] = b_phi;
388 geo[ind++] = b_sigma;
389 z = t.z;
390 Xmax = Ymax = t.R;
391 } else {
392 // AZ CbmStsStation *st = StsDigi.GetStation(ist - NMvdStations);
393 CbmStsStation* st = StsDigi->GetStation(ist - NMvdStations);
394
395 geo[ind++] = st->GetZ();
396 geo[ind++] = st->GetD();
397 geo[ind++] = st->GetRmin();
398 geo[ind++] = st->GetRmax();
399 geo[ind++] = st->GetRadLength();
400
401 CbmStsSector* sector = st->GetSector(0);
402 fscal f_phi, f_sigma, b_phi, b_sigma; // angle and sigma front/back side
403 f_phi = sector->GetRotation();
404 b_phi = sector->GetRotation();
405 if (sector->GetType() == 1) {
406 b_phi += 3.14159265358 / 2.;
407 b_sigma = sector->GetDy() / sqrt(12.);
408 f_sigma = b_sigma; // CHECKME
409 } else {
410 f_phi += sector->GetStereoF();
411 b_phi += sector->GetStereoB();
412 f_sigma = sector->GetDx() / TMath::Sqrt(12);
413 // AZ - use "real" errors
414 if (f_sigma < 0.01)
415 f_sigma = 0.01; // Si
416 else
417 f_sigma = 0.015; // GEMs
418 // AZ
419 b_sigma = f_sigma;
420 }
421 // AZ f_sigma *= cos(f_phi); // TODO: think about this
422 // AZ b_sigma *= cos(b_phi);
423 geo[ind++] = f_phi;
424 geo[ind++] = f_sigma;
425 geo[ind++] = b_phi;
426 geo[ind++] = b_sigma;
427 z = st->GetZ();
428
429 Xmax = -100;
430 Ymax = -100;
431
432 double x, y;
433 for (int isec = 0; isec < st->GetNSectors(); isec++) {
434 CbmStsSector* sect = L1_DYNAMIC_CAST<CbmStsSector*>(st->GetSector(isec));
435 for (int isen = 0; isen < sect->GetNSensors(); isen++) {
436 x = sect->GetSensor(isen)->GetX0() + sect->GetSensor(isen)->GetLx() / 2.;
437 y = sect->GetSensor(isen)->GetY0() + sect->GetSensor(isen)->GetLy() / 2.;
438 if (x > Xmax)
439 Xmax = x;
440 if (y > Ymax)
441 Ymax = y;
442 }
443 }
444 // cout << "Station "<< ist << ", Xmax " << Xmax<<", Ymax" << Ymax<<endl;
445 if (fVerbose)
446 cout << " AZ - station " << ist << " " << f_phi * TMath::RadToDeg() << " " << b_phi * TMath::RadToDeg()
447 << " " << Xmax << " " << Ymax << endl;
448 Xmax = TMath::Abs(Xmax); // AZ-040623
449 Ymax = TMath::Abs(Ymax); // AZ-040623
450 }
451
452 double dx = 1.; // step for the field approximation
453 double dy = 1.;
454
455 if (dx > Xmax / N / 2)
456 dx = Xmax / N / 4.;
457 if (dy > Ymax / N / 2)
458 dy = Ymax / N / 4.;
459
460 for (int i = 0; i < 3; i++)
461 for (int k = 0; k < N; k++)
462 C[i][k] = 0;
463 TMatrixD A(N, N);
464 TVectorD b0(N), b1(N), b2(N);
465 for (int i = 0; i < N; i++) {
466 for (int j = 0; j < N; j++)
467 A(i, j) = 0.;
468 b0(i) = b1(i) = b2(i) = 0.;
469 }
470 for (double x = -Xmax; x <= Xmax; x += dx)
471 for (double y = -Ymax; y <= Ymax; y += dy) {
472 double r = sqrt(fabs(x * x / Xmax / Xmax + y / Ymax * y / Ymax));
473 if (r > 1.)
474 continue;
475 Double_t w = 1. / (r * r + 1);
476 Double_t p[3] = {x, y, z};
477 Double_t B[3] = {0., 0., 0.};
478 CbmKF::Instance()->GetMagneticField()->GetFieldValue(p, B);
479 TVectorD m(N);
480 m(0) = 1;
481 for (int i = 1; i <= M; i++) {
482 int k = (i - 1) * (i) / 2;
483 int l = i * (i + 1) / 2;
484 for (int j = 0; j < i; j++)
485 m(l + j) = x * m(k + j);
486 m(l + i) = y * m(k + i - 1);
487 }
488
489 TVectorD mt = m;
490 for (int i = 0; i < N; i++) {
491 for (int j = 0; j < N; j++)
492 A(i, j) += w * m(i) * m(j);
493 b0(i) += w * B[0] * m(i);
494 b1(i) += w * B[1] * m(i);
495 b2(i) += w * B[2] * m(i);
496 }
497 }
498 double det;
499 A.Invert(&det);
500 TVectorD c0 = A * b0, c1 = A * b1, c2 = A * b2;
501 for (int i = 0; i < N; i++) {
502 C[0][i] = c0(i);
503 C[1][i] = c1(i);
504 C[2][i] = c2(i);
505 }
506
507 geo[ind++] = N;
508 for (int k = 0; k < 3; k++) {
509 for (int j = 0; j < N; j++)
510 geo[ind++] = C[k][j];
511 }
512 }
513
514 geo[ind++] = fTrackingLevel;
515 geo[ind++] = fMomentumCutOff;
516 geo[ind++] = fGhostSuppression;
517
518 {
519 if (fSTAPDataMode % 2 == 1) { // 1,3
520 WriteSTAPGeoData(static_cast<void*>(geo), ind);
521 };
522 if (fSTAPDataMode >= 2) { // 2,3
523 int ind2;
524 ReadSTAPGeoData(static_cast<void*>(geo), ind2);
525 if (ind2 != ind)
526 cout << "-E- CbmL1: Read geometry from file " << fSTAPDataDir + "geo_algo.txt was NOT successful."
527 << endl;
528 };
529 }
530
531 algo->Init(geo);
532
533 algo->fRadThick.resize(algo->NStations);
534 // MVD does not use map
535 for (int iSta = 0; iSta < algo->NMvdStations; iSta++) {
536 algo->fRadThick[iSta].SetBins(1, 100); // mvd should be in +-100 cm square
537 algo->fRadThick[iSta].table.resize(1);
538 algo->fRadThick[iSta].table[0].resize(1);
539 algo->fRadThick[iSta].table[0][0] = algo->vStations[iSta].materialInfo.RadThick[0];
540 }
541
542 // Read STS Radiation Thickness table
543 if (fMatBudgetFileName != "") {
544 TFile* oldfile = gFile;
545 TFile* rlFile = new TFile(fMatBudgetFileName);
546
547 cout << "STS Material budget file is " << fMatBudgetFileName << "." << endl;
548 for (int j = 0, iSta = algo->NMvdStations; iSta < algo->NStations; iSta++, j++) {
549 TString name = "Radiation Thickness [%]";
550 name += ", Station";
551 name += j + 1;
552
553 TProfile2D* hStaRadLen = (TProfile2D*)rlFile->Get(name);
554 if (!hStaRadLen) {
555 cout << "L1: incorrect " << fMatBudgetFileName << " file. No " << name << endl;
556 exit(1);
557 }
558
559 const int NBins = hStaRadLen->GetNbinsX(); // should be same in Y
560 const float RMax = hStaRadLen->GetXaxis()->GetXmax(); // should be same as min
561 algo->fRadThick[iSta].SetBins(NBins, RMax); // TODO
562 algo->fRadThick[iSta].table.resize(NBins);
563
564 for (int iB = 0; iB < NBins; iB++) {
565 algo->fRadThick[iSta].table[iB].resize(NBins);
566 for (int iB2 = 0; iB2 < NBins; iB2++) {
567 algo->fRadThick[iSta].table[iB][iB2] =
568 0.01 * hStaRadLen->GetBinContent(iB, iB2); // 0.0034;//0.003209;//
569 if (algo->fRadThick[iSta].table[iB][iB2] < algo->vStations[iSta].materialInfo.RadThick[0])
570 algo->fRadThick[iSta].table[iB][iB2] = algo->vStations[iSta].materialInfo.RadThick[0];
571 // cout << iB << " " << iB2 << " " << algo->fRadThick[iSta].table[iB][iB2] << endl; // dbg
572 }
573 }
574 }
575 rlFile->Close();
576 rlFile->Delete();
577 gFile = oldfile;
578 } else {
579 if (fVerbose)
580 cout << "No material budget file is found. Homogenious budget will be used" << endl;
581 for (int iSta = algo->NMvdStations; iSta < algo->NStations; iSta++) {
582 // cout << iSta << endl;
583 algo->fRadThick[iSta].SetBins(1, 100); // mvd should be in +-100 cm square
584 algo->fRadThick[iSta].table.resize(1);
585 algo->fRadThick[iSta].table[0].resize(1);
586 algo->fRadThick[iSta].table[0][0] = algo->vStations[iSta].materialInfo.RadThick[0];
587 }
588 }
589
590 return kSUCCESS;
591}
592
593void CbmL1::Exec(Option_t* option) {}
594
595void CbmL1::Reconstruct()
596{
597 static int nevent = 0;
598 if (fVerbose > 1)
599 cout << endl << "CbmL1::Exec event " << ++nevent << " ..." << endl << endl;
600
601 // repack data
602 ReadEvent();
603
604 { // save data for standalone
605 if (fSTAPDataMode % 2 == 1) { // 1,3
606 WriteSTAPAlgoData();
607 WriteSTAPPerfData();
608 };
609 if (fSTAPDataMode >= 2) { // 2,3
610 ReadSTAPAlgoData();
611 ReadSTAPPerfData();
612 };
613 }
614
615 if (0) { // correct hits on MC // dbg
616 TRandom3 random;
617 vector<int> sF, sB, zP;
618 sF.clear();
619 sB.clear();
620 zP.clear();
621 for (unsigned int iH = 0; iH < algo->vStsHits.size(); ++iH) {
622 // for( unsigned int iH = algo->vStsHits.size() - 1; iH >= 0; --iH ) {
623 if (vStsHits[iH].mcPointIds.size() == 0)
624 continue;
625 L1StsHit& h = algo->vStsHits[iH];
626 const CbmL1MCPoint& mcp = vMCPoints[vStsHits[iH].mcPointIds[0]];
627 const int ista = algo->vSFlag[h.f] / 4;
628 const L1Station& sta = algo->vStations[ista];
629 if (std::find(sF.begin(), sF.end(), h.f) != sF.end()) { // separate strips
630 algo->vSFlag.push_back(algo->vSFlag[h.f]);
631 h.f = algo->vStsStrips.size();
632 algo->vStsStrips.push_back(L1Strip());
633 }
634 sF.push_back(h.f);
635 if (std::find(sB.begin(), sB.end(), h.b) != sB.end()) {
636 algo->vSFlagB.push_back(algo->vSFlagB[h.b]);
637 h.b = algo->vStsStripsB.size();
638 algo->vStsStripsB.push_back(L1Strip());
639 }
640 sB.push_back(h.b);
641 if (std::find(zP.begin(), zP.end(), h.iz) != zP.end()) { // TODO why do we need it??gives prob=0
642 h.iz = algo->vStsZPos.size();
643 algo->vStsZPos.push_back(0);
644 }
645 zP.push_back(h.iz);
646
647 const fscal idet = 1 / (sta.xInfo.sin_phi * sta.yInfo.sin_phi - sta.xInfo.cos_phi * sta.yInfo.cos_phi)[0];
648#if 1 // GAUSS
649 algo->vStsStrips[h.f] = idet * (+sta.yInfo.sin_phi[0] * mcp.x - sta.xInfo.cos_phi[0] * mcp.y)
650 + random.Gaus(0, sqrt(sta.frontInfo.sigma2)[0]);
651 algo->vStsStripsB[h.b] = idet * (-sta.yInfo.cos_phi[0] * mcp.x + sta.xInfo.sin_phi[0] * mcp.y)
652 + random.Gaus(0, sqrt(sta.backInfo.sigma2)[0]);
653#else // UNIFORM
654 algo->vStsStrips[h.f] =
655 idet * (+sta.yInfo.sin_phi[0] * mcp.x - sta.xInfo.cos_phi[0] * mcp.y)
656 + random.Uniform(-sqrt(sta.frontInfo.sigma2)[0] * 3.5, sqrt(sta.frontInfo.sigma2)[0] * 3.5);
657 algo->vStsStripsB[h.b] =
658 idet * (-sta.yInfo.cos_phi[0] * mcp.x + sta.xInfo.sin_phi[0] * mcp.y)
659 + random.Uniform(-sqrt(sta.backInfo.sigma2)[0] * 3.5, sqrt(sta.backInfo.sigma2)[0] * 3.5);
660#endif
661 algo->vStsZPos[h.iz] = mcp.z;
662 }
663 }
664
665 // input performance
666 if (fPerformance) {
667 InputPerformance();
668 }
669 // FieldApproxCheck();
670 // FieldIntegralCheck();
671
672 if (fVerbose > 1)
673 cout << "L1 Track finder..." << endl;
674 algo->CATrackFinder();
675 // IdealTrackFinder();
676 if (fVerbose > 1)
677 cout << "L1 Track finder ok" << endl;
678 algo->L1KFTrackFitter(fExtrapolateToTheEndOfSTS);
679 // algo->KFTrackFitter_simple();
680 if (fVerbose > 1)
681 cout << "L1 Track fitter ok" << endl;
682
683 // save recontstructed tracks
684 vRTracks.clear();
685 int start_hit = 0;
686 for (vector<L1Track>::iterator it = algo->vTracks.begin(); it != algo->vTracks.end(); it++) {
687 CbmL1Track t;
688 for (int i = 0; i < 6; i++)
689 t.T[i] = it->TFirst[i];
690 for (int i = 0; i < 15; i++)
691 t.C[i] = it->CFirst[i];
692 for (int i = 0; i < 6; i++)
693 t.TLast[i] = it->TLast[i];
694 for (int i = 0; i < 15; i++)
695 t.CLast[i] = it->CLast[i];
696 t.chi2 = it->chi2;
697 t.NDF = it->NDF;
698 // t.T[4] = it->Momentum;
699 for (int i = 0; i < it->NHits; i++) {
700 t.StsHits.push_back(algo->vRecoHits[start_hit++]);
701 }
702 t.mass = 0.1395679; // pion mass
703 t.is_electron = 0;
704
705 t.SetId(vRTracks.size());
706 vRTracks.push_back(t);
707 }
708
709 // Find Primary vertex, Ks, Lambdas,...
710 // static CbmL1ParticlesFinder PF;
711
712 // output performance
713 if (fPerformance) {
714 if (fVerbose > 1)
715 cout << "Performance..." << endl;
716 TrackMatch();
717 }
718
719 if (fFindParticlesMode) {
720 // Find Primary vertex, Ks, Lambdas,...
721 PF->FindParticles(vRTracks);
722 vRParticles = PF->GetParticles();
723 }
724
725 if (fPerformance) {
726 EfficienciesPerformance();
727 HistoPerformance();
728 TrackFitPerformance();
729
730 if (fFindParticlesMode) {
731 GetMCParticles();
732 FindReconstructableMCParticles();
733 MatchParticles();
734 PartEffPerformance();
735 PartHistoPerformance();
736 }
738 }
739 if (fVerbose > 1)
740 cout << "End of L1" << endl;
741
742 static bool ask = 0;
743 char symbol;
744 if (ask) {
745 std::cout << std::endl << "L1 run (any/r/q) > ";
746 do {
747 std::cin.get(symbol);
748 if (symbol == 'r')
749 ask = false;
750 if ((symbol == 'e') || (symbol == 'q'))
751 exit(0);
752 } while (symbol != '\n');
753 }
754}
755
756// ----- Finish CbmStsFitPerformanceTask task -----------------------------
757void CbmL1::Finish()
758{
759 TDirectory* curr = gDirectory;
760 TFile* currentFile = gFile;
761 // Open output file and write histograms
762 TFile* outfile = new TFile("L1_histo.root", "RECREATE");
763 outfile->cd();
764 writedir2current(histodir);
765 outfile->Close();
766 outfile->Delete();
767 gFile = currentFile;
768 gDirectory = curr;
769}
770
771void CbmL1::writedir2current(TObject* obj)
772{
773 if (!obj->IsFolder())
774 obj->Write();
775 else {
776 TDirectory* cur = gDirectory;
777 TDirectory* sub = cur->mkdir(obj->GetName());
778 sub->cd();
779 TList* listSub = (L1_DYNAMIC_CAST<TDirectory*>(obj))->GetList();
780 TIter it(listSub);
781 while (TObject* obj1 = it())
782 writedir2current(obj1);
783 cur->cd();
784 }
785}
786
788
789void CbmL1::IdealTrackFinder()
790{
791 algo->vTracks.clear();
792 algo->vRecoHits.clear();
793
794 for (vector<CbmL1MCTrack>::iterator i = vMCTracks.begin(); i != vMCTracks.end(); ++i) {
795 CbmL1MCTrack& MC = *i;
796
797 if (!MC.IsReconstructable())
798 continue;
799 if (!(MC.ID >= 0))
800 continue;
801
802 if (MC.StsHits.size() < 4)
803 continue;
804 L1Track algoTr;
805 algoTr.NHits = 0;
806 int lastStation = -1;
807 for (unsigned int iH = 0; iH < MC.StsHits.size(); iH++) {
808 const int hitI = MC.StsHits[iH];
809 const CbmL1StsHit& hit = vStsHits[hitI];
810 if (vMCPoints[hit.mcPointIds[0]].iStation <= lastStation) { // one hit per station
811 continue;
812 }
813 lastStation = vMCPoints[hit.mcPointIds[0]].iStation;
814 algo->vRecoHits.push_back(hitI);
815 algoTr.NHits++;
816 }
817 algoTr.Momentum = MC.p;
818
819 algo->vTracks.push_back(algoTr);
820 }
821
822}; // void CbmL1::IdealTrackFinder()
823
825
826void CbmL1::WriteSTAPGeoData(void* geo_, int size)
827{
828 fscal* geo = static_cast<fscal*>(geo_);
829 // write geo in file
830 TString fgeo_name = fSTAPDataDir + "geo_algo.txt";
831 ofstream fgeo(fgeo_name);
832 fgeo.setf(ios::scientific, ios::floatfield);
833 fgeo.precision(20);
834 for (int i = 0; i < size; i++) {
835 fgeo << geo[i] << endl;
836 };
837 fgeo.close();
838 cout << "-I- CbmL1: Geometry data has been written in " << fgeo_name << endl;
839} // void CbmL1::WriteSTAPGeoData(void* geo_, int size)
840
841void CbmL1::WriteSTAPAlgoData() // must be called after ReadEvent
842{
843 // write algo data in file
844 static int vNEvent = 1;
845 fstream fadata;
846
847 TString fadata_name = fSTAPDataDir + "data_algo.txt";
848 // if ( vNEvent <= maxNEvent ) {
849 if (1) {
850
851 if (vNEvent == 1)
852 fadata.open(fadata_name, fstream::out); // begin new file
853 else
854 fadata.open(fadata_name, fstream::out | fstream::app);
855
856 fadata << "Event:" << " ";
857 fadata << vNEvent << endl;
858 // write vStsStrips
859 int n = algo->vStsStrips.size(); // number of elements
860 fadata << n << endl;
861 for (int i = 0; i < n; i++) {
862 fadata << algo->vStsStrips[i] << endl;
863 };
864 if (fVerbose >= 4)
865 cout << "vStsStrips[" << n << "]" << " have been written." << endl;
866 // write vStsStripsB
867 n = algo->vStsStripsB.size();
868 fadata << n << endl;
869 for (int i = 0; i < n; i++) {
870 fadata << algo->vStsStripsB[i] << endl;
871 };
872 if (fVerbose >= 4)
873 cout << "vStsStripsB[" << n << "]" << " have been written." << endl;
874 // write vStsZPos
875 n = algo->vStsZPos.size();
876 fadata << n << endl;
877 for (int i = 0; i < n; i++) {
878 fadata << algo->vStsZPos[i] << endl;
879 };
880 if (fVerbose >= 4)
881 cout << "vStsZPos[" << n << "]" << " have been written." << endl;
882 // write vSFlag
883 n = algo->vSFlag.size();
884 fadata << n << endl;
885 unsigned char element;
886 for (int i = 0; i < n; i++) {
887 element = algo->vSFlag[i];
888 fadata << static_cast<int>(element) << endl;
889 };
890 if (fVerbose >= 4)
891 cout << "vSFlag[" << n << "]" << " have been written." << endl;
892 // write vSFlagB
893 n = algo->vSFlagB.size();
894 fadata << n << endl;
895 for (int i = 0; i < n; i++) {
896 element = algo->vSFlagB[i];
897 fadata << static_cast<int>(element) << endl;
898 };
899 if (fVerbose >= 4)
900 cout << "vSFlagB[" << n << "]" << " have been written." << endl;
901 // write vStsHits
902 n = algo->vStsHits.size();
903 fadata << n << endl;
904 for (int i = 0; i < n; i++) {
905 fadata << static_cast<int>(algo->vStsHits[i].f) << " ";
906 fadata << static_cast<int>(algo->vStsHits[i].b) << " ";
907 fadata << static_cast<int>(algo->vStsHits[i].iz) << endl;
908 };
909 if (fVerbose >= 4)
910 cout << "vStsHits[" << n << "]" << " have been written." << endl;
911 // write StsHitsStartIndex and StsHitsStopIndex
912 n = 20;
913 for (int i = 0; i < n; i++) {
914 if (algo->MaxNStations + 1 > i)
915 fadata << algo->StsHitsStartIndex[i] << endl;
916 else
917 fadata << 0 << endl;
918 };
919 for (int i = 0; i < n; i++) {
920 if (algo->MaxNStations + 1 > i)
921 fadata << algo->StsHitsStopIndex[i] << endl;
922 else
923 fadata << 0 << endl;
924 };
925
926 fadata.close();
927 }
928 cout << "-I- CbmL1: CATrackFinder data for event number " << vNEvent << " have been written in file " << fadata_name
929 << endl;
930 vNEvent++;
931} // void CbmL1::WriteSTAPAlgoData()
932
933void CbmL1::WriteSTAPPerfData() // must be called after ReadEvent
934{
935 fstream fpdata;
936 fpdata << setprecision(8);
937
938 static int vNEvent = 1;
939
940 TString fpdata_name = fSTAPDataDir + "data_perfo.txt";
941 // write data for performance in file
942 // if ( vNEvent <= maxNEvent ) {
943 if (1) {
944
945 if (vNEvent == 1)
946 fpdata.open(fpdata_name, fstream::out); // begin new file
947 else
948 fpdata.open(fpdata_name, fstream::out | fstream::app);
949
950 fpdata << "Event: ";
951 fpdata << vNEvent << endl;
952 // write vMCPoints
953 Int_t n = vMCPoints.size(); // number of elements
954 fpdata << n << endl;
955 for (int i = 0; i < n; i++) {
956 fpdata << vMCPoints[i].xIn << " ";
957 fpdata << vMCPoints[i].yIn << " ";
958 fpdata << vMCPoints[i].zIn << " ";
959 fpdata << vMCPoints[i].pxIn << " ";
960 fpdata << vMCPoints[i].pyIn << " ";
961 fpdata << vMCPoints[i].pzIn << " " << endl;
962 fpdata << vMCPoints[i].xOut << " ";
963 fpdata << vMCPoints[i].yOut << " ";
964 fpdata << vMCPoints[i].zOut << " ";
965 fpdata << vMCPoints[i].pxOut << " ";
966 fpdata << vMCPoints[i].pyOut << " ";
967 fpdata << vMCPoints[i].pzOut << " " << endl;
968
969 fpdata << vMCPoints[i].p << " ";
970 fpdata << vMCPoints[i].q << " ";
971 fpdata << vMCPoints[i].mass << " ";
972
973 fpdata << vMCPoints[i].pdg << " ";
974 fpdata << vMCPoints[i].ID << " ";
975 fpdata << vMCPoints[i].mother_ID << " ";
976 fpdata << vMCPoints[i].iStation << endl;
977
978 const int nhits = vMCPoints[i].hitIds.size();
979 fpdata << nhits << endl << " ";
980 for (int k = 0; k < nhits; k++) {
981 fpdata << vMCPoints[i].hitIds[k] << " ";
982 };
983 fpdata << endl;
984 };
985 if (fVerbose >= 4)
986 cout << "vMCPoints[" << n << "]" << " have been written." << endl;
987
988 // write vMCTracks . without Points
989 n = vMCTracks.size(); // number of elements
990 fpdata << n << endl;
991 for (int i = 0; i < n; i++) {
992 fpdata << vMCTracks[i].x << " ";
993 fpdata << vMCTracks[i].y << " ";
994 fpdata << vMCTracks[i].z << " ";
995 fpdata << vMCTracks[i].px << " ";
996 fpdata << vMCTracks[i].py << " ";
997 fpdata << vMCTracks[i].pz << " ";
998 fpdata << vMCTracks[i].p << " ";
999 fpdata << vMCTracks[i].q << " ";
1000 fpdata << vMCTracks[i].mass << " ";
1001
1002 fpdata << vMCTracks[i].pdg << " ";
1003 fpdata << vMCTracks[i].ID << " ";
1004 fpdata << vMCTracks[i].mother_ID << endl;
1005
1006 int nhits = vMCTracks[i].StsHits.size();
1007 fpdata << " " << nhits << endl << " ";
1008 for (int k = 0; k < nhits; k++) {
1009 fpdata << vMCTracks[i].StsHits[k] << " ";
1010 };
1011 fpdata << endl;
1012
1013 const int nPoints = vMCTracks[i].Points.size();
1014 fpdata << nPoints << endl << " ";
1015 for (int k = 0; k < nPoints; k++) {
1016 fpdata << vMCTracks[i].Points[k] << " ";
1017 };
1018 fpdata << endl;
1019
1020 fpdata << vMCTracks[i].nMCContStations << " ";
1021 fpdata << vMCTracks[i].nHitContStations << " ";
1022 fpdata << vMCTracks[i].maxNStaMC << " ";
1023 fpdata << vMCTracks[i].maxNSensorMC << " ";
1024 fpdata << vMCTracks[i].maxNStaHits << " ";
1025 fpdata << vMCTracks[i].nStations << endl;
1026 };
1027 if (fVerbose >= 4)
1028 cout << "vMCTracks[" << n << "]" << " have been written." << endl;
1029
1030 // write vHitMCRef
1031 n = vHitMCRef.size(); // number of elements
1032 fpdata << n << endl;
1033 for (int i = 0; i < n; i++) {
1034 fpdata << vHitMCRef[i] << endl;
1035 };
1036 if (fVerbose >= 4)
1037 cout << "vHitMCRef[" << n << "]" << " have been written." << endl;
1038
1039 // write vHitStore
1040 n = vHitStore.size(); // number of elements
1041 fpdata << n << endl;
1042 for (int i = 0; i < n; i++) {
1043 fpdata << vHitStore[i].ExtIndex << " ";
1044 fpdata << vHitStore[i].iStation << " ";
1045
1046 fpdata << vHitStore[i].x << " ";
1047 fpdata << vHitStore[i].y << endl;
1048 };
1049 if (fVerbose >= 4)
1050 cout << "vHitStore[" << n << "]" << " have been written." << endl;
1051
1052 // write vStsHits
1053 n = vStsHits.size(); // number of elements
1054 fpdata << n << endl;
1055 for (int i = 0; i < n; i++) {
1056 fpdata << vStsHits[i].hitId << " ";
1057 fpdata << vStsHits[i].extIndex << endl;
1058
1059 const int nPoints = vStsHits[i].mcPointIds.size();
1060 fpdata << nPoints << endl << " ";
1061 for (int k = 0; k < nPoints; k++) {
1062 fpdata << vStsHits[i].mcPointIds[k] << " ";
1063 };
1064 fpdata << endl;
1065 };
1066 if (fVerbose >= 4)
1067 cout << "vStsHits[" << n << "]" << " have been written." << endl;
1068
1069 fpdata.close();
1070 }
1071 cout << "-I- CbmL1: Data for performance of event number " << vNEvent << " have been written in file "
1072 << fpdata_name << endl;
1073 vNEvent++;
1074} // void CbmL1::WriteSTAPPerfData()
1075
1076istream& CbmL1::eatwhite(istream& is) // skip spaces
1077{
1078 char c;
1079 while (is.get(c)) {
1080 if (isspace(c) == 0) {
1081 is.putback(c);
1082 break;
1083 }
1084 }
1085 return is;
1086}
1087
1088void CbmL1::ReadSTAPGeoData(void* geo_, int& size)
1089{
1090 fscal* geo = static_cast<fscal*>(geo_);
1091 TString fgeo_name = fSTAPDataDir + "geo_algo.txt";
1092 ifstream fgeo(fgeo_name);
1093
1094 cout << "-I- CbmL1: Read geometry from file " << fgeo_name << endl;
1095 int i;
1096 for (i = 0; !fgeo.eof(); i++) {
1097 fscal tmp;
1098 fgeo >> tmp >> eatwhite;
1099 geo[i] = tmp;
1100 };
1101 size = i;
1102 fgeo.close();
1103} // void CbmL1::ReadSTAPGeoData(void* geo_, int &size)
1104
1105void CbmL1::ReadSTAPAlgoData()
1106{
1107 static int nEvent = 1;
1108 static fstream fadata;
1109 TString fadata_name = fSTAPDataDir + "data_algo.txt";
1110 // if (nEvent <= maxNEvent){
1111 if (1) {
1112 if (nEvent == 1) {
1113 fadata.open(fadata_name, fstream::in);
1114 };
1115
1116 algo->vStsHits.clear();
1117 algo->vStsStrips.clear();
1118 algo->vStsStripsB.clear();
1119 algo->vStsZPos.clear();
1120 algo->vSFlag.clear();
1121 algo->vSFlagB.clear();
1122
1123 // check correct position in file
1124 char s[] = "Event: ";
1125 int nEv;
1126 fadata >> s;
1127 fadata >> nEv;
1128 if (nEv != nEvent)
1129 cout << "-E- CbmL1: Can't read event number " << nEvent << " from file " << fadata_name << endl;
1130
1131 int n; // number of elements
1132 // read algo->vStsStrips
1133 fadata >> n;
1134 for (int i = 0; i < n; i++) {
1135 fscal element;
1136 fadata >> element;
1137 algo->vStsStrips.push_back(element);
1138 };
1139 if (fVerbose >= 4)
1140 cout << "vStsStrips[" << n << "]" << " have been read." << endl;
1141 // read algo->vStsStripsB
1142 fadata >> n;
1143 for (int i = 0; i < n; i++) {
1144 fscal element;
1145 fadata >> element;
1146 algo->vStsStripsB.push_back(element);
1147 };
1148 if (fVerbose >= 4)
1149 cout << "vStsStripsB[" << n << "]" << " have been read." << endl;
1150 // read algo->vStsZPos
1151 fadata >> n;
1152 for (int i = 0; i < n; i++) {
1153 fscal element;
1154 fadata >> element;
1155 algo->vStsZPos.push_back(element);
1156 };
1157 if (fVerbose >= 4)
1158 cout << "vStsZPos[" << n << "]" << " have been read." << endl;
1159 // read algo->vSFlag
1160 fadata >> n;
1161 for (int i = 0; i < n; i++) {
1162 int element;
1163 fadata >> element;
1164 algo->vSFlag.push_back(static_cast<unsigned char>(element));
1165 };
1166 if (fVerbose >= 4)
1167 cout << "vSFlag[" << n << "]" << " have been read." << endl;
1168 // read algo->vSFlagB
1169 fadata >> n;
1170 for (int i = 0; i < n; i++) {
1171 int element;
1172 fadata >> element;
1173 algo->vSFlagB.push_back(static_cast<unsigned char>(element));
1174 };
1175 if (fVerbose >= 4)
1176 cout << "vSFlagB[" << n << "]" << " have been read." << endl;
1177 // read algo->vStsHits
1178 fadata >> n;
1179 int element_f; // for convert
1180 int element_b;
1181 int element_iz;
1182 for (int i = 0; i < n; i++) {
1183 L1StsHit element;
1184 fadata >> element_f >> element_b >> element_iz;
1185 element.f = static_cast<THitI>(element_f);
1186 element.b = static_cast<THitI>(element_b);
1187 element.iz = static_cast<TZPosI>(element_iz);
1188 algo->vStsHits.push_back(element);
1189 };
1190 if (fVerbose >= 4)
1191 cout << "vStsHits[" << n << "]" << " have been read." << endl;
1192 // read StsHitsStartIndex and StsHitsStopIndex
1193 n = 20;
1194 for (int i = 0; i < n; i++) {
1195 int tmp;
1196 fadata >> tmp;
1197 if (algo->MaxNStations + 1 > i)
1198 algo->StsHitsStartIndex[i] = tmp;
1199 };
1200 for (int i = 0; i < n; i++) {
1201 int tmp;
1202 fadata >> tmp;
1203 if (algo->MaxNStations + 1 > i)
1204 algo->StsHitsStopIndex[i] = tmp;
1205 };
1206
1207 cout << "-I- CbmL1: CATrackFinder data for event " << nEvent << " has been read from file " << fadata_name
1208 << " successfully." << endl;
1209 // if (nEvent == maxNEvent) fadata.close();
1210 }
1211 nEvent++;
1212} // void CbmL1::ReadSTAPAlgoData()
1213
1214void CbmL1::ReadSTAPPerfData()
1215{
1216 static int nEvent = 1;
1217 static fstream fpdata;
1218 TString fpdata_name = fSTAPDataDir + "data_perfo.txt";
1219 // if (nEvent <= maxNEvent){
1220 if (1) {
1221 if (nEvent == 1) {
1222 fpdata.open(fpdata_name, fstream::in);
1223 };
1224
1225 vMCPoints.clear();
1226 vMCTracks.clear();
1227 vHitMCRef.clear();
1228 vHitStore.clear();
1229 vStsHits.clear();
1230
1231 // check if it is right position in file
1232 char s[] = "EVENT: "; // buffer
1233 int nEv = 0; // event number
1234 fpdata >> s;
1235 fpdata >> nEv;
1236
1237 if (nEv != nEvent)
1238 cout << "-E- CbmL1: Performance: can't read event number " << nEvent << " from file " << "data_perfo.txt"
1239 << endl;
1240 // vMCPoints
1241 int n; // number of elements
1242 fpdata >> n;
1243 for (int i = 0; i < n; i++) {
1244 CbmL1MCPoint element;
1245
1246 fpdata >> element.xIn;
1247 fpdata >> element.yIn;
1248 fpdata >> element.zIn;
1249 fpdata >> element.pxIn;
1250 fpdata >> element.pyIn;
1251 fpdata >> element.pzIn;
1252
1253 fpdata >> element.xOut;
1254 fpdata >> element.yOut;
1255 fpdata >> element.zOut;
1256 fpdata >> element.pxOut;
1257 fpdata >> element.pyOut;
1258 fpdata >> element.pzOut;
1259
1260 fpdata >> element.p;
1261 fpdata >> element.q;
1262 fpdata >> element.mass;
1263
1264 fpdata >> element.pdg;
1265 fpdata >> element.ID;
1266 fpdata >> element.mother_ID;
1267 fpdata >> element.iStation;
1268
1269 int nhits;
1270 fpdata >> nhits;
1271 for (int k = 0; k < nhits; k++) {
1272 int helement;
1273 fpdata >> helement;
1274 element.hitIds.push_back(helement);
1275 };
1276
1277 vMCPoints.push_back(element);
1278 };
1279 if (fVerbose >= 4)
1280 cout << "vMCPoints[" << n << "]" << " have been read." << endl;
1281
1282 // vMCTracks . without Points
1283 fpdata >> n;
1284 for (int i = 0; i < n; i++) {
1285 CbmL1MCTrack element;
1286
1287 fpdata >> element.x;
1288 fpdata >> element.y;
1289 fpdata >> element.z;
1290 fpdata >> element.px;
1291 fpdata >> element.py;
1292 fpdata >> element.pz;
1293 fpdata >> element.p;
1294 fpdata >> element.q;
1295 fpdata >> element.mass;
1296
1297 fpdata >> element.pdg;
1298 fpdata >> element.ID;
1299 fpdata >> element.mother_ID;
1300
1301 int nhits;
1302 fpdata >> nhits;
1303 for (int k = 0; k < nhits; k++) {
1304 int helement;
1305 fpdata >> helement;
1306 element.StsHits.push_back(helement);
1307 };
1308 fpdata >> nhits;
1309 for (int k = 0; k < nhits; k++) {
1310 int helement;
1311 fpdata >> helement;
1312 element.Points.push_back(helement);
1313 };
1314
1315 fpdata >> element.nMCContStations;
1316 fpdata >> element.nHitContStations;
1317 fpdata >> element.maxNStaMC;
1318 fpdata >> element.maxNSensorMC;
1319 fpdata >> element.maxNStaHits;
1320 fpdata >> element.nStations;
1321
1322 element.CalculateIsReconstructable();
1323 vMCTracks.push_back(element);
1324 };
1325 if (fVerbose >= 4)
1326 cout << "vMCTracks[" << n << "]" << " have been read." << endl;
1327
1328 // vHitMCRef
1329 fpdata >> n;
1330 for (int i = 0; i < n; i++) {
1331 int element;
1332 fpdata >> element;
1333 vHitMCRef.push_back(element);
1334 };
1335 if (fVerbose >= 4)
1336 cout << "vHitMCRef[" << n << "]" << " have been read." << endl;
1337
1338 // vHitStore
1339 fpdata >> n;
1340 for (int i = 0; i < n; i++) {
1341 CbmL1HitStore element;
1342 fpdata >> element.ExtIndex;
1343 fpdata >> element.iStation;
1344
1345 fpdata >> element.x;
1346 fpdata >> element.y;
1347 vHitStore.push_back(element);
1348 };
1349 if (fVerbose >= 4)
1350 cout << "vHitStore[" << n << "]" << " have been read." << endl;
1351
1352 // vStsHits
1353 fpdata >> n;
1354 for (int i = 0; i < n; i++) {
1355 CbmL1StsHit element;
1356 fpdata >> element.hitId;
1357 fpdata >> element.extIndex;
1358
1359 int nPoints;
1360 fpdata >> nPoints;
1361 for (int k = 0; k < nPoints; k++) {
1362 int id;
1363 fpdata >> id;
1364 element.mcPointIds.push_back(id);
1365 };
1366 vStsHits.push_back(element);
1367 };
1368 if (fVerbose >= 4)
1369 cout << "vStsHits[" << n << "]" << " have been read." << endl;
1370
1371 // if (nEvent == maxNEvent) { // file open on begin of all work class and close at end
1372 // fpdata.close();
1373 // cout << " -I- Performance: data read from file " << "data_perfo.txt" << " successfully"<< endl;
1374 // }
1375 cout << "-I- CbmL1: L1Performance data for event " << nEvent << " has been read from file " << fpdata_name
1376 << " successfully." << endl;
1377
1378 } // if (nEvent <= maxNEvent)
1379 nEvent++;
1380} // void CbmL1::ReadSTAPPerfData()
1381
1382void CbmL1::WriteSIMDKFData()
1383{
1384 static bool first = 1;
1385
1387 if (first) {
1388 FairField* dMF = CbmKF::Instance()->GetMagneticField();
1389
1390 fstream FileGeo;
1391 FileGeo.open("geo.dat", ios::out);
1392
1393 fstream FieldCheck;
1394 FieldCheck.open("field.dat", ios::out);
1395
1396 Double_t bfg[3], rfg[3];
1397
1398 rfg[0] = 0.;
1399 rfg[1] = 0.;
1400 rfg[2] = 0.;
1401 dMF->GetFieldValue(rfg, bfg);
1402 FileGeo << rfg[2] << " " << bfg[0] << " " << bfg[1] << " " << bfg[2] << " " << endl;
1403
1404 rfg[0] = 0.;
1405 rfg[1] = 0.;
1406 rfg[2] = 2.5;
1407 dMF->GetFieldValue(rfg, bfg);
1408 FileGeo << rfg[2] << " " << bfg[0] << " " << bfg[1] << " " << bfg[2] << " " << endl;
1409
1410 rfg[0] = 0.;
1411 rfg[1] = 0.;
1412 rfg[2] = 5.0;
1413 dMF->GetFieldValue(rfg, bfg);
1414 FileGeo << rfg[2] << " " << bfg[0] << " " << bfg[1] << " " << bfg[2] << " " << endl << endl;
1415 FileGeo << NStation << endl;
1416
1417 const int M = 5; // polinom order
1418 const int N = (M + 1) * (M + 2) / 2;
1419
1420 for (Int_t ist = 0; ist < NStation; ist++) {
1421 fscal f_phi, f_sigma, b_phi, b_sigma;
1422
1423 double C[3][N];
1424 double z = 0;
1425 double Xmax, Ymax;
1426 if (ist < NMvdStations) {
1427 CbmKFTube& t = CbmKF::Instance()->vMvdMaterial[ist];
1428 f_phi = 0;
1429 f_sigma = 5.e-4;
1430 b_phi = 3.14159265358 / 2.;
1431 b_sigma = 5.e-4;
1432 z = t.z;
1433 Xmax = Ymax = t.R;
1434 } else {
1435 // AZ CbmStsStation *st = StsDigi.GetStation(ist - NMvdStations);
1436 CbmStsStation* st = StsDigi->GetStation(ist - NMvdStations);
1437 CbmStsSector* sector = st->GetSector(0);
1438 f_phi = sector->GetRotation();
1439 b_phi = sector->GetRotation();
1440 if (sector->GetType() == 1) {
1441 b_phi += 3.14159265358 / 2.;
1442 b_sigma = sector->GetDy() / sqrt(12.);
1443 f_sigma = b_sigma; // CHECKME
1444 } else {
1445 f_phi += sector->GetStereoF();
1446 b_phi += sector->GetStereoB();
1447 f_sigma = sector->GetDx() / TMath::Sqrt(12);
1448 b_sigma = f_sigma;
1449 }
1450 // AZ f_sigma *= cos(f_phi); // TODO: think about this
1451 // AZ b_sigma *= cos(b_phi);
1452 z = st->GetZ();
1453
1454 Xmax = -100;
1455 Ymax = -100;
1456
1457 double x, y;
1458 for (int isec = 0; isec < st->GetNSectors(); isec++) {
1459 CbmStsSector* sect = L1_DYNAMIC_CAST<CbmStsSector*>(st->GetSector(isec));
1460 for (int isen = 0; isen < sect->GetNSensors(); isen++) {
1461 x = sect->GetSensor(isen)->GetX0() + sect->GetSensor(isen)->GetLx() / 2.;
1462 y = sect->GetSensor(isen)->GetY0() + sect->GetSensor(isen)->GetLy() / 2.;
1463 if (x > Xmax)
1464 Xmax = x;
1465 if (y > Ymax)
1466 Ymax = y;
1467 }
1468 }
1469 }
1470
1471 double dx = 1.; // step for the field approximation
1472 double dy = 1.;
1473
1474 if (dx > Xmax / N / 2)
1475 dx = Xmax / N / 4.;
1476 if (dy > Ymax / N / 2)
1477 dy = Ymax / N / 4.;
1478
1479 for (int i = 0; i < 3; i++)
1480 for (int k = 0; k < N; k++)
1481 C[i][k] = 0;
1482 TMatrixD A(N, N);
1483 TVectorD b0(N), b1(N), b2(N);
1484 for (int i = 0; i < N; i++) {
1485 for (int j = 0; j < N; j++)
1486 A(i, j) = 0.;
1487 b0(i) = b1(i) = b2(i) = 0.;
1488 }
1489 for (double x = -Xmax; x <= Xmax; x += dx)
1490 for (double y = -Ymax; y <= Ymax; y += dy) {
1491 double r = sqrt(fabs(x * x / Xmax / Xmax + y / Ymax * y / Ymax));
1492 if (r > 1.)
1493 continue;
1494 Double_t w = 1. / (r * r + 1);
1495 Double_t p[3] = {x, y, z};
1496 Double_t B[3] = {0., 0., 0.};
1497 CbmKF::Instance()->GetMagneticField()->GetFieldValue(p, B);
1498 TVectorD m(N);
1499 m(0) = 1;
1500 for (int i = 1; i <= M; i++) {
1501 int k = (i - 1) * (i) / 2;
1502 int l = i * (i + 1) / 2;
1503 for (int j = 0; j < i; j++)
1504 m(l + j) = x * m(k + j);
1505 m(l + i) = y * m(k + i - 1);
1506 }
1507
1508 TVectorD mt = m;
1509 for (int i = 0; i < N; i++) {
1510 for (int j = 0; j < N; j++)
1511 A(i, j) += w * m(i) * m(j);
1512 b0(i) += w * B[0] * m(i);
1513 b1(i) += w * B[1] * m(i);
1514 b2(i) += w * B[2] * m(i);
1515 }
1516 }
1517 double det;
1518 A.Invert(&det);
1519 TVectorD c0 = A * b0, c1 = A * b1, c2 = A * b2;
1520 for (int i = 0; i < N; i++) {
1521 C[0][i] = c0(i);
1522 C[1][i] = c1(i);
1523 C[2][i] = c2(i);
1524 }
1525
1526 double c_f = cos(f_phi);
1527 double s_f = sin(f_phi);
1528 double c_b = cos(b_phi);
1529 double s_b = sin(b_phi);
1530
1531 double det_m = c_f * s_b - s_f * c_b;
1532 det_m *= det_m;
1533 // double C00 = ( s_b*s_b*f_sigma*f_sigma + s_f*s_f*b_sigma*b_sigma )/det_m;
1534 // double C10 =-( s_b*c_b*f_sigma*f_sigma + s_f*c_f*b_sigma*b_sigma )/det_m;
1535 // double C11 = ( c_b*c_b*f_sigma*f_sigma + c_f*c_f*b_sigma*b_sigma )/det_m;
1536
1537 // float delta_x = sqrt(C00);
1538 // float delta_y = sqrt(C11);
1539 FileGeo << " " << ist << " ";
1540 if (ist < NMvdStations) {
1541 CbmKFTube& t = CbmKF::Instance()->vMvdMaterial[ist];
1542 FileGeo << t.z << " ";
1543 FileGeo << t.dz << " ";
1544 FileGeo << t.RadLength << " ";
1545 } else if (ist < (NStsStations + NMvdStations)) {
1546 // AZ CbmStsStation *st = StsDigi.GetStation(ist - NMvdStations);
1547 CbmStsStation* st = StsDigi->GetStation(ist - NMvdStations);
1548 FileGeo << st->GetZ() << " ";
1549 FileGeo << st->GetD() << " ";
1550 FileGeo << st->GetRadLength() << " ";
1551 }
1552 FileGeo << f_sigma << " ";
1553 FileGeo << b_sigma << " ";
1554 FileGeo << f_phi << " ";
1555 FileGeo << b_phi << endl;
1556 FileGeo << " " << N << endl;
1557 FileGeo << " ";
1558 for (int ik = 0; ik < N; ik++)
1559 FileGeo << C[0][ik] << " ";
1560 FileGeo << endl;
1561 FileGeo << " ";
1562 for (int ik = 0; ik < N; ik++)
1563 FileGeo << C[1][ik] << " ";
1564 FileGeo << endl;
1565 FileGeo << " ";
1566 for (int ik = 0; ik < N; ik++)
1567 FileGeo << C[2][ik] << " ";
1568 FileGeo << endl;
1569 }
1570 FileGeo.close();
1571 }
1572
1574
1575 static int TrNumber = 0;
1576 fstream Tracks, McTracksCentr, McTracksIn, McTracksOut;
1577 if (first) {
1578 Tracks.open("tracks.dat", fstream::out);
1579 McTracksCentr.open("mctrackscentr.dat", fstream::out);
1580 McTracksIn.open("mctracksin.dat", fstream::out);
1581 McTracksOut.open("mctracksout.dat", fstream::out);
1582 first = 0;
1583 } else {
1584 Tracks.open("tracks.dat", fstream::out | fstream::app);
1585 McTracksCentr.open("mctrackscentr.dat", fstream::out | fstream::app);
1586 McTracksIn.open("mctracksin.dat", fstream::out | fstream::app);
1587 McTracksOut.open("mctracksout.dat", fstream::out | fstream::app);
1588 }
1589
1590 for (vector<CbmL1Track>::iterator RecTrack = vRTracks.begin(); RecTrack != vRTracks.end(); ++RecTrack) {
1591 if (RecTrack->IsGhost())
1592 continue;
1593
1594 CbmL1MCTrack* MCTrack = RecTrack->GetMCTrack();
1595 if (!(MCTrack->IsPrimary()))
1596 continue;
1597
1598 int NHits = (RecTrack->StsHits).size();
1599 float x[20], y[20], z[20];
1600 int st[20];
1601 int jHit = 0;
1602 for (int iHit = 0; iHit < NHits; iHit++) {
1603 CbmL1HitStore& h = vHitStore[RecTrack->StsHits[iHit]];
1604 st[jHit] = h.iStation;
1605 if (h.ExtIndex < 0) {
1606 CbmMvdHit* MvdH = (CbmMvdHit*)listMvdHits->At(-h.ExtIndex - 1);
1607 x[jHit] = MvdH->GetX();
1608 y[jHit] = MvdH->GetY();
1609 z[jHit] = MvdH->GetZ();
1610 jHit++;
1611 } else {
1612 CbmStsHit* StsH = (CbmStsHit*)listStsHits->At(h.ExtIndex);
1613 x[jHit] = StsH->GetX();
1614 y[jHit] = StsH->GetY();
1615 z[jHit] = StsH->GetZ();
1616 jHit++;
1617 }
1618 }
1619
1620 Tracks << TrNumber << " " << MCTrack->x << " " << MCTrack->y << " " << MCTrack->z << " " << MCTrack->px << " "
1621 << MCTrack->py << " " << MCTrack->pz << " " << MCTrack->q << " " << NHits << endl;
1622
1623 float AngleXAxis = 0, AngleYAxis = 0;
1624 for (int i = 0; i < NHits; i++)
1625 Tracks << " " << st[i] << " " << x[i] << " " << y[i] << " " << z[i] << " " << AngleXAxis << " "
1626 << AngleYAxis << endl;
1627 Tracks << endl;
1628
1629 int NMCPoints = (MCTrack->Points).size();
1630
1631 McTracksIn << TrNumber << " " << MCTrack->x << " " << MCTrack->y << " " << MCTrack->z << " " << MCTrack->px
1632 << " " << MCTrack->py << " " << MCTrack->pz << " " << MCTrack->q << " " << NMCPoints << endl;
1633 McTracksOut << TrNumber << " " << MCTrack->x << " " << MCTrack->y << " " << MCTrack->z << " " << MCTrack->px
1634 << " " << MCTrack->py << " " << MCTrack->pz << " " << MCTrack->q << " " << NMCPoints << endl;
1635 McTracksCentr << TrNumber << " " << MCTrack->x << " " << MCTrack->y << " " << MCTrack->z << " " << MCTrack->px
1636 << " " << MCTrack->py << " " << MCTrack->pz << " " << MCTrack->q << " " << NMCPoints << endl;
1637
1638 for (int iPoint = 0; iPoint < NMCPoints; iPoint++) {
1639 CbmL1MCPoint& MCPoint = vMCPoints[MCTrack->Points[iPoint]];
1640 McTracksIn << " " << MCPoint.iStation << " " << MCPoint.xIn << " " << MCPoint.yIn << " " << MCPoint.zIn
1641 << " " << MCPoint.pxIn << " " << MCPoint.pyIn << " " << MCPoint.pzIn << endl;
1642 McTracksOut << " " << MCPoint.iStation << " " << MCPoint.xOut << " " << MCPoint.yOut << " "
1643 << MCPoint.zOut << " " << MCPoint.pxOut << " " << MCPoint.pyOut << " " << MCPoint.pzOut << endl;
1644 McTracksCentr << " " << MCPoint.iStation << " " << MCPoint.x << " " << MCPoint.y << " " << MCPoint.z
1645 << " " << MCPoint.px << " " << MCPoint.py << " " << MCPoint.pz << endl;
1646 }
1647 McTracksIn << endl;
1648 McTracksOut << endl;
1649 McTracksCentr << endl;
1650
1651 TrNumber++;
1652 }
1653}
W
vector< Double_t > W(vector< Double_t > dist, Double_t sig)
Definition BmnMath.cxx:878
NBins
const int NBins
Definition BmnMonProf.cxx:47
TZPosI
unsigned short int TZPosI
Definition L1StsHit.h:7
THitI
unsigned int THitI
Definition L1StsHit.h:6
fscal
float fscal
Definition P4_F32vec4.h:232
sin
friend F32vec4 sin(const F32vec4 &a)
Definition P4_F32vec4.h:124
sqrt
friend F32vec4 sqrt(const F32vec4 &a)
Definition P4_F32vec4.h:34
fabs
friend F32vec4 fabs(const F32vec4 &a)
Definition P4_F32vec4.h:52
cos
friend F32vec4 cos(const F32vec4 &a)
Definition P4_F32vec4.h:125
_fvecalignment
#define _fvecalignment
Definition P4_F32vec4.h:236
i
int i
Definition P4_F32vec4.h:22
m
__m128 m
Definition P4_F32vec4.h:27
BmnTask
BmnTask.
Definition BmnTask.h:13
CbmKFMaterial::RadLength
Double_t RadLength
Definition CbmKFMaterial.h:34
CbmKFTube::R
Double_t R
Definition CbmKFMaterial.h:86
CbmKFTube::r
Double_t r
Definition CbmKFMaterial.h:86
CbmKFTube::dz
Double_t dz
Definition CbmKFMaterial.h:85
CbmKFTube::z
Double_t z
Definition CbmKFMaterial.h:85
CbmKF::Instance
static CbmKF * Instance()
Definition CbmKF.h:35
CbmKF::GetMagneticField
BmnNewFieldMap * GetMagneticField()
Definition CbmKF.h:73
CbmKF::vMvdMaterial
std::vector< CbmKFTube > vMvdMaterial
Definition CbmKF.h:60
CbmL1HitStore::y
double y
Definition CbmL1.h:45
CbmL1HitStore::ExtIndex
int ExtIndex
Definition CbmL1.h:43
CbmL1HitStore::x
double x
Definition CbmL1.h:45
CbmL1HitStore::iStation
int iStation
Definition CbmL1.h:44
CbmL1MCTrack::IsReconstructable
bool IsReconstructable() const
Definition CbmL1MCTrack.h:50
CbmL1MCTrack::StsHits
vector< int > StsHits
Definition CbmL1MCTrack.h:35
CbmL1MCTrack::IsPrimary
bool IsPrimary() const
Definition CbmL1MCTrack.h:49
CbmL1MCTrack::Points
vector< int > Points
Definition CbmL1MCTrack.h:34
CbmL1ParticlesFinder::GetParticles
vector< CbmKFParticle > & GetParticles()
Definition CbmL1ParticlesFinder.h:35
CbmL1ParticlesFinder::FindParticles
void FindParticles(vector< CbmL1Track > &vRTracks)
Definition CbmL1ParticlesFinder.cxx:49
CbmL1StsHit::mcPointIds
vector< int > mcPointIds
Definition CbmL1StsHit.h:17
CbmL1Track::TLast
double TLast[6]
Definition CbmL1Track.h:53
CbmL1Track::StsHits
vector< int > StsHits
Definition CbmL1Track.h:54
CbmL1Track::CLast
double CLast[15]
Definition CbmL1Track.h:53
CbmL1
Definition CbmL1.h:49
CbmL1::algo
L1Algo * algo
Definition CbmL1.h:55
CbmL1::vHitStore
vector< CbmL1HitStore > vHitStore
Definition CbmL1.h:87
CbmL1::Init
virtual InitStatus Init()
Definition CbmL1.cxx:202
CbmL1::CbmL1
CbmL1()
Definition CbmL1.cxx:51
CbmL1::vRTracks
vector< CbmL1Track > vRTracks
Definition CbmL1.h:58
CbmL1::Finish
void Finish()
Definition CbmL1.cxx:757
CbmL1::Exec
void Exec(Option_t *option)
Definition CbmL1.cxx:593
CbmL1::PF
CbmL1ParticlesFinder * PF
for access to L1 Algorithm from L1::Instance
Definition CbmL1.h:56
CbmL1::Reconstruct
void Reconstruct()
Definition CbmL1.cxx:595
CbmL1::SetParContainers
void SetParContainers()
Definition CbmL1.cxx:182
CbmL1::~CbmL1
~CbmL1()
Definition CbmL1.cxx:175
CbmL1::ReInit
virtual InitStatus ReInit()
Definition CbmL1.cxx:194
CbmStsDigiScheme::Instance
static CbmStsDigiScheme * Instance(int version=1)
Definition CbmStsDigiScheme.cxx:65
CbmStsDigiScheme::GetStation
CbmStsStation * GetStation(Int_t iStation)
Definition CbmStsDigiScheme.cxx:184
CbmStsSector::GetType
Int_t GetType() const
Definition CbmStsSector.h:93
CbmStsSector::GetDx
Double_t GetDx() const
Definition CbmStsSector.h:95
CbmStsSector::GetDy
Double_t GetDy() const
Definition CbmStsSector.h:96
CbmStsSector::GetStereoF
Double_t GetStereoF() const
Definition CbmStsSector.h:97
CbmStsSector::GetRotation
Double_t GetRotation() const
Definition CbmStsSector.h:94
CbmStsSector::GetSensor
CbmStsSensor * GetSensor(Int_t iSensor)
Definition CbmStsSector.h:85
CbmStsSector::GetStereoB
Double_t GetStereoB() const
Definition CbmStsSector.h:98
CbmStsSector::GetNSensors
Int_t GetNSensors() const
Definition CbmStsSector.h:82
CbmStsSensor::GetX0
Double_t GetX0() const
Definition CbmStsSensor.h:94
CbmStsSensor::GetLy
Double_t GetLy() const
Definition CbmStsSensor.h:99
CbmStsSensor::GetY0
Double_t GetY0() const
Definition CbmStsSensor.h:95
CbmStsSensor::GetLx
Double_t GetLx() const
Definition CbmStsSensor.h:98
CbmStsStation::GetRadLength
Double_t GetRadLength() const
Definition CbmStsStation.h:59
CbmStsStation::GetNSectors
Int_t GetNSectors() const
Definition CbmStsStation.h:62
CbmStsStation::GetD
Double_t GetD() const
Definition CbmStsStation.h:58
CbmStsStation::GetRmin
Double_t GetRmin() const
Definition CbmStsStation.h:60
CbmStsStation::GetZ
Double_t GetZ(Int_t it=0)
Definition CbmStsStation.cxx:152
CbmStsStation::GetRmax
Double_t GetRmax() const
Definition CbmStsStation.h:61
CbmStsStation::GetSector
CbmStsSector * GetSector(Int_t iSector)
Definition CbmStsStation.h:65
L1Algo::vSFlag
vector< unsigned char > vSFlag
Definition L1Algo.h:134
L1Algo::NStations
int NStations
Definition L1Algo.h:123
L1Algo::StsHitsStopIndex
THitI StsHitsStopIndex[MaxNStations+1]
Definition L1Algo.h:136
L1Algo::vStsStrips
vector< L1Strip > vStsStrips
Definition L1Algo.h:129
L1Algo::MaxNStations
@ MaxNStations
Definition L1Algo.h:122
L1Algo::vTracks
vector< L1Track > vTracks
--— Output data --—
Definition L1Algo.h:151
L1Algo::vSFlagB
vector< unsigned char > vSFlagB
Definition L1Algo.h:135
L1Algo::StsHitsStartIndex
THitI StsHitsStartIndex[MaxNStations+1]
Definition L1Algo.h:136
L1Algo::vRecoHits
vector< THitI > vRecoHits
Definition L1Algo.h:152
L1Algo::NMvdStations
int NMvdStations
Definition L1Algo.h:124
L1Algo::vStsHits
vector< L1StsHit > vStsHits
Definition L1Algo.h:132
L1Algo::vStsStripsB
vector< L1Strip > vStsStripsB
Definition L1Algo.h:130
L1Algo::CATrackFinder
void CATrackFinder()
The main procedure - find tracks.
Definition L1CATrackFinder.cxx:1332
L1Algo::vStsZPos
vector< fscal > vStsZPos
Definition L1Algo.h:131
L1Algo::L1KFTrackFitter
void L1KFTrackFitter(bool extrapolateToTheEndOfSTS=false)
Definition L1TrackFitter.cxx:329
L1Algo::Init
void Init(const fscal geo[])
Definition L1Algo.cxx:3
L1Algo::fRadThick
vector< L1Material > fRadThick
Definition L1Algo.h:127
L1Station::backInfo
L1UMeasurementInfo backInfo
Definition L1Station.h:22
L1Station::xInfo
L1UMeasurementInfo xInfo
Definition L1Station.h:22
L1Station::yInfo
L1UMeasurementInfo yInfo
Definition L1Station.h:22
L1Station::frontInfo
L1UMeasurementInfo frontInfo
Definition L1Station.h:22
L1StsHit::iz
TZPosI iz
Definition L1StsHit.h:13
L1StsHit::f
TStripI f
Definition L1StsHit.h:12
L1StsHit::b
TStripI b
Definition L1StsHit.h:12
L1Track::Momentum
float Momentum
Definition L1Track.h:24
L1Track::NHits
unsigned char NHits
Definition L1Track.h:23
check_success.s
s
Definition check_success.py:47
std
STL namespace.
CbmL1MCPoint::hitIds
vector< int > hitIds
Definition CbmL1MCPoint.h:52
CbmL1TrackPar::C
double C[15]
Definition CbmL1TrackPar.h:18