sundar_ima - 4 months ago 124x

Python Question

I have a netcdf file with high resolution grided data (.1 x .1 deg) and plotting on the basemap using matplotlib. The plot I am trying to do is contour lines. For a specific reason I would like to plot the data at the interval of 1 x 1 deg resolution. For which I have used the following code sample example from here Regridding regular netcdf data.

See the update 1 for link to actual data.

For the sake of clearity following is the code I am trying to execute for regridding to lower resolution:-

`from mpl_toolkits import basemap`

from netCDF4 import Dataset

filename = 'path/to/netcdf/file.nc'

with Dataset(filename, mode='r') as fh:

lons = fh.variables['LON'][:]

lats = fh.variables['LAT'][:]

data = fh.variables['Data'][:].squeeze()

lons_sub, lats_sub = np.meshgrid(lons[::4], lats[::4], sparse=True)

data_coarse = basemap.interp(data, lons, lats, lons_sub, lats_sub, order=1)

The code seems to be correct. But when I execute the code I get the following error in the line

`data_coarse = basemap.interp(data, lons, lats, lons_sub, lats_sub, order=3)`

`/__init__.py", line 4897, in interp`

if xin[-1]-xin[0] < 0 or yin[-1]-yin[0] < 0:

ValueError: The truth value of an array with more than one element is ambiguous. Use a.any() or a.all()

I have not understood where the problem lies and how to solve this?

Any help is appreciated.

link to the actual data https://www.dropbox.com/s/ddlpbw5vvj5kz5w/mslp.txt?dl=0

link to lats https://www.dropbox.com/s/lpwjavxtwtt3r13/xlat.txt?dl=0

and the link to lons https://www.dropbox.com/s/1a0q49drfcd2o9h/xlon.txt?dl=0

Answer

In your code snippet, `xin = lons`

and `yin = lats`

. From the `basemap.interp`

docs,

xin, yin: rank-1 arrays containing x and y of datain grid in increasing order.

meaning `xin`

(`lons`

) and `yin`

(`lats`

) have to be continuously increasing. From your error traceback, it looks like that is not the case since either `xin[-1]-xin[0] < 0`

and/or `yin[-1]-yin[0] < 0`

.

Without knowing what `lons`

or `lats`

are exactly (since your code does not run with the data you link to), it is tough to elaborate any further.

I think the main problem you were running into was the shape of your latitude and longitude data. Since `lons`

and `lats`

were essentially stacks of the same arrays, I converted them both to 1d arrays as you'll see below.

```
import numpy as np
from mpl_toolkits import basemap
import matplotlib.pyplot as plt
# load data
data = np.loadtxt('mslp.txt', delimiter=',')
lons = np.loadtxt('xlon.txt', delimiter=',')
lats = np.loadtxt('xlat.txt', delimiter=',')
# take 1 dimensional slice of lons and lats
longitude = lons[0] # array([31.1, 31.4, 31.7,...,119.1, 119.3, 119.6])
latitude = lats[:,0] # array([-2.3, -2.1, -1.8,...,39.4, 39.6, 39.8])
# subdivide longitude and latitude into ~1deg grid data
lons_sub, lats_sub = np.meshgrid(longitude[::4], latitude[::4])
# implement basemap.interp to interpolate along ~1deg grid data
data_course = basemap.interp(datain=data, xin=longitude, yin=latitude,
xout=lons_sub, yout=lats_sub, order=3)
m = basemap.Basemap(llcrnrlon=longitude[0], llcrnrlat=latitude[0],
urcrnrlon=longitude[-1], urcrnrlat=latitude[-1],
projection='merc', resolution='c')
fig = plt.figure()
# ax1 -- contour plot of resampled data
ax1 = fig.add_subplot(211)
m.contour(lons_sub, lats_sub, data_course, linewidths=2.5, latlon=True)
m.fillcontinents(color=[0.7, 0.9, 0.8], lake_color=[0.3, 0.7, 0.9])
m.drawmapboundary(fill_color=[0.3, 0.7, 0.9])
m.drawcoastlines()
# ax2 -- imshow plot of resampled data
ax2 = fig.add_subplot(212)
m.drawcoastlines()
m.imshow(data_course, interpolation='nearest',
extent=[lons_sub[0], lons_sub[-1], lats[0], lats[-1]])
```

Source (Stackoverflow)

Comments