code: def mergebydepth(args):
      
      p = OptionParser(mergebydepth.__doc__)
      p.add_option("--mindepth", default=3, type="int",
                   help="Minimum depth required")
      opts, args = p.parse_args(args)
  
      if len(args) != 2:
          sys.exit(not p.print_help())
  
      bedfile, fastafile = args
      mindepth = opts.mindepth
      bedgraph = make_bedgraph(bedfile)
  
      bedgraphfiltered = bedgraph + ".d{0}".format(mindepth)
      if need_update(bedgraph, bedgraphfiltered):
          filter([bedgraph, "--minaccn={0}".format(mindepth),
                  "--outfile={0}".format(bedgraphfiltered)])
  
      merged = bedgraphfiltered + ".merge.fasta"
      if need_update(bedgraphfiltered, merged):
          mergeBed(bedgraphfiltered, sorted=True)
docs: %prog mergebydepth reads.bed genome.fasta
  
      Similar to mergeBed, but only returns regions beyond certain depth.