Merge branch 'develop' of github.com:samtools/bcftools into develop

Author: pd3

NEWS

@@ -23,6 +23,13 @@ Changes affecting specific commands:
       when the site is monoallelic in both query and genotype file. The new option --keep-refs
       allows to always include monoallelic sites.
 
+    - Fix an error in parsing -i/-e command line options where the `qry:` and `gt:` prefix was
+      not stripped (#2432)
+
+* bcftools mpileup
+
+    - Make `-d, --max-depth 0` set the depth to unlimited (#2435)
+
 * bcftools norm
 
     - Make the -i/-e filtering option work for all options, such as line merging and

doc/bcftools.txt

@@ -2308,7 +2308,8 @@ those scenarios.
     This behavior was problematic when working with a combination of
     single- and multi-sample bams, therefore in *bcftools mpileup* the user
     is given the full control (and responsibility), and an informative message
-    is printed instead [250]
+    is printed instead. Passing zero for this option sets it to the highest possible
+    value, effectively removing the depth limit [250]
 
 *-E, --redo-BAQ*::
     Recalculate BAQ on the fly, ignore existing BQ tags

mpileup.c

@@ -930,11 +930,19 @@ static int mpileup(mplp_conf_t *conf)
     // init mpileup
     conf->iter = bam_mplp_init(conf->nfiles, mplp_func, (void**)conf->mplp_data);
     if ( conf->flag & MPLP_SMART_OVERLAPS ) bam_mplp_init_overlaps(conf->iter);
-    fprintf(stderr, "[%s] maximum number of reads per input file set to -d %d\n",  __func__, conf->max_depth);
-    if ( (double)conf->max_depth * conf->nfiles > 1<<20)
-        fprintf(stderr, "Warning: Potential memory hog, up to %.0fM reads in the pileup!\n", (double)conf->max_depth*conf->nfiles);
-    if ( (double)conf->max_depth * conf->nfiles / nsmpl < 250 )
-        fprintf(stderr, "Note: The maximum per-sample depth with -d %d is %.1fx\n", conf->max_depth,(double)conf->max_depth * conf->nfiles / nsmpl);
+    if ( !conf->max_depth )
+    {
+        conf->max_depth = INT_MAX;
+        fprintf(stderr, "[%s] Max depth set to maximum value (%d)\n", __func__, INT_MAX);
+    }
+    else
+    {
+        fprintf(stderr, "[%s] maximum number of reads per input file set to -d %d\n",  __func__, conf->max_depth);
+        if ( (double)conf->max_depth * conf->nfiles > 1<<20)
+            fprintf(stderr, "Warning: Potential memory hog, up to %.0fM reads in the pileup!\n", (double)conf->max_depth*conf->nfiles);
+        if ( (double)conf->max_depth * conf->nfiles / nsmpl < 250 )
+            fprintf(stderr, "Note: The maximum per-sample depth with -d %d is %.1fx\n", conf->max_depth,(double)conf->max_depth * conf->nfiles / nsmpl);
+    }
     bam_mplp_set_maxcnt(conf->iter, conf->max_depth);
     conf->max_indel_depth = conf->max_indel_depth * nsmpl;
     conf->bcf_rec = bcf_init1();

test/norm.vcf

@@ -35,6 +35,7 @@
 ##INFO=<ID=ISTR,Number=1,Type=String,Description="Test String in INFO">
 ##INFO=<ID=END,Number=1,Type=Integer,Description="End position">
 #CHROM	POS	ID	REF	ALT	QUAL	FILTER	INFO	FORMAT	XY00001	XY00002
+1	60	.	A	T	999	PASS	AN=4;AC=4	GT	0/0	0/1
 1	105	.	TAAACCCTAAA	TAA,TAACCCTAAA	999	PASS	INDEL;AN=4;AC=2,2;DP=19;ISTR=SomeString;XRF=1e+06,2e+06,500000;XRI=1111,2222,5555;XRS=AAA,BBB,DDD;XAF=1e+06,500000;XAI=1111,5555;XAS=AAA,DDD;XGF=1e+06,2e+06,3e+06,500000,.,9e+09;XGI=1111,2222,3333,5555,.,9999;XGS=A,B,C,E,.,F	GT:PL:DP:FRF:FRI:FRS:FAF:FAI:FAS:FGF:FGI:FGS	1/2:1,2,3,4,5,6:1:1e+06,2e+06,500000:1111,2222,5555:AAAA,BBB,CC:1e+06,500000:1111,5555:A,BB:1e+06,2e+06,3e+06,500000,.,9e+09:1111,2222,3333,5555,.,9999:A,BB,CCC,EEEE,.,FFFFF	1/2:1,2,3,4,5,6:1:1e+06,2e+06,500000:1111,2222,5555:AAAA,BBB,CC:1e+06,500000:1111,5555:A,BB:1e+06,2e+06,3e+06,500000,.,9e+09:1111,2222,3333,5555,.,9999:A,BB,CCC,EEEE,.,FFFFF
 2	1	.	GGGCGTCTCATAGCTGGAGCAATGGCGAGCGCCTGGACAAGGGAGGGGAAGGGGTTCTTATTACTGACGCGGGTAGCCCCTACTGCTGTGTGGTTCCCCTATTTTTTTTTTTTTCTTTTTGAGACGGAGTCTCGCTCTGTCACCCAGGCTGGAGTGCAGTGGCACAATCTCGGCTCACTGCAAGCTCCACCT	ACGT	999	PASS	INDEL;AN=4;AC=2;END=192	GT:DP	1/0:1	1/0:1
 2	101	.	ATTTTTTTTTTTTT	ATTTTTTTTTTTTTTT	999	PASS	INDEL;AN=4;AC=4;END=114	GT:DP	1/1:1	1/1:1

vcfgtcheck.c

@@ -1324,13 +1324,13 @@ int main_vcfgtcheck(int argc, char *argv[])
                 if ( !strncasecmp("gt:",optarg,3) )
                 {
                     if ( args->gt_filter_str ) error("Error: only one -i or -e expression can be given for gt:, and they cannot be combined\n");
-                    args->gt_filter_str = optarg;
+                    args->gt_filter_str = optarg+3;
                     args->gt_filter_logic |= FLT_EXCLUDE;
                 }
                 else if ( !strncasecmp("qry:",optarg,4) )
                 {
                     if ( args->qry_filter_str ) error("Error: only one -i or -e expression can be given for qry:, and they cannot be combined\n");
-                    args->qry_filter_str = optarg;
+                    args->qry_filter_str = optarg+4;
                     args->qry_filter_logic |= FLT_EXCLUDE;
                 }
                 else
@@ -1355,13 +1355,13 @@ int main_vcfgtcheck(int argc, char *argv[])
                 if ( !strncasecmp("gt:",optarg,3) )
                 {
                     if ( args->gt_filter_str ) error("Error: only one -i or -e expression can be given for gt:, and they cannot be combined\n");
-                    args->gt_filter_str = optarg;
+                    args->gt_filter_str = optarg+3;
                     args->gt_filter_logic |= FLT_INCLUDE;
                 }
                 else if ( !strncasecmp("qry:",optarg,4) )
                 {
                     if ( args->qry_filter_str ) error("Error: only one -i or -e expression can be given for qry:, and they cannot be combined\n");
-                    args->qry_filter_str = optarg;
+                    args->qry_filter_str = optarg+4;
                     args->qry_filter_logic |= FLT_INCLUDE;
                 }
                 else

vcfnorm.c

@@ -2328,8 +2328,8 @@ static void destroy_data(args_t *args)
     if ( args->mseq ) free(args->seq);
 }
 
-
-static void normalize_line(args_t *args, bcf1_t *line)
+// return 0 on success, -1 if line was skipped (due to ref mismatch)
+static int normalize_line(args_t *args, bcf1_t *line)
 {
     if ( args->fai )
     {
@@ -2346,7 +2346,7 @@ static void normalize_line(args_t *args, bcf1_t *line)
             if ( ret==ERR_REF_MISMATCH && args->check_ref & CHECK_REF_SKIP )
             {
                 args->nskipped++;
-                return;
+                return -1;
             }
             if ( ret==ERR_DUP_ALLELE )
             {
@@ -2388,9 +2388,10 @@ static void normalize_line(args_t *args, bcf1_t *line)
         }
         if ( !args->filter_pass || args->atomize!=SPLIT ) break;
     }
+    return 0;
 }
 
-// return 0 on success, 1 when done
+// return 0 on success, 1 when done, -1 if line skipped
 static int split_and_normalize(args_t *args)
 {
     if ( !bcf_sr_next_line(args->files) ) return 1;
@@ -2408,8 +2409,7 @@ static int split_and_normalize(args_t *args)
     if ( args->mrows_op!=MROWS_SPLIT || line->n_allele<=2 || !args->filter_pass )
     {
         // normal operation, no splitting
-        normalize_line(args, line);
-        return 0;
+        return normalize_line(args, line);
     }
 
     // any restrictions on variant types to split?
@@ -2418,8 +2418,7 @@ static int split_and_normalize(args_t *args)
         int type = args->mrows_collapse==COLLAPSE_SNPS ? VCF_SNP : VCF_INDEL;
         if ( !(bcf_get_variant_types(line) & type) )
         {
-            normalize_line(args, line);
-            return 0;
+            return normalize_line(args, line);
         }
     }
 
@@ -2428,7 +2427,11 @@ static int split_and_normalize(args_t *args)
 
     int j;
     for (j=0; j<args->ntmp_lines; j++)
-        normalize_line(args, args->tmp_lines[j]);
+    {
+        int ret = normalize_line(args, args->tmp_lines[j]);
+        if (ret)
+            return ret;
+    }
 
     return 0;
 }
@@ -2453,7 +2456,10 @@ static void normalize_vcf(args_t *args)
         int done = 0;
         while (1)
         {
-            done = split_and_normalize(args);
+            do
+            {
+                done = split_and_normalize(args);
+            } while (done < 0); // Skipped line
             if ( done ) break;      // no more lines available
             int i = args->rbuf.f;
             int j = rbuf_last(&args->rbuf);