CbmSttTrackFinderIdeal.cxx 21 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834
  1. // -------------------------------------------------------------------------
  2. // ----- CbmSttTrackFinderIdeal source file -----
  3. // ----- Created 28/03/06 by V. Friese -----
  4. // -------------------------------------------------------------------------
  5. // C++ includes
  6. #include <iostream>
  7. #include <map>
  8. // ROOT includes
  9. #include "TClonesArray.h"
  10. #include "TCanvas.h"
  11. #include "TH2F.h"
  12. #include "TArc.h"
  13. #include "TDatabasePDG.h"
  14. #include "TRandom.h"
  15. // CBM includes
  16. #include "FairMCPoint.h"
  17. #include "FairRootManager.h"
  18. #include "CbmSttHit.h"
  19. #include "CbmSttTrack.h"
  20. #include "CbmSttTrackFinderIdeal.h"
  21. //#include "CbmSttHoughAccumulatorNew.h"
  22. #include "CbmSttHoughDefines.h"
  23. // ----- Default constructor -------------------------------------------
  24. CbmSttTrackFinderIdeal::CbmSttTrackFinderIdeal()
  25. {
  26. rootoutput = kFALSE;
  27. fMCTrackArray = NULL;
  28. fVerbose = 1;
  29. }
  30. // -------------------------------------------------------------------------
  31. // ----- Standard constructor ------------------------------------------
  32. CbmSttTrackFinderIdeal::CbmSttTrackFinderIdeal(Int_t verbose)
  33. {
  34. fMCTrackArray = NULL;
  35. fVerbose = verbose;
  36. if (verbose > 2) rootoutput = kTRUE;
  37. else rootoutput = kFALSE; // stt1
  38. }
  39. // -------------------------------------------------------------------------
  40. // ----- Destructor ----------------------------------------------------
  41. CbmSttTrackFinderIdeal::~CbmSttTrackFinderIdeal()
  42. {
  43. }
  44. // -------------------------------------------------------------------------
  45. // ----- Public method Init --------------------------------------------
  46. void CbmSttTrackFinderIdeal::Init()
  47. {
  48. // Get and check FairRootManager
  49. FairRootManager* ioman = FairRootManager::Instance();
  50. if (!ioman)
  51. {
  52. cout << "-E- CbmSttTrackFinderIdeal::Init: "
  53. << "RootManager not instantised!" << endl;
  54. return;
  55. }
  56. // Get MCTrack array
  57. // fMCTrackArray = (TClonesArray*) ioman->ActivateBranch("MCTrack");
  58. fMCTrackArray = (TClonesArray*) ioman->GetObject("MCTrack"); // EL
  59. if ( ! fMCTrackArray)
  60. {
  61. cout << "-E- CbmSttTrackFinderIdeal::Init: No MCTrack array!"
  62. << endl;
  63. return;
  64. }
  65. }
  66. // -------------------------------------------------------------------------
  67. // ----- Public method DoFind ------------------------------------------
  68. Int_t CbmSttTrackFinderIdeal::DoFind(TClonesArray* trackArray)
  69. {
  70. // Check pointers
  71. if ( !fMCTrackArray )
  72. {
  73. cout << "-E- CbmSttTrackFinderIdeal::DoFind: "
  74. << "MCTrack array missing! " << endl;
  75. return -1;
  76. }
  77. if (fHitCollectionList.GetEntries() == 0)
  78. {
  79. cout << "-E- CbmSttTrackFinderIdeal::DoFind: "
  80. << "No hit arrays present, call AddHitCollection() first (at least once)! " << endl;
  81. return -1;
  82. }
  83. if (fPointCollectionList.GetEntries() == 0)
  84. {
  85. cout << "-E- CbmSttTrackFinderIdeal::DoFind: "
  86. << "No point arrays present, call AddHitCollection() first (at least once)! " << endl;
  87. return -1;
  88. }
  89. if ( !trackArray )
  90. {
  91. cout << "-E- CbmSttTrackFinderIdeal::DoFind: "
  92. << "Track array missing! " << endl;
  93. return -1;
  94. }
  95. TCanvas
  96. *finderCanvas;
  97. if (rootoutput)
  98. {
  99. TH2F
  100. myRange("range", "range", 100, -50., 50., 100, -50., 50);
  101. finderCanvas = new TCanvas("findercanvas", "findercanvas", 800, 600);
  102. myRange.DrawCopy();
  103. plotAllStraws();
  104. }
  105. // Initialise control counters
  106. Int_t nNoMCTrack = 0;
  107. Int_t nNoTrack = 0;
  108. Int_t nNoSttPoint = 0;
  109. Int_t nNoSttHit = 0;
  110. // Create pointers to hit and SttPoint
  111. CbmSttHit* pMhit = NULL;
  112. FairMCPoint* pMCpt = NULL;
  113. FairMCTrack* pMCtr = NULL;
  114. CbmSttTrack* pTrck = NULL;
  115. // Number of STT hits
  116. Int_t
  117. nHits = 0;
  118. for (Int_t hitListCounter = 0; hitListCounter < fHitCollectionList.GetEntries(); hitListCounter++)
  119. {
  120. nHits += ((TClonesArray *)fHitCollectionList.At(hitListCounter))->GetEntriesFast();
  121. }
  122. // Declare some variables outside the loops
  123. Int_t ptIndex = 0; // MCPoint index
  124. Int_t mcTrackIndex = 0; // MCTrack index
  125. Int_t trackIndex = 0; // STTTrack index
  126. // Create STL map from MCtrack index to number of valid SttHits
  127. map<Int_t, map<Int_t, Int_t> >
  128. hitMap;
  129. // Loop over hits
  130. for (Int_t iHit = 0; iHit < nHits; iHit++)
  131. {
  132. pMhit = GetHitFromCollections(iHit);
  133. if ( ! pMhit )
  134. continue;
  135. ptIndex = pMhit->GetRefIndex();
  136. if (rootoutput)
  137. {
  138. TArc
  139. myArc;
  140. myArc.SetFillStyle(0);
  141. myArc.SetLineColor(1);
  142. if (pMhit->GetIsochrone() == 0)
  143. myArc.DrawArc(pMhit->GetX(), pMhit->GetY(), 1.0);
  144. else
  145. myArc.DrawArc(pMhit->GetX(), pMhit->GetY(), pMhit->GetIsochrone());
  146. }
  147. if (ptIndex < 0)
  148. continue; // fake or background hit
  149. pMCpt = GetPointFromCollections(iHit);
  150. if ( ! pMCpt )
  151. continue;
  152. mcTrackIndex = pMCpt->GetTrackID();
  153. Double_t
  154. wireX = pMhit->GetX(),
  155. wireY = pMhit->GetY();
  156. (hitMap[mcTrackIndex])[(Int_t)(wireX * wireX + wireY * wireY)]++;
  157. }
  158. // Create STL map from MCTrack index to SttTrack index
  159. map<Int_t, Int_t>
  160. correlationMap,
  161. trackMap;
  162. // Create STTTracks for reconstructable MCTracks
  163. Int_t nMCacc = 0; // accepted MCTracks (more than 3 points)
  164. Int_t nTracks = 0; // reconstructable MCTracks
  165. Int_t nMCTracks = fMCTrackArray->GetEntriesFast();
  166. for (Int_t iMCTrack = 0; iMCTrack < nMCTracks; iMCTrack++)
  167. {
  168. pMCtr = (FairMCTrack*) fMCTrackArray->At(iMCTrack);
  169. if ( ! pMCtr )
  170. continue;
  171. if (hitMap[iMCTrack].size() < 3)
  172. continue;
  173. // temporary hack, fix this in monte carlo ....
  174. if (TDatabasePDG::Instance()->GetParticle(pMCtr->GetPdgCode()) == NULL)
  175. continue;
  176. if (!(fabs(TDatabasePDG::Instance()->GetParticle(pMCtr->GetPdgCode())->Charge() / 3.0) > 0.))
  177. continue;
  178. nMCacc++;
  179. new((*trackArray)[nTracks]) CbmSttTrack();
  180. if (fVerbose) cout << "-I- CbmSttTrackFinderIdeal: STTTrack "
  181. << nTracks << " created from MCTrack "
  182. << iMCTrack << " (" << pMCtr->GetNPoints(kTOF)
  183. << " STTPoints)" << endl;
  184. correlationMap[nTracks] = iMCTrack;
  185. trackMap[iMCTrack] = nTracks++;
  186. }
  187. // Loop over hits. Get corresponding MCPoint and MCTrack index
  188. for (Int_t iHit = 0; iHit < nHits; iHit++)
  189. {
  190. pMhit = GetHitFromCollections(iHit);
  191. if ( ! pMhit )
  192. {
  193. cout << "-E- CbmSttTrackFinderIdeal::DoFind: Empty slot "
  194. << "in HitArray at position " << iHit << endl;
  195. nNoSttHit++;
  196. continue;
  197. }
  198. if (pMhit->GetIsochrone() == 0.)
  199. {
  200. cout << "mvd hit: " << pMhit->GetX() << " " << pMhit->GetY() << " " << pMhit->GetZ() << endl;
  201. }
  202. ptIndex = pMhit->GetRefIndex();
  203. if (ptIndex < 0)
  204. continue; // fake or background hit
  205. pMCpt = GetPointFromCollections(iHit);
  206. if ( ! pMCpt )
  207. {
  208. nNoSttPoint++;
  209. continue;
  210. }
  211. mcTrackIndex = pMCpt->GetTrackID();
  212. if (mcTrackIndex < 0 || mcTrackIndex > nMCTracks)
  213. {
  214. cout << "-E- CbmSttTrackFinderIdeal::DoFind: "
  215. << "MCTrack index out of range. " << mcTrackIndex << " "
  216. << nMCTracks << endl;
  217. nNoMCTrack++;
  218. continue;
  219. }
  220. if (trackMap.find(mcTrackIndex) == trackMap.end())
  221. continue;
  222. trackIndex = trackMap[mcTrackIndex];
  223. pTrck = (CbmSttTrack*) trackArray->At(trackIndex);
  224. if ( ! pTrck )
  225. {
  226. cout << "-E- CbmSttTrackFinderIdeal::DoFind: "
  227. << "No SttTrack pointer. " << iHit << " " << ptIndex
  228. << " " << mcTrackIndex << " " << trackIndex << endl;
  229. nNoTrack++;
  230. continue;
  231. }
  232. pTrck->AddHit(iHit, pMhit);
  233. if (rootoutput)
  234. {
  235. TArc
  236. myArc;
  237. myArc.SetFillStyle(0);
  238. myArc.SetLineColor(trackIndex + 2);
  239. if (pMhit->GetIsochrone() == 0)
  240. myArc.DrawArc(pMhit->GetX(), pMhit->GetY(), 1.0);
  241. else
  242. myArc.DrawArc(pMhit->GetX(), pMhit->GetY(), pMhit->GetIsochrone());
  243. }
  244. if (fVerbose > 2) cout << "Hit " << iHit << " from STTPoint "
  245. << ptIndex << " (MCTrack "
  246. << mcTrackIndex << ") added to STTTrack "
  247. << trackIndex << endl;
  248. }
  249. for (Int_t trackTeller = 0; trackTeller < nTracks; trackTeller++)
  250. {
  251. // loop over
  252. pTrck = (CbmSttTrack*) trackArray->At(trackTeller);
  253. if (pTrck != NULL)
  254. {
  255. Double_t
  256. dSeed,
  257. rSeed,
  258. phiSeed;
  259. GetTrack(dSeed, phiSeed, rSeed, correlationMap[trackTeller]);
  260. Double_t
  261. xSeed = (dSeed + rSeed) * cos(phiSeed),
  262. ySeed = (dSeed + rSeed) * sin(phiSeed);
  263. Double_t
  264. xSeed_old = xSeed,
  265. ySeed_old = ySeed,
  266. rSeed_old = rSeed;
  267. // ZoomTrack(dSeed, phiSeed, rSeed, pTrck); // CHECK!!!
  268. xSeed = (dSeed + rSeed) * cos(phiSeed);
  269. ySeed = (dSeed + rSeed) * sin(phiSeed);
  270. if (rootoutput)
  271. {
  272. TArc
  273. myArc;
  274. myArc.SetFillStyle(0);
  275. //myArc.SetLineColor(trackTeller + 2);
  276. //myArc.DrawArc(xSeed_old, ySeed_old, rSeed_old);
  277. myArc.SetLineColor(trackTeller + 2);
  278. myArc.DrawArc(xSeed, ySeed, rSeed);
  279. }
  280. pTrck->GetParamLast()->SetX(dSeed);
  281. pTrck->GetParamLast()->SetY(phiSeed);
  282. pTrck->GetParamLast()->SetZ(0.);
  283. pTrck->GetParamLast()->SetTx(rSeed);
  284. pTrck->GetParamLast()->SetTy(0.);
  285. pTrck->GetParamLast()->SetQp(0.);
  286. }
  287. }
  288. if (rootoutput)
  289. {
  290. finderCanvas->Update();
  291. finderCanvas->Show();
  292. char
  293. waitchar;
  294. cout << "press any key to continue." << endl;
  295. cin >> waitchar;
  296. delete finderCanvas;
  297. }
  298. if (fVerbose)
  299. {
  300. cout << endl;
  301. cout << "-------------------------------------------------------"
  302. << endl;
  303. cout << "-I- Ideal STT track finding -I-"
  304. << endl;
  305. cout << "Hits: " << nHits << endl;
  306. cout << "MCTracks: total " << nMCTracks << ", accepted " << nMCacc
  307. << ", reconstructable: " << nTracks << endl;
  308. cout << "SttHits not found : "
  309. << nNoSttHit << endl;
  310. cout << "SttPoints not found : "
  311. << nNoSttPoint << endl;
  312. cout << "MCTracks not found : "
  313. << nNoMCTrack << endl;
  314. cout << "SttTracks not found : "
  315. << nNoTrack << endl;
  316. cout << "-------------------------------------------------------"
  317. << endl;
  318. }
  319. else cout << "-I- CbmSttTrackFinderIdeal: all " << nMCTracks
  320. << ", acc. " << nMCacc << ", rec. " << nTracks << endl;
  321. return nTracks;
  322. }
  323. void CbmSttTrackFinderIdeal::GetTrackletCircular(Double_t firstX, Double_t firstY, Double_t firstR,
  324. Double_t secondX, Double_t secondY, Double_t secondR,
  325. Double_t thirdX, Double_t thirdY, Double_t thirdR,
  326. Double_t *circleRadii, Double_t *circleCentersX,
  327. Double_t *circleCentersY) const
  328. {
  329. Int_t
  330. trackletCounter = 0;
  331. Double_t
  332. small_limit = 0.0001;
  333. for (Int_t sign1 = -1; sign1 < 3; sign1 += 2)
  334. {
  335. for (Int_t sign2 = -1; sign2 < 3; sign2 += 2)
  336. {
  337. for (Int_t sign3 = -1; sign3 < 3; sign3 += 2)
  338. {
  339. // if three points are collinear, shift middle point by 50 microns
  340. if ((firstX - secondX == 0) && (secondX - thirdX == 0))
  341. {
  342. secondX += 0.005;
  343. secondY += 0.005;
  344. }
  345. else if (!((fabs(firstX - secondX) < small_limit) || (fabs(secondX - thirdX) < small_limit)))
  346. {
  347. if (fabs((firstY - secondY) / (firstX - secondX) - (secondY - thirdY) / (secondX - thirdX)) < small_limit)
  348. {
  349. secondX += 0.005;
  350. secondY += 0.005;
  351. }
  352. }
  353. Double_t
  354. a = -2. * (firstX - secondX),
  355. b = -2. * (firstY - secondY),
  356. c = 2. * (sign1 * firstR - sign2 * secondR),
  357. d = ((firstX * firstX + firstY * firstY - firstR * firstR) -
  358. (secondX * secondX + secondY * secondY - secondR * secondR)),
  359. e = -2. * (firstX - thirdX),
  360. f = -2. * (firstY - thirdY),
  361. g = 2. * (sign1 * firstR - sign3 * thirdR),
  362. h = ((firstX * firstX + firstY * firstY - firstR * firstR) -
  363. (thirdX * thirdX + thirdY * thirdY - thirdR * thirdR));
  364. Double_t
  365. A = -1. * (f * d - b * h) / (a * f - b * e),
  366. B = -1. * (f * c - b * g) / (a * f - b * e),
  367. C = (e * d - a * h) / (a * f - e * b),
  368. D = (e * c - a * g) / (a * f - e * b),
  369. I = B * B + D * D - 1.,
  370. II = 2. * A * B - 2. * B * firstX + 2. * C * D - 2. * D * firstY + sign1 * 2. * firstR,
  371. III = A * A - 2. * A * firstX + firstX * firstX + C * C - 2. * C * firstY + firstY * firstY - firstR * firstR;
  372. if ((fabs(I) > small_limit) && ((II * II - 4. * I * III) > 0.))
  373. {
  374. Double_t
  375. r = (-1. * II - sqrt(II * II - 4. * I * III)) / (2. * I),
  376. x = A + B * r,
  377. y = C + D * r;
  378. circleRadii[trackletCounter] = sqrt(r * r);
  379. circleCentersX[trackletCounter] = x;
  380. circleCentersY[trackletCounter] = y;
  381. trackletCounter++;
  382. }
  383. else
  384. {
  385. circleRadii[trackletCounter] = -1.;
  386. circleCentersX[trackletCounter] = 0.;
  387. circleCentersY[trackletCounter] = 0.;
  388. trackletCounter++;
  389. }
  390. }
  391. }
  392. }
  393. }
  394. // void CbmSttTrackFinderIdeal::ZoomTrack(Double_t &dSeed, Double_t &phiSeed,
  395. // Double_t &rSeed, CbmSttTrack *track)
  396. // {
  397. // track->SortHits();
  398. // if (track->GetNofHits() < 3)
  399. // {
  400. // cout << "-W- CbmSttTrackFinderHough::GetSeed2Circular: less than three hits found, not possible to fit!"
  401. // << endl;
  402. // return;
  403. // }
  404. // else
  405. // {
  406. // CbmSttHoughAccumulatorNew
  407. // *accumulator = new CbmSttHoughAccumulatorNew(2, kTRUE,
  408. // 150, dSeed - 5., dSeed + 5.,
  409. // 150, phiSeed - 0.5, phiSeed + 0.5,
  410. // 100, rSeed - 200., rSeed + 200.);
  411. // Int_t
  412. // segmentLength1 = track->GetNofHits() / 3,
  413. // segmentLength2 = track->GetNofHits() / 3,
  414. // segmentLength3 = track->GetNofHits() - segmentLength2 - segmentLength1;
  415. // Int_t
  416. // segmentStart1 = 0,
  417. // segmentStart2 = segmentLength1,
  418. // segmentStart3 = segmentLength1 + segmentLength2;
  419. // // select n random combinations
  420. // for (Int_t teller = 0; teller < 100; teller++)
  421. // {
  422. // Int_t
  423. // firstSample = gRandom->Integer(segmentLength1),
  424. // secondSample = gRandom->Integer(segmentLength2),
  425. // thirdSample = gRandom->Integer(segmentLength3);
  426. // Int_t
  427. // i = segmentStart1 + firstSample,
  428. // ii = segmentStart2 + secondSample,
  429. // iii = segmentStart3 + thirdSample;
  430. // CbmSttHit
  431. // *iFirstHit = GetHitFromCollections(track->GetHitIndex(i)),
  432. // *iSecondHit = GetHitFromCollections(track->GetHitIndex(ii)),
  433. // *iThirdHit = GetHitFromCollections(track->GetHitIndex(iii));
  434. // if ((!iFirstHit) || (!iSecondHit) || (!iThirdHit))
  435. // continue;
  436. // Double_t
  437. // firstX = iFirstHit->GetX(),
  438. // firstY = iFirstHit->GetY(),
  439. // firstR = iFirstHit->GetIsochrone(),
  440. // secondX = iSecondHit->GetX(),
  441. // secondY = iSecondHit->GetY(),
  442. // secondR = iSecondHit->GetIsochrone(),
  443. // thirdX = iThirdHit->GetX(),
  444. // thirdY = iThirdHit->GetY(),
  445. // thirdR = iThirdHit->GetIsochrone();
  446. // Double_t
  447. // circleRadii[8],
  448. // circleCentersX[8],
  449. // circleCentersY[8];
  450. // GetTrackletCircular(firstX, firstY, firstR,
  451. // secondX, secondY, secondR,
  452. // thirdX, thirdY, thirdR,
  453. // circleRadii, circleCentersX, circleCentersY);
  454. // // loop over all tracklets between this and next
  455. // for (Int_t trackletteller = 0; trackletteller < 8; trackletteller++)
  456. // {
  457. // if (circleRadii[trackletteller] > 0.)
  458. // {
  459. // Double_t
  460. // phi = atan(circleCentersY[trackletteller] / circleCentersX[trackletteller]),
  461. // dist = sqrt(circleCentersX[trackletteller] * circleCentersX[trackletteller] +
  462. // circleCentersY[trackletteller] * circleCentersY[trackletteller])
  463. // - circleRadii[trackletteller];
  464. // if (circleCentersX[trackletteller] < 0.)
  465. // {
  466. // if (circleCentersY[trackletteller] < 0.)
  467. // phi -= dPi;
  468. // else
  469. // phi += dPi;
  470. // }
  471. // accumulator->AddHit(dist, i,
  472. // phi, ii,
  473. // circleRadii[trackletteller], iii);
  474. // }
  475. // }
  476. // }
  477. // accumulator->Cluster();
  478. // dSeed = accumulator->GetHighestPeakX();
  479. // phiSeed = accumulator->GetHighestPeakY();
  480. // rSeed = accumulator->GetHighestPeakZ();
  481. // delete accumulator;
  482. // }
  483. // }
  484. void CbmSttTrackFinderIdeal::GetTrack(Double_t &dSeed, Double_t &phiSeed,
  485. Double_t &rSeed, Int_t mcTrackNo)
  486. {
  487. FairMCTrack
  488. *mcTrack = (FairMCTrack*) fMCTrackArray->At(mcTrackNo);
  489. TVector3 vectorMomentum;
  490. mcTrack->GetMomentum(vectorMomentum);
  491. TVector3 vectorTrack;
  492. mcTrack->GetStartVertex(vectorTrack);
  493. // TODO: read field from container
  494. rSeed = sqrt(vectorMomentum.X() * vectorMomentum.X() +
  495. vectorMomentum.Y() * vectorMomentum.Y()) / 0.006;
  496. Double_t phiStart = atan(vectorMomentum.Y() / vectorMomentum.X()),
  497. xStart = vectorTrack.X(),
  498. yStart = vectorTrack.Y();
  499. if (vectorMomentum.X() < 0.)
  500. {
  501. if (vectorMomentum.Y() < 0.)
  502. phiStart -= dPi;
  503. else
  504. phiStart += dPi;
  505. }
  506. Double_t
  507. sign = 1.;
  508. if (TDatabasePDG::Instance()->GetParticle(mcTrack->GetPdgCode()) != NULL)
  509. {
  510. if (TDatabasePDG::Instance()->GetParticle(mcTrack->GetPdgCode())->Charge() < 0.)
  511. {
  512. sign = -1.;
  513. }
  514. }
  515. Double_t
  516. xCircleCenter = xStart + sign * (rSeed * sin(phiStart)),
  517. yCircleCenter = yStart - sign * (rSeed * cos(phiStart));
  518. dSeed = sqrt(xCircleCenter * xCircleCenter + yCircleCenter * yCircleCenter) - rSeed;
  519. phiSeed = atan(yCircleCenter / xCircleCenter);
  520. if (xCircleCenter < 0.)
  521. {
  522. if (yCircleCenter < 0.)
  523. phiSeed -= dPi;
  524. else
  525. phiSeed += dPi;
  526. }
  527. }
  528. // -------------------------------------------------------------------------
  529. Bool_t CbmSttTrackFinderIdeal::putStraw(Double_t xpos, Double_t ypos, Double_t radius)
  530. {
  531. Double_t
  532. pipeDiam = 0.42,
  533. tubeOuterDiam = 1.032,
  534. CoverThickness = 0.1,
  535. outerRadius = 42.,
  536. innerRadius = 15.;
  537. if ((sqrt(xpos * xpos + ypos * ypos) < outerRadius - CoverThickness - (tubeOuterDiam / 2.)) &&
  538. (sqrt(xpos * xpos + ypos * ypos) > innerRadius + CoverThickness + (tubeOuterDiam / 2.)) &&
  539. (sqrt(ypos * ypos) > (pipeDiam + tubeOuterDiam / 2.)))
  540. {
  541. TArc
  542. myArc;
  543. myArc.SetFillStyle(0);
  544. myArc.SetLineColor(17);
  545. myArc.DrawArc(xpos, ypos, radius);
  546. }
  547. else
  548. {
  549. return false;
  550. }
  551. return true;
  552. }
  553. void CbmSttTrackFinderIdeal::plotAllStraws()
  554. {
  555. Double_t
  556. tubeOuterDiam = 1.032,
  557. sttCenterX = 0.,
  558. sttCenterY = 0.;
  559. Int_t
  560. ringmax = -1;
  561. Bool_t
  562. started = kFALSE,
  563. goOn = kTRUE;
  564. Int_t
  565. ringteller = 18;
  566. while (goOn)
  567. {
  568. goOn = kFALSE;
  569. Double_t
  570. sqrt3 = sqrt(3.),
  571. radius = tubeOuterDiam / 2.;
  572. Double_t
  573. zpos = radius;
  574. // place first
  575. Double_t
  576. xpos = sttCenterX + ((ringteller) * 2 * radius),
  577. ypos = sttCenterY;
  578. for(Int_t i = 0; i < ringteller; i++)
  579. {
  580. xpos -= radius;
  581. ypos += sqrt3 * radius;
  582. if (putStraw(xpos, ypos, zpos))
  583. goOn = kTRUE;
  584. }
  585. for(Int_t i = 0; i < ringteller; i++)
  586. {
  587. xpos -= 2 * radius;
  588. if (putStraw(xpos, ypos, zpos))
  589. goOn = kTRUE;
  590. }
  591. for(Int_t i = 0; i < ringteller; i++)
  592. {
  593. xpos -= radius;
  594. ypos -= sqrt3 * radius;
  595. if (putStraw(xpos, ypos, zpos))
  596. goOn = kTRUE;
  597. }
  598. for(Int_t i = 0; i < ringteller; i++)
  599. {
  600. xpos += radius;
  601. ypos -= sqrt3 * radius;
  602. if (putStraw(xpos, ypos, zpos))
  603. goOn = kTRUE;
  604. }
  605. for(Int_t i = 0; i < ringteller; i++)
  606. {
  607. xpos += 2 * radius;
  608. if (putStraw(xpos, ypos, zpos))
  609. goOn = kTRUE;
  610. }
  611. for(Int_t i = 0; i < ringteller; i++)
  612. {
  613. xpos += radius;
  614. ypos += sqrt3 * radius;
  615. if (putStraw(xpos, ypos, zpos))
  616. goOn = kTRUE;
  617. }
  618. if (goOn)
  619. started = kTRUE;
  620. if (!started)
  621. goOn = kTRUE;
  622. ringteller++;
  623. if ((ringmax > 0) && (ringteller == ringmax))
  624. goOn = kFALSE;
  625. }
  626. }
  627. CbmSttHit* CbmSttTrackFinderIdeal::GetHitFromCollections(Int_t hitCounter)
  628. {
  629. CbmSttHit
  630. *retval = NULL;
  631. Int_t
  632. relativeCounter = hitCounter;
  633. for (Int_t collectionCounter = 0; collectionCounter < fHitCollectionList.GetEntries(); collectionCounter++)
  634. {
  635. Int_t
  636. size = ((TClonesArray *)fHitCollectionList.At(collectionCounter))->GetEntriesFast();
  637. if (relativeCounter < size)
  638. {
  639. retval = (CbmSttHit*) ((TClonesArray *)fHitCollectionList.At(collectionCounter))->At(relativeCounter);
  640. break;
  641. }
  642. else
  643. {
  644. relativeCounter -= size;
  645. }
  646. }
  647. return retval;
  648. }
  649. FairMCPoint* CbmSttTrackFinderIdeal::GetPointFromCollections(Int_t hitCounter)
  650. {
  651. FairMCPoint
  652. *retval = NULL;
  653. Int_t
  654. relativeCounter = hitCounter;
  655. for (Int_t collectionCounter = 0; collectionCounter < fHitCollectionList.GetEntries(); collectionCounter++)
  656. {
  657. Int_t
  658. size = ((TClonesArray *)fHitCollectionList.At(collectionCounter))->GetEntriesFast();
  659. if (relativeCounter < size)
  660. {
  661. Int_t
  662. tmpHit = ((CbmSttHit*) ((TClonesArray *)fHitCollectionList.At(collectionCounter))->At(relativeCounter))->GetRefIndex();
  663. retval = (FairMCPoint*) ((TClonesArray *)fPointCollectionList.At(collectionCounter))->At(tmpHit);
  664. break;
  665. }
  666. else
  667. {
  668. relativeCounter -= size;
  669. }
  670. }
  671. return retval;
  672. }
  673. ClassImp(CbmSttTrackFinderIdeal)