/home/bes3soft/bes3soft/Boss/7.0.2/dist/7.0.2/Reconstruction/TrkReco/TrkReco-00-08-59-patch4-slc6tag/src/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    : yoshihito.iwasaki@kek.jp
//-----------------------------------------------------------------------------
// 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;
}

---