examples/calDepth.c

   1 #include <stdlib.h>
   2 #include <string.h>
   3 #include <stdio.h>
   4 #include <unistd.h>
   5 #include "sam.h"
   6
   7 typedef struct {
   8         int beg, end;
   9         samfile_t *in;
  10 } tmpstruct_t;
  11
  12 // callback for bam_fetch()
  13 static int fetch_func(const bam1_t *b, void *data)
  14 {
  15         bam_plbuf_t *buf = (bam_plbuf_t*)data;
  16         bam_plbuf_push(b, buf);
  17         return 0;
  18 }
  19 // callback for bam_plbuf_init()
  20 static int pileup_func(uint32_t tid, uint32_t pos, int n, const bam_pileup1_t *pl, void *data)
  21 {
  22         tmpstruct_t *tmp = (tmpstruct_t*)data;
  23         if ((int)pos >= tmp->beg && (int)pos < tmp->end)
  24                 printf("%s\t%d\t%d\n", tmp->in->header->target_name[tid], pos + 1, n);
  25         return 0;
  26 }
  27
  28 int main(int argc, char *argv[])
  29 {
  30         tmpstruct_t tmp;
  31         if (argc == 1) {
  32                 fprintf(stderr, "Usage: calDepth <in.bam> [region]\n");
  33                 return 1;
  34         }
  35         tmp.beg = 0; tmp.end = 0x7fffffff;
  36         tmp.in = samopen(argv[1], "rb", 0);
  37         if (tmp.in == 0) {
  38                 fprintf(stderr, "Fail to open BAM file %s\n", argv[1]);
  39                 return 1;
  40         }
  41         if (argc == 2) {
  42                 sampileup(tmp.in, -1, pileup_func, &tmp);
  43         } else {
  44                 int ref;
  45                 bam_index_t *idx;
  46                 bam_plbuf_t *buf;
  47                 idx = bam_index_load(argv[1]);
  48                 if (idx == 0) {
  49                         fprintf(stderr, "BAM indexing file is not available.\n");
  50                         return 1;
  51                 }
  52                 bam_parse_region(tmp.in->header, argv[2], &ref, &tmp.beg, &tmp.end);
  53                 if (ref < 0) {
  54                         fprintf(stderr, "Invalid region %s\n", argv[2]);
  55                         return 1;
  56                 }
  57                 buf = bam_plbuf_init(pileup_func, &tmp);
  58                 bam_fetch(tmp.in->x.bam, idx, ref, tmp.beg, tmp.end, buf, fetch_func);
  59                 bam_plbuf_push(0, buf);
  60                 bam_index_destroy(idx);
  61                 bam_plbuf_destroy(buf);
  62         }
  63         samclose(tmp.in);
  64         return 0;
  65 }