Recipe: Python program to read data from downloaded file

From GEOS-5
Jump to: navigation, search

Back to G5NR Data Access Guide.

Problem

We want to read a downloaded data file using Python.

Solution

For the purpose of this example, we assume that we have already downloaded the file c1440_NR.inst01hr_3d_T_Cv.20060918_0900z.nc4 from the ftp server. For more information about file naming conventions, and how to download a file from the ftp server, please follow the links in the #See Also section.

Code (using netcdf4-python)

Here we use the netcdf4-python package to read the data file. If the command

> python -c "import netCDF4"

does not return any error, we have the package installed.

The code below reads the global temperature data from file, computes the maximum and minimum temperatures and plots (using the matplotlib package) the surface (level=71) temperature.

 1 #!/usr/bin/env python
 2 
 3 import sys
 4 import numpy as np
 5 import netCDF4 as nc4
 6 import matplotlib.pyplot as plt
 7 
 8 from mpl_toolkits.basemap import Basemap
 9 
10 rootgrp = nc4.Dataset('c1440_NR.inst01hr_3d_T_Cv.20060918_0900z.nc4', 'r')
11 print "rootgrp.variables['T'].shape", rootgrp.variables['T'].shape
12 
13 # read global air temperature for all levels
14 print 'Reading T...',; sys.stdout.flush()
15 T = rootgrp.variables['T'][0,:,:,:]
16 print 'done.'; sys.stdout.flush()
17 print 'T.shape:', T.shape
18 
19 # min/max
20 print 'min(T): %.4f' % np.min(T)
21 print 'max(T): %.4f' % np.max(T)
22 
23 # set up cylindrical map
24 m = Basemap(
25     projection='cyl',
26     llcrnrlat=-90, urcrnrlat=90,
27     llcrnrlon=-180, urcrnrlon=180,
28     resolution='c'
29     )
30 m.drawcoastlines(linewidth=0.5)
31 m.drawmapboundary()
32 
33 # plot contour
34 level = 71 # surface
35 X = np.arange(-180.0, 180.0, .5)
36 Y = np.arange(-90.0, 90.1, .5) # 90 is the last element
37 cp = plt.contour(X, Y, T[level,:,:], 20, zorder=2)
38 plt.clabel(cp, inline=1, fontsize=9)
39 plt.title('Air temperature at the surface')
40 plt.show()

Output

Running this python script, we get the text output

rootgrp.variables['T'].shape (1, 72, 361, 720)
Reading T... done.
T.shape: (72, 361, 720)
min(T): 180.3667
max(T): 315.6512

and the plot

SurfaceTemp.png

Discussions

Modifications to read a subset of the data

The above python code can be easily modified to read a subset of the data instead. For example, to read temperature data inside the box bounded by latitudes 25oN, 50oN and longitudes -130oW, -65oW, we will replace line 15 in the code above by

T = rootgrp.variables['T'][0,:,229:280,99:230]

The conversion of lat/lon co-ordinates to indices is left as an exercise for the user.

Running the modified script, produces the text output

rootgrp.variables['T'].shape (1, 72, 361, 720)
Reading T... done.
T.shape: (72, 51, 131)
min(T): 191.6956
max(T): 305.6048

See Also

  1. File Spec: File:G5NR-Ganymed-7km FileSpec-ON6-V1.0.pdf
  2. Recipe: Retrieve (global) data from FTP server

No Warranty

Copyright