BmnRoot
Loading...
Searching...
No Matches
BmnGemStripHitMaker.cxx
Go to the documentation of this file.
1#include "BmnGemStripHitMaker.h"
2
3#include "BmnEventHeader.h"
4#include "FairLinkManager.h"
5#include "FairRunAna.h"
6#include "TSystem.h"
7#include "UniDetectorParameter.h"
8
9#include <TChain.h>
10#include <TRandom.h>
11#include <TStopwatch.h>
12#include <zlib.h>
13
14static Double_t workTime = 0.0;
15
16BmnGemStripHitMaker::BmnGemStripHitMaker()
17 : fHitMatching(kTRUE)
18{
19
20 fInputPointsBranchName = "StsPoint";
21 fInputDigitsBranchName = "BmnGemStripDigit";
22 fInputDigitMatchesBranchName = "BmnGemStripDigitMatch";
23
24 fOutputHitsBranchName = "BmnGemStripHit";
25
26 fField = NULL;
27 fFieldScale = 0.0;
28
29 fCurrentConfig = BmnGemStripConfiguration::None;
30 StationSet = nullptr;
31 TransfSet = nullptr;
32
33 fBmnGemStripHitsArray = nullptr;
34 fBmnGemUpperClustersArray = nullptr;
35 fBmnGemLowerClustersArray = nullptr;
36}
37
38BmnGemStripHitMaker::BmnGemStripHitMaker(Int_t run_period, Int_t run_number, Bool_t isExp, Bool_t isSrc)
39 : fHitMatching(kTRUE)
40{
41
42 fInputPointsBranchName = "StsPoint";
43 fInputDigitsBranchName = (!isExp) ? "BmnGemStripDigit" : "GEM";
44
45 fRunPeriod = run_period;
46 fRunNumber = run_number;
47
48 fIsExp = isExp;
49 fIsSrc = isSrc;
50
51 fSignalLow = 0.;
52 fSignalUp = DBL_MAX;
53
54 fInputDigitMatchesBranchName = "BmnGemStripDigitMatch";
55
56 fOutputHitsBranchName = "BmnGemStripHit";
57
58 fField = NULL;
59 fFieldScale = 0.0;
60
61 fCurrentConfig = BmnGemStripConfiguration::None;
62 StationSet = nullptr;
63 TransfSet = nullptr;
64
65 fBmnGemStripHitsArray = nullptr;
66 fBmnGemUpperClustersArray = nullptr;
67 fBmnGemLowerClustersArray = nullptr;
68
69 switch (fRunPeriod) {
70 case 6: // BM@N RUN-6
71 fCurrentConfig = BmnGemStripConfiguration::RunSpring2017;
72 break;
73 case 7: // BM@N RUN-7 (and SRC)
74 if (fIsSrc) {
75 fCurrentConfig = BmnGemStripConfiguration::RunSRCSpring2018;
76 } else {
77 fCurrentConfig = BmnGemStripConfiguration::RunSpring2018;
78 }
79 break;
80 case 8: // BM@N RUN-8
81 if (fIsSrc) {
82 fCurrentConfig = BmnGemStripConfiguration::RunSRC2021;
83 } else {
84 fCurrentConfig = BmnGemStripConfiguration::Run8;
85 }
86 break;
87 case 9: // BM@N RUN-9
88 fCurrentConfig = BmnGemStripConfiguration::Run9;
89 break;
90
91 default:
92 fCurrentConfig = BmnGemStripConfiguration::None;
93 }
94}
95
96void BmnGemStripHitMaker::createGemDetector()
97{
98 TString gPathGemConfig = gSystem->Getenv("VMCWORKDIR");
99 gPathGemConfig += "/parameters/gem/XMLConfigs/";
100
101 switch (fCurrentConfig) {
102 case BmnGemStripConfiguration::RunSpring2017:
103 XMLConfigFile = gPathGemConfig + "GemRunSpring2017.xml";
104 if (fVerbose)
105 cout << " Current GEM Configuration : RunSpring2017"
106 << "\n";
107 break;
108
109 case BmnGemStripConfiguration::RunSpring2018:
110 XMLConfigFile = gPathGemConfig + "GemRunSpring2018.xml";
111 if (fVerbose)
112 cout << " Current GEM Configuration : RunSpring2018"
113 << "\n";
114 break;
115
116 case BmnGemStripConfiguration::RunSRCSpring2018:
117
118 XMLConfigFile = gPathGemConfig + "GemRunSRCSpring2018.xml";
119 if (fVerbose)
120 cout << " Current GEM Configuration : GemRunSRCSpring2018"
121 << "\n";
122 break;
123
124 case BmnGemStripConfiguration::RunSRC2021:
125 XMLConfigFile = gPathGemConfig + "GemRunSRC2021.xml";
126 if (fVerbose)
127 cout << " Current GEM Configuration : RunSRC2021"
128 << "\n";
129 break;
130
131 case BmnGemStripConfiguration::Run8:
132 XMLConfigFile = gPathGemConfig + "GemRun8.xml";
133 if (fVerbose)
134 cout << " Current GEM Configuration : Run8"
135 << "\n";
136 break;
137
138 case BmnGemStripConfiguration::Run9:
139 XMLConfigFile = gPathGemConfig + "GemRun9.xml";
140 if (fVerbose)
141 cout << " Current GEM Configuration : Run9"
142 << "\n";
143 break;
144
145 default:;
146 }
147
148 if (!gSystem->AccessPathName(XMLConfigFile)) { // returns FALSE if the file exists
149 if (fVerbose)
150 cout << " XMLConfigFile : " << XMLConfigFile << "\n";
151 StationSet = new BmnGemStripStationSet(XMLConfigFile);
152 TransfSet = new BmnGemStripTransform();
153 TransfSet->LoadFromXMLFile(XMLConfigFile);
154 } else {
155 Fatal("BmnGemStripHitMaker::createGemDetector()", " !!! Current configuration is not set !!! ");
156 }
157
158 const Int_t nStat = StationSet->GetNStations();
159 fLorCor = new Double_t*[nStat];
160 for (Int_t iStat = 0; iStat < nStat; iStat++) {
161 const Int_t nParams = 3;
162 fLorCor[iStat] = new Double_t[nParams];
163 for (Int_t iPar = 0; iPar < nParams; iPar++) {
164 fLorCor[iStat][iPar] = 0.0;
165 }
166 }
167
168 if (fRunPeriod == 7) {
169 if (fIsExp) {
170 UniDetectorParameter* coeffLorCorrs =
171 UniDetectorParameter::GetDetectorParameter("GEM", "lorentz_shift", fRunPeriod, fRunNumber);
172 vector<UniValue*> shifts;
173 if (coeffLorCorrs)
174 coeffLorCorrs->GetValue(shifts);
175 for (Int_t iStat = 0; iStat < nStat; iStat++) {
176 const Int_t nParams = 3;
177 for (Int_t iPar = 0; iPar < nParams; iPar++) {
178 fLorCor[iStat][iPar] = (coeffLorCorrs) ? ((LorentzShiftValue*)shifts[iStat])->ls[iPar] : 0.0;
179 // printf("fLorCor[%d][%d] = %f\n", iStat, iPar, fLorCor[iStat][iPar]);
180 }
181 }
182 }
183 // AZ-110625 } else if (fRunPeriod == 8) {
184 } else if (fRunPeriod >= 8) { // AZ-110625
185 for (Int_t iStat = 0; iStat < nStat; iStat++) {
186 // Pol2 approximation of the next configuration:
187 // BmnGemStripMediumConfiguration::ARC4H10_80_20_E_1720_2240_3230_3730_B_0_8T
188 fLorCor[iStat][0] = -0.01710;
189 fLorCor[iStat][1] = 0.01730;
190 fLorCor[iStat][2] = 0.00024;
191 }
192 }
193
194 if (fIsExp) {
195 UniDetectorParameter* coeffAlignCorrs =
196 UniDetectorParameter::GetDetectorParameter("GEM", "alignment_shift", fRunPeriod, fRunNumber);
197 vector<UniValue*> algnShifts;
198 if (coeffAlignCorrs)
199 coeffAlignCorrs->GetValue(algnShifts);
200 fAlignCor = new Double_t**[nStat];
201 for (Int_t iStat = 0; iStat < nStat; iStat++) {
202 const Int_t nParams = 3;
203 Int_t nModul = StationSet->GetGemStation(iStat)->GetNModules();
204 fAlignCor[iStat] = new Double_t*[nModul];
205 for (Int_t iMod = 0; iMod < nModul; iMod++) {
206 fAlignCor[iStat][iMod] = new Double_t[nParams];
207 for (Int_t iPar = 0; iPar < nParams; iPar++) {
208 // cout << iStat << " " << iMod << ": " << ((AlignmentValue*)algnShifts[iStat * nModul +
209 // iMod])->value[iPar] << endl;
210 fAlignCor[iStat][iMod][iPar] =
211 (coeffAlignCorrs) ? ((AlignmentValue*)algnShifts[iStat * nModul + iMod])->value[iPar] : 0.0;
212 }
213 }
214 }
215 }
216}
217
218BmnGemStripHitMaker::~BmnGemStripHitMaker()
219{
220 if (fIsExp) {
221 for (Int_t iStat = 0; iStat < StationSet->GetNStations(); iStat++) {
222 for (Int_t iMod = 0; iMod < StationSet->GetGemStation(iStat)->GetNModules(); iMod++) {
223 delete[] fAlignCor[iStat][iMod];
224 }
225 delete[] fAlignCor[iStat];
226 }
227 delete[] fAlignCor;
228 }
229 for (Int_t iStat = 0; iStat < StationSet->GetNStations(); iStat++) {
230 delete[] fLorCor[iStat];
231 }
232 delete[] fLorCor;
233
234 if (StationSet) {
235 delete StationSet;
236 }
237
238 if (TransfSet) {
239 delete TransfSet;
240 }
241
242 if (fBmnGemStripHitsArray) {
243 fBmnGemStripHitsArray->Delete();
244 delete fBmnGemStripHitsArray;
245 }
246 if (fBmnGemUpperClustersArray) {
247 fBmnGemUpperClustersArray->Delete();
248 delete fBmnGemUpperClustersArray;
249 }
250 if (fBmnGemLowerClustersArray) {
251 fBmnGemLowerClustersArray->Delete();
252 delete fBmnGemLowerClustersArray;
253 }
254}
255
256InitStatus BmnGemStripHitMaker::Init()
257{
258
259 if (fVerbose > 1)
260 cout << "=================== BmnGemStripHitMaker::Init() started ===============" << endl;
261
262 createGemDetector();
263
264 FairRootManager* ioman = FairRootManager::Instance();
265
266 fBmnGemStripDigitsArray = (TClonesArray*)ioman->GetObject(fInputDigitsBranchName);
267 if (!fBmnGemStripDigitsArray) {
268 cout << "BmnGemStripHitMaker::Init(): branch " << fInputDigitsBranchName
269 << " not found! Task will be deactivated" << endl;
270 SetActive(kFALSE);
271 return kERROR;
272 }
273
274 fBmnGemStripDigitMatchesArray = (TClonesArray*)ioman->GetObject(fInputDigitMatchesBranchName);
275
276 if (fVerbose > 1) {
277 if (fBmnGemStripDigitMatchesArray)
278 cout << " Strip matching information exists!\n";
279 else
280 cout << " Strip matching information doesn`t exist!\n";
281 }
282
283 fBmnGemStripHitsArray = new TClonesArray(fOutputHitsBranchName);
284 ioman->Register(fOutputHitsBranchName, "GEM", fBmnGemStripHitsArray, kTRUE);
285 fBmnGemUpperClustersArray = new TClonesArray("StripCluster");
286 ioman->Register("BmnGemUpperCluster", "GEM", fBmnGemUpperClustersArray, kTRUE);
287 fBmnGemLowerClustersArray = new TClonesArray("StripCluster");
288 ioman->Register("BmnGemLowerCluster", "GEM", fBmnGemLowerClustersArray, kTRUE);
289
290 fField = FairRunAna::Instance()->GetField();
291 if (!fField)
292 // Fatal("Init", "No Magnetic Field found");
293 Warning("Init", "No Magnetic Field found");
294
295 //--------------------------------------------------------------------------
296
297 if (fVerbose > 1)
298 cout << "=================== BmnGemStripHitMaker::Init() finished ==============" << endl;
299
300 return kSUCCESS;
301}
302
303void BmnGemStripHitMaker::Exec(Option_t* opt)
304{
305
306 TStopwatch sw;
307 sw.Start();
308
309 if (!IsActive())
310 return;
311
312 fBmnGemStripHitsArray->Delete();
313 fBmnGemUpperClustersArray->Delete();
314 fBmnGemLowerClustersArray->Delete();
315
316 BmnGemStripLayer::SetLowerUniqueID(0);
317 BmnGemStripLayer::SetUpperUniqueID(0);
318
319 if (fVerbose > 1)
320 cout << "=================== BmnGemStripHitMaker::Exec() started ===============" << endl;
321
322 if (fVerbose > 1)
323 cout << " BmnGemStripHitMaker::Exec(), Number of BmnGemStripDigits = "
324 << fBmnGemStripDigitsArray->GetEntriesFast() << "\n";
325
326 ProcessDigits();
327
328 if (fVerbose > 1)
329 cout << "=================== BmnGemStripHitMaker::Exec() finished ==============" << endl;
330
331 sw.Stop();
332 workTime += sw.RealTime();
333}
334
335void BmnGemStripHitMaker::ProcessDigits()
336{
337
338 // FairMCPoint* MCPoint;
339 BmnGemStripDigit* digit;
340 BmnMatch* strip_match; // MC-information for a strip
341
342 // Loading digits ---------------------------------------------------------------
343 Int_t AddedDigits = 0;
344 Int_t AddedStripDigitMatches = 0;
345
346 for (Int_t idigit = 0; idigit < fBmnGemStripDigitsArray->GetEntriesFast(); idigit++) {
347 digit = (BmnGemStripDigit*)fBmnGemStripDigitsArray->At(idigit);
348 if (!digit->IsGoodDigit())
349 continue;
350
351 if (digit->GetStripSignal() < fSignalLow || digit->GetStripSignal() > fSignalUp)
352 continue;
353
354 BmnGemStripStation* station = StationSet->GetGemStation(digit->GetStation());
355 BmnGemStripModule* module = station->GetModule(digit->GetModule());
356
357 if (module->SetStripSignalInLayer(digit->GetStripLayer(), digit->GetStripNumber(), digit->GetStripSignal()))
358 AddedDigits++;
359
360 // Add a MC-match to the current strip if this MC-match array exists
361 if (fBmnGemStripDigitMatchesArray) {
362 strip_match = (BmnMatch*)fBmnGemStripDigitMatchesArray->At(idigit);
363 if (module->SetStripMatchInLayer(digit->GetStripLayer(), digit->GetStripNumber(), *strip_match))
364 AddedStripDigitMatches++;
365 }
366
367 // Add a digit number match to the current strip
368 BmnMatch stripDigitNumberMatch; // digit number information for the current strip
369 stripDigitNumberMatch.AddLink(1.0, idigit);
370 module->SetStripDigitNumberMatchInLayer(digit->GetStripLayer(), digit->GetStripNumber(), stripDigitNumberMatch);
371 }
372
373 if (fVerbose > 1)
374 cout << " Processed strip digits : " << AddedDigits << "\n";
375 if (fVerbose > 1 && fBmnGemStripDigitMatchesArray)
376 cout << " Added strip digit MC-matches : " << AddedStripDigitMatches << "\n";
377 //------------------------------------------------------------------------------
378
379 // Processing digits
380 StationSet->ProcessPointsInDetector();
381
382 Int_t NCalculatedPoints = StationSet->CountNProcessedPointsInDetector();
383 if (fVerbose > 1)
384 cout << " Calculated points : " << NCalculatedPoints << "\n";
385 if (fVerbose == 1)
386 cout << "BmnGemStripHitMaker: " << NCalculatedPoints << " hits\n";
387
388 Int_t clear_matched_points_cnt = 0; // points with the only one match-index
389
390 map<Int_t, StripCluster> UniqueUpperClusters;
391 map<Int_t, StripCluster> UniqueLowerClusters;
392
393 int fieldFlag = 0; // AZ-120823
394 Double_t Bx = 0.;
395 Double_t By = 0.;
396 if (fField) {
397 if (Abs(fField->GetBy(0., 0., 0.)) > FLT_EPSILON)
398 fieldFlag = 1; // AZ-120823
399 Bx = Abs(fField->GetBx(0.0, 0.0, 135.0)); // AZ-090823 - fixed value (for more stability)
400 By = Abs(fField->GetBy(0.0, 0.0, 135.0)); // AZ-090823
401 }
402
403 for (Int_t iStation = 0; iStation < StationSet->GetNStations(); ++iStation) {
404 BmnGemStripStation* station = StationSet->GetGemStation(iStation);
405
406 for (Int_t iModule = 0; iModule < station->GetNModules(); ++iModule) {
407 BmnGemStripModule* module = station->GetModule(iModule);
408
409 for (Int_t iLayer = 0; iLayer < module->GetNStripLayers(); ++iLayer) {
410 BmnGemStripLayer layer = module->GetStripLayer(iLayer);
411
412 for (size_t iCluster = 0; iCluster < layer.GetStripClusters().size(); ++iCluster) {
413 StripCluster cluster = layer.GetStripClusters()[iCluster];
414 cluster.SetModule(iModule);
415 cluster.SetStation(iStation);
416
417 // Saving all strip clusters into branches (each cluster have its own unique id in the current
418 // event)
419 if (cluster.GetType() == UpperStripLayer) {
420 new ((*fBmnGemUpperClustersArray)[fBmnGemUpperClustersArray->GetEntriesFast()])
421 StripCluster(cluster);
422 } else {
423 new ((*fBmnGemLowerClustersArray)[fBmnGemLowerClustersArray->GetEntriesFast()])
424 StripCluster(cluster);
425 }
426 }
427 }
428
429 Int_t NIntersectionPointsInModule = module->GetNIntersectionPoints();
430
431 // Int_t NPseudoIntersectionsInModule = module->GetNPseudoIntersections();
432 Int_t NPseudoIntersectionsInModule = 0; // !!! IGNORE PSEUDO INTERSECTIONS !!!
433
434 for (Int_t iPoint = 0; iPoint < NIntersectionPointsInModule + NPseudoIntersectionsInModule; ++iPoint) {
435
436 // Double_t threshold = 2000000; // ADC
437
438 // Double_t sigL;
439 // Double_t sigU;
440 Double_t x, y, z;
441 Double_t x_err, y_err, z_err;
442
443 if (iPoint < NIntersectionPointsInModule) {
444 // sigL = module->GetIntersectionPoint_LowerLayerSripTotalSignal(iPoint);
445 // sigU = module->GetIntersectionPoint_UpperLayerSripTotalSignal(iPoint);
446 x = module->GetIntersectionPointX(iPoint);
447 y = module->GetIntersectionPointY(iPoint);
448 z = module->GetZPositionRegistered();
449 x_err = module->GetIntersectionPointXError(iPoint);
450 y_err = module->GetIntersectionPointYError(iPoint);
451 z_err = 0.0;
452 } else {
453 // sigL =
454 // module->GetPseudoIntersection_LowerLayerSripTotalSignal(iPoint -
455 // NIntersectionPointsInModule);
456 // sigU =
457 // module->GetPseudoIntersection_UpperLayerSripTotalSignal(iPoint -
458 // NIntersectionPointsInModule);
459 x = module->GetPseudoIntersectionX(iPoint - NIntersectionPointsInModule);
460 y = module->GetPseudoIntersectionY(iPoint - NIntersectionPointsInModule);
461 z = module->GetZPositionRegistered();
462 x_err = module->GetPseudoIntersectionXError(iPoint - NIntersectionPointsInModule);
463 y_err = module->GetPseudoIntersectionYError(iPoint - NIntersectionPointsInModule);
464 z_err = 0.0;
465 }
466
467 if (fIsExp && !fIsSrc)
468 z += fAlignCor[iStation][iModule][2]; // alignment implementation
469
470 // Transform hit coordinates from local coordinate system of GEM-planes to global
471 if (TransfSet) {
472 Plane3D::Point glob_point = TransfSet->ApplyTransforms(Plane3D::Point(-x, y, z), iStation, iModule);
473 x = -glob_point.X();
474 y = glob_point.Y();
475 z = glob_point.Z();
476 }
477
478 Int_t RefMCIndex = -1;
479
480 // MC-matching for the current hit (define RefMCIndex)) ---------
481 BmnMatch mc_match_hit;
482 if (iPoint < NIntersectionPointsInModule) {
483 mc_match_hit = module->GetIntersectionPointMatch(iPoint);
484 } else {
485 mc_match_hit = module->GetPseudoIntersectionMatch(iPoint - NIntersectionPointsInModule);
486 }
487
488 Int_t most_probably_index = -1;
489 Double_t max_weight = 0;
490
491 Int_t n_links = mc_match_hit.GetNofLinks();
492 if (n_links == 1)
493 clear_matched_points_cnt++;
494 for (Int_t ilink = 0; ilink < n_links; ilink++) {
495 Int_t index = mc_match_hit.GetLink(ilink).GetIndex();
496 Double_t weight = mc_match_hit.GetLink(ilink).GetWeight();
497 if (weight > max_weight) {
498 max_weight = weight;
499 most_probably_index = index;
500 }
501 }
502
503 RefMCIndex = most_probably_index;
504 //--------------------------------------------------------------
505
506 // Add hit ------------------------------------------------------
507 x *= -1; // invert to global X
508
509 if (fIsExp) {
510 x += fAlignCor[iStation][iModule][0];
511 y += fAlignCor[iStation][iModule][1];
512 }
513
514 if (fieldFlag) { // AZ-120823
515 // AZ-120823 if (Abs(fField->GetBy(0., 0., 0.)) > FLT_EPSILON) {
516 if (!fIsSrc) { // For SRC lorentz corrections included into fAlignCor
517 // AZ-090823 Double_t Bx = Abs(fField->GetBx(x, y, z));
518 // AZ-090823 Double_t By = Abs(fField->GetBy(x, y, z));
519 // Bx = Abs(fField->GetBx(x, y, z)); //AZ-130823
520 // By = Abs(fField->GetBy(x, y, z)); //AZ-130823
521 Double_t yCor =
522 fLorCor[iStation][0] + fLorCor[iStation][1] * Bx + fLorCor[iStation][2] * Bx * Bx;
523 Double_t xCor =
524 fLorCor[iStation][0] + fLorCor[iStation][1] * By + fLorCor[iStation][2] * By * By;
525 if (!fIsExp)
526 yCor += 0.014; // AZ-191223
527 if (int(TMath::Abs(By) + 0.5) == 4) {
528 // AZ-180822 - 0.4T
529 xCor += 0.0160;
530 yCor += 0.0140;
531 } else if (int(TMath::Abs(By) + 0.5) == 6) {
532 // AZ-180822 - 0.6T
533 xCor += 0.0180;
534 yCor += 0.0140;
535 } else {
536 // AZ-030823 - new field
537 // xCor += 0.0175;
538 // yCor += 0.0150;
539 }
540
541 Int_t sign = (module->GetElectronDriftDirection() == ForwardZAxisEDrift) ? +1 : -1;
542 x += xCor * sign;
543 y += yCor * sign;
544 }
545 }
546
547 new ((*fBmnGemStripHitsArray)[fBmnGemStripHitsArray->GetEntriesFast()])
548 BmnGemStripHit(0, TVector3(x, y, z), TVector3(x_err, y_err, z_err), RefMCIndex);
549
550 BmnGemStripHit* hit =
551 (BmnGemStripHit*)fBmnGemStripHitsArray->At(fBmnGemStripHitsArray->GetEntriesFast() - 1);
552 hit->SetStation(iStation);
553 hit->SetModule(iModule);
554 hit->SetIndex(fBmnGemStripHitsArray->GetEntriesFast() - 1);
555
556 if (iPoint < NIntersectionPointsInModule) {
557 hit->SetDigitNumberMatch(
558 module->GetIntersectionPointDigitNumberMatch(iPoint)); // digit number match for the hit
559 hit->SetUpperClusterIndex(module->GetUpperCluster(iPoint).GetUniqueID());
560 hit->SetLowerClusterIndex(module->GetLowerCluster(iPoint).GetUniqueID());
561 hit->SetPseudo(false);
562 } else {
563 hit->SetDigitNumberMatch(module->GetPseudoIntersectionDigitNumberMatch(
564 iPoint - NIntersectionPointsInModule)); // digit number match for the hit
565 hit->SetUpperClusterIndex(
566 module->GetUpperCluster_PseudoIntersections(iPoint - NIntersectionPointsInModule)
567 .GetUniqueID());
568 hit->SetLowerClusterIndex(
569 module->GetLowerCluster_PseudoIntersections(iPoint - NIntersectionPointsInModule)
570 .GetUniqueID());
571 hit->SetPseudo(true);
572 }
573 /*
574 StripCluster ucls = module->GetUpperCluster(iPoint);
575 StripCluster lcls = module->GetLowerCluster(iPoint);
576 ucls.SetModule(iModule);
577 lcls.SetModule(iModule);
578 ucls.SetStation(iStation);
579 lcls.SetStation(iStation);
580 UniqueUpperClusters[ucls.GetUniqueID()] = ucls;
581 UniqueLowerClusters[lcls.GetUniqueID()] = lcls;
582 hit->SetUpperClusterIndex(ucls.GetUniqueID());
583 hit->SetLowerClusterIndex(lcls.GetUniqueID());
584 */
585 if (fHitMatching) {
586 BmnMatch digiMatch;
587 Int_t idx0, idx1;
588 if (iPoint < NIntersectionPointsInModule) {
589 digiMatch = module->GetIntersectionPointDigitNumberMatch(iPoint);
590 idx0 = digiMatch.GetLink(0).GetIndex();
591 idx1 = digiMatch.GetLink(1).GetIndex();
592 } else {
593 digiMatch = module->GetPseudoIntersectionDigitNumberMatch(iPoint - NIntersectionPointsInModule);
594 idx0 = digiMatch.GetLink(0).GetIndex();
595 idx1 = digiMatch.GetLink(0).GetIndex();
596 }
597 BmnMatch* digiMatch0 = (BmnMatch*)fBmnGemStripDigitMatchesArray->At(idx0);
598 BmnMatch* digiMatch1 = (BmnMatch*)fBmnGemStripDigitMatchesArray->At(idx1);
599
600 Bool_t hitOk = kFALSE;
601 for (Int_t ilink = 0; ilink < digiMatch0->GetNofLinks(); ilink++) {
602 Int_t iindex = digiMatch0->GetLink(ilink).GetIndex();
603 for (Int_t jlink = 0; jlink < digiMatch1->GetNofLinks(); jlink++) {
604 Int_t jindex = digiMatch1->GetLink(jlink).GetIndex();
605 if (iindex == jindex) {
606 hitOk = kTRUE;
607 break;
608 }
609 }
610 if (hitOk)
611 break;
612 }
613
614 hit->SetType(hitOk);
615 if (!hitOk)
616 hit->SetRefIndex(-1);
617 //--------------------------------------------------------------
618
619 // hit MC-matching ----------------------------------------------
620 FairRootManager::Instance()->SetUseFairLinks(kTRUE);
621 BmnMatch hitMatch;
622 if (iPoint < NIntersectionPointsInModule) {
623 hitMatch = module->GetIntersectionPointMatch(iPoint);
624 } else {
625 hitMatch = module->GetPseudoIntersectionMatch(iPoint - NIntersectionPointsInModule);
626 }
627
628 for (BmnLink lnk : hitMatch.GetLinks())
629 hit->AddLink(FairLink(-1, lnk.GetIndex(), lnk.GetWeight()));
630 FairRootManager::Instance()->SetUseFairLinks(kFALSE);
631 }
632 //--------------------------------------------------------------
633 }
634 }
635 }
636 /*
637 for (auto it : UniqueUpperClusters) {
638 for (Int_t i = 0; i < fBmnGemStripHitsArray->GetEntriesFast(); i++) {
639 BmnGemStripHit* hit = (BmnGemStripHit*)fBmnGemStripHitsArray->At(i);
640 if (hit->GetUpperClusterIndex() != it.first) continue;
641 hit->SetUpperClusterIndex(fBmnGemUpperClustersArray->GetEntriesFast());
642 }
643 it.second.SetUniqueID(fBmnGemUpperClustersArray->GetEntriesFast());
644 new ((*fBmnGemUpperClustersArray)[fBmnGemUpperClustersArray->GetEntriesFast()]) StripCluster(it.second);
645 }
646 for (auto it : UniqueLowerClusters) {
647 for (Int_t i = 0; i < fBmnGemStripHitsArray->GetEntriesFast(); i++) {
648 BmnGemStripHit* hit = (BmnGemStripHit*)fBmnGemStripHitsArray->At(i);
649 if (hit->GetLowerClusterIndex() != it.first) continue;
650 hit->SetLowerClusterIndex(fBmnGemLowerClustersArray->GetEntriesFast());
651 }
652 it.second.SetUniqueID(fBmnGemLowerClustersArray->GetEntriesFast());
653 new ((*fBmnGemLowerClustersArray)[fBmnGemLowerClustersArray->GetEntriesFast()]) StripCluster(it.second);
654 }
655 */
656 if (fVerbose > 1)
657 cout << " N clear matches with MC-points = " << clear_matched_points_cnt << "\n";
658 //------------------------------------------------------------------------------
659 StationSet->Reset();
660}
661
662void BmnGemStripHitMaker::OnlineWrite(const std::unique_ptr<TTree>& resultTree)
663{
664 if (!IsActive())
665 return;
666
667 resultTree->Branch(fOutputHitsBranchName, &fBmnGemStripHitsArray);
668 resultTree->Branch("BmnGemUpperCluster", &fBmnGemUpperClustersArray);
669 resultTree->Branch("BmnGemLowerCluster", &fBmnGemLowerClustersArray);
670 resultTree->Fill();
671}
672
673void BmnGemStripHitMaker::Finish()
674{
675 printf("Work time of BmnGemStripHitMaker: %4.2f sec.\n", workTime);
676}
UpperStripLayer
@ UpperStripLayer
Definition BmnStripData.h:13
BmnGemStripHitMaker::Init
virtual InitStatus Init()
Definition BmnGemStripHitMaker.cxx:256
BmnGemStripHitMaker::Exec
virtual void Exec(Option_t *opt)
Definition BmnGemStripHitMaker.cxx:303
BmnGemStripHitMaker::OnlineWrite
virtual void OnlineWrite(const std::unique_ptr< TTree > &resultTree)
Write task resul to tree.
Definition BmnGemStripHitMaker.cxx:662
BmnGemStripLayer::SetUpperUniqueID
static void SetUpperUniqueID(Int_t id)
Definition BmnGemStripLayer.h:112
BmnGemStripLayer::GetStripClusters
vector< StripCluster > GetStripClusters()
Definition BmnGemStripLayer.h:110
BmnGemStripLayer::SetLowerUniqueID
static void SetLowerUniqueID(Int_t id)
Definition BmnGemStripLayer.h:111
BmnGemStripStationSet::GetGemStation
BmnGemStripStation * GetGemStation(Int_t station_num)
Definition BmnGemStripStationSet.cxx:150
BmnGemStripTransform::LoadFromXMLFile
Bool_t LoadFromXMLFile(TString xml_config_file)
Definition BmnGemStripTransform.cxx:91
BmnGemStripTransform::ApplyTransforms
Plane3D::Point ApplyTransforms(Plane3D::Point point, Int_t station, Int_t module)
Definition BmnGemStripTransform.cxx:223
BmnHit::SetModule
void SetModule(Int_t mod)
Definition BmnHit.h:73
BmnHit::SetType
void SetType(Int_t type)
Definition BmnHit.h:81
BmnHit::SetUpperClusterIndex
void SetUpperClusterIndex(Int_t idx)
Definition BmnHit.h:133
BmnHit::SetDigitNumberMatch
void SetDigitNumberMatch(BmnMatch match)
Definition BmnHit.h:121
BmnHit::SetIndex
void SetIndex(Int_t id)
Definition BmnHit.h:57
BmnHit::SetPseudo
void SetPseudo(Bool_t pseudo)
Definition BmnHit.h:141
BmnHit::SetStation
void SetStation(Short_t st)
Definition BmnHit.h:69
BmnHit::SetLowerClusterIndex
void SetLowerClusterIndex(Int_t idx)
Definition BmnHit.h:137
BmnMatch::GetLink
const BmnLink & GetLink(Int_t i) const
Definition BmnMatch.h:37
BmnMatch::AddLink
void AddLink(const BmnMatch &match)
Definition BmnMatch.cxx:43
BmnMatch::GetLinks
const vector< BmnLink > & GetLinks() const
Definition BmnMatch.h:38
BmnMatch::GetNofLinks
Int_t GetNofLinks() const
Definition BmnMatch.h:40
BmnStripDigit::IsGoodDigit
Bool_t IsGoodDigit()
Definition BmnStripDigit.h:32
BmnStripDigit::GetStripNumber
Int_t GetStripNumber()
Definition BmnStripDigit.h:29
BmnStripDigit::GetStripLayer
Int_t GetStripLayer()
Definition BmnStripDigit.h:28
BmnStripDigit::GetStation
Int_t GetStation()
Definition BmnStripDigit.h:26
BmnStripDigit::GetStripSignal
Double_t GetStripSignal()
Definition BmnStripDigit.h:30
StripCluster::SetStation
void SetStation(Int_t st)
Definition BmnStripData.h:215
StripCluster::GetType
Int_t GetType()
Definition BmnStripData.h:211
StripCluster::SetModule
void SetModule(Int_t m)
Definition BmnStripData.h:217
UniDetectorParameter::GetValue
int GetValue(vector< UniValue * > &parameter_value)
get value of detector parameter presented by an array
Definition UniDetectorParameter.cxx:1107
UniDetectorParameter::GetDetectorParameter
static UniDetectorParameter * GetDetectorParameter(int value_id)
get detector parameter from the database
Definition UniDetectorParameter.cxx:163