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