ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/gclib/tophat_cpp/bam_dump.cpp
Revision: 29
Committed: Tue Aug 2 21:24:54 2011 UTC (8 years, 10 months ago) by gpertea
File size: 1386 byte(s)
Log Message:
adding tophat source work

Line User Rev File contents
1 gpertea 29 #include <cstdlib>
2     #include <cstdio>
3    
4     #include "bam/bam.h"
5     #include "bam/sam.h"
6    
7     int main(int argc, char *argv[]) {
8    
9     samfile_t *fp_in = NULL;
10     bam1_t *b=NULL;
11    
12     fp_in = samopen(argv[1], "rb", 0);
13    
14     if(NULL == fp_in) {
15     printf("Could not open file\n");
16     }
17    
18     b = bam_init1();
19     int pos=0;
20     int lastpos=0;
21    
22     while(samread(fp_in, b) > 0) {
23     lastpos = pos;
24     pos = b->core.pos;
25    
26     if(pos != lastpos) {
27     printf("tid : %d\n",b->core.tid);
28     printf("pos : %d\n",b->core.pos);
29     char *name = bam1_qname(b);
30     char *qual = (char*)bam1_qual(b);
31    
32     int n=0;
33     char *qseq = (char *) malloc(b->core.l_qseq+1);
34     char *s = (char*)bam1_seq(b);
35     for(n=0;n<(b->core.l_qseq);n++) {
36     char v = bam1_seqi(s,n);
37     qseq[n] = bam_nt16_rev_table[v];
38     }
39     qseq[n] = 0;
40    
41     printf("name : %s\n",name);
42     //print SAM text format of this BAM record:
43     char * samtext=bam_format1(fp_in->header, b);
44     printf("%s\n",samtext);
45    
46     printf("qseq : %s\n",qseq);
47     //printf("s cigar: %s\n",cigar);
48    
49     printf("qual : ");
50     for(n=0;n<(b->core.l_qseq);n++) {
51     //printf(" %d",qual[n]);
52     qseq[n]=qual[n]+33;
53     }
54     qseq[n]=0;
55     printf("%s\n",qseq);
56     }
57     bam_destroy1(b);
58     b = bam_init1();
59     }
60     bam_destroy1(b);
61    
62     samclose(fp_in);
63     return 0;
64     }