[SciPy-User] Memory error with scipy.interpolate

ashwinD12 . winash12 at gmail.com
Sat Jul 30 03:54:48 EDT 2016


I was able to get this to work in the following way but I would really
appreciate a code review of my approach.When I used linear I got a qHull
error. Then I swtiched to nearest and the code seems to run. But the
details of the why and what would be deeply appreciated.

grb = grbs.select(name='Geopotential',level=500)[0]
data = grb.values

hgt = [x/10 for x in data]

lat,lon = grb.latlons()

hgt = np.array(hgt)

lona,lata =
np.linspace(lon.min(),lon.max(),100),np.linspace(lat.min(),lat.max(),100)
hgt1 = np.linspace(hgt.min(),hgt.max(),100)

hgt1 = np.resize(hgt,(100,100))

longr,latgr = np.meshgrid(lona, lata)

x,y = m(longr,latgr)


zi =
scipy.interpolate.griddata((lona,lata),hgt1,(longr,latgr),method='nearest')


On Sat, Jul 30, 2016 at 12:59 PM, ashwinD12 . <winash12 at gmail.com> wrote:

> Hello David,
>                     Thanks for your response. I am not sure whether I was
> clear in the example I gave you. I gave gridded data.
>
> The data array is a 2D array and has shape 561 x 401. Over on github I was
> told by Pauli Virtanen my data array  is too big for Rbf, Since I have grid
> data I was asked to use a method corresponding to that.
>
> I tried using scipy.interpolate.griddata and I am just really confused how
> to proceed further. I have a 2-dimensional array of heights data and I want
> to make a contour plot with x axis longitude and y axis longitude. I chose
> the interpolation method griddata. Following is my code -
>
> Originally hgt is array of shape 561 x 401 and lat and lon have the same
> dimensions.
>
> grb = grbs.select(name='Geopotential',level=500)[0]
> data = grb.values
>
> hgt = [x/10 for x in data]
>
> lat,lon = grb.latlons()
>
> hgt = np.array(hgt)
>
>
> lona,lata =
> np.linspace(lon.min(),lon.max(),100),np.linspace(lat.min(),lat.max(),100)
> hgt1 = np.linspace(hgt.min(),hgt.max(),100)
>
> hgt1 = np.resize(hgt,(100,100))
>
> longr,latgr = np.meshgrid(lona, lata)
>
>
> print(x.shape,y.shape,hgt.shape)
>
> gd =
> scipy.interpolate.griddata((lon,lat),hgt,(longr,latgr),method='linear')
>
>
> When I do this I get this error message -
>
> Traceback (most recent call last):
>   File "contour.py", line 44, in <module>
>     rbf =
> scipy.interpolate.griddata((lon,lat),hgt,(longr,latgr),method='linear')
>   File
> "/usr/local/lib/python3.4/dist-packages/scipy/interpolate/ndgriddata.py",
> line 217, in griddata
>     rescale=rescale)
>   File "scipy/interpolate/interpnd.pyx", line 243, in
> scipy.interpolate.interpnd.LinearNDInterpolator.__init__
> (scipy/interpolate/interpnd.c:4934)
>   File "scipy/interpolate/interpnd.pyx", line 78, in
> scipy.interpolate.interpnd.NDInterpolatorBase.__init__
> (scipy/interpolate/interpnd.c:2386)
>   File "scipy/interpolate/interpnd.pyx", line 191, in
> scipy.interpolate.interpnd._check_init_shape
> (scipy/interpolate/interpnd.c:4593)
> ValueError: invalid shape for input data points
>
>
> On Wed, Jul 20, 2016 at 9:03 PM, ashwinD12 . <winash12 at gmail.com> wrote:
>
>> From this old thread -
>> https://mail.scipy.org/pipermail/scipy-user/2009-April/020607.html
>>
>> it appears that Scipy.Rbf has a limitation with number of points. If so
>> what is the alternative to doing an interpolation so I can do a contour
>> plot ?
>>
>>
>> On Wed, Jul 20, 2016 at 12:13 PM, ashwinD12 . <winash12 at gmail.com> wrote:
>>
>>> Hello,
>>>          I get a memory error when I interpolate a data set of 567
>>> points as shown in this stack trace
>>>
>>> Traceback (most recent call last):
>>>  File "contour.py", line 34, in <module>
>>>  rbf = scipy.interpolate.Rbf(x, y, data, function='linear')
>>>  File "/usr/local/lib/python3.4/dist-packages/scipy/interpolate/rbf.py", line 200, in __init__
>>>  r = self._call_norm(self.xi, self.xi)
>>>  File "/usr/local/lib/python3.4/dist-packages/scipy/interpolate/rbf.py", line 231, in _call_norm
>>>  return self.norm(x1, x2)
>>>  File "/usr/local/lib/python3.4/dist-packages/scipy/interpolate/rbf.py", line 118, in _euclidean_norm
>>>  return np.sqrt(((x1 - x2)**2).sum(axis=0))
>>>  MemoryError
>>>
>>> With numpy earlier I had gotten the same error and I had moved past it
>>>
>>> by this option -
>>> lon,lat = np.meshgrid(lon, lat,sparse=True,copy=False)
>>>
>>> This is the minimum code snippet that should reproduce the problem
>>> grbs = pygrib.open('20000')
>>>
>>>
>>> grb = grbs.select(name='Geopotential',level=500)[0]
>>> data = grb.values
>>>
>>> hgt = [x/10 for x in data]
>>> lat,lon = grb.latlons()
>>>
>>> lon,lat = np.meshgrid(lon, lat,sparse=True,copy=False)
>>>
>>> m=Basemap(projection='mill',lat_ts=10,llcrnrlon=lon.min(), \
>>> urcrnrlon=lon.max(),llcrnrlat=lat.min(),urcrnrlat=lat.max(), \
>>>       resolution=None)
>>>
>>>
>>> x,y = m(lon,lat)
>>>
>>> rbf = scipy.interpolate.Rbf(x, y, data, function='linear')
>>> zi = rbf(x, y)
>>>
>>> Regards,
>>>
>>> Ashwin.
>>>
>>>
>>>
>>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.scipy.org/pipermail/scipy-user/attachments/20160730/9fcc4719/attachment.html>


More information about the SciPy-User mailing list