BOSS 7.0.6
BESIII Offline Software System
Loading...
Searching...
No Matches
TConformalFinder0.cxx
Go to the documentation of this file.
1//-----------------------------------------------------------------------------
2// $Id: TConformalFinder0.cxx,v 1.14 2010/03/31 09:58:59 liucy Exp $
3//-----------------------------------------------------------------------------
4// Filename : TConformalFinder0.cc
5// Section : Tracking
6// Owner : Yoshi Iwasaki
7// Email : [email protected]
8//-----------------------------------------------------------------------------
9// Description : A class to find tracks with the conformal method.
10// See http://bsunsrv1.kek.jp/~yiwasaki/tracking/
11//-----------------------------------------------------------------------------
12// TConformalLink removed. TSegment added.
13//-----------------------------------------------------------------------------
14
15#include "TrkReco/TMDCUtil.h"
16#include "TrkReco/TMDCWire.h"
17#include "TrkReco/TMDCWireHit.h"
18#include "TrkReco/TMDCWireHitMC.h"
19#include "TrkReco/TConformalFinder0.h"
20#include "TrkReco/TMLink.h"
21#include "TrkReco/THistogram.h"
22#include "TrkReco/TCircle.h"
23#include "TrkReco/TTrack.h"
24#include "TrkReco/TSegment0.h"
25//#include "cdc/Range.h"
26#include "TrkReco/Range.h"
27
28std::string
29TConformalFinder0::version(void) const {
30 return "1.63";
31}
32
33TConformalFinder0::TConformalFinder0(float maxSigma,
34 float fraction,
35 float stereoZ3,
36 float stereoZ4,
37 float stereoChisq3,
38 float stereoChisq4,
39 float stereoMaxSigma,
40 unsigned fittingCorrections,
41 float salvageLevel,
42 bool cosmic)
43: TFinderBase(),
44 _builder(0),
45 _doStereo(true),
46 _doSalvage(true), //liuqg
47 _fraction(fraction) {
48
49 //...Parameters for a track...
50 _trackSelector.nLinks(4);
51 _trackSelector.nSuperLayers(2);
52 _trackSelector.minPt(0.05);
53 _trackSelector.maxImpact(100.);
54 _trackSelector.maxSigma(maxSigma);
55 _trackSelector.nLinksStereo(3);
56 _trackSelector.maxDistance(30.);
57
58 //...Make a builder...
59 if (cosmic) _builder = new TBuilderCosmic("cosmic builder", salvageLevel);
60 else _builder = new TBuilder0("conformal builder",
61 stereoZ3,
62 stereoZ4,
63 stereoChisq3,
64 stereoChisq4,
65 stereoMaxSigma,
66 fittingCorrections,
67 salvageLevel);
68
69 //...Set up TBuilder...
70 _builder->trackSelector(_trackSelector);
71}
72
73TConformalFinder0::~TConformalFinder0() {
74 delete _builder;
75}
76
77void
78TConformalFinder0::dump(const std::string & msg, const std::string & pre) const {
79 std::cout << pre;
80 TFinderBase::dump(msg);
81 std::cout << pre;
82 if (msg.find("state") != std::string::npos) {
83 std::cout << "#axialConfPos=" << _axialConfLinks.length();
84 std::cout << ",#stereoConfPos=" << _stereoConfLinks.length();
85 }
86}
87
88void
89TConformalFinder0::clear(void) {
90 HepAListDeleteAll(_axialConfLinks);
91 HepAListDeleteAll(_stereoConfLinks);
92 _unusedAxialConfLinks.removeAll();
93 _unusedStereoConfLinks.removeAll();
94 _goodAxialConfLinks.removeAll();
95 HepAListDeleteAll(_circles);
96 _tracks.removeAll();
97}
98
99void
100TConformalFinder0::conformalTransformation(const HepPoint3D & center,
101 const AList<TMDCWireHit> & hits,
102 AList<TMLink> & links) {
103
104 unsigned nHits = hits.length();
105 if (center == ORIGIN) {
106 for (unsigned i = 0; i < nHits; i++) {
107 TMDCWireHit * h = hits[i];
108 const HepPoint3D & p = h->xyPosition();
109
110//zsl HepPoint3D cp(2. * p.x() / p.mag2(), 2. * p.y() / p.mag2());
111 HepPoint3D cp(2. * p.x() / p.mag2(), 2. * p.y() / p.mag2(), 0.);
112// std::cout<<" xypo "<<p<<" pos "<<cp<<std::endl;
113 links.append(new TMLink(0, h, cp));
114 }
115 }
116 else {
117 for (unsigned i = 0; i < nHits; i++) {
118 TMDCWireHit * h = hits[i];
119 HepPoint3D p(h->xyPosition() - center);
120 HepPoint3D cp(2. * p.x() / p.mag2(), 2. * p.y() / p.mag2(), 0.);
121 links.append(new TMLink(0, h, cp));
122 }
123 }
124}
125
126void
127TConformalFinder0::conformalTransformationDriftCircle(const HepPoint3D & center,
128 const AList<TMDCWireHit> & hits,
129 AList<TMLink> & links) { //added by Liuqg for Tsf
130 unsigned nHits = hits.length();
131 if (center == ORIGIN) {
132 for (unsigned i = 0; i < nHits; i++) {
133 TMDCWireHit * h = hits[i];
134 const HepPoint3D & p = h->xyPosition();
135
136 const double r = 0.5*(h->drift(0) + h->drift(1));
137 HepPoint3D cp(2. * p.x() / p.mag2(), 2. * p.y() / p.mag2(), 0.);
138// HepPoint3D cp2(100 * 2. * p.x() / (p.mag2() - r*r), 100 * 2. * p.y() / (p.mag2() - r*r), 0.);
139 HepPoint3D cp2(2. * p.x() / (p.mag2() - r*r), 2. * p.y() / (p.mag2() - r*r), 0.);
140// float cDrift = 100 * 2. * r / (p.mag2() - r*r);
141 double cDrift = 2. * r / (p.mag2() - r*r);
142
143 links.append(new TMLink(0, h, cp, cp2, cDrift));
144 }
145 }
146 else {
147 for (unsigned i = 0; i < nHits; i++) {
148 TMDCWireHit * h = hits[i];
149 HepPoint3D p(h->xyPosition() - center);
150
151 const double r = 0.5*(h->drift(0) + h->drift(1));
152 HepPoint3D cp(2. * p.x() / p.mag2(), 2. * p.y() / p.mag2(), 0.);
153 //unit of the following is km-1, origin is cm.
154// HepPoint3D cp2(100 * 2. * p.x() / (p.mag2() - r*r), 100 * 2. * p.y() / (p.mag2() - r*r), 0.);
155 HepPoint3D cp2(2. * p.x() / (p.mag2() - r*r), 2. * p.y() / (p.mag2() - r*r), 0.);
156// float cDrift = 100 * 2. * r / (p.mag2() - r*r);
157 double cDrift = 2. * r / (p.mag2() - r*r);
158
159 links.append(new TMLink(0, h, cp, cp2, cDrift));
160 }
161 }
162}
163
164void
165TConformalFinder0::conformalTransformationRphi(const HepPoint3D & center,
166 const AList<TMDCWireHit> & hits,
167 AList<TMLink> & links) {
168
169 unsigned nHits = hits.length();
170 if (center == ORIGIN) {
171 for (unsigned i = 0; i < nHits; i++) {
172 TMDCWireHit * h = hits[i];
173 const HepPoint3D & p = h->xyPosition();
174 HepPoint3D cp(2. * p.x() / p.mag2(), 2. * p.y() / p.mag2(), 0.);
175 double r = log(cp.mag()) + 4.;
176 double phi = atan2(cp.y(), cp.x()) + M_PI;
177 HepPoint3D cpt(phi, r, 0.);
178 links.append(new TMLink(0, h, cpt));
179 }
180 }
181 else {
182 for (unsigned i = 0; i < nHits; i++) {
183 TMDCWireHit * h = hits[i];
184 HepPoint3D p(h->xyPosition() - center);
185 HepPoint3D cp(2. * p.x() / p.mag2(), 2. * p.y() / p.mag2(), 0.);
186 double r = log(cp.mag()) + 4.;
187 double phi = atan2(cp.y(), cp.x()) + M_PI;
188 HepPoint3D cpt(phi, r, 0.);
189 links.append(new TMLink(0, h, cpt));
190 }
191 }
192}
193
194AList<TSegment0>
195TConformalFinder0::findClusters(const THistogram & hist) const {
196
197 //...Obtain raw clusters...
198 AList<TSegment0> list = hist.clusters0();
199 unsigned n = list.length();
200 if (n == 0) return list;
201
202#ifdef TRKRECO_DEBUG_DETAIL
203 // static TChecker chk0("clusters before splitting");
204 // chk0.check(list);
205 // chk0.dump("detail", " ");
206#endif
207
208 //...Examine each cluster...
209 AList<TSegment0> splitted;
210 for (unsigned i = 0; i < n; i++) {
211 TSegment0 * c = list[i];
212
213 AList<TSegment0> newClusters = c->split();
214 if (newClusters.length() == 0) {
215 c->solveDualHits();
216 continue;
217 }
218
219 list.append(newClusters);
220 splitted.append(c);
221#ifdef TRKRECO_DEBUG_DETAIL
222 c->dump("hits", " ");
223 std::cout << " ... splitted as" << std::endl;
224 for (unsigned j = 0; j < newClusters.length(); j++) {
225 std::cout << " " << j << " : ";
226 newClusters[j]->dump("hits");
227 }
228#endif
229 }
230 list.remove(splitted);
231 HepAListDeleteAll(splitted);
232
233#ifdef TRKRECO_DEBUG_DETAIL
234 // static TChecker chk1("clusters after splitting");
235 // chk1.check(list);
236 // chk1.dump("detail", " ");
237#endif
238
239 return list;
240}
241
242AList<TSegment0>
243TConformalFinder0::findClusters2(const THistogram & hist) const {
244
245 //...Obtain raw clusters...
246 AList<TSegment0> list = hist.clusters0();
247 unsigned n = list.length();
248 if (n == 0) return list;
249
250#ifdef TRKRECO_DEBUG_DETAIL
251 // static TChecker chk0("clusters before splitting (2)");
252 // chk0.check(list);
253 // chk0.dump("detail", " ");
254#endif
255
256 //...Examine each cluster...
257 for (unsigned i = 0; i < n; i++) {
258 TSegment0 * c = list[i];
259 unsigned type = c->clusterType();
260
261 if ((type == 1) || (type == 2)) {
262 c->dump("hits mc", " ");
263 c->solveDualHits();
264 }
265 }
266
267#ifdef TRKRECO_DEBUG_DETAIL
268 // static TChecker chk1("clusters after splitting (2)");
269 // chk1.check(list);
270 // chk1.dump("detail", " ");
271#endif
272
273 return list;
274}
275
276AList<TMLink>
277TConformalFinder0::findCloseHits(const AList<TMLink> & links,
278 const TTrack & track) const {
279 //
280 // Coded by J.Suzuki
281 //
282 AList<TMLink> list;
283
284 //...Check condition...
285 if (track.links().length() == 0) {
286#ifdef TRKRECO_DEBUG_DETAIL
287 std::cout << "TConformalFinder0::findCloseHits !!! ";
288 std::cout << " no links found in a track : This should not be happened";
289 std::cout << std::endl;
290#endif
291
292 return list;
293 }
294
295 //...Parameters...
296 // float dRcut[11] = {0, 3.5, 0., 5.5, 0., 6.5, 0., 7.5, 0., 9.5, 0.};
297// float dRcut[11] = {0., 4.3, 0., 6.5, 0., 7.5, 0., 8.0, 0., 9.5, 0.};
298 float dRcut[11] = {4.3, 6.5, 0., 0., 0., 7.5, 8.0, 9.5, 11.0, 0., 0.}; //Liuqg, Bes
299
300 //...Select Stereo hits associated to the current r-phi curve...
301 double R0 = track.helix().curv();
302 double xInnerWire = track.links()[0]->wire()->xyPosition().x();
303 double yInnerWire = track.links()[0]->wire()->xyPosition().y();
304 unsigned nall = links.length();
305 for (unsigned j = 0; j < nall; j++) {
306 TMLink & t = * links[j];
307 const TMDCWire & w = * t.wire();
308 HepVector3D X = w.xyPosition() - track.helix().center();
309 double Rmag2 = X.mag2();
310 double DR = fabs(sqrt(Rmag2) - fabs(R0));
311 t.zStatus(-10);
312 t.zPair(0);
313 if (DR < dRcut[w.superLayerId()] &&
314 (xInnerWire*w.xyPosition().x()+yInnerWire*w.xyPosition().y())>0.){
315 list.append(t);
316 }
317 }
318
319 //cout<<"stereo hits: "<<list.length()<<endl;
320 return list;
321}
322
323TSegment0 *
324TConformalFinder0::findBestLink(const TSegment0 & base,
325 const AList<TSegment0> & candidates) const {
326 //...Parameters...
327 double minAngle = 0.80;
328 double maxDistance = 0.3;
329// double minAngle = 0.40;
330// double maxDistance = 0.5;
331
332 //...Candidate loop...
333 unsigned n = candidates.length();
334 double minDistance = 999.;
335 TSegment0 * best = NULL;
336 for (unsigned j = 0; j < n; j++) {
337 TSegment0 * current = candidates[j];
338 if (current->nLinks() < 2) continue;
339
340 float angle = base.direction().dot(current->direction());
341#ifdef TRKRECO_DEBUG_DETAIL
342 current->dump("vector hits mc", " ");
343 std::cout << " angle=" << angle;
344 if (angle < minAngle) std::cout << std::endl;
345#endif
346// cout<<" angle = "<<angle<<endl;
347// if(angle < 0.3) cout<<" base.direction(): "<<base.direction()
348// <<" current->direction(): "<<current->direction()<<endl;
349 if (angle < minAngle) continue;
350
351 float distance = base.distance(* current);
352 if (distance < minDistance) {
353 minDistance = distance;
354 best = current;
355 }
356#ifdef TRKRECO_DEBUG_DETAIL
357 std::cout << ",dist=" << distance << std::endl;
358#endif
359 }
360
361// if (minDistance < maxDistance) return best;
362 return best;
363// return NULL; //Liuqg.....
364}
365
366TSegment0 *
367TConformalFinder0::appendCluster(TTrack & t, AList<TSegment0> & list) const {
368
369 //...Candidate loop...
370 unsigned n = list.length();
371 TSegment0 * best = NULL;
372 unsigned nBest = 0;
373 for (unsigned j = 0; j < n; j++) {
374 TSegment0 * c = list[j];
375
376 unsigned nOk = t.testByApproach(c->links(),
377 _trackSelector.maxSigma());
378 if (nOk > nBest) {
379 nBest = nOk;
380 best = c;
381 }
382 }
383
384 //...Try to append...
385 if (best) {
386#ifdef TRKRECO_DEBUG_DETAIL
387 std::cout << " ... appending a cluster" << std::endl;
388 best->dump("hits mc", " ");
389#endif
390 AList<TMLink> links(best->links());
391 t.appendByApproach(links, _trackSelector.maxSigma());
392 return best;
393 }
394
395 return NULL;
396}
397
398AList<TSegment0>
399TConformalFinder0::findClusterLink(TSegment0 & base,
400 const AList<TSegment0> * const list) const {
401
402#ifdef TRKRECO_DEBUG_DETAIL
403 std::cout << name() << " ... finding cluster linkage" << std::endl;
404 if (base.links().length() == 0)
405 std::cout << name() << " !!! base doesn't have any TMLink." << std::endl;
406 std::cout << "... base cluster" << std::endl;
407 base.dump("cluster hits mc", " ->");
408#endif
409
410 //...Preparation of return value...
411 AList<TSegment0> seeds;
412 seeds.append(base);
413
414 //...Which super layer?...
415// unsigned outerMost = (base.links())[0]->wire()->superLayerId() / 2;
416 unsigned outerMost = (base.links())[0]->wire()->axialStereoLayerId() / 4;
417
418 //...Inner super layer loop...
419 int next = outerMost;
420 TSegment0 * last = & base;
421 while (next) {
422 --next;
423 const AList<TSegment0> & candidates = list[next];
424 if (candidates.length() == 0) continue;
425
426#ifdef TRKRECO_DEBUG_DETAIL
427 std::cout << "... clusters in super layer " << next << std::endl;
428#endif
429
430 //...Find best match...
431 TSegment0 * best = findBestLink(* last, candidates);
432 if (best != NULL) {
433 seeds.append(best);
434 last = best;
435#ifdef TRKRECO_DEBUG_DETAIL
436 std::cout << " ->Best is ";
437 std::cout << best->position() << " ";
438 best->dump("hits mc");
439#endif
440 }
441 }
442
443 return seeds;
444}
445
446TTrack *
447TConformalFinder0::makeTrack(const AList<TSegment0> & list) const {
448 AList<TMLink> links;
449 for (unsigned i = 0; i < list.length(); i++) {
450 const AList<TMLink> & tmp = list[i]->links();
451 unsigned n = tmp.length();
452 for (unsigned j = 0; j < n; j++) {
453 if (tmp[j]->hit()->track()) continue;
454 links.append(tmp[j]);
455 }
456 }
457
458 TTrack * t = _builder->buildRphi(links);
459
460 return t;
461}
462
463AList<TSegment0>
464TConformalFinder0::findCloseClusters(const TTrack & track,
465 const AList<TSegment0> & list,
466 double maxDistance) const {
467
468 //...Cal. direction of rotation of track...
469 double radius = fabs(track.helix().radius());
470 const HepPoint3D & center = track.helix().center();
471 int rotation =
472 (center.cross(track.links()[0]->xyPosition()).z() > 0.) ? 1 : -1;
473
474#ifdef TRKRECO_DEBUG_DETAIL
475 std::cout << name() << " ... finding close clusters:maxDistance=";
476 std::cout << maxDistance << std::endl;
477 std::cout << " radius,center,rotation=" << radius << ",";
478 std::cout << center << "," << rotation << std::endl;
479#endif
480
481 //...Cluster loop...
482 AList<TSegment0> close;
483 unsigned n = list.length();
484 for (unsigned i = 0; i < n; i++) {
485 TSegment0 & c = * list[i];
486
487 //...Cal. position of cluster in the real plane...
488 HepPoint3D position(ORIGIN);
489 unsigned m = c.links().length();
490 for (unsigned j = 0; j < m; j++) {
491 position += c.links()[j]->xyPosition();
492 }
493 position *= 1. / double(m);
494
495#ifdef TRKRECO_DEBUG_DETAIL
496 c.dump("cluster hits mc", " ");
497 std::cout << " position=" << position;
498 std::cout << ",diff=" << (position - center).mag() - radius << std::endl;
499#endif
500 //...Cal. distance to a track...
501 HepVector3D diff = position - center;
502 if ((diff.mag() - radius) < maxDistance) {
503
504 //...Same side?...
505 int direction = (center.cross(position).z() > 0.) ? 1 : -1;
506 if (direction == rotation) close.append(c);
507 }
508 }
509
510#ifdef TRKRECO_DEBUG_DETAIL
511 std::cout << " found clusters" << std::endl;
512 for (unsigned i = 0; i < close.length(); i++) {
513 close[i]->dump("hits mc", " ");
514 }
515#endif
516 return close;
517}
518
519void
520TConformalFinder0::appendClusters2(TTrack & track,
521 AList<TSegment0> & list) const {
522
523#ifdef TRKRECO_DEBUG_DETAIL
524 std::cout << name() << " ... appending clusters remained" << std::endl;
525 std::cout << " clusters to be tested : " << std::endl;
526 for (unsigned i = 0; i < list.length(); i++) {
527 list[i]->dump("cluster hits mc", " ");
528 }
529#endif
530
531 unsigned n = list.length();
532 if (n == 0) return;
533
534 AList<TSegment0> closer;
535 closer.append(findCloseClusters(track, list, 1.));
536
537#ifdef TRKRECO_DEBUG_DETAIL
538 std::cout << " found clusters" << std::endl;
539 for (unsigned i = 0; i < closer.length(); i++) {
540 closer[i]->dump("cluster hits mc", " ");
541 }
542#endif
543 n = closer.length();
544 if (closer.length() == 0) return;
545
546 //...Append them...
547 AList<TMLink> candidates;
548 for (unsigned i = 0; i < n; i++)
549 candidates.append(closer[i]->links());
550 _builder->appendClusters(track, candidates);
551
552 //...Remove TMLinks from clusters...
553 for (unsigned i = 0; i < n; i++) {
554 closer[i]->TTrackBase::remove(track.links());
555 if (closer[i]->nLinks() == 0) list.remove(closer[i]);
556 }
557}
558
559int
560TConformalFinder0::doit(const AList<TMDCWireHit> & axialHits,
561 const AList<TMDCWireHit> & stereoHits,
562 AList<TTrack> & tracks,
563 AList<TTrack> &) {
564
565 //...For debug...
566 if (debugLevel()) {
567 std::cout << name() << " ... processing" << std::endl;
568 std::cout << " axialHits=" << axialHits.length();
569 std::cout << ",stereoHits=" << stereoHits.length();
570 std::cout << ",tracks=" << tracks.length();
571 std::cout << std::endl;
572
573 if (debugLevel() > 1)
574 std::cout << name() << " ... conformal transformation0" << std::endl;
575 }
576
577 //...Conformal transformation with IP constraint...
578 conformalTransformationRphi(ORIGIN, axialHits, _axialConfLinks);
579 conformalTransformationRphi(ORIGIN, stereoHits, _stereoConfLinks);
580 _unusedAxialConfLinks.append(_axialConfLinks);
581 _unusedStereoConfLinks.append(_stereoConfLinks);
582 AList<TMLink> unusedConfLinks;
583 if (_doSalvage) {
584 unusedConfLinks.append(_axialConfLinks);
585 unusedConfLinks.append(_stereoConfLinks);
586 }
587
588 //...For debug...
589 if (debugLevel() > 1)
590 std::cout << name() << " ... selecting good hits" << std::endl;
591
592 //...Select good axial hits...
593 AList<TMLink> goodHits;
594 int nLinks = _axialConfLinks.length();
595 for (unsigned i = 0; i < nLinks; i++) {
596 TMLink * l = _axialConfLinks[i];
597 const TMDCWireHit & h = * l->hit();
598//liuqg if ((h.state() & WireHitIsolated) &&
599// (h.state() & WireHitContinuous))
600 goodHits.append(l);
601 }
602 //...Main algorithm...
603 standardFinding(goodHits, unusedConfLinks, _fraction);
604
605 //...Main algorithm for second trial...
606 specialFinding(goodHits, unusedConfLinks, _fraction);
607
608 //...For debug...
609 if (debugLevel()) {
610 std::cout << name() << " ... processed : ";
611 std::cout << "good hits=" << goodHits.length();
612 std::cout << ",tracks=" << _tracks.length();
613 std::cout << std::endl;
614 }
615
616 tracks.append(_tracks);
617 return 0;
618}
619
620void
621TConformalFinder0::standardFinding(AList<TMLink> & list,
622 AList<TMLink> & unusedLinks,
623 double fraction) {
624#ifdef TRKRECO_DEBUG_DETAIL
625 std::cout << name() << " ... standard finding with salvage : given hits :";
626 Dump(list, "sort");
627#endif
628
629 //...Find segments...
630 AList< AList<TSegment0> > segments = findSegments(list);
631 AList<TSegment0> segmentList[5];
632 AList<TSegment0> original[5];
633 for (unsigned i = 0; i < 5; i++) {
634 segmentList[i] = * segments[i];
635 original[i] = * segments[i];
636// cout<< " segment: " << i << " : " << segmentList[i].length() << endl;
637// for (unsigned j = 0; j < segmentList[i].length(); j++)
638// cout<<"links in seg: "<<segmentList[i][j]->nLinks()<<endl;
639 }
640
641#ifdef TRKRECO_DEBUG_DETAIL
642 for (unsigned i = 0; i < 5; i++) {
643 std::cout << "... clusters in super layer " << i << std::endl;
644 for (unsigned j = 0; j < segmentList[i].length(); j++) {
645 segmentList[i][j]->dump("", " ");
646 segmentList[i][j]->dump("hits mc", " ");
647 }
648 }
649#endif
650
651 //...Main loop...
652 AList<TSegment0> retryList;
653 AList<TTrack> salvageList;
654 unsigned outerMost = 4;
655 while (outerMost) {
656
657 while (TSegment0 * base = segmentList[outerMost][0]) {
658
659 //...Get linked clusters...
660 AList<TSegment0> clusters = findClusterLink(* base, segmentList);
661
662 //...Make a track...
663 TTrack * t = makeTrack(clusters); // don't change 'clusters'
664 if (t == NULL) {
665 retryList.append(base);
666 segmentList[outerMost].remove(base);
667 continue;
668 }
669
670 //...Check track quality...
671 // double f = float(t->nLinks()) / float(nTLinks(clusters));
672 double f = float(t->nCores()) / float(NCoreLinks(clusters));
673 if (f < fraction) {
674#ifdef TRKRECO_DEBUG_DETAIL
675 std::cout << "... fraction too low:" << f << std::endl;
676 std::cout << " used cores=" << t->nCores();
677 std::cout << ", candidate cores=" << NCoreLinks(clusters);
678 std::cout << " ... retry later" << std::endl;
679#endif
680 retryList.append(base);
681 segmentList[outerMost].remove(base);
682 delete t;
683 continue;
684 }
685
686 //...Append other hits...
687 appendClusters2(* t, retryList);
688
689#ifdef TRKRECO_DEBUG_DETAIL
690 std::cout << name() << " ... 2D result :" << std::endl;
691 t->dump("detail", " ");
692#endif
693
694 //...Make it 3D...
695 TTrack * ts = t;
696 if (_doStereo) {
697 ts = _builder->buildStereo(* t,
698 findCloseHits(_unusedStereoConfLinks, * t));
699 }
700
701 //...Check track quality...
702 if (! ts) {
703#ifdef TRKRECO_DEBUG_DETAIL
704 std::cout << "... failed to make a track 3D" << std::endl;
705#endif
706 retryList.append(base);
707 segmentList[outerMost].remove(base);
708 delete t;
709 continue;
710 }
711
712 //...Salvaging...
713// _builder->salvage(* t, unusedLinks);
714
715 //...OK...
716 t->assign(WireHitConformalFinder);
717 t->finder(TrackOldConformalFinder);
718// TrackOldConformalFinder | TrackValid | Track3D);
719 _tracks.append(t);
720
721 //...Remove used links...
722 const AList<TMLink> & usedLinks = t->links();
723 list.remove(usedLinks);
724 unusedLinks.remove(usedLinks);
725 _unusedStereoConfLinks.remove(usedLinks);
726 for (unsigned i = 0; i <= outerMost; i++)
727 segmentList[i].remove(clusters);
728
729 //...For debug...
730 if (debugLevel() > 1) {
731 std::cout << name() << " ... track # " << _tracks.length() - 1;
732 std::cout << " found" << std::endl;
733 t->dump("detail", " ");
734 }
735 }
736
737 //...Loop termination...
738 --outerMost;
739 }
740
741 //...Termination...
742 for (unsigned i = 0; i < 5; i++) {
743 HepAListDeleteAll(original[i]);
744 delete segments[i];
745 }
746}
747
748void
749TConformalFinder0::specialFinding(AList<TMLink> & list,
750 AList<TMLink> & unusedLinks,
751 double fraction) {
752#ifdef TRKRECO_DEBUG_DETAIL
753 std::cout << name() << " ... standard finding with salvage : given hits :";
754 Dump(list, "sort");
755#endif
756
757 //...Find segments...
758 AList< AList<TSegment0> > segments = findSegments(list);
759 AList<TSegment0> segmentList[5];
760 AList<TSegment0> original[5];
761 for (unsigned i = 0; i < 5; i++) {
762 segmentList[i] = * segments[i];
763 original[i] = * segments[i];
764 }
765
766#ifdef TRKRECO_DEBUG_DETAIL
767 for (unsigned i = 0; i < 5; i++) {
768 std::cout << "... clusters in super layer " << i << std::endl;
769 for (unsigned j = 0; j < segmentList[i].length(); j++) {
770 segmentList[i][j]->dump("", " ");
771 segmentList[i][j]->dump("hits mc", " ");
772 }
773 }
774#endif
775
776 //...Main loop...
777 AList<TSegment0> retryList;
778 unsigned outerMost = 4;
779 while (outerMost) {
780
781 while (TSegment0 * base = segmentList[outerMost][0]) {
782
783 //...Get linked clusters...
784 AList<TSegment0> clusters = findClusterLink(* base, segmentList);
785
786 again:;
787
788 //...Check # of clusters...
789 if (clusters.length() < 2) {
790 segmentList[outerMost].remove(base);
791 continue;
792 }
793
794 //...Make a track...
795 TTrack * t = makeTrack(clusters); // don't change 'clusters'
796 if (t == NULL) {
797 clusters.remove(clusters.last());
798 goto again;
799 }
800
801 //...Check track quality...
802 // double f = float(t->nLinks()) / float(nTLinks(clusters));
803 double f = float(t->nCores()) / float(NCoreLinks(clusters));
804 if (f < fraction) {
805#ifdef TRKRECO_DEBUG_DETAIL
806 std::cout << "... fraction too low:" << f << std::endl;
807 std::cout << " retry later" << std::endl;
808#endif
809 delete t;
810 clusters.remove(clusters.last());
811 goto again;
812 }
813
814 //...Append other hits...
815 appendClusters2(* t, retryList);
816
817 //...Make it 3D...
818
819 TTrack * ts = t;
820 if (_doStereo) {
821 ts = _builder->buildStereo(* t,
822 findCloseHits(_unusedStereoConfLinks, * t));
823 }
824
825 //...Check track quality...
826 if (! ts) {
827 clusters.remove(clusters.last());
828 delete t;
829 goto again;
830 }
831
832 //...Salvaging...
833// _builder->salvage(* t, unusedLinks);
834
835 //...OK...
836 t->assign(WireHitConformalFinder);
837 t->finder(TrackOldConformalFinder);
838// TrackOldConformalFinder | TrackValid | Track3D);
839 _tracks.append(t);
840
841 //...Remove used links...
842 const AList<TMLink> & usedLinks = t->links();
843 list.remove(usedLinks);
844 unusedLinks.remove(usedLinks);
845 _unusedStereoConfLinks.remove(usedLinks);
846 for (unsigned i = 0; i <= outerMost; i++)
847 segmentList[i].remove(clusters);
848
849 //...For debug...
850 if (debugLevel() > 1) {
851 std::cout << name() << " ... track # " << _tracks.length() - 1;
852 std::cout << " found" << std::endl;
853 t->dump("detail", " ");
854 }
855 }
856
857 //...Loop termination...
858 --outerMost;
859 }
860
861 //...Termination...
862 for (unsigned i = 0; i < 5; i++) {
863 HepAListDeleteAll(original[i]);
864 delete segments[i];
865 }
866}
867
868AList< AList<TSegment0> >
869TConformalFinder0::findSegments(const AList<TMLink> & in) const {
870 AList< AList<TSegment0> > a;
871
872#ifdef TRKRECO_DEBUG_DETAIL
873 std::cout << name() << " ... finding segments : given hits =" << std::endl;
874 Dump(in, "sort");
875#endif
876
877 //...Create lists of links for each super layer...
878 AList<TMLink> links[5];
879 unsigned n = in.length();
880 for (unsigned i = 0; i < n; i++) {
881 TMLink & l = * in[i];
882// links[l.wire()->superLayerId() / 2].append(l);
883 links[l.wire()->axialStereoLayerId()/4].append(l);
884 }
885
886 //...Create phi hists and clusters for each super layer...
887 THistogram * hist[5];
888 hist[0] = new THistogram(76);
889 hist[1] = new THistogram(100);
890 hist[2] = new THistogram(128);
891 hist[3] = new THistogram(256);
892 hist[4] = new THistogram(288);
893 for (unsigned i = 0; i < 5; i++) {
894 hist[i]->fillX(links[i]);
895 AList<TSegment0> * b = new AList<TSegment0>();
896 a.append(b);
897 b->append(findClusters(* hist[i]));
898 delete hist[i];
899 }
900
901 return a;
902}
w
double w
Definition: EvtPhokharaDef.hh:26
M_PI
#define M_PI
Definition: TConstant.h:4
ORIGIN
const HepPoint3D ORIGIN
Constants.
Definition: TMDCUtil.cxx:47
WireHitConformalFinder
#define WireHitConformalFinder
Definition: TMDCWireHit.h:51
NCoreLinks
unsigned NCoreLinks(const CAList< TSegment > &list)
returns # of core links in segments.
Definition: TSegment.cxx:947
TrackOldConformalFinder
#define TrackOldConformalFinder
Definition: TTrack.h:23
t
TTree * t
Definition: binning.cxx:23
Helix::center
const HepPoint3D & center(void) const
returns position of helix center(z = 0.);
Definition: Reconstruction/TrackUtil/TrackUtil-00-00-08/TrackUtil/Helix.h:203
Helix::radius
double radius(void) const
returns radious of helix.
Definition: Reconstruction/TrackUtil/TrackUtil-00-00-08/TrackUtil/Helix.h:215
TBuilder0
A class to build a track.
Definition: TBuilder0.h:35
TBuilder0::buildRphi
TTrack * buildRphi(const AList< TMLink > &) const
builds a r/phi track from TMLinks or from Segments.
Definition: TBuilder0.cxx:83
TBuilder0::trackSelector
const TMSelector & trackSelector(void) const
returns a track selector.
Definition: TBuilder0.h:131
TBuilder0::buildStereo
virtual TTrack * buildStereo(TTrack &track, const AList< TMLink > &) const
appends stereo hits to a track.
Definition: TBuilder0.cxx:535
TBuilder0::appendClusters
void appendClusters(TTrack &track, const AList< TMLink > &) const
appends TMLinks in a list.
Definition: TBuilder0.cxx:1417
TBuilderCosmic
A class to build a cosmic track.
Definition: TBuilderCosmic.h:36
TConformalFinder0::clear
void clear(void)
clear internal information.
Definition: TConformalFinder0.cxx:89
TConformalFinder0::findSegments
AList< AList< TSegment0 > > findSegments(const AList< TMLink > &in) const
finds segments.
Definition: TConformalFinder0.cxx:869
TConformalFinder0::findClusters2
AList< TSegment0 > findClusters2(const THistogram &) const
Definition: TConformalFinder0.cxx:243
TConformalFinder0::conformalTransformationRphi
static void conformalTransformationRphi(const HepPoint3D &center, const AList< TMDCWireHit > &hits, AList< TMLink > &links)
transforms hits into a conformal plane. 'center' is a center of the transformation....
Definition: TConformalFinder0.cxx:165
TConformalFinder0::conformalTransformation
static void conformalTransformation(const HepPoint3D &center, const AList< TMDCWireHit > &hits, AList< TMLink > &links)
transforms hits into a conformal plane. 'center' is a center of the transformation....
Definition: TConformalFinder0.cxx:100
TConformalFinder0::dump
void dump(const std::string &message=std::string(""), const std::string &prefix=std::string("")) const
dumps debug information.
Definition: TConformalFinder0.cxx:78
TConformalFinder0::name
std::string name(void) const
returns name.
Definition: TConformalFinder0.h:214
TConformalFinder0::~TConformalFinder0
virtual ~TConformalFinder0()
Destructor.
Definition: TConformalFinder0.cxx:73
TConformalFinder0::doit
int doit(const AList< TMDCWireHit > &axialHits, const AList< TMDCWireHit > &stereoHits, AList< TTrack > &tracks, AList< TTrack > &tracks3D)
finds tracks.
Definition: TConformalFinder0.cxx:560
TConformalFinder0::version
std::string version(void) const
returns version.
Definition: TConformalFinder0.cxx:29
TConformalFinder0::findClusters
AList< TSegment0 > findClusters(const THistogram &) const
finds segments. (obsolete functions)
Definition: TConformalFinder0.cxx:195
TConformalFinder0::TConformalFinder0
TConformalFinder0(float maxSigma, float fraction, float stereoZ3, float stereoZ4, float stereoChisq3, float stereoChisq4, float stereoMaxSigma, unsigned fittingCorrections, float salvageLevel, bool cosmic)
Constructor.
Definition: TConformalFinder0.cxx:33
TConformalFinder0::conformalTransformationDriftCircle
static void conformalTransformationDriftCircle(const HepPoint3D &center, const AList< TMDCWireHit > &hits, AList< TMLink > &links)
transforms drift circle of hits into a conformal plane. transformed positions( x0,...
Definition: TConformalFinder0.cxx:127
TFinderBase
A virtual class for a track finder in tracking.
Definition: TFinderBase.h:31
TFinderBase::debugLevel
virtual int debugLevel(void) const
returns debug level.
Definition: TFinderBase.h:90
TFinderBase::dump
virtual void dump(const std::string &message=std::string(""), const std::string &prefix=std::string("")) const
dumps debug information.
Definition: TFinderBase.cxx:23
THistogram
A class for a histogram used in tracking.
Definition: THistogram.h:34
THistogram::clusters0
AList< TSegment0 > clusters0(void) const
returns an AList<TSegment0> of clusters.
Definition: THistogram.cxx:162
THistogram::fillX
void fillX(const AList< TMLink > &links)
fills with hits.
Definition: THistogram.cxx:76
TMDCWireHit::drift
float drift(unsigned) const
returns drift distance.
Definition: TMDCWireHit.h:236
TMDCWireHit::xyPosition
const HepPoint3D & xyPosition(void) const
returns drift time
Definition: TMDCWireHit.h:262
TMDCWire
A class to represent a wire in MDC.
Definition: TMDCWire.h:55
TMDCWire::axialStereoLayerId
unsigned axialStereoLayerId(void) const
returns id of axial or stereo id.
Definition: TMDCWire.h:360
TMSelector::maxDistance
double maxDistance(void) const
returns max. distance required for stereo hits.
Definition: TMSelector.h:235
TMSelector::nSuperLayers
unsigned nSuperLayers(void) const
returns min. # of super layers required.
Definition: TMSelector.h:150
TMSelector::maxImpact
double maxImpact(void) const
returns max. impact(2D) required.
Definition: TMSelector.h:184
TMSelector::nLinks
unsigned nLinks(void) const
returns min. # of hits(TMLinks) requried.
Definition: TMSelector.h:133
TMSelector::maxSigma
double maxSigma(void) const
returns max. sigma for each TMLink.
Definition: TMSelector.h:201
TMSelector::nLinksStereo
unsigned nLinksStereo(void) const
returns min. # of stereo hits(TMLinks) requried.
Definition: TMSelector.h:218
TMSelector::minPt
double minPt(void) const
returns min. pt required.
Definition: TMSelector.h:167
TSegment0
A class to relate TMDCWireHit and TTrack objects.
Definition: TSegment0.h:41
TSegment0::dump
void dump(const std::string &message=std::string(""), const std::string &prefix=std::string("")) const
dumps debug information.
Definition: TSegment0.cxx:51
TSegment0::solveDualHits
int solveDualHits(void)
Definition: TSegment0.cxx:624
TSegment0::distance
double distance(const TSegment0 &) const
calculates distance between two clusters. Smaller value indicates closer.
Definition: TSegment0.cxx:131
TSegment0::split
AList< TSegment0 > split(void) const
returns a list of sub TSegments in this cluster. If cluster type is 1, 2, or 7, no cluster is returne...
Definition: TSegment0.cxx:281
TSegment0::direction
const HepVector3D & direction(void) const
returns direction.
Definition: TSegment0.h:179
TSegment0::position
const HepPoint3D & position(void) const
returns position.
Definition: TSegment0.h:172
TSegment0::clusterType
unsigned clusterType(void) const
returns cluster type. 0:empty, 1:short line, 2:long line, 3:V shage(from outside),...
Definition: TSegment0.h:214
TTrackBase::links
const AList< TMLink > & links(unsigned mask=0) const
returns a list of masked TMLinks assigned to this track. 'mask' will be applied if mask is not 0.
Definition: TTrackBase.cxx:297
TTrackBase::nLinks
unsigned nLinks(unsigned mask=0) const
returns # of masked TMLinks assigned to this track object.
Definition: TTrackBase.cxx:305
TTrack
A class to represent a track in tracking.
Definition: TTrack.h:129
TTrack::helix
const Helix & helix(void) const
returns helix parameter.
Definition: TTrack.h:477
f
TFile f("ana_bhabha660a_dqa_mcPat_zy_old.root")
close
file close()
EvtCyclic3::next
Index next(Index i)
Definition: EvtCyclic3.cc:107
b
const double b
Definition: slope.cxx:9