add KHASH to accumulate PASS site

htsutils.c

@@ -14,6 +14,10 @@ int parseBAM(Settings *props)
     bam_hdr_t *header = sam_hdr_read(fh_in);
     bam1_t *bam = bam_init1();
 
+    int absent;
+    khiter_t k;
+    khash_t(pass_hash) *hptr = kh_init(pass_hash);
+
     uint32_t lines = 0;
     uint32_t poly = 0;
 
@@ -47,23 +51,48 @@ int parseBAM(Settings *props)
         if ((polyTailSize >= props->sizepoly) && (polyReadSize <= (bam->core.l_qseq - props->sizepoly)))
         {
             poly++;
+
+            // add polyBasePosition to hash
+            k = kh_put(pass_hash, hptr, polyBasePosition, &absent);
+            if (absent)
+                kh_value(hptr, k) = 1;
+            else
+                kh_value(hptr, k)++;
+
         }
 
+
+
         //printf("Tail: %d|%d\tHead: %d|%d\tSeq:: ", polyTailA,polyTailT,polyHeadA,polyHeadT);
         //dumpBamSeq(bam);
         //printf("\n");
 
         lines++;
-        //if (lines == 10)
+        //if (lines == 1)
         //    break;
 
     }
 
+    fprintf(stderr, "BAM Lines: %d / %d\n", poly, lines);
+
+    // loop through hash and print values
+    fprintf(stderr, "HASH Size: %d\n", kh_size(hptr));
+    for (k = kh_begin(hptr); k != kh_end(hptr); ++k)
+    {
+        if (kh_exist(hptr, k))
+        {
+            int key = kh_key(hptr, k);
+            int val = kh_val(hptr, k);
+            printf("%d\t%d\n", key, val);
+        }
+    }
+
+    kh_destroy(pass_hash, hptr);
     bam_destroy1(bam);
     bam_hdr_destroy(header);
     sam_close(fh_in);
 
-    printf("BAM Lines: %d / %d\n", poly, lines);
+
 
     return 0;
 }

htsutils.h

@@ -12,12 +12,17 @@
 #include <stdio.h>
 #include <stdlib.h>
 #include <ctype.h>
+
 #include <htslib/sam.h>
+#include <htslib/khash.h>
+
 #include "config.h"
 
 #define POLY_HEAD 0
 #define POLY_TAIL 1
 
+KHASH_MAP_INIT_INT(pass_hash, int);
+
 int parseBAM(Settings *);
 void dumpBamSeq(bam1_t *);
 

main.c

@@ -31,7 +31,7 @@ int main(int argc, const char * argv[])
 
     double time_spent = (double)(end - begin) / CLOCKS_PER_SEC;
 
-    printf("Time: %.4f sec\n", time_spent);
+    fprintf(stderr, "Time: %.4f sec\n", time_spent);
 
     return 0;
 }