netcdf4 extract for subset of lat lon

Solution 1:

Well this is pretty easy, you have to find the index for the upper and lower bound in latitude and longitude. You can do it by finding the value that is closest to the ones you're looking for.

latbounds = [ 40 , 43 ]
lonbounds = [ -96 , -89 ] # degrees east ? 
lats = f.variables['latitude'][:] 
lons = f.variables['longitude'][:]

# latitude lower and upper index
latli = np.argmin( np.abs( lats - latbounds[0] ) )
latui = np.argmin( np.abs( lats - latbounds[1] ) ) 

# longitude lower and upper index
lonli = np.argmin( np.abs( lons - lonbounds[0] ) )
lonui = np.argmin( np.abs( lons - lonbounds[1] ) )

Then just subset the variable array.

# Air (time, latitude, longitude) 
airSubset = f.variables['air'][ : , latli:latui , lonli:lonui ]

Note, i'm assuming the longitude dimension variable is in degrees east, and the air variable has time, latitude, longitude dimensions.

Solution 2:

Favo's answer works (I assume; haven't checked). A more direct and idiomatic way is to use numpy's where function to find the necessary indices.

lats = f.variables['latitude'][:] 
lons = f.variables['longitude'][:]
lat_bnds, lon_bnds = [40, 43], [-96, -89]

lat_inds = np.where((lats > lat_bnds[0]) & (lats < lat_bnds[1]))
lon_inds = np.where((lons > lon_bnds[0]) & (lons < lon_bnds[1]))

air_subset = f.variables['air'][:,lat_inds,lon_inds]

Solution 3:

If you like pandas, then you should think about checking out xarray.

import xarray as xr

ds = xr.open_dataset('http://geoport.whoi.edu/thredds/dodsC/usgs/data2/rsignell/models/ncep/narr/air.2m.1980.nc',
                     decode_cf=False)
lat_bnds, lon_bnds = [40, 43], [-96, -89]
ds.sel(lat=slice(*lat_bnds), lon=slice(*lon_bnds))

Solution 4:

Note that this can be accomplished even quicker on the command line using NCO's ncks.

ncks -v air -d latitude,40.,43. -d longitude,-89.,-96. infile.nc -O subset_infile.nc

Solution 5:

To mirror the response from N1B4, you can also do it on one line with climate data operators (cdo):

cdo sellonlatbox,-96.5,-89,40,43 in.nc out.nc

Thus to loop over a set of file, I would do this in a BASH script, using cdo to process each file and then calling your python script:

#!/bin/bash

# pick up a list of files (I'm presuming the loop is over the years)
files=`ls /usgs/data2/rsignell/models/ncep/narr/air.2m.*.nc`

for file in $files ; do 
   # extract the location, I haven't used your exact lat/lons
   cdo sellonlatbox,-96.5,-89,40,43 $file iowa.nc

   # Call your python or R script here to process file iowa.nc
   python script
done

I always try and do my file processing "offline" as I find it less prone to error. cdo is an alternative to ncks, I'm not saying it is better, I just find it easier to remember the commands. nco in general is more powerful and I resort to it when cdo can't perform the task I wish to carry out.