gridf.py

   1 import multiprocessing
   2
   3 import Image, ImageDraw, ImageFilter
   4
   5 from geometry import V
   6 from manual import lines as g_grid, l2ad, intersection, line as g_line
   7 from intrsc import intersections_from_angl_dist
   8 from linef import line_from_angl_dist
   9 import pcf
  10
  11 class GridFittingFailedError(Exception):
  12     pass
  13
  14 class MyGaussianBlur(ImageFilter.Filter):
  15     name = "GaussianBlur"
  16
  17     def __init__(self, radius=2):
  18         self.radius = radius
  19     def filter(self, image):
  20         return image.gaussian_blur(self.radius)
  21
  22 def projection(point, line, vector):
  23     return V(*intersection(g_line(point, point + vector.normal), g_line(*line)))
  24
  25 def job_br1(args):
  26     X, Y, im_l, a, b, c, d, s, v1, v2, k, hough, size = args
  27     return [(distance(im_l,
  28                      get_grid(a + X[y] * s * v1,
  29                               b + Y[y] * s * v1,
  30                               c, d, hough, size),
  31                      size), a + X[y] * s * v1, b + Y[y] * s * v1) for y in range(2 *k)]
  32
  33 def job_br2(args):
  34     X, Y, im_l, a, b, c, d, s, v1, v2, k, hough, size = args
  35     return [(distance(im_l,
  36                      get_grid(a, b, c + X[y] * s * v2,
  37                               d + Y[y] * s * v2,
  38                               hough, size),
  39                      size), c + X[y] * s * v2, d + Y[y] * s * v2) for y in range(2 *k)]
  40
  41 def find(lines, size, l1, l2, bounds, hough, do_something, im_h):
  42     a, b, c, d = [V(*a) for a in bounds]
  43     l1 = line_from_angl_dist(l1, size)
  44     l2 = line_from_angl_dist(l2, size)
  45     v1 = V(*l1[0]) - V(*l1[1])
  46     v2 = V(*l2[0]) - V(*l2[1])
  47     a = projection(a, l1, v1)
  48     b = projection(b, l1, v1)
  49     c = projection(c, l2, v2)
  50     d = projection(d, l2, v2)
  51
  52     im_l = Image.new('L', size)
  53     dr_l = ImageDraw.Draw(im_l)
  54     for line in sum(lines, []):
  55         dr_l.line(line_from_angl_dist(line, size), width=1, fill=255)
  56     #im_l = im_l.filter(MyGaussianBlur(radius=30))
  57     #GaussianBlur is undocumented class, may not work in future versions of PIL
  58     #im_l = im_l.tostring()
  59     im_l = im_h.tostring() # hocus pocus
  60
  61     #from gridf_analyzer import error_surface
  62     #error_surface(im_l, a, b, c, d, hough, size, v1 ,v2)
  63
  64     grid = get_grid(a, b, c, d, hough, size)
  65     dist = distance(im_l, grid, size)
  66
  67     #let's try the bruteforce aproach:
  68     s = 0.001
  69     k = 50
  70     X, Y = [], []
  71     for i in range(-k, k):
  72         X.append(range(-k, k))
  73         Y.append(2*k*[i])
  74
  75     pool = multiprocessing.Pool(None)
  76
  77     tasks = [(X[x], Y[x], im_l, a, b, c, d, s, v1, v2, k, hough, size) for x in xrange(0, 2 * k)]
  78
  79     #start = time.time()
  80     opt_ab = pool.map(job_br1, tasks, 1)
  81     opt_cd = pool.map(job_br2, tasks, 1)
  82     an, bn, cn, dn = 4 * [0]
  83     d1 = 0
  84     for lst in opt_ab:
  85         for tpl in lst:
  86             if tpl[0] > d1:
  87                 d1 = tpl[0]
  88                 a, b = tpl[1], tpl[2]
  89     d1 = 0
  90     for lst in opt_cd:
  91         for tpl in lst:
  92             if tpl[0] > d1:
  93                 d1 = tpl[0]
  94                 c, d = tpl[1], tpl[2]
  95     #print time.time() - start
  96     grid = get_grid(a, b, c, d, hough, size)
  97     grid_lines = [[l2ad(l, size) for l in grid[0]], [l2ad(l, size) for l in grid[1]]]
  98     return grid, grid_lines
  99
 100     #old optimization experiments:
 101     print dist
 102
 103     path = [(0,0)] #MNTR
 104     s = 0.01
 105     for _ in range(10):
 106         ts1 = [(s, 0), (0, s), (-s, 0), (0, -s)]
 107         grids = [(get_grid(a + t[0] * v1, b + t[1] * v1,
 108                            c, d, hough, size), t) for t in ts1]
 109         distances = [distance(im_l, grid, size) for (grid, t) in grids]
 110         gradient = [(di - dist) for di in distances]
 111         gradient = [gradient[0] - gradient[2], gradient[1] - gradient[3]]
 112         norm = (gradient[0] ** 2 + gradient[1] ** 2) ** 0.5
 113         gradient = [g / (100 * norm) for g in gradient]
 114         path.append(gradient)
 115         a, b = a + gradient[0] * v1, b + gradient[1] * v1
 116         dist = distance(im_l, grid, size)
 117         print dist
 118
 119     ###MNTR
 120     import matplotlib.pyplot as plt
 121     from matplotlib import cm
 122     import pickle
 123
 124     X, Y, Z = pickle.load(open('surface250'))
 125
 126     plt.imshow(Z, cmap=cm.jet, interpolation='none',
 127                origin='upper', extent=(-0.250, 0.250, -0.250, 0.250), aspect='equal')
 128     plt.colorbar()
 129     plt.plot([y for (x, y) in path], [x for (x, y) in path], 'go-')
 130
 131     plt.show()
 132     ###MNTR
 133
 134     print "---"
 135
 136     s = 0.02
 137     while True:
 138         ts1 = [(s, 0), (-s, 0), (s, s), (-s, -s), (-s, s), (s, -s), (0, s),  (0, -s)]
 139         grids = [(get_grid(a, b,
 140                            c + t[0] * v2, d + t[1] * v2, hough, size), t) for t in ts1]
 141         distances = [(distance(im_l, grid, size),
 142                       grid, t) for grid, t in grids]
 143         distances.sort(reverse=True)
 144         if distances[0][0] > dist:
 145             dist = distances[0][0]
 146             grid = distances[0][1]
 147             t = distances[0][2]
 148             c, d = c + t[0] * v2, d + t[1] * v2
 149             print dist
 150             s *= 0.75
 151         else:
 152             break
 153
 154     grid_lines = [[l2ad(l, size) for l in grid[0]], [l2ad(l, size) for l in grid[1]]]
 155     return grid, grid_lines
 156
 157 def get_grid(a, b, c, d, hough, size):
 158     l1 = hough.lines_from_list([a, b])
 159     l2 = hough.lines_from_list([c, d])
 160     c = intersections_from_angl_dist([l1, l2], size, get_all=True)
 161     #TODO do something when a corner is outside the image
 162     corners = (c[0] + c[1])
 163     if len(corners) < 4:
 164         print l1, l2, c
 165         raise GridFittingFailedError
 166     grid = g_grid(corners)
 167     return grid
 168
 169 def line_out(line, size):
 170     for p in line:
 171         if p[0] < 0 or p[0] > size[0] or p[1] < 0 or p[1] > size[1]:
 172             return True
 173     else:
 174         return False
 175
 176 def distance(im_l, grid, size):
 177     im_g = Image.new('L', size)
 178     dr_g = ImageDraw.Draw(im_g)
 179     for line in grid[0] + grid[1]:
 180         dr_g.line(line, width=1, fill=255)
 181         if line_out(line, size):
 182             return 0
 183     #im_g = im_g.filter(MyGaussianBlur(radius=3))
 184     #GaussianBlur is undocumented class, may not work in future versions of PIL
 185     #im_d, distance = combine(im_l, im_g)
 186     distance = pcf.combine(im_l, im_g.tostring())
 187     return distance
 188
 189 def combine(bg, fg):
 190     bg_l = bg.load()
 191     fg_l = fg.load()
 192     #res = Image.new('L', fg.size)
 193     #res_l = res.load()
 194
 195     score = 0
 196     area = 0
 197
 198     for x in xrange(fg.size[0]):
 199         for y in xrange(fg.size[1]):
 200             if fg_l[x, y]:
 201                 #res_l[x, y] = bg_l[x, y] * fg_l[x, y]
 202                 score +=  bg_l[x, y]
 203                 area += 1
 204
 205     #return res, float(score)/area
 206     return None, float(score)/area