#======================================================================= # Additional Examples for Testing Function curl_2d in Module curl_2d. # # Copyright (c) 2004 by Johnny Lin. For licensing, distribution # conditions, contact information, and additional documentation see # the URL http://www.johnny-lin.com/py_pkgs/gemath/doc/. #======================================================================= # --------------------------- Import Modules --------------------------- try: execfile('/home/jlin/.pythonrc.py') except: print "test_curl_2d_add.py: Ignore Johnny Lin's path settings." import cdms, os, MV import Numeric as N import gemath from gemath.curl_2d import curl_2d import IaGraph from IaGraph import * # --------------------- Examples From CCM3 Dataset --------------------- #- Data: From CCM3 output. u and v are in m/s and x and y are in # deg (except for the final examples, when we take x and y as being # in m): f = cdms.open('uv_501hPa_ccm3.nc') u = f('U') v = f('V') x = MV.filled( u.getAxis(-1)[:] ) y = MV.filled( u.getAxis(-2)[:] ) u = MV.filled( u[:,:] ) v = MV.filled( v[:,:] ) #- Set system variables for lon/lat axes tick labels: sysvar = IaGraph.Sysvar() sysvar.__class__.x_tickvalues = \ {0:'0', 90:'90E', 180:'180E', 270:'270E', 360:'360E'} sysvar.__class__.y_tickvalues = \ {-90:'90S', -45:'45S', 0:'EQ', 45:'45N', 90:'90N'} #- Add more missing values at the edges, use "default_spherical" # algorithm. The order 1 spherical algorithm will be chosen: u[30:41,20:37] = 1e+20 u[0:9,50:67] = 1e+20 curl = curl_2d(x, y, u, v, algorithm='default_spherical') base_fn = 'test_curl_2d_ccm3_defsph_edge' contour(MA.masked_values(curl, 1e+20), x, y, continents=1, title=base_fn) active2gif(base_fn+'.gif') os.system('convert -verbose -trim -border 8x8 -bordercolor white ' \ +'-mattecolor black -frame 2x2 ' \ +base_fn+'.gif gem_'+base_fn+'.jpg') #- Use order 1 Cartesian algorithm: If x and y are in meters, and u # and v are in m/s, the curl is in 1/s. (I reset x and y tickvalues # for this and the following plot): sysvar.__class__.x_tickvalues = sysvar.__class__.y_tickvalues = None u[30:41,20:37] = 1e+20 u[0:9,50:67] = 1e+20 curl = curl_2d(x, y, u, v, algorithm='order1_cartesian') base_fn = 'test_curl_2d_ccm3_o1cart' contour(MA.masked_values(curl, 1e+20), x, y, title=base_fn) active2gif(base_fn+'.gif') os.system('convert -verbose -trim -border 8x8 -bordercolor white ' \ +'-mattecolor black -frame 2x2 ' \ +base_fn+'.gif gem_'+base_fn+'.jpg') #- Use "default" algorithm, which for this domain chooses the order 1 # Cartesian algorithm on the dataset with above missing values: u[30:41,20:37] = 1e+20 u[0:9,50:67] = 1e+20 curl = curl_2d(x, y, u, v, algorithm='default') base_fn = 'test_curl_2d_ccm3_def' contour(MA.masked_values(curl, 1e+20), x, y, title=base_fn) active2gif(base_fn+'.gif') os.system('convert -verbose -trim -border 8x8 -bordercolor white ' \ +'-mattecolor black -frame 2x2 ' \ +base_fn+'.gif gem_'+base_fn+'.jpg') #- Close data file: f.close() # ------------------ Examples From Synthetic Dataset ------------------- #- Data: Zonal velocity u is non-zero, except poleward of 60N/S, in # which zonal velocity is zero. Likewise for meridional velocity v: nx = 180 ny = 91 x = N.arange(nx) * 2.0 - 180.0 y = N.arange(ny) * 2.0 - 90.0 y_2d = N.reshape( N.repeat(y, nx), (ny, nx) ) x_2d = N.reshape( N.repeat(N.reshape(x,(1,nx)), ny), (ny, nx) ) u = N.sin(x_2d*N.pi/180.)*N.sin((y_2d-90.)*N.pi/30.) v = N.sin((x_2d-30.)*N.pi/180.)*N.sin((y_2d-90.)*N.pi/30.) N.putmask( u, N.where(N.absolute(u) < 1e-10, 1, 0), 0. ) N.putmask( v, N.where(N.absolute(v) < 1e-10, 1, 0), 0. ) N.putmask( u, N.where(N.absolute(y_2d) > 60., 1, 0), 0. ) N.putmask( v, N.where(N.absolute(y_2d) > 60., 1, 0), 0. ) #- Set system variables for lon/lat axes tick labels, which all the # rest of the examples assume x and y are (in deg), respectively: sysvar = IaGraph.Sysvar() sysvar.__class__.x_tickvalues = \ {-180:'-180', -90:'-90', 0:'0', 90:'90', 180:'180'} sysvar.__class__.y_tickvalues = \ {-90:'-90', -45:'-45', 0:'0', 45:'45', 90:'90'} #- Plot u and v data: base_fn = 'test_curl_2d_syn_udata' contour(u, x, y, title=base_fn) active2gif(base_fn+'.gif') os.system('convert -verbose -trim -border 8x8 -bordercolor white ' \ +'-mattecolor black -frame 2x2 ' \ +base_fn+'.gif gem_'+base_fn+'.jpg') base_fn = 'test_curl_2d_syn_vdata' contour(v, x, y, title=base_fn) active2gif(base_fn+'.gif') os.system('convert -verbose -trim -border 8x8 -bordercolor white ' \ +'-mattecolor black -frame 2x2 ' \ +base_fn+'.gif gem_'+base_fn+'.jpg') #- Use order 1 spherical algorithm: x and y are in deg and u and v # are in m/s. The curl is in 1/s. Note that these results should # be nearly identical to those from SPHEREPACK. Assume we're on # the earth: curl = curl_2d(x, y, u, v, algorithm='order1_spherical') base_fn = 'test_curl_2d_syn_o1sph' contour(curl, x, y, title=base_fn) active2gif(base_fn+'.gif') os.system('convert -verbose -trim -border 8x8 -bordercolor white ' \ +'-mattecolor black -frame 2x2 ' \ +base_fn+'.gif gem_'+base_fn+'.jpg') #- Use spectral method from SPHEREPACK (need to first remove the # repeating dateline point from the domain). Here we take x and y in # deg, while u and v are in m/s. If so the curl is in 1/s. Assume # we're on the earth: curl = curl_2d(x, y, u, v, algorithm='spherepack') base_fn = 'test_curl_2d_syn_sphpk' contour(curl, x, y, title=base_fn) active2gif(base_fn+'.gif') os.system('convert -verbose -trim -border 8x8 -bordercolor white ' \ +'-mattecolor black -frame 2x2 ' \ +base_fn+'.gif gem_'+base_fn+'.jpg') #- Use order 1 spherical algorithm, but assume we're on Mars # (radius here is the volumetric mean radius): curl = curl_2d(x, y, u, v, algorithm='order1_spherical', R_sphere=3390e+3) base_fn = 'test_curl_2d_syn_o1sph_mars' contour(curl, x, y, title=base_fn) active2gif(base_fn+'.gif') os.system('convert -verbose -trim -border 8x8 -bordercolor white ' \ +'-mattecolor black -frame 2x2 ' \ +base_fn+'.gif gem_'+base_fn+'.jpg') #- Use SPHEREPACK, but assume we're on Mars (radius here is the # volumetric mean radius): curl = curl_2d(x, y, u, v, algorithm='spherepack', R_sphere=3390e+3) base_fn = 'test_curl_2d_syn_sphpk_mars' contour(curl, x, y, title=base_fn) active2gif(base_fn+'.gif') os.system('convert -verbose -trim -border 8x8 -bordercolor white ' \ +'-mattecolor black -frame 2x2 ' \ +base_fn+'.gif gem_'+base_fn+'.jpg') #- Use order 1 Cartesian: curl = curl_2d(x, y, u, v, algorithm='order1_cartesian') base_fn = 'test_curl_2d_syn_o1cart' contour(curl, x, y, title=base_fn) active2gif(base_fn+'.gif') os.system('convert -verbose -trim -border 8x8 -bordercolor white ' \ +'-mattecolor black -frame 2x2 ' \ +base_fn+'.gif gem_'+base_fn+'.jpg') #- Use order 1 spherical algorithm, but assume we're on Mars, and add # some missing values: u[30:46,90:114] = 1e+20 curl = curl_2d(x, y, u, v, algorithm='order1_spherical', R_sphere=3390e+3) base_fn = 'test_curl_2d_syn_o1sph_miss_mars' contour(MA.masked_values(curl, 1e+20), x, y, title=base_fn) active2gif(base_fn+'.gif') os.system('convert -verbose -trim -border 8x8 -bordercolor white ' \ +'-mattecolor black -frame 2x2 ' \ +base_fn+'.gif gem_'+base_fn+'.jpg') # ===== end file =====