

def coverage(args):
    '\n    %prog coverage fastafile ctg bedfile1 bedfile2 ..\n\n    Plot coverage from a set of BED files that contain the read mappings. The\n    paired read span will be converted to a new bedfile that contain the happy\n    mates. ctg is the chr/scf/ctg that you want to plot the histogram on.\n\n    If the bedfiles already contain the clone spans, turn on --spans.\n    '
    from jcvi.formats.bed import mates, bedpe
    p = OptionParser(coverage.__doc__)
    p.add_option('--ymax', default=None, type='int', help='Limit ymax [default: %default]')
    p.add_option('--spans', default=False, action='store_true', help='BED files already contain clone spans [default: %default]')
    (opts, args, iopts) = p.set_image_options(args, figsize='8x5')
    if (len(args) < 3):
        sys.exit((not p.print_help()))
    (fastafile, ctg) = args[0:2]
    bedfiles = args[2:]
    sizes = Sizes(fastafile)
    size = sizes.mapping[ctg]
    plt.figure(1, (iopts.w, iopts.h))
    ax = plt.gca()
    bins = 100
    lines = []
    legends = []
    not_covered = []
    yy = 0.9
    for (bedfile, c) in zip(bedfiles, 'rgbcky'):
        if (not opts.spans):
            pf = bedfile.rsplit('.', 1)[0]
            matesfile = (pf + '.mates')
            if need_update(bedfile, matesfile):
                (matesfile, matesbedfile) = mates([bedfile, '--lib'])
            bedspanfile = (pf + '.spans.bed')
            if need_update(matesfile, bedspanfile):
                (bedpefile, bedspanfile) = bedpe([bedfile, '--span', '--mates={0}'.format(matesfile)])
            bedfile = bedspanfile
        bedsum = Bed(bedfile).sum(seqid=ctg)
        notcoveredbases = (size - bedsum)
        legend = bedfile.split('.')[0]
        msg = '{0}: {1} bp not covered'.format(legend, thousands(notcoveredbases))
        not_covered.append(msg)
        ((print >> sys.stderr), msg)
        ax.text(0.1, yy, msg, color=c, size=9, transform=ax.transAxes)
        yy -= 0.08
        cov = Coverage(bedfile, sizes.filename)
        (x, y) = cov.get_plot_data(ctg, bins=bins)
        (line,) = ax.plot(x, y, '-', color=c, lw=2, alpha=0.5)
        lines.append(line)
        legends.append(legend)
    leg = ax.legend(lines, legends, shadow=True, fancybox=True)
    leg.get_frame().set_alpha(0.5)
    ylabel = 'Average depth per {0}Kb'.format(((size / bins) / 1000))
    ax.set_xlim(0, size)
    ax.set_ylim(0, opts.ymax)
    ax.set_xlabel(ctg)
    ax.set_ylabel(ylabel)
    set_human_base_axis(ax)
    figname = '{0}.{1}.pdf'.format(fastafile, ctg)
    savefig(figname, dpi=iopts.dpi, iopts=iopts)

dataset/ETHPy150Open tanghaibao/jcvi/graphics/assembly.py coverage/original