#include "sigma_hz.h"
#include <math.h>

// Cross section for e+e- -> HZ from Barger et al., Phys Rev D
// 49 (1994) 79, equation (2).

float sigma_hz(float Ecm, float mh, float mz){

  const float pi = 3.14159265;
  const float hbarc = 0.1973;           // GeV-fm
  const float sw2 = 0.23;
  const float G_mu =  1.16639e-5;

  float sigma;
  float ae = 1.;
  float ve = 1. - 4.*sw2;
  float s = pow(Ecm,2);
  float mz2 = pow(mz,2);
  float mh2 = pow(mh,2);
  float beta2 = ( 1 - pow((mh+mz),2)/s ) *
                ( 1 - pow((mh-mz),2)/s );
  float beta;
  if ( beta2 <= 0. ) {
    sigma = 0.;
  }
  else {
    sigma = ( ( pow(G_mu,2)*pow(mz,4) ) / ( 96.* pi * s ) ) *  
      ( pow(ae,2) + pow(ve,2) ) * sqrt(beta2) * 
      ( beta2 + 12.*mz2/s ) / pow((1.-mz2/s), 2);;
  }
  sigma = sigma * pow(hbarc,2) * 1.e10;    // convert to pb
  return sigma;
}