

def qc(args):
    '\n    %prog qc prefix\n\n    Expects data files including:\n    1. `prefix.bedpe` draws Bezier curve between paired reads\n    2. `prefix.sizes` draws length of the contig/scaffold\n    3. `prefix.gaps.bed` mark the position of the gaps in sequence\n    4. `prefix.bed.coverage` plots the base coverage\n    5. `prefix.pairs.bed.coverage` plots the clone coverage\n\n    See assembly.coverage.posmap() for the generation of these files.\n    '
    from jcvi.graphics.glyph import Bezier
    p = OptionParser(qc.__doc__)
    (opts, args) = p.parse_args(args)
    if (len(args) != 1):
        sys.exit(p.print_help())
    (prefix,) = args
    scf = prefix
    bedpefile = (prefix + '.bedpe')
    fastafile = (prefix + '.fasta')
    sizesfile = (prefix + '.sizes')
    gapsbedfile = (prefix + '.gaps.bed')
    bedfile = (prefix + '.bed')
    bedpefile = (prefix + '.bedpe')
    pairsbedfile = (prefix + '.pairs.bed')
    sizes = Sizes(fastafile).mapping
    size = sizes[scf]
    fig = plt.figure(1, (8, 5))
    root = fig.add_axes([0, 0, 1, 1])
    root.add_patch(Rectangle((0.1, 0.15), 0.8, 0.03, fc='k'))
    ax = fig.add_axes([0.1, 0.45, 0.8, 0.45])
    bins = 200
    basecoverage = Coverage(bedfile, sizesfile)
    matecoverage = Coverage(pairsbedfile, sizesfile)
    (x, y) = basecoverage.get_plot_data(scf, bins=bins)
    (baseline,) = ax.plot(x, y, 'g-')
    (x, y) = matecoverage.get_plot_data(scf, bins=bins)
    (mateline,) = ax.plot(x, y, 'r-')
    legends = ('Base coverage', 'Mate coverage')
    leg = ax.legend((baseline, mateline), legends, shadow=True, fancybox=True)
    leg.get_frame().set_alpha(0.5)
    ax.set_xlim(0, size)
    fp = open(bedpefile)
    pairs = []
    for row in fp:
        (scf, astart, aend, scf, bstart, bend, clonename) = row.split()
        (astart, bstart) = (int(astart), int(bstart))
        (aend, bend) = (int(aend), int(bend))
        start = (min(astart, bstart) + 1)
        end = max(aend, bend)
        pairs.append((start, end))
    bpratio = (0.8 / size)
    cutoff = 1000
    pos = (lambda x: (0.1 + (x * bpratio)))
    ypos = (0.15 + 0.03)
    for (start, end) in pairs:
        dist = (end - start)
        if (dist < cutoff):
            continue
        dist = min(dist, 10000)
        height = ((0.25 * dist) / 10000)
        xstart = pos(start)
        xend = pos(end)
        p0 = (xstart, ypos)
        p1 = (xstart, (ypos + height))
        p2 = (xend, (ypos + height))
        p3 = (xend, ypos)
        Bezier(root, p0, p1, p2, p3)
    fp = open(gapsbedfile)
    for row in fp:
        b = BedLine(row)
        (start, end) = (b.start, b.end)
        xstart = pos(start)
        xend = pos(end)
        root.add_patch(Rectangle((xstart, 0.15), (xend - xstart), 0.03, fc='w'))
    root.text(0.5, 0.1, scf, color='b', ha='center')
    warn_msg = 'Only the inserts > {0}bp are shown'.format(cutoff)
    root.text(0.5, 0.1, scf, color='b', ha='center')
    root.text(0.5, 0.05, warn_msg, color='gray', ha='center')
    set_human_base_axis(ax)
    root.set_xlim(0, 1)
    root.set_ylim(0, 1)
    root.set_axis_off()
    figname = (prefix + '.pdf')
    savefig(figname, dpi=300)

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