|
4
|
1 #!/usr/bin/env python
|
|
|
2
|
|
|
3
|
|
|
4 """
|
|
|
5 interpolation
|
|
|
6 """
|
|
|
7
|
|
|
8 __all__ = ['neighbors', 'linear_interpolation']
|
|
|
9
|
|
|
10 def neighbors(start, finish):
|
|
|
11 """
|
|
|
12 returns the neighbors in finish from start
|
|
|
13 assumes both are sorted
|
|
|
14 """
|
|
|
15 assert finish
|
|
|
16 index = 0
|
|
|
17 retval = []
|
|
|
18 for x in start:
|
|
|
19 # endpoints
|
|
|
20 if x < finish[0]:
|
|
|
21 retval.append((None, 0))
|
|
|
22 continue
|
|
|
23 if x > finish[-1]:
|
|
|
24 retval.append((len(finish)-1, None))
|
|
|
25 continue
|
|
|
26 # traverse
|
|
|
27 try:
|
|
|
28 while True:
|
|
|
29 if x < finish[index] or x > finish[index+1]:
|
|
|
30 index += 1
|
|
|
31 continue
|
|
|
32 else:
|
|
|
33 break
|
|
|
34 retval.append((index, index+1))
|
|
|
35 except IndexError:
|
|
|
36 retval.append((len(finish)-2, len(finish)-1))
|
|
|
37
|
|
|
38 return retval
|
|
|
39
|
|
|
40
|
|
|
41 def linear_interpolation(data, points):
|
|
|
42 """
|
|
|
43 linearly interpolate data to points
|
|
|
44
|
|
|
45 data -- iterable of 2-tuples (or equivalent) of `x,y`
|
|
|
46 points -- `x`-values to interpolate to
|
|
|
47 """
|
|
|
48
|
|
|
49 # ensure we are sorted
|
|
|
50 data = sorted(data, key=lambda x: x[0])
|
|
|
51 points = sorted(points)
|
|
|
52
|
|
|
53 # get the neighbors
|
|
|
54 x = [value[0] for value in data]
|
|
|
55 nearest_neighbors = neighbors(points, x)
|
|
|
56
|
|
|
57 # we don't support endpoints yet; this is interpolation, not extrapolation
|
|
|
58 assert not any([(neighbor[0] is None or neighbor[1] is None)
|
|
|
59 for neighbor in nearest_neighbors])
|
|
|
60
|
|
|
61 retval = []
|
|
|
62 for index, (left, right) in enumerate(nearest_neighbors):
|
|
|
63 # linearly interpolate
|
|
|
64 ratio = (points[index] - data[left][0])/float(data[right][0] - data[left][0])
|
|
|
65 retval.append(ratio*data[right][1] + (1.-ratio)*data[left][1])
|
|
|
66 return retval
|
|
|
67
|
|
|
68
|
