Goodman, Alexander (398J-Affiliate)

June 24, 2013 4:30 PM

Just something I

thought I should add to this, but would the following code seem like a reasonable workaround? Basically, since map_coordinates() cannot handle mask arrays, I thought that perhaps passing in the actual mask itself to

interp() with order=3 may produce close to reasonable results, eg:

dataout = interp(datain, xin, yin, xout, yout, masked=True, order=3)

maskin = datain.mask.astype(int)

maskout = interp(datain.mask, xin, yin, xout, yout, order=3)

dataout.mask = dataout.mask | maskout

Thanks,

Alex

Alex: The basemap.interp docstring includes a note describing a trick I often use with masked arrays:

Note

If datain is a masked array and order=1 (bilinear interpolation) is

used, elements of dataout will be masked if any of the four surrounding

points in datain are masked. To avoid this, do the interpolation in two

passes, first with order=1 (producing dataout1), then with order=0

(producing dataout2). Then replace all the masked values in dataout1

with the corresponding elements in dataout2 (using numpy.where).

This effectively uses nearest neighbor interpolation if any of the

four surrounding points in datain are masked, and bilinear interpolation

otherwise.

I suppose the same trick might work with order=3, but I have never tried it.

-Jeff