Difference between revisions of "Python/PlotMetData"

From mn/geo/geoit
Jump to: navigation, search
(Plotting Meteorological Data with Python)
Line 16: Line 16:
 
<source lang='py'>
 
<source lang='py'>
  
    import matplotlib.pyplot as plt
+
import matplotlib.pyplot as plt
    import pygrib
+
import pygrib
    from netCDF4 import Dataset as NetCDFFile
+
from netCDF4 import Dataset as NetCDFFile
 
</source>
 
</source>
  
Line 24: Line 24:
 
Now, we'll write some python code to plot the netcdf variables::
 
Now, we'll write some python code to plot the netcdf variables::
 
<source lang='py'>
 
<source lang='py'>
    nci = NetCDFFile('geopot_april_01.nc')
+
nci = NetCDFFile('geopot_april_01.nc')
    nci.variables #will print a listing of the variables (note: nci.variables is a dictionary)
+
nci.variables #will print a listing of the variables (note: nci.variables is a dictionary)
    z = nci.variables['z']
+
z = nci.variables['z']
    z.shape  #what shape is the 'z' array?
+
z.shape  #what shape is the 'z' array?
  
    for i in range(3):
+
for i in range(3):
        plt.figure()
+
    plt.figure()
        plt.imshow(z[0,i,:,:]) #we take the i'th slab
+
    plt.imshow(z[0,i,:,:]) #we take the i'th slab
        plt.title('netcdf geopot {0}'.format(i))
+
    plt.title('netcdf geopot {0}'.format(i))
        plt.colorbar()
+
    plt.colorbar()
        plt.savefig('nc_geopot_{}'.format(i))
+
    plt.savefig('nc_geopot_{}'.format(i))
 
      
 
      
 
</source>
 
</source>
 
Now, let's take a look at the grib data::
 
Now, let's take a look at the grib data::
 
<source lang='py'>
 
<source lang='py'>
    grbs = pygrib.open('geop_april_01_2011.grib')
+
grbs = pygrib.open('geop_april_01_2011.grib')
    for grb in grbs:
+
for grb in grbs:
        grb
+
    grb
 
</source>
 
</source>
 
This will print a listing of all the gribs in the grib file. Make sure we reset the grbs back to the beginning::
 
This will print a listing of all the gribs in the grib file. Make sure we reset the grbs back to the beginning::
 
<source lang='py'>
 
<source lang='py'>
    grbs.seek(0)
+
grbs.seek(0)
    ghs = grbs.select(name='Geopotential')
+
ghs = grbs.select(name='Geopotential')
    for i,grb in enumerate(ghs):
+
 
        plt.figure()
+
for i,grb in enumerate(ghs):
        g = grb
+
    plt.figure()
        plt.imshow(g.values)
+
    g = grb
        title = g.name + ' at level ' + str(g.level) + ' [' + g.units + ']'
+
    plt.imshow(g.values)
        plt.title(title)
+
    title = g.name + ' at level ' + str(g.level) + ' [' + g.units + ']'
        plt.colorbar()
+
    plt.title(title)
        plt.savefig(g.name + str(i) +  '.png')
+
    plt.colorbar()
 +
    plt.savefig(g.name + str(i) +  '.png')
 
</source>
 
</source>
 
That's it! In your folder you should have some test plots.
 
That's it! In your folder you should have some test plots.

Revision as of 14:03, 14 September 2011

Plotting Meteorological Data with Python

This is a very brief introduction to using pygrib and NetCDFFile in order to plot met data in Python. Note that you should read the `Setting Python Paths<PythonPaths>` section in order to be sure that you'll be able to import the appropriate modules.

Let's get started with some test data::

   mkdir test
   cd test
   cp ~sec/kleinproject/FUKU/geop_april* .
   

Now, start ipython in your shell, or create a python file in the test directory::

   ipython
   

Next, let's take care of some basic imports::

import matplotlib.pyplot as plt
import pygrib
from netCDF4 import Dataset as NetCDFFile


Now, we'll write some python code to plot the netcdf variables::

nci = NetCDFFile('geopot_april_01.nc')
nci.variables #will print a listing of the variables (note: nci.variables is a dictionary)
z = nci.variables['z']
z.shape  #what shape is the 'z' array?

for i in range(3):
    plt.figure()
    plt.imshow(z[0,i,:,:]) #we take the i'th slab
    plt.title('netcdf geopot {0}'.format(i))
    plt.colorbar()
    plt.savefig('nc_geopot_{}'.format(i))

Now, let's take a look at the grib data::

grbs = pygrib.open('geop_april_01_2011.grib')
for grb in grbs:
    grb

This will print a listing of all the gribs in the grib file. Make sure we reset the grbs back to the beginning::

grbs.seek(0)
ghs = grbs.select(name='Geopotential')

for i,grb in enumerate(ghs):
    plt.figure()
    g = grb
    plt.imshow(g.values)
    title = g.name + ' at level ' + str(g.level) + ' [' + g.units + ']'
    plt.title(title)
    plt.colorbar()
    plt.savefig(g.name + str(i) +  '.png')

That's it! In your folder you should have some test plots.