[GRASS-user] Rotated pole coordinate system (a howto)

[Sebastian Schubert]

Hi all,

the following tries to describe how to use a rotated pole coordinate
system in grass (developed with the help of Hermann Pfeiffer). Such a
system is common in regional climate and weather prediction systems.
Here mainly because of numerical issues, wgs84 coordinates, which are
assumed to lie on a sphere, are modified such that the north pole is put
on another place on the sphere (to have the equator in the region of
interest). The new coordinates of the pole are called o_lat_p and
o_lon_p here.

The proj4 library handles that case with the "ob_tran" projection but
with opposite direction! To get normal lat lon coordinated from rotated
ones, use

$ proj -v -f "%.6f" -m 57.2957795130823 +proj=ob_tran +o_proj=latlon
+o_lon_p=LON_POLE +o_lat_p=LAT_POLE +lon_0=180

and for the opposite direction

$ invproj -v -f "%.6f" -m 57.2957795130823 +proj=ob_tran +o_proj=latlon
+o_lon_p=LON_POLE +o_lat_p=LAT_POLE +lon_0=180

The factor 57.2957795130823 = 180/pi ensures the output of degree
instead of radient.

OGR does not support ob_tran:
http://trac.osgeo.org/gdal/ticket/4285

Furthermore, I was not able to create a ob_tran location with the grass
gui, I always got proj4 parse errors! Because of that and the opposite
direction of the transformation, the transformation in grass is a bit
tricky.

Unrotated to rotated system
---------------------------

* read in the data in a wgs84 latlon location

* change the PERMANENT/PROJ_INFO of that location to

name: General Oblique Transformation
datum: wgs84
towgs84: 0.000,0.000,0.000
proj: ob_tran
o_proj: latlon
ellps: wgs84
a: 6378137.0000000000
es: 0.0066943800
f: 298.2572235630
lat_0: 0.0000000000
lon_0: 180.0000000000
o_lat_p: LAT_POLE
o_lon_p: LON_POLE

and PERMANENT/PROJ_UNITS to

unit: degree
units: degrees
meters: .0174532925

with .0174532925 = pi/180. v.in.ogr did not work for me after changing
to that projection, I guess because of ogr does not support it.

* create another wgs84 latlon location (this will be the rotated system)

* project the data from first system into the second


Rotated to unrotated system
---------------------------

* create a wgs84 latlon system and read in the data with coordinate
system check disables

* create another location with the PERMANENT/PROJ_INFO and
PERMANENT/PROJ_UNITS from above

* project from wgs84 to ob_tran system


HTH,
Sebastian

-------------- next part --------------
A non-text attachment was scrubbed...
Name: signature.asc
Type: application/pgp-signature
Size: 262 bytes
Desc: OpenPGP digital signature
Url : http://lists.osgeo.org/pipermail/grass-user/attachments/20111013/0ea7aa18/signature-0001.bin