|
|
| Line 25: |
Line 25: |
|
| |
|
| == Python == | | == Python == |
| To use python scripts to visualize GEOS-5 cubed-sphere products :
| | For the new format of the output, users can use the python codes here to transform the data to lat-lon grid and visualize the products: |
| #Install python 3 (or above) and the modules numpy, netcdf4, matplotlib and mpl_toolkits on your computer. | | #Install python 3 (or above) and the modules numpy, netcdf4, matplotlib and mpl_toolkits on your computer. |
| #Download source codes CSnative.py, example1.py, example2.py and the example data file HU_getc180forweyium_c180.geosgcm_prog.20000414_2200z.nc4 to the same folder. | | #Download source codes CSnative.py, example1.py, example2.py and the example data file HU_getc180forweyium_c180.geosgcm_prog.20000414_2200z.nc4 to the same folder. |
| Line 31: |
Line 31: |
| #To show example 2, $python example2.py | | #To show example 2, $python example2.py |
|
| |
|
| | ==== [[Recipe: python program reads cubed-sphere data| CSnative.py]]==== |
| | ==== [[Recipe: python program example1 | example1.py]] ==== |
| | ==== [[Recipe: python program example2| example2.py]] ==== |
|
| |
|
| For the new format of the output, users can use the python code here to transform the data to lat-lon grid and visualize it.
| |
|
| |
|
| <syntaxhighlight lang="python" line> | | <syntaxhighlight lang="python" > |
| import numpy as np
| |
| import matplotlib.pyplot as plt
| |
| from mpl_toolkits.axes_grid1 import make_axes_locatable
| |
| import CSgrid
| |
|
| |
|
| def main() :
| |
|
| |
| Time=0
| |
| Level=0
| |
| Nlat=500
| |
| Nlon=1000
| |
| filename = 'example_prog.20141118_2130z.nc4'
| |
| #filename = 'MAX_Discharge.geosgcm_surf.20141118_2130z.nc4'
| |
| var='T'
| |
| #var = 'PHIS'
| |
| # convert data trough interpolation
| |
| grid_T=CSgrid.Convert2LatLon(filename,var,Nlat,Nlon,Time,Level)
| |
| #grid_T=CSgrid.Convert2LatLon(filename,var,Nlat,Nlon,Time)
| |
| # color plot
| |
| plt.figure()
| |
| ax=plt.gca()
| |
| plt.xlabel('lon')
| |
| plt.ylabel('lat')
| |
| im=plt.imshow(grid_T.T,extent=(0,360,-90,90),origin='lower')
| |
| divider = make_axes_locatable(ax)
| |
| cax = divider.append_axes("right", size="5%", pad=0.05)
| |
| plt.colorbar(im, cax=cax)
| |
| plt.title(var)
| |
| plt.show()
| |
| if __name__ == '__main__' :
| |
| main()
| |
|
| |
| # save below as CSgrid.py
| |
|
| |
| import sys
| |
| from netCDF4 import Dataset
| |
| import numpy as np
| |
| from scipy.interpolate import griddata
| |
|
| |
| def Convert2LatLon(filename,var,Nlat,Nlon,Time,Level=None):
| |
| fh=Dataset(filename,mode='r')
| |
|
| |
| Times=fh.dimensions['time'].size
| |
| #Levs=fh.dimensions['lev'].size
| |
| nf=fh.dimensions['nf'].size
| |
| Ydim=fh.dimensions['y'].size
| |
| Xdim=fh.dimensions['x'].size
| |
| if Level is None:
| |
| T=fh.variables[var][Time,:,:,:]
| |
| else :
| |
| T = fh.variables[var][Time,Level, :, :, :]
| |
| lon3d=fh.variables['lons'][:,:,:]
| |
| lat3d=fh.variables['lats'][:,:,:]
| |
| fh.close()
| |
|
| |
| lat3d=np.reshape(lat3d,(nf*Ydim*Xdim,))
| |
| lon3d=np.reshape(lon3d,(nf*Ydim*Xdim,))
| |
| T=np.reshape(T,(nf*Ydim*Xdim,))
| |
| Tm=T
| |
| # get rid of non values
| |
| if(hasattr(T,'mask')):
| |
| Tm=T[~T.mask]
| |
| lat3d=lat3d[~T.mask]
| |
| lon3d=lon3d[~T.mask]
| |
|
| |
| lat=np.linspace(min(lat3d),max(lat3d),Nlat)
| |
| lon=np.linspace(min(lon3d),max(lon3d),Nlon)
| |
|
| |
| grid_x,grid_y=np.meshgrid(lat,lon)
| |
| points=(lat3d,lon3d)
| |
| grid_T=griddata(points,Tm,(grid_x,grid_y),'linear')
| |
|
| |
| # filled in extrapolated value with nearest
| |
|
| |
| lat3d =np.reshape(grid_x,(Nlon*Nlat,))
| |
| lon3d =np.reshape(grid_y,(Nlon*Nlat,))
| |
| Tm=np.reshape(grid_T,(Nlon*Nlat,))
| |
| lat3d=lat3d[~np.isnan(Tm)]
| |
| lon3d=lon3d[~np.isnan(Tm)]
| |
| Tm=Tm[~np.isnan(Tm)]
| |
| points=(lat3d,lon3d)
| |
|
| |
| grid_T=griddata(points,Tm,(grid_x,grid_y),'nearest')
| |
|
| |
| return grid_T
| |
| </syntaxhighlight> | | </syntaxhighlight> |
|
| |
|