|
2 | 2 |
|
3 | 3 | /* The MIT License |
4 | 4 |
|
5 | | - Copyright (c) 2014-2024 Genome Research Ltd. |
| 5 | + Copyright (c) 2014-2025 Genome Research Ltd. |
6 | 6 |
|
7 | 7 | Author: Petr Danecek <[email protected]> |
8 | 8 |
|
@@ -230,24 +230,24 @@ static void init_data(args_t *args) |
230 | 230 | if ( !bcf_sr_add_reader(args->files,args->fname) ) error("Failed to read from %s: %s\n", !strcmp("-",args->fname)?"standard input":args->fname, bcf_sr_strerror(args->files->errnum)); |
231 | 231 | args->hdr = args->files->readers[0].header; |
232 | 232 | args->isample = -1; |
233 | | - if ( !args->sample ) |
| 233 | + if ( args->sample_fname ) |
234 | 234 | { |
235 | | - args->smpl = smpl_ilist_init(args->hdr,NULL,0,SMPL_NONE|SMPL_VERBOSE); |
236 | | - if ( !args->smpl->n ) |
237 | | - { |
238 | | - smpl_ilist_destroy(args->smpl); |
239 | | - args->smpl = NULL; |
240 | | - } |
| 235 | + args->smpl = smpl_ilist_init(args->hdr,args->sample_fname,1,SMPL_NONE|SMPL_VERBOSE); |
| 236 | + if ( args->smpl && !args->smpl->n ) error("No matching sample found\n"); |
241 | 237 | } |
242 | 238 | else if ( args->sample && strcmp("-",args->sample) ) |
243 | 239 | { |
244 | 240 | args->smpl = smpl_ilist_init(args->hdr,args->sample,0,SMPL_NONE|SMPL_VERBOSE); |
245 | 241 | if ( args->smpl && !args->smpl->n ) error("No matching sample found\n"); |
246 | 242 | } |
247 | | - else if ( args->sample_fname ) |
| 243 | + else if ( !args->sample ) |
248 | 244 | { |
249 | | - args->smpl = smpl_ilist_init(args->hdr,args->sample_fname,1,SMPL_NONE|SMPL_VERBOSE); |
250 | | - if ( args->smpl && !args->smpl->n ) error("No matching sample found\n"); |
| 245 | + args->smpl = smpl_ilist_init(args->hdr,NULL,0,SMPL_NONE|SMPL_VERBOSE); |
| 246 | + if ( !args->smpl->n ) |
| 247 | + { |
| 248 | + smpl_ilist_destroy(args->smpl); |
| 249 | + args->smpl = NULL; |
| 250 | + } |
251 | 251 | } |
252 | 252 | if ( args->smpl ) |
253 | 253 | { |
@@ -770,12 +770,26 @@ static void apply_variant(args_t *args, bcf1_t *rec) |
770 | 770 | } |
771 | 771 | if ( ialt==-1 ) |
772 | 772 | { |
773 | | - char alleles[4]; |
774 | | - alleles[0] = rec->d.allele[0][0]; |
775 | | - alleles[1] = ','; |
776 | | - alleles[2] = args->missing_allele; |
777 | | - alleles[3] = 0; |
778 | | - bcf_update_alleles_str(args->hdr, rec, alleles); |
| 773 | + // missing allele, it can be a single position or an entire gvcf block |
| 774 | + if ( rec->rlen>1 && bcf_has_variant_types(rec,VCF_REF,bcf_match_exact)>0 ) |
| 775 | + { |
| 776 | + kstring_t str = {0,0,0}; |
| 777 | + int idx = rec->pos - args->fa_ori_pos + args->fa_mod_off; // position of the variant within the modified fasta sequence |
| 778 | + kputsn(args->fa_buf.s+idx,rec->rlen, &str); |
| 779 | + kputc(',', &str); |
| 780 | + for (i=0; i<rec->rlen; i++) kputc(args->missing_allele, &str); |
| 781 | + bcf_update_alleles_str(args->hdr, rec, str.s); |
| 782 | + free(str.s); |
| 783 | + } |
| 784 | + else |
| 785 | + { |
| 786 | + char alleles[4]; |
| 787 | + alleles[0] = rec->d.allele[0][0]; |
| 788 | + alleles[1] = ','; |
| 789 | + alleles[2] = args->missing_allele; |
| 790 | + alleles[3] = 0; |
| 791 | + bcf_update_alleles_str(args->hdr, rec, alleles); |
| 792 | + } |
779 | 793 | ialt = 1; |
780 | 794 | } |
781 | 795 |
|
@@ -1205,6 +1219,7 @@ static void usage(args_t *args) |
1205 | 1219 | fprintf(bcftools_stderr, " --regions-overlap 0|1|2 Include if POS in the region (0), record overlaps (1), variant overlaps (2) [1]\n"); |
1206 | 1220 | fprintf(bcftools_stderr, " -s, --samples LIST Comma-separated list of samples to include, \"-\" to ignore samples and use REF,ALT\n"); |
1207 | 1221 | fprintf(bcftools_stderr, " -S, --samples-file FILE File of samples to include\n"); |
| 1222 | + fprintf(bcftools_stderr, " -v, --verbosity INT Verbosity level\n"); |
1208 | 1223 | fprintf(bcftools_stderr, "Examples:\n"); |
1209 | 1224 | fprintf(bcftools_stderr, " # Get the consensus for one region. The fasta header lines are then expected\n"); |
1210 | 1225 | fprintf(bcftools_stderr, " # in the form \">chr:from-to\".\n"); |
@@ -1242,13 +1257,17 @@ int main_consensus(int argc, char *argv[]) |
1242 | 1257 | {"chain",1,0,'c'}, |
1243 | 1258 | {"prefix",required_argument,0,'p'}, |
1244 | 1259 | {"regions-overlap",required_argument,0,5}, |
| 1260 | + {"verbosity",required_argument,NULL,'v'}, |
1245 | 1261 | {0,0,0,0} |
1246 | 1262 | }; |
1247 | 1263 | int c; |
1248 | | - while ((c = getopt_long(argc, argv, "h?s:S:1Ii:e:H:f:o:m:c:M:p:a:",loptions,NULL)) >= 0) |
| 1264 | + while ((c = getopt_long(argc, argv, "h?s:S:1Ii:e:H:f:o:m:c:M:p:a:v:",loptions,NULL)) >= 0) |
1249 | 1265 | { |
1250 | 1266 | switch (c) |
1251 | 1267 | { |
| 1268 | + case 'v': |
| 1269 | + if ( apply_verbosity(optarg) < 0 ) error("Could not parse argument: --verbosity %s\n", optarg); |
| 1270 | + break; |
1252 | 1271 | case 1 : args->mark_del = optarg[0]; break; |
1253 | 1272 | case 2 : |
1254 | 1273 | if ( !strcasecmp(optarg,"uc") ) args->mark_ins = TO_UPPER; |
|
0 commit comments