On Sat, Jul 7, 2012 at 7:30 PM, Amit Aronovitch <aro...@gm...>wrote: > > The current implementation of Delaunay interpolator returns NaN for grids > whose x or y has dimension 1 (i.e. when you try to interpolate along a > horizontal/vertical line, or in a single point). See example below. > > By looking at the code, it seems that this can be fixed by simple > rearrangement of calculations. > Suggested patch provided here: > https://github.com/AmitAronovitch/matplotlib/commit/f312d864da9c72681eb3db3b5920ae64793c713e(let me know if you want a pull request). > The suggested implementation is almost identical. It might actually > perform faster in some cases (there is one less multiplication op in the > inner loop). There might be some differences in accuracy, but I believe > they should only become observable in cases where the grid size is very > large (which would probably cause memory problems anyway). > > Example (before suggested patch): > > >>> from matplotlib.delaunay import Triangulation > >>> tri = Triangulation([0,10,10,0],[0,0,10,10]) > >>> lin = tri.linear_interpolator([1,10,5,2.0]) > >>> # 2x2 grid works fine > >>> lin[3:6:2j,1:4:2j] > array([[ 1.6, 3.1], > [ 1.9, 2.8]]) > >>> # but not when 1x2, 2x1, 1x1: > >>> lin[3:6:2j,1:1:1j] > array([[ nan], > [ nan]]) > >>> lin[3:3:1j,1:1:1j] > array([[ nan]]) > >>> > > After suggested patch: > > >>> from matplotlib.delaunay import Triangulation > >>> tri = Triangulation([0,10,10,0],[0,0,10,10]) > >>> lin = tri.linear_interpolator([1,10,5,2.0]) > >>> # 2x2 grid: same same > >>> lin[3:6:2j,1:4:2j] > array([[ 1.6, 3.1], > [ 1.9, 2.8]]) > >>> # but these work now > >>> lin[3:6:2j,1:1:1j] > array([[ 1.6], > [ 1.9]]) > >>> lin[3:3:1j,1:1:1j] > array([[ 1.6]]) > >>> > > I am always a fan of people who test and design their methods against edge cases like these, so my hat is off to you. I would suggest putting together a pull request so that we can properly test the potential impact such a change could have. Thanks! Ben Root