TRobustLineFitter Class Reference

A class to fit a TTrackBase object to a line. More...

#include <TRobustLineFitter.h>

Inheritance diagram for TRobustLineFitter:

List of all members.

Public Member Functions
	TRobustLineFitter (const std::string &name)
	Constructor.
virtual	~TRobustLineFitter ()
	Destructor.
void	dump (const std::string &message=std::string(""), const std::string &prefix=std::string("")) const
	dumps debug information.
double	a (void) const
double	b (void) const
double	det (void) const
virtual int	fit (TTrackBase &) const
const std::string &	name (void) const
	returns name.
Protected Member Functions
void	fitDone (TTrackBase &) const
	sets the fitted flag. (Bad implementation)
Private Member Functions
double	rofunc (const TTrackBase &, double) const
double	select (unsigned k, unsigned n, double *) const
Private Attributes
double	_a
double	_b
double	_det
unsigned	_n

---

Detailed Description

A class to fit a TTrackBase object to a line.

Definition at line 28 of file TRobustLineFitter.h.

---

Constructor & Destructor Documentation

TRobustLineFitter::TRobustLineFitter

(

const std::string &

name

)

Constructor.

Definition at line 37 of file TRobustLineFitter.cxx.

: TLineFitter(name), _a(0.), _b(0.), _det(0.) {
}

TRobustLineFitter::~TRobustLineFitter

(

)

[virtual]

Destructor.

Definition at line 41 of file TRobustLineFitter.cxx.

                                      {
}

---

Member Function Documentation

double TRobustLineFitter::a

(

void

)

const [inline]

Reimplemented from TLineFitter.

Definition at line 72 of file TRobustLineFitter.h.

References _a.

Referenced by select().

                               {
    return _a;
}

double TRobustLineFitter::b

(

void

)

const [inline]

Reimplemented from TLineFitter.

Definition at line 78 of file TRobustLineFitter.h.

References _b.

                               {
    return _b;
}

double TRobustLineFitter::det

(

void

)

const [inline]

Reimplemented from TLineFitter.

Definition at line 84 of file TRobustLineFitter.h.

References _det.

                                 {
    return _det;
}

void TRobustLineFitter::dump	(	const std::string &	message = std::string(""),
		const std::string &	prefix = std::string("")
	)			const

dumps debug information.

Reimplemented from TLineFitter.

int TRobustLineFitter::fit

(

TTrackBase &

)

const [virtual]

Reimplemented from TLineFitter.

Definition at line 45 of file TRobustLineFitter.cxx.

References _a, _b, _det, _n, TLineFitter::a(), TLineFitter::b(), TLineFitter::det(), showlog::err, f1, f2, TLineFitter::fit(), TMFitter::fitDone(), genRecEmupikp::i, Line, rofunc(), and t().

Referenced by TBuilderCurl::buildStereo(), TBuilder::buildStereo(), TBuilder::initialLine2(), and TBuilder::searchLine().

                                           {
    //...Initial guess...
    int err = TLineFitter::fit(t);
    if (err) return err;

    //...Check # of hits...
    const AList<TMLink> & links = t.links();
    _n = links.length();
    if (_n < 3) return 0;

    //...Standard deviation...
    _a = TLineFitter::a();
    _b = TLineFitter::b();
    _det = TLineFitter::det();
    double chisq = 0.;
    for (unsigned i = 0; i < _n; i++) {
        const HepPoint3D & p = links[i]->position();
        double tmp = p.y() - (_a * p.x() + _b);
        chisq += tmp * tmp;
    }
    double siga = sqrt(chisq / _det);

    //...Decide iteration step...
    double a1 = _a;
    double f1 = rofunc(t, a1);
    double a2 = _a + copysign(3.0 * siga, f1);
    double f2 = rofunc(t, a2);

    // if initial value f2 >f1, change the search direction
    if( f1 * f2 > 0. && fabs(f2) > fabs(f1) ){
        a2 = _a - copysign(3.0 * siga, f1);
        f2 = rofunc(t, a2);
    }

    int backwardSearch=0;
    while (f1 * f2 > 0. ) {
      _a = 2.0 * a2 - a1;
      a1 = a2;
      f1 = f2;
      a2 = _a;
      f2 = rofunc(t, a2);

      if( f1 * f2 > 0. && ( fabs(f2) > fabs(f1) || fabs(f2-f1)<0.01 ) ){
        backwardSearch++;
        if(backwardSearch==2){
          break;
        }

        double f=f2;
        a2 = a1;
        f2 = f1;
        a1 = _a;
        f1 = f;
      }
    }

    if(backwardSearch==2){
      // for case of  no zero cross
      //search minimun fabs(f)

      double siga1 = 0.01 * siga;
      double a21(fabs(a2-a1));
      while (a21 > siga1) {
        _a = 0.5 * (a1 + a2);
        if (_a == a1 || _a == a2) break;
        double f = rofunc(t, _a);

        if( f * f1 <0 ){
          f1=f;
          a1=_a;
          backwardSearch--;
          break;
        }

        if ( fabs(f) <= fabs(f1) && fabs(f) <= fabs(f2)) {
          if( fabs(f-f1) > fabs(f-f2) ){
            f1 = f;
            a1 = _a;
          }else{
            f2 = f;
            a2 = _a;
          }
        }else if( fabs(f) <= fabs(f1) ){
          f1 = f;
          a1 = _a;
        }else if ( fabs(f) <= fabs(f2) ){
          f2 = f;
          a2 = _a;
        }else{
          if( fabs(f2) > fabs(f1) ){
            f1=f;
            a1 = _a;
          }else{
            f2 = f;
            a2 = _a;
          }
        }
        if (fabs(a2-a1) >= a21) break;
        a21 = fabs(a2-a1);
      }
    }

    if(backwardSearch<=1){

      //search zero cross
      siga = 0.01 * siga;

      double a21(fabs(a2-a1));
      while (a21 > siga) {
        _a = 0.5 * (a1 + a2);
        if (_a == a1 || _a == a2) break;
        double f = rofunc(t, _a);
        if (f * f1 >= 0.) {
          f1 = f;
          a1 = _a;
        }
        else {
          f2 = f;
          a2 = _a;
        }
        if (fabs(a2-a1) >= a21) break;
        a21 = fabs(a2-a1);
      }
    }

    _det = _det / double(_n);

    if (t.objectType() == Line)
        ((TMLine &) t).property(_a, _b, _det);
    fitDone(t);
    return 0;
}

void TMFitter::fitDone

(

TTrackBase &

)

const [protected, inherited]

sets the fitted flag. (Bad implementation)

Definition at line 24 of file TMFitter.cxx.

References t().

Referenced by fit(), TLineFitter::fit(), and TCircleFitter::fit().

                                      {
    t._fitted = true;
}

const std::string & TMFitter::name

(

void

)

const [inline, inherited]

returns name.

Definition at line 73 of file TMFitter.h.

References TMFitter::_name.

                         {
    return _name;
}

double TRobustLineFitter::rofunc	(	const TTrackBase &	,
		double
	)			const [private]

Definition at line 179 of file TRobustLineFitter.cxx.

References _b, _det, _n, genRecEmupikp::i, ganga-rec::j, select(), t(), and x.

Referenced by fit().

                                                              {
    double * arr = (double *) malloc(_n * sizeof(double));
    for (unsigned i = 0; i < _n; i++)
        arr[i] = t.links()[i]->position().y()
            - a * t.links()[i]->position().x();
    if (_n & 1) {
        _b = select((_n + 1) >> 1, _n, arr);
    }
    else {
        unsigned j = _n >> 1;
        _b = 0.5 * (select(j, _n, arr) + select(j + 1, _n, arr));
    }

    _det = 0.;
    double sum = 0.;
    for (unsigned i = 0; i < _n; i++) {
        double x = t.links()[i]->position().x();
        double y = t.links()[i]->position().y();
        double d = y - (a * x + _b);
        _det += fabs(d);
        if (y != 0.) d /= fabs(y);
        if (fabs(d) > 1.0e-7) sum += d > 0.0 ? x : - x;
    }
    free(arr);
    return sum;
}

double TRobustLineFitter::select	(	unsigned	k,
		unsigned	n,
		double *
	)			const [private]

Definition at line 209 of file TRobustLineFitter.cxx.

References a(), genRecEmupikp::i, ir, ganga-rec::j, and SWAP.

Referenced by rofunc().

                                                                    {
    --k;
    double tmp;
    unsigned l = 0;
    unsigned ir = n - 1;
    while (1) {
        if (ir <= l + 1) {
            if (ir == l + 1 && arr[ir] < arr[l]) {
                SWAP(arr[l],arr[ir]);
            }

//          std::cout << "k = " << k << std::endl;
//          for (unsigned i = 0; i < _n; i++)
//              std::cout << i << " : " << arr[i] << std::endl;

            return arr[k];
        }
        else {
            unsigned mid = (l + ir) >> 1;
            SWAP(arr[mid],arr[l+1]);
            if (arr[l + 1] > arr[ir]) {
                SWAP(arr[l+1],arr[ir]);
            }
            if (arr[l] > arr[ir]) {
                SWAP(arr[l],arr[ir]);
            }
            if (arr[l + 1] > arr[l]) {
                SWAP(arr[l+1],arr[l]);
            }
            unsigned i = l + 1;
            unsigned j = ir;
            double a = arr[l];
            while (1) {
//              do i++; while (arr[i] < a);
//              do j--; while (arr[j] > a);
//              while (arr[i] < a) ++i;
//              while (arr[j] > a) --j;
                while (arr[++i] < a);
                while (arr[--j] > a);
                if (j < i) break;
                SWAP(arr[i],arr[j]);
            }
            arr[l] = arr[j];
            arr[j] = a;
            if (j >= k) ir = j - 1;
            if (j <= k) l = i;
        }
    }
}

---

Member Data Documentation

double TRobustLineFitter::_a [mutable, private]

Reimplemented from TLineFitter.

Definition at line 53 of file TRobustLineFitter.h.

Referenced by a(), and fit().

double TRobustLineFitter::_b [mutable, private]

Reimplemented from TLineFitter.

Definition at line 54 of file TRobustLineFitter.h.

Referenced by b(), fit(), and rofunc().

double TRobustLineFitter::_det [mutable, private]

Reimplemented from TLineFitter.

Definition at line 55 of file TRobustLineFitter.h.

Referenced by det(), fit(), and rofunc().

unsigned TRobustLineFitter::_n [mutable, private]

Definition at line 56 of file TRobustLineFitter.h.

Referenced by fit(), and rofunc().

---