+
+/***********
+ * mpileup *
+ ***********/
+
+static int mpileup(int n, char **fn)
+{
+ bamFile *fp;
+ int i, tid, pos, *n_plp;
+ const bam_pileup1_t **plp;
+ bam_mplp_t iter;
+ bam_header_t *h = 0;
+ fp = calloc(n, sizeof(void*));
+ plp = calloc(n, sizeof(void*));
+ n_plp = calloc(n, sizeof(int*));
+ for (i = 0; i < n; ++i) {
+ bam_header_t *h_tmp;
+ fp[i] = bam_open(fn[i], "r");
+ h_tmp = bam_header_read(fp[i]);
+ if (i == 0) h = h_tmp;
+ else {
+ // FIXME: to check consistency
+ bam_header_destroy(h_tmp);
+ }
+ }
+ iter = bam_mplp_init(n, bam_read1, fp);
+ while (bam_mplp_auto(iter, &tid, &pos, n_plp, plp) > 0) {
+ printf("%s\t%d\tN", h->target_name[tid], pos + 1);
+ for (i = 0; i < n; ++i) {
+ int j;
+ printf("\t%d\t", n_plp[i]);
+ if (n_plp[i] == 0) printf("*\t*");
+ else {
+ for (j = 0; j < n_plp[i]; ++j)
+ pileup_seq(plp[i] + j, pos, 0, 0);
+ putchar('\t');
+ for (j = 0; j < n_plp[i]; ++j) {
+ const bam_pileup1_t *p = plp[i] + j;
+ int c = bam1_qual(p->b)[p->qpos] + 33;
+ if (c > 126) c = 126;
+ putchar(c);
+ }
+ }
+ }
+ putchar('\n');
+ }
+ bam_mplp_destroy(iter);
+ bam_header_destroy(h);
+ for (i = 0; i < n; ++i) bam_close(fp[i]);
+ free(fp); free(plp);
+ return 0;
+}
+
+int bam_mpileup(int argc, char *argv[])
+{
+ if (argc == 1) {
+ fprintf(stderr, "Usage: samtools mpileup <in1.bam> [<in2.bam> [...]]\n");
+ return 1;
+ }
+ mpileup(argc - 1, argv + 1);
+ return 0;
+}