TBuilderCosmic.cxx

Go to the documentation of this file.

//-----------------------------------------------------------------------------
// $Id: TBuilderCosmic.cxx,v 1.6 2010/03/31 09:58:59 liucy Exp $
//-----------------------------------------------------------------------------
// Filename : TBuilderCosmic.cc
// Section  : Tracking
// Owner    : Yoshi Iwasaki
// Email    : [email protected]
//-----------------------------------------------------------------------------
// Description : A class to build a cosmic track.
//               See http://bsunsrv1.kek.jp/~yiwasaki/tracking/
//-----------------------------------------------------------------------------
 
#ifdef TRKRECO_DEBUG_DETAIL
#ifndef TRKRECO_DEBUG
#define TRKRECO_DEBUG
#endif
#endif
#include "TrkReco/TBuilderCosmic.h"
#include "TrkReco/TMLink.h"
#include "TrkReco/TTrack.h"
#include "TrkReco/TLine0.h"
#include "TrkReco/T3DLine.h"
 
TBuilderCosmic::TBuilderCosmic(const std::string & name, float salvageLevel)
: TBuilder0(name, salvageLevel), _fitter("TBuilderCosmic Fitter") {
}
 
TBuilderCosmic::~TBuilderCosmic() {
}
 
TTrack *
TBuilderCosmic::buildStereo(TTrack & track, const AList<TMLink> & list) const {
#ifdef TRKRECO_DEBUG_DETAIL
    std::cout << name() << "(stereo) ... dump of stereo candidate hits" << std::endl;
    AList<TMLink> tmp = list;
    tmp.sort(SortByWireId);
    std::cout << "    ";
    for (unsigned i = 0; i < tmp.length(); i++) {
    TMLink * l = tmp[i];
    std::cout << l->wire()->layerId() << "-";
    std::cout << l->wire()->localId() << ",";
    }
    std::cout << std::endl;
#endif
 
    //...Check # of links...
    if (list.length() < _lineSelector.nLinksStereo()) {
#ifdef TRKRECO_DEBUG_DETAIL
    std::cout << name() << "(stereo) ... rejected by nLinks(";
    std::cout << list.length() << ") < ";
    std::cout << _lineSelector.nLinks() << std::endl;
#endif  
    return NULL;
    }
 
    //...Calculate s and z for every links...
    unsigned n = list.length();
    AList<TMLink> forLine;
    for (unsigned i = 0; i < n; i++) {
    TMLink * l = list[i];
 
        //... Require Fitting vaildation 
        if(!(l->hit()->state()& WireHitFittingValid)) continue;
 
    TMLink * t = new TMLink(* l);
 
    //...Assuming wire position...
    t->leftRight(2);
    int err = track.szPosition(* t);
    if (err) {
        delete t;
        continue;
    }
 
    //...Store the sz link...
    t->link(l);
    forLine.append(t);
    }
 
#ifdef TRKRECO_DEBUG_DETAIL
    std::cout << name() << "(stereo) ... dump of sz links" << std::endl;
    std::cout << "    ";
    tmp = forLine;
    tmp.sort(SortByWireId);
    for (unsigned i = 0; i < tmp.length(); i++) {
    TMLink * l = tmp[i];
    std::cout << l->wire()->layerId() << "-";
    std::cout << l->wire()->localId() << ",";
    }
    std::cout << std::endl;
#endif  
    
    //...Check # of sz links...
    if (forLine.length() < _lineSelector.nLinksStereo()) {
#ifdef TRKRECO_DEBUG_DETAIL
    std::cout << name() << "(stereo) ... rejected by sz nLinks(";
    std::cout << forLine.length() << ") < ";
    std::cout << _lineSelector.nLinks() << std::endl;
#endif  
    HepAListDeleteAll(forLine);
    return NULL;
    }
 
    //...Make a line...
    unsigned nLine = forLine.length();
    TLine0 line(forLine);
    int err = line.fit();
 
    //...Linear fit...
    if (err < 0) {
#ifdef TRKRECO_DEBUG_DETAIL
    std::cout << name() << "(stereo) ... linear fit failure. nLinks(";
    std::cout << forLine.length() << ")" << std::endl;
#endif  
    HepAListDeleteAll(forLine);
    return NULL;
    }
 
#ifdef TRKRECO_DEBUG_DETAIL
    std::cout << name() << "(stereo) ... dump of left-right" << std::endl;
#endif
 
    //...Decide Left or Right...
    AList<TMLink> forNewLine;
    for (unsigned i = 0; i < nLine; i++) {
    TMLink * t = forLine[i];
    TMLink * tl = new TMLink(* t);
    TMLink * tr = new TMLink(* t);
 
    tl->leftRight(WireHitLeft);
    tr->leftRight(WireHitRight);
 
    int err = track.szPosition(* tl);
    if (err) {
        delete tl;
        tl = NULL;
    }
    err = track.szPosition(* tr);
    if (err) {
        delete tr;
        tr = NULL;
    }
    if ((tl == NULL) && (tr == NULL)) continue;
 
    TMLink * best;
    if (tl == NULL) best = tr;
    else if (tr == NULL) best = tl;
    else {
        if (line.distance(* tl) < line.distance(* tr)) {
        best = tl;
        delete tr;
        }
        else {
        best = tr;
        delete tl;
        }
    }
 
#ifdef TRKRECO_DEBUG_DETAIL
    std::cout << "    ";
    std::cout << t->wire()->layerId() << "-";
    std::cout << t->wire()->localId();
    if (tl != NULL)
        std::cout << ",left " << tl->position() << "," << line.distance(* tl);
    if (tr != NULL)
        std::cout << ",right " << tr->position() << "," << line.distance(* tr);
    std::cout << std::endl;
#endif
 
    best->link(t->link());
    forNewLine.append(best);
    }
 
    //...Check # of sz links...
    if (forNewLine.length() < _lineSelector.nLinksStereo()) {
#ifdef TRKRECO_DEBUG_DETAIL
    std::cout << name() << "(stereo) ... rejected by lr nLinks(";
    std::cout << forNewLine.length() << ") < ";
    std::cout << _lineSelector.nLinks() << std::endl;
#endif  
    HepAListDeleteAll(forLine);
    HepAListDeleteAll(forNewLine);
    return NULL;
    }
 
    //...Create new line...
#ifdef TRKRECO_DEBUG_DETAIL
    std::cout << name() << "(stereo) ... creating a new line" << std::endl;
#endif  
    unsigned nNewLine = forNewLine.length();
    TLine0 newLine(forNewLine);
 
    //...Make a seed track again
    err = newLine.fit();
 
    //...Linear fit...
    if (err < 0) {
#ifdef TRKRECO_DEBUG_DETAIL
    std::cout << name() << "(stereo) ... 2nd linear fit failure. nLinks(";
    std::cout << forNewLine.length() << ")" << std::endl;
#endif  
    HepAListDeleteAll(forLine);
    HepAListDeleteAll(forNewLine);
    return NULL;
    }
 
    //...Remove bad points...
    AList<TMLink> bad;
//    newLine.refine(bad, 40.);  //Liuqg, meaningless while without magnetic field
//    err = newLine.fit();
//    newLine.refine(bad, 20.);
//    err = newLine.fit();
//    newLine.refine(bad, 10.);
//    err = newLine.fit();
 
    //...Linear fit again...
    if (err < 0) {
    HepAListDeleteAll(forLine);
    HepAListDeleteAll(forNewLine);
#ifdef TRKRECO_DEBUG_DETAIL
    std::cout << "    appendStereo cut ... new line 2nd linear fit failure. ";
    std::cout << "# of links = " << n << "," << nLine;
    std::cout << "," << nNewLine << std::endl;
#endif  
    return NULL;
    }
 
    //...3D fit...
    const AList<TMLink> & good = newLine.links();
    unsigned nn = good.length();
    for (unsigned i = 0; i < nn; i++) {
    track.append(* good[i]->link());
    }
    Vector a(5);
    a = track.helix().a();
    a[3] = newLine.b();
    a[4] = track.charge() * newLine.a();
    track._helix->a(a);
 
#ifdef LINE_COSMIC
    T3DLine * Ltrack = new T3DLine(track);
    
    //...Refine...
    err = _fitter.fit(*Ltrack);
    Ltrack->refine(bad, _trackSelector.maxSigma() * 30.);
    err = _fitter.fit(*Ltrack);
    Ltrack->refine(bad, _trackSelector.maxSigma() * 3.);
    err = _fitter.fit(*Ltrack);
//    Ltrack->refine(bad, _trackSelector.maxSigma() * 0.21);  //liuqg, for prelimilary test
    track.refine(bad, _trackSelector.maxSigma()); 
    err = _fitter.fit(*Ltrack);
 
    //...Test it...
/*    if (! _trackSelector.select(*Ltrack)) {
    HepAListDeleteAll(forLine);
    HepAListDeleteAll(forNewLine);
    delete Ltrack;
    return NULL;
    }
*/
    //...Termination...
    HepAListDeleteAll(forLine);
    HepAListDeleteAll(forNewLine);
    track.removeLinks();
    track.append(Ltrack->links());
    Vector a1(5);
    a1 = Ltrack->helix().a();
    track._helix->a(a1);
    delete Ltrack;
    return & track;
#endif
 
    //...Refine...
    err = _fitter.fit(track);
    track.refine(bad, _trackSelector.maxSigma() * 30.);
    err = _fitter.fit(track);
    track.refine(bad, _trackSelector.maxSigma() * 3.);
    err = _fitter.fit(track);
    track.refine(bad, _trackSelector.maxSigma());  //liuqg, for prelimilary test
//    track.refine(bad, _trackSelector.maxSigma() * 0.21); 
    err = _fitter.fit(track);
    //...Test it...
    if (! _trackSelector.select(track)) {
    HepAListDeleteAll(forLine);
    HepAListDeleteAll(forNewLine);
    return NULL;
    }
 
    //...Termination...
    HepAListDeleteAll(forLine);
    HepAListDeleteAll(forNewLine);
    return & track;
}