How to use C++ htslib to read VCF contig name and size?

Question

A typical VCF file has:

##contig=<ID=chr1,length=248956422>
##contig=<ID=chr2,length=242193529>

I would like to use htslib in C++ to read it. My attempt:

htsFile *fp = bcf_open("my.vcf", "r");
bcf_hdr_t *hdr = bcf_hdr_read(fp);

In https://github.com/samtools/htslib/blob/develop/htslib/vcf.h, I'm not able to find a function that can do that for me.

How to read chr1 and 248956422 in C++?

As far as I can tell, the (internal) function bcf_hdr_parse_line does that. It returns an object x of type bcf_hrec_t*, for which strcmp(x->key, "contig") == 0) for your lines. You can then access the x->keys and x->vals. Untested. — Konrad Rudolph, Feb 28 '18 at 10:49
@KonradRudolph Can I use it externally properly? Without hacking? — SmallChess, Feb 28 '18 at 11:03
Depends on what you mean by “properly”. Since the API isn’t properly documented anyway, the distinction between external and internal functions is somewhat fluid. It might break with future versions of htslib, if that’s what you’re concerned about. — Konrad Rudolph, Feb 28 '18 at 11:13
@KonradRudolph Yeah. Breaking with future versions is my concern. If I see don't any other alternative, I might go with that or just parse the file myself. — SmallChess, Feb 28 '18 at 11:14

score 4 · Accepted Answer · answered Feb 28 '18 at 14:10

#include "htslib/vcf.h"

int main(int argc, char *argv[])
{
    htsFile *fp;
    bcf_hdr_t *hdr;
    bcf_idpair_t *ctg;
    int i;
    if (argc == 1) {
        fprintf(stderr, "Usage: print-ctg <in.vcf>\n");
        return 1;
    }
    fp = vcf_open(argv[1], "r");
    hdr = vcf_hdr_read(fp);
    ctg = hdr->id[BCF_DT_CTG];
    for (i = 0; i < hdr->n[BCF_DT_CTG]; ++i)
        printf("%s\t%d\n", ctg[i].key, ctg[i].val->info[0]);
    bcf_hdr_destroy(hdr);
    vcf_close(fp);
    return 0;
}

On stability: this use has been in htslib forever. In general, functions/structs/variables in the public headers are meant to be stable. However, there is no guarantee that future versions will always keep the same APIs.

How to use C++ htslib to read VCF contig name and size?

1 Answers1