Browse code

add support for the noncovered argument; tells bam_tally to report even those positions that have no coverage

git-svn-id: file:///home/git/hedgehog.fhcrc.org/bioconductor/trunk/madman/Rpacks/gmapR@89724 bc3139a8-67e5-0310-9ffc-ced21a209358

Michael Lawrence authored on 30/04/2014 17:05:27
Showing 5 changed files

... ...
@@ -18,7 +18,8 @@ setClass("BamTallyParam",
18 18
                         variant_strand = "integer",
19 19
                         ignore_query_Ns = "logical",
20 20
                         indels = "logical",
21
-                        include_soft_clips = "integer"))
21
+                        include_soft_clips = "integer",
22
+                        noncovered = "logical"))
22 23
 
23 24
 ### - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
24 25
 ### Constructor
... ...
@@ -39,7 +40,8 @@ BamTallyParam <- function(genome, which = GRanges(),
39 40
                           primary_only = FALSE, ignore_duplicates = FALSE,
40 41
                           min_depth = 0L, variant_strand = 0L,
41 42
                           ignore_query_Ns = FALSE,
42
-                          indels = FALSE, include_soft_clips = 0L)
43
+                          indels = FALSE, include_soft_clips = 0L,
44
+                          noncovered = FALSE)
43 45
 {
44 46
   if (!is.null(desired_read_group) && !isSingleString(desired_read_group))
45 47
     stop("'desired_read_group' must be NULL or a single, non-NA string")
... ...
@@ -70,7 +72,7 @@ BamTallyParam <- function(genome, which = GRanges(),
70 72
   integer_params <- c("minimum_mapq", "min_depth", "variant_strand",
71 73
                       "include_soft_clips")
72 74
   params[integer_params] <- lapply(params[integer_params], as.integer)
73
-  do.call(new, c("BamTallyParam", params))  
75
+  do.call(new, c("BamTallyParam", params))
74 76
 }
75 77
 
76 78
 ### - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
... ...
@@ -195,7 +195,7 @@ normArgTRUEorFALSE <- function(x) {
195 195
                          ignore_query_Ns = FALSE,
196 196
                          indels = FALSE,
197 197
                          blocksize = 1000L, verbosep = FALSE,
198
-                         include_soft_clips = 0L)
198
+                         include_soft_clips = 0L, noncovered = FALSE)
199 199
 {
200 200
   if (!is(bamreader, "GmapBamReader"))
201 201
     stop("'bamreader' must be a GmapBamReader")
... ...
@@ -217,7 +217,7 @@ normArgTRUEorFALSE <- function(x) {
217 217
     if (is.unsorted(read_pos_breaks))
218 218
       stop("'read_pos_breaks' must be sorted")
219 219
   }
220
-  .Call(R_Bamtally_iit, bamreader@.extptr, genome_dir, db, which,
220
+  ptr <- .Call(R_Bamtally_iit, bamreader@.extptr, genome_dir, db, which,
221 221
         desired_read_group,
222 222
         normArgSingleInteger(alloclength),
223 223
         normArgSingleInteger(minimum_mapq),
... ...
@@ -233,7 +233,8 @@ normArgTRUEorFALSE <- function(x) {
233 233
         normArgTRUEorFALSE(indels),
234 234
         normArgSingleInteger(blocksize),
235 235
         normArgTRUEorFALSE(verbosep),
236
-        normArgSingleInteger(include_soft_clips))
236
+        normArgSingleInteger(include_soft_clips),
237
+        normArgTRUEorFALSE(noncovered))
237 238
 }
238 239
 
239 240
 ### - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
... ...
@@ -7,7 +7,7 @@
7 7
 static const R_CallMethodDef callMethods[] = {
8 8
 
9 9
   /* bamtally.c */
10
-  CALLMETHOD_DEF(R_Bamtally_iit, 20),
10
+  CALLMETHOD_DEF(R_Bamtally_iit, 21),
11 11
   CALLMETHOD_DEF(R_tally_iit_parse, 5),
12 12
   
13 13
   /* bamreader.c */
... ...
@@ -23,7 +23,7 @@ R_Bamtally_iit (SEXP bamreader_R, SEXP genome_dir_R, SEXP db_R,
23 23
                 SEXP ignore_query_Ns_p_R,
24 24
                 SEXP print_indels_p_R,
25 25
                 SEXP blocksize_R, 
26
-                SEXP verbosep_R, SEXP max_softclip_R)
26
+                SEXP verbosep_R, SEXP max_softclip_R, SEXP noncovered_R)
27 27
 {
28 28
   Bamreader_T bamreader = (Bamreader_T) R_ExternalPtrAddr(bamreader_R);
29 29
   const char *genome_dir =
... ...
@@ -47,6 +47,7 @@ R_Bamtally_iit (SEXP bamreader_R, SEXP genome_dir_R, SEXP db_R,
47 47
   int blocksize = asInteger(blocksize_R);
48 48
   int verbosep = asLogical(verbosep_R);
49 49
   int max_softclip = asInteger(max_softclip_R);
50
+  int noncovered = asLogical(noncovered_R);
50 51
 
51 52
   Genome_T genome = createGenome(genome_dir, db);
52 53
   IIT_T chromosome_iit = readChromosomeIIT(genome_dir, db);
... ...
@@ -72,7 +73,7 @@ R_Bamtally_iit (SEXP bamreader_R, SEXP genome_dir_R, SEXP db_R,
72 73
                                  ignore_duplicates_p,
73 74
                                  min_depth, variant_strands, ignore_query_Ns_p,
74 75
                                  print_indels_p, blocksize, verbosep,
75
-                                 /*readlevel_p*/false, max_softclip);
76
+                                 /*readlevel_p*/false, max_softclip, noncovered);
76 77
   IIT_free(&chromosome_iit);
77 78
   Genome_free(&genome);
78 79
 
... ...
@@ -20,7 +20,7 @@ R_Bamtally_iit (SEXP bamreader_R, SEXP genome_dir_R, SEXP db_R,
20 20
                 SEXP ignore_query_Ns_p,
21 21
                 SEXP print_indels_p_R,
22 22
                 SEXP blocksize_R, 
23
-                SEXP verbosep_R, SEXP max_softclip_R);
23
+                SEXP verbosep_R, SEXP max_softclip_R, SEXP noncovered_R);
24 24
 
25 25
 SEXP
26 26
 R_tally_iit_parse(SEXP tally_iit_R, SEXP cycle_breaks_R,