13#include "TrkReco/TLine0.h"
14#include "TrkReco/TMDCUtil.h"
15#include "TrkReco/TMDCWire.h"
16#include "TrkReco/TMDCWireHit.h"
17#include "TrkReco/TMDCWireHitMC.h"
18#include "TrkReco/TTrackHEP.h"
27 , _fittedUpdated( false )
29 , _reducedChi2( 0. ) {
32 fitter( &TLine0::_fitter );
40 , _fittedUpdated( false )
42 , _reducedChi2( 0. ) {
45 fitter( &TLine0::_fitter );
50void TLine0::dump(
const std::string& msg,
const std::string& pre )
const {
52 if ( msg ==
"" ) def =
true;
54 if ( def || msg.find(
"line" ) != std::string::npos ||
55 msg.find(
"detail" ) != std::string::npos )
58 std::cout <<
"#links=" <<
_links.length();
59 std::cout <<
",a=" << _a;
60 std::cout <<
",b=" << _b;
61 std::cout <<
",det=" << _det;
62 std::cout << std::endl;
112#ifdef TRKRECO_DEBUG_DETAIL
113 if ( !
_fitted ) std::cout <<
"TLine0::chi2 !!! fit not performed" << std::endl;
116 if ( _fittedUpdated )
return _chi2;
119 for (
unsigned i = 0; i <
n; i++ )
125 double c = y - _a * x - _b;
128 _fittedUpdated =
true;
135 for (
unsigned i = 0; i <
n; i++ )
139 if ( dist > maxSigma ) bad.append( l );
142#ifdef TRKRECO_DEBUG_DETAIL
143 std::cout <<
" TLine0::refine ... rejected hits:max distance=" << maxSigma;
144 std::cout << std::endl;
146 for (
unsigned i = 0; i < bad.length(); i++ )
156 else std::cout <<
"0";
160 if (
distance( l ) > maxSigma ) std::cout <<
" X";
161 std::cout << std::endl;
170 _fittedUpdated =
false;
178 int mask[100] = { 0 };
179 int nsl[11] = { 64, 80, 96, 128, 144, 160, 176, 208, 240, 256, 288 };
180 int npos = 0,
nneg = 0;
181 for (
unsigned i = 0; i <
n - 1; i++ )
184 for (
unsigned j = i + 1; j <
n; j++ )
190 if ( i > 0 && ( mask[i - 1] == 1 && mask[j] == 1 ) )
193 if ( l.
hit()->
wire()->
layerId() ==
t.hit()->wire()->layerId() ) { mask[i] = 1; }
197 int jlocal =
s.hit()->wire()->localId();
198 if ( ilocal > 0 && ilocal < ilast )
200 if (
abs( jlocal - ilocal ) > 1 )
206 else if ( ilocal == 0 )
208 if ( jlocal > 1 && jlocal < ilast )
214 else if ( ilocal == ilast )
216 if ( jlocal > 0 && jlocal < ilast - 1 )
227 if ( l.
position().y() >= 0 ) npos += 1;
233 double sumX = 0., sumY = 0., sumX2 = 0., sumXY = 0., sumY2 = 0.;
236 for (
unsigned i = 0; i <
n; i++ )
239 if ( mask[i] == 1 )
continue;
247 if ( npos >
nneg && y < 0 )
continue;
248 if ( npos < nneg && y > 0 )
continue;
259 if ( nused < 2 || ( nused == 2 && lyid[0] == lyid[1] ) ) {
return -2; }
260 double sum = double( nused );
261 _det = sum * sumX2 - sumX * sumX;
262 if ( _det == 0. ) {
return -1; }
263 _a = ( sumXY * sum - sumX * sumY ) / _det;
264 _b = ( sumX2 * sumY - sumX * sumXY ) / _det;
274 int mask[100] = { 0 };
275 int npos = 0,
nneg = 0;
276 for (
unsigned i = 0; i <
n - 1; i++ )
279 for (
unsigned j = i + 1; j <
n; j++ )
291 if ( l.
position().y() >= 0 ) npos += 1;
297 double sumX = 0., sumY = 0., sumX2 = 0., sumXY = 0., sumY2 = 0.;
300 for (
unsigned i = 0; i <
n; i++ )
303 if ( mask[i] == 1 )
continue;
309 if ( npos >
nneg && y < 0 )
continue;
310 if ( npos < nneg && y > 0 )
continue;
320 if ( nused < 4 ) {
return -2; }
321 double sum = double( nused );
322 _det = sum * sumX2 - sumX * sumX;
323 if ( _det == 0. ) {
return -1; }
324 _a = ( sumXY * sum - sumX * sumY ) / _det;
325 _b = ( sumX2 * sumY - sumX * sumXY ) / _det;
334 int mask[100] = { 0 };
337 double Crad = 180. / 3.141592;
338 for (
unsigned i = 0; i <
n - 1; i++ )
341 for (
unsigned j = i + 1; j <
n; j++ )
360 if ( mask[
n - 1] != 1 )
371 phi_max = phi_ave /
n + 40;
372 phi_min = phi_ave /
n - 40;
376 double sumX = 0., sumY = 0., sumX2 = 0., sumXY = 0., sumY2 = 0.;
379 for (
unsigned i = 0; i <
n; i++ )
382 if ( mask[i] == 1 )
continue;
389 if ( phi > phi_max && phi < phi_min )
continue;
399 if ( nused < 4 ) {
return -2; }
400 double sum = double( nused );
401 _det = sum * sumX2 - sumX * sumX;
402 if ( _det == 0. ) {
return -1; }
403 _a = ( sumXY * sum - sumX * sumY ) / _det;
404 _b = ( sumX2 * sumY - sumX * sumXY ) / _det;
414 int mask[100] = { 0 };
415 int nsl[11] = { 64, 80, 96, 128, 144, 160, 176, 208, 240, 256, 288 };
418 double Crad = 180. / 3.141592;
419 for (
unsigned i = 0; i <
n - 1; i++ )
422 for (
unsigned j = i + 1; j <
n; j++ )
428 if ( i > 0 && ( mask[i - 1] == 1 && mask[j] == 1 ) )
431 if ( l.
hit()->
wire()->
layerId() ==
t.hit()->wire()->layerId() ) { mask[i] = 1; }
435 int jlocal =
s.hit()->wire()->localId();
436 if ( ilocal > 0 && ilocal < ilast )
438 if (
abs( jlocal - ilocal ) > 1 )
444 else if ( ilocal == 0 )
446 if ( jlocal > 1 && jlocal < ilast )
452 else if ( ilocal == ilast )
454 if ( jlocal > 0 && jlocal < ilast - 1 )
473 if ( mask[
n - 1] != 1 )
484 phi_max = phi_ave /
n + 40;
485 phi_min = phi_ave /
n - 40;
489 double sumX = 0., sumY = 0., sumX2 = 0., sumXY = 0., sumY2 = 0.;
492 for (
unsigned i = 0; i <
n; i++ )
495 if ( mask[i] == 1 )
continue;
502 if ( phi > phi_max && phi < phi_min )
continue;
513 if ( nused < 2 || ( nused == 2 && lyid[0] == lyid[1] ) ) {
return -2; }
514 double sum = double( nused );
515 _det = sum * sumX2 - sumX * sumX;
516 if ( _det == 0. ) {
return -1; }
517 _a = ( sumXY * sum - sumX * sumY ) / _det;
518 _b = ( sumX2 * sumY - sumX * sumXY ) / _det;
527 int nlyr[50] = { 0 };
528 int nneg = 0, npos = 0;
529 for (
unsigned i = 0; i <
n - 1; i++ )
539 for (
unsigned i = 0; i <
n; i++ )
553 if ( npos >
nneg && l.
position().y() < 0 ) bad.append( l );
554 if ( npos <
nneg && l.
position().y() > 0 ) bad.append( l );
558 if ( bad.length() > 0 && bad.length() <
n ) {
_links.remove( bad ); }
562 _fittedUpdated =
false;
568 _fittedUpdated =
false;
574 _fittedUpdated =
false;
580 unsigned nb =
_links.length();
583 unsigned n = list.length();
584 for (
unsigned i = 0; i <
n; i++ )
590 if ( dist < maxSigma ) {
_links.append( l ); }
594 unsigned na =
_links.length();
599 for (
unsigned i = 0; i < na; i++ )
612 if ( bad.length() > 0 )
_links.remove( bad );
614 _fittedUpdated =
false;
619#ifdef TRKRECO_DEBUG_DETAIL
620 if ( !
_fitted ) std::cout <<
"TLine0::reducedChi2 !!! fit not performed" << std::endl;
623 if ( _fittedUpdated )
return _reducedChi2;
627 for (
unsigned i = 0; i <
n; i++ )
633 double c = y - _a * x - _b;
635 chi2 += c * c / err / err;
638 _reducedChi2 =
chi2 / (
n - 2 );
639 _fittedUpdated =
true;
int SortByWireId(const void *a, const void *b)
Sorter.
void appendSLY(AList< TMLink > &list)
double chi2(void) const
returns chi2.
double distance(const TMLink &) const
returns distance to a position of TMLink itself. (not to a wire)
void appendByszdistance(AList< TMLink > &list, unsigned isl, float maxSigma)
double reducedChi2(void) const
returns reduced-chi2.
int fit2()
fits itself. Error was happened if return value is not zero.
virtual ~TLine0()
Destructor.
void removeChits()
remove extremly bad points.
void removeSLY(AList< TMLink > &list)
int fit2p()
fits itself using isolated hits. Error was happened if return value is not zero.
void refine(AList< TMLink > &list, float maxSigma)
double a(void) const
returns coefficient a.
void dump(const std::string &message=std::string(""), const std::string &prefix=std::string("")) const
dumps debug information.
A class to fit a TTrackBase object to a line.
const TTrackHEP *const hep(void) const
returns a pointer to a GEN_HEPEVT.
float dDrift(unsigned) const
returns drift distance error.
const TMDCWire *const wire(void) const
returns a pointer to a TMDCWire.
const TMDCWireHitMC *const mc(void) const
returns a pointer to TMDCWireHitMC.
unsigned localId(void) const
returns local id in a wire layer.
unsigned layerId(void) const
returns layer id.
unsigned superLayerId(void) const
returns super layer id.
A class to relate TMDCWireHit and TTrack objects.
const TMDCWireHit * hit(void) const
returns a pointer to a hit.
const HepPoint3D & position(void) const
returns position.
const TMDCWire *const wire(void) const
returns a pointer to a wire.
virtual void dump(const std::string &message=std::string(""), const std::string &prefix=std::string("")) const
dumps debug information.
const TMFitter *const fitter(void) const
returns a pointer to a default fitter.
unsigned id(void) const
returns an id started from 0.