X-Git-Url: https://git.donarmstrong.com/?a=blobdiff_plain;f=bam_plcmd.c;h=002297a1229678a256e9045f420b89d5f273b1ea;hb=82797ab5f5ba03634967e19822e11a491c8e6198;hp=7f13eee4ed51a987a6386796eb4eeae2bad65891;hpb=4c8c9dfc1e3b3b066a62a703fd3ba04db6ad5a45;p=samtools.git diff --git a/bam_plcmd.c b/bam_plcmd.c index 7f13eee..002297a 100644 --- a/bam_plcmd.c +++ b/bam_plcmd.c @@ -563,9 +563,10 @@ static int mplp_func(void *data, bam1_t *b) mplp_aux_t *ma = (mplp_aux_t*)data; int ret, skip = 0; do { - int has_ref = (ma->ref && ma->ref_id == b->core.tid)? 1 : 0; + int has_ref; ret = ma->iter? bam_iter_read(ma->fp, ma->iter, b) : bam_read1(ma->fp, b); if (ret < 0) break; + has_ref = (ma->ref && ma->ref_id == b->core.tid)? 1 : 0; skip = 0; if (has_ref && (ma->flag&MPLP_REALN)) bam_prob_realn_core(b, ma->ref, 1); if (has_ref && ma->capQ_thres > 10) { @@ -738,31 +739,16 @@ static int mpileup(mplp_conf_t *conf, int n, char **fn) (conf->flag&MPLP_FMT_SP), 0, 0); bcf_write(bp, bh, b); bcf_destroy(b); - if (!(conf->flag&MPLP_NO_INDEL)) { - // call MNPs - if (bcf_call_mnp_prep(gplp.n, gplp.n_plp, gplp.plp, pos, bca, ref) >= 0) { - for (i = 0; i < gplp.n; ++i) - bcf_call_glfgen(gplp.n_plp[i], gplp.plp[i], B2B_REF_MNP, bca, bcr + i); - if (bcf_call_combine(gplp.n, bcr, B2B_REF_MNP, &bc) >= 0) { - b = calloc(1, sizeof(bcf1_t)); - bcf_call2bcf(tid, pos, &bc, b, (conf->flag&(MPLP_FMT_DP|MPLP_FMT_SP))? bcr : 0, - (conf->flag&MPLP_FMT_SP), bca, ref); - bcf_write(bp, bh, b); - bcf_destroy(b); - bca->last_mnp_pos = pos; - } - } - // call indels - if (bcf_call_gap_prep(gplp.n, gplp.n_plp, gplp.plp, pos, bca, ref, rghash) >= 0) { - for (i = 0; i < gplp.n; ++i) - bcf_call_glfgen(gplp.n_plp[i], gplp.plp[i], B2B_REF_INDEL, bca, bcr + i); - if (bcf_call_combine(gplp.n, bcr, B2B_REF_INDEL, &bc) >= 0) { - b = calloc(1, sizeof(bcf1_t)); - bcf_call2bcf(tid, pos, &bc, b, (conf->flag&(MPLP_FMT_DP|MPLP_FMT_SP))? bcr : 0, - (conf->flag&MPLP_FMT_SP), bca, ref); - bcf_write(bp, bh, b); - bcf_destroy(b); - } + // call indels + if (!(conf->flag&MPLP_NO_INDEL) && bcf_call_gap_prep(gplp.n, gplp.n_plp, gplp.plp, pos, bca, ref, rghash) >= 0) { + for (i = 0; i < gplp.n; ++i) + bcf_call_glfgen(gplp.n_plp[i], gplp.plp[i], -1, bca, bcr + i); + if (bcf_call_combine(gplp.n, bcr, -1, &bc) >= 0) { + b = calloc(1, sizeof(bcf1_t)); + bcf_call2bcf(tid, pos, &bc, b, (conf->flag&(MPLP_FMT_DP|MPLP_FMT_SP))? bcr : 0, + (conf->flag&MPLP_FMT_SP), bca, ref); + bcf_write(bp, bh, b); + bcf_destroy(b); } } } else {