GsfBetheHeitlerUpdator.h

Go to the documentation of this file.

#ifndef GsfBetheHeitlerUpdator_h_
#define GsfBetheHeitlerUpdator_h_

#include "TrackingTools/GsfTracking/interface/GsfMaterialEffectsUpdator.h"
#include "TrackingTools/GsfTracking/interface/Triplet.h"

#include <vector>
#include <iosfwd>
#include <string>

class GsfBetheHeitlerUpdator : public GsfMaterialEffectsUpdator {

private:
  class Polynomial {
  public:
    Polynomial () {}
    Polynomial (const std::vector<double>& coefficients) :
      theCoeffs(coefficients) {}
    double operator() (const double& x) const {
      double sum(0.);
      for ( std::vector<double>::const_iterator i=theCoeffs.begin();
        i!=theCoeffs.end(); i++ )  sum = x*sum + *i;
      return sum;
    }
  private:
    std::vector<double> theCoeffs;
  };

public:
  enum CorrectionFlag { NoCorrection=0, MeanCorrection=1, FullCorrection=2 };

public:
  virtual GsfBetheHeitlerUpdator* clone() const
  {
    return new GsfBetheHeitlerUpdator(*this);
  }
  
public:
  GsfBetheHeitlerUpdator (const std::string fileName, const int correctionFlag);
//   GsfBetheHeitlerUpdator (const std::string fileName, const CorrectionFlag correctionFlag);

private:
  typedef std::vector< Triplet<double,double,double> > GSContainer;

private:
  virtual void compute (const TrajectoryStateOnSurface&, const PropagationDirection) const;
  void readParameters (const std::string);
  Polynomial readPolynomial (std::ifstream&,const int);
  inline double logisticFunction (const double x) const {return 1./(1.+exp(-x));}
  inline double BetheHeitlerMean (const double rl) const
  {
    return exp(-rl);
  }
  inline double BetheHeitlerVariance (const double rl) const
  {
    return exp(-rl*log(3.)/log(2.)) - exp(-2*rl);
  }
  void getMixtureParameters (const double, GSContainer&) const;
  void correctWeights (GSContainer&) const;
  double correctedFirstMean (const double, const GSContainer&) const;
  double correctedFirstVar (const double,const GSContainer&) const;
  
protected:
  virtual bool newArguments (const TrajectoryStateOnSurface&, const PropagationDirection) const;
  virtual void storeArguments (const TrajectoryStateOnSurface&, const PropagationDirection) const;

private:
  int theNrComponents;                  
  int theTransformationCode;            
  std::vector<Polynomial> thePolyWeights;    
  std::vector<Polynomial> thePolyMeans;      
  std::vector<Polynomial> thePolyVars;       

  int theCorrectionFlag;                

  mutable float theLastDz;
  mutable float theLastP;
  mutable PropagationDirection theLastPropDir;
  mutable float theLastRadLength;
};

#endif