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