WIP: adding new tests and fixing bugs in pg.interpolateArray

This commit is contained in:
Luke Campagnola 2015-02-16 11:04:58 -05:00
parent 1d7cbca64d
commit accafcce36
3 changed files with 54 additions and 27 deletions

View File

@ -460,7 +460,7 @@ def affineSlice(data, shape, origin, vectors, axes, order=1, returnCoords=False,
ind = (Ellipsis,) + inds ind = (Ellipsis,) + inds
output[ind] = scipy.ndimage.map_coordinates(data[ind], x, order=order, **kargs) output[ind] = scipy.ndimage.map_coordinates(data[ind], x, order=order, **kargs)
else: else:
# map_coordinates expects the indexes as the first axis, whereas # map_coordinates expects the indexes as the first axis, whereas
# interpolateArray expects indexes at the last axis. # interpolateArray expects indexes at the last axis.
tr = tuple(range(1,x.ndim)) + (0,) tr = tuple(range(1,x.ndim)) + (0,)
output = interpolateArray(data, x.transpose(tr)) output = interpolateArray(data, x.transpose(tr))
@ -483,7 +483,7 @@ def affineSlice(data, shape, origin, vectors, axes, order=1, returnCoords=False,
def interpolateArray(data, x, default=0.0): def interpolateArray(data, x, default=0.0):
""" """
N-dimensional interpolation similar scipy.ndimage.map_coordinates. N-dimensional interpolation similar to scipy.ndimage.map_coordinates.
This function returns linearly-interpolated values sampled from a regular This function returns linearly-interpolated values sampled from a regular
grid of data. grid of data.
@ -492,7 +492,7 @@ def interpolateArray(data, x, default=0.0):
*x* is an array with (shape[-1] <= data.ndim) containing the locations *x* is an array with (shape[-1] <= data.ndim) containing the locations
within *data* to interpolate. within *data* to interpolate.
Returns array of shape (x.shape[:-1] + data.shape) Returns array of shape (x.shape[:-1] + data.shape[x.shape[-1]:])
For example, assume we have the following 2D image data:: For example, assume we have the following 2D image data::
@ -535,7 +535,7 @@ def interpolateArray(data, x, default=0.0):
This is useful for interpolating from arrays of colors, vertexes, etc. This is useful for interpolating from arrays of colors, vertexes, etc.
""" """
print "x:\n", x.shape, x
prof = debug.Profiler() prof = debug.Profiler()
nd = data.ndim nd = data.ndim
@ -552,7 +552,7 @@ def interpolateArray(data, x, default=0.0):
for ax in range(md): for ax in range(md):
mask = (xmin[...,ax] >= 0) & (x[...,ax] <= data.shape[ax]-1) mask = (xmin[...,ax] >= 0) & (x[...,ax] <= data.shape[ax]-1)
# keep track of points that need to be set to default # keep track of points that need to be set to default
totalMask &= mask totalMask &= mask
# ..and keep track of indexes that are out of bounds # ..and keep track of indexes that are out of bounds
# (note that when x[...,ax] == data.shape[ax], then xmax[...,ax] will be out # (note that when x[...,ax] == data.shape[ax], then xmax[...,ax] will be out
@ -564,7 +564,7 @@ def interpolateArray(data, x, default=0.0):
axisIndex[axisIndex >= data.shape[ax]] = 0 axisIndex[axisIndex >= data.shape[ax]] = 0
fieldInds.append(axisIndex) fieldInds.append(axisIndex)
prof() prof()
print "fieldInds:\n", fieldInds
# Get data values surrounding each requested point # Get data values surrounding each requested point
# fieldData[..., i] contains all 2**nd values needed to interpolate x[i] # fieldData[..., i] contains all 2**nd values needed to interpolate x[i]
fieldData = data[tuple(fieldInds)] fieldData = data[tuple(fieldInds)]
@ -585,8 +585,11 @@ def interpolateArray(data, x, default=0.0):
result = result.sum(axis=0) result = result.sum(axis=0)
prof() prof()
totalMask.shape = totalMask.shape + (1,) * (nd - md) #totalMask.shape = totalMask.shape + (1,) * (nd - md)
print "mask:\n", totalMask.shape, totalMask
print "result:\n", result.shape, result
result[~totalMask] = default result[~totalMask] = default
print "masked:\n", result.shape, result
prof() prof()
return result return result

View File

@ -1061,8 +1061,8 @@ class ROI(GraphicsObject):
=================== ==================================================== =================== ====================================================
This method uses :func:`affineSlice <pyqtgraph.affineSlice>` to generate This method uses :func:`affineSlice <pyqtgraph.affineSlice>` to generate
the slice from *data* and uses :func:`getAffineSliceParams <pyqtgraph.ROI.getAffineSliceParams>` to determine the parameters to the slice from *data* and uses :func:`getAffineSliceParams <pyqtgraph.ROI.getAffineSliceParams>`
pass to :func:`affineSlice <pyqtgraph.affineSlice>`. to determine the parameters to pass to :func:`affineSlice <pyqtgraph.affineSlice>`.
If *returnMappedCoords* is True, then the method returns a tuple (result, coords) If *returnMappedCoords* is True, then the method returns a tuple (result, coords)
such that coords is the set of coordinates used to interpolate values from the original such that coords is the set of coordinates used to interpolate values from the original
@ -1079,24 +1079,16 @@ class ROI(GraphicsObject):
else: else:
kwds['returnCoords'] = True kwds['returnCoords'] = True
result, coords = fn.affineSlice(data, shape=shape, vectors=vectors, origin=origin, axes=axes, **kwds) result, coords = fn.affineSlice(data, shape=shape, vectors=vectors, origin=origin, axes=axes, **kwds)
#tr = fn.transformToArray(img.transform())[:2] ## remove perspective transform values
### separate translation from scale/rotate
#translate = tr[:,2]
#tr = tr[:,:2]
#tr = tr.reshape((2,2) + (1,)*(coords.ndim-1))
#coords = coords[np.newaxis, ...]
### map coordinates and return ### map coordinates and return
#mapped = (tr*coords).sum(axis=0) ## apply scale/rotate
#mapped += translate.reshape((2,1,1))
mapped = fn.transformCoordinates(img.transform(), coords) mapped = fn.transformCoordinates(img.transform(), coords)
return result, mapped return result, mapped
def getAffineSliceParams(self, data, img, axes=(0,1)): def getAffineSliceParams(self, data, img, axes=(0,1)):
""" """
Returns the parameters needed to use :func:`affineSlice <pyqtgraph.affineSlice>` to Returns the parameters needed to use :func:`affineSlice <pyqtgraph.affineSlice>`
extract a subset of *data* using this ROI and *img* to specify the subset. (shape, vectors, origin) to extract a subset of *data* using this ROI
and *img* to specify the subset.
See :func:`getArrayRegion <pyqtgraph.ROI.getArrayRegion>` for more information. See :func:`getArrayRegion <pyqtgraph.ROI.getArrayRegion>` for more information.
""" """
@ -1138,8 +1130,6 @@ class ROI(GraphicsObject):
relativeTo['scale'] = relativeTo['size'] relativeTo['scale'] = relativeTo['size']
st['scale'] = st['size'] st['scale'] = st['size']
t1 = SRTTransform(relativeTo) t1 = SRTTransform(relativeTo)
t2 = SRTTransform(st) t2 = SRTTransform(st)
return t2/t1 return t2/t1

View File

@ -22,18 +22,39 @@ def testSolve3D():
def test_interpolateArray(): def test_interpolateArray():
def interpolateArray(data, x):
result = pg.interpolateArray(data, x)
assert result.shape == x.shape[:-1] + data.shape[x.shape[-1]:]
return result
data = np.array([[ 1., 2., 4. ], data = np.array([[ 1., 2., 4. ],
[ 10., 20., 40. ], [ 10., 20., 40. ],
[ 100., 200., 400.]]) [ 100., 200., 400.]])
# test various x shapes
interpolateArray(data, np.ones((1,)))
interpolateArray(data, np.ones((2,)))
interpolateArray(data, np.ones((1, 1)))
interpolateArray(data, np.ones((1, 2)))
interpolateArray(data, np.ones((5, 1)))
interpolateArray(data, np.ones((5, 2)))
interpolateArray(data, np.ones((5, 5, 1)))
interpolateArray(data, np.ones((5, 5, 2)))
with pytest.raises(TypeError):
interpolateArray(data, np.ones((3,)))
with pytest.raises(TypeError):
interpolateArray(data, np.ones((1, 3,)))
with pytest.raises(TypeError):
interpolateArray(data, np.ones((5, 5, 3,)))
x = np.array([[ 0.3, 0.6], x = np.array([[ 0.3, 0.6],
[ 1. , 1. ], [ 1. , 1. ],
[ 0.5, 1. ], [ 0.5, 1. ],
[ 0.5, 2.5], [ 0.5, 2.5],
[ 10. , 10. ]]) [ 10. , 10. ]])
result = pg.interpolateArray(data, x) result = interpolateArray(data, x)
#import scipy.ndimage #import scipy.ndimage
#spresult = scipy.ndimage.map_coordinates(data, x.T, order=1) #spresult = scipy.ndimage.map_coordinates(data, x.T, order=1)
spresult = np.array([ 5.92, 20. , 11. , 0. , 0. ]) # generated with the above line spresult = np.array([ 5.92, 20. , 11. , 0. , 0. ]) # generated with the above line
@ -44,9 +65,10 @@ def test_interpolateArray():
x = np.array([[ 0.3, 0], x = np.array([[ 0.3, 0],
[ 0.3, 1], [ 0.3, 1],
[ 0.3, 2]]) [ 0.3, 2]])
r1 = interpolateArray(data, x)
x = np.array([0.3]) # should broadcast across axis 1
r2 = interpolateArray(data, x)
r1 = pg.interpolateArray(data, x)
r2 = pg.interpolateArray(data, x[0,:1])
assert_array_almost_equal(r1, r2) assert_array_almost_equal(r1, r2)
@ -54,13 +76,25 @@ def test_interpolateArray():
x = np.array([[[0.5, 0.5], [0.5, 1.0], [0.5, 1.5]], x = np.array([[[0.5, 0.5], [0.5, 1.0], [0.5, 1.5]],
[[1.5, 0.5], [1.5, 1.0], [1.5, 1.5]]]) [[1.5, 0.5], [1.5, 1.0], [1.5, 1.5]]])
r1 = pg.interpolateArray(data, x) r1 = interpolateArray(data, x)
#r2 = scipy.ndimage.map_coordinates(data, x.transpose(2,0,1), order=1) #r2 = scipy.ndimage.map_coordinates(data, x.transpose(2,0,1), order=1)
r2 = np.array([[ 8.25, 11. , 16.5 ], # generated with the above line r2 = np.array([[ 8.25, 11. , 16.5 ], # generated with the above line
[ 82.5 , 110. , 165. ]]) [ 82.5 , 110. , 165. ]])
assert_array_almost_equal(r1, r2) assert_array_almost_equal(r1, r2)
# test interpolate where data.ndim > x.shape[1]
data = np.array([[[1, 2, 3], [4, 5, 6]], [[7, 8, 9], [10, 11, 12]]]) # 2x2x3
x = np.array([[1, 1], [0, 0.5], [5, 5]])
r1 = interpolateArray(data, x)
assert np.all(r1[0] == data[1, 1])
assert np.all(r1[1] == 0.5 * (data[0, 0] + data[0, 1]))
assert np.all(r1[2] == 0)
def test_subArray(): def test_subArray():
a = np.array([0, 0, 111, 112, 113, 0, 121, 122, 123, 0, 0, 0, 211, 212, 213, 0, 221, 222, 223, 0, 0, 0, 0]) a = np.array([0, 0, 111, 112, 113, 0, 121, 122, 123, 0, 0, 0, 211, 212, 213, 0, 221, 222, 223, 0, 0, 0, 0])
b = pg.subArray(a, offset=2, shape=(2,2,3), stride=(10,4,1)) b = pg.subArray(a, offset=2, shape=(2,2,3), stride=(10,4,1))