Skip to content

Commits on Source 5

Expand all Show whitespace changes
Inline Side-by-side
LICENSE
View file @ 18118746
The MIT/Expat License
Copyright (C) 2008-2018 Genome Research Ltd.
Copyright (C) 2008-2019 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
......
Makefile
View file @ 18118746
# Makefile for samtools, utilities for the Sequence Alignment/Map format.
#
# Copyright (C) 2008-2017 Genome Research Ltd.
# Copyright (C) 2008-2019 Genome Research Ltd.
# Portions copyright (C) 2010-2012 Broad Institute.
#
# Permission is hereby granted, free of charge, to any person obtaining a copy
......@@ -35,13 +35,13 @@ LZ4_LDFLAGS = -L$(LZ4DIR)
LOBJS= bam_aux.o bam.o bam_import.o sam.o \
sam_header.o bam_plbuf.o
AOBJS= bam_index.o bam_plcmd.o sam_view.o \
LOBJS= bam_aux.o bam.o sam.o \
bam_plbuf.o
AOBJS= bam_index.o bam_plcmd.o sam_view.o bam_fastq.o \
bam_cat.o bam_md.o bam_reheader.o bam_sort.o bedidx.o \
bam_rmdup.o bam_rmdupse.o bam_mate.o bam_stat.o bam_color.o \
bamtk.o bam2bcf.o bam2bcf_indel.o sample.o \
cut_target.o phase.o bam2depth.o padding.o bedcov.o bamshuf.o \
cut_target.o phase.o bam2depth.o coverage.o padding.o bedcov.o bamshuf.o \
faidx.o dict.o stats.o stats_isize.o bam_flags.o bam_split.o \
bam_tview.o bam_tview_curses.o bam_tview_html.o bam_lpileup.o \
bam_quickcheck.o bam_addrprg.o bam_markdup.o tmp_file.o
......@@ -77,8 +77,7 @@ MISC_SCRIPTS = \
misc/interpolate_sam.pl misc/novo2sam.pl \
misc/plot-bamstats misc/psl2sam.pl \
misc/sam2vcf.pl misc/samtools.pl misc/seq_cache_populate.pl \
misc/soap2sam.pl \
misc/varfilter.py misc/wgsim_eval.pl misc/zoom2sam.pl
misc/soap2sam.pl misc/wgsim_eval.pl misc/zoom2sam.pl
TEST_PROGRAMS = \
test/merge/test_bam_translate \
......@@ -158,55 +157,62 @@ bedidx_h = bedidx.h $(htslib_hts_h)
sam_h = sam.h $(htslib_sam_h) $(bam_h)
sam_opts_h = sam_opts.h $(htslib_hts_h)
sample_h = sample.h $(htslib_kstring_h)
samtools_h = samtools.h $(htslib_hts_defs_h) $(htslib_sam_h)
stats_isize_h = stats_isize.h $(htslib_khash_h)
tmp_file_h = tmp_file.h $(htslib_sam_h) $(LZ4DIR)/lz4.h
bam.o: bam.c config.h $(bam_h) $(htslib_kstring_h) sam_header.h
bam.o: bam.c config.h $(bam_h) $(htslib_kstring_h)
bam2bcf.o: bam2bcf.c config.h $(htslib_hts_h) $(htslib_sam_h) $(htslib_kstring_h) $(htslib_kfunc_h) $(bam2bcf_h)
bam2bcf_indel.o: bam2bcf_indel.c config.h $(htslib_hts_h) $(htslib_sam_h) $(bam2bcf_h) $(htslib_khash_h) $(htslib_ksort_h)
bam2depth.o: bam2depth.c config.h $(htslib_sam_h) samtools.h $(sam_opts_h)
bam_addrprg.o: bam_addrprg.c config.h $(htslib_sam_h) $(htslib_kstring_h) samtools.h $(htslib_thread_pool_h) $(sam_opts_h)
bam2depth.o: bam2depth.c config.h $(htslib_sam_h) $(samtools_h) $(bedidx_h) $(sam_opts_h)
coverage.o: coverage.c config.h $(htslib_sam_h) $(htslib_hts_h) $(samtools_h) $(sam_opts_h)
bam_addrprg.o: bam_addrprg.c config.h $(htslib_sam_h) $(htslib_kstring_h) $(samtools_h) $(htslib_thread_pool_h) $(sam_opts_h)
bam_aux.o: bam_aux.c config.h $(bam_h)
bam_cat.o: bam_cat.c config.h $(htslib_bgzf_h) $(htslib_sam_h) $(htslib_cram_h) $(htslib_khash_h) samtools.h
bam_cat.o: bam_cat.c config.h $(htslib_bgzf_h) $(htslib_sam_h) $(htslib_cram_h) $(htslib_kstring_h) $(samtools_h) $(sam_opts_h)
bam_color.o: bam_color.c config.h $(bam_h)
bam_import.o: bam_import.c config.h $(htslib_kstring_h) $(bam_h) $(htslib_kseq_h)
bam_index.o: bam_index.c config.h $(htslib_hts_h) $(htslib_sam_h) $(htslib_khash_h) samtools.h $(sam_opts_h)
bam_fastq.o: bam_fastq.c config.h $(htslib_sam_h) $(htslib_klist_h) $(htslib_kstring_h) $(htslib_bgzf_h) $(htslib_thread_pool_h) $(samtools_h) $(sam_opts_h)
bam_index.o: bam_index.c config.h $(htslib_hts_h) $(htslib_sam_h) $(htslib_khash_h) $(samtools_h) $(sam_opts_h)
bam_lpileup.o: bam_lpileup.c config.h $(bam_plbuf_h) $(bam_lpileup_h) $(htslib_ksort_h)
bam_mate.o: bam_mate.c config.h $(htslib_thread_pool_h) $(sam_opts_h) $(htslib_kstring_h) $(htslib_sam_h) samtools.h
bam_md.o: bam_md.c config.h $(htslib_faidx_h) $(htslib_sam_h) $(htslib_kstring_h) $(htslib_thread_pool_h) $(sam_opts_h) samtools.h
bam_mate.o: bam_mate.c config.h $(htslib_thread_pool_h) $(sam_opts_h) $(htslib_kstring_h) $(htslib_sam_h) $(samtools_h)
bam_md.o: bam_md.c config.h $(htslib_faidx_h) $(htslib_sam_h) $(htslib_kstring_h) $(htslib_thread_pool_h) $(sam_opts_h) $(samtools_h)
bam_plbuf.o: bam_plbuf.c config.h $(htslib_hts_h) $(htslib_sam_h) $(bam_plbuf_h)
bam_plcmd.o: bam_plcmd.c config.h $(htslib_sam_h) $(htslib_faidx_h) $(htslib_kstring_h) $(htslib_khash_str2int_h) sam_header.h samtools.h $(sam_opts_h) $(bam2bcf_h) $(sample_h)
bam_plcmd.o: bam_plcmd.c config.h $(htslib_sam_h) $(htslib_faidx_h) $(htslib_kstring_h) $(htslib_klist_h) $(htslib_khash_str2int_h) $(samtools_h) $(bedidx_h) $(sam_opts_h) $(bam2bcf_h) $(sample_h)
bam_quickcheck.o: bam_quickcheck.c config.h $(htslib_hts_h) $(htslib_sam_h)
bam_reheader.o: bam_reheader.c config.h $(htslib_bgzf_h) $(htslib_sam_h) $(htslib_hfile_h) $(htslib_cram_h) samtools.h
bam_rmdup.o: bam_rmdup.c config.h $(htslib_sam_h) $(sam_opts_h) samtools.h $(bam_h) $(htslib_khash_h)
bam_rmdupse.o: bam_rmdupse.c config.h $(bam_h) $(htslib_sam_h) $(htslib_khash_h) $(htslib_klist_h) samtools.h
bam_sort.o: bam_sort.c config.h $(htslib_ksort_h) $(htslib_hts_os_h) $(htslib_khash_h) $(htslib_klist_h) $(htslib_kstring_h) $(htslib_sam_h) $(sam_opts_h) samtools.h
bam_split.o: bam_split.c config.h $(htslib_sam_h) $(htslib_khash_h) $(htslib_kstring_h) $(htslib_cram_h) $(htslib_thread_pool_h) $(sam_opts_h) samtools.h
bam_stat.o: bam_stat.c config.h $(htslib_sam_h) samtools.h $(sam_opts_h)
bam_tview.o: bam_tview.c config.h $(bam_tview_h) $(htslib_faidx_h) $(htslib_sam_h) $(htslib_bgzf_h) samtools.h $(sam_opts_h)
bam_reheader.o: bam_reheader.c config.h $(htslib_bgzf_h) $(htslib_sam_h) $(htslib_hfile_h) $(htslib_cram_h) $(samtools_h)
bam_rmdup.o: bam_rmdup.c config.h $(htslib_sam_h) $(sam_opts_h) $(samtools_h) $(bam_h) $(htslib_khash_h)
bam_rmdupse.o: bam_rmdupse.c config.h $(bam_h) $(htslib_sam_h) $(htslib_khash_h) $(htslib_klist_h) $(samtools_h)
bam_sort.o: bam_sort.c config.h $(htslib_ksort_h) $(htslib_hts_os_h) $(htslib_khash_h) $(htslib_klist_h) $(htslib_kstring_h) $(htslib_sam_h) $(htslib_hts_endian_h) $(sam_opts_h) $(samtools_h)
bam_split.o: bam_split.c config.h $(htslib_sam_h) $(htslib_khash_h) $(htslib_kstring_h) $(htslib_cram_h) $(htslib_thread_pool_h) $(sam_opts_h) $(samtools_h)
bam_stat.o: bam_stat.c config.h $(htslib_sam_h) $(samtools_h) $(sam_opts_h)
bam_tview.o: bam_tview.c config.h $(bam_tview_h) $(htslib_faidx_h) $(htslib_sam_h) $(htslib_bgzf_h) $(samtools_h) $(sam_opts_h)
bam_tview_curses.o: bam_tview_curses.c config.h $(bam_tview_h)
bam_tview_html.o: bam_tview_html.c config.h $(bam_tview_h)
bam_flags.o: bam_flags.c config.h $(htslib_sam_h)
bamshuf.o: bamshuf.c config.h $(htslib_sam_h) $(htslib_hts_h) $(htslib_ksort_h) samtools.h $(htslib_thread_pool_h) $(sam_opts_h) $(htslib_khash_h)
bamtk.o: bamtk.c config.h $(htslib_hts_h) samtools.h version.h
bedcov.o: bedcov.c config.h $(htslib_kstring_h) $(htslib_sam_h) $(htslib_thread_pool_h) $(sam_opts_h) $(htslib_kseq_h)
bamshuf.o: bamshuf.c config.h $(htslib_sam_h) $(htslib_hts_h) $(htslib_ksort_h) $(samtools_h) $(htslib_thread_pool_h) $(sam_opts_h) $(htslib_khash_h)
bamtk.o: bamtk.c config.h $(htslib_hts_h) $(samtools_h) version.h
bedcov.o: bedcov.c config.h $(htslib_kstring_h) $(htslib_sam_h) $(htslib_thread_pool_h) $(samtools_h) $(sam_opts_h) $(htslib_kseq_h)
bedidx.o: bedidx.c config.h $(bedidx_h) $(htslib_ksort_h) $(htslib_kseq_h) $(htslib_khash_h)
cut_target.o: cut_target.c config.h $(htslib_hts_h) $(htslib_sam_h) $(htslib_faidx_h) samtools.h $(sam_opts_h)
cut_target.o: cut_target.c config.h $(htslib_hts_h) $(htslib_sam_h) $(htslib_faidx_h) $(samtools_h) $(sam_opts_h)
dict.o: dict.c config.h $(htslib_kseq_h) $(htslib_hts_h)
faidx.o: faidx.c config.h $(htslib_faidx_h) $(htslib_hts_h) $(htslib_hfile_h) $(htslib_kstring_h) samtools.h
padding.o: padding.c config.h $(htslib_kstring_h) $(htslib_sam_h) $(htslib_faidx_h) sam_header.h $(sam_opts_h) samtools.h
phase.o: phase.c config.h $(htslib_hts_h) $(htslib_sam_h) $(htslib_kstring_h) $(sam_opts_h) samtools.h $(htslib_hts_os_h) $(htslib_kseq_h) $(htslib_khash_h) $(htslib_ksort_h)
faidx.o: faidx.c config.h $(htslib_faidx_h) $(htslib_hts_h) $(htslib_hfile_h) $(htslib_kstring_h) $(samtools_h)
padding.o: padding.c config.h $(htslib_kstring_h) $(htslib_sam_h) $(htslib_faidx_h) $(sam_opts_h) $(samtools_h)
phase.o: phase.c config.h $(htslib_hts_h) $(htslib_sam_h) $(htslib_kstring_h) $(sam_opts_h) $(samtools_h) $(htslib_hts_os_h) $(htslib_kseq_h) $(htslib_khash_h) $(htslib_ksort_h)
sam.o: sam.c config.h $(htslib_faidx_h) $(sam_h)
sam_header.o: sam_header.c config.h sam_header.h $(htslib_khash_h)
sam_opts.o: sam_opts.c config.h $(sam_opts_h)
sam_utils.o: sam_utils.c config.h samtools.h
sam_view.o: sam_view.c config.h $(htslib_sam_h) $(htslib_faidx_h) $(htslib_kstring_h) $(htslib_khash_h) $(htslib_klist_h) $(htslib_thread_pool_h) $(htslib_bgzf_h) samtools.h $(sam_opts_h) $(bedidx_h)
sam_utils.o: sam_utils.c config.h $(samtools_h)
sam_view.o: sam_view.c config.h $(htslib_sam_h) $(htslib_faidx_h) $(htslib_khash_h) $(htslib_thread_pool_h) $(samtools_h) $(sam_opts_h) $(bedidx_h)
sample.o: sample.c config.h $(sample_h) $(htslib_khash_h)
stats_isize.o: stats_isize.c config.h $(stats_isize_h) $(htslib_khash_h)
stats.o: stats.c config.h $(htslib_faidx_h) $(htslib_sam_h) $(htslib_hts_h) sam_header.h $(htslib_khash_str2int_h) samtools.h $(htslib_khash_h) $(htslib_kstring_h) $(stats_isize_h) $(sam_opts_h) $(bedidx_h)
bam_markdup.o: bam_markdup.c config.h $(htslib_thread_pool_h) $(htslib_sam_h) $(sam_opts_h) samtools.h $(htslib_khash_h) $(htslib_klist_h) $(htslib_kstring_h) $(tmp_file_h)
stats.o: stats.c config.h $(htslib_faidx_h) $(htslib_sam_h) $(htslib_hts_h) $(htslib_hts_defs_h) $(htslib_khash_str2int_h) $(samtools_h) $(htslib_khash_h) $(htslib_kstring_h) $(stats_isize_h) $(sam_opts_h) $(bedidx_h)
bam_markdup.o: bam_markdup.c config.h $(htslib_thread_pool_h) $(htslib_sam_h) $(sam_opts_h) $(samtools_h) $(htslib_khash_h) $(htslib_klist_h) $(htslib_kstring_h) $(tmp_file_h)
tmp_file.o: tmp_file.c config.h $(tmp_file_h) $(htslib_sam_h)
# Maintainer source code checks
# - copyright boilerplate presence
# - tab and trailing space detection
maintainer-check:
test/maintainer/check_copyright.pl .
test/maintainer/check_spaces.pl .
# test programs
......@@ -217,26 +223,26 @@ tmp_file.o: tmp_file.c config.h $(tmp_file_h) $(htslib_sam_h)
# If using MSYS, avoid poor shell expansion via:
# MSYS2_ARG_CONV_EXCL="*" make check
check test: samtools $(BGZIP) $(TEST_PROGRAMS)
test/split/test_count_rg
test/split/test_expand_format_string
test/split/test_filter_header_rg
test/split/test_parse_args
REF_PATH=: test/test.pl --exec bgzip=$(BGZIP) $${TEST_OPTS:-}
test/merge/test_bam_translate test/merge/test_bam_translate.tmp
test/merge/test_rtrans_build
test/merge/test_trans_tbl_init
cd test/mpileup && ./regression.sh mpileup.reg
cd test/mpileup && ./regression.sh depth.reg
test/split/test_count_rg
test/split/test_expand_format_string
test/split/test_filter_header_rg
test/split/test_parse_args
test/merge/test_bam_translate: test/merge/test_bam_translate.o test/test.o libst.a $(HTSLIB)
$(CC) $(ALL_LDFLAGS) -o $@ test/merge/test_bam_translate.o test/test.o libst.a $(HTSLIB_LIB) $(ALL_LIBS) -lpthread
test/merge/test_rtrans_build: test/merge/test_rtrans_build.o libst.a $(HTSLIB)
$(CC) $(ALL_LDFLAGS) -o $@ test/merge/test_rtrans_build.o libst.a $(HTSLIB_LIB) $(ALL_LIBS) -lpthread
test/merge/test_rtrans_build: test/merge/test_rtrans_build.o test/test.o libst.a $(HTSLIB)
$(CC) $(ALL_LDFLAGS) -o $@ test/merge/test_rtrans_build.o test/test.o libst.a $(HTSLIB_LIB) $(ALL_LIBS) -lpthread
test/merge/test_trans_tbl_init: test/merge/test_trans_tbl_init.o libst.a $(HTSLIB)
$(CC) $(ALL_LDFLAGS) -o $@ test/merge/test_trans_tbl_init.o libst.a $(HTSLIB_LIB) $(ALL_LIBS) -lpthread
test/merge/test_trans_tbl_init: test/merge/test_trans_tbl_init.o test/test.o libst.a $(HTSLIB)
$(CC) $(ALL_LDFLAGS) -o $@ test/merge/test_trans_tbl_init.o test/test.o libst.a $(HTSLIB_LIB) $(ALL_LIBS) -lpthread
test/split/test_count_rg: test/split/test_count_rg.o test/test.o libst.a $(HTSLIB)
$(CC) $(ALL_LDFLAGS) -o $@ test/split/test_count_rg.o test/test.o libst.a $(HTSLIB_LIB) $(ALL_LIBS) -lpthread
......@@ -260,7 +266,7 @@ test/merge/test_rtrans_build.o: test/merge/test_rtrans_build.c config.h bam_sort
test/merge/test_trans_tbl_init.o: test/merge/test_trans_tbl_init.c config.h bam_sort.o
test/split/test_count_rg.o: test/split/test_count_rg.c config.h bam_split.o $(test_test_h)
test/split/test_expand_format_string.o: test/split/test_expand_format_string.c config.h bam_split.o $(test_test_h)
test/split/test_filter_header_rg.o: test/split/test_filter_header_rg.c config.h bam_split.o $(test_test_h)
test/split/test_filter_header_rg.o: test/split/test_filter_header_rg.c config.h $(test_test_h) $(samtools_h) $(htslib_kstring_h)
test/split/test_parse_args.o: test/split/test_parse_args.c config.h bam_split.o $(test_test_h)
test/test.o: test/test.c config.h $(htslib_sam_h) $(test_test_h)
test/vcf-miniview.o: test/vcf-miniview.c config.h $(htslib_vcf_h)
......@@ -303,7 +309,7 @@ install: $(PROGRAMS) $(MISC_PROGRAMS)
$(INSTALL_PROGRAM) $(PROGRAMS) $(DESTDIR)$(bindir)
$(INSTALL_PROGRAM) $(MISC_PROGRAMS) $(DESTDIR)$(misc_bindir)
$(INSTALL_SCRIPT) $(MISC_SCRIPTS) $(DESTDIR)$(misc_bindir)
$(INSTALL_MAN) samtools.1 misc/wgsim.1 $(DESTDIR)$(man1dir)
$(INSTALL_MAN) doc/samtools*.1 misc/wgsim.1 $(DESTDIR)$(man1dir)
testclean:
......
Makefile.mingw deleted 100644 → 0
View file @ e3705033
CC= gcc.exe
AR= ar.exe
CFLAGS= -g -Wall -O2
DFLAGS= -D_USE_KNETFILE -D_CURSES_LIB=2
KNETFILE_O= knetfile.o
LOBJS= bgzf.o kstring.o bam_aux.o bam.o bam_import.o sam.o bam_index.o \
bam_pileup.o bam_lpileup.o bam_md.o razf.o faidx.o \
$(KNETFILE_O) bam_sort.o sam_header.o bam_reheader.o kprobaln.o bedidx.o
AOBJS= bam_tview.o bam_plcmd.o sam_view.o \
bam_rmdup.o bam_rmdupse.o bam_mate.o bam_stat.o bam_color.o \
bamtk.o kaln.o bam2bcf.o bam2bcf_indel.o errmod.o sample.o \
cut_target.o phase.o bam_cat.o bam2depth.o
BCFOBJS= bcftools/bcf.o bcftools/fet.o bcftools/bcf2qcall.o bcftools/bcfutils.o \
bcftools/call1.o bcftools/index.o bcftools/kfunc.o bcftools/em.o \
bcftools/kmin.o bcftools/prob1.o bcftools/vcf.o bcftools/mut.o
PROG= samtools.exe bcftools.exe
INCLUDES= -I. -Iwin32
SUBDIRS= .
LIBPATH=
.SUFFIXES:.c .o
.c.o:
$(CC) -c $(CFLAGS) $(DFLAGS) $(INCLUDES) $< -o $@
all:$(PROG)
.PHONY:all lib clean cleanlocal
.PHONY:all-recur lib-recur clean-recur cleanlocal-recur install-recur
lib:libbam.a
libbam.a:$(LOBJS)
$(AR) -cru $@ $(LOBJS)
samtools.exe:$(AOBJS) libbam.a $(BCFOBJS)
$(CC) $(CFLAGS) -o $@ $(AOBJS) $(BCFOBJS) $(LIBPATH) -lm -L. -lbam -Lwin32 -lz -lcurses -lws2_32
bcftools.exe:$(BCFOBJS) bcftools/main.o kstring.o bgzf.o knetfile.o bedidx.o
$(CC) $(CFLAGS) -o $@ $(BCFOBJS) bcftools/main.o kstring.o bgzf.o knetfile.o bedidx.o -lm -Lwin32 -lz -lws2_32
razip.o:razf.h
bam.o:bam.h razf.h bam_endian.h kstring.h sam_header.h
sam.o:sam.h bam.h
bam_import.o:bam.h kseq.h khash.h razf.h
bam_pileup.o:bam.h razf.h ksort.h
bam_plcmd.o:bam.h faidx.h bcftools/bcf.h bam2bcf.h
bam_index.o:bam.h khash.h ksort.h razf.h bam_endian.h
bam_lpileup.o:bam.h ksort.h
bam_tview.o:bam.h faidx.h
bam_sort.o:bam.h ksort.h razf.h
bam_md.o:bam.h faidx.h
sam_header.o:sam_header.h khash.h
bcf.o:bcftools/bcf.h
bam2bcf.o:bam2bcf.h errmod.h bcftools/bcf.h
bam2bcf_indel.o:bam2bcf.h
errmod.o:errmod.h
faidx.o:faidx.h razf.h khash.h
faidx_main.o:faidx.h razf.h
clean:
rm -fr gmon.out *.o a.out *.exe *.dSYM razip bgzip $(PROG) *~ *.a *.so.* *.so *.dylib
NEWS
View file @ 18118746
Release 1.10 (6th December 2019)
--------------------------------
Changes affecting the whole of samtools, or multiple sub-commands:
* Samtools now uses the new HTSlib header API. As this adds more checks
for invalid headers, it is possible that some illegal files will now
be rejected when they would have been allowed by earlier versions. (#998)
Examples of problems that will now be rejected include @SQ lines with
no SN: tag, and @RG or @PG lines with no ID: tag.
* samtools sub-commands will now add '@PG' header lines to output sam/bam/cram
files. To disable this, use the '--no-PG' option. (#1087; #1097)
* samtools now supports alignment records with reference positions greater
than 2 gigabases. This allows samtools to process alignments for
species which have large chromosomes, like axolotl and lungfish. Note that
due to file format limitations, data with large reference positions
must use the SAM format. (#1107; #1117)
* Improved the efficiency of reading and writing SAM format data by 2 fold
(single thread). This is further improved by the ability to use multiple
threads, as previously done with BAM and CRAM.
* samtools can now write BGZF-compressed SAM format. To enable this,
either save files with a '.sam.gz' suffix, or use '--output-fmt sam.gz'.
* samtools can now index BGZF-compressed SAM files.
* The region parsing code has been improved to handle colons in reference
names. Strings can be disambiguated by the use of braces, so for
example when reference sequences called "chr1" and "chr1:100-200"
are both present, the regions "{chr1}:100-200" and "{chr1:100-200}"
unambiguously indicate which reference is being used. (#864)
* samtools flags, flagstats, idxstats and stats now have aliases
flag, flagstat, idxstat and stat. (#934)
* A new global '--write-index' option has been added. This allows output
sam.gz/bam/cram files to be indexed while they are being written out.
This should work with addreplacerg, depad, markdup, merge, sort, split,
and view. (#1062)
* A global '--verbosity' option has been added to enable/disable
debugging output. (#1124, thanks to John Marshall)
* It is now possible to have data and index files stored in different
locations. There are two ways to tell samtools where to find the
index:
1. Samtools bedcov, depth, merge, mpileup, stats, tview, and view
accept a new option (-X). When this is used, each input sam/bam/cram
listed on the command line should have a corresponding index file.
Note that all the data files should be listed first, followed by all
the index files. (#978, thanks to Mingfei Shao)
2. A delimiter '##idx##' can be appended to the data file name followed
by the index file name. This can be used both for input files and
outputs when indexing on-the-fly.
* HTSlib (and therefore SAMtools) now uses version 4 signatures by default
for its s3:// plug-in. It can also write to S3 buckets, as long as
version 4 signatures are in use. See HTSlib's NEWS file and
htslib-s3-plugin manual page for more information.
* HTSlib (and therefore SAMtools) no longer considers a zero-length
file to be a valid SAM file. This has been changed so that pipelines such
as `somecmd | samtools ...` with `somecmd` aborting before outputting
anything will now propagate the error to the second command.
* The samtools manual page has been split up into one for each
sub-command. The main samtools.1 manual page now lists the sub-commands
and describes the common global options. (#894)
* The meaning of decode_md, store_md and store_nm in the fmt-option section
of the samtools.1 man page has been clarified. (#898, thanks to Evan Benn)
* Fixed numerous memory leaks. (#892)
* Fixed incorrect macro definition on Windows. (#950)
* bedcov, phase, misc/ace2sam and misc/wgsim now check for failure to open
files. (#1013, thanks to Julie Blommaert and John Marshall)
Changes affecting specific sub-commands:
* A new "coverage" sub-command has been added. This prints a tabular format
of the average coverage and percent coverage for each reference sequence,
as well as number of aligned reads, average mapping quality and base
quality. It can also (with the '-m' option) plot a histogram of
coverage across the genome. (#992, thanks to Florian Breitwieser)
* samtools calmd:
- Reference bases in MD: tags are now converted to upper case. (#981, #988)
* samtools depth:
- Add new options to write a header to the output (-H) and to direct
the output to a file (-o). (#937, thanks to Pierre Lindenbaum)
- New options '-g' and '-G' can be used to filter reads. (#953)
- Fix memory leak when failing to set CRAM options. (#985, thanks
to Florian Breitwieser)
- Fix bug when using region filters where the '-a' option did not
work for regions with no coverage. (#1113; #1112 reported by
Paweł Sztromwasser)
* samtools fasta and fastq:
- '-1 FILE -2 FILE' with the same filename now works properly. (#1042)
- '-o FILE' is added as a synonym for '-1 FILE -2 FILE'. (#1042)
- The '-F' option now defaults to 0x900 (SECONDARY,SUPPLEMENTARY).
Previously secondary and supplementary records were filtered internally
in a way that could not be turned off. (#1042; #939 reported
by @finswimmer)
- Allow reading from a pipe without an explicit '-' on the command line.
(#1042; #874 reported by John Marshall)
- Turn on multi-threading for bgzf compressed output files. (#908)
- Fixed bug where the samtools fastq -i would output incorrect information
in the Casava tags for dual-index reads. It also now prints the tags
for dual indices in the same way as bcl2fastq, using a '+' sign between
the two parts of the index. (#1059; #1047 reported by Denis Loginov)
* samtools flagstat:
- Samtools flagstat can now optionally write its output in JSON format or
as a tab-separated values file. (#1106, thanks to Vivek Rai).
* samtools markdup:
- It can optionally tag optical duplicates (reads following Illumina
naming conventions only). The is enabled with the '-d' option,
which sets the distance for duplicates to be considered as optical.
(#1091; #1103; #1121; #1128; #1134)
- The report stats (-s) option now outputs counts for optical and
non-primary (supplementary / secondary) duplicates. It also reports
the Picard "estimate library size" statistic. A new '-f' option
can be used to save the statistics in a given file. (#1091)
- The rules for calling duplicates can be changed using the new --mode
option. This mainly changes the position associated with each read in
a pair. '--mode t' (the default) is the existing behaviour where the
position used is that of the outermost template base associated with the
read. Alternatively '--mode s' always uses the first unclipped sequence
base. In practice, this only makes a difference for read pairs where the
two reads are aligned in the same direction. (#1091)
- A new '-c' option can be used to clear any existing duplicate tags.
(#1091)
- A new '--include-fails' option makes markdup include QC-failed reads.
(#1091)
- Fixed buffer overflow in temporary file writer when writing a mixture
of long and short alignment records. (#911; #909)
* samtools mpileup:
- mpileup can now process alignments including CIGAR P (pad) operators
correctly. They will now also produce the correct output for alignments
where insertions are immediately followed by deletions, or deletions by
insertions. Note that due to limitations in HTSlib, they are still
unable to output sequences that have been inserted before the first
aligned base of a read. (#847; #842 reported by Tiffany Delhomme.
See also htslib issue #59 and pull request #699).
- In samtools mpileup, a deletion or pad on the reverse strand is now
marked with a different character ('#') than the one used on a forward
strand ('*'), if the '--reverse-del' option is used. (#1070)
- New option '--output-extra' can be used to add columns for user
selected alignment fields or aux tags. (#1073)
- Fixed double-counting of overlapping bases in alignment records with
deletions or reference skips longer than twice the insert size.
(#989; #987 reported by @dariomel)
- Improved manual page with documentation about what each output column
means. (#1055, thanks to John Marshall)
* samtools quickcheck:
- Add unmapped (-u) option, which disables the check for @SQ lines in
the header. (#920, thanks to Shane McCarthy)
* samtools reheader:
- A new option '-c' allows the input header to be passed to a given
command. Samtools then takes the output of this command and uses it
as the replacement header. (#1007)
- Make it clear in help message that reheader --in-place only works on
CRAM files. (#921, thanks to Julian Gehring)
- Refuse to in-place reheader BAM files, instead of unexpectedly writing
a BAM file to stdout. (#935)
* samtools split:
- In samtools split, the '-u' option no longer accepts an extra file name
from which a replacement header was read. The two file names were
separated using a colon, which caused problems on Windows and prevented
the use of URLs. A new '-h' option has been added to allow the replacement
header file to be specified in its own option. (#961)
- Fixed bug where samtools split would crash if it read a SAM header that
contained an @RG line with no ID tag. (#954, reported by @blue-bird1)
* samtools stats:
- stats will now compute base compositions for BC, CR, OX and RX tags,
and quality histograms for QT, CY, BZ and QX tags. (#904)
- New stats FTC and LTC showing total number of nucleotides for first and
last fragments. (#946)
- The rules for classifying reads as "first" or "last" fragment have been
tightened up. (#949)
- Fixed bug where stats could over-estimate coverage when using the
target-regions option or when a region was specified on the command-line.
(#1027; #1025, reported by Miguel Machado; #1029, reported by Jody Phelan).
- Fixed error in stats GCD percentile depth calculation when the depth to be
reported fell between two bins. It would report the depth entirely from
the lower bin instead of taking a weighted average of the two. (#1048)
- Better catching and reporting of out of memory conditions. (#984;
#982, reported by Jukka Matilainen)
- Improved manual page. (#927)
* samtools tview:
- tview can now display alignments including CIGAR P operators, D followed
by I and I followed by D correctly. See mpileup above for more
information. (#847; #734, reported by Ryan Lorig-Roach)
- The "go to position" text entry box has been made wider. (#968, thanks
to John Marshall)
- Fixed samtools tview -s option which was not filtering reads correctly.
It now only shows reads from the requested sample or read group. (#1089)
* samtools view:
- New options '-d' and '-D' to only output alignments which have a tag
with a given type and value. (#1001, thanks to Gert Hulselmans)
* misc/plot-bamstats script:
- Fixed merge (-m) option. (#923, #924 both thanks to Marcus D Sherman)
- Made the quality heatmap work with gnuplot version 5.2.7 and later.
(#1068; #1065 reported by Martin Mokrejš)
- Fixed --do-ref-stats bug where fasta header lines would be counted
as part of the sequence when the --targets option was used. (#1120,
thanks to Neil Goodgame)
* Removed the misc/varfilter.py Python script, as it takes consensus-pileup
as input, which was removed from samtools in release 0.1.17 in 2011. (#1125)
Release 1.9 (18th July 2018)
----------------------------
......@@ -91,8 +364,8 @@ Release 1.9 (18th July 2018)
- Fixed compilation problems when using C_INCLUDE_PATH. (#870; #817 reported
by Robert Boissy)
- Fixed --version when built from a Git repository. (#844, thanks to John
Marshall)
- Fixed --version when built from a Git repository. (#844, thanks to
John Marshall)
- Use noenhanced mode for title in plot-bamstats. Prevents unwanted
interpretation of characters like underscore in gnuplot version 5. (#829,
......
README
View file @ 18118746
......@@ -9,7 +9,7 @@ Building samtools
The typical simple case of building Samtools using the HTSlib bundled within
this Samtools release tarball is done as follows:
cd .../samtools-1.9 # Within the unpacked release directory
cd .../samtools-1.10 # Within the unpacked release directory
./configure
make
......@@ -21,7 +21,7 @@ install samtools etc properly into a directory of your choosing. Building for
installation using the HTSlib bundled within this Samtools release tarball,
and building the various HTSlib utilities such as bgzip is done as follows:
cd .../samtools-1.9 # Within the unpacked release directory
cd .../samtools-1.10 # Within the unpacked release directory
./configure --prefix=/path/to/location
make all all-htslib
make install install-htslib
......@@ -48,7 +48,7 @@ There are two advantages to this:
To build with plug-ins, you need to use the --enable-plugins configure option
as follows:
cd .../samtools-1.9 # Within the unpacked release directory
cd .../samtools-1.10 # Within the unpacked release directory
./configure --enable-plugins --prefix=/path/to/location
make all all-htslib
make install install-htslib
......@@ -66,8 +66,8 @@ Setting --with-plugin-path is useful if you want to run directly from
the source distribution instead of installing the package. In that case
you can use:
cd .../samtools-1.9 # Within the unpacked release directory
./configure --enable-plugins --with-plugin-path=$PWD/htslib-1.9
cd .../samtools-1.10 # Within the unpacked release directory
./configure --enable-plugins --with-plugin-path=$PWD/htslib-1.10
make all all-htslib
It is possible to override the built-in search path using the HTS_PATH
......
bam.c
View file @ 18118746
/* bam.c -- BAM format.
Copyright (C) 2008-2013, 2015 Genome Research Ltd.
Copyright (C) 2008-2013, 2015, 2019 Genome Research Ltd.
Portions copyright (C) 2009-2012 Broad Institute.
Author: Heng Li <lh3@sanger.ac.uk>
......@@ -30,7 +30,6 @@ DEALINGS IN THE SOFTWARE. */
#include <errno.h>
#include "bam.h"
#include "htslib/kstring.h"
#include "sam_header.h"
char *bam_format1(const bam_header_t *header, const bam1_t *b)
{
......@@ -59,7 +58,7 @@ int bam_validate1(const bam_header_t *header, const bam1_t *b)
char *s;
if (b->core.tid < -1 || b->core.mtid < -1) return 0;
if (header && (b->core.tid >= header->n_targets || b->core.mtid >= header->n_targets)) return 0;
if (header && (b->core.tid >= sam_hdr_nref(header) || b->core.mtid >= sam_hdr_nref(header))) return 0;
if (b->data_len < b->core.l_qname) return 0;
s = memchr(bam1_qname(b), '\0', b->core.l_qname);
......@@ -77,9 +76,8 @@ int bam_validate1(const bam_header_t *header, const bam1_t *b)
// FIXME: we should also check the LB tag associated with each alignment
const char *bam_get_library(bam_header_t *h, const bam1_t *b)
{
// Slow and inefficient. Rewrite once we get a proper header API.
const char *rg;
char *cp = h->text;
kstring_t lib = { 0, 0, NULL };
rg = (char *)bam_aux_get(b, "RG");
if (!rg)
......@@ -87,52 +85,20 @@ const char *bam_get_library(bam_header_t *h, const bam1_t *b)
else
rg++;
// Header is guaranteed to be nul terminated, so this is valid.
while (*cp) {
char *ID, *LB;
char last = '\t';
// Find a @RG line
if (strncmp(cp, "@RG", 3) != 0) {
while (*cp && *cp != '\n') cp++; // skip line
if (*cp) cp++;
continue;
}
// Find ID: and LB: keys
cp += 4;
ID = LB = NULL;
while (*cp && *cp != '\n') {
if (last == '\t') {
if (strncmp(cp, "LB:", 3) == 0)
LB = cp+3;
else if (strncmp(cp, "ID:", 3) == 0)
ID = cp+3;
}
last = *cp++;
}
if (sam_hdr_find_tag_id(h, "RG", "ID", rg, "LB", &lib) < 0)
return NULL;
if (!ID || !LB)
continue;
static char LB_text[1024];
int len = lib.l < sizeof(LB_text) - 1 ? lib.l : sizeof(LB_text) - 1;
// Check it's the correct ID
if (strncmp(rg, ID, strlen(rg)) != 0 || ID[strlen(rg)] != '\t')
continue;
memcpy(LB_text, lib.s, len);
LB_text[len] = 0;
// Valid until next query
static char LB_text[1024];
for (cp = LB; *cp && *cp != '\t' && *cp != '\n'; cp++)
;
strncpy(LB_text, LB, MIN(cp-LB, 1023));
LB_text[MIN(cp-LB, 1023)] = 0;
free(lib.s);
// Return it; valid until the next query.
return LB_text;
}
return NULL;
}
int bam_fetch(bamFile fp, const bam_index_t *idx, int tid, int beg, int end, void *data, bam_fetch_f func)
{
int ret;
......
bam.h
View file @ 18118746
/* bam.h -- BAM API.
Copyright (C) 2008-2014 Genome Research Ltd.
Copyright (C) 2008-2014, 2019 Genome Research Ltd.
Portions copyright (C) 2010-2012 Broad Institute.
Author: Heng Li <lh3@sanger.ac.uk>
......@@ -38,7 +38,7 @@ DEALINGS IN THE SOFTWARE. */
@copyright Genome Research Ltd.
*/
#define BAM_VERSION "1.9"
#define BAM_VERSION "1.10"
#include <stdint.h>
#include <stdlib.h>
......@@ -223,16 +223,6 @@ extern "C" {
*/
// int sam_read1(tamFile fp, bam_header_t *header, bam1_t *b);
/*!
@abstract Read header information from a TAB-delimited list file.
@param fn_list file name for the list
@return a pointer to the header structure
@discussion Each line in this file consists of chromosome name and
the length of chromosome.
*/
bam_header_t *sam_header_read2(const char *fn_list);
/*!
@abstract Read header from a SAM file (if present)
@param fp SAM file handler
......@@ -252,13 +242,13 @@ extern "C" {
@abstract Initialize a header structure.
@return the pointer to the header structure
*/
static inline bam_header_t *bam_header_init(void) { return bam_hdr_init(); }
static inline bam_header_t *bam_header_init(void) { return sam_hdr_init(); }
/*!
@abstract Destroy a header structure.
@param header pointer to the header
*/
static inline void bam_header_destroy(bam_header_t *header) { bam_hdr_destroy(header); }
static inline void bam_header_destroy(bam_header_t *header) { sam_hdr_destroy(header); }
/*!
@abstract Read a header structure from BAM.
......@@ -277,7 +267,7 @@ extern "C" {
@param header pointer to the header structure
@return always 0 currently
*/
static inline int bam_header_write(bamFile fp, const bam_header_t *header) { return bam_hdr_write(fp, header); }
static inline int bam_header_write(bamFile fp, bam_header_t *header) { return bam_hdr_write(fp, header); }
/*!
@abstract Read an alignment from BAM.
......
bam2bcf.c
View file @ 18118746
/* bam2bcf.c -- variant calling.
Copyright (C) 2010-2012 Broad Institute.
Copyright (C) 2012-2014 Genome Research Ltd.
Copyright (C) 2012-2015 Genome Research Ltd.
Author: Heng Li <lh3@sanger.ac.uk>
......
bam2bcf.h
View file @ 18118746
/* bam2bcf.h -- variant calling.
Copyright (C) 2010-2012 Broad Institute.
Copyright (C) 2012-2014 Genome Research Ltd.
Copyright (C) 2012-2014, 2019 Genome Research Ltd.
Author: Heng Li <lh3@sanger.ac.uk>
......@@ -99,7 +99,8 @@ typedef struct {
} bcf_callret1_t;
typedef struct {
int tid, pos;
int tid;
hts_pos_t pos;
bcf_hdr_t *bcf_hdr;
int a[5]; // alleles: ref, alt, alt2, alt3
float qsum[5]; // for the QS tag
......@@ -128,7 +129,7 @@ extern "C" {
int bcf_call_combine(int n, const bcf_callret1_t *calls, bcf_callaux_t *bca, int ref_base /*4-bit*/, bcf_call_t *call);
int bcf_call2bcf(bcf_call_t *bc, bcf1_t *b, bcf_callret1_t *bcr, int fmt_flag,
const bcf_callaux_t *bca, const char *ref);
int bcf_call_gap_prep(int n, int *n_plp, bam_pileup1_t **plp, int pos, bcf_callaux_t *bca, const char *ref,
int bcf_call_gap_prep(int n, int *n_plp, bam_pileup1_t **plp, hts_pos_t pos, bcf_callaux_t *bca, const char *ref,
const void *rghash);
void bcf_callaux_clean(bcf_callaux_t *bca, bcf_call_t *call);
......
bam2bcf_indel.c
View file @ 18118746
/* bam2bcf_indel.c -- indel caller.
Copyright (C) 2010, 2011 Broad Institute.
Copyright (C) 2012-2014 Genome Research Ltd.
Copyright (C) 2012-2014, 2019 Genome Research Ltd.
Author: Heng Li <lh3@sanger.ac.uk>
......@@ -87,9 +87,10 @@ void bcf_call_del_rghash(void *_hash)
kh_destroy(rg, hash);
}
static int tpos2qpos(const bam1_core_t *c, const uint32_t *cigar, int32_t tpos, int is_left, int32_t *_tpos)
static int tpos2qpos(const bam1_core_t *c, const uint32_t *cigar, hts_pos_t tpos, hts_pos_t is_left, hts_pos_t *_tpos)
{
int k, x = c->pos, y = 0, last_y = 0;
int k, y = 0, last_y = 0;
hts_pos_t x = c->pos;
*_tpos = c->pos;
for (k = 0; k < c->n_cigar; ++k) {
int op = cigar[k] & BAM_CIGAR_MASK;
......@@ -124,9 +125,10 @@ static inline int est_seqQ(const bcf_callaux_t *bca, int l, int l_run)
return q < qh? q : qh;
}
static inline int est_indelreg(int pos, const char *ref, int l, char *ins4)
static inline int est_indelreg(hts_pos_t pos, const char *ref, int l, char *ins4)
{
int i, j, max = 0, max_i = pos, score = 0;
int j, max = 0, score = 0;
hts_pos_t i, max_i = pos;
l = abs(l);
for (i = pos + 1, j = 0; ref[i]; ++i, ++j) {
if (ins4) score += (toupper(ref[i]) != "ACGTN"[(int)ins4[j%l]])? -10 : 1;
......@@ -146,11 +148,12 @@ static inline int est_indelreg(int pos, const char *ref, int l, char *ins4)
- 8: estimated sequence quality .. (aux>>8)&0xff
- 8: indel quality .. aux&0xff
*/
int bcf_call_gap_prep(int n, int *n_plp, bam_pileup1_t **plp, int pos, bcf_callaux_t *bca, const char *ref,
int bcf_call_gap_prep(int n, int *n_plp, bam_pileup1_t **plp, hts_pos_t pos, bcf_callaux_t *bca, const char *ref,
const void *rghash)
{
int i, s, j, k, t, n_types, *types, max_rd_len, left, right, max_ins, *score1, *score2, max_ref2;
int s, k, t, n_types, *types, max_rd_len, max_ins, *score1, *score2, max_ref2;
int N, K, l_run, ref_type, n_alt;
hts_pos_t i, j, left, right;
char *inscns = 0, *ref2, *query, **ref_sample;
khash_t(rg) *hash = (khash_t(rg)*)rghash;
if (ref == 0 || bca == 0) return -1;
......@@ -225,7 +228,7 @@ int bcf_call_gap_prep(int n, int *n_plp, bam_pileup1_t **plp, int pos, bcf_calla
free(aux);
// TODO revisit how/whether to control printing this warning
if (hts_verbose >= 2)
fprintf(stderr, "[%s] excessive INDEL alleles at position %d. Skip the position.\n", __func__, pos + 1);
fprintf(stderr, "[%s] excessive INDEL alleles at position %"PRIhts_pos". Skip the position.\n", __func__, pos + 1);
return -1;
}
types = (int*)calloc(n_types, sizeof(int));
......@@ -274,7 +277,7 @@ int bcf_call_gap_prep(int n, int *n_plp, bam_pileup1_t **plp, int pos, bcf_calla
bam1_t *b = p->b;
uint32_t *cigar = bam_get_cigar(b);
uint8_t *seq = bam_get_seq(b);
int x = b->core.pos, y = 0;
hts_pos_t x = b->core.pos, y = 0;
for (k = 0; k < b->core.n_cigar; ++k) {
int op = cigar[k]&0xf;
int j, l = cigar[k]>>4;
......@@ -382,7 +385,8 @@ int bcf_call_gap_prep(int n, int *n_plp, bam_pileup1_t **plp, int pos, bcf_calla
// align each read to ref2
for (i = 0; i < n_plp[s]; ++i, ++K) {
bam_pileup1_t *p = plp[s] + i;
int qbeg, qend, tbeg, tend, sc, kk;
int qbeg, qend, sc, kk;
hts_pos_t tbeg, tend;
uint8_t *seq = bam_get_seq(p->b);
uint32_t *cigar = bam_get_cigar(p->b);
if (p->b->core.flag&4) continue; // unmapped reads
......
bam2depth.c
View file @ 18118746
/* bam2depth.c -- depth subcommand.
Copyright (C) 2011, 2012 Broad Institute.
Copyright (C) 2012-2014 Genome Research Ltd.
Copyright (C) 2012-2016, 2018, 2019 Genome Research Ltd.
Author: Heng Li <lh3@sanger.ac.uk>
......@@ -39,20 +39,19 @@ DEALINGS IN THE SOFTWARE. */
#include <unistd.h>
#include "htslib/sam.h"
#include "samtools.h"
#include "bedidx.h"
#include "sam_opts.h"
#define BAM_FMAX ((BAM_FSUPPLEMENTARY << 1) - 1)
typedef struct { // auxiliary data structure
samFile *fp; // the file handle
bam_hdr_t *hdr; // the file header
sam_hdr_t *hdr; // the file header
hts_itr_t *iter; // NULL if a region not specified
int min_mapQ, min_len; // mapQ filter; length filter
uint32_t flags; // read filtering flags
} aux_t;
void *bed_read(const char *fn); // read a BED or position list file
void bed_destroy(void *_h); // destroy the BED data structure
int bed_overlap(const void *_h, const char *chr, int beg, int end); // test if chr:beg-end overlaps
int bed_query(const void *_h, const char *chr, int pos, int *beg, int *end);
// This function reads a BAM alignment from one BAM file.
static int read_bam(void *data, bam1_t *b) // read level filters better go here to avoid pileup
{
......@@ -62,7 +61,7 @@ static int read_bam(void *data, bam1_t *b) // read level filters better go here
{
ret = aux->iter? sam_itr_next(aux->fp, aux->iter, b) : sam_read1(aux->fp, aux->hdr, b);
if ( ret<0 ) break;
if ( b->core.flag & (BAM_FUNMAP | BAM_FSECONDARY | BAM_FQCFAIL | BAM_FDUP) ) continue;
if ( b->core.flag & aux->flags) continue;
if ( (int)b->core.qual < aux->min_mapQ ) continue;
if ( aux->min_len && bam_cigar2qlen(b->core.n_cigar, bam_get_cigar(b)) < aux->min_len ) continue;
break;
......@@ -79,15 +78,21 @@ static int usage() {
fprintf(stderr, " -a output all positions (including zero depth)\n");
fprintf(stderr, " -a -a (or -aa) output absolutely all positions, including unused ref. sequences\n");
fprintf(stderr, " -b <bed> list of positions or regions\n");
fprintf(stderr, " -X use customized index files\n");
fprintf(stderr, " -f <list> list of input BAM filenames, one per line [null]\n");
fprintf(stderr, " -H print a file header\n");
fprintf(stderr, " -l <int> read length threshold (ignore reads shorter than <int>) [0]\n");
fprintf(stderr, " -d/-m <int> maximum coverage depth [8000]. If 0, depth is set to the maximum\n"
" integer value, effectively removing any depth limit.\n"); // the htslib's default
fprintf(stderr, " -o FILE where to write output to [stdout]\n");
fprintf(stderr, " -q <int> base quality threshold [0]\n");
fprintf(stderr, " -Q <int> mapping quality threshold [0]\n");
fprintf(stderr, " -r <chr:from-to> region\n");
fprintf(stderr, " -g <flags> include reads that have any of the specified flags set [0]\n");
fprintf(stderr, " -G <flags> filter out reads that have any of the specified flags set"
" [UNMAP,SECONDARY,QCFAIL,DUP]\n");
sam_global_opt_help(stderr, "-.--.-");
sam_global_opt_help(stderr, "-.--.--.");
fprintf(stderr, "\n");
fprintf(stderr, "The output is a simple tab-separated table with three columns: reference name,\n");
......@@ -95,21 +100,27 @@ static int usage() {
fprintf(stderr, "omitted by default; see the -a option.\n");
fprintf(stderr, "\n");
return 1;
return EXIT_FAILURE;
}
int main_depth(int argc, char *argv[])
{
int i, n, tid, reg_tid, beg, end, pos, *n_plp, baseQ = 0, mapQ = 0, min_len = 0;
int i, n, tid, reg_tid, *n_plp, baseQ = 0, mapQ = 0, min_len = 0, has_index_file = 0;
hts_pos_t beg, end, pos, last_pos = -1;
int all = 0, status = EXIT_SUCCESS, nfiles, max_depth = -1;
const bam_pileup1_t **plp;
char *reg = 0; // specified region
void *bed = 0; // BED data structure
char *file_list = NULL, **fn = NULL;
bam_hdr_t *h = NULL; // BAM header of the 1st input
sam_hdr_t *h = NULL; // BAM header of the 1st input
aux_t **data;
bam_mplp_t mplp;
int last_pos = -1, last_tid = -1, ret;
int last_tid = -1, ret;
int print_header = 0;
char *output_file = NULL;
FILE *file_out = stdout;
uint32_t flags = (BAM_FUNMAP | BAM_FSECONDARY | BAM_FQCFAIL | BAM_FDUP);
int tflags = 0;
sam_global_args ga = SAM_GLOBAL_ARGS_INIT;
static const struct option lopts[] = {
......@@ -118,19 +129,41 @@ int main_depth(int argc, char *argv[])
};
// parse the command line
while ((n = getopt_long(argc, argv, "r:b:q:Q:l:f:am:d:", lopts, NULL)) >= 0) {
while ((n = getopt_long(argc, argv, "r:b:Xq:Q:l:f:am:d:Ho:g:G:", lopts, NULL)) >= 0) {
switch (n) {
case 'l': min_len = atoi(optarg); break; // minimum query length
case 'r': reg = strdup(optarg); break; // parsing a region requires a BAM header
case 'b':
bed = bed_read(optarg); // BED or position list file can be parsed now
if (!bed) { print_error_errno("depth", "Could not read file \"%s\"", optarg); return 1; }
if (!bed) {
print_error_errno("depth", "Could not read file \"%s\"", optarg);
return EXIT_FAILURE;
}
break;
case 'X': has_index_file = 1; break;
case 'q': baseQ = atoi(optarg); break; // base quality threshold
case 'Q': mapQ = atoi(optarg); break; // mapping quality threshold
case 'f': file_list = optarg; break;
case 'a': all++; break;
case 'd': case 'm': max_depth = atoi(optarg); break; // maximum coverage depth
case 'H': print_header = 1; break;
case 'o': output_file = optarg; break;
case 'g':
tflags = bam_str2flag(optarg);
if (tflags < 0 || tflags > BAM_FMAX) {
print_error_errno("depth", "Flag value \"%s\" is not supported", optarg);
return 1;
}
flags &= ~tflags;
break;
case 'G':
tflags = bam_str2flag(optarg);
if (tflags < 0 || tflags > BAM_FMAX) {
print_error_errno("depth", "Flag value \"%s\" is not supported", optarg);
return 1;
}
flags |= tflags;
break;
default: if (parse_sam_global_opt(n, optarg, lopts, &ga) == 0) break;
/* else fall-through */
case '?': return usage();
......@@ -139,18 +172,40 @@ int main_depth(int argc, char *argv[])
if (optind == argc && !file_list)
return usage();
/* output file provided by user */
if (output_file != NULL && strcmp(output_file,"-")!=0) {
file_out = fopen( output_file, "w" );
if (file_out == NULL) {
print_error_errno("depth", "Cannot open \"%s\" for writing.", output_file);
return EXIT_FAILURE;
}
}
// initialize the auxiliary data structures
if (file_list)
{
if ( read_file_list(file_list,&nfiles,&fn) ) return 1;
if (has_index_file) {
print_error("depth", "The -f option cannot be combined with -X");
return 1;
}
if ( read_file_list(file_list,&nfiles,&fn) ) return EXIT_FAILURE;
n = nfiles;
argv = fn;
optind = 0;
}
else
n = argc - optind; // the number of BAMs on the command line
else if (has_index_file) { // Calculate # of input BAM files
if ((argc - optind) % 2 != 0) {
fprintf(stderr, "Error: Odd number of filenames detected! Each BAM file should have an index file\n");
return 1;
}
n = (argc - optind) / 2;
} else {
n = argc - optind;
}
data = calloc(n, sizeof(aux_t*)); // data[i] for the i-th input
reg_tid = 0; beg = 0; end = INT_MAX; // set the default region
reg_tid = 0; beg = 0; end = HTS_POS_MAX; // set the default region
for (i = 0; i < n; ++i) {
int rf;
data[i] = calloc(1, sizeof(aux_t));
......@@ -163,24 +218,32 @@ int main_depth(int argc, char *argv[])
rf = SAM_FLAG | SAM_RNAME | SAM_POS | SAM_MAPQ | SAM_CIGAR | SAM_SEQ;
if (baseQ) rf |= SAM_QUAL;
if (hts_set_opt(data[i]->fp, CRAM_OPT_REQUIRED_FIELDS, rf)) {
fprintf(stderr, "Failed to set CRAM_OPT_REQUIRED_FIELDS value\n");
return 1;
print_error_errno("depth", "Failed to set CRAM_OPT_REQUIRED_FIELDS value");
status = EXIT_FAILURE;
goto depth_end;
}
if (hts_set_opt(data[i]->fp, CRAM_OPT_DECODE_MD, 0)) {
fprintf(stderr, "Failed to set CRAM_OPT_DECODE_MD value\n");
return 1;
print_error_errno("depth", "Failed to set CRAM_OPT_DECODE_MD value");
status = EXIT_FAILURE;
goto depth_end;
}
data[i]->min_mapQ = mapQ; // set the mapQ filter
data[i]->min_len = min_len; // set the qlen filter
data[i]->hdr = sam_hdr_read(data[i]->fp); // read the BAM header
if (data[i]->hdr == NULL) {
fprintf(stderr, "Couldn't read header for \"%s\"\n",
print_error_errno("depth", "Couldn't read header for \"%s\"",
argv[optind+i]);
status = EXIT_FAILURE;
goto depth_end;
}
if (reg) { // if a region is specified
hts_idx_t *idx = sam_index_load(data[i]->fp, argv[optind+i]); // load the index
hts_idx_t *idx = NULL;
// If index filename has not been specfied, look in BAM folder
if (has_index_file) {
idx = sam_index_load2(data[i]->fp, argv[optind+i], argv[optind+i+n]); // load the index
} else {
idx = sam_index_load(data[i]->fp, argv[optind+i]);
}
if (idx == NULL) {
print_error("depth", "can't load index for \"%s\"", argv[optind+i]);
status = EXIT_FAILURE;
......@@ -194,8 +257,16 @@ int main_depth(int argc, char *argv[])
goto depth_end;
}
}
data[i]->flags = flags;
}
if (print_header) {
fputs("#CHROM\tPOS", file_out);
for (i = 0; i < n; ++i) {
fputc('\t', file_out);
fputs(argv[optind+i], file_out);
}
fputc('\n', file_out);
}
h = data[0]->hdr; // easy access to the header of the 1st BAM
if (reg) {
beg = data[0]->iter->beg; // and to the parsed region coordinates
......@@ -211,21 +282,22 @@ int main_depth(int argc, char *argv[])
bam_mplp_set_maxcnt(mplp,INT_MAX);
n_plp = calloc(n, sizeof(int)); // n_plp[i] is the number of covering reads from the i-th BAM
plp = calloc(n, sizeof(bam_pileup1_t*)); // plp[i] points to the array of covering reads (internal in mplp)
while ((ret=bam_mplp_auto(mplp, &tid, &pos, n_plp, plp)) > 0) { // come to the next covered position
while ((ret=bam_mplp64_auto(mplp, &tid, &pos, n_plp, plp)) > 0) { // come to the next covered position
if (pos < beg || pos >= end) continue; // out of range; skip
if (tid >= h->n_targets) continue; // diff number of @SQ lines per file?
if (tid >= sam_hdr_nref(h)) continue; // diff number of @SQ lines per file?
if (all) {
while (tid > last_tid) {
if (last_tid >= 0 && !reg) {
// Deal with remainder or entirety of last tid.
while (++last_pos < h->target_len[last_tid]) {
while (++last_pos < sam_hdr_tid2len(h, last_tid)) {
// Horribly inefficient, but the bed API is an obfuscated black box.
if (bed && bed_overlap(bed, h->target_name[last_tid], last_pos, last_pos + 1) == 0)
if (bed && bed_overlap(bed, sam_hdr_tid2name(h, last_tid), last_pos, last_pos + 1) == 0)
continue;
fputs(h->target_name[last_tid], stdout); printf("\t%d", last_pos+1);
fputs(sam_hdr_tid2name(h, last_tid), file_out);
fprintf(file_out, "\t%"PRIhts_pos, last_pos+1);
for (i = 0; i < n; i++)
putchar('\t'), putchar('0');
putchar('\n');
fputc('\t', file_out), fputc('0', file_out);
fputc('\n', file_out);
}
}
last_tid++;
......@@ -237,19 +309,21 @@ int main_depth(int argc, char *argv[])
// Deal with missing portion of current tid
while (++last_pos < pos) {
if (last_pos < beg) continue; // out of range; skip
if (bed && bed_overlap(bed, h->target_name[tid], last_pos, last_pos + 1) == 0)
if (bed && bed_overlap(bed, sam_hdr_tid2name(h, tid), last_pos, last_pos + 1) == 0)
continue;
fputs(h->target_name[tid], stdout); printf("\t%d", last_pos+1);
fputs(sam_hdr_tid2name(h, tid), file_out);
fprintf(file_out, "\t%"PRIhts_pos, last_pos+1);
for (i = 0; i < n; i++)
putchar('\t'), putchar('0');
putchar('\n');
fputc('\t', file_out), fputc('0', file_out);
fputc('\n', file_out);
}
last_tid = tid;
last_pos = pos;
}
if (bed && bed_overlap(bed, h->target_name[tid], pos, pos + 1) == 0) continue;
fputs(h->target_name[tid], stdout); printf("\t%d", pos+1); // a customized printf() would be faster
if (bed && bed_overlap(bed, sam_hdr_tid2name(h, tid), pos, pos + 1) == 0) continue;
fputs(sam_hdr_tid2name(h, tid), file_out);
fprintf(file_out, "\t%"PRIhts_pos, pos+1); // a customized printf() would be faster
for (i = 0; i < n; ++i) { // base level filters have to go here
int j, m = 0;
for (j = 0; j < n_plp[i]; ++j) {
......@@ -258,9 +332,9 @@ int main_depth(int argc, char *argv[])
else if (p->qpos < p->b->core.l_qseq &&
bam_get_qual(p->b)[p->qpos] < baseQ) ++m; // low base quality
}
printf("\t%d", n_plp[i] - m); // this the depth to output
fprintf(file_out, "\t%d", n_plp[i] - m); // this the depth to output
}
putchar('\n');
fputc('\n', file_out);
}
if (ret < 0) status = EXIT_FAILURE;
free(n_plp); free(plp);
......@@ -268,19 +342,20 @@ int main_depth(int argc, char *argv[])
if (all) {
// Handle terminating region
if (last_tid < 0 && reg && all > 1) {
if (last_tid < 0 && reg) {
last_tid = reg_tid;
last_pos = beg-1;
}
while (last_tid >= 0 && last_tid < h->n_targets) {
while (++last_pos < h->target_len[last_tid]) {
while (last_tid >= 0 && last_tid < sam_hdr_nref(h)) {
while (++last_pos < sam_hdr_tid2len(h, last_tid)) {
if (last_pos >= end) break;
if (bed && bed_overlap(bed, h->target_name[last_tid], last_pos, last_pos + 1) == 0)
if (bed && bed_overlap(bed, sam_hdr_tid2name(h, last_tid), last_pos, last_pos + 1) == 0)
continue;
fputs(h->target_name[last_tid], stdout); printf("\t%d", last_pos+1);
fputs(sam_hdr_tid2name(h, last_tid), file_out);
fprintf(file_out, "\t%"PRIhts_pos, last_pos+1);
for (i = 0; i < n; i++)
putchar('\t'), putchar('0');
putchar('\n');
fputc('\t', file_out), fputc('0', file_out);
fputc('\n', file_out);
}
last_tid++;
last_pos = -1;
......@@ -290,8 +365,17 @@ int main_depth(int argc, char *argv[])
}
depth_end:
if (fclose(file_out) != 0) {
if (status == EXIT_SUCCESS) {
print_error_errno("depth", "error on closing \"%s\"",
(output_file && strcmp(output_file, "-") != 0
? output_file : "stdout"));
status = EXIT_FAILURE;
}
}
for (i = 0; i < n && data[i]; ++i) {
bam_hdr_destroy(data[i]->hdr);
sam_hdr_destroy(data[i]->hdr);
if (data[i]->fp) sam_close(data[i]->fp);
hts_itr_destroy(data[i]->iter);
free(data[i]);
......
bam_addrprg.c
View file @ 18118746
/* bam_addrprg.c -- samtools command to add or replace readgroups.
Copyright (c) 2013, 2015, 2016 Genome Research Limited.
Copyright (c) 2013, 2015-2017, 2019 Genome Research Limited.
Author: Martin O. Pollard <mp15@sanger.ac.uk>
......@@ -47,6 +47,7 @@ struct parsed_opts {
char* output_name;
char* rg_id;
char* rg_line;
int no_pg;
rg_mode mode;
sam_global_args ga;
htsThreadPool p;
......@@ -58,9 +59,9 @@ typedef struct state state_t;
struct state {
samFile* input_file;
bam_hdr_t* input_header;
sam_hdr_t* input_header;
samFile* output_file;
bam_hdr_t* output_header;
sam_hdr_t* output_header;
char* rg_id;
void (*mode_func)(const state_t*, bam1_t*);
};
......@@ -71,6 +72,7 @@ static void cleanup_opts(parsed_opts_t* opts)
free(opts->rg_id);
free(opts->output_name);
free(opts->input_name);
free(opts->rg_line);
if (opts->p.pool) hts_tpool_destroy(opts->p.pool);
sam_global_args_free(&opts->ga);
free(opts);
......@@ -81,9 +83,9 @@ static void cleanup_state(state_t* state)
if (!state) return;
free(state->rg_id);
if (state->output_file) sam_close(state->output_file);
bam_hdr_destroy(state->output_header);
sam_hdr_destroy(state->output_header);
if (state->input_file) sam_close(state->input_file);
bam_hdr_destroy(state->input_header);
sam_hdr_destroy(state->input_header);
free(state);
}
......@@ -147,20 +149,6 @@ static char *dup_substring(const char *s, const char *slim, size_t *lenp)
return ns;
}
// These are to be replaced by samtools header parser
// Extracts the first @RG line from a string.
static char* get_rg_line(const char* text, size_t* last)
{
const char* rg = text;
if (rg[0] != '@' || rg[1] != 'R' || rg[2] != 'G' ) {
if ((rg = (const char*)strstr(text,"\n@RG")) == NULL) {
return NULL;
}
rg++;//skip initial \n
}
// duplicate the line for return
return dup_substring(rg, strchr(rg, '\n'), last);
}
// Given a @RG line return the id
static char* get_rg_id(const char *line)
......@@ -172,44 +160,6 @@ static char* get_rg_id(const char *line)
return dup_substring(id, strchr(id, '\t'), NULL);
}
// Confirms the existance of an RG line with a given ID in a bam header
static bool confirm_rg( const bam_hdr_t *hdr, const char* rgid )
{
assert( hdr != NULL && rgid != NULL );
const char *ptr = hdr->text;
bool found = false;
while (ptr != NULL && *ptr != '\0' && found == false ) {
size_t end = 0;
char* line = get_rg_line(ptr, &end);
if (line == NULL) break; // No more @RG
char* id;
if (((id = get_rg_id(line)) != NULL) && !strcmp(id, rgid)) {
found = true;
}
free(id);
free(line);
ptr += end;
}
return found;
}
static char* get_first_rgid( const bam_hdr_t *hdr )
{
assert( hdr != NULL );
const char *ptr = hdr->text;
char* found = NULL;
while (ptr != NULL && *ptr != '\0' && found == NULL ) {
size_t end = 0;
char* line = get_rg_line(ptr, &end);
if ( line ) {
found = get_rg_id(line);
} else break;
free(line);
ptr += end;
}
return found;
}
static void usage(FILE *fp)
{
......@@ -221,8 +171,9 @@ static void usage(FILE *fp)
" -o FILE Where to write output to [stdout]\n"
" -r STRING @RG line text\n"
" -R STRING ID of @RG line in existing header to use\n"
" --no-PG Do not add a PG line\n"
);
sam_global_opt_help(fp, "..O..@");
sam_global_opt_help(fp, "..O..@..");
}
static bool parse_args(int argc, char** argv, parsed_opts_t** opts)
......@@ -242,6 +193,7 @@ static bool parse_args(int argc, char** argv, parsed_opts_t** opts)
sam_global_args_init(&retval->ga);
static const struct option lopts[] = {
SAM_OPT_GLOBAL_OPTIONS(0, 0, 'O', 0, 0, '@'),
{"no-PG", no_argument, NULL, 1},
{ NULL, 0, NULL, 0 }
};
kstring_t rg_line = {0,0,NULL};
......@@ -280,6 +232,9 @@ static bool parse_args(int argc, char** argv, parsed_opts_t** opts)
usage(stdout);
free(retval);
return true;
case 1:
retval->no_pg = 1;
break;
case '?':
usage(stderr);
free(retval);
......@@ -316,6 +271,7 @@ static bool parse_args(int argc, char** argv, parsed_opts_t** opts)
cleanup_opts(retval);
return false;
}
free(retval->rg_line);
retval->rg_line = tmp;
}
retval->input_name = strdup(argv[optind+0]);
......@@ -375,7 +331,7 @@ static bool init(const parsed_opts_t* opts, state_t** state_out) {
}
retval->input_header = sam_hdr_read(retval->input_file);
retval->output_header = bam_hdr_dup(retval->input_header);
retval->output_header = sam_hdr_dup(retval->input_header);
if (opts->output_name) // File format auto-detection
sam_open_mode(output_mode + 1, opts->output_name, NULL);
retval->output_file = sam_open_format(opts->output_name == NULL?"-":opts->output_name, output_mode, &opts->ga.out);
......@@ -393,34 +349,39 @@ static bool init(const parsed_opts_t* opts, state_t** state_out) {
if (opts->rg_line) {
// Append new RG line to header.
// Check does not already exist
if ( confirm_rg(retval->output_header, opts->rg_id) ) {
kstring_t hdr_line = { 0, 0, NULL };
if (sam_hdr_find_line_id(retval->output_header, "RG", "ID", opts->rg_id, &hdr_line) == 0) {
fprintf(stderr, "[init] ID of new RG line specified conflicts with that of an existing header RG line. Overwrite not yet implemented.\n");
free(hdr_line.s);
return false;
}
retval->rg_id = strdup(opts->rg_id);
size_t new_len = strlen( retval->output_header->text ) + strlen( opts->rg_line ) + 2;
char* new_header = malloc(new_len);
if (!new_header) {
fprintf(stderr, "[init] Out of memory whilst writing new header.\n");
if (-1 == sam_hdr_add_lines(retval->output_header, opts->rg_line, strlen(opts->rg_line))) {
fprintf(stderr, "[init] Error adding RG line with ID:%s to the output header.\n", opts->rg_id);
return false;
}
if (opts->mode == overwrite_all &&
-1 == sam_hdr_remove_except(retval->output_header, "RG", "ID", opts->rg_id)) {
fprintf(stderr, "[init] Error removing the old RG lines from the output header.\n");
return false;
}
sprintf(new_header,"%s%s\n", retval->output_header->text, opts->rg_line);
free(retval->output_header->text);
retval->output_header->text = new_header;
retval->output_header->l_text = (int)new_len - 1;
retval->rg_id = strdup(opts->rg_id);
} else {
if (opts->rg_id) {
// Confirm what has been supplied exists
if ( !confirm_rg(retval->output_header, opts->rg_id) ) {
kstring_t hdr_line = { 0, 0, NULL };
if (sam_hdr_find_line_id(retval->output_header, "RG", "ID", opts->rg_id, &hdr_line) < 0) {
fprintf(stderr, "RG ID supplied does not exist in header. Supply full @RG line with -r instead?\n");
return false;
}
retval->rg_id = strdup(opts->rg_id);
free(hdr_line.s);
} else {
if ((retval->rg_id = get_first_rgid(retval->output_header)) == NULL ) {
kstring_t rg_id = { 0, 0, NULL };
if (sam_hdr_find_tag_id(retval->output_header, "RG", NULL, NULL, "ID", &rg_id) < 0) {
fprintf(stderr, "No RG specified on command line or in existing header.\n");
return false;
}
retval->rg_id = ks_release(&rg_id);
}
}
......@@ -436,12 +397,24 @@ static bool init(const parsed_opts_t* opts, state_t** state_out) {
return true;
}
static bool readgroupise(state_t* state)
static bool readgroupise(parsed_opts_t *opts, state_t* state, char *arg_list)
{
if (!opts->no_pg && sam_hdr_add_pg(state->output_header, "samtools",
"VN", samtools_version(),
arg_list ? "CL": NULL,
arg_list ? arg_list : NULL,
NULL))
return false;
if (sam_hdr_write(state->output_file, state->output_header) != 0) {
print_error_errno("addreplacerg", "[%s] Could not write header to output file", __func__);
return false;
}
char *idx_fn = NULL;
if (opts->ga.write_index) {
if (!(idx_fn = auto_index(state->output_file, opts->output_name, state->output_header)))
return false;
}
bam1_t* file_read = bam_init1();
int ret;
......@@ -451,14 +424,25 @@ static bool readgroupise(state_t* state)
if (sam_write1(state->output_file, state->output_header, file_read) < 0) {
print_error_errno("addreplacerg", "[%s] Could not write read to output file", __func__);
bam_destroy1(file_read);
free(idx_fn);
return false;
}
}
bam_destroy1(file_read);
if (ret != -1) {
print_error_errno("addreplacerg", "[%s] Error reading from input file", __func__);
free(idx_fn);
return false;
} else {
if (opts->ga.write_index) {
if (sam_idx_save(state->output_file) < 0) {
print_error_errno("addreplacerg", "[%s] Writing index failed", __func__);
free(idx_fn);
return false;
}
}
free(idx_fn);
return true;
}
}
......@@ -467,20 +451,25 @@ int main_addreplacerg(int argc, char** argv)
{
parsed_opts_t* opts = NULL;
state_t* state = NULL;
char *arg_list = stringify_argv(argc+1, argv-1);
if (!arg_list)
return EXIT_FAILURE;
if (!parse_args(argc, argv, &opts)) goto error;
if (opts == NULL) return EXIT_SUCCESS; // Not an error but user doesn't want us to proceed
if (!opts || !init(opts, &state)) goto error;
if (!readgroupise(state)) goto error;
if (opts) { // Not an error but user doesn't want us to proceed
if (!init(opts, &state) || !readgroupise(opts, state, arg_list))
goto error;
}
cleanup_state(state);
cleanup_opts(opts);
free(arg_list);
return EXIT_SUCCESS;
error:
cleanup_state(state);
cleanup_opts(opts);
free(arg_list);
return EXIT_FAILURE;
}
bam_aux.c
View file @ 18118746
/* bam_aux.c -- remaining aux field handling.
Copyright (C) 2008-2010, 2013 Genome Research Ltd.
Copyright (C) 2008-2010, 2013, 2015, 2019 Genome Research Ltd.
Portions copyright (C) 2011 Broad Institute.
Author: Heng Li <lh3@sanger.ac.uk>
......@@ -61,21 +61,15 @@ int bam_aux_drop_other(bam1_t *b, uint8_t *s)
return 0;
}
// Only here due to libbam.a being used by some applications.
int bam_parse_region(bam_header_t *header, const char *str, int *ref_id, int *beg, int *end)
{
const char *name_lim = hts_parse_reg(str, beg, end);
if (name_lim) {
char *name = malloc(name_lim - str + 1);
memcpy(name, str, name_lim - str);
name[name_lim - str] = '\0';
*ref_id = bam_name2id(header, name);
free(name);
}
else {
// not parsable as a region, but possibly a sequence named "foo:a"
*ref_id = bam_name2id(header, str);
*beg = 0; *end = INT_MAX;
}
if (*ref_id == -1) return -1;
return *beg <= *end? 0 : -1;
hts_pos_t beg64, end64;
int r;
r = sam_parse_region(header, str, ref_id, &beg64, &end64, 0) ? 0 : -1;
if (beg64 > INT_MAX || end64 > INT_MAX)
return -1;
*beg = beg64;
*end = end64;
return r;
}
bam_cat.c
View file @ 18118746
/* bam_cat.c -- efficiently concatenates bam files.
Copyright (C) 2008-2009, 2011-2013, 2015-2016 Genome Research Ltd.
Copyright (C) 2008-2009, 2011-2013, 2015-2017, 2019 Genome Research Ltd.
Modified SAMtools work copyright (C) 2010 Illumina, Inc.
Permission is hereby granted, free of charge, to any person obtaining a copy
......@@ -45,162 +45,43 @@ Illumina.
#include "htslib/bgzf.h"
#include "htslib/sam.h"
#include "htslib/cram.h"
#include "htslib/khash.h"
#include "htslib/kstring.h"
#include "samtools.h"
KHASH_MAP_INIT_STR(s2i, int)
// Bi-directional lookup.
// We can go from name to ID or ID to name.
typedef struct khash_s2i {
khash_t(s2i) *h;
int n_id, a_id;
const char **id; // map Nth entry back to key
const char **line;
} khash_s2i;
static int hash_s2i_inc(khash_s2i *hash, const char *str, const char *line, int *added) {
// loosly based on khash_str2int_inc
khint_t k;
int n;
if ( !hash ) return -1;
// inefficient, but works
char *my_str = strdup(str);
k = kh_put(s2i, hash->h, my_str, added);
if (*added == 0) {
free(my_str);
return kh_val(hash->h, k);
}
n = hash->n_id++;
kh_val(hash->h, k) = n;
if (hash->a_id <= n) {
const char **id;
hash->a_id = (n+1)*2;
if (!(id = realloc(hash->id, hash->a_id*sizeof(*hash->id))))
return -1;
hash->id = id;
if (!(id = realloc(hash->line, hash->a_id*sizeof(*hash->line))))
return -1;
hash->line = id;
}
hash->id[n] = my_str; // reverse map
if (line)
hash->line[n] = line;
return n;
}
khash_s2i *hash_s2i_create(void) {
khash_s2i *h = calloc(1, sizeof(*h));
if (!h)
return NULL;
h->h = kh_init(s2i);
if (!h->h) {
free(h);
return NULL;
}
return h;
}
static void hash_s2i_free(khash_s2i *hash) {
// based on khash_str2int_destroy_free
khint_t k;
if (!hash) return;
if (hash->h) {
for (k = 0; k < kh_end(hash->h); ++k)
if (kh_exist(hash->h, k)) free((char*)kh_key(hash->h, k));
kh_destroy(s2i, hash->h);
}
if (hash->id)
free(hash->id);
if (hash->line)
free(hash->line);
free(hash);
}
static khash_s2i *hash_rg(const bam_hdr_t *h) {
khash_s2i *rg2id = hash_s2i_create();
char *cp, *line;
int j, l;
if (!h)
return rg2id;
if (!rg2id)
return NULL;
cp = h->text;
for (l = 0; l+3 < h->l_text; l++) {
line = &cp[l];
if (!(cp[l] == '@' && cp[l+1] == 'R' && cp[l+2] == 'G')) {
while (l < h->l_text && cp[l] != '\n')
l++;
continue;
}
// Found an @RG line; add to hash
while (cp[l] != '\n') {
while (l < h->l_text && cp[l] != '\n' && cp[l] != '\t')
l++;
if (l+4 < h->l_text && cp[l+1] == 'I' && cp[l+2] == 'D')
break;
}
if (cp[l] == '\n')
continue;
l = (j = l+4);
while (l < h->l_text && cp[l] != '\n' && cp[l] != '\t')
l++;
// To do: save id and keep realloc as needed, as hash_s2i_inc strdups.
char *id = malloc(l-j+1);
strncpy(id, &cp[j], l-j);
id[l-j] = 0;
int added;
hash_s2i_inc(rg2id, id, line, &added);
free(id);
while (l < h->l_text && cp[l] != '\n')
l++;
}
return rg2id;
}
#include "sam_opts.h"
/*
* Check the files are consistent and capable of being concatenated.
* Also fills out the rg2id read-group hash and the version numbers
* and produces a new bam_hdr_t structure with merged RG lines.
* Note it is only a simple merge, as we lack the niceties of a proper
* header API.
* Also fills out the version numbers and produces a new sam_hdr_t
* structure with merged RG lines.
* Note it is only a simple merge.
*
* Returns updated header on success;
* NULL on failure.
*/
static bam_hdr_t *cram_cat_check_hdr(int nfn, char * const *fn, const bam_hdr_t *h,
khash_s2i **rg2id, int *vers_maj_p, int *vers_min_p) {
static sam_hdr_t *cram_cat_check_hdr(int nfn, char * const *fn, const sam_hdr_t *h,
int *vers_maj_p, int *vers_min_p) {
int i, vers_maj = -1, vers_min = -1;
bam_hdr_t *new_h = NULL;
sam_hdr_t *new_h = NULL, *old_h = NULL;
samFile *in = NULL;
kstring_t ks = KS_INITIALIZE;
if (h) {
new_h = bam_hdr_dup(h);
*rg2id = hash_rg(new_h);
new_h = sam_hdr_dup(h);
if (!new_h) {
fprintf(stderr, "[%s] ERROR: header duplication failed.\n",
__func__);
goto fail;
}
}
for (i = 0; i < nfn; ++i) {
samFile *in;
cram_fd *in_c;
khint_t ki;
int new_rg = -1;
int ki;
in = sam_open(fn[i], "rc");
if (in == 0) {
print_error_errno("cat", "fail to open file '%s'", fn[i]);
return NULL;
goto fail;
}
in_c = in->fp.cram;
......@@ -210,55 +91,81 @@ static bam_hdr_t *cram_cat_check_hdr(int nfn, char * const *fn, const bam_hdr_t
(vers_min != -1 && vers_min != vmin)) {
fprintf(stderr, "[%s] ERROR: input files have differing version numbers.\n",
__func__);
return NULL;
goto fail;
}
vers_maj = vmaj;
vers_min = vmin;
bam_hdr_t *old = sam_hdr_read(in);
khash_s2i *rg2id_in = hash_rg(old);
old_h = sam_hdr_read(in);
if (!old_h) {
fprintf(stderr, "[%s] ERROR: header reading for file '%s' filed.\n",
__func__, fn[i]);
goto fail;
}
if (!new_h) {
new_h = bam_hdr_dup(old);
*rg2id = hash_rg(new_h);
new_h = sam_hdr_dup(old_h);
if (!new_h) {
fprintf(stderr, "[%s] ERROR: header duplication for file '%s' failed.\n",
__func__, fn[i]);
goto fail;
}
sam_hdr_destroy(old_h);
sam_close(in);
continue;
}
// Add any existing @RG entries to our global @RG hash.
for (ki = 0; ki < rg2id_in->n_id; ki++) {
int added;
new_rg = hash_s2i_inc(*rg2id, rg2id_in->id[ki], rg2id_in->line[ki], &added);
//fprintf(stderr, "RG %s: #%d -> #%d\n",
// rg2id_in->id[ki], ki, new_rg);
if (added) {
// Also add to new_h
const char *line = rg2id_in->line[ki];
const char *line_end = line;
while (*line && *line_end++ != '\n')
;
new_h->l_text += line_end - line;
new_h->text = realloc(new_h->text, new_h->l_text+1);
strncat(&new_h->text[new_h->l_text - (line_end - line)],
line, line_end - line);
int old_count = sam_hdr_count_lines(old_h, "RG");
for (ki = 0; ki < old_count; ki++) {
const char *old_name = sam_hdr_line_name(old_h, "RG", ki);
if (old_name) {
int new_i = sam_hdr_line_index(new_h, "RG", old_name);
if (-1 == new_i) { // line does not exist in the new header
if (sam_hdr_find_line_pos(old_h, "RG", ki, &ks) ||
!ks.s || sam_hdr_add_lines(new_h, ks.s, ks.l)) {
fprintf(stderr, "[%s] ERROR: failed to add @RG line 'ID:%s' from file '%s'\n",
__func__, old_name, fn[i]);
goto fail;
}
ks_free(&ks);
}
} else {
fprintf(stderr, "[%s] ERROR: failed to read %d @RG line from file '%s'\n",
__func__, ki, fn[i]);
goto fail;
}
}
if (new_rg != ki && rg2id_in->n_id > 1) {
if (old_count > 1 && sam_hdr_count_lines(new_h, "RG") == old_count) {
for (ki = 0; ki < old_count; ki++) {
const char *old_name = sam_hdr_line_name(old_h, "RG", ki);
const char *new_name = sam_hdr_line_name(new_h, "RG", ki);
if (!old_name || !new_name || strcmp(old_name, new_name)) {
fprintf(stderr, "[%s] ERROR: Same size @RG lists but differing order / contents\n",
__func__);
return NULL;
goto fail;
}
}
}
hash_s2i_free(rg2id_in);
bam_hdr_destroy(old);
sam_hdr_destroy(old_h);
sam_close(in);
}
ks_free(&ks);
*vers_maj_p = vers_maj;
*vers_min_p = vers_min;
return new_h;
fail:
ks_free(&ks);
if (old_h) sam_hdr_destroy(old_h);
if (new_h) sam_hdr_destroy(new_h);
if (in) sam_close(in);
return NULL;
}
......@@ -289,22 +196,21 @@ static bam_hdr_t *cram_cat_check_hdr(int nfn, char * const *fn, const bam_hdr_t
* huffman code. In this situation we can change the meta-data in the
* compression header to renumber an RG value..
*/
int cram_cat(int nfn, char * const *fn, const bam_hdr_t *h, const char* outcram)
int cram_cat(int nfn, char * const *fn, const sam_hdr_t *h, const char* outcram, sam_global_args *ga, char *arg_list, int no_pg)
{
samFile *out;
cram_fd *out_c;
int i, vers_maj, vers_min;
khash_s2i *rg2id = NULL;
bam_hdr_t *new_h = NULL;
sam_hdr_t *new_h = NULL;
/* Check consistent versioning and compatible headers */
if (!(new_h = cram_cat_check_hdr(nfn, fn, h, &rg2id, &vers_maj, &vers_min)))
if (!(new_h = cram_cat_check_hdr(nfn, fn, h, &vers_maj, &vers_min)))
return -1;
/* Open the file with cram_vers */
char vers[100];
sprintf(vers, "%d.%d", vers_maj, vers_min);
out = sam_open(outcram, "wc");
out = sam_open_format(outcram, "wc", &ga->out);
if (out == 0) {
print_error_errno("cat", "fail to open output file '%s'", outcram);
return -1;
......@@ -313,7 +219,13 @@ int cram_cat(int nfn, char * const *fn, const bam_hdr_t *h, const char* outcram)
cram_set_option(out_c, CRAM_OPT_VERSION, vers);
//fprintf(stderr, "Creating cram vers %s\n", vers);
cram_fd_set_header(out_c, sam_hdr_parse_(new_h->text, new_h->l_text)); // needed?
if (!no_pg && sam_hdr_add_pg(new_h, "samtools",
"VN", samtools_version(),
arg_list ? "CL": NULL,
arg_list ? arg_list : NULL,
NULL))
return -1;
if (sam_hdr_write(out, new_h) < 0) {
print_error_errno("cat", "Couldn't write header");
return -1;
......@@ -323,7 +235,7 @@ int cram_cat(int nfn, char * const *fn, const bam_hdr_t *h, const char* outcram)
samFile *in;
cram_fd *in_c;
cram_container *c;
bam_hdr_t *old;
sam_hdr_t *old_h;
int new_rg = -1;
in = sam_open(fn[i], "rc");
......@@ -333,20 +245,29 @@ int cram_cat(int nfn, char * const *fn, const bam_hdr_t *h, const char* outcram)
}
in_c = in->fp.cram;
old = sam_hdr_read(in);
khash_s2i *rg2id_in = hash_rg(old);
old_h = sam_hdr_read(in);
if (!old_h) {
print_error("cat", "fail to read the header of file '%s'", fn[i]);
return -1;
}
// Compute RG mapping if suitable for changing.
if (rg2id_in->n_id == 1) {
int _;
new_rg = hash_s2i_inc(rg2id, rg2id_in->id[0], NULL, &_);
if (sam_hdr_count_lines(old_h, "RG") == 1) {
const char *old_name = sam_hdr_line_name(old_h, "RG", 0);
if (old_name) {
new_rg = sam_hdr_line_index(new_h, "RG", old_name);
if (new_rg < 0) {
print_error("cat", "fail to find @RG line '%s' in the new header", old_name);
return -1;
}
} else {
print_error("cat", "fail to find @RG line in file '%s'", fn[i]);
return -1;
}
} else {
new_rg = 0;
}
hash_s2i_free(rg2id_in);
// Copy contains and blocks within them
while ((c = cram_read_container(in_c))) {
cram_block *blk;
......@@ -400,13 +321,11 @@ int cram_cat(int nfn, char * const *fn, const bam_hdr_t *h, const char* outcram)
cram_free_container(c);
}
bam_hdr_destroy(old);
sam_hdr_destroy(old_h);
sam_close(in);
}
sam_close(out);
hash_s2i_free(rg2id);
bam_hdr_destroy(new_h);
sam_hdr_destroy(new_h);
return 0;
}
......@@ -419,7 +338,7 @@ int cram_cat(int nfn, char * const *fn, const bam_hdr_t *h, const char* outcram)
#define BGZF_EMPTY_BLOCK_SIZE 28
int bam_cat(int nfn, char * const *fn, const bam_hdr_t *h, const char* outbam)
int bam_cat(int nfn, char * const *fn, sam_hdr_t *h, const char* outbam, char *arg_list, int no_pg)
{
BGZF *fp, *in = NULL;
uint8_t *buf = NULL;
......@@ -433,6 +352,13 @@ int bam_cat(int nfn, char * const *fn, const bam_hdr_t *h, const char* outbam)
return -1;
}
if (h) {
if (!no_pg && sam_hdr_add_pg(h, "samtools",
"VN", samtools_version(),
arg_list ? "CL": NULL,
arg_list ? arg_list : NULL,
NULL))
goto fail;
if (bam_hdr_write(fp, h) < 0) {
print_error_errno("cat", "Couldn't write header");
goto fail;
......@@ -445,7 +371,7 @@ int bam_cat(int nfn, char * const *fn, const bam_hdr_t *h, const char* outbam)
goto fail;
}
for(i = 0; i < nfn; ++i){
bam_hdr_t *old;
sam_hdr_t *old;
int len,j;
in = strcmp(fn[i], "-")? bgzf_open(fn[i], "r") : bgzf_fdopen(fileno(stdin), "r");
......@@ -462,6 +388,13 @@ int bam_cat(int nfn, char * const *fn, const bam_hdr_t *h, const char* outbam)
goto fail;
}
if (h == 0 && i == 0) {
if (!no_pg && sam_hdr_add_pg(old, "samtools",
"VN", samtools_version(),
arg_list ? "CL": NULL,
arg_list ? arg_list : NULL,
NULL))
goto fail;
if (bam_hdr_write(fp, old) < 0) {
print_error_errno("cat", "Couldn't write header");
goto fail;
......@@ -507,7 +440,7 @@ int bam_cat(int nfn, char * const *fn, const bam_hdr_t *h, const char* outbam)
if (bgzf_raw_write(fp, ebuf, es) < 0) goto write_fail;
}
}
bam_hdr_destroy(old);
sam_hdr_destroy(old);
bgzf_close(in);
in = NULL;
}
......@@ -530,14 +463,25 @@ int bam_cat(int nfn, char * const *fn, const bam_hdr_t *h, const char* outbam)
int main_cat(int argc, char *argv[])
{
bam_hdr_t *h = 0;
sam_hdr_t *h = 0;
char *outfn = 0;
char **infns = NULL; // files to concatenate
int infns_size = 0;
int c, ret = 0;
int c, ret = 0, no_pg = 0;
samFile *in;
sam_global_args ga;
static const struct option lopts[] = {
SAM_OPT_GLOBAL_OPTIONS('-', '-', '-', 0, '-', '@'),
{"no-PG", no_argument, NULL, 1},
{ NULL, 0, NULL, 0 }
};
while ((c = getopt(argc, argv, "h:o:b:")) >= 0) {
char *arg_list = NULL;
sam_global_args_init(&ga);
while ((c = getopt_long(argc, argv, "h:o:b:", lopts, NULL)) >= 0) {
switch (c) {
case 'h': {
samFile *fph = sam_open(optarg, "r");
......@@ -573,9 +517,19 @@ int main_cat(int argc, char *argv[])
}
break;
}
case 1:
no_pg = 1;
break;
default:
if (parse_sam_global_opt(c, optarg, lopts, &ga) == 0) break;
}
}
if (!no_pg && !(arg_list = stringify_argv(argc+1, argv-1))) {
print_error("cat", "failed to create arg_list");
return 1;
}
// Append files specified in argv to the list.
int nargv_fns = argc - optind;
if (nargv_fns > 0) {
......@@ -592,6 +546,8 @@ int main_cat(int argc, char *argv[])
fprintf(stderr, "Options: -b FILE list of input BAM/CRAM file names, one per line\n");
fprintf(stderr, " -h FILE copy the header from FILE [default is 1st input file]\n");
fprintf(stderr, " -o FILE output BAM/CRAM\n");
fprintf(stderr, " --no-PG do not add a PG line\n");
sam_global_opt_help(stderr, "--..-@-.");
return 1;
}
......@@ -604,13 +560,13 @@ int main_cat(int argc, char *argv[])
switch (hts_get_format(in)->format) {
case bam:
sam_close(in);
if (bam_cat(infns_size+nargv_fns, infns, h, outfn? outfn : "-") < 0)
if (bam_cat(infns_size+nargv_fns, infns, h, outfn? outfn : "-", arg_list, no_pg) < 0)
ret = 1;
break;
case cram:
sam_close(in);
if (cram_cat(infns_size+nargv_fns, infns, h, outfn? outfn : "-") < 0)
if (cram_cat(infns_size+nargv_fns, infns, h, outfn? outfn : "-", &ga, arg_list, no_pg) < 0)
ret = 1;
break;
......@@ -629,9 +585,9 @@ int main_cat(int argc, char *argv[])
free(outfn);
free(infns);
free(arg_list);
if (h)
bam_hdr_destroy(h);
sam_hdr_destroy(h);
return ret;
}
bam_fastq.c 0 → 100644
View file @ 18118746
This diff is collapsed.
bam_import.c deleted 100644 → 0
View file @ e3705033
/* bam_import.c -- SAM format parsing.
Copyright (C) 2008-2013 Genome Research Ltd.
Author: Heng Li <lh3@sanger.ac.uk>
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 <zlib.h>
#include <stdio.h>
#include <string.h>
#include <unistd.h>
#include "htslib/kstring.h"
#include "bam.h"
#include "htslib/kseq.h"
KSTREAM_INIT(gzFile, gzread, 16384)
bam_header_t *sam_header_read2(const char *fn)
{
bam_header_t *header;
int c, dret, n_targets = 0;
gzFile fp;
kstream_t *ks;
kstring_t *str;
kstring_t samstr = { 0, 0, NULL };
if (fn == 0) return 0;
fp = (strcmp(fn, "-") == 0)? gzdopen(fileno(stdin), "r") : gzopen(fn, "r");
if (fp == 0) return 0;
ks = ks_init(fp);
str = (kstring_t*)calloc(1, sizeof(kstring_t));
while (ks_getuntil(ks, 0, str, &dret) > 0) {
ksprintf(&samstr, "@SQ\tSN:%s", str->s);
ks_getuntil(ks, 0, str, &dret);
ksprintf(&samstr, "\tLN:%d\n", atoi(str->s));
n_targets++;
if (dret != '\n')
while ((c = ks_getc(ks)) != '\n' && c != -1);
}
ks_destroy(ks);
gzclose(fp);
free(str->s); free(str);
header = sam_hdr_parse(samstr.l, samstr.s? samstr.s : "");
free(samstr.s);
fprintf(stderr, "[sam_header_read2] %d sequences loaded.\n", n_targets);
return header;
}
bam_index.c
View file @ 18118746
/* bam_index.c -- index and idxstats subcommands.
Copyright (C) 2008-2011, 2013, 2014 Genome Research Ltd.
Copyright (C) 2008-2011, 2013-2016, 2018, 2019 Genome Research Ltd.
Portions copyright (C) 2010 Broad Institute.
Portions copyright (C) 2013 Peter Cock, The James Hutton Institute.
......@@ -114,20 +114,20 @@ int bam_index(int argc, char *argv[])
* Returns 0 on success,
* -1 on failure.
*/
int slow_idxstats(samFile *fp, bam_hdr_t *header) {
int slow_idxstats(samFile *fp, sam_hdr_t *header) {
int ret, last_tid = -2;
bam1_t *b = bam_init1();
if (hts_set_opt(fp, CRAM_OPT_REQUIRED_FIELDS, SAM_RNAME | SAM_FLAG))
return -1;
uint64_t (*count0)[2] = calloc(header->n_targets+1, sizeof(*count0));
uint64_t (*count0)[2] = calloc(sam_hdr_nref(header)+1, sizeof(*count0));
uint64_t (*counts)[2] = count0+1;
if (!count0)
return -1;
while ((ret = sam_read1(fp, header, b)) >= 0) {
if (b->core.tid >= header->n_targets || b->core.tid < -1) {
if (b->core.tid >= sam_hdr_nref(header) || b->core.tid < -1) {
free(count0);
return -1;
}
......@@ -148,10 +148,10 @@ int slow_idxstats(samFile *fp, bam_hdr_t *header) {
if (ret == -1) {
int i;
for (i = 0; i < header->n_targets; i++) {
printf("%s\t%d\t%"PRIu64"\t%"PRIu64"\n",
header->target_name[i],
header->target_len[i],
for (i = 0; i < sam_hdr_nref(header); i++) {
printf("%s\t%"PRId64"\t%"PRIu64"\t%"PRIu64"\n",
sam_hdr_tid2name(header, i),
(int64_t) sam_hdr_tid2len(header, i),
counts[i][0], counts[i][1]);
}
printf("*\t0\t%"PRIu64"\t%"PRIu64"\n", counts[-1][0], counts[-1][1]);
......@@ -167,14 +167,14 @@ int slow_idxstats(samFile *fp, bam_hdr_t *header) {
static void usage_exit(FILE *fp, int exit_status)
{
fprintf(fp, "Usage: samtools idxstats [options] <in.bam>\n");
sam_global_opt_help(fp, "-.---@");
sam_global_opt_help(fp, "-.---@-.");
exit(exit_status);
}
int bam_idxstats(int argc, char *argv[])
{
hts_idx_t* idx;
bam_hdr_t* header;
sam_hdr_t* header;
samFile* fp;
int c;
......@@ -227,9 +227,9 @@ int bam_idxstats(int argc, char *argv[])
}
int i;
for (i = 0; i < header->n_targets; ++i) {
for (i = 0; i < sam_hdr_nref(header); ++i) {
// Print out contig name and length
printf("%s\t%d", header->target_name[i], header->target_len[i]);
printf("%s\t%"PRId64, sam_hdr_tid2name(header, i), (int64_t) sam_hdr_tid2len(header, i));
// Now fetch info about it from the meta bin
uint64_t u, v;
hts_idx_get_stat(idx, i, &u, &v);
......@@ -240,7 +240,7 @@ int bam_idxstats(int argc, char *argv[])
hts_idx_destroy(idx);
}
bam_hdr_destroy(header);
sam_hdr_destroy(header);
sam_close(fp);
return 0;
}
bam_lpileup.c
View file @ 18118746
......@@ -100,7 +100,7 @@ void bam_lplbuf_reset(bam_lplbuf_t *buf)
buf->n_nodes = 0;
}
static int tview_func(uint32_t tid, uint32_t pos, int n, const bam_pileup1_t *pl, void *data)
static int tview_func(uint32_t tid, hts_pos_t pos, int n, const bam_pileup1_t *pl, void *data)
{
bam_lplbuf_t *tv = (bam_lplbuf_t*)data;
freenode_t *p;
......
bam_lpileup.h
View file @ 18118746
......@@ -33,7 +33,7 @@ typedef struct __bam_lplbuf_t bam_lplbuf_t;
#ifndef BAM_PILEUP_F_DEFINED
#define BAM_PILEUP_F_DEFINED
typedef int (*bam_pileup_f)(uint32_t tid, uint32_t pos, int n, const bam_pileup1_t *pl, void *data);
typedef int (*bam_pileup_f)(uint32_t tid, hts_pos_t pos, int n, const bam_pileup1_t *pl, void *data);
#endif //BAM_PILEUP_F_DEFINED
......
bam_markdup.c
View file @ 18118746
This diff is collapsed.