-static double mc_ref_prob(const mc_aux_t *ma, double *_PD, double *f_exp)
-{
- int k, i;
- long double PD = 0., Pref = 0., Ef = 0.;
- for (k = 0; k <= ma->M; ++k) {
- long double x = 1., y = 0.;
- double *bk = ma->beta + k * 3;
- for (i = 0; i < ma->n; ++i) {
- double *pdg = ma->pdg + i * 3;
- double z = pdg[0] * bk[0] + pdg[1] * bk[1] + pdg[2] * bk[2];
- x *= z;
- y += (pdg[1] * bk[1] + 2. * pdg[2] * bk[2]) / z;
- }
- PD += x * ma->alpha[k];
- Ef += x * y * ma->alpha[k];
- }
- for (k = 0; k <= ma->n * 2; ++k) {
- long double x = 1.0;
- for (i = 0; i < ma->n; ++i)
- x *= ma->pdg[i * 3 + 2] * ma->beta[k * 3 + 2];
- Pref += x * ma->alpha[k];
- }
- *f_exp = (double)(Ef / PD / ma->M);
- *_PD = PD;
- return Pref / PD;
-}
-