test_snps <- function() {
  ## FIXME: rewrite this test to use the output of bam_tally, once we
  ## can convert those results to a VCF.
  library(gmapR)
  genome <- GmapGenome("hg19_IGIS21", create = TRUE, k = 12)
  dir <- GmapSnpDirectory(genome)
  gr <- as(seqinfo(genome), "GenomicRanges")
  gr <- GRanges("3", IRanges(3e6,4e6))
  options(error=recover)
  snps(dir, genome = genome, which = gr) <- "~/share/tracks/00-All.vcf.gz"
  snps(genome, which = gr, name = "snps2") <-
    "~/share/tracks/00-All.vcf.gz"
}