[Science]   [BAS home]   [Met home]   [WMC home]	Antarctic Meteorology

Differential stuff on the sphere using SPHEREPACK

Note for non-pp users: all the IDL progs here assume the data is given wrapped up in a pp-field. For just about all purposes, all you need to know is that this is a structure with various components, including:

1. pp.data - nxm float array of data
2. pp.lbrow - number of rows (m; m=73 implies 2.5 step in latitude)
3. pp.lbnpt - number of columns (n; n=96 implies 3.75 step in longitude)
4. pp.bdx - step in x-direction (longitude)
5. pp.bdy - step in y-direction (latitude)

The main purpose of this was the Rossby Wave Source term stuff: RWS: see below. But the routines here would also be useful if you just want gradients, laplacians etc of PP fields on the sphere. If you use this, read the cautions...

SPHEREPACK page: http://www.scd.ucar.edu/css/software/spherepack/

NOTE (this caused me so much pain...) that spherepack uses a slightly different wind-component-sign-convention: the NS one is the negative of what I/UKMO normally use. Ie, the first component is easterly (as usual): a +ve wind is blowing from W to E; but the second component is such that a +ve wind is blowing from N to S. REMEMBER THIS! In practice, this means that all NS components need to be negated before being handed to spherepack, and then (if the result is a wind) the NS component negated on return.

Cautions, notes

1. My strategy was to take spherepack routines, write wrappers for the various functions, the IDL wrappers, then use. This means that these programs will only work on architectures for which the fortan was compiled - currently alpha (bscomp, bsaicde) or linux
2. These progs probably won't run direct on ERA data (they are too big, hmmm). What you need to do is regrid them smaller, probably to the UKMO grid, by:

   	new_field=pp_regrid(era_field,get_orog(),/pole)
   

3. Note: you must include the /pole switch. If you don't, then the regridded field will contain junk in the polar row and terrible things will happen.
4. Note: the IDL pros have names like fpoisson_solver, indicating that they are spawning to fortran.
5. Note: they write data files, *spawn* to fortran, then read in data files. This might be thought slow... but its not too bad. A few seconds.
6. Note: all this is on the UNIT SPHERE. To get into geophysical units, its often necessary to multiply/divide as appropriate by !met_values.radius_earth.
7. And finally... I've checked this stuff, see below. But be careful.

Laplacian

Poissons Equation

Stream function and velocity potential

My code: sfvp.f compiled by compile1_sfvp and then used by IDL prog fsfvp.

test_sfvp1. Uses a real ERA 850 hPA field. The result is:

which seems to show we've got it right. It also checks that div(rot_wind) is nearly zero (yes) and curl(div_wind) too (yes).

ALSO: since v_rot = k ^ grad(stream function), one can test this by taking grad (using fgrad.pro), then using kcross.pro, then computing fsfvp on the result, and you should be back to where you started. And you are.

Rossby Wave Source (rws)

See: rws for more.

Curl

My code: curl.f compiled by compile1_curl and then used by IDL prog fcurl.

Testing

Poisson solver

@test_poisson
z=get_field('test-1',nf=3)
gettwogifs,out='test-1'