scripts/last-dotplot
changeset 639 047baa399047
parent 482 f5b731ee703e
child 640 374bc3ee1115
     1.1 --- a/scripts/last-dotplot	Mon Sep 29 07:37:32 2014 +0900
     1.2 +++ b/scripts/last-dotplot	Mon Oct 19 10:53:05 2015 +0900
     1.3 @@ -205,6 +205,30 @@
     1.4  seq_pix2, seq_starts2, height = get_pix_info(seq_sizes2, margin2)
     1.5  seq_start_dic1 = dict(zip(seq_names1, seq_starts1))
     1.6  seq_start_dic2 = dict(zip(seq_names2, seq_starts2))
     1.7 +
     1.8 +def drawLineForward(hits, width, bp_per_pix, origin, beg1, beg2, size):
     1.9 +    while True:
    1.10 +        q1, r1 = divmod(beg1, bp_per_pix)
    1.11 +        q2, r2 = divmod(beg2, bp_per_pix)
    1.12 +        hits[origin + q2 * width + q1] |= 1
    1.13 +        next_pix = min(bp_per_pix - r1, bp_per_pix - r2)
    1.14 +        if next_pix >= size: break
    1.15 +        beg1 += next_pix
    1.16 +        beg2 += next_pix
    1.17 +        size -= next_pix
    1.18 +
    1.19 +def drawLineReverse(hits, width, bp_per_pix, origin, beg1, beg2, size):
    1.20 +    beg2 = -1 - beg2
    1.21 +    while True:
    1.22 +        q1, r1 = divmod(beg1, bp_per_pix)
    1.23 +        q2, r2 = divmod(beg2, bp_per_pix)
    1.24 +        hits[origin + q2 * width + q1] |= 2
    1.25 +        next_pix = min(bp_per_pix - r1, r2 + 1)
    1.26 +        if next_pix >= size: break
    1.27 +        beg1 += next_pix
    1.28 +        beg2 -= next_pix
    1.29 +        size -= next_pix
    1.30 +
    1.31  hits = [0] * (width * height)  # the image data
    1.32  
    1.33  sys.stderr.write(my_name + ": processing alignments...\n")
    1.34 @@ -213,18 +237,15 @@
    1.35      seq_start2 = seq_start_dic2[seq2]
    1.36      my_start = seq_start2 * width + seq_start1
    1.37      for beg1, beg2, size in blocks:
    1.38 -        end1 = beg1 + size
    1.39 -        end2 = beg2 + size
    1.40 -        if beg1 >= 0: j = xrange(beg1, end1)
    1.41 -        else:         j = xrange(-1 - beg1, -1 - end1, -1)
    1.42 -        if beg2 >= 0: k = xrange(beg2, end2)
    1.43 -        else:         k = xrange(-1 - beg2, -1 - end2, -1)
    1.44 -        if beg2 >= 0: store_value = 1
    1.45 -        else:         store_value = 2
    1.46 -        for real_pos1, real_pos2 in itertools.izip(j, k):
    1.47 -            pix1 = real_pos1 // bp_per_pix
    1.48 -            pix2 = real_pos2 // bp_per_pix
    1.49 -            hits[my_start + pix2 * width + pix1] |= store_value
    1.50 +        if beg1 < 0:
    1.51 +            beg1 = -(beg1 + size)
    1.52 +            beg2 = -(beg2 + size)
    1.53 +        if beg2 >= 0:
    1.54 +            drawLineForward(hits, width, bp_per_pix, my_start,
    1.55 +                            beg1, beg2, size)
    1.56 +        else:
    1.57 +            drawLineReverse(hits, width, bp_per_pix, my_start,
    1.58 +                            beg1, beg2, size)
    1.59  sys.stderr.write(my_name + ": done\n")
    1.60  
    1.61  def make_label(text, text_size, range_start, range_size):