Re: [CF-metadata] Rotated Pole projection possibly wrong
The Z, Y', X'' rotation scheme, intrinsic or extrinsic, sounds like an application of Euler (sounds like Oiler) angles [1,2]. Tom Kunicki Center for Integrated Data Analytics U.S. Geological Survey 8505 Research Way Middleton, WI 53562 [1] http://en.wikipedia.org/wiki/Euler_angles [2] http://mathworld.wolfram.com/EulerAngles.html On Mar 6, 2013, at 2:17 PM, Seth McGinnis mcgin...@ucar.edu wrote: So there are an infinite number of ways to decompose a 3-D orientation into three orthogonal rotations. Specifying the projection in terms of the final location of the north pole and an optional rotation around that point means that CF is actually agnostic about exactly how you get there. That said, given that it's being specified in that way, the decomposition that maps most straightforwardly to those parameters is a Z-Y-X rotation sequence using an extrinsic (space-fixed) coordinate system: First rotate by angle phi around the north pole in space-fixed coordinates. Then rotate by angle theta around 0N 90E in space-fixed coordinates. Finally, rotate by angle psi around 0N 0E in space-fixed coordinates. The location of the grid's north pole in the original lat-lon coordinate system is now at longitude = phi, colatitude = minus theta. (I'm not sure how psi relates to north_pole_grid_longitude mathematically; it defaults to zero and I've never seen it in use, so I think it's included mostly for completeness. I agree that we may well want to change it to something more obvious like an azimuth angle.) Switching to the GRIB2 convention where you use the south pole instead of the north should be just a matter of adding 180 degrees to phi and negating theta. Cheers, --Seth On Wed, 6 Mar 2013 18:29:47 + David Hassell d.c.hass...@reading.ac.uk wrote: Hello, It is a bit of mess! As I understand it, the full rotation described is a sequence of rotations about three different axes: Z = 90S - 90N Y' = 90W' - 90E' X = 0E - 180E where it is understood that the definitions of the N-S and W-E axes change after each rotation (hence the primes and double primes). Therefore the order in which they are done matters. I suspect that the usual and assumed order is Z, Y', X? From the GRIB-2 stuff John posted the north_pole_grid_longitude gives the rotation about the X axis. Something like John's Angle of rotation seems right to me. The apparent lack of consistency between the parameter names irks me. Perhaps one solution could be to: i) Add some text to the conventions state the order (Z, Y', X, say) and direction of rotations. ii) add three new, consistent, self-describing parameters (e.g.) angle_of_rotation_z angle_of_rotation_y angle_of_rotation_x iii) allow aliases for backwards compatibility grid_north_pole_longitude = angle_of_rotation_z grid_north_pole_latitude = angle_of_rotation_y north_pole_grid_longitude = angle_of_rotation_x Or would that just obfuscate things even more? All the best, David P.S. If you have a copy available, there are some nice descriptions in Coordinate Systems and Map Projections by D. H. Maling -- David Hassell National Centre for Atmospheric Science (NCAS) Department of Meteorology, University of Reading, Earley Gate, PO Box 243, Reading RG6 6BB, U.K. Tel : 0118 3785613 Fax : 0118 3788316 E-mail: d.c.hass...@reading.ac.uk ___ CF-metadata mailing list CF-metadata@cgd.ucar.edu http://mailman.cgd.ucar.edu/mailman/listinfo/cf-metadata ___ CF-metadata mailing list CF-metadata@cgd.ucar.edu http://mailman.cgd.ucar.edu/mailman/listinfo/cf-metadata ___ CF-metadata mailing list CF-metadata@cgd.ucar.edu http://mailman.cgd.ucar.edu/mailman/listinfo/cf-metadata
Re: [CF-metadata] Rotated Pole projection possibly wrong
Hi Heiko: So grid_north_pole_longitude = normalizeLongitude180(longitudeOfSouthernPoleInDegrees) grid_north_pole_latitude = -1 * longitudeOfSouthernPoleInDegrees ? where does one find documentation on proj4's ob_tran routine? thanks, John On 3/6/2013 8:03 AM, Heiko Klein wrote: Hi John, all our models use rotated_latitute_longitude. And we have a consistent way of translating it between grib, proj and netcdf-java, i.e. publicly available at: http://thredds.met.no/thredds/catalog/metno/proff4km/default/catalog.html grid_mapping_name: rotated_latitude_longitude grid_north_pole_longitude: 156.0 grid_north_pole_latitude: 23.5 earth_radius: 6367470.0 to proj: proj4: +proj=ob_tran +o_proj=longlat +lon_0=-24 +o_lat_p=23.5 +a=6367470.0 +no_defs from grib: lonRot = longitudeOfSouthernPoleInDegrees latRot = latitudeOfSouthernPoleInDegrees +proj=ob_tran +o_proj=longlat +lon_0= normalizeLongitude180(lonRot) +o_lat_p= (-1 * latRot); This is the set of parameters which have been useful. I'm not sure if they are named in a useful way. proj4's ob_tran even allows for a o_lon_p parameter, but this doesn't seem to be used by grib or CF. Best regards, Heiko On 2013-03-06 15:18, John Caron wrote: The Rotated Pole projection here: http://cf-pcmdi.llnl.gov/documents/cf-conventions/1.6/cf-conventions.html#appendix-grid-mappings Rotated pole grid_mapping_name = rotated_latitude_longitude Map parameters: * grid_north_pole_latitude * grid_north_pole_longitude * north_pole_grid_longitude - This parameter is option (default is 0). Map coordinates: The rotated latitude and longitude coordinates are identified by the standard_name attribute values grid_latitude and grid_longitude respectively. is problematic. It has no reference to Proj.4, and is not in Snyder. Its likely the same as the GRIB-2 GDS Template 3.1: Grid Definition Template 3.1: Rotated Latitude/longitude (or equidistant cylindrical, or Plate Carrée) Octet No. Contents 15-72 Same as Grid Definition Template 3.0 (see Note 1) 73-76 Latitude of the southern pole of projection 77-80 Longitude of the southern pole of projection 81-84 Angle of rotation of projection 85-nnList of number of points along each meridian or parallel (These octets are only present for quasi-regular grids as described in Note 3) Notes: (2) Three parameters define a general latitude/longitude coordinate system, formed by a general rotation of the sphere. One choice for these parameters is: (a) The geographic latitude in degrees of the southern pole of the coordinate system, θp for example. (b) The geographic longitude in degrees of the southern pole of the coordinate system, λp for example. (c) The angle of rotation in degrees about the new polar axis (measured clockwise when looking from the southern to the northern pole) of the coordinate system, assuming the new axis to have been obtained by first rotating the sphere through λp degrees about the geographic polar axis, and then rotating through (90 + θp) degrees so that the southern pole moved along the (previously rotated) Greenwich meridian. Except note that GRIB-2 uses southern pole of projection. Not sure how that affects the algorithm. The CF parameter north_pole_grid_longitude likely should be Angle of rotation. Anyone using this or have some insight on it ?? ___ CF-metadata mailing list CF-metadata@cgd.ucar.edu http://mailman.cgd.ucar.edu/mailman/listinfo/cf-metadata ___ CF-metadata mailing list CF-metadata@cgd.ucar.edu http://mailman.cgd.ucar.edu/mailman/listinfo/cf-metadata
Re: [CF-metadata] Rotated Pole projection possibly wrong
Hello, It is a bit of mess! As I understand it, the full rotation described is a sequence of rotations about three different axes: Z = 90S - 90N Y' = 90W' - 90E' X = 0E - 180E where it is understood that the definitions of the N-S and W-E axes change after each rotation (hence the primes and double primes). Therefore the order in which they are done matters. I suspect that the usual and assumed order is Z, Y', X? From the GRIB-2 stuff John posted the north_pole_grid_longitude gives the rotation about the X axis. Something like John's Angle of rotation seems right to me. The apparent lack of consistency between the parameter names irks me. Perhaps one solution could be to: i) Add some text to the conventions state the order (Z, Y', X, say) and direction of rotations. ii) add three new, consistent, self-describing parameters (e.g.) angle_of_rotation_z angle_of_rotation_y angle_of_rotation_x iii) allow aliases for backwards compatibility grid_north_pole_longitude = angle_of_rotation_z grid_north_pole_latitude = angle_of_rotation_y north_pole_grid_longitude = angle_of_rotation_x Or would that just obfuscate things even more? All the best, David P.S. If you have a copy available, there are some nice descriptions in Coordinate Systems and Map Projections by D. H. Maling -- David Hassell National Centre for Atmospheric Science (NCAS) Department of Meteorology, University of Reading, Earley Gate, PO Box 243, Reading RG6 6BB, U.K. Tel : 0118 3785613 Fax : 0118 3788316 E-mail: d.c.hass...@reading.ac.uk ___ CF-metadata mailing list CF-metadata@cgd.ucar.edu http://mailman.cgd.ucar.edu/mailman/listinfo/cf-metadata
Re: [CF-metadata] Rotated Pole projection possibly wrong
Hi Burkhardt: Of course we cant change the existing conventions as there are files that already use it. But we need to document it so that we know when we are doing it right. Im cleaning up my code and noticed that i have two implementations, one for grib and one for CF. So Id like to merge the two. John On 3/6/2013 12:32 PM, Burkhardt Rockel wrote: John, the rotated pole definition was requested by me in 2002 (see the thread projections/default cell_methods). In this thread the present description was chosen. Please note: the rotated pole grid mapping is NOT a PROJECTION! it is a TRANSFORMATION (therefore the GRIB description is misleading). This is probably the reason why you did not find it in Proj4. The a, b, c in the GRIB definition description is fine. Whether the rotated South Pole (grib) or North Pole (CF-Conventions) is used in the definition does not matter, the values can be easily converted to each other. The rotated pole grid mapping as defined in the CF-Conventions document is already used widely as netCDF output standard for regional climate models in several international projects. It is also the standard in the IPCC-AR5 CORDEX model output for regional climate models. Therefore I am not in favor of changing the definition. However, a more precise description in the CF-Conventions would certainly be helpful for those not familiar with the rotated pole transformation. Regards Burkhardt --- Dr. Burkhardt Rockel Helmholtz-Zentrum Geesthacht Institute of Coastal Research / Group Regional Atmospheric Modeling Max-Planck-Strasse 1 D-21502 Geesthacht Germany Phone: +49 4152 87 1803 Fax: +49 4152 87 4 1803 Email: Burkhardt.Rockel (at) hzg.de www: http://rockel.staff.coast.hzg.de coordinates: 53.40579 N, 10.428647 E --- Am 06.03.2013 um 17:09 schrieb John Caron ca...@unidata.ucar.edu: Hi Heiko: So grid_north_pole_longitude = normalizeLongitude180(longitudeOfSouthernPoleInDegrees) grid_north_pole_latitude = -1 * longitudeOfSouthernPoleInDegrees ? where does one find documentation on proj4's ob_tran routine? thanks, John On 3/6/2013 8:03 AM, Heiko Klein wrote: Hi John, all our models use rotated_latitute_longitude. And we have a consistent way of translating it between grib, proj and netcdf-java, i.e. publicly available at: http://thredds.met.no/thredds/catalog/metno/proff4km/default/catalog.html grid_mapping_name: rotated_latitude_longitude grid_north_pole_longitude: 156.0 grid_north_pole_latitude: 23.5 earth_radius: 6367470.0 to proj: proj4: +proj=ob_tran +o_proj=longlat +lon_0=-24 +o_lat_p=23.5 +a=6367470.0 +no_defs from grib: lonRot = longitudeOfSouthernPoleInDegrees latRot = latitudeOfSouthernPoleInDegrees +proj=ob_tran +o_proj=longlat +lon_0= normalizeLongitude180(lonRot) +o_lat_p= (-1 * latRot); This is the set of parameters which have been useful. I'm not sure if they are named in a useful way. proj4's ob_tran even allows for a o_lon_p parameter, but this doesn't seem to be used by grib or CF. Best regards, Heiko On 2013-03-06 15:18, John Caron wrote: The Rotated Pole projection here: http://cf-pcmdi.llnl.gov/documents/cf-conventions/1.6/cf-conventions.html#appendix-grid-mappings Rotated pole grid_mapping_name = rotated_latitude_longitude Map parameters: * grid_north_pole_latitude * grid_north_pole_longitude * north_pole_grid_longitude - This parameter is option (default is 0). Map coordinates: The rotated latitude and longitude coordinates are identified by the standard_name attribute values grid_latitude and grid_longitude respectively. is problematic. It has no reference to Proj.4, and is not in Snyder. Its likely the same as the GRIB-2 GDS Template 3.1: Grid Definition Template 3.1: Rotated Latitude/longitude (or equidistant cylindrical, or Plate Carrée) Octet No. Contents 15-72 Same as Grid Definition Template 3.0 (see Note 1) 73-76 Latitude of the southern pole of projection 77-80 Longitude of the southern pole of projection 81-84 Angle of rotation of projection 85-nnList of number of points along each meridian or parallel (These octets are only present for quasi-regular grids as described in Note 3) Notes: (2) Three parameters define a general latitude/longitude coordinate system, formed by a general rotation of the sphere. One choice for these parameters is: (a) The geographic latitude in degrees of the southern pole of the coordinate system, θp for example. (b) The geographic longitude in degrees of the southern pole of the coordinate system, λp for example. (c) The angle of rotation in degrees about the new polar axis (measured clockwise when looking from the southern to the northern pole) of the coordinate system, assuming the new axis to have been obtained by first rotating the sphere
Re: [CF-metadata] Rotated Pole projection possibly wrong
Hi Seth: Thanks for the description, its really helpful. Im just wondering how we document this in the CF convention, so that implementers have something to check against. Is there a reference implementation we can use ? John On 3/6/2013 1:17 PM, Seth McGinnis wrote: So there are an infinite number of ways to decompose a 3-D orientation into three orthogonal rotations. Specifying the projection in terms of the final location of the north pole and an optional rotation around that point means that CF is actually agnostic about exactly how you get there. That said, given that it's being specified in that way, the decomposition that maps most straightforwardly to those parameters is a Z-Y-X rotation sequence using an extrinsic (space-fixed) coordinate system: First rotate by angle phi around the north pole in space-fixed coordinates. Then rotate by angle theta around 0N 90E in space-fixed coordinates. Finally, rotate by angle psi around 0N 0E in space-fixed coordinates. The location of the grid's north pole in the original lat-lon coordinate system is now at longitude = phi, colatitude = minus theta. (I'm not sure how psi relates to north_pole_grid_longitude mathematically; it defaults to zero and I've never seen it in use, so I think it's included mostly for completeness. I agree that we may well want to change it to something more obvious like an azimuth angle.) Switching to the GRIB2 convention where you use the south pole instead of the north should be just a matter of adding 180 degrees to phi and negating theta. Cheers, --Seth On Wed, 6 Mar 2013 18:29:47 + David Hassell d.c.hass...@reading.ac.uk wrote: Hello, It is a bit of mess! As I understand it, the full rotation described is a sequence of rotations about three different axes: Z = 90S - 90N Y' = 90W' - 90E' X = 0E - 180E where it is understood that the definitions of the N-S and W-E axes change after each rotation (hence the primes and double primes). Therefore the order in which they are done matters. I suspect that the usual and assumed order is Z, Y', X? From the GRIB-2 stuff John posted the north_pole_grid_longitude gives the rotation about the X axis. Something like John's Angle of rotation seems right to me. The apparent lack of consistency between the parameter names irks me. Perhaps one solution could be to: i) Add some text to the conventions state the order (Z, Y', X, say) and direction of rotations. ii) add three new, consistent, self-describing parameters (e.g.) angle_of_rotation_z angle_of_rotation_y angle_of_rotation_x iii) allow aliases for backwards compatibility grid_north_pole_longitude = angle_of_rotation_z grid_north_pole_latitude = angle_of_rotation_y north_pole_grid_longitude = angle_of_rotation_x Or would that just obfuscate things even more? All the best, David P.S. If you have a copy available, there are some nice descriptions in Coordinate Systems and Map Projections by D. H. Maling -- David Hassell National Centre for Atmospheric Science (NCAS) Department of Meteorology, University of Reading, Earley Gate, PO Box 243, Reading RG6 6BB, U.K. Tel : 0118 3785613 Fax : 0118 3788316 E-mail: d.c.hass...@reading.ac.uk ___ CF-metadata mailing list CF-metadata@cgd.ucar.edu http://mailman.cgd.ucar.edu/mailman/listinfo/cf-metadata ___ CF-metadata mailing list CF-metadata@cgd.ucar.edu http://mailman.cgd.ucar.edu/mailman/listinfo/cf-metadata
Re: [CF-metadata] Rotated Pole projection possibly wrong
Hi John, grid_north_pole_latitude = -1 * longitudeOfSouthernPoleInDegrees is ok The normalize function is just used because grib is 0-360, while proj is -180 to 180. But you need to flip longitude when going from north to south: grid_north_pole_longitude = norm180(longitudeOfSouthernPoleInDegrees + 180) It's documented in the proj documentation ftp://ftp.remotesensing.org/proj/proj.4.3.I2.pdf Chapter 4 on page 19. Heiko On 2013-03-06 17:09, John Caron wrote: Hi Heiko: So grid_north_pole_longitude = normalizeLongitude180(longitudeOfSouthernPoleInDegrees) grid_north_pole_latitude = -1 * longitudeOfSouthernPoleInDegrees ? where does one find documentation on proj4's ob_tran routine? thanks, John On 3/6/2013 8:03 AM, Heiko Klein wrote: Hi John, all our models use rotated_latitute_longitude. And we have a consistent way of translating it between grib, proj and netcdf-java, i.e. publicly available at: http://thredds.met.no/thredds/catalog/metno/proff4km/default/catalog.html grid_mapping_name: rotated_latitude_longitude grid_north_pole_longitude: 156.0 grid_north_pole_latitude: 23.5 earth_radius: 6367470.0 to proj: proj4: +proj=ob_tran +o_proj=longlat +lon_0=-24 +o_lat_p=23.5 +a=6367470.0 +no_defs from grib: lonRot = longitudeOfSouthernPoleInDegrees latRot = latitudeOfSouthernPoleInDegrees +proj=ob_tran +o_proj=longlat +lon_0= normalizeLongitude180(lonRot) +o_lat_p= (-1 * latRot); This is the set of parameters which have been useful. I'm not sure if they are named in a useful way. proj4's ob_tran even allows for a o_lon_p parameter, but this doesn't seem to be used by grib or CF. Best regards, Heiko On 2013-03-06 15:18, John Caron wrote: The Rotated Pole projection here: http://cf-pcmdi.llnl.gov/documents/cf-conventions/1.6/cf-conventions.html#appendix-grid-mappings Rotated pole grid_mapping_name = rotated_latitude_longitude Map parameters: * grid_north_pole_latitude * grid_north_pole_longitude * north_pole_grid_longitude - This parameter is option (default is 0). Map coordinates: The rotated latitude and longitude coordinates are identified by the standard_name attribute values grid_latitude and grid_longitude respectively. is problematic. It has no reference to Proj.4, and is not in Snyder. Its likely the same as the GRIB-2 GDS Template 3.1: Grid Definition Template 3.1: Rotated Latitude/longitude (or equidistant cylindrical, or Plate Carrée) Octet No. Contents 15-72 Same as Grid Definition Template 3.0 (see Note 1) 73-76 Latitude of the southern pole of projection 77-80 Longitude of the southern pole of projection 81-84 Angle of rotation of projection 85-nnList of number of points along each meridian or parallel (These octets are only present for quasi-regular grids as described in Note 3) Notes: (2) Three parameters define a general latitude/longitude coordinate system, formed by a general rotation of the sphere. One choice for these parameters is: (a) The geographic latitude in degrees of the southern pole of the coordinate system, θp for example. (b) The geographic longitude in degrees of the southern pole of the coordinate system, λp for example. (c) The angle of rotation in degrees about the new polar axis (measured clockwise when looking from the southern to the northern pole) of the coordinate system, assuming the new axis to have been obtained by first rotating the sphere through λp degrees about the geographic polar axis, and then rotating through (90 + θp) degrees so that the southern pole moved along the (previously rotated) Greenwich meridian. Except note that GRIB-2 uses southern pole of projection. Not sure how that affects the algorithm. The CF parameter north_pole_grid_longitude likely should be Angle of rotation. Anyone using this or have some insight on it ?? ___ CF-metadata mailing list CF-metadata@cgd.ucar.edu http://mailman.cgd.ucar.edu/mailman/listinfo/cf-metadata -- Dr. Heiko Klein Tel. + 47 22 96 32 58 Development Section / IT Department Fax. + 47 22 69 63 55 Norwegian Meteorological Institute http://www.met.no P.O. Box 43 Blindern 0313 Oslo NORWAY ___ CF-metadata mailing list CF-metadata@cgd.ucar.edu http://mailman.cgd.ucar.edu/mailman/listinfo/cf-metadata