7df6e18265341f87a69fba808aa1f92f8ebca841
markd
  Wed Apr 15 13:39:42 2026 -0700
move copy of htslib

diff --git src/htslib/bgzf.c src/htslib/bgzf.c
deleted file mode 100644
index 3086c1167d1..00000000000
--- src/htslib/bgzf.c
+++ /dev/null
@@ -1,1127 +0,0 @@
-/* The MIT License
-
-   Copyright (c) 2008 Broad Institute / Massachusetts Institute of Technology
-                 2011, 2012 Attractive Chaos <attractor@live.co.uk>
-   Copyright (C) 2009, 2013, 2014 Genome Research Ltd
-
-   Permission is hereby granted, free of charge, to any person obtaining a copy
-   of this software and associated documentation files (the "Software"), to deal
-   in the Software without restriction, including without limitation the rights
-   to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
-   copies of the Software, and to permit persons to whom the Software is
-   furnished to do so, subject to the following conditions:
-
-   The above copyright notice and this permission notice shall be included in
-   all copies or substantial portions of the Software.
-
-   THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
-   IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
-   FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
-   AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
-   LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
-   OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
-   THE SOFTWARE.
-*/
-
-#include <config.h>
-
-#include <stdio.h>
-#include <stdlib.h>
-#include <string.h>
-#include <errno.h>
-#include <unistd.h>
-#include <assert.h>
-#include <pthread.h>
-#include <sys/types.h>
-#include <inttypes.h>
-
-#include "htslib/hts.h"
-#include "htslib/bgzf.h"
-#include "htslib/hfile.h"
-
-#define BGZF_CACHE
-#define BGZF_MT
-
-#define BLOCK_HEADER_LENGTH 18
-#define BLOCK_FOOTER_LENGTH 8
-
-
-/* BGZF/GZIP header (speciallized from RFC 1952; little endian):
- +---+---+---+---+---+---+---+---+---+---+---+---+---+---+---+---+---+---+
- | 31|139|  8|  4|              0|  0|255|      6| 66| 67|      2|BLK_LEN|
- +---+---+---+---+---+---+---+---+---+---+---+---+---+---+---+---+---+---+
-  BGZF extension:
-                ^                              ^   ^   ^
-                |                              |   |   |
-               FLG.EXTRA                     XLEN  B   C
-
-  BGZF format is compatible with GZIP. It limits the size of each compressed
-  block to 2^16 bytes and adds and an extra "BC" field in the gzip header which
-  records the size.
-
-*/
-static const uint8_t g_magic[19] = "\037\213\010\4\0\0\0\0\0\377\6\0\102\103\2\0\0\0";
-
-#ifdef BGZF_CACHE
-typedef struct {
-    int size;
-    uint8_t *block;
-    int64_t end_offset;
-} cache_t;
-#include "htslib/khash.h"
-KHASH_MAP_INIT_INT64(cache, cache_t)
-#endif
-
-typedef struct
-{
-    uint64_t uaddr;  // offset w.r.t. uncompressed data
-    uint64_t caddr;  // offset w.r.t. compressed data
-}
-bgzidx1_t;
-
-struct __bgzidx_t
-{
-    int noffs, moffs;       // the size of the index, n:used, m:allocated
-    bgzidx1_t *offs;        // offsets
-    uint64_t ublock_addr;   // offset of the current block (uncompressed data)
-};
-
-void bgzf_index_destroy(BGZF *fp);
-int bgzf_index_add_block(BGZF *fp);
-
-static inline void packInt16(uint8_t *buffer, uint16_t value)
-{
-    buffer[0] = value;
-    buffer[1] = value >> 8;
-}
-
-static inline int unpackInt16(const uint8_t *buffer)
-{
-    return buffer[0] | buffer[1] << 8;
-}
-
-static inline void packInt32(uint8_t *buffer, uint32_t value)
-{
-    buffer[0] = value;
-    buffer[1] = value >> 8;
-    buffer[2] = value >> 16;
-    buffer[3] = value >> 24;
-}
-
-static BGZF *bgzf_read_init(hFILE *hfpr)
-{
-    BGZF *fp;
-    uint8_t magic[18];
-    ssize_t n = hpeek(hfpr, magic, 18);
-    if (n < 0) return NULL;
-
-    fp = (BGZF*)calloc(1, sizeof(BGZF));
-    if (fp == NULL) return NULL;
-
-    fp->is_write = 0;
-    fp->is_compressed = (n==2 && magic[0]==0x1f && magic[1]==0x8b);
-    fp->uncompressed_block = malloc(BGZF_MAX_BLOCK_SIZE);
-    fp->compressed_block = malloc(BGZF_MAX_BLOCK_SIZE);
-    fp->is_compressed = (n==18 && magic[0]==0x1f && magic[1]==0x8b) ? 1 : 0;
-    fp->is_gzip = ( !fp->is_compressed || ((magic[3]&4) && memcmp(&magic[12], "BC\2\0",4)==0) ) ? 0 : 1;
-#ifdef BGZF_CACHE
-    fp->cache = kh_init(cache);
-#endif
-    return fp;
-}
-
-// get the compress level from the mode string: compress_level==-1 for the default level, -2 plain uncompressed
-static int mode2level(const char *__restrict mode)
-{
-    int i, compress_level = -1;
-    for (i = 0; mode[i]; ++i)
-        if (mode[i] >= '0' && mode[i] <= '9') break;
-    if (mode[i]) compress_level = (int)mode[i] - '0';
-    if (strchr(mode, 'u')) compress_level = -2;
-    return compress_level;
-}
-static BGZF *bgzf_write_init(const char *mode)
-{
-    BGZF *fp;
-    fp = (BGZF*)calloc(1, sizeof(BGZF));
-    fp->is_write = 1;
-    int compress_level = mode2level(mode);
-    if ( compress_level==-2 )
-    {
-        fp->is_compressed = 0;
-        return fp;
-    }
-    fp->is_compressed = 1;
-    fp->uncompressed_block = malloc(BGZF_MAX_BLOCK_SIZE);
-    fp->compressed_block = malloc(BGZF_MAX_BLOCK_SIZE);
-    fp->compress_level = compress_level < 0? Z_DEFAULT_COMPRESSION : compress_level; // Z_DEFAULT_COMPRESSION==-1
-    if (fp->compress_level > 9) fp->compress_level = Z_DEFAULT_COMPRESSION;
-    if ( strchr(mode,'g') )
-    {
-        // gzip output
-        fp->is_gzip = 1;
-        fp->gz_stream = (z_stream*)calloc(1,sizeof(z_stream));
-        fp->gz_stream->zalloc = NULL;
-        fp->gz_stream->zfree  = NULL;
-        if ( deflateInit2(fp->gz_stream, fp->compress_level, Z_DEFLATED, 15|16, 8, Z_DEFAULT_STRATEGY)!=Z_OK ) return NULL;
-    }
-    return fp;
-}
-
-BGZF *bgzf_open(const char *path, const char *mode)
-{
-    BGZF *fp = 0;
-    assert(compressBound(BGZF_BLOCK_SIZE) < BGZF_MAX_BLOCK_SIZE);
-    if (strchr(mode, 'r')) {
-        hFILE *fpr;
-        if ((fpr = hopen(path, mode)) == 0) return 0;
-        fp = bgzf_read_init(fpr);
-        if (fp == 0) { hclose_abruptly(fpr); return NULL; }
-        fp->fp = fpr;
-    } else if (strchr(mode, 'w') || strchr(mode, 'a')) {
-        hFILE *fpw;
-        if ((fpw = hopen(path, mode)) == 0) return 0;
-        fp = bgzf_write_init(mode);
-        fp->fp = fpw;
-    }
-    else { errno = EINVAL; return 0; }
-
-    fp->is_be = ed_is_big();
-    return fp;
-}
-
-BGZF *bgzf_dopen(int fd, const char *mode)
-{
-    BGZF *fp = 0;
-    assert(compressBound(BGZF_BLOCK_SIZE) < BGZF_MAX_BLOCK_SIZE);
-    if (strchr(mode, 'r')) {
-        hFILE *fpr;
-        if ((fpr = hdopen(fd, mode)) == 0) return 0;
-        fp = bgzf_read_init(fpr);
-        if (fp == 0) { hclose_abruptly(fpr); return NULL; } // FIXME this closes fd
-        fp->fp = fpr;
-    } else if (strchr(mode, 'w') || strchr(mode, 'a')) {
-        hFILE *fpw;
-        if ((fpw = hdopen(fd, mode)) == 0) return 0;
-        fp = bgzf_write_init(mode);
-        fp->fp = fpw;
-    }
-    else { errno = EINVAL; return 0; }
-
-    fp->is_be = ed_is_big();
-    return fp;
-}
-
-BGZF *bgzf_hopen(hFILE *hfp, const char *mode)
-{
-    BGZF *fp = NULL;
-    assert(compressBound(BGZF_BLOCK_SIZE) < BGZF_MAX_BLOCK_SIZE);
-    if (strchr(mode, 'r')) {
-        fp = bgzf_read_init(hfp);
-        if (fp == NULL) return NULL;
-    } else if (strchr(mode, 'w') || strchr(mode, 'a')) {
-        fp = bgzf_write_init(mode);
-    }
-    else { errno = EINVAL; return 0; }
-
-    fp->fp = hfp;
-    fp->is_be = ed_is_big();
-    return fp;
-}
-
-int bgzf_compress(void *_dst, size_t *dlen, const void *src, size_t slen, int level)
-{
-    uint32_t crc;
-    z_stream zs;
-    uint8_t *dst = (uint8_t*)_dst;
-
-    // compress the body
-    zs.zalloc = NULL; zs.zfree = NULL;
-    zs.next_in  = (Bytef*)src;
-    zs.avail_in = slen;
-    zs.next_out = dst + BLOCK_HEADER_LENGTH;
-    zs.avail_out = *dlen - BLOCK_HEADER_LENGTH - BLOCK_FOOTER_LENGTH;
-    if (deflateInit2(&zs, level, Z_DEFLATED, -15, 8, Z_DEFAULT_STRATEGY) != Z_OK) return -1; // -15 to disable zlib header/footer
-    if (deflate(&zs, Z_FINISH) != Z_STREAM_END) return -1;
-    if (deflateEnd(&zs) != Z_OK) return -1;
-    *dlen = zs.total_out + BLOCK_HEADER_LENGTH + BLOCK_FOOTER_LENGTH;
-    // write the header
-    memcpy(dst, g_magic, BLOCK_HEADER_LENGTH); // the last two bytes are a place holder for the length of the block
-    packInt16(&dst[16], *dlen - 1); // write the compressed length; -1 to fit 2 bytes
-    // write the footer
-    crc = crc32(crc32(0L, NULL, 0L), (Bytef*)src, slen);
-    packInt32((uint8_t*)&dst[*dlen - 8], crc);
-    packInt32((uint8_t*)&dst[*dlen - 4], slen);
-    return 0;
-}
-
-static int bgzf_gzip_compress(BGZF *fp, void *_dst, size_t *dlen, const void *src, size_t slen, int level)
-{
-    uint8_t *dst = (uint8_t*)_dst;
-    z_stream *zs = fp->gz_stream;
-    int flush = slen ? Z_NO_FLUSH : Z_FINISH;
-    zs->next_in   = (Bytef*)src;
-    zs->avail_in  = slen;
-    zs->next_out  = dst;
-    zs->avail_out = *dlen;
-    if ( deflate(zs, flush) == Z_STREAM_ERROR ) return -1;
-    *dlen = *dlen - zs->avail_out;
-    return 0;
-}
-
-// Deflate the block in fp->uncompressed_block into fp->compressed_block. Also adds an extra field that stores the compressed block length.
-static int deflate_block(BGZF *fp, int block_length)
-{
-    size_t comp_size = BGZF_MAX_BLOCK_SIZE;
-    int ret;
-    if ( !fp->is_gzip )
-        ret = bgzf_compress(fp->compressed_block, &comp_size, fp->uncompressed_block, block_length, fp->compress_level);
-    else
-        ret = bgzf_gzip_compress(fp, fp->compressed_block, &comp_size, fp->uncompressed_block, block_length, fp->compress_level);
-
-    if ( ret != 0 )
-    {
-        fp->errcode |= BGZF_ERR_ZLIB;
-        return -1;
-    }
-    fp->block_offset = 0;
-    return comp_size;
-}
-
-// Inflate the block in fp->compressed_block into fp->uncompressed_block
-static int inflate_block(BGZF* fp, int block_length)
-{
-    z_stream zs;
-    zs.zalloc = NULL;
-    zs.zfree = NULL;
-    zs.next_in = (Bytef*)fp->compressed_block + 18;
-    zs.avail_in = block_length - 16;
-    zs.next_out = (Bytef*)fp->uncompressed_block;
-    zs.avail_out = BGZF_MAX_BLOCK_SIZE;
-
-    if (inflateInit2(&zs, -15) != Z_OK) {
-        fp->errcode |= BGZF_ERR_ZLIB;
-        return -1;
-    }
-    if (inflate(&zs, Z_FINISH) != Z_STREAM_END) {
-        inflateEnd(&zs);
-        fp->errcode |= BGZF_ERR_ZLIB;
-        return -1;
-    }
-    if (inflateEnd(&zs) != Z_OK) {
-        fp->errcode |= BGZF_ERR_ZLIB;
-        return -1;
-    }
-    return zs.total_out;
-}
-
-static int inflate_gzip_block(BGZF *fp, int cached)
-{
-    int ret = Z_OK;
-    do
-    {
-        if ( !cached && fp->gz_stream->avail_out!=0 )
-        {
-            fp->gz_stream->avail_in = hread(fp->fp, fp->compressed_block, BGZF_BLOCK_SIZE);
-            if ( fp->gz_stream->avail_in<=0 ) return fp->gz_stream->avail_in;
-            if ( fp->gz_stream->avail_in==0 ) break;
-            fp->gz_stream->next_in = fp->compressed_block;
-        }
-        else cached = 0;
-        do
-        {
-            fp->gz_stream->next_out = (Bytef*)fp->uncompressed_block + fp->block_offset;
-            fp->gz_stream->avail_out = BGZF_MAX_BLOCK_SIZE - fp->block_offset;
-            ret = inflate(fp->gz_stream, Z_NO_FLUSH);
-            if ( ret==Z_BUF_ERROR ) continue;   // non-critical error
-            if ( ret<0 ) return -1;
-            unsigned int have = BGZF_MAX_BLOCK_SIZE - fp->gz_stream->avail_out;
-            if ( have ) return have;
-        }
-        while ( fp->gz_stream->avail_out == 0 );
-    }
-    while (ret != Z_STREAM_END);
-    return BGZF_MAX_BLOCK_SIZE - fp->gz_stream->avail_out;
-}
-
-// Returns: 0 on success (BGZF header); -1 on non-BGZF GZIP header; -2 on error
-static int check_header(const uint8_t *header)
-{
-    if ( header[0] != 31 || header[1] != 139 || header[2] != 8 ) return -2;
-    return ((header[3] & 4) != 0
-            && unpackInt16((uint8_t*)&header[10]) == 6
-            && header[12] == 'B' && header[13] == 'C'
-            && unpackInt16((uint8_t*)&header[14]) == 2) ? 0 : -1;
-}
-
-#ifdef BGZF_CACHE
-static void free_cache(BGZF *fp)
-{
-    khint_t k;
-    khash_t(cache) *h = (khash_t(cache)*)fp->cache;
-    if (fp->is_write) return;
-    for (k = kh_begin(h); k < kh_end(h); ++k)
-        if (kh_exist(h, k)) free(kh_val(h, k).block);
-    kh_destroy(cache, h);
-}
-
-static int load_block_from_cache(BGZF *fp, int64_t block_address)
-{
-    khint_t k;
-    cache_t *p;
-    khash_t(cache) *h = (khash_t(cache)*)fp->cache;
-    k = kh_get(cache, h, block_address);
-    if (k == kh_end(h)) return 0;
-    p = &kh_val(h, k);
-    if (fp->block_length != 0) fp->block_offset = 0;
-    fp->block_address = block_address;
-    fp->block_length = p->size;
-    memcpy(fp->uncompressed_block, p->block, BGZF_MAX_BLOCK_SIZE);
-    if ( hseek(fp->fp, p->end_offset, SEEK_SET) < 0 )
-    {
-        // todo: move the error up
-        fprintf(stderr,"Could not hseek to %"PRId64"\n", p->end_offset);
-        exit(1);
-    }
-    return p->size;
-}
-
-static void cache_block(BGZF *fp, int size)
-{
-    int ret;
-    khint_t k;
-    cache_t *p;
-    khash_t(cache) *h = (khash_t(cache)*)fp->cache;
-    if (BGZF_MAX_BLOCK_SIZE >= fp->cache_size) return;
-    if ((kh_size(h) + 1) * BGZF_MAX_BLOCK_SIZE > (uint32_t)fp->cache_size) {
-        /* A better way would be to remove the oldest block in the
-         * cache, but here we remove a random one for simplicity. This
-         * should not have a big impact on performance. */
-        for (k = kh_begin(h); k < kh_end(h); ++k)
-            if (kh_exist(h, k)) break;
-        if (k < kh_end(h)) {
-            free(kh_val(h, k).block);
-            kh_del(cache, h, k);
-        }
-    }
-    k = kh_put(cache, h, fp->block_address, &ret);
-    if (ret == 0) return; // if this happens, a bug!
-    p = &kh_val(h, k);
-    p->size = fp->block_length;
-    p->end_offset = fp->block_address + size;
-    p->block = (uint8_t*)malloc(BGZF_MAX_BLOCK_SIZE);
-    memcpy(kh_val(h, k).block, fp->uncompressed_block, BGZF_MAX_BLOCK_SIZE);
-}
-#else
-static void free_cache(BGZF *fp) {}
-static int load_block_from_cache(BGZF *fp, int64_t block_address) {return 0;}
-static void cache_block(BGZF *fp, int size) {}
-#endif
-
-int bgzf_read_block(BGZF *fp)
-{
-    uint8_t header[BLOCK_HEADER_LENGTH], *compressed_block;
-    int count, size = 0, block_length, remaining;
-
-    // Reading an uncompressed file
-    if ( !fp->is_compressed )
-    {
-        count = hread(fp->fp, fp->uncompressed_block, BGZF_MAX_BLOCK_SIZE);
-        if ( count==0 )
-        {
-            fp->block_length = 0;
-            return 0;
-        }
-        if (fp->block_length != 0) fp->block_offset = 0;
-        fp->block_address += count;
-        fp->block_length = count;
-        return 0;
-    }
-
-    // Reading compressed file
-    int64_t block_address;
-    block_address = htell(fp->fp);
-    if ( fp->is_gzip && fp->gz_stream ) // is this is a initialized gzip stream?
-    {
-        count = inflate_gzip_block(fp, 0);
-        if ( count<0 )
-        {
-            fp->errcode |= BGZF_ERR_ZLIB;
-            return -1;
-        }
-        fp->block_length = count;
-        fp->block_address = block_address;
-        return 0;
-    }
-    if (fp->cache_size && load_block_from_cache(fp, block_address)) return 0;
-    count = hread(fp->fp, header, sizeof(header));
-    if (count == 0) { // no data read
-        fp->block_length = 0;
-        return 0;
-    }
-    int ret;
-    if ( count != sizeof(header) || (ret=check_header(header))==-2 )
-    {
-        fp->errcode |= BGZF_ERR_HEADER;
-        return -1;
-    }
-    if ( ret==-1 )
-    {
-        // GZIP, not BGZF
-        uint8_t *cblock = (uint8_t*)fp->compressed_block;
-        memcpy(cblock, header, sizeof(header));
-        count = hread(fp->fp, cblock+sizeof(header), BGZF_BLOCK_SIZE - sizeof(header)) + sizeof(header);
-        int nskip = 10;
-
-        // Check optional fields to skip: FLG.FNAME,FLG.FCOMMENT,FLG.FHCRC,FLG.FEXTRA
-        // Note: Some of these fields are untested, I did not have appropriate data available
-        if ( header[3] & 0x4 ) // FLG.FEXTRA
-        {
-            nskip += unpackInt16(&cblock[nskip]) + 2;
-        }
-        if ( header[3] & 0x8 ) // FLG.FNAME
-        {
-            while ( nskip<count && cblock[nskip] ) nskip++;
-            nskip++;
-        }
-        if ( header[3] & 0x10 ) // FLG.FCOMMENT
-        {
-            while ( nskip<count && cblock[nskip] ) nskip++;
-            nskip++;
-        }
-        if ( header[3] & 0x2 ) nskip += 2;  //  FLG.FHCRC
-
-        /* FIXME: Should handle this better.  There's no reason why
-           someone shouldn't include a massively long comment in their
-           gzip stream. */
-        if ( nskip >= count )
-        {
-            fp->errcode |= BGZF_ERR_HEADER;
-            return -1;
-        }
-
-        fp->is_gzip = 1;
-        fp->gz_stream = (z_stream*) calloc(1,sizeof(z_stream));
-        int ret = inflateInit2(fp->gz_stream, -15);
-        if (ret != Z_OK)
-        {
-            fp->errcode |= BGZF_ERR_ZLIB;
-            return -1;
-        }
-        fp->gz_stream->avail_in = count - nskip;
-        fp->gz_stream->next_in  = cblock + nskip;
-        count = inflate_gzip_block(fp, 1);
-        if ( count<0 )
-        {
-            fp->errcode |= BGZF_ERR_ZLIB;
-            return -1;
-        }
-        fp->block_length = count;
-        fp->block_address = block_address;
-        if ( fp->idx_build_otf ) return -1; // cannot build index for gzip
-        return 0;
-    }
-    size = count;
-    block_length = unpackInt16((uint8_t*)&header[16]) + 1; // +1 because when writing this number, we used "-1"
-    compressed_block = (uint8_t*)fp->compressed_block;
-    memcpy(compressed_block, header, BLOCK_HEADER_LENGTH);
-    remaining = block_length - BLOCK_HEADER_LENGTH;
-    count = hread(fp->fp, &compressed_block[BLOCK_HEADER_LENGTH], remaining);
-    if (count != remaining) {
-        fp->errcode |= BGZF_ERR_IO;
-        return -1;
-    }
-    size += count;
-    if ((count = inflate_block(fp, block_length)) < 0) return -1;
-    if (fp->block_length != 0) fp->block_offset = 0; // Do not reset offset if this read follows a seek.
-    fp->block_address = block_address;
-    fp->block_length = count;
-    if ( fp->idx_build_otf )
-    {
-        bgzf_index_add_block(fp);
-        fp->idx->ublock_addr += count;
-    }
-    cache_block(fp, size);
-    return 0;
-}
-
-ssize_t bgzf_read(BGZF *fp, void *data, size_t length)
-{
-    ssize_t bytes_read = 0;
-    uint8_t *output = (uint8_t*)data;
-    if (length <= 0) return 0;
-    assert(fp->is_write == 0);
-    while (bytes_read < length) {
-        int copy_length, available = fp->block_length - fp->block_offset;
-        uint8_t *buffer;
-        if (available <= 0) {
-            if (bgzf_read_block(fp) != 0) return -1;
-            available = fp->block_length - fp->block_offset;
-            if (available <= 0) break;
-        }
-        copy_length = length - bytes_read < available? length - bytes_read : available;
-        buffer = (uint8_t*)fp->uncompressed_block;
-        memcpy(output, buffer + fp->block_offset, copy_length);
-        fp->block_offset += copy_length;
-        output += copy_length;
-        bytes_read += copy_length;
-    }
-    if (fp->block_offset == fp->block_length) {
-        fp->block_address = htell(fp->fp);
-        fp->block_offset = fp->block_length = 0;
-    }
-    fp->uncompressed_address += bytes_read;
-    return bytes_read;
-}
-
-ssize_t bgzf_raw_read(BGZF *fp, void *data, size_t length)
-{
-    return hread(fp->fp, data, length);
-}
-
-#ifdef BGZF_MT
-
-typedef struct {
-    struct bgzf_mtaux_t *mt;
-    void *buf;
-    int i, errcode, toproc, compress_level;
-} worker_t;
-
-typedef struct bgzf_mtaux_t {
-    int n_threads, n_blks, curr, done;
-    volatile int proc_cnt;
-    void **blk;
-    int *len;
-    worker_t *w;
-    pthread_t *tid;
-    pthread_mutex_t lock;
-    pthread_cond_t cv;
-} mtaux_t;
-
-static int worker_aux(worker_t *w)
-{
-    int i, stop = 0;
-    // wait for condition: to process or all done
-    pthread_mutex_lock(&w->mt->lock);
-    while (!w->toproc && !w->mt->done)
-        pthread_cond_wait(&w->mt->cv, &w->mt->lock);
-    if (w->mt->done) stop = 1;
-    w->toproc = 0;
-    pthread_mutex_unlock(&w->mt->lock);
-    if (stop) return 1; // to quit the thread
-    w->errcode = 0;
-    for (i = w->i; i < w->mt->curr; i += w->mt->n_threads) {
-        size_t clen = BGZF_MAX_BLOCK_SIZE;
-        if (bgzf_compress(w->buf, &clen, w->mt->blk[i], w->mt->len[i], w->compress_level) != 0)
-            w->errcode |= BGZF_ERR_ZLIB;
-        memcpy(w->mt->blk[i], w->buf, clen);
-        w->mt->len[i] = clen;
-    }
-    __sync_fetch_and_add(&w->mt->proc_cnt, 1);
-    return 0;
-}
-
-static void *mt_worker(void *data)
-{
-    while (worker_aux((worker_t*)data) == 0);
-    return 0;
-}
-
-int bgzf_mt(BGZF *fp, int n_threads, int n_sub_blks)
-{
-    int i;
-    mtaux_t *mt;
-    pthread_attr_t attr;
-    if (!fp->is_write || fp->mt || n_threads <= 1) return -1;
-    mt = (mtaux_t*)calloc(1, sizeof(mtaux_t));
-    mt->n_threads = n_threads;
-    mt->n_blks = n_threads * n_sub_blks;
-    mt->len = (int*)calloc(mt->n_blks, sizeof(int));
-    mt->blk = (void**)calloc(mt->n_blks, sizeof(void*));
-    for (i = 0; i < mt->n_blks; ++i)
-        mt->blk[i] = malloc(BGZF_MAX_BLOCK_SIZE);
-    mt->tid = (pthread_t*)calloc(mt->n_threads, sizeof(pthread_t)); // tid[0] is not used, as the worker 0 is launched by the master
-    mt->w = (worker_t*)calloc(mt->n_threads, sizeof(worker_t));
-    for (i = 0; i < mt->n_threads; ++i) {
-        mt->w[i].i = i;
-        mt->w[i].mt = mt;
-        mt->w[i].compress_level = fp->compress_level;
-        mt->w[i].buf = malloc(BGZF_MAX_BLOCK_SIZE);
-    }
-    pthread_attr_init(&attr);
-    pthread_attr_setdetachstate(&attr, PTHREAD_CREATE_JOINABLE);
-    pthread_mutex_init(&mt->lock, 0);
-    pthread_cond_init(&mt->cv, 0);
-    for (i = 1; i < mt->n_threads; ++i) // worker 0 is effectively launched by the master thread
-        pthread_create(&mt->tid[i], &attr, mt_worker, &mt->w[i]);
-    fp->mt = mt;
-    return 0;
-}
-
-static void mt_destroy(mtaux_t *mt)
-{
-    int i;
-    // signal all workers to quit
-    pthread_mutex_lock(&mt->lock);
-    mt->done = 1; mt->proc_cnt = 0;
-    pthread_cond_broadcast(&mt->cv);
-    pthread_mutex_unlock(&mt->lock);
-    for (i = 1; i < mt->n_threads; ++i) pthread_join(mt->tid[i], 0); // worker 0 is effectively launched by the master thread
-    // free other data allocated on heap
-    for (i = 0; i < mt->n_blks; ++i) free(mt->blk[i]);
-    for (i = 0; i < mt->n_threads; ++i) free(mt->w[i].buf);
-    free(mt->blk); free(mt->len); free(mt->w); free(mt->tid);
-    pthread_cond_destroy(&mt->cv);
-    pthread_mutex_destroy(&mt->lock);
-    free(mt);
-}
-
-static void mt_queue(BGZF *fp)
-{
-    mtaux_t *mt = fp->mt;
-    assert(mt->curr < mt->n_blks); // guaranteed by the caller
-    memcpy(mt->blk[mt->curr], fp->uncompressed_block, fp->block_offset);
-    mt->len[mt->curr] = fp->block_offset;
-    fp->block_offset = 0;
-    ++mt->curr;
-}
-
-static int mt_flush_queue(BGZF *fp)
-{
-    int i;
-    mtaux_t *mt = fp->mt;
-    // signal all the workers to compress
-    pthread_mutex_lock(&mt->lock);
-    for (i = 0; i < mt->n_threads; ++i) mt->w[i].toproc = 1;
-    mt->proc_cnt = 0;
-    pthread_cond_broadcast(&mt->cv);
-    pthread_mutex_unlock(&mt->lock);
-    // worker 0 is doing things here
-    worker_aux(&mt->w[0]);
-    // wait for all the threads to complete
-    while (mt->proc_cnt < mt->n_threads);
-    // dump data to disk
-    for (i = 0; i < mt->n_threads; ++i) fp->errcode |= mt->w[i].errcode;
-    for (i = 0; i < mt->curr; ++i)
-        if (hwrite(fp->fp, mt->blk[i], mt->len[i]) != mt->len[i]) {
-            fp->errcode |= BGZF_ERR_IO;
-            break;
-        }
-    mt->curr = 0;
-    return (fp->errcode == 0)? 0 : -1;
-}
-
-static int lazy_flush(BGZF *fp)
-{
-    if (fp->mt) {
-        if (fp->block_offset) mt_queue(fp);
-        return (fp->mt->curr < fp->mt->n_blks)? 0 : mt_flush_queue(fp);
-    }
-    else return bgzf_flush(fp);
-}
-
-#else  // ~ #ifdef BGZF_MT
-
-int bgzf_mt(BGZF *fp, int n_threads, int n_sub_blks)
-{
-    return 0;
-}
-
-static inline int lazy_flush(BGZF *fp)
-{
-    return bgzf_flush(fp);
-}
-
-#endif // ~ #ifdef BGZF_MT
-
-int bgzf_flush(BGZF *fp)
-{
-    if (!fp->is_write) return 0;
-#ifdef BGZF_MT
-    if (fp->mt) {
-        if (fp->block_offset) mt_queue(fp); // guaranteed that assertion does not fail
-        return mt_flush_queue(fp);
-    }
-#endif
-    while (fp->block_offset > 0) {
-        if ( fp->idx_build_otf )
-        {
-            bgzf_index_add_block(fp);
-            fp->idx->ublock_addr += fp->block_offset;
-        }
-        int block_length = deflate_block(fp, fp->block_offset);
-        if (block_length < 0) return -1;
-        if (hwrite(fp->fp, fp->compressed_block, block_length) != block_length) {
-            fp->errcode |= BGZF_ERR_IO; // possibly truncated file
-            return -1;
-        }
-        fp->block_address += block_length;
-    }
-    return 0;
-}
-
-int bgzf_flush_try(BGZF *fp, ssize_t size)
-{
-    if (fp->block_offset + size > BGZF_BLOCK_SIZE) return lazy_flush(fp);
-    return 0;
-}
-
-ssize_t bgzf_write(BGZF *fp, const void *data, size_t length)
-{
-    if ( !fp->is_compressed )
-        return hwrite(fp->fp, data, length);
-
-    const uint8_t *input = (const uint8_t*)data;
-    ssize_t remaining = length;
-    assert(fp->is_write);
-    while (remaining > 0) {
-        uint8_t* buffer = (uint8_t*)fp->uncompressed_block;
-        int copy_length = BGZF_BLOCK_SIZE - fp->block_offset;
-        if (copy_length > remaining) copy_length = remaining;
-        memcpy(buffer + fp->block_offset, input, copy_length);
-        fp->block_offset += copy_length;
-        input += copy_length;
-        remaining -= copy_length;
-        if (fp->block_offset == BGZF_BLOCK_SIZE) {
-            if (lazy_flush(fp) != 0) return -1;
-        }
-    }
-    return length - remaining;
-}
-
-ssize_t bgzf_raw_write(BGZF *fp, const void *data, size_t length)
-{
-    return hwrite(fp->fp, data, length);
-}
-
-int bgzf_close(BGZF* fp)
-{
-    int ret, block_length;
-    if (fp == 0) return -1;
-    if (fp->is_write && fp->is_compressed) {
-        if (bgzf_flush(fp) != 0) return -1;
-        fp->compress_level = -1;
-        block_length = deflate_block(fp, 0); // write an empty block
-        if (hwrite(fp->fp, fp->compressed_block, block_length) < 0
-            || hflush(fp->fp) != 0) {
-            fp->errcode |= BGZF_ERR_IO;
-            return -1;
-        }
-#ifdef BGZF_MT
-        if (fp->mt) mt_destroy(fp->mt);
-#endif
-    }
-    if ( fp->is_gzip )
-    {
-        if (!fp->is_write) (void)inflateEnd(fp->gz_stream);
-        else (void)deflateEnd(fp->gz_stream);
-        free(fp->gz_stream);
-    }
-    ret = hclose(fp->fp);
-    if (ret != 0) return -1;
-    bgzf_index_destroy(fp);
-    free(fp->uncompressed_block);
-    free(fp->compressed_block);
-    free_cache(fp);
-    free(fp);
-    return 0;
-}
-
-void bgzf_set_cache_size(BGZF *fp, int cache_size)
-{
-    if (fp) fp->cache_size = cache_size;
-}
-
-#ifndef UCSC_CRAM
-int bgzf_check_EOF(BGZF *fp)
-{
-    uint8_t buf[28];
-    off_t offset = htell(fp->fp);
-    if (hseek(fp->fp, -28, SEEK_END) < 0) {
-        if (errno == ESPIPE) { hclearerr(fp->fp); return 2; }
-        else return -1;
-    }
-    if ( hread(fp->fp, buf, 28) != 28 ) return -1;
-    if ( hseek(fp->fp, offset, SEEK_SET) < 0 ) return -1;
-    return (memcmp("\037\213\010\4\0\0\0\0\0\377\6\0\102\103\2\0\033\0\3\0\0\0\0\0\0\0\0\0", buf, 28) == 0)? 1 : 0;
-}
-#endif
-
-int64_t bgzf_seek(BGZF* fp, int64_t pos, int where)
-{
-    int block_offset;
-    int64_t block_address;
-
-    if (fp->is_write || where != SEEK_SET) {
-        fp->errcode |= BGZF_ERR_MISUSE;
-        return -1;
-    }
-    block_offset = pos & 0xFFFF;
-    block_address = pos >> 16;
-    if (hseek(fp->fp, block_address, SEEK_SET) < 0) {
-        fp->errcode |= BGZF_ERR_IO;
-        return -1;
-    }
-    fp->block_length = 0;  // indicates current block has not been loaded
-    fp->block_address = block_address;
-    fp->block_offset = block_offset;
-    return 0;
-}
-
-int bgzf_is_bgzf(const char *fn)
-{
-    uint8_t buf[16];
-    int n;
-    hFILE *fp;
-    if ((fp = hopen(fn, "r")) == 0) return 0;
-    n = hread(fp, buf, 16);
-    if ( hclose(fp) < 0 ) return -1;
-    if (n != 16) return 0;
-    return memcmp(g_magic, buf, 16) == 0? 1 : 0;
-}
-
-int bgzf_getc(BGZF *fp)
-{
-    int c;
-    if (fp->block_offset >= fp->block_length) {
-        if (bgzf_read_block(fp) != 0) return -2; /* error */
-        if (fp->block_length == 0) return -1; /* end-of-file */
-    }
-    c = ((unsigned char*)fp->uncompressed_block)[fp->block_offset++];
-    if (fp->block_offset == fp->block_length) {
-        fp->block_address = htell(fp->fp);
-        fp->block_offset = 0;
-        fp->block_length = 0;
-    }
-    fp->uncompressed_address++;
-    return c;
-}
-
-#ifndef kroundup32
-#define kroundup32(x) (--(x), (x)|=(x)>>1, (x)|=(x)>>2, (x)|=(x)>>4, (x)|=(x)>>8, (x)|=(x)>>16, ++(x))
-#endif
-
-int bgzf_getline(BGZF *fp, int delim, kstring_t *str)
-{
-    int l, state = 0;
-    unsigned char *buf = (unsigned char*)fp->uncompressed_block;
-    str->l = 0;
-    do {
-        if (fp->block_offset >= fp->block_length) {
-            if (bgzf_read_block(fp) != 0) { state = -2; break; }
-            if (fp->block_length == 0) { state = -1; break; }
-        }
-        for (l = fp->block_offset; l < fp->block_length && buf[l] != delim; ++l);
-        if (l < fp->block_length) state = 1;
-        l -= fp->block_offset;
-        if (str->l + l + 1 >= str->m) {
-            str->m = str->l + l + 2;
-            kroundup32(str->m);
-            str->s = (char*)realloc(str->s, str->m);
-        }
-        memcpy(str->s + str->l, buf + fp->block_offset, l);
-        str->l += l;
-        fp->block_offset += l + 1;
-        if (fp->block_offset >= fp->block_length) {
-            fp->block_address = htell(fp->fp);
-            fp->block_offset = 0;
-            fp->block_length = 0;
-        }
-    } while (state == 0);
-    if (str->l == 0 && state < 0) return state;
-    fp->uncompressed_address += str->l;
-    if ( delim=='\n' && str->l>0 && str->s[str->l-1]=='\r' ) str->l--;
-    str->s[str->l] = 0;
-    return str->l;
-}
-
-void bgzf_index_destroy(BGZF *fp)
-{
-    if ( !fp->idx ) return;
-    free(fp->idx->offs);
-    free(fp->idx);
-    fp->idx = NULL;
-    fp->idx_build_otf = 0;
-}
-
-int bgzf_index_build_init(BGZF *fp)
-{
-    bgzf_index_destroy(fp);
-    fp->idx = (bgzidx_t*) calloc(1,sizeof(bgzidx_t));
-    if ( !fp->idx ) return -1;
-    fp->idx_build_otf = 1;  // build index on the fly
-    return 0;
-}
-
-int bgzf_index_add_block(BGZF *fp)
-{
-    fp->idx->noffs++;
-    if ( fp->idx->noffs > fp->idx->moffs )
-    {
-        fp->idx->moffs = fp->idx->noffs;
-        kroundup32(fp->idx->moffs);
-        fp->idx->offs = (bgzidx1_t*) realloc(fp->idx->offs, fp->idx->moffs*sizeof(bgzidx1_t));
-        if ( !fp->idx->offs ) return -1;
-    }
-    fp->idx->offs[ fp->idx->noffs-1 ].uaddr = fp->idx->ublock_addr;
-    fp->idx->offs[ fp->idx->noffs-1 ].caddr = fp->block_address;
-    return 0;
-}
-
-int bgzf_index_dump(BGZF *fp, const char *bname, const char *suffix)
-{
-    if (bgzf_flush(fp) != 0) return -1;
-
-    assert(fp->idx);
-    char *tmp = NULL;
-    if ( suffix )
-    {
-        int blen = strlen(bname);
-        int slen = strlen(suffix);
-        tmp = (char*) malloc(blen + slen + 1);
-        if ( !tmp ) return -1;
-        memcpy(tmp,bname,blen);
-        memcpy(tmp+blen,suffix,slen+1);
-    }
-
-    FILE *idx = fopen(tmp?tmp:bname,"wb");
-    if ( tmp ) free(tmp);
-    if ( !idx ) return -1;
-
-    // Note that the index contains one extra record when indexing files opened
-    // for reading. The terminating record is not present when opened for writing.
-    // This is not a bug.
-
-    int i;
-    if ( fp->is_be )
-    {
-        uint64_t x = fp->idx->noffs - 1;
-        fwrite(ed_swap_8p(&x), 1, sizeof(x), idx);
-        for (i=1; i<fp->idx->noffs; i++)
-        {
-            x = fp->idx->offs[i].caddr; fwrite(ed_swap_8p(&x), 1, sizeof(x), idx);
-            x = fp->idx->offs[i].uaddr; fwrite(ed_swap_8p(&x), 1, sizeof(x), idx);
-        }
-    }
-    else
-    {
-        uint64_t x = fp->idx->noffs - 1;
-        fwrite(&x, 1, sizeof(x), idx);
-        for (i=1; i<fp->idx->noffs; i++)
-        {
-            fwrite(&fp->idx->offs[i].caddr, 1, sizeof(fp->idx->offs[i].caddr), idx);
-            fwrite(&fp->idx->offs[i].uaddr, 1, sizeof(fp->idx->offs[i].uaddr), idx);
-        }
-    }
-    fclose(idx);
-    return 0;
-}
-
-
-int bgzf_index_load(BGZF *fp, const char *bname, const char *suffix)
-{
-    char *tmp = NULL;
-    if ( suffix )
-    {
-        int blen = strlen(bname);
-        int slen = strlen(suffix);
-        tmp = (char*) malloc(blen + slen + 1);
-        if ( !tmp ) return -1;
-        memcpy(tmp,bname,blen);
-        memcpy(tmp+blen,suffix,slen+1);
-    }
-
-    FILE *idx = fopen(tmp?tmp:bname,"rb");
-    if ( tmp ) free(tmp);
-    if ( !idx ) return -1;
-
-    fp->idx = (bgzidx_t*) calloc(1,sizeof(bgzidx_t));
-    uint64_t x;
-    if ( fread(&x, 1, sizeof(x), idx) != sizeof(x) ) return -1;
-
-    fp->idx->noffs = fp->idx->moffs = 1 + (fp->is_be ? ed_swap_8(x) : x);
-    fp->idx->offs  = (bgzidx1_t*) malloc(fp->idx->moffs*sizeof(bgzidx1_t));
-    fp->idx->offs[0].caddr = fp->idx->offs[0].uaddr = 0;
-
-    int i;
-    if ( fp->is_be )
-    {
-        int ret = 0;
-        for (i=1; i<fp->idx->noffs; i++)
-        {
-            ret += fread(&x, 1, sizeof(x), idx); fp->idx->offs[i].caddr = ed_swap_8(x);
-            ret += fread(&x, 1, sizeof(x), idx); fp->idx->offs[i].uaddr = ed_swap_8(x);
-        }
-        if ( ret != sizeof(x)*2*(fp->idx->noffs-1) ) return -1;
-    }
-    else
-    {
-        int ret = 0;
-        for (i=1; i<fp->idx->noffs; i++)
-        {
-            ret += fread(&x, 1, sizeof(x), idx); fp->idx->offs[i].caddr = x;
-            ret += fread(&x, 1, sizeof(x), idx); fp->idx->offs[i].uaddr = x;
-        }
-        if ( ret != sizeof(x)*2*(fp->idx->noffs-1) ) return -1;
-    }
-    fclose(idx);
-    return 0;
-
-}
-
-int bgzf_useek(BGZF *fp, long uoffset, int where)
-{
-    if ( !fp->is_compressed )
-    {
-        if (hseek(fp->fp, uoffset, SEEK_SET) < 0)
-        {
-            fp->errcode |= BGZF_ERR_IO;
-            return -1;
-        }
-        fp->block_length = 0;  // indicates current block has not been loaded
-        fp->block_address = uoffset;
-        fp->block_offset = 0;
-        bgzf_read_block(fp);
-        fp->uncompressed_address = uoffset;
-        return 0;
-    }
-
-    if ( !fp->idx )
-    {
-        fp->errcode |= BGZF_ERR_IO;
-        return -1;
-    }
-
-    // binary search
-    int ilo = 0, ihi = fp->idx->noffs - 1;
-    while ( ilo<=ihi )
-    {
-        int i = (ilo+ihi)*0.5;
-        if ( uoffset < fp->idx->offs[i].uaddr ) ihi = i - 1;
-        else if ( uoffset >= fp->idx->offs[i].uaddr ) ilo = i + 1;
-        else break;
-    }
-    int i = ilo-1;
-    if (hseek(fp->fp, fp->idx->offs[i].caddr, SEEK_SET) < 0)
-    {
-        fp->errcode |= BGZF_ERR_IO;
-        return -1;
-    }
-    fp->block_length = 0;  // indicates current block has not been loaded
-    fp->block_address = fp->idx->offs[i].caddr;
-    fp->block_offset = 0;
-    if ( bgzf_read_block(fp) < 0 ) return -1;
-    if ( uoffset - fp->idx->offs[i].uaddr > 0 )
-    {
-        fp->block_offset = uoffset - fp->idx->offs[i].uaddr;
-        assert( fp->block_offset <= fp->block_length );     // todo: skipped, unindexed, blocks
-    }
-    fp->uncompressed_address = uoffset;
-    return 0;
-}
-
-long bgzf_utell(BGZF *fp)
-{
-    return fp->uncompressed_address;    // currently maintained only when reading
-}
-