fix bug in gridf2
[imago.git] / src / geometry.py
1 """Imago geometry module."""
2
3 from math import sin, cos, atan, pi
4
5 class V(object):
6     """Class for vector manipulation."""
7
8     def __init__(self, x, y):
9         self.x = x
10         self.y = y
11     
12     def __add__(self, other):
13         return V(self.x + other.x, self.y + other.y)
14
15     def __sub__(self, other):
16         return V(self.x - other.x, self.y - other.y)
17
18     def __rmul__(self, other):
19         return V(other * self.x, other * self.y)
20
21     def __len__(self):
22         return 2
23
24     def __getitem__(self, key):
25         if key == 0:
26             return self.x
27         elif key == 1:
28             return self.y
29         elif type(key) != int:
30             raise TypeError("V indices must be integers") 
31         else:
32             raise KeyError("V index ({}) out of range".format(key))
33
34     def __iter__(self):
35         yield self.x
36         yield self.y
37
38     @property
39     def normal(self):
40         return V(-self.y, self.x)
41
42 def projection(p, l, v):
43     #TODO what is this?
44     return V(*intersection(line(p, p + v.normal), line(*l)))
45
46 def l2ad((a, b), size):
47     """Represent line as (angle, distance).
48     
49     Take a line (represented by two points) and image size.
50     Return the line represented by its angle and distance
51     from the center of the image.
52
53     """
54     if (a[0] - b[0]) == 0:
55         angle = pi / 2
56     else:
57         q = float(a[1] - b[1]) / (a[0] - b[0])
58         angle = atan(q)
59
60     if angle < 0:
61         angle += pi
62     if angle > pi:
63         angle -= pi
64
65     distance = (((a[0] - (size[0] / 2)) * sin(angle)) + 
66                 ((a[1] - (size[1] / 2)) * - cos(angle)))
67     return (angle, distance)
68
69 def line(x, y):
70     """Return parametric representation of line."""
71     a = x[1] - y[1]
72     b = y[0] - x[0]
73     c = a * y[0] + b * y[1]
74     return (a, b, c)
75
76 def intersection(p, q):
77     """Return intersection of two lines."""
78     det = p[0] * q[1] - p[1] * q[0]
79     if det == 0:
80         return None
81     return (int(round(float(q[1] * p[2] - p[1] * q[2]) / det)), 
82             int(round(float(p[0] * q[2] - q[0] * p[2]) / det)))