+ if (ma->n1 > 0 && ma->n1 < ma->n) {
+ for (k = 0; k < 2; ++k) {
+ flast = rst->f_em;
+ for (i = 0; i < MC_MAX_EM_ITER; ++i) {
+ rst->f_em2[k] = k? mc_freq_iter(flast, ma, ma->n1, ma->n) : mc_freq_iter(flast, ma, 0, ma->n1);
+ if (fabs(rst->f_em2[k] - flast) < MC_EM_EPS) break;
+ flast = rst->f_em2[k];
+ }
+ }
+ }
+ }
+ { // compute g[3]
+ rst->g[0] = (1. - rst->f_em) * (1. - rst->f_em);
+ rst->g[1] = 2. * rst->f_em * (1. - rst->f_em);
+ rst->g[2] = rst->f_em * rst->f_em;
+ for (i = 0; i < MC_MAX_EM_ITER; ++i)
+ if (mc_gtfreq_iter(rst->g, ma, 0, ma->n) < MC_EM_EPS) break;