BmnRoot 26.10.0-18-g4e3414f9
Software framework of the BM@N experiment
Loading...
Searching...
No Matches
CbmStsIdealFindHits.cxx
Go to the documentation of this file.
1//* $Id: */
2
3// -------------------------------------------------------------------------
4// ----- CbmStsIdealFindHits source file -----
5// ----- Created 11/09/06 by V. Friese -----
6// -------------------------------------------------------------------------
7
8#include "CbmStsIdealFindHits.h"
9
10#include "CbmGeoStsPar.h"
11#include "CbmStsDigi.h"
12#include "CbmStsDigiMatch.h" //AZ
13#include "CbmStsDigiPar.h"
14#include "CbmStsDigiScheme.h"
15#include "CbmStsHit.h"
16#include "CbmStsPoint.h" //AZ
17#include "CbmStsSector.h"
18#include "CbmStsStation.h"
19#include "FairRootManager.h"
20#include "FairRunAna.h"
21#include "FairRuntimeDb.h"
22#include "TClonesArray.h"
23#include "TMath.h"
24
25#include <iomanip>
26
27using std::cerr;
28using std::cout;
29using std::endl;
30using std::fixed;
31using std::flush;
32using std::left;
33using std::map;
34using std::right;
35using std::set;
36using std::setprecision;
37using std::setw;
38
39// ----- Default constructor ------------------------------------------
40CbmStsIdealFindHits::CbmStsIdealFindHits()
41 : FairTask("STS Ideal Hit Finder", 1)
42 , fGeoPar(NULL)
43 , fDigiPar(NULL)
44 , fDigis(NULL)
45 , fDigiMatch(NULL)
46 , // AZ
47 fPoints(NULL)
48 , // AZ
49 fHits(NULL)
50 , fDigiScheme(new CbmStsDigiScheme())
51 , fNStations(0)
52 , fNEvents(0)
53 , fDigiMapF()
54 , fDigiMapB()
55 , fTimer()
56 , fTime1(0.)
57{}
58// -------------------------------------------------------------------------
59
60// ----- Standard constructor ------------------------------------------
61CbmStsIdealFindHits::CbmStsIdealFindHits(Int_t iVerbose)
62 : FairTask("STSIdealFindHits", iVerbose)
63 , fGeoPar(NULL)
64 , fDigiPar(NULL)
65 , fDigis(NULL)
66 , fDigiMatch(NULL)
67 , // AZ
68 fPoints(NULL)
69 , // AZ
70 fHits(NULL)
71 , fDigiScheme(new CbmStsDigiScheme())
72 , fNStations(0)
73 , fNEvents(0)
74 , fDigiMapF()
75 , fDigiMapB()
76 , fTimer()
77 , fTime1(0.)
78{}
79// -------------------------------------------------------------------------
80
81// ----- Constructor with name -----------------------------------------
82CbmStsIdealFindHits::CbmStsIdealFindHits(const char* name, Int_t iVerbose)
83 : FairTask(name, iVerbose)
84 , fGeoPar(NULL)
85 , fDigiPar(NULL)
86 , fDigis(NULL)
87 , fDigiMatch(NULL)
88 , // AZ
89 fPoints(NULL)
90 , // AZ
91 fHits(NULL)
92 , fDigiScheme(new CbmStsDigiScheme())
93 , fNStations(0)
94 , fNEvents(0)
95 , fDigiMapF()
96 , fDigiMapB()
97 , fTimer()
98 , fTime1(0.)
99{}
100// -------------------------------------------------------------------------
101
102// ----- Destructor ----------------------------------------------------
103CbmStsIdealFindHits::~CbmStsIdealFindHits()
104{
105 if (fHits) {
106 fHits->Delete();
107 delete fHits;
108 }
109}
110// -------------------------------------------------------------------------
111
112// ----- Public method Exec --------------------------------------------
113void CbmStsIdealFindHits::Exec(Option_t* opt)
114{
115 fTimer.Start();
116
117 Bool_t warn = kFALSE;
118
119 // Check for digi scheme
120 if (!fDigiScheme) {
121 cerr << "-E- " << fName << "::Exec: No digi scheme!" << endl;
122 return;
123 }
124
125 // Clear output array
126 // fHits->Clear();
127 fHits->Delete();
128
129 // Sort STS digis with respect to sectors
130 SortDigis();
131
132 // Find hits in sectors
133 Int_t nDigisF = 0;
134 Int_t nDigisB = 0;
135 Int_t nHits = 0;
136 Int_t nStations = fDigiScheme->GetNStations();
137 CbmStsStation* station = NULL;
138 for (Int_t iStation = 0; iStation < nStations; iStation++) {
139 station = fDigiScheme->GetStation(iStation);
140
141 Int_t nDigisFInStation = 0;
142 Int_t nDigisBInStation = 0;
143 Int_t nHitsInStation = 0;
144 Int_t nSectors = station->GetNSectors();
145 for (Int_t iSector = 0; iSector < nSectors; iSector++) {
146 CbmStsSector* sector = station->GetSector(iSector);
147 set<Int_t> fSet, bSet;
148 if (fDigiMapF.find(sector) == fDigiMapF.end()) {
149 cout << "-E- " << fName << "::Exec: sector " << sector->GetSectorNr() << " of station "
150 << station->GetStationNr() << "not found in front map!" << endl;
151 warn = kTRUE;
152 continue;
153 }
154 fSet = fDigiMapF[sector];
155 if (sector->GetType() == 2 || sector->GetType() == 3) {
156 if (fDigiMapB.find(sector) == fDigiMapB.end()) {
157 cout << "-E- " << fName << "::Exec: sector " << sector->GetSectorNr() << " of station "
158 << station->GetStationNr() << "not found in back map!" << endl;
159 warn = kTRUE;
160 continue;
161 }
162 }
163 bSet = fDigiMapB[sector];
164 Int_t nDigisFInSector = fSet.size();
165 Int_t nDigisBInSector = bSet.size();
166 Int_t nHitsInSector = FindHits(station, sector, fSet, bSet);
167 if (fVerbose > 2)
168 cout << "Sector " << sector->GetSectorNr() << ", Digis front " << nDigisFInSector << ", Digis Back "
169 << nDigisBInSector << ", Hits " << nHitsInSector << endl;
170 nHitsInStation += nHitsInSector;
171 nDigisFInStation += nDigisFInSector;
172 nDigisBInStation += nDigisBInSector;
173 } // Sector loop
174
175 if (fVerbose > 1)
176 cout << "Total for station " << station->GetStationNr() << ": Digis front " << nDigisFInStation
177 << ", digis back " << nDigisBInStation << ", hits " << nHitsInStation << endl;
178 nDigisB += nDigisBInStation;
179 nDigisF += nDigisFInStation;
180 nHits += nHitsInStation;
181
182 } // Station loop
183
184 fTimer.Stop();
185 if (fVerbose > 1) {
186 cout << endl;
187 cout << "-I- " << fName << ":Event summary" << endl;
188 cout << " Active channels front side: " << nDigisF << endl;
189 cout << " Active channels back side : " << nDigisB << endl;
190 cout << " Hits created : " << nHits << endl;
191 cout << " Real time : " << fTimer.RealTime() << endl;
192 }
193 if (fVerbose == 1) {
194 if (warn)
195 cout << "- ";
196 else
197 cout << "+ ";
198 cout << setw(15) << left << fName << ": " << setprecision(4) << setw(8) << fixed << right << fTimer.RealTime()
199 << " s, digis " << nDigisF << " / " << nDigisB << ", hits: " << nHits << endl;
200 }
201}
202// -------------------------------------------------------------------------
203
204// ----- Private method SetParContainers -------------------------------
205void CbmStsIdealFindHits::SetParContainers()
206{
207
208 // Get run and runtime database
209 FairRunAna* run = FairRunAna::Instance();
210 if (!run)
211 Fatal("SetParContainers", "No analysis run");
212
213 FairRuntimeDb* db = run->GetRuntimeDb();
214 if (!db)
215 Fatal("SetParContainers", "No runtime database");
216
217 // Get STS geometry parameter container
218 fGeoPar = (CbmGeoStsPar*)db->getContainer("CbmGeoStsPar");
219
220 // Get STS digitisation parameter container
221 fDigiPar = (CbmStsDigiPar*)db->getContainer("CbmStsDigiPar");
222}
223// -------------------------------------------------------------------------
224
225// ----- Private method Init -------------------------------------------
226InitStatus CbmStsIdealFindHits::Init()
227{
228
229 // Get input array
230 FairRootManager* ioman = FairRootManager::Instance();
231 if (!ioman)
232 Fatal("Init", "No FairRootManager");
233 fDigis = (TClonesArray*)ioman->GetObject("StsDigi");
234 fDigiMatch = (TClonesArray*)ioman->GetObject("StsDigiMatch"); // AZ
235 fPoints = (TClonesArray*)ioman->GetObject("StsPoint"); // AZ
236
237 // Register output array
238 fHits = new TClonesArray("CbmStsHit", 1000);
239 ioman->Register("StsHit", "Hit in STS", fHits, kTRUE);
240
241 // Build digitisation scheme
242 Bool_t success = fDigiScheme->Init(fGeoPar, fDigiPar);
243 if (!success)
244 return kERROR;
245
246 // Create sectorwise digi sets
247 MakeSets();
248
249 // Control output
250 if (fVerbose == 1 || fVerbose == 2)
251 fDigiScheme->Print(kFALSE);
252 else if (fVerbose > 2)
253 fDigiScheme->Print(kTRUE);
254 cout << "-I- " << fName << "::Init: "
255 << "STS digitisation scheme succesfully initialised" << endl;
256 cout << " Stations: " << fDigiScheme->GetNStations() << ", Sectors: " << fDigiScheme->GetNSectors()
257 << ", Channels: " << fDigiScheme->GetNChannels() << endl;
258
259 return kSUCCESS;
260}
261// -------------------------------------------------------------------------
262
263// ----- Private method ReInit -----------------------------------------
264InitStatus CbmStsIdealFindHits::ReInit()
265{
266
267 // Clear digitisation scheme
268 fDigiScheme->Clear();
269
270 // Build new digitisation scheme
271 Bool_t success = fDigiScheme->Init(fGeoPar, fDigiPar);
272 if (!success)
273 return kERROR;
274
275 // Create sectorwise digi sets
276 MakeSets();
277
278 // Control output
279 if (fVerbose == 1 || fVerbose == 2)
280 fDigiScheme->Print(kFALSE);
281 else if (fVerbose > 2)
282 fDigiScheme->Print(kTRUE);
283 cout << "-I- " << fName << "::Init: "
284 << "STS digitisation scheme succesfully reinitialised" << endl;
285 cout << " Stations: " << fDigiScheme->GetNStations() << ", Sectors: " << fDigiScheme->GetNSectors()
286 << ", Channels: " << fDigiScheme->GetNChannels() << endl;
287
288 return kSUCCESS;
289}
290// -------------------------------------------------------------------------
291
292// ----- Private method MakeSets ---------------------------------------
293void CbmStsIdealFindHits::MakeSets()
294{
295
296 fDigiMapF.clear();
297 fDigiMapB.clear();
298 Int_t nStations = fDigiScheme->GetNStations();
299 for (Int_t iStation = 0; iStation < nStations; iStation++) {
300 CbmStsStation* station = fDigiScheme->GetStation(iStation);
301 Int_t nSectors = station->GetNSectors();
302 for (Int_t iSector = 0; iSector < nSectors; iSector++) {
303 CbmStsSector* sector = station->GetSector(iSector);
304 set<Int_t> a;
305 fDigiMapF[sector] = a;
306 if (sector->GetType() == 2 || sector->GetType() == 3) {
307 set<Int_t> b;
308 fDigiMapB[sector] = b;
309 }
310 }
311 }
312}
313// -------------------------------------------------------------------------
314
315// ----- Private method SortDigis --------------------------------------
316void CbmStsIdealFindHits::SortDigis()
317{
318
319 // Check input array
320 if (!fDigis) {
321 cout << "-E- " << fName << "::SortDigis: No input array!" << endl;
322 return;
323 }
324
325 // Clear sector digi sets
326 map<CbmStsSector*, set<Int_t>>::iterator mapIt;
327 for (mapIt = fDigiMapF.begin(); mapIt != fDigiMapF.end(); mapIt++)
328 ((*mapIt).second).clear();
329 for (mapIt = fDigiMapB.begin(); mapIt != fDigiMapB.end(); mapIt++)
330 ((*mapIt).second).clear();
331
332 // Fill digis into sets
333 CbmStsDigi* digi = NULL;
334 CbmStsSector* sector = NULL;
335 Int_t stationNr = -1;
336 Int_t sectorNr = -1;
337 Int_t iSide = -1;
338 Int_t nDigis = fDigis->GetEntriesFast();
339 for (Int_t iDigi = 0; iDigi < nDigis; iDigi++) {
340 digi = (CbmStsDigi*)fDigis->At(iDigi);
341 stationNr = digi->GetStationNr();
342 sectorNr = digi->GetSectorNr();
343 iSide = digi->GetSide();
344 sector = fDigiScheme->GetSector(stationNr, sectorNr);
345 if (iSide == 0) {
346 if (fDigiMapF.find(sector) == fDigiMapF.end()) {
347 cerr << "-E- " << fName << "::SortDigits:: sector " << sectorNr << " of station " << stationNr
348 << " not found in digi scheme (F)!" << endl;
349 continue;
350 }
351 fDigiMapF[sector].insert(iDigi);
352 } else if (iSide == 1) {
353 if (fDigiMapB.find(sector) == fDigiMapB.end()) {
354 cerr << "-E- " << fName << "::SortDigits:: sector " << sectorNr << " of station " << stationNr
355 << " not found in digi scheme (B)!" << endl;
356 continue;
357 }
358 fDigiMapB[sector].insert(iDigi);
359 }
360 }
361}
362// -------------------------------------------------------------------------
363
364// ----- Private method FindHits ---------------------------------------
365Int_t CbmStsIdealFindHits::FindHits(CbmStsStation* station, CbmStsSector* sector, set<Int_t>& fSet, set<Int_t>& bSet)
366{
367
368 // Counter
369 Int_t nNew = 0;
370
371 // Get sector parameters
372 Int_t detId = sector->GetDetectorId();
373 Int_t iType = sector->GetType();
374
375 Double_t rot = sector->GetRotation();
376 Double_t dx = sector->GetDx();
377 Double_t dy = sector->GetDy();
378 Double_t stereoF = sector->GetStereoF();
379 Double_t stereoB = sector->GetStereoB();
380
381 // Double_t z = station->GetZ();
382 Int_t stationNr = station->GetStationNr();
383 Int_t sectorNr = sector->GetSectorNr();
384
385 // Some auxiliary values
386 Double_t sinrot = TMath::Sin(rot);
387 Double_t cosrot = TMath::Cos(rot);
388 Double_t tanstrB = TMath::Tan(stereoB);
389 Double_t tanstrF = TMath::Tan(stereoF);
390
391 // Calculate error matrix in sector system
392 Double_t vX, vY, vXY;
393 // AZ if ( iType == 1 ) {
394 if (iType == 1
395 || TMath::Abs(stereoF) < 0.01
396 && (TMath::Abs(stereoB * TMath::RadToDeg() - 90) < 0.01
397 || TMath::Abs(stereoB * TMath::RadToDeg() + 90) < 0.01))
398 { // AZ
399 vX = dx / TMath::Sqrt(12.);
400 vY = dy / TMath::Sqrt(12.);
401 vXY = 0.;
402 } else if (iType == 2 || iType == 3) {
403 if (stereoF == 0.) {
404 vX = dx / TMath::Sqrt(12.);
405 vY = dx / TMath::Sqrt(6.) / TMath::Abs(tanstrB);
406 vXY = -1. * dx * dx / 12. / tanstrB;
407 } else {
408 vX = dx / TMath::Sqrt(24.);
409 vY = dx / TMath::Sqrt(24.) / TMath::Abs(tanstrB);
410 vXY = 0.;
411 }
412 } else {
413 cerr << "-E- " << fName << "::FindHits: Illegal sector type " << iType << endl;
414 return 0;
415 }
416
417 // Transform variances into global c.s.
418 Double_t wX = vX * vX * cosrot * cosrot - 2. * vXY * cosrot * sinrot + vY * vY * sinrot * sinrot;
419 Double_t wY = vX * vX * sinrot * sinrot + 2. * vXY * cosrot * sinrot + vY * vY * cosrot * cosrot;
420 Double_t wXY = (vX * vX - vY * vY) * cosrot * sinrot + vXY * (cosrot * cosrot - sinrot * sinrot);
421 Double_t sigmaX = TMath::Sqrt(wX);
422 Double_t sigmaY = TMath::Sqrt(wY);
423
424 // Now perform the loop over active channels
425 set<Int_t>::iterator it1;
426 set<Int_t>::iterator it2;
427
428 // ----- Type 1 : Pixel sector ---------------------------------------
429 if (iType == 1) {
430 Fatal("FindHits", "Sorry, not implemented yet");
431 } // Pixel sensor
432 // ---------------------------------------------------------------------
433
434 // ----- Type 2: Strip sector OSU -----------------------------------
435 else if (iType == 2)
436 {
437 Fatal("FindHits", "Sorry, not implemented yet");
438 } // Strip OSU
439 // ---------------------------------------------------------------------
440
441 // ----- Type 3: Strip sector GSI -----------------------------------
442 else if (iType == 3)
443 {
444 Int_t iDigiF = -1;
445 Int_t iDigiB = -1;
446 Int_t iChanF = -1;
447 Int_t iChanB = -1;
448 Int_t nHits = fHits->GetEntriesFast();
449 Double_t xHit;
450 Double_t yHit;
451 Double_t zHit;
452 TVector3 pos, dpos;
453 CbmStsDigi* digiF = NULL;
454 CbmStsDigi* digiB = NULL;
455 CbmStsDigiMatch *digiFM = NULL, *digiBM = NULL; // AZ
456
457 for (it1 = fSet.begin(); it1 != fSet.end(); it1++) {
458 iDigiF = (*it1);
459 digiF = (CbmStsDigi*)fDigis->At(iDigiF);
460 if (!digiF) {
461 cout << "-W- " << GetName() << "::FindHits: Invalid digi index " << iDigiF << " in front set of sector "
462 << sector->GetSectorNr() << ", station " << station->GetStationNr() << endl;
463 continue;
464 }
465 iChanF = digiF->GetChannelNr();
466 digiFM = (CbmStsDigiMatch*)fDigiMatch->UncheckedAt(iDigiF); // AZ
467
468 for (it2 = bSet.begin(); it2 != bSet.end(); it2++) {
469 iDigiB = (*it2);
470 digiB = (CbmStsDigi*)fDigis->At(iDigiB);
471 if (!digiB) {
472 cout << "-W- " << GetName() << "::FindHits: Invalid digi index " << iDigiB
473 << " in back set of sector " << sector->GetSectorNr() << ", station "
474 << station->GetStationNr() << endl;
475 continue;
476 }
477 // AZ - pixel simulation: intersect only digits from the same point
478 digiBM = (CbmStsDigiMatch*)fDigiMatch->UncheckedAt(iDigiB); // AZ
479 /*
480 if (TMath::Abs(TMath::Abs(stereoF-stereoB)*TMath::RadToDeg()-90) < 0.1) {
481 Bool_t ok = kFALSE;
482 for (Int_t ii = 0; ii < 3; ++ii) {
483 if (digiFM->GetRefIndex(ii) < 0) break;
484 for (Int_t jj = 0; jj < 3; ++jj) {
485 if (digiBM->GetRefIndex(jj) < 0) break;
486 if (digiFM->GetRefIndex(ii) == digiBM->GetRefIndex(jj)) {
487 ok = kTRUE;
488 // Take MCPoint coordinates
489 CbmStsPoint *point = (CbmStsPoint*) fPoints->UncheckedAt(digiFM->GetRefIndex(ii));
490 xHit = point->GetXIn();
491 yHit = point->GetYIn();
492 break;
493 }
494 }
495 }
496 if (!ok) continue;
497 }
498 */
499 // AZ
500
501 iChanB = digiB->GetChannelNr();
502
503 Int_t sensorDetId = sector->Intersect(iChanF, iChanB, xHit, yHit, zHit);
504
505 if (sensorDetId == -1)
506 continue;
507
508 pos.SetXYZ(xHit, yHit, zHit);
509 dpos.SetXYZ(sigmaX, sigmaY, 0.);
510
511 Int_t statLayer = -1;
512 for (Int_t istatL = station->GetNofZ(); istatL > 0; istatL--)
513 if (TMath::Abs(zHit - station->GetZ(istatL - 1)) < 0.00001) {
514 statLayer = istatL - 1;
515 break;
516 }
517
518 if (statLayer == -1)
519 cout << "unknown layer for hit at z = " << zHit << endl;
520
521 new ((*fHits)[nHits++])
522 CbmStsHit(sensorDetId, pos, dpos, wXY, iDigiF, iDigiB, iChanF, iChanB, statLayer);
523
524 nNew++;
525 if (fVerbose > 3)
526 cout << "New StsHit at (" << xHit << ", " << yHit << ", " << zHit << "), station " << stationNr
527 << ", sector " << sectorNr << ", channel " << iChanF << " / " << iChanB << endl;
528
529 } // back side strip loop
530 } // front side strip loop
531
532 } // strip GSI
533 // ---------------------------------------------------------------------
534
535 return nNew;
536}
537// -------------------------------------------------------------------------
538
539/*
540// ----- Private method FindHits ---------------------------------------
541Int_t CbmStsIdealFindHits::FindHits(CbmStsStation* station,
542 CbmStsSector* sector,
543 set<Int_t>& fSet, set<Int_t>& bSet) {
544
545 // Counter
546 Int_t nNew = 0;
547
548 // Get sector parameters
549 Int_t detId = sector->GetDetectorId();
550 Int_t iType = sector->GetType();
551 Double_t xc = sector->GetX0();
552 Double_t yc = sector->GetY0();
553 Double_t rot = sector->GetRotation();
554 Double_t lx = sector->GetLx();
555 Double_t ly = sector->GetLy();
556 Double_t dx = sector->GetDx();
557 Double_t dy = sector->GetDy();
558 Double_t stereoB = sector->GetStereoB();
559 Double_t z = sector->GetZ0();
560 // Double_t z = station->GetZ();
561 Int_t stationNr = station->GetStationNr();
562 Int_t sectorNr = sector->GetSectorNr();
563
564 Int_t statLayer = 0;
565 for ( Int_t isz = 1 ; isz < 11 ; isz++ )
566 if ( TMath::Abs(station->GetZ(isz)-z) < 0.001 ) {
567 statLayer = isz;
568 break;
569 }
570
571 // Some auxiliary values
572 Double_t sinrot = TMath::Sin(rot);
573 Double_t cosrot = TMath::Cos(rot);
574 Double_t tanstr = TMath::Tan(stereoB);
575
576 // Calculate error matrix in sector system
577 Double_t vX, vY, vXY;
578 if ( iType == 1 ) {
579 vX = dx / TMath::Sqrt(12.);
580 vY = dy / TMath::Sqrt(12.);
581 vXY = 0.;
582 }
583 else if ( iType == 2 || iType == 3 ) {
584 vX = dx / TMath::Sqrt(12.);
585 vY = dx / TMath::Sqrt(6.) / TMath::Abs(tanstr);
586 vXY = -1. * dx * dx / 12. / tanstr;
587 }
588 else {
589 cerr << "-E- " << fName << "::FindHits: Illegal sector type "
590 << iType << endl;
591 return 0;
592 }
593
594 // Transform variances into global c.s.
595 Double_t wX = vX * vX * cosrot * cosrot
596 - 2. * vXY * cosrot * sinrot
597 + vY * vY * sinrot * sinrot;
598 Double_t wY = vX * vX * sinrot * sinrot
599 + 2. * vXY * cosrot * sinrot
600 + vY * vY * cosrot * cosrot;
601 Double_t wXY = (vX*vX - vY*vY) * cosrot * sinrot
602 + vXY * ( cosrot*cosrot - sinrot*sinrot );
603 Double_t sigmaX = TMath::Sqrt(wX);
604 Double_t sigmaY = TMath::Sqrt(wY);
605
606 // Now perform the loop over active channels
607 set<Int_t>::iterator it1;
608 set<Int_t>::iterator it2;
609
610 // ----- Type 1 : Pixel sector ---------------------------------------
611 if ( iType == 1 ) {
612 Int_t nColumns = Int_t(TMath::Ceil ( lx / dx ));
613 Int_t iDigi = -1;
614 Int_t iChan = -1;
615 Int_t iCol = -1;
616 Int_t iRow = -1;
617 Int_t nHits = fHits->GetEntriesFast();
618 Double_t xint, yint;
619 Double_t x, y;
620 TVector3 pos, dpos;
621 CbmStsDigi* digi = NULL;
622 for (it1=fSet.begin(); it1!=fSet.end(); it1++) {
623 iDigi = (*it1);
624 digi = (CbmStsDigi*) fDigis->At(iDigi);
625 if ( ! digi ) {
626 cout << "-W- " << GetName() << "::FindHits: Invalid digi index "
627 << iDigi << " in front set of sector "
628 << sector->GetSectorNr() << ", station "
629 << station->GetStationNr() << endl;
630 continue;
631 }
632 iChan = digi->GetChannelNr();
633 iRow = Int_t( iChan / nColumns );
634 iCol = iChan - iRow * nColumns;
635 xint = ( Double_t(iCol) + 0.5 ) * dx;
636 yint = ( Double_t(iRow) + 0.5 ) * dy;
637
638 // Translation to centre of sector
639 xint = xint - lx/2.;
640 yint = yint - ly/2.;
641
642 // Rotation around sector centre
643 x = xint * cosrot - yint * sinrot;
644 y = xint * sinrot + yint * cosrot;
645
646 // Translation into global c.s.
647 x = x + xc;
648 y = y + yc;
649
650 // Make new hit
651 pos.SetXYZ(x, y, z);
652 dpos.SetXYZ(sigmaX, sigmaY, 0.);
653 new ((*fHits)[nHits++]) CbmStsHit(detId, pos, dpos, wXY, iDigi, -1);
654 nNew++;
655 if ( fVerbose > 3 ) cout << "New StsHit at (" << x << ", " << y
656 << ", " << z << "), station " << stationNr
657 << ", sector " << sectorNr << ", channel "
658 << iChan << endl;
659 } // Loop over digi set
660
661 } // Pixel sensor
662 // ---------------------------------------------------------------------
663
664
665
666 // ----- Type 2: Strip sector OSU -----------------------------------
667 else if ( iType == 2 ) {
668 Int_t iDigiF = -1;
669 Int_t iDigiB = -1;
670 Int_t iChanF = -1;
671 Int_t iChanB = -1;
672 Int_t nHits = fHits->GetEntriesFast();
673 Double_t x0 = 0.;
674 Double_t xint, yint;
675 Double_t xtemp, ytemp;
676 Double_t x, y;
677 TVector3 pos, dpos;
678 CbmStsDigi* digiF = NULL;
679 CbmStsDigi* digiB = NULL;
680 for (it1=fSet.begin(); it1!=fSet.end(); it1++) {
681 iDigiF = (*it1);
682 digiF = (CbmStsDigi*) fDigis->At(iDigiF);
683 if ( ! digiF ) {
684 cout << "-W- " << GetName() << "::FindHits: Invalid digi index "
685 << iDigiF << " in front set of sector "
686 << sector->GetSectorNr() << ", station "
687 << station->GetStationNr() << endl;
688 continue;
689 }
690 iChanF = digiF->GetChannelNr();
691 xint = ( Double_t(iChanF) + 0.5 ) * dx;
692 for (it2=bSet.begin(); it2!=bSet.end(); it2++ ) {
693 iDigiB = (*it2);
694 digiB = (CbmStsDigi*) fDigis->At(iDigiB);
695 if ( ! digiB ) {
696 cout << "-W- " << GetName() << "::FindHits: Invalid digi index "
697 << iDigiB << " in back set of sector "
698 << sector->GetSectorNr() << ", station "
699 << station->GetStationNr() << endl;
700 continue;
701 }
702 iChanB = digiB->GetChannelNr();;
703 x0 = ( Double_t(iChanB) + 0.5 ) * dx;
704 yint = (x0 - xint) / tanstr;
705 if ( ! ( yint>0. && yint<ly) ) continue;
706
707 // Translation to centre of sector
708 xtemp = xint - lx/2.;
709 ytemp = yint - ly/2.;
710
711 // Rotation around sector centre
712 x = xtemp * cosrot - ytemp * sinrot;
713 y = xtemp * sinrot + ytemp * cosrot;
714
715 // Translation into global c.s.
716 x = x + xc;
717 y = y + yc;
718
719 // Make new hit
720 pos.SetXYZ(x, y, z);
721 dpos.SetXYZ(sigmaX, sigmaY, 0.);
722 new ((*fHits)[nHits++]) CbmStsHit(detId, pos, dpos, wXY,
723 iDigiF, iDigiB);
724 nNew++;
725 if ( fVerbose > 3 ) cout << "New StsHit at (" << x << ", " << y
726 << ", " << z << "), station " << stationNr
727 << ", sector " << sectorNr << ", channel "
728 << iChanF << " / " << iChanB << endl;
729 } // Back side set loop
730 } // Front side set sloop
731
732 } // Strip OSU
733 // ---------------------------------------------------------------------
734
735
736
737 // ----- Type 3: Strip sector GSI -----------------------------------
738 else if (iType == 3 ) {
739 Int_t iDigiF = -1;
740 Int_t iDigiB = -1;
741 Int_t iChanF = -1;
742 Int_t iChanB = -1;
743 Int_t nHits = fHits->GetEntriesFast();
744 Int_t nStripMax = ( stereoB<0. ? 0 : Int_t(ly*tanstr/lx)+1 ); // max. number of strips
745 Int_t nStripBeg = ( stereoB>0. ? 0 : -Int_t(ly*tanstr/lx)-1 );
746 Double_t x0 = 0.;
747 Double_t xint, yint;
748 Double_t xtemp, ytemp;
749 Double_t x,y;
750 TVector3 pos, dpos;
751 CbmStsDigi* digiF = NULL;
752 CbmStsDigi* digiB = NULL;
753 for (it1=fSet.begin(); it1!=fSet.end(); it1++) {
754 iDigiF = (*it1);
755 digiF = (CbmStsDigi*) fDigis->At(iDigiF);
756 if ( ! digiF ) {
757 cout << "-W- " << GetName() << "::FindHits: Invalid digi index "
758 << iDigiF << " in front set of sector "
759 << sector->GetSectorNr() << ", station "
760 << station->GetStationNr() << endl;
761 continue;
762 }
763 iChanF = digiF->GetChannelNr();
764 xint = ( Double_t(iChanF) + 0.5 ) * dx;
765 for (it2=bSet.begin(); it2!=bSet.end(); it2++ ) {
766 iDigiB = (*it2);
767 digiB = (CbmStsDigi*) fDigis->At(iDigiB);
768 if ( ! digiB ) {
769 cout << "-W- " << GetName() << "::FindHits: Invalid digi index "
770 << iDigiB << " in back set of sector "
771 << sector->GetSectorNr() << ", station "
772 << station->GetStationNr() << endl;
773 continue;
774 }
775 iChanB = digiB->GetChannelNr();;
776 x0 = ( Double_t(iChanB) + 0.5 ) * dx;
777 for (Int_t iStrip=nStripBeg; iStrip<=nStripMax; iStrip++) {
778 yint = ( x0 - xint + Double_t(iStrip) * lx ) / tanstr;
779
780 if ( ! ( yint>0. && yint<ly ) ) continue;
781
782 Int_t discreteNetYPos = (Int_t)(0.5+yint*tanstr/dx);
783
784 // Translation to centre of sector
785 xtemp = xint - lx/2.;
786 ytemp = yint - ly/2.;
787
788 // Rotation around sector centre
789 x = xtemp * cosrot - ytemp * sinrot;
790 y = xtemp * sinrot + ytemp * cosrot;
791
792 // Translation into global c.s.
793 x = x + xc;
794 y = y + yc;
795
796 // Make new hit
797 pos.SetXYZ(x, y, z);
798 dpos.SetXYZ(sigmaX, sigmaY, 0.);
799 new ((*fHits)[nHits++]) CbmStsHit(detId, pos, dpos, wXY,
800 iDigiF, iDigiB, iChanF, discreteNetYPos, statLayer);
801 nNew++;
802 if ( fVerbose > 3 ) cout << "New StsHit at (" << x << ", " << y
803 << ", " << z << "), station "
804 << stationNr << ", sector " << sectorNr
805 << ", channel " << iChanF << " / "
806 << iChanB << ", segment " << iStrip
807 << endl;
808
809 } // strip segment loop
810 } // back side strip loop
811 } // front side strip loop
812
813 } // strip GSI
814 // ---------------------------------------------------------------------
815
816
817 return nNew;
818}
819// -------------------------------------------------------------------------
820
821*/
Bool_t
bool Bool_t
Int_t
int Int_t
kFALSE
constexpr Bool_t kFALSE
Double_t
double Double_t
kTRUE
constexpr Bool_t kTRUE
Option_t
const char Option_t
Fatal
void Fatal(const char *location, const char *msgfmt,...)
CbmStsDigiMatch
CbmStsDigiMatch.h.
Definition CbmStsDigiMatch.h:27
CbmStsDigi
CbmStsDigi.h.
Definition CbmStsDigi.h:32
CbmStsDigi::GetSide
Int_t GetSide() const
Side ( 0 = front, 1 = back ).
Definition CbmStsDigi.h:73
CbmStsDigi::GetSectorNr
Int_t GetSectorNr() const
Sector number.
Definition CbmStsDigi.h:70
CbmStsDigi::GetChannelNr
Int_t GetChannelNr() const
Channel number.
Definition CbmStsDigi.h:76
CbmStsDigi::GetStationNr
Int_t GetStationNr() const
Station number.
Definition CbmStsDigi.h:67
CbmStsIdealFindHits::fDigiPar
CbmStsDigiPar * fDigiPar
Geometry parameters.
Definition CbmStsIdealFindHits.h:55
CbmStsIdealFindHits::Exec
virtual void Exec(Option_t *opt)
Execution.
Definition CbmStsIdealFindHits.cxx:113
CbmStsIdealFindHits::~CbmStsIdealFindHits
virtual ~CbmStsIdealFindHits()
Destructor.
Definition CbmStsIdealFindHits.cxx:103
CbmStsIdealFindHits::fDigiScheme
CbmStsDigiScheme * fDigiScheme
Digitisation parameters.
Definition CbmStsIdealFindHits.h:56
CbmStsIdealFindHits::FindHits
Int_t FindHits(CbmStsStation *station, CbmStsSector *sector, std::set< Int_t > &fSet, std::set< Int_t > &bSet)
Find hits in one sector.
Definition CbmStsIdealFindHits.cxx:365
CbmStsIdealFindHits::CbmStsIdealFindHits
CbmStsIdealFindHits()
Default constructor.
Definition CbmStsIdealFindHits.cxx:40
CbmStsIdealFindHits::fTimer
TStopwatch fTimer
sector digis (back)
Definition CbmStsIdealFindHits.h:63
CbmStsIdealFindHits::fTime1
Double_t fTime1
Number of events with success.
Definition CbmStsIdealFindHits.h:92
CbmStsIdealFindHits::Init
virtual InitStatus Init()
Intialisation.
Definition CbmStsIdealFindHits.cxx:226
CbmStsIdealFindHits::ReInit
virtual InitStatus ReInit()
Reinitialisation.
Definition CbmStsIdealFindHits.cxx:264
CbmStsIdealFindHits::fPoints
TClonesArray * fPoints
AZ: Input array of CbmStsDigiMatch.
Definition CbmStsIdealFindHits.h:59
CbmStsIdealFindHits::MakeSets
void MakeSets()
Make sectorwise sets for sigis.
Definition CbmStsIdealFindHits.cxx:293
CbmStsIdealFindHits::fHits
TClonesArray * fHits
AZ: Input array of CbmStsPoint.
Definition CbmStsIdealFindHits.h:60
CbmStsIdealFindHits::fDigiMatch
TClonesArray * fDigiMatch
Input array of CbmStsDigi.
Definition CbmStsIdealFindHits.h:58
CbmStsIdealFindHits::SortDigis
void SortDigis()
Sort digis sectorwise.
Definition CbmStsIdealFindHits.cxx:316
CbmStsIdealFindHits::fDigis
TClonesArray * fDigis
Digitisation scheme.
Definition CbmStsIdealFindHits.h:57
CbmStsIdealFindHits::fDigiMapB
std::map< CbmStsSector *, std::set< Int_t > > fDigiMapB
sector digis (front)
Definition CbmStsIdealFindHits.h:62
CbmStsIdealFindHits::fDigiMapF
std::map< CbmStsSector *, std::set< Int_t > > fDigiMapF
Output array of CbmStsHit.
Definition CbmStsIdealFindHits.h:61
CbmStsIdealFindHits::SetParContainers
virtual void SetParContainers()
Get parameter containers.
Definition CbmStsIdealFindHits.cxx:205
CbmStsSector::GetType
Int_t GetType() const
Definition CbmStsSector.h:81
CbmStsSector::GetSectorNr
Int_t GetSectorNr() const
Definition CbmStsSector.h:63
CbmStsSector::GetDx
Double_t GetDx() const
Definition CbmStsSector.h:83
CbmStsSector::GetDy
Double_t GetDy() const
Definition CbmStsSector.h:84
CbmStsSector::GetStereoF
Double_t GetStereoF() const
Definition CbmStsSector.h:85
CbmStsSector::GetDetectorId
Int_t GetDetectorId() const
Definition CbmStsSector.h:60
CbmStsSector::GetRotation
Double_t GetRotation() const
Definition CbmStsSector.h:82
CbmStsSector::Intersect
Int_t Intersect(Int_t iFStrip, Int_t iBStrip, Double_t &xCross, Double_t &yCross, Double_t &zCross)
Definition CbmStsSector.cxx:272
CbmStsSector::GetStereoB
Double_t GetStereoB() const
Definition CbmStsSector.h:86
CbmStsStation::GetNSectors
Int_t GetNSectors() const
Definition CbmStsStation.h:62
CbmStsStation::GetZ
Double_t GetZ(Int_t it=0)
Definition CbmStsStation.cxx:151
CbmStsStation::GetStationNr
Int_t GetStationNr() const
Definition CbmStsStation.h:56
CbmStsStation::GetSector
CbmStsSector * GetSector(Int_t iSector)
Definition CbmStsStation.h:65
TVector3::SetXYZ
void SetXYZ(Double_t x, Double_t y, Double_t z)
TMath::Sqrt
Double_t Sqrt(Double_t x)
TMath::Cos
Double_t Cos(Double_t)
TMath::Abs
Double_t Abs(Double_t d)
TMath::Sin
Double_t Sin(Double_t)
TMath::Tan
Double_t Tan(Double_t)
TMath::RadToDeg
constexpr Double_t RadToDeg()