MultiscatteringModel.h

Go to the documentation of this file.

#ifndef MULTISCATTERINGMODEL_H
#define MULTISCATTERINGMODEL_H
 
//#include <cmath>
#include "TMath.h"
 
//#include <omp.h>
 
 
class MultiscatteringModel {
  std::vector<double> xlag32, wlag32;
 
public:
  MultiscatteringModel();
 
  virtual ~MultiscatteringModel(void) {}
 
  double GetSigt(double apt) { return TMath::Sqrt(2. / TMath::Pi()) * apt; }
 
  double GetSigz(double apz2, double apz4) {
    double ret = 3. / 2. * apz2 * apz2 - apz4 / 2.;
    if (ret < 0.)
      return TMath::Sqrt(apz2);
    else
      return TMath::Sqrt(apz2 - TMath::Sqrt(ret));
  }
 
  double GetQz(double apz2, double apz4) {
    double ret = 3. / 2. * apz2 * apz2 - apz4 / 2.;
    if (ret < 0.) return 0.;
    return TMath::Power(ret, 1. / 4.);
  }
 
  double fpt(double pt, double sigt) {
    return pt / sigt / sigt * TMath::Exp(-pt * pt / 2. / sigt / sigt);
  }
 
  double fpz(double pz, double sigz, double qz) {
    return 1. / 2. / TMath::Sqrt(2. * TMath::Pi()) / sigz
           * (TMath::Exp(-(pz - qz) * (pz - qz) / 2. / sigz / sigz)
              + TMath::Exp(-(pz + qz) * (pz + qz) / 2. / sigz / sigz));
  }
 
  double dndy(double y, double m, double sigt, double sigz, double qz) {
    double ret = 0.;
    for (int i = 0; i < 32; i++) {
      double mt = TMath::Sqrt(m * m + xlag32[i] * xlag32[i]);
      ret += wlag32[i] * mt * fpt(xlag32[i], sigt)
             * fpz(mt * TMath::SinH(y), sigz, qz);
    }
    return ret * TMath::CosH(y);
  }
 
  double
  fypt(double y, double pt, double m, double sigt, double sigz, double qz) {
    double mt = TMath::Sqrt(m * m + pt * pt);
    return mt * TMath::CosH(y) * fpt(pt, sigt)
           * fpz(mt * TMath::SinH(y), sigz, qz);
  }
};
 
#endif