GFAbsRecoHit.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 GFABSRECOHIT_H
#define GFABSRECOHIT_H

#include <iostream>
#include <stdexcept> // std::logic_error

#include "TMatrixT.h"
#include "TObject.h"

#include "larreco/Genfit/GFDetPlane.h"
#include "larreco/Genfit/GFException.h" // PrintROOTobject()
#include <cmath>

namespace genf {
  class GFAbsTrackRep;
}

namespace genf {

  class GFAbsRecoHit : public TObject {
  protected:
    TMatrixT<Double_t> fHitCoord;

    TMatrixT<Double_t> fHitCov;

  private:
    int fNparHit;

  public:
    virtual ~GFAbsRecoHit();

    GFAbsRecoHit(int NparHit);

    GFAbsRecoHit();

    virtual TMatrixT<Double_t> getHMatrix(const GFAbsTrackRep* stateVector) = 0;
    virtual TMatrixT<Double_t> getHMatrix(const GFAbsTrackRep* stateVector,
                                          const Double_t&,
                                          const Double_t&) = 0;

    virtual TMatrixT<Double_t> residualVector(const GFAbsTrackRep* stateVector,
                                              const TMatrixT<Double_t>& state,
                                              const GFDetPlane& d)
    {
      std::cout << "GFAbsRecoHit::residualVector(3args) Not correctly Using theta -- multiple "
                   "scattering -- information !!! Fix this if you really want to use getChi2Hit"
                << std::endl;
      TMatrixT<Double_t> H = getHMatrix(stateVector);
      return (getHitCoord(d) - (H * state));
    }

    virtual TMatrixT<Double_t> residualVector(const GFAbsTrackRep* stateVector,
                                              const TMatrixT<Double_t>& state,
                                              const GFDetPlane& d,
                                              const GFDetPlane& dPrev,
                                              const double& mass)
    {
      Double_t dist = (d.getO() - dPrev.getO()).Mag();
      //    Double_t mass = 0.104; // close enough for muons, pions, I think.
      Double_t mom = fabs(1.0 / state[0][0]);
      Double_t beta = mom / sqrt(mass * mass + mom * mom);
      if (std::isnan(dist) || dist < 0.2) dist = 0.2; // don't allow 0s here.
      if (std::isnan(beta) || beta < 0.04) beta = 0.04;

      TMatrixT<Double_t> H = getHMatrix(stateVector, beta, dist);
      return (getHitCoord(d, dPrev) - (H * state));
    }

    TMatrixT<Double_t> getRawHitCov() const { return fHitCov; }

    TMatrixT<Double_t> getRawHitCoord() const { return fHitCoord; }

    virtual TMatrixT<Double_t> getHitCov(const GFDetPlane&) = 0;
    virtual TMatrixT<Double_t> getHitCov(const GFDetPlane&,
                                         const GFDetPlane&,
                                         const TMatrixT<Double_t>&,
                                         const Double_t&) = 0;

    virtual TMatrixT<Double_t> getHitCoord(const GFDetPlane&, const GFDetPlane&) = 0;
    virtual TMatrixT<Double_t> getHitCoord(const GFDetPlane&) = 0;

    virtual const GFDetPlane& getDetPlane(GFAbsTrackRep*) = 0;

    virtual GFAbsRecoHit* clone() = 0;

    virtual void Print(std::ostream& out = std::cout) const { PrintROOTobject(out, fHitCoord); }

    virtual const std::string& getPolicyName();

    int getNparHit() { return fNparHit; }

    // public:
    //ClassDef(GFAbsRecoHit,3)

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

#endif //FITTER_ABSHIT_H