TUtil.h

Go to the documentation of this file.

#include "TLorentzVector.h"
#include "TConstant.h"
#include <complex>
 
typedef std::complex<double> complex_t;
 
inline complex_t operator+(const double &x, const complex_t &y) {
  return y + x;
}
 
inline complex_t operator-(const double &x, const complex_t &y) {
  return -(y - x);
}
 
namespace rb{
  inline double Gamma(const double Q2, 
              const double M2, const double G,
              const double mkpp2, const double mkmp2){
    const double mkmp22 = mkmp2*mkmp2;
    double beta2_q = (Q2*Q2 - 2*Q2*mkpp2 + mkmp22)/pow(Q2 - mkpp2,2);
    double beta2_M = (M2*M2 - 2*M2*mkpp2 + mkmp22)/pow(M2 - mkpp2,2);
    double beta2 = beta2_q/beta2_M;
    double beta = sqrt(beta2);
    return G*M2/Q2*beta2*beta;
  }
  
  inline complex_t R(double Q2, 
             double M2, double G, 
             double Mp2, double Mm2){
    if(Q2<Mp2) return 0;
    const complex_t I(0,1);
    double Q = sqrt(Q2);
    complex_t prop = M2/(Q2 - M2 + I*(Gamma(Q2,M2,G,Mp2,Mm2)*Q));
    return prop;
  }  
 
  inline TLorentzVector PseudoScalars3(const TLorentzVector &qp,
                       const TLorentzVector &qm,
                       const TLorentzVector &q0){
    // The hadronic current
    // here for 3pi J_\mu=\epsilon_{a b g \mu}q+_{a}q-_{b}q0_{g}
 
    double Jx = qp.Y()*qm.Z()*q0.T()
      -    qp.Y()*qm.T()*q0.Z()
      +    qp.T()*qm.Y()*q0.Z()
      -    qp.T()*qm.Z()*q0.Y()
      +    qp.Z()*qm.T()*q0.Y()
      -    qp.Z()*qm.Y()*q0.T();
    
    double Jy = qp.X()*qm.Z()*q0.T()
      -    qp.X()*qm.T()*q0.Z()
      +    qp.T()*qm.X()*q0.Z()
      -    qp.T()*qm.Z()*q0.X()
      +    qp.Z()*qm.T()*q0.X()
      -    qp.Z()*qm.X()*q0.T();
    
    double Jz = qp.Y()*qm.X()*q0.T()
      -    qp.Y()*qm.T()*q0.X()
      +    qp.T()*qm.Y()*q0.X()
      -    qp.T()*qm.X()*q0.Y()
      +    qp.X()*qm.T()*q0.Y()
      -    qp.X()*qm.Y()*q0.T();
    
    double Jt = qp.Y()*qm.Z()*q0.X()
      -    qp.Y()*qm.X()*q0.Z()
      +    qp.X()*qm.Y()*q0.Z()
      -    qp.X()*qm.Z()*q0.Y()
      +    qp.Z()*qm.X()*q0.Y()
      -    qp.Z()*qm.Y()*q0.X();
    
    TLorentzVector J(Jx,Jy,Jz,Jt);
    return J;
  }
 
}