GFTrack.h

Go to the documentation of this file.

/* Copyright 2008-2010, Technische Universitaet Muenchen,
   Authors: Christian Hoeppner & Sebastian Neubert

   This file is part of GENFIT.

   GENFIT is free software: you can redistribute it and/or modify
   it under the terms of the GNU Lesser General Public License as published
   by the Free Software Foundation, either version 3 of the License, or
   (at your option) any later version.

   GENFIT is distributed in the hope that it will be useful,
   but WITHOUT ANY WARRANTY; without even the implied warranty of
   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
   GNU Lesser General Public License for more details.

   You should have received a copy of the GNU Lesser General Public License
   along with GENFIT.  If not, see <http://www.gnu.org/licenses/>.
*/
#ifndef GFTRACK_H
#define GFTRACK_H

#include <vector>

#include "larreco/Genfit/GFAbsTrackRep.h"
#include "larreco/Genfit/GFException.h"

#include "TObjArray.h"

#include "larreco/Genfit/GFBookkeeping.h"
#include "larreco/Genfit/GFDetPlane.h"
#include "larreco/Genfit/GFTrackCand.h"

class TVirtualGeoTrack;

namespace genf {
  class GFAbsRecoHit;
}

namespace genf {

  class GFTrack : public TObject {
  private:
    TObjArray* fTrackReps; //->

    std::vector<GFAbsRecoHit*> fHits;

    std::vector<GFBookkeeping*> fBookkeeping;

    std::vector<int> fRepAtHit;

    GFTrackCand fCand; // list of hits
    // Cov and Updates at each step (on last iteration!) as calculated
    // in GFKalman.cxx.
    std::vector<TMatrixT<Double_t>> fHitMeasCov;
    std::vector<TMatrixT<Double_t>> fHitUpdate;
    std::vector<TMatrixT<Double_t>> fHitState;
    std::vector<Double_t> fHitChi2;
    std::vector<TMatrixT<Double_t>> fHitCov7x7;
    std::vector<TMatrixT<Double_t>> fHitCov;
    std::vector<TVector3> fHitPlaneXYZ;
    std::vector<TVector3> fHitPlaneUxUyUz;
    std::vector<TVector3> fHitPlaneU;
    std::vector<TVector3> fHitPlaneV;

    int fPDG;

    static const int fDefNumTrackReps = 10; 
    unsigned int fCardinal_rep;             // THE selected rep, default=0;

    unsigned int fNextHitToFit;

  public:
    GFTrack();

    GFTrack(const GFTrack&);

    GFTrack& operator=(const GFTrack&);

    GFTrack(GFAbsTrackRep*);

    virtual ~GFTrack();

    // -----------------------
    // Accessors
    // -----------------------

    void reset(); // deletes the RecoHits!

    int getFailedHits(int repId = -1)
    {
      int theRep;
      if (repId == -1)
        theRep = fCardinal_rep;
      else
        theRep = repId;
      return fBookkeeping.at(theRep)->getNumFailed();
    }

    std::vector<GFAbsRecoHit*> getHits() { return fHits; }

    const GFTrackCand& getCand() const { return fCand; }

    GFAbsRecoHit* getHit(int id) const { return fHits.at(id); }

    unsigned int getNumHits() const { return fHits.size(); }

    void mergeHits(GFTrack* trk);

    void releaseHits() { fHits.clear(); }

    unsigned int getNextHitToFit() const { return fNextHitToFit; }

    void setNextHitToFit(unsigned int i) { fNextHitToFit = i; }

    GFAbsTrackRep* getTrackRep(int id) const
    {
      return reinterpret_cast<GFAbsTrackRep*>(fTrackReps->At(id));
    }

    unsigned int getNumReps() const { return fTrackReps->GetEntriesFast(); }

    GFAbsTrackRep* getCardinalRep() const
    {
      return ((GFAbsTrackRep*)fTrackReps->At(fCardinal_rep));
    }

    TVector3 getMom() const { return getCardinalRep()->getMom(); }

    TVector3 getMom(const GFDetPlane& pl) const { return getCardinalRep()->getMom(pl); }

    TVector3 getPos() const { return getCardinalRep()->getPos(); }

    TVector3 getPos(const GFDetPlane& pl) const { return getCardinalRep()->getPos(pl); }

    void getPosMomCov(TVector3& pos, TVector3& mom, TMatrixT<Double_t>& cov)
    {
      getCardinalRep()->getPosMomCov(pos, mom, cov);
    }

    void getPosMomCov(const GFDetPlane& pl, TVector3& pos, TVector3& mom, TMatrixT<Double_t>& cov)
    {
      getCardinalRep()->getPosMomCov(pl, pos, mom, cov);
    }

    double getChiSqu() const { return getCardinalRep()->getChiSqu(); }

    unsigned int getNDF() const { return getCardinalRep()->getNDF(); }

    double getRedChiSqu() const { return getCardinalRep()->getRedChiSqu(); }

    double getCharge() const { return getCardinalRep()->getCharge(); }

    void fillGeoTrack(TVirtualGeoTrack* tr) const { fillGeoTrack(tr, fCardinal_rep); }

    void fillGeoTrack(TVirtualGeoTrack* tr, unsigned int repid) const;

    // ---------------------
    // Modifiers
    // ---------------------

    void addFailedHit(unsigned int irep, unsigned int id)
    {
      fBookkeeping.at(irep)->addFailedHit(id);
    }

    inline void addHit(GFAbsRecoHit* theHit) { fHits.push_back(theHit); }

    void addHit(GFAbsRecoHit* theHit,
                unsigned int detId,
                unsigned int hitId,
                double rho = 0.,
                unsigned int planeId = 0)
    {
      fHits.push_back(theHit);
      fCand.addHit(detId, hitId, rho, planeId);
    }

    void addHitVector(std::vector<GFAbsRecoHit*> hits) { fHits = hits; }

    void addTrackRep(GFAbsTrackRep* theTrackRep)
    {
      if (fTrackReps == NULL) fTrackReps = new TObjArray(fDefNumTrackReps);
      fTrackReps->Add(theTrackRep);
      fBookkeeping.push_back(new GFBookkeeping());
      fRepAtHit.push_back(-1);
    }

    GFBookkeeping* getBK(int index = -1)
    {
      if (index == -1) return fBookkeeping.at(fCardinal_rep);
      if ((unsigned int)index >= getNumReps())
        throw GFException("genf::GFTrack::getBK(): index out of range", __LINE__, __FILE__)
          .setFatal();
      return fBookkeeping.at(index);
    }

    void setCandidate(const GFTrackCand& cand, bool doreset = false);

    void setCardinalRep(unsigned int r)
    {
      if ((int)r < fTrackReps->GetEntriesFast()) fCardinal_rep = r;
    }

    void setHitMeasuredCov(TMatrixT<Double_t> mat) { fHitMeasCov.push_back(mat); }
    void setHitUpdate(TMatrixT<Double_t> mat) { fHitUpdate.push_back(mat); }
    void setHitChi2(Double_t mat) { fHitChi2.push_back(mat); }
    void setHitState(TMatrixT<Double_t> mat) { fHitState.push_back(mat); }
    void setHitCov(TMatrixT<Double_t> mat) { fHitCov.push_back(mat); }
    void setHitCov7x7(TMatrixT<Double_t> mat) { fHitCov7x7.push_back(mat); }
    void setHitPlaneXYZ(TVector3 pl) { fHitPlaneXYZ.push_back(pl); }
    void setHitPlaneUxUyUz(TVector3 pl) { fHitPlaneUxUyUz.push_back(pl); }
    void setHitPlaneU(TVector3 pl) { fHitPlaneU.push_back(pl); }
    void setHitPlaneV(TVector3 pl) { fHitPlaneV.push_back(pl); }
    void setPDG(int pdgt) { fPDG = pdgt; }
    std::vector<TMatrixT<Double_t>> getHitMeasuredCov() { return fHitMeasCov; }
    std::vector<TMatrixT<Double_t>> getHitUpdate() { return fHitUpdate; }
    std::vector<Double_t> getHitChi2() { return fHitChi2; }
    std::vector<TMatrixT<Double_t>> getHitState() { return fHitState; }
    std::vector<TMatrixT<Double_t>> getHitCov() { return fHitCov; }
    std::vector<TMatrixT<Double_t>> getHitCov7x7() { return fHitCov; }
    std::vector<TVector3> getHitPlaneXYZ() { return fHitPlaneXYZ; }
    std::vector<TVector3> getHitPlaneUxUyUz() { return fHitPlaneUxUyUz; }
    std::vector<TVector3> getHitPlaneU() { return fHitPlaneU; }
    std::vector<TVector3> getHitPlaneV() { return fHitPlaneV; }

    int getPDG() { return fPDG; }
    void getResiduals(unsigned int detId, // which detector?
                      unsigned int dim,   // which projection?
                      unsigned int rep,   // which trackrep ?
                      std::vector<double>& result);

    void setRepAtHit(unsigned int irep, int ihit)
    {
      if (irep >= getNumReps())
        throw GFException("genf::GFTrack::setRepAtHit(): index out of range", __LINE__, __FILE__)
          .setFatal();
      fRepAtHit.at(irep) = ihit;
    }

    int getRepAtHit(unsigned int irep)
    {
      if (irep >= getNumReps())
        throw GFException("genf::GFTrack::getRepAtHit(): index out of range", __LINE__, __FILE__)
          .setFatal();
      return fRepAtHit.at(irep);
    }

    void clearRepAtHit()
    {
      for (unsigned int i = 0; i < getNumReps(); ++i) {
        fRepAtHit.at(i) = -1;
      }
    }

    void printBookkeeping(std::ostream& out = std::cout) const;

    void Print(std::ostream& out = std::cout) const;

    void clearBookkeeping()
    {
      for (unsigned int i = 0; i < getNumReps(); ++i) {
        fBookkeeping.at(i)->clearAll();
      }
    }

    void clearFailedHits()
    {
      for (unsigned int i = 0; i < getNumReps(); ++i) {
        fBookkeeping.at(i)->clearFailedHits();
      }
    }

    void getHitsByPlane(std::vector<std::vector<int>*>& retVal);

  private:
    virtual void Print(Option_t*) const
    {
      throw std::logic_error(std::string(__func__) + "::Print(Option_t*) not available");
    }

    //public:
    //ClassDef(GFTrack,1)
  };
} // namespace genf

#endif