CbmRoot
riplet/Lx.cxx
Go to the documentation of this file.
1 #include "Lx.h"
2 #include "CbmKFParticle.h"
3 #include "CbmKFTrack.h"
4 #include "CbmMCTrack.h"
5 #include "CbmMuchCluster.h"
6 #include "CbmMuchDigiMatch.h"
7 #include "CbmMuchGeoScheme.h"
8 #include "CbmMuchPoint.h"
9 #include "CbmMuchTrack.h"
10 #include "CbmStsAddress.h"
11 #include "CbmStsPoint.h"
12 #include "LxDraw.h"
13 #include "TDatabasePDG.h"
14 #include "TGeoManager.h"
15 #include <iostream>
16 #include <sys/time.h>
17 
18 ClassImp(LxFinderTriplet)
19 
20  using namespace std;
21 
22 LxFinderTriplet* LxFinderTriplet::fInstance = 0;
23 TH1F* LxFinderTriplet::massHisto = 0;
24 //TH1F* LxFinderTriplet::backgroundMassHisto = 0;
25 CbmStsTrack* LxFinderTriplet::superEventData = 0;
26 TTree* LxFinderTriplet::superEventTracks = 0;
27 TH1F* LxFinderTriplet::signalChi2Histo = 0;
28 TH1F* LxFinderTriplet::backgroundChi2Histo = 0;
29 TH1F* LxFinderTriplet::unknownChi2Histo = 0;
30 TProfile* LxFinderTriplet::effByMomentumProfile = 0;
31 TH1F* LxFinderTriplet::signalInterTracksDistance = 0;
32 TH1F* LxFinderTriplet::bgrInterTracksDistance = 0;
33 TH1F* LxFinderTriplet::signalSignDefect = 0;
34 TH1F* LxFinderTriplet::bgrSignDefect = 0;
35 TH1F* LxFinderTriplet::signalYAtZ0 = 0;
36 TH1F* LxFinderTriplet::bgrYAtZ0 = 0;
37 TH1F* LxFinderTriplet::signalYAtZ0_2 = 0;
38 TH1F* LxFinderTriplet::bgrYAtZ0_2 = 0;
39 
40 #ifdef MAKE_DISPERSE_2D_HISTOS
41 TProfile2D* disperseLHistos[LXSTATIONS];
42 TProfile2D* disperseRHistos[LXSTATIONS];
43 TProfile2D* disperseDHistos[LXSTATIONS];
44 #endif //MAKE_DISPERSE_2D_HISTOS
45 
46 static TH1F* stsTrackChi2 = 0;
47 static TH1F* stsTrackX = 0;
48 static TH1F* stsTrackY = 0;
49 
50 LxFinderTriplet* LxFinderTriplet::Instance() { return fInstance; }
51 
52 std::list<CbmStsTrack> positiveTracks;
53 std::list<CbmStsTrack> negativeTracks;
54 bool generateInvMass;
55 bool generateBackground;
56 bool generateChi2;
57 bool linkWithSts;
58 bool useMCPInsteadOfHits;
59 bool calcMiddlePoints;
60 Double_t cutCoeff;
61 bool pPtCut;
62 std::vector<LxMCPoint>
63  MCPoints; // Points should lay here in the same order as in listMuchPts.
64 std::vector<LxMCTrack>
65  MCTracks; // Tracks should lay here in the same order as in listMCTracks.
66 std::list<LxStsMCPoint> MCStsPoints;
67 std::list<LxStsMCPoint*> MCStsPointsByStations[8];
68 #ifdef MAKE_DISPERSE_2D_HISTOS
69 std::list<LxMCPoint*> MCPointsByStations[LXSTATIONS][LXLAYERS];
70 Double_t zCoordsByStations[LXSTATIONS][LXLAYERS];
71 #endif //MAKE_DISPERSE_2D_HISTOS
72 LxSpace caSpace;
73 std::map<Int_t, std::map<Int_t, int>> particleCounts;
74 #ifdef MAKE_EFF_CALC
75 std::ofstream
76  incomplete_events; // Events where not all tracks are reconstructed.
77 #endif //MAKE_EFF_CALC
78 Int_t eventNumber;
79 #ifdef CALC_MUCH_DETECTORS_EFF
80 Int_t mcPointsCount;
81 Int_t mcPointsTriggered;
82 #endif //CALC_MUCH_DETECTORS_EFF
83 
84 
85 LxFinderTriplet::LxFinderTriplet()
86  : muchPixelHits(0)
87  , listMCTracks(0)
88  , listMuchPts(0)
89  , listMuchClusters(0)
90  , listMuchPixelDigiMatches(0)
91  , listStsTracks(0)
92  , listStsMatches(0)
93  , listStsPts(0)
94  , listRecoTracks(0)
95  , effCounter(*this)
96  , extFitter()
97  , fPrimVtx(0)
98  , positiveTracks()
99  , negativeTracks()
100  , generateInvMass(false)
101  , generateBackground(false)
102  , generateChi2(false)
103  , linkWithSts(true)
104  , useMCPInsteadOfHits(false)
105  , calcMiddlePoints(true)
106  , cutCoeff(4.0)
107  , pPtCut(true)
108  , MCPoints()
109  , MCTracks()
110  , MCStsPoints()
111  , MCStsPointsByStations()
112 #ifdef MAKE_DISPERSE_2D_HISTOS
113  , MCPointsByStations()
114  , zCoordsByStations()
115 #endif //MAKE_DISPERSE_2D_HISTOS
116  , caSpace()
117  , particleCounts()
118 #ifdef MAKE_EFF_CALC
119  , incomplete_events()
120 #endif //MAKE_EFF_CALC
121  , eventNumber(0)
122 #ifdef CALC_MUCH_DETECTORS_EFF
123  , mcPointsCount(0)
124  , mcPointsTriggered(0)
125 #endif //CALC_MUCH_DETECTORS_EFF
126 {
127  fInstance = this;
128 }
129 
130 LxFinderTriplet::~LxFinderTriplet() {}
131 
132 TString lxFinderParticleType = "jpsi";
133 
134 static bool
135 GetHistoRMS(const char* histoNameBase, Int_t histoNumber, Double_t& retVal) {
136  char name[256];
137  char dir_name[256];
138  sprintf(dir_name, "configuration.%s", lxFinderParticleType.Data());
139  bool result = false;
140  sprintf(name, "%s/%s_%d.root", dir_name, histoNameBase, histoNumber);
141  TFile* curFile = TFile::CurrentFile();
142  TFile* f = new TFile(name);
143 
144  if (!f->IsZombie()) {
145  sprintf(name, "%s_%d", histoNameBase, histoNumber);
146  TH1F* h = static_cast<TH1F*>(f->Get(name));
147  retVal = h->GetRMS();
148  result = true;
149  }
150 
151  delete f;
152  TFile::CurrentFile() = curFile;
153  return result;
154 }
155 
156 static bool GetHistoCOV(const char* histoNameBase,
157  Int_t histoNumber,
158  Int_t axis1,
159  Int_t axis2,
160  Double_t& retVal) {
161  char name[256];
162  bool result = false;
163  sprintf(name, "configuration/%s_%d.root", histoNameBase, histoNumber);
164  TFile* curFile = TFile::CurrentFile();
165  TFile* f = new TFile(name);
166 
167  if (!f->IsZombie()) {
168  sprintf(name, "%s_%d", histoNameBase, histoNumber);
169  TH2F* h = static_cast<TH2F*>(f->Get(name));
170  retVal = h->GetCovariance(axis1, axis2);
171  result = true;
172  }
173 
174  delete f;
175  TFile::CurrentFile() = curFile;
176  return result;
177 }
178 
179 InitStatus LxFinderTriplet::Init() {
180  static Int_t nTimes = 1;
181  cout << "LxFinderTriplet::Init() called at " << nTimes++ << " time" << endl;
182 
183  superEventData = new CbmStsTrack;
184 
185  FairRootManager* fManager = FairRootManager::Instance();
186 
187  muchPixelHits =
188  LX_DYNAMIC_CAST<TClonesArray*>(fManager->GetObject("MuchPixelHit"));
189  listMCTracks = LX_DYNAMIC_CAST<TClonesArray*>(fManager->GetObject("MCTrack"));
190  listMuchPts =
191  LX_DYNAMIC_CAST<TClonesArray*>(fManager->GetObject("MuchPoint"));
192  listMuchClusters =
193  LX_DYNAMIC_CAST<TClonesArray*>(fManager->GetObject("MuchCluster"));
194  listMuchPixelDigiMatches =
195  LX_DYNAMIC_CAST<TClonesArray*>(fManager->GetObject("MuchDigiMatch"));
196  listStsTracks =
197  LX_DYNAMIC_CAST<TClonesArray*>(fManager->GetObject("StsTrack"));
198  listStsMatches =
199  LX_DYNAMIC_CAST<TClonesArray*>(fManager->GetObject("StsTrackMatch"));
200  listStsPts = LX_DYNAMIC_CAST<TClonesArray*>(fManager->GetObject("StsPoint"));
201  //fPrimVtx = LX_DYNAMIC_CAST<CbmVertex*> (fManager->GetObject("PrimaryVertex"));
202  /*
203  // Get pointer to PrimaryVertex object from IOManager if it exists
204  // The old name for the object is "PrimaryVertex" the new one
205  // "PrimaryVertex." Check first for the new name
206  fPrimVtx = LX_DYNAMIC_CAST<CbmVertex*>(fManager->GetObject("PrimaryVertex."));
207  if (nullptr == fPrimVtx) {
208  fPrimVtx = LX_DYNAMIC_CAST<CbmVertex*>(fManager->GetObject("PrimaryVertex"));
209  }
210  if (nullptr == fPrimVtx) {
211  Error("CbmL1SttTrackFinder::ReInit","vertex not found!");
212  return kERROR;
213  }
214  */
215  fPrimVtx = new CbmVertex;
216 
217  // Read Z-positions of MUCH station layers and save them in LxLayer objects.
218 
219  for (Int_t i = 0; i < LXSTATIONS; ++i) {
220  LxStation* aStation = caSpace.stations[i];
221 
222  for (Int_t j = 0; j < LXLAYERS; ++j) {
223  TString muchStLrPath =
224  Form("/cave_1/much_0/muchstation%02i_0/muchstation%02ilayer%01i_0",
225  i + 1,
226  i + 1,
227  j + 1);
228  gGeoManager->cd(muchStLrPath.Data());
229  Double_t localCoords[3] = {0., 0., 0.};
230  Double_t globalCoords[3];
231  gGeoManager->LocalToMaster(localCoords, globalCoords);
232  aStation->layers[j]->zCoord = globalCoords[2];
233 #ifdef MAKE_DISPERSE_2D_HISTOS
234  zCoordsByStations[i][j] = globalCoords[2];
235 #endif //MAKE_DISPERSE_2D_HISTOS
236  }
237 
238  aStation->zCoord = aStation->layers[LXMIDDLE]->zCoord;
239 
240  bool readHistoResult = true;
241 
242  readHistoResult &=
243  GetHistoRMS("muchInStationXDispLeft", i, aStation->xDispLeft);
244  aStation->xDispLeft2 = aStation->xDispLeft * aStation->xDispLeft;
245  readHistoResult &=
246  GetHistoRMS("muchInStationYDispLeft", i, aStation->yDispLeft);
247  aStation->yDispLeft2 = aStation->yDispLeft * aStation->yDispLeft;
248  readHistoResult &=
249  GetHistoRMS("muchInStationXDispRight", i, aStation->xDispRight);
250  aStation->xDispRight2 = aStation->xDispRight * aStation->xDispRight;
251  readHistoResult &=
252  GetHistoRMS("muchInStationYDispRight", i, aStation->yDispRight);
253  aStation->yDispRight2 = aStation->yDispRight * aStation->yDispRight;
254  readHistoResult &=
255  GetHistoRMS("muchInStationXDispRL", i, aStation->xDispRL);
256  aStation->xDispRL2 = aStation->xDispRL * aStation->xDispRL;
257  readHistoResult &=
258  GetHistoRMS("muchInStationYDispRL", i, aStation->yDispRL);
259  aStation->yDispRL2 = aStation->yDispRL * aStation->yDispRL;
260 
261  if (i > 0) {
262 #ifdef OUT_DISP_BY_TRIPLET_DIR
263  readHistoResult &= GetHistoRMS(
264  "muchOutStationXDispByTriplet", i, aStation->xOutDispTriplet);
265  aStation->xOutDispTriplet2 =
266  aStation->xOutDispTriplet * aStation->xOutDispTriplet;
267  readHistoResult &= GetHistoRMS(
268  "muchOutStationYDispByTriplet", i, aStation->yOutDispTriplet);
269  aStation->yOutDispTriplet2 =
270  aStation->yOutDispTriplet * aStation->yOutDispTriplet;
271 #else //OUT_DISP_BY_TRIPLET_DIR
272  readHistoResult &=
273  GetHistoRMS("muchOutStationXDispByVertex", i, aStation->xOutDispVertex);
274  aStation->xOutDispVertex2 =
275  aStation->xOutDispVertex * aStation->xOutDispVertex;
276  readHistoResult &=
277  GetHistoRMS("muchOutStationYDispByVertex", i, aStation->yOutDispVertex);
278  aStation->yOutDispVertex2 =
279  aStation->yOutDispVertex * aStation->yOutDispVertex;
280 #endif //OUT_DISP_BY_TRIPLET_DIR
281 
282 #ifdef USE_SEGMENTS
283  LxStation* anLStation = caSpace.stations[i - 1];
284  readHistoResult &=
285  GetHistoRMS("muchOutStationTxBreakLeft", i, anLStation->txBreakRight);
286  anLStation->txBreakRight2 =
287  anLStation->txBreakRight * anLStation->txBreakRight;
288  readHistoResult &=
289  GetHistoRMS("muchOutStationTyBreakLeft", i, anLStation->tyBreakRight);
290  anLStation->tyBreakRight2 =
291  anLStation->tyBreakRight * anLStation->tyBreakRight;
292  readHistoResult &=
293  GetHistoRMS("muchOutStationTxBreakRight", i, aStation->txBreakLeft);
294  aStation->txBreakLeft2 = aStation->txBreakLeft * aStation->txBreakLeft;
295  readHistoResult &=
296  GetHistoRMS("muchOutStationTyBreakRight", i, aStation->tyBreakLeft);
297  aStation->tyBreakLeft2 = aStation->tyBreakLeft * aStation->tyBreakLeft;
298 #else //USE_SEGMENTS
299  readHistoResult &=
300  GetHistoRMS("muchTripletTxBreak", i, aStation->txInterTripletBreak);
301  aStation->txInterTripletBreak2 =
302  aStation->txInterTripletBreak * aStation->txInterTripletBreak;
303  readHistoResult &=
304  GetHistoRMS("muchTripletTyBreak", i, aStation->tyInterTripletBreak);
305  aStation->tyInterTripletBreak2 =
306  aStation->tyInterTripletBreak * aStation->tyInterTripletBreak;
307 #endif //USE_SEGMENTS
308 
309 #ifdef CLUSTER_MODE
310  if (!GetHistoRMS("muchClusterXDispHisto", i, rms)) return kFATAL;
311 
312  aStation->clusterXLimit = rms;
313  aStation->clusterXLimit2 = rms * rms;
314 
315  if (!GetHistoRMS("muchClusterYDispHisto", i, rms)) return kFATAL;
316 
317  aStation->clusterYLimit = rms;
318  aStation->clusterYLimit2 = rms * rms;
319 
320  if (!GetHistoRMS("muchClusterTxDispHisto", i, rms)) return kFATAL;
321 
322  aStation->clusterTxLimit = rms;
323  aStation->clusterTxLimit2 = rms * rms;
324 
325  if (!GetHistoRMS("muchClusterTyDispHisto", i, rms)) return kFATAL;
326 
327  aStation->clusterTyLimit = rms;
328  aStation->clusterTyLimit2 = rms * rms;
329 #endif //CLUSTER_MODE
330 
331 #ifdef USE_SEGMENTS
332  if (i < LXSTATIONS - 1) {
333  readHistoResult &= GetHistoRMS(
334  "muchSegmentTxBreakHisto", i, aStation->txInterStationBreak);
335  aStation->txInterStationBreak2 =
336  aStation->txInterStationBreak * aStation->txInterStationBreak;
337  readHistoResult &= GetHistoRMS(
338  "muchSegmentTyBreakHisto", i, aStation->tyInterStationBreak);
339  aStation->tyInterStationBreak2 =
340  aStation->tyInterStationBreak * aStation->tyInterStationBreak;
341  }
342 #endif //USE_SEGMENTS
343  }
344 
345  if (!readHistoResult) return kFATAL;
346 
347 #ifdef USE_KALMAN_FIT
348  if (i < LXSTATIONS - 1) {
349  for (Int_t j = 0; j <= 1; ++j) {
350  for (Int_t k = 0; k <= 1; ++k) {
351  //if (!GetHistoCOV("muchXTxCovHisto", i, j + 1, k + 1, aStation->MSNoise[0][j][k]))
352  //return kFATAL;
353 
354  //if (!GetHistoCOV("muchYTyCovHisto", i, j + 1, k + 1, aStation->MSNoise[1][j][k]))
355  //return kFATAL;
356  }
357  }
358  }
359 #endif //USE_KALMAN_FIT
360  }
361 
362  // Create an output array.
363  listRecoTracks = new TClonesArray("CbmMuchTrack", 100);
364  fManager->Register(
365  "MuchTrack", "Much", listRecoTracks, IsOutputBranchPersistent("MuchTrack"));
366 
367  if (generateInvMass)
368  massHisto = new TH1F("jpsi_mass", "jpsi_mass", 100, 2., 4.);
369 
370  if (generateBackground) {
371  //backgroundMassHisto = new TH1F("background_mass", "background_mass", 400, 2., 4.);
372  superEventTracks =
373  new TTree("SuperEventTracks", "Tracks for building a super event");
374  superEventTracks->Branch("tracks", "CbmStsTrack", &superEventData);
375  }
376 
377  //if (generateChi2)
378  //{
379  signalChi2Histo = new TH1F("signal_chi2", "signal_chi2", 200, 0., 20.);
380  signalChi2Histo->StatOverflows();
381  backgroundChi2Histo =
382  new TH1F("background_chi2", "background_chi2", 200, 0., 20.);
383  backgroundChi2Histo->StatOverflows();
384  unknownChi2Histo = new TH1F("unknown_chi2", "unknown_chi2", 200, 0., 20.);
385  unknownChi2Histo->StatOverflows();
386  //}
387 
388 #ifdef MAKE_HISTOS
389  char histoName[128];
390 
391  for (int i = 0; i < 6; ++i) {
392 #ifdef MAKE_DISPERSE_2D_HISTOS
393  sprintf(histoName, "disperseL_%d", i);
394  disperseLHistos[i] =
395  new TProfile2D(histoName, histoName, 30, -3., 3., 30, -3., 3.);
396  disperseLHistos[i]->StatOverflows();
397  sprintf(histoName, "disperseR_%d", i);
398  disperseRHistos[i] =
399  new TProfile2D(histoName, histoName, 30, -3., 3., 30, -3., 3.);
400  disperseRHistos[i]->StatOverflows();
401  sprintf(histoName, "disperseD_%d", i);
402  disperseDHistos[i] =
403  new TProfile2D(histoName, histoName, 30, -3., 3., 30, -3., 3.);
404  disperseDHistos[i]->StatOverflows();
405 #endif //MAKE_DISPERSE_2D_HISTOS
406  }
407 
408 #endif //MAKE_HISTOS
409 
410  effByMomentumProfile = new TProfile("Reconstruction efficiency by momentum",
411  "Reconstruction efficiency by momentum",
412  100,
413  0.0,
414  50.0,
415  0.0,
416  100.0);
417 
418  stsTrackChi2 = new TH1F("stsTrackChi2", "stsTrackChi2", 200, 0., 10.);
419  stsTrackChi2->StatOverflows();
420  stsTrackX = new TH1F("stsTrackX", "stsTrackX", 200, -0.2, 0.2);
421  stsTrackX->StatOverflows();
422  stsTrackY = new TH1F("stsTrackY", "stsTrackY", 200, -0.2, 0.2);
423  stsTrackY->StatOverflows();
424 
425  signalInterTracksDistance = new TH1F("signalInterTracksDistance",
426  "Distance between signal tracks",
427  200,
428  0.0,
429  40.0);
430  signalInterTracksDistance->StatOverflows();
431  bgrInterTracksDistance = new TH1F("bgrInterTracksDistance",
432  "Distance between background tracks",
433  200,
434  0.0,
435  40.0);
436  bgrInterTracksDistance->StatOverflows();
437  signalSignDefect =
438  new TH1F("signalSignDefect", "signalSignDefect", 200, -0.15, 0.15);
439  signalSignDefect->StatOverflows();
440  bgrSignDefect = new TH1F("bgrSignDefect", "bgrSignDefect", 200, -0.15, 0.15);
441  bgrSignDefect->StatOverflows();
442  signalYAtZ0 =
443  new TH1F("signalYAtZ0", "Signal track Y at Z=0", 200, -40.0, 40.0);
444  signalYAtZ0->StatOverflows();
445  bgrYAtZ0 =
446  new TH1F("bgrYAtZ0", "Background track Y at Z=0", 200, -40.0, 40.0);
447  bgrYAtZ0->StatOverflows();
448  signalYAtZ0_2 =
449  new TH1F("signalYAtZ0_2", "Signal track Y at Z=0 (2)", 200, -40.0, 40.0);
450  signalYAtZ0_2->StatOverflows();
451  bgrYAtZ0_2 =
452  new TH1F("bgrYAtZ0_2", "Background track Y at Z=0 (2)", 200, -40.0, 40.0);
453  bgrYAtZ0_2->StatOverflows();
454 
455 #ifdef MAKE_EFF_CALC
456  incomplete_events.open("incomplete_events.txt");
457 #endif //MAKE_EFF_CALC
458 
459  return kSUCCESS;
460 }
461 
462 static Int_t nTimes = 1;
463 
464 void LxFinderTriplet::Exec(Option_t*) {
465  cout << "LxFinderTriplet::Exec() called at " << nTimes++ << " time" << endl;
466  timeval bTime, eTime;
467  int exeDuration;
468  gettimeofday(&bTime, 0);
469 
470  // Clean all previously created data.
471  MCTracks.clear();
472  MCPoints.clear();
473  MCStsPoints.clear();
474 
475  for (int i = 0; i < 8; ++i)
476  MCStsPointsByStations[i].clear();
477 
478 #ifdef MAKE_DISPERSE_2D_HISTOS
479  for (int i = 0; i < LXSTATIONS; ++i) {
480  for (int j = 0; j < LXLAYERS; ++j)
481  MCPointsByStations[i][j].clear();
482  }
483 #endif //MAKE_DISPERSE_2D_HISTOS
484 
485  caSpace.Clear();
486  listRecoTracks->Clear();
487 
488  // Read MC tracks
489  Int_t nEnt = listMCTracks->GetEntries();
490  cout << "There are: " << nEnt << " of MC tracks" << endl;
491  LxMCTrack mcTrack;
492 
493  MCTracks.reserve(nEnt);
494  Int_t* root2lxmctrackmap = new Int_t[nEnt];
495  Int_t mapCnt = 0;
496 
497  for (int i = 0; i < nEnt; ++i) {
498  CbmMCTrack* mct = LX_DYNAMIC_CAST<CbmMCTrack*>(listMCTracks->At(i));
499  mcTrack.p = mct->GetP();
500  Int_t pdgCode = mct->GetPdgCode();
501 
502  if (abs(pdgCode) >= 10000) {
503  root2lxmctrackmap[i] = -1;
504  continue;
505  }
506 
507  root2lxmctrackmap[i] = mapCnt++;
508 
509  mcTrack.q = TDatabasePDG::Instance()->GetParticle(pdgCode)->Charge() / 3.0;
510  mcTrack.x = mct->GetStartX();
511  mcTrack.y = mct->GetStartY();
512  mcTrack.z = mct->GetStartZ();
513  mcTrack.px = mct->GetPx();
514  mcTrack.py = mct->GetPy();
515  mcTrack.pz = mct->GetPz();
516  mcTrack.mother_ID = mct->GetMotherId();
517  mcTrack.fUniqueID = mct->GetUniqueID();
518  mcTrack.pdg = pdgCode;
519  mcTrack.externalTrack = 0;
520  MCTracks.push_back(mcTrack);
521  }
522 
523  nEnt = listStsPts->GetEntries();
524  cout << "There are: " << nEnt << " of STS MC points" << endl;
525 
526  for (int i = 0; i < nEnt; ++i) {
527  TVector3 xyzI, PI, xyzO, PO;
528  CbmStsPoint* stsPt = LX_DYNAMIC_CAST<CbmStsPoint*>(listStsPts->At(i));
529 
530  if (0 == stsPt) continue;
531 
532  LxStsMCPoint stsMCPoint;
533  stsPt->Position(xyzI);
534  stsPt->Momentum(PI);
535  stsPt->PositionOut(xyzO);
536  stsPt->MomentumOut(PO);
537  TVector3 xyz = .5 * (xyzI + xyzO);
538  TVector3 P = .5 * (PI + PO);
539  stsMCPoint.x = xyz.X();
540  stsMCPoint.y = xyz.Y();
541  stsMCPoint.z = xyz.Z();
542  stsMCPoint.px = P.X();
543  stsMCPoint.py = P.Y();
544  stsMCPoint.pz = P.Z();
545  stsMCPoint.p =
546  sqrt(fabs(stsMCPoint.px * stsMCPoint.px + stsMCPoint.py * stsMCPoint.py
547  + stsMCPoint.pz * stsMCPoint.pz));
548  stsMCPoint.stationNumber =
549  CbmStsAddress::GetElementId(stsPt->GetDetectorID(), kStsStation);
550  Int_t trackId = root2lxmctrackmap[stsPt->GetTrackID()];
551 
552  if (-1 != trackId) {
553  stsMCPoint.mcTrack = &MCTracks[trackId];
554  MCStsPoints.push_back(stsMCPoint);
555  MCTracks[trackId].stsPoints[stsMCPoint.stationNumber].push_back(
556  &MCStsPoints.back());
557  MCStsPointsByStations[stsMCPoint.stationNumber].push_back(
558  &MCStsPoints.back());
559  }
560  }
561 
562  nEnt = listMuchPts->GetEntries();
563  cout << "There are: " << nEnt << " of MUCH MC points" << endl;
564  LxMCPoint mcPoint;
565 
566  MCPoints.reserve(nEnt);
567  Int_t* root2lxmcpointmap = new Int_t
568  [nEnt]; // Unfortunately we have to use this map because in the loop
569  // below some iterations can not to produce a point.
570  mapCnt = 0;
571  // Int_t mcPtsCount = nEnt;
572  Int_t maxReferencedPtsIndex = 0;
573 
574  for (int i = 0; i < nEnt; ++i) {
575  TVector3 xyzI, PI, xyzO, PO;
576  CbmMuchPoint* pt = LX_DYNAMIC_CAST<CbmMuchPoint*>(listMuchPts->At(i));
577 
578  if (!pt) {
579  root2lxmcpointmap[i] = -1;
580  continue;
581  }
582 
583  Int_t trackId = root2lxmctrackmap[pt->GetTrackID()];
584 
585  if (-1 == trackId) {
586  root2lxmcpointmap[i] = -1;
587  continue;
588  }
589 
590  mcPoint.trackId = trackId;
591  mcPoint.track = &MCTracks[trackId];
592 
593  root2lxmcpointmap[i] = mapCnt++;
594 
595  pt->Position(xyzI);
596  pt->Momentum(PI);
597  pt->PositionOut(xyzO);
598  pt->MomentumOut(PO);
599  TVector3 xyz = .5 * (xyzI + xyzO);
600  TVector3 P = .5 * (PI + PO);
601  mcPoint.x = xyz.X();
602  mcPoint.y = xyz.Y();
603  mcPoint.z = xyz.Z();
604  mcPoint.px = P.X();
605  mcPoint.py = P.Y();
606  mcPoint.pz = P.Z();
607  mcPoint.p = sqrt(fabs(mcPoint.px * mcPoint.px + mcPoint.py * mcPoint.py
608  + mcPoint.pz * mcPoint.pz));
609  mcPoint.stationNumber =
610  CbmMuchGeoScheme::GetStationIndex(pt->GetDetectorId());
611  mcPoint.layerNumber = CbmMuchGeoScheme::GetLayerIndex(pt->GetDetectorId());
612  MCPoints.push_back(mcPoint);
613  Int_t ptId = root2lxmcpointmap[i];
614 
615  /*
616 #ifdef MAKE_HISTOS
617  Double_t trackPt2 = MCTracks[trackId].px * MCTracks[trackId].px + MCTracks[trackId].py * MCTracks[trackId].py;
618 #endif//MAKE_HISTOS
619 */
620 
621  MCTracks[trackId].Points.push_back(&MCPoints[ptId]);
622 #ifdef MAKE_DISPERSE_2D_HISTOS
623  MCPointsByStations[mcPoint.stationNumber][mcPoint.layerNumber].push_back(
624  &MCPoints[ptId]);
625 #endif //MAKE_DISPERSE_2D_HISTOS
626  }
627 
628 #ifdef MAKE_HISTOS
629  // Build angle (tangent) breaks distribution <
630  for (vector<LxMCTrack>::iterator i = MCTracks.begin(); i != MCTracks.end();
631  ++i) {
632  LxMCTrack& track = *i;
633 
634  if ((13 != track.pdg && -13 == track.pdg) || track.mother_ID >= 0) continue;
635 
636  if (track.p < 3) continue;
637 
638  Double_t pt2 = track.px * track.px + track.py * track.py;
639 
640  if (pt2 < 1) continue;
641 
642  LxMCPoint* points[LXSTATIONS][LXLAYERS];
643  memset(points, 0, sizeof(points));
644 
645  for (vector<LxMCPoint*>::iterator j = track.Points.begin();
646  j != track.Points.end();
647  ++j) {
648  LxMCPoint* point = *j;
649  points[point->stationNumber][point->layerNumber] = point;
650  }
651 
652  for (int j = 0; j < LXSTATIONS; ++j) {
653 #ifdef MAKE_DISPERSE_2D_HISTOS
654  if (0 != points[j][1]) {
655  LxMCPoint* mPoint = points[j][1];
656  Double_t tx = mPoint->x / mPoint->z;
657  Double_t ty = mPoint->y / mPoint->z;
658  Double_t diffZ = zCoordsByStations[j][0] - mPoint->z;
659  Double_t x = mPoint->x + tx * diffZ;
660  Double_t y = mPoint->y + ty * diffZ;
661  list<LxMCPoint*>& lPoints = MCPointsByStations[j][0];
662 
663  for (list<LxMCPoint*>::iterator k = lPoints.begin(); k != lPoints.end();
664  ++k) {
665  LxMCPoint* lPoint = *k;
666 
667  if (lPoint->trackId == mPoint->trackId)
668  disperseLHistos[j]->Fill(lPoint->x - x, lPoint->y - y, 100.);
669  else
670  disperseLHistos[j]->Fill(lPoint->x - x, lPoint->y - y, 0.);
671  }
672 
673  diffZ = zCoordsByStations[j][2] - mPoint->z;
674  x = mPoint->x + tx * diffZ;
675  y = mPoint->y + ty * diffZ;
676  list<LxMCPoint*>& rPoints = MCPointsByStations[j][2];
677 
678  for (list<LxMCPoint*>::iterator k = rPoints.begin(); k != rPoints.end();
679  ++k) {
680  LxMCPoint* rPoint = *k;
681 
682  if (rPoint->trackId == mPoint->trackId)
683  disperseRHistos[j]->Fill(rPoint->x - x, rPoint->y - y, 100.);
684  else
685  disperseRHistos[j]->Fill(rPoint->x - x, rPoint->y - y, 0.);
686  }
687  }
688 
689  if (0 != points[j][2]) {
690  LxMCPoint* rPoint = points[j][2];
691  Double_t tx = rPoint->x / rPoint->z;
692  Double_t ty = rPoint->y / rPoint->z;
693  Double_t diffZ = zCoordsByStations[j][0] - rPoint->z;
694  Double_t x = rPoint->x + tx * diffZ;
695  Double_t y = rPoint->y + ty * diffZ;
696  list<LxMCPoint*>& lPoints = MCPointsByStations[j][0];
697 
698  for (list<LxMCPoint*>::iterator k = lPoints.begin(); k != lPoints.end();
699  ++k) {
700  LxMCPoint* lPoint = *k;
701 
702  if (lPoint->trackId == rPoint->trackId)
703  disperseDHistos[j]->Fill(lPoint->x - x, lPoint->y - y, 100.);
704  else
705  disperseDHistos[j]->Fill(lPoint->x - x, lPoint->y - y, 0.);
706  }
707  }
708 #endif //MAKE_DISPERSE_2D_HISTOS
709  } //for (int j = 0; j < LXSTATIONS; ++j)
710  } //for (vector<LxMCTrack>::iterator i = MCTracks.begin(); i != MCTracks.end(); ++i)
711  // > angle (tangent) breaks distribution
712 #endif //MAKE_HISTOS
713 
714  // Read MUCH pixel hits.
715  nEnt = muchPixelHits->GetEntries();
716  cout << "There are: " << nEnt << " of MUCH pixel hits" << endl;
717 
718  Double_t minXErr = 1000;
719  Double_t minYErr = 1000;
720 
721  for (Int_t i = 0; i < nEnt; ++i) {
722  CbmMuchPixelHit* mh =
723  LX_DYNAMIC_CAST<CbmMuchPixelHit*>(muchPixelHits->At(i));
724 
725  Int_t stationNumber = CbmMuchGeoScheme::GetStationIndex(mh->GetAddress());
726  Int_t layerNumber = CbmMuchGeoScheme::GetLayerIndex(mh->GetAddress());
727  TVector3 pos, err;
728  mh->Position(pos);
729  mh->PositionError(err);
730  Double_t x = pos.X();
731  Double_t y = pos.Y();
732  Double_t z = pos.Z();
733 
734  if (x != x || y != y || z != z) // Test for NaN
735  continue;
736 
737  if (minXErr > err.X()) minXErr = err.X();
738 
739  if (minYErr > err.Y()) minYErr = err.Y();
740 
741  LxPoint* lxPoint = 0;
742 
743  if (!useMCPInsteadOfHits)
744  lxPoint = caSpace.AddPoint(
745  stationNumber, layerNumber, i, x, y, z, err.X(), err.Y(), err.Z());
746 
747  Int_t clusterId = mh->GetRefId();
748  CbmMuchCluster* cluster =
749  LX_DYNAMIC_CAST<CbmMuchCluster*>(listMuchClusters->At(clusterId));
750  Int_t nDigis = cluster->GetNofDigis();
751 
752  for (Int_t j = 0; j < nDigis; ++j) {
753  CbmMuchDigiMatch* digiMatch = LX_DYNAMIC_CAST<CbmMuchDigiMatch*>(
754  listMuchPixelDigiMatches->At(cluster->GetDigi(j)));
755  Int_t nMCs = digiMatch->GetNofLinks();
756 
757  for (Int_t k = 0; k < nMCs; ++k) {
758  const CbmLink& lnk = digiMatch->GetLink(k);
759  Int_t mcIndex = lnk.GetIndex();
760 
761  if (mcIndex > maxReferencedPtsIndex) maxReferencedPtsIndex = mcIndex;
762 
763  Int_t mcIndexMapped = root2lxmcpointmap[mcIndex];
764 
765  if (-1 == mcIndexMapped) continue;
766 
767  if (!useMCPInsteadOfHits) {
768  MCPoints[mcIndexMapped].lxPoints.push_back(lxPoint);
769  lxPoint->mcPoints.push_back(&MCPoints[mcIndexMapped]);
770  } else if (MCPoints[mcIndexMapped].lxPoints.empty()) {
771  LxMCPoint& mcp = MCPoints[mcIndexMapped];
772  lxPoint = caSpace.AddPoint(stationNumber,
773  layerNumber,
774  i,
775  mcp.x,
776  mcp.y,
777  mcp.z,
778  err.X(),
779  err.Y(),
780  err.Z());
781  MCPoints[mcIndexMapped].lxPoints.push_back(lxPoint);
782  lxPoint->mcPoints.push_back(&MCPoints[mcIndexMapped]);
783  }
784  }
785  }
786  } // for (Int_t i = 0; i < nEnt; ++i)
787 
788  cout << "minXErr = " << minXErr << " ; minYErr = " << minYErr << endl << endl;
789 
790  // At last start invocation of CA procedures.
791  gettimeofday(&eTime, 0);
792  exeDuration =
793  (eTime.tv_sec - bTime.tv_sec) * 1000000 + eTime.tv_usec - bTime.tv_usec;
794  cout << "Execution duration 1 was: " << exeDuration << endl;
795 
796  //if (calcMiddlePoints)
797  //caSpace.RestoreMiddlePoints();
798 
799  bTime.tv_sec = eTime.tv_sec;
800  bTime.tv_usec = eTime.tv_usec;
801  gettimeofday(&eTime, 0);
802  exeDuration =
803  (eTime.tv_sec - bTime.tv_sec) * 1000000 + eTime.tv_usec - bTime.tv_usec;
804  cout << "Execution duration 2 was: " << exeDuration << endl;
805 #ifdef CLUSTER_MODE
806  caSpace.BuildRays2();
807 #else //CLUSTER_MODE
808  //caSpace.BuildRays();
809 #endif //CLUSTER_MODE
810  bTime.tv_sec = eTime.tv_sec;
811  bTime.tv_usec = eTime.tv_usec;
812  gettimeofday(&eTime, 0);
813  exeDuration =
814  (eTime.tv_sec - bTime.tv_sec) * 1000000 + eTime.tv_usec - bTime.tv_usec;
815  cout << "Execution duration 3 was: " << exeDuration << endl;
816 #ifdef CLUSTER_MODE
817  caSpace.ConnectNeighbours2();
818 #else //CLUSTER_MODE
819  //caSpace.ConnectNeighbours();
820 #endif //CLUSTER_MODE
821  bTime.tv_sec = eTime.tv_sec;
822  bTime.tv_usec = eTime.tv_usec;
823  gettimeofday(&eTime, 0);
824  exeDuration =
825  (eTime.tv_sec - bTime.tv_sec) * 1000000 + eTime.tv_usec - bTime.tv_usec;
826  cout << "Execution duration 4 was: " << exeDuration << endl;
827 #ifdef CLUSTER_MODE
828  caSpace.Reconstruct2();
829 #else //CLUSTER_MODE
830  caSpace.Reconstruct();
831 #endif //CLUSTER_MODE
832  //caSpace.FitTracks();
833  bTime.tv_sec = eTime.tv_sec;
834  bTime.tv_usec = eTime.tv_usec;
835  gettimeofday(&eTime, 0);
836  exeDuration =
837  (eTime.tv_sec - bTime.tv_sec) * 1000000 + eTime.tv_usec - bTime.tv_usec;
838  cout << "Execution duration 5 was: " << exeDuration << endl;
839  caSpace.FitTracks();
840  //cout << "maxReferencedPtsIndex=" << maxReferencedPtsIndex << " mcPtsCount=" << mcPtsCount << " GEF=" << listMuchPts->GetEntriesFast() << endl;
841  //cout << "Hits=" << muchPixelHits->GetEntries() << " clusters=" << listMuchClusters->GetEntries() << " digi ms=" << listMuchPixelDigiMatches->GetEntries() << endl;
842  //sleep(3);
843 
845  // Join MUCH-tracks with STS-tracks.
847 
848  // STS-tracks fulfilling the limitations peculiar for muons are picked up in the loop below.
849  nEnt = listStsTracks->GetEntries();
850 
851  for (int i = 0; i < nEnt; ++i) {
852  CbmStsTrack* stsTrack = LX_DYNAMIC_CAST<CbmStsTrack*>(listStsTracks->At(i));
853 
854  Double_t lpa[5] = {stsTrack->GetParamLast()->GetX(),
855  stsTrack->GetParamLast()->GetY(),
856  stsTrack->GetParamLast()->GetTx(),
857  stsTrack->GetParamLast()->GetTy(),
858  stsTrack->GetParamLast()->GetQp()};
859 
860  if (lpa[0] != lpa[0] || lpa[1] != lpa[1] || lpa[2] != lpa[2]
861  || lpa[3] != lpa[3] || lpa[4] != lpa[4])
862  continue;
863 
864  CbmStsTrack aTrack = *stsTrack;
865 
866  Int_t pdgHypo = 13;
867 
868  if (stsTrack->GetParamLast()->GetTx() > stsTrack->GetParamFirst()->GetTx())
869  pdgHypo = -13;
870 
871  extFitter.DoFit(&aTrack, pdgHypo);
872  Double_t chi2Prim = extFitter.GetChiToVertex(&aTrack, fPrimVtx);
873  FairTrackParam params;
874  extFitter.Extrapolate(&aTrack, fPrimVtx->GetZ(), &params);
875  Double_t p = 1 / params.GetQp();
876  Double_t p2 = p * p;
877 
878  if (p2 < 3.0 * 3.0) continue;
879 
880  /*
881  Double_t xDelta = 0.05;//5.0 * sqrt(params.GetCovariance(0, 0));
882  Double_t yDelta = 0.05;//5.0 * sqrt(params.GetCovariance(1, 1));
883 
884  if (params.GetX() < -xDelta || params.GetX() > xDelta || params.GetY() < -yDelta || params.GetY() > yDelta)
885  continue;
886 */
887  Double_t tx2 = params.GetTx() * params.GetTx();
888  Double_t ty2 = params.GetTy() * params.GetTy();
889  Double_t pt2 = p2 * (tx2 + ty2) / (1 + tx2 + ty2);
890 
891  if (pt2 < 1.0) continue;
892 
893  LxExtTrack extTrack;
894  extTrack.track = stsTrack;
895  extTrack.extId = i;
896 
897  CbmTrackMatch* match =
898  LX_DYNAMIC_CAST<CbmTrackMatch*>(listStsMatches->At(i));
899  Int_t numberOfHits = match->GetNofTrueHits() + match->GetNofWrongHits()
900  + match->GetNofFakeHits();
901 
902  if (match->GetNofTrueHits() >= 0.7 * numberOfHits) {
903  Int_t mcTrackId = match->GetMCTrackId();
904  Int_t mappedId = root2lxmctrackmap[mcTrackId];
905 
906  if (-1 != mappedId) {
907  MCTracks[mappedId].externalTrack = stsTrack;
908  extTrack.mcTrack = &MCTracks[mappedId];
909 
910  if (13 == MCTracks[mappedId].pdg || -13 == MCTracks[mappedId].pdg) {
911  stsTrackChi2->Fill(chi2Prim);
912  stsTrackX->Fill(params.GetX());
913  stsTrackY->Fill(params.GetY());
914  }
915  }
916  }
917 
918  caSpace.extTracks.push_back(extTrack);
919  }
920 
921  cout << "External tracks are read" << endl;
922 
923  // Join reconstructed tracks with external tracks.
924  caSpace.JoinExtTracks();
925  cout << "External tracks are connected" << endl;
926 
927  // Measure a reconstruction efficiency.
928  effCounter.CalcRecoEff(linkWithSts);
929  MatchRecoToMC();
930 
931  // Draw some result.
932  //static LxDraw drawer;
933 
934  //drawer.ClearView();
935  //drawer.DrawMuch();
936  //drawer.SaveCanvas("Geometry_");
937 
938  //drawer.ClearView();
939  //drawer.DrawInputHits();
940  //drawer.SaveCanvas("Hits_");
941  //drawer.Ask();
942 
943  //drawer.ClearView();
944  //drawer.DrawInputHits();
945  //drawer.DrawMCTracks();
946  //drawer.DrawRecoTracks();
947  //drawer.DrawExtTracks();
948  //drawer.SaveCanvas("MC_");
949  //gPad->WaitPrimitive();
950 
951  /*drawer.ClearView();
952  drawer.DrawInputHits();
953  drawer.DrawMCTracks();
954  drawer.DrawRays();
955  gPad->WaitPrimitive();*/
956 
957  //drawer.ClearView();
958  //drawer.DrawInputHits();
959  //drawer.DrawRecoTracks();
960  //drawer.SaveCanvas("Reco_");
961  //drawer.Ask();
962 
963  //drawer.ClearView();
964  //drawer.DrawInputHits();
965  //drawer.DrawRecoTracks();
966  //drawer.DrawExtTracks();
967  //drawer.Ask();
968 
969  if (generateInvMass) CalcInvMass();
970 
971  if (generateBackground) SaveEventTracks();
972 
973  // Finally save reconstructed tracks.
974  //SaveRecoTracks();
975 
976  // Hopefully we don't need these two maps any more.
977  delete[] root2lxmctrackmap;
978  delete[] root2lxmcpointmap;
979  ++eventNumber;
980 }
981 
982 void LxFinderTriplet::SaveRecoTracks() {
983  Int_t trackNo = listRecoTracks->GetEntriesFast();
984 
985  for (list<LxTrack*>::iterator i = caSpace.tracks.begin();
986  i != caSpace.tracks.end();
987  ++i) {
988  LxTrack* recoTrack = *i;
989  LxExtTrack* stsTrack = recoTrack->externalTrack;
990 
991  if (0 == stsTrack) continue;
992 
993  CbmMuchTrack muchTrack;
994 
995  for (int j = 0; j < LXSTATIONS * LXLAYERS; ++j) {
996  LxPoint* pPoint = recoTrack->points[j];
997 
998  if (0 != pPoint) muchTrack.AddHit(pPoint->hitId, kMUCHPIXELHIT);
999  }
1000 
1001  muchTrack.SetChiSq(recoTrack->chi2);
1002  muchTrack.SetNDF(
1003  4 * LXSTATIONS); // Probably need to calculate it more accurately.
1004  muchTrack.SetPreviousTrackId(stsTrack->extId);
1005 
1006  FairTrackParam parFirst;
1007  parFirst = *stsTrack->track->GetParamLast();
1008  muchTrack.SetFlag(0); // Good track.
1009 
1010  FairTrackParam parLast(parFirst);
1011  muchTrack.SetParamFirst(&parFirst);
1012  muchTrack.SetParamLast(&parLast);
1013 
1014  new ((*listRecoTracks)[trackNo++]) CbmMuchTrack(muchTrack);
1015  }
1016 }
1017 
1018 void LxFinderTriplet::CalcInvMass() {
1019  for (list<LxTrack*>::iterator i = caSpace.tracks.begin();
1020  i != caSpace.tracks.end();
1021  ++i) {
1022  LxTrack* firstTrack = *i;
1023  LxMCTrack* mcTrack1 = firstTrack->mcTrack;
1024 
1025  if (0 == firstTrack->externalTrack) continue;
1026 
1027  LxMCTrack* mcTrack2 = firstTrack->externalTrack->mcTrack;
1028 
1029  if (mcTrack1) {
1030  Int_t pdg1 = mcTrack1 ? mcTrack1->pdg : 10000;
1031  map<Int_t, map<Int_t, int>>::iterator j = particleCounts.find(pdg1);
1032 
1033  if (j == particleCounts.end())
1034  j = particleCounts
1035  .insert(pair<Int_t, map<Int_t, int>>(pdg1, map<Int_t, int>()))
1036  .first;
1037 
1038  Int_t pdg2 = mcTrack2 ? mcTrack2->pdg : 10000;
1039 
1040  map<Int_t, int>::iterator k = j->second.find(pdg2);
1041 
1042  if (k != j->second.end())
1043  ++(k->second);
1044  else
1045  j->second.insert(pair<Int_t, int>(pdg2, 1));
1046  }
1047 
1048  CbmStsTrack t1 = *firstTrack->externalTrack->track;
1049 
1050  // Add entry to the respective chi2 histogram if needed. This can be used for placing chi2 cuts determining J/psi.
1051  // pdg == 13 | -13 -- it is a muon.
1052  // mother_ID < 0 -- it is a primary particle. URQMD doesn't produce primary muons. So it from PLUTO <=> is J/psi.
1053 
1054  if (generateChi2) {
1055  // Double_t normalizedChi2 = firstTrack->chi2 / (LXSTATIONS * 4);// length * 4 == NDF.
1056 
1057  //if (mcTrack2 && (mcTrack2->pdg == 13 || mcTrack2->pdg == -13) && mcTrack2->mother_ID < 0)
1058  //signalChi2Histo->Fill(normalizedChi2);
1059  //else
1060  //backgroundChi2Histo->Fill(normalizedChi2);
1061  }
1062 
1063  if (!mcTrack2 || (mcTrack2->pdg != 13 && mcTrack2->pdg != -13)
1064  || mcTrack2->mother_ID >= 0)
1065  continue;
1066 
1067  extFitter.DoFit(&t1, 13);
1068  Double_t chi2Prim = extFitter.GetChiToVertex(&t1, fPrimVtx);
1069  FairTrackParam t1param;
1070  extFitter.Extrapolate(&t1, fPrimVtx->GetZ(), &t1param);
1071 
1072  if (t1param.GetQp() <= 0) continue;
1073 
1074  Double_t p1 = 1 / t1param.GetQp();
1075  Double_t tx12 = t1param.GetTx() * t1param.GetTx();
1076  Double_t ty12 = t1param.GetTy() * t1param.GetTy();
1077  Double_t pt12 = p1 * p1 * (tx12 + ty12) / (1 + tx12 + ty12);
1078 
1079  if (pt12 < 1) continue;
1080 
1081  CbmKFTrack muPlus(t1);
1082 
1083  for (list<LxTrack*>::iterator j = caSpace.tracks.begin();
1084  j != caSpace.tracks.end();
1085  ++j) {
1086  LxTrack* secondTrack = *j;
1087 
1088  if (0 == secondTrack->externalTrack) continue;
1089 
1090  LxMCTrack* mcSecondTrack = secondTrack->externalTrack->mcTrack;
1091 
1092  if (!mcSecondTrack
1093  || (mcSecondTrack->pdg != 13 && mcSecondTrack->pdg != -13)
1094  || mcSecondTrack->mother_ID >= 0)
1095  continue;
1096 
1097  CbmStsTrack t2 = *secondTrack->externalTrack->track;
1098  extFitter.DoFit(&t2, 13);
1099  chi2Prim = extFitter.GetChiToVertex(&t2, fPrimVtx);
1100  FairTrackParam t2param;
1101  extFitter.Extrapolate(&t2, fPrimVtx->GetZ(), &t2param);
1102 
1103  if (t2param.GetQp() >= 0) continue;
1104 
1105  Double_t p2 = 1 / t2param.GetQp();
1106  Double_t tx22 = t2param.GetTx() * t2param.GetTx();
1107  Double_t ty22 = t2param.GetTy() * t2param.GetTy();
1108  Double_t pt22 = p2 * p2 * (tx22 + ty22) / (1 + tx22 + ty22);
1109 
1110  if (pt22 < 1) continue;
1111 
1112  CbmKFTrack muMinus(t2);
1113  vector<CbmKFTrackInterface*> kfData;
1114  kfData.push_back(&muPlus);
1115  kfData.push_back(&muMinus);
1116  /*CbmKFParticle DiMu;
1117  DiMu.Construct(kfData, 0);
1118  DiMu.TransportToDecayVertex();
1119  Double_t m, merr;
1120  DiMu.GetMass(m, merr);
1121 
1122  massHisto->Fill(m);*/
1123  }
1124  }
1125 }
1126 
1127 //#ifdef MAKE_HISTOS
1128 // It also deletes the histogram.
1129 static void SaveHisto(TH1* histo, const char* name) {
1130  TFile fh(name, "RECREATE");
1131  histo->Write();
1132  fh.Close();
1133  delete histo;
1134 }
1135 //#endif//MAKE_HISTOS
1136 
1137 void LxFinderTriplet::FinishTask() {
1138  if (generateInvMass) SaveInvMass();
1139 
1140  if (generateBackground) SaveBackground();
1141 
1142  /*if (generateChi2)
1143  {
1144  SaveSignalChi2();
1145  SaveBackgroundChi2();
1146  }
1147 
1148  cout << "Statistics on particles found:" << endl;
1149 
1150  for (map<Int_t, map<Int_t, int> >::iterator i = particleCounts.begin(); i != particleCounts.end(); ++i)
1151  {
1152  string muchParticle = (10000 <= i->first || -10000 >= i->first) ? "Undefined MUCH particle" : TDatabasePDG::Instance()->GetParticle(i->first)->GetName();
1153  cout << "For MUCH particle [ " << muchParticle << " ] ( PDG = " << i->first << " ) found STS joints to:" << endl;
1154 
1155  for (map<Int_t, int>::iterator j = i->second.begin(); j != i->second.end(); ++j)
1156  {
1157  string stsParticle = (10000 <= j->first || -10000 >= j->first) ? "Undefined STS particle" : TDatabasePDG::Instance()->GetParticle(j->first)->GetName();
1158  cout << "\t" << j->second << " [ " << stsParticle << " ] ( PDG = " << j->first << " )" << endl;
1159  }
1160  }*/
1161  TFile* curFile = TFile::CurrentFile();
1162 #ifdef MAKE_HISTOS
1163  char histoFileName[128];
1164 
1165  for (int i = 0; i < 6; ++i) {
1166 #ifdef MAKE_DISPERSE_2D_HISTOS
1167  sprintf(histoFileName, "disperseL_histo_%d.root", i);
1168  SaveHisto(disperseLHistos[i], histoFileName);
1169  sprintf(histoFileName, "disperseR_histo_%d.root", i);
1170  SaveHisto(disperseRHistos[i], histoFileName);
1171  sprintf(histoFileName, "disperseD_histo_%d.root", i);
1172  SaveHisto(disperseDHistos[i], histoFileName);
1173 #endif //MAKE_DISPERSE_2D_HISTOS
1174  }
1175 #endif //MAKE_HISTOS
1176 
1177  {
1178  TFile fh("effByMomentumProfile.root", "RECREATE");
1179  effByMomentumProfile->Write();
1180  fh.Close();
1181  delete effByMomentumProfile;
1182  }
1183 
1184  SaveHisto(stsTrackChi2, "stsTrackChi2.root");
1185  SaveHisto(stsTrackX, "stsTrackX.root");
1186  SaveHisto(stsTrackY, "stsTrackY.root");
1187 
1188  SaveHisto(signalInterTracksDistance, "signalInterTracksDistance.root");
1189  SaveHisto(bgrInterTracksDistance, "bgrInterTracksDistance.root");
1190  SaveHisto(signalChi2Histo, "signalChi2Histo.root");
1191  SaveHisto(backgroundChi2Histo, "backgroundChi2Histo.root");
1192  SaveHisto(unknownChi2Histo, "unknownChi2Histo.root");
1193  SaveHisto(signalSignDefect, "signalSignDefectHisto.root");
1194  SaveHisto(bgrSignDefect, "bgrSignDefectHisto.root");
1195  SaveHisto(signalYAtZ0, "signalYAtZ0.root");
1196  SaveHisto(bgrYAtZ0, "bgrYAtZ0.root");
1197  SaveHisto(signalYAtZ0_2, "signalYAtZ0_2.root");
1198  SaveHisto(bgrYAtZ0_2, "bgrYAtZ0_2.root");
1199 
1200  TFile::CurrentFile() = curFile;
1201 
1202 #ifdef MAKE_EFF_CALC
1203  incomplete_events.close();
1204 #endif //MAKE_EFF_CALC
1205 
1206  delete superEventData;
1207 
1208  FairTask::FinishTask();
1209 }
1210 
1211 void LxFinderTriplet::SaveInvMass() {
1212  TFile* curFile = TFile::CurrentFile();
1213 
1214  TFile fh("jpsi_inv_mass_histo.root", "RECREATE");
1215  massHisto->Write();
1216  fh.Close();
1217  delete massHisto;
1218 
1219  TFile::CurrentFile() = curFile;
1220 }
1221 
1222 void LxFinderTriplet::SaveBackground() {
1223  for (list<CbmStsTrack>::iterator i = positiveTracks.begin();
1224  i != positiveTracks.end();
1225  ++i) {
1226  CbmStsTrack& t1 = *i;
1227  CbmKFTrack muPlus(t1);
1228 
1229  for (list<CbmStsTrack>::iterator j = negativeTracks.begin();
1230  j != negativeTracks.end();
1231  ++j) {
1232  CbmStsTrack& t2 = *j;
1233  CbmKFTrack muMinus(t2);
1234  vector<CbmKFTrackInterface*> kfData;
1235  kfData.push_back(&muPlus);
1236  kfData.push_back(&muMinus);
1237  /*CbmKFParticle DiMu;
1238  DiMu.Construct(kfData, 0);
1239  DiMu.TransportToDecayVertex();
1240  Double_t m, merr;
1241  DiMu.GetMass(m, merr);*/
1242 
1243  //for (int k = 0; k < 1000; ++k)
1244  //backgroundMassHisto->Fill(m);
1245  }
1246  }
1247 
1248  TFile* curFile = TFile::CurrentFile();
1249 
1250  TFile fh("tracks_tree.root", "RECREATE");
1251  superEventTracks->Write();
1252  fh.Close();
1253  delete superEventTracks;
1254 
1255  TFile::CurrentFile() = curFile;
1256 }
LxFinderTriplet::MCPoints
std::vector< LxMCPoint > MCPoints
Definition: riplet/Lx.h:115
LxFinderTriplet::superEventData
static CbmStsTrack * superEventData
Definition: riplet/Lx.h:87
LxTrack
Definition: LxCA.h:268
LxMCPoint
Definition: Simple/LxMC.h:16
CbmMCTrack::GetMotherId
Int_t GetMotherId() const
Definition: CbmMCTrack.h:71
LxStation::layers
std::vector< LxLayer * > layers
Definition: LxCA.h:190
disperseLHistos
TProfile2D * disperseLHistos[LXSTATIONS]
Definition: riplet/Lx.cxx:41
LxFinderTriplet::MCStsPointsByStations
std::list< LxStsMCPoint * > MCStsPointsByStations[8]
Definition: riplet/Lx.h:119
CbmMCTrack::GetStartX
Double_t GetStartX() const
Definition: CbmMCTrack.h:75
h
Generates beam ions for transport simulation.
Definition: CbmBeamGenerator.h:17
zCoordsByStations
Double_t zCoordsByStations[LXSTATIONS][LXLAYERS]
Definition: riplet/Lx.cxx:70
CbmTrack::GetParamLast
const FairTrackParam * GetParamLast() const
Definition: CbmTrack.h:62
LxFinderTriplet::CalcInvMass
void CalcInvMass()
Definition: riplet/Lx.cxx:1018
LxFinderTriplet::extFitter
CbmStsKFTrackFitter extFitter
Definition: riplet/Lx.h:101
MCStsPointsByStations
std::list< LxStsMCPoint * > MCStsPointsByStations[8]
Definition: riplet/Lx.cxx:67
LxStation::yDispLeft
Double_t yDispLeft
Definition: LxCATriplets.h:148
LxFinderTriplet::signalYAtZ0_2
static TH1F * signalYAtZ0_2
Definition: riplet/Lx.h:99
CbmMuchPoint
Definition: CbmMuchPoint.h:21
f
float f
Definition: L1/vectors/P4_F32vec4.h:24
LxStation::yDispLeft2
Double_t yDispLeft2
Definition: LxCATriplets.h:152
LxFinderTriplet::eventNumber
Int_t eventNumber
Definition: riplet/Lx.h:130
LxStsMCPoint
Definition: Simple/Lx.h:67
CbmPixelHit::Position
void Position(TVector3 &pos) const
Copies hit position to pos.
Definition: CbmPixelHit.cxx:64
LxFinderTriplet::bgrSignDefect
static TH1F * bgrSignDefect
Definition: riplet/Lx.h:96
LxFinderTriplet::fInstance
static LxFinderTriplet * fInstance
Definition: riplet/Lx.h:73
stsTrackX
static TH1F * stsTrackX
Definition: riplet/Lx.cxx:47
CbmMatch::GetLink
const CbmLink & GetLink(Int_t i) const
Definition: CbmMatch.h:35
LxMCTrack::fUniqueID
UInt_t fUniqueID
Definition: Simple/LxMC.h:35
CbmMatch::GetNofLinks
Int_t GetNofLinks() const
Definition: CbmMatch.h:38
stsTrackChi2
static TH1F * stsTrackChi2
Definition: riplet/Lx.cxx:46
LxFinderTriplet
Definition: riplet/Lx.h:28
CbmStsKFTrackFitter::Extrapolate
void Extrapolate(CbmStsTrack *track, Double_t z, FairTrackParam *e_track)
Definition: CbmStsKFTrackFitter.cxx:134
LxMCTrack::q
scaltype q
Definition: Simple/LxMC.h:28
LxMCPoint::p
scaltype p
Definition: Simple/LxMC.h:17
LxFinderTriplet::Instance
static LxFinderTriplet * Instance()
Definition: riplet/Lx.cxx:50
memset
void memset(T *dest, T i, size_t num)
Definition: L1Grid.cxx:21
LxMCTrack::pdg
Int_t pdg
Definition: Simple/LxMC.h:30
LxFinderTriplet::effByMomentumProfile
static TProfile * effByMomentumProfile
Definition: riplet/Lx.h:92
sqrt
friend F32vec4 sqrt(const F32vec4 &a)
Definition: L1/vectors/P4_F32vec4.h:41
LxPoint::hitId
Int_t hitId
Definition: LxCA.h:61
LxFinderTriplet::listStsMatches
TClonesArray * listStsMatches
Definition: riplet/Lx.h:80
LxStsMCPoint::stationNumber
Int_t stationNumber
Definition: Simple/Lx.h:69
LxMCPoint::x
scaltype x
Definition: Simple/LxMC.h:17
CbmMuchCluster
Data container for MUCH clusters.
Definition: CbmMuchCluster.h:20
LxStation::yDispRight
Double_t yDispRight
Definition: LxCATriplets.h:150
LxFinderTriplet::muchPixelHits
TClonesArray * muchPixelHits
Definition: riplet/Lx.h:74
LxFinderTriplet::listRecoTracks
TClonesArray * listRecoTracks
Definition: riplet/Lx.h:83
CbmStsKFTrackFitter::GetChiToVertex
Double_t GetChiToVertex(CbmStsTrack *track, CbmVertex *vtx=0)
Definition: CbmStsKFTrackFitter.cxx:164
incomplete_events
std::ofstream incomplete_events
Definition: riplet/Lx.cxx:76
MCPoints
std::vector< LxMCPoint > MCPoints
Definition: riplet/Lx.cxx:63
CbmMCTrack::GetPdgCode
Int_t GetPdgCode() const
Definition: CbmMCTrack.h:70
LxSpace::tracks
std::list< LxTrack * > tracks
Definition: LxCA.h:326
disperseRHistos
TProfile2D * disperseRHistos[LXSTATIONS]
Definition: riplet/Lx.cxx:42
i
int i
Definition: L1/vectors/P4_F32vec4.h:25
LxFinderTriplet::positiveTracks
std::list< CbmStsTrack > positiveTracks
Definition: riplet/Lx.h:104
LxStation::yOutDispVertex2
Double_t yOutDispVertex2
Definition: LxCATriplets.h:170
CbmTrackMatch::GetNofFakeHits
Int_t GetNofFakeHits() const
Definition: CbmTrackMatch.h:53
LxMCTrack::p
scaltype p
Definition: Simple/LxMC.h:28
LxMCPoint::z
scaltype z
Definition: Simple/LxMC.h:17
CbmTrackMatch
Definition: CbmTrackMatch.h:18
CbmMuchPoint::GetDetectorId
Int_t GetDetectorId() const
Definition: CbmMuchPoint.h:69
LxStsMCPoint::y
scaltype y
Definition: Simple/Lx.h:68
LxMCPoint::py
scaltype py
Definition: Simple/LxMC.h:17
CbmTrack::SetPreviousTrackId
void SetPreviousTrackId(Int_t previousTrackId)
Definition: CbmTrack.h:72
LxMCTrack::Points
std::vector< LxMCPoint * > Points
Definition: Simple/LxMC.h:31
LxSpace::extTracks
std::list< LxExtTrack > extTracks
Definition: LxCA.h:327
LxFinderTriplet::effCounter
LxEff effCounter
Definition: riplet/Lx.h:84
LxMCPoint::layerNumber
Int_t layerNumber
Definition: Simple/LxMC.h:18
CbmHit::GetRefId
Int_t GetRefId() const
Definition: CbmHit.h:72
LXLAYERS
#define LXLAYERS
Definition: Simple/LxSettings.h:8
LxFinderTriplet::generateBackground
bool generateBackground
Definition: riplet/Lx.h:107
caSpace
LxSpace caSpace
Definition: riplet/Lx.cxx:72
LxFinderTriplet::listStsPts
TClonesArray * listStsPts
Definition: riplet/Lx.h:82
CbmMCTrack::GetPx
Double_t GetPx() const
Definition: CbmMCTrack.h:72
LxFinderTriplet::LxFinderTriplet
LxFinderTriplet()
Definition: riplet/Lx.cxx:85
SaveHisto
static void SaveHisto(TH1 *histo, const char *name)
Definition: riplet/Lx.cxx:1129
LxMCPoint::px
scaltype px
Definition: Simple/LxMC.h:17
CbmTrackMatch::GetNofTrueHits
Int_t GetNofTrueHits() const
Definition: CbmTrackMatch.h:47
CbmMCTrack::GetPy
Double_t GetPy() const
Definition: CbmMCTrack.h:73
LxStation::xOutDispVertex
Double_t xOutDispVertex
Definition: LxCATriplets.h:167
LxSpace
Definition: LxCA.h:309
LxFinderTriplet::MatchRecoToMC
void MatchRecoToMC()
Definition: riplet/LxEff.cxx:1052
CbmTrack::SetParamLast
void SetParamLast(const FairTrackParam *par)
Definition: CbmTrack.h:76
CbmStsPoint
Definition: CbmStsPoint.h:27
CbmTrack::SetNDF
void SetNDF(Int_t ndf)
Definition: CbmTrack.h:71
LxTrack::chi2
scaltype chi2
Definition: LxCA.h:283
LXMIDDLE
#define LXMIDDLE
Definition: Simple/LxSettings.h:10
LxMCTrack::y
scaltype y
Definition: Simple/LxMC.h:28
CbmStsKFTrackFitter::DoFit
Int_t DoFit(CbmStsTrack *track, Int_t pidHypo=211)
Definition: CbmStsKFTrackFitter.cxx:79
LxStation::xDispLeft2
Double_t xDispLeft2
Definition: LxCATriplets.h:151
CbmMuchTrack
Definition: CbmMuchTrack.h:16
lxFinderParticleType
TString lxFinderParticleType
Definition: riplet/Lx.cxx:132
LxMCPoint::trackId
Int_t trackId
Definition: Simple/LxMC.h:21
CbmMuchDigiMatch
Definition: CbmMuchDigiMatch.h:17
stsTrackY
static TH1F * stsTrackY
Definition: riplet/Lx.cxx:48
LxPoint
Definition: LxCA.h:52
LxStation::tyInterTripletBreak2
Double_t tyInterTripletBreak2
Definition: LxCATriplets.h:191
CbmTrack::SetFlag
void SetFlag(Int_t flag)
Definition: CbmTrack.h:69
CbmCluster::GetNofDigis
Int_t GetNofDigis() const
Number of digis in cluster.
Definition: CbmCluster.h:69
linkWithSts
bool linkWithSts
Definition: riplet/Lx.cxx:57
CbmMuchCluster.h
Data container for MUCH clusters.
MAKE_EFF_CALC
#define MAKE_EFF_CALC
Definition: Simple/LxSettings.h:32
LxStsMCPoint::x
scaltype x
Definition: Simple/Lx.h:68
CbmMCTrack::GetStartZ
Double_t GetStartZ() const
Definition: CbmMCTrack.h:77
LxFinderTriplet::Init
InitStatus Init()
Definition: riplet/Lx.cxx:179
LxExtTrack
Definition: LxCA.h:249
CbmStsPoint::PositionOut
void PositionOut(TVector3 &pos)
Definition: CbmStsPoint.h:96
LxFinderTriplet::superEventTracks
static TTree * superEventTracks
Definition: riplet/Lx.h:88
LxMCTrack::externalTrack
CbmStsTrack * externalTrack
Definition: Simple/LxMC.h:36
LxFinderTriplet::negativeTracks
std::list< CbmStsTrack > negativeTracks
Definition: riplet/Lx.h:105
PI
#define PI
Definition: HagedornSpectrum.h:5
LxFinderTriplet::SaveEventTracks
void SaveEventTracks()
Definition: riplet/LxUtils.cxx:6
LxFinderTriplet::massHisto
static TH1F * massHisto
Definition: riplet/Lx.h:85
CbmMuchGeoScheme::GetLayerIndex
static Int_t GetLayerIndex(Int_t address)
Definition: CbmMuchGeoScheme.h:71
LxExtTrack::mcTrack
LxMCTrack * mcTrack
Definition: LxCA.h:252
CbmMuchGeoScheme::GetStationIndex
static Int_t GetStationIndex(Int_t address)
Definition: CbmMuchGeoScheme.h:68
LxFinderTriplet::FinishTask
void FinishTask()
Definition: riplet/Lx.cxx:1137
LxMCPoint::stationNumber
Int_t stationNumber
Definition: Simple/LxMC.h:18
LxFinderTriplet::incomplete_events
std::ofstream incomplete_events
Definition: riplet/Lx.h:128
CbmMuchPoint.h
LxStation::tyInterTripletBreak
Double_t tyInterTripletBreak
Definition: LxCATriplets.h:189
LxStsMCPoint::py
scaltype py
Definition: Simple/Lx.h:68
CbmMuchTrack.h
LxStsMCPoint::mcTrack
LxMCTrack * mcTrack
Definition: Simple/Lx.h:70
positiveTracks
std::list< CbmStsTrack > positiveTracks
Definition: riplet/Lx.cxx:52
LxStation::txInterTripletBreak2
Double_t txInterTripletBreak2
Definition: LxCATriplets.h:190
LxFinderTriplet::bgrYAtZ0
static TH1F * bgrYAtZ0
Definition: riplet/Lx.h:98
LxFinderTriplet::signalSignDefect
static TH1F * signalSignDefect
Definition: riplet/Lx.h:95
LxStation::yDispRight2
Double_t yDispRight2
Definition: LxCATriplets.h:154
CbmStsAddress::GetElementId
UInt_t GetElementId(Int_t address, Int_t level)
Get the index of an element.
Definition: CbmStsAddress.cxx:180
CbmHit::GetAddress
Int_t GetAddress() const
Definition: CbmHit.h:73
CbmVertex
Definition: CbmVertex.h:26
LxSpace::stations
std::vector< LxStation * > stations
Definition: LxCA.h:325
CbmStsPoint::MomentumOut
void MomentumOut(TVector3 &mom)
Definition: CbmStsPoint.h:97
LxDraw.h
LxStation::yOutDispVertex
Double_t yOutDispVertex
Definition: LxCATriplets.h:168
LxFinderTriplet::listMuchClusters
TClonesArray * listMuchClusters
Definition: riplet/Lx.h:77
LxStation::xDispRight
Double_t xDispRight
Definition: LxCATriplets.h:149
LxFinderTriplet::zCoordsByStations
Double_t zCoordsByStations[LXSTATIONS][LXLAYERS]
Definition: riplet/Lx.h:122
LxFinderTriplet::SaveRecoTracks
void SaveRecoTracks()
Definition: riplet/Lx.cxx:982
GetHistoRMS
static bool GetHistoRMS(const char *histoNameBase, Int_t histoNumber, Double_t &retVal)
Definition: riplet/Lx.cxx:135
Lx.h
LxFinderTriplet::generateChi2
bool generateChi2
Definition: riplet/Lx.h:108
CbmTrack::SetChiSq
void SetChiSq(Double_t chiSq)
Definition: CbmTrack.h:70
LxFinderTriplet::MCStsPoints
std::list< LxStsMCPoint > MCStsPoints
Definition: riplet/Lx.h:118
LxMCPoint::track
LxMCTrack * track
Definition: Simple/LxMC.h:22
LxSpace::AddPoint
LxPoint * AddPoint(int stationNumber, int layerNumber, int hitId, scaltype x, scaltype y, scaltype z, scaltype dx, scaltype dy, scaltype dz)
Definition: LxCA.h:339
CbmMuchPoint::PositionOut
void PositionOut(TVector3 &pos) const
Definition: CbmMuchPoint.h:80
nTimes
static Int_t nTimes
Definition: riplet/Lx.cxx:462
CbmVertex::GetZ
Double_t GetZ() const
Definition: CbmVertex.h:70
LxFinderTriplet::listStsTracks
TClonesArray * listStsTracks
Definition: riplet/Lx.h:79
CbmKFTrack.h
LxFinderTriplet::signalYAtZ0
static TH1F * signalYAtZ0
Definition: riplet/Lx.h:97
LxSpace::Reconstruct
void Reconstruct()
Definition: LxCA.cxx:2107
LxFinderTriplet::generateInvMass
bool generateInvMass
Definition: riplet/Lx.h:106
LxTrack::mcTrack
LxMCTrack * mcTrack
Definition: LxCA.h:276
LxFinderTriplet::Exec
void Exec(Option_t *opt)
Definition: riplet/Lx.cxx:464
LxMCTrack::px
scaltype px
Definition: Simple/LxMC.h:28
LxPoint::mcPoints
std::list< LxMCPoint * > mcPoints
Definition: LxCA.h:63
CbmMuchPoint::MomentumOut
void MomentumOut(TVector3 &mom) const
Definition: CbmMuchPoint.h:81
LxFinderTriplet::~LxFinderTriplet
~LxFinderTriplet()
Definition: riplet/Lx.cxx:130
LxStation::txInterTripletBreak
Double_t txInterTripletBreak
Definition: LxCATriplets.h:188
cutCoeff
Double_t cutCoeff
Definition: riplet/Lx.cxx:60
CbmTrack::GetParamFirst
const FairTrackParam * GetParamFirst() const
Definition: CbmTrack.h:61
LxStsMCPoint::px
scaltype px
Definition: Simple/Lx.h:68
LxStation::xDispRight2
Double_t xDispRight2
Definition: LxCATriplets.h:153
LxEff::CalcRecoEff
void CalcRecoEff(bool joinExt)
Definition: riplet/LxEff.cxx:647
LxMCTrack::z
scaltype z
Definition: Simple/LxMC.h:28
generateChi2
bool generateChi2
Definition: riplet/Lx.cxx:56
CbmKFParticle.h
LxStsMCPoint::p
scaltype p
Definition: Simple/Lx.h:68
GetHistoCOV
static bool GetHistoCOV(const char *histoNameBase, Int_t histoNumber, Int_t axis1, Int_t axis2, Double_t &retVal)
Definition: riplet/Lx.cxx:156
LxExtTrack::extId
Int_t extId
Definition: LxCA.h:251
ClassImp
ClassImp(LxFinderTriplet) using namespace std
CbmTrackMatch::GetNofWrongHits
Int_t GetNofWrongHits() const
Definition: CbmTrackMatch.h:50
LxExtTrack::track
CbmStsTrack * track
Definition: LxCA.h:250
LxFinderTriplet::listMuchPixelDigiMatches
TClonesArray * listMuchPixelDigiMatches
Definition: riplet/Lx.h:78
LxFinderTriplet::SaveBackground
void SaveBackground()
Definition: riplet/Lx.cxx:1222
LxStation::xDispLeft
Double_t xDispLeft
Definition: LxCATriplets.h:147
negativeTracks
std::list< CbmStsTrack > negativeTracks
Definition: riplet/Lx.cxx:53
LxStation::xOutDispVertex2
Double_t xOutDispVertex2
Definition: LxCATriplets.h:169
LxMCTrack::x
scaltype x
Definition: Simple/LxMC.h:28
LxSpace::FitTracks
void FitTracks()
Definition: LxCA.cxx:2227
LxStation::xDispRL
Double_t xDispRL
Definition: LxCATriplets.h:156
LxFinderTriplet::signalInterTracksDistance
static TH1F * signalInterTracksDistance
Definition: riplet/Lx.h:93
LxSpace::JoinExtTracks
void JoinExtTracks()
Definition: LxCA.cxx:2232
LxFinderTriplet::signalChi2Histo
static TH1F * signalChi2Histo
Definition: riplet/Lx.h:89
CbmMCTrack::GetStartY
Double_t GetStartY() const
Definition: CbmMCTrack.h:76
LxMCTrack::pz
scaltype pz
Definition: Simple/LxMC.h:28
CbmStsPoint.h
points
TClonesArray * points
Definition: Analyze_matching.h:18
CbmMCTrack.h
LxFinderTriplet::MCPointsByStations
std::list< LxMCPoint * > MCPointsByStations[LXSTATIONS][LXLAYERS]
Definition: riplet/Lx.h:121
LxStation::zCoord
scaltype zCoord
Definition: LxCA.h:207
LxTrack::points
LxPoint * points[LXSTATIONS *LXLAYERS]
Definition: LxCA.h:282
LxStation
Definition: LxCA.h:189
LxFinderTriplet::bgrYAtZ0_2
static TH1F * bgrYAtZ0_2
Definition: riplet/Lx.h:100
LxFinderTriplet::SaveInvMass
void SaveInvMass()
Definition: riplet/Lx.cxx:1211
CbmMCTrack
Definition: CbmMCTrack.h:34
x
Double_t x
Definition: CbmMvdSensorDigiToHitTask.cxx:68
fabs
friend F32vec4 fabs(const F32vec4 &a)
Definition: L1/vectors/P4_F32vec4.h:60
LxStation::xDispRL2
Double_t xDispRL2
Definition: LxCATriplets.h:158
LxFinderTriplet::unknownChi2Histo
static TH1F * unknownChi2Histo
Definition: riplet/Lx.h:91
calcMiddlePoints
bool calcMiddlePoints
Definition: riplet/Lx.cxx:59
CbmTrackMatch::GetMCTrackId
Int_t GetMCTrackId() const
Definition: CbmTrackMatch.h:44
y
Double_t y
Definition: CbmMvdSensorDigiToHitTask.cxx:68
pPtCut
bool pPtCut
Definition: riplet/Lx.cxx:61
CbmTrack::SetParamFirst
void SetParamFirst(const FairTrackParam *par)
Definition: CbmTrack.h:75
LxStation::yDispRL2
Double_t yDispRL2
Definition: LxCATriplets.h:159
pos
TVector3 pos
Definition: CbmMvdSensorDigiToHitTask.cxx:60
LxFinderTriplet::fPrimVtx
CbmVertex * fPrimVtx
Definition: riplet/Lx.h:102
CbmStsAddress.h
generateInvMass
bool generateInvMass
Definition: riplet/Lx.cxx:54
kMUCHPIXELHIT
@ kMUCHPIXELHIT
Definition: CbmHit.h:23
LxMCPoint::pz
scaltype pz
Definition: Simple/LxMC.h:17
CbmMuchDigiMatch.h
LxFinderTriplet::MCTracks
std::vector< LxMCTrack > MCTracks
Definition: riplet/Lx.h:117
LxFinderTriplet::backgroundChi2Histo
static TH1F * backgroundChi2Histo
Definition: riplet/Lx.h:90
CbmMuchPixelHit
Definition: CbmMuchPixelHit.h:17
LxMCTrack
Definition: Simple/LxMC.h:27
LxFinderTriplet::listMCTracks
TClonesArray * listMCTracks
Definition: riplet/Lx.h:75
LxFinderTriplet::particleCounts
std::map< Int_t, std::map< Int_t, int > > particleCounts
Definition: riplet/Lx.h:125
LxFinderTriplet::linkWithSts
bool linkWithSts
Definition: riplet/Lx.h:109
LxFinderTriplet::caSpace
LxSpace caSpace
Definition: riplet/Lx.h:124
CbmMuchGeoScheme.h
generateBackground
bool generateBackground
Definition: riplet/Lx.cxx:55
MAKE_DISPERSE_2D_HISTOS
#define MAKE_DISPERSE_2D_HISTOS
Definition: riplet/LxSettings.h:18
LxFinderTriplet::useMCPInsteadOfHits
bool useMCPInsteadOfHits
Definition: riplet/Lx.h:110
MCPointsByStations
std::list< LxMCPoint * > MCPointsByStations[LXSTATIONS][LXLAYERS]
Definition: riplet/Lx.cxx:69
CbmStsTrack
Definition: CbmStsTrack.h:37
CbmKFTrack
Definition: CbmKFTrack.h:21
LxSpace::Clear
void Clear()
Definition: LxCA.cxx:1834
LxStation::yDispRL
Double_t yDispRL
Definition: LxCATriplets.h:157
LXSTATIONS
#define LXSTATIONS
Definition: Simple/LxSettings.h:9
CbmCluster::GetDigi
Int_t GetDigi(Int_t index) const
Get digi at position index.
Definition: CbmCluster.h:76
CbmTrack::AddHit
void AddHit(Int_t index, HitType type)
Definition: CbmTrack.cxx:75
CbmPixelHit::PositionError
void PositionError(TVector3 &dpos) const
Copies hit position error to pos.
Definition: CbmPixelHit.cxx:68
LxMCPoint::y
scaltype y
Definition: Simple/LxMC.h:17
eventNumber
Int_t eventNumber
Definition: riplet/Lx.cxx:78
CbmMCTrack::GetP
Double_t GetP() const
Definition: CbmMCTrack.h:100
LxFinderTriplet::bgrInterTracksDistance
static TH1F * bgrInterTracksDistance
Definition: riplet/Lx.h:94
useMCPInsteadOfHits
bool useMCPInsteadOfHits
Definition: riplet/Lx.cxx:58
disperseDHistos
TProfile2D * disperseDHistos[LXSTATIONS]
Definition: riplet/Lx.cxx:43
LxStsMCPoint::z
scaltype z
Definition: Simple/Lx.h:68
CbmMCTrack::GetPz
Double_t GetPz() const
Definition: CbmMCTrack.h:74
LxMCTrack::mother_ID
Int_t mother_ID
Definition: Simple/LxMC.h:29
LxFinderTriplet::listMuchPts
TClonesArray * listMuchPts
Definition: riplet/Lx.h:76
LxMCTrack::py
scaltype py
Definition: Simple/LxMC.h:28
LxStsMCPoint::pz
scaltype pz
Definition: Simple/Lx.h:68
particleCounts
std::map< Int_t, std::map< Int_t, int > > particleCounts
Definition: riplet/Lx.cxx:73
LxTrack::externalTrack
LxExtTrack * externalTrack
Definition: LxCA.h:269
MCTracks
std::vector< LxMCTrack > MCTracks
Definition: riplet/Lx.cxx:65
MCStsPoints
std::list< LxStsMCPoint > MCStsPoints
Definition: riplet/Lx.cxx:66