code: def chimera(args):
      
      import pysam
      from jcvi.utils.natsort import natsorted
  
      p = OptionParser(chimera.__doc__)
      p.set_verbose()
      opts, args = p.parse_args(args)
      if len(args) != 1:
          sys.exit(not p.print_help())
  
      samfile, = args
      samfile = pysam.AlignmentFile(samfile)
      rstore = defaultdict(list)
      hstore = defaultdict(int)
      for r in samfile.fetch():
          rstore[r.query_name] += list(breakpoint(r))
          hstore[r.query_name] += 1
          if opts.verbose:
              print(r.query_name, "+-"[r.is_reverse], \
                          sum(l for o, l in r.cigartuples), r.cigarstring, list(breakpoint(r)), file=sys.stderr)
  
      for rn, bps in natsorted(rstore.items()):
          bps = "|".join(str(x) for x in sorted(bps)) if bps else "na"
          print("\t".join((rn, str(hstore[rn]), bps)))
docs: %prog chimera bamfile
  
      Parse BAM file from `bwasw` and list multi-hit reads and breakpoints.