NOTICE: This version of the NSF Unidata web site (archive.unidata.ucar.edu) is no longer being updated.
Current content can be found at unidata.ucar.edu.
To learn about what's going on, see About the Archive Site.
Tor,Ok, this make sense to me. I implemented the following and the projection worked correctly now. Do all Rotated Lat/Lon grids follow this scheme?
Thanks for your help, this was my first projection and I went off in the wrong direction some times.
Robb... On Mon, 29 Dec 2008, Tor Christian Bekkvik wrote:
Hi, Since the projection code is very likely correct, I guess the error is in interpretation of the reference coordinates and deltas given in the grib-files; - In the unrotated image, the grid usually cannot be expressed in linear first degree polynomials (longitude = lonRef + DX* i_index, latitude = latRef + DY*j_index). But in the rotated lon/lat reference system (where the south pole could be on equator), the grid is a regular, linear grid where parameters as lower left corner and grid deltas makes sense:Rotated LatLon Grid (468x378) lon: -46.5 to 46.9 (dx 9.99990041077413E24) lat: -36.5 to 38.9 (dy 9.99990041077413E24)I interpret the dx,dy parameters above as undefined since they are too large to give meaning (1 degree north is approx 111 km). Dividing delta lon by 467 (468-1 cells in transformed east-direction) gives 0.2 degrees grid resolution east. Similiar calculation gives 0.2 degrees north: | diff| N-1 | diff/(N-1) lon -46,5 46,9 | 93,4 | 467 | 0,2 deg lat -36,5 38,9 | 75,4 | 377 | 0,2 deg PS. These lon, lat coordinates are not real, geographical east/north coordinates, but angular positions relative to the rotated reference system. The inverse rotated latlon (lon,lat) projection will give the real geographical coordinates of any positions in the grid. Hope this helps.. Tor C ----- Original Message ----- From: Robb Kambic <address@hidden> Sent: Wed, 24-12-2008 01:03 To: Tor Christian Bekkvik <address@hidden> Cc: decoders <address@hidden> Subject: RE: [netcdf-java] Adding a new projection to GRIB netcdf decoder Tor, I changed the Lower Left corner of the unrotated image to be approximately 28 N and 47W and the image almost lines up with the map. Also it appears that the image needs to be rotated some. the rotated coordinates I used are: La1 = -20.25 Lo1 = -43.5 yours are: inv [-46.5, -36.5] -> [-37.030363026843844, 14.484143589394591] since we have verified that the projection code is the same, what's the next step. Robb... On Tue, 23 Dec 2008, Robb Kambic wrote:Tor, I was on vacation for a while, so no work was done on the projection. I tested my code implemetation against the code you sent, I got all the same values. Currently, my results are south and east of the desired map location. Here's my thoughts, the Dx and Dy values in the file are not correct, I get Dx = Dy = .2 degrees. These values are in millidegrees according to the documentation. The values that you have are grid_dx = 9.99990041077413E24 that I assume are meters. I think that if the Dx and Dy values were correct then the image would line up correctly. What do you think? RObb... C:\data\rotatedlatlon.grb index_version = 6.0 grid_edition = 1 location = C:\data\rotatedlatlon.grb length = 265454 created = 2008-12-23T21:44:09Z center = 96 sub_center = 0 table_version = 1 tiles = 1 thin = false ensemble = false -------------------------------------------------------------------- 1 0 -1 7 1 105 0.0 255 0.0 2003-12-15T18:00:00Z 0 2891248751 36 265454 0 false 96 0 1 -------------------------------------------------------------------- GDSkey = 2891248751 grid_type = 10 grid_name = Rotated latitude/longitude grid grid_shape_code = 0 grid_shape = spherical grid_radius_spherical_earth = 6367.47 Nx = 468 Ny = 378 La1 = -36.5 Lo1 = -46.5 ResCompFlag = 136 Winds = Relative La2 = 38.9 Lo2 = 46.9 Dx = 0.2 Dy = 0.2 ScanningMode = 64 SpLat = -25.0 SpLon = 0.0 Angle = 0.0 On Mon, 24 Nov 2008, Tor Christian Bekkvik wrote:Hi, I did som debugging on our JGRIB application(rotatedlatlon.grb), hope this helps: 1) Grib GDS Section content (JGRIB interpretation): EARTH_RADIUS = 6367470.0 length = 46 nv = 1 octetPos = 43 grid_type = 10 grid_type_string = {java.lang.String@616}"Rotated LatLon Grid" grid_nx = 468 grid_ny = 378 grid_lat1 = -36.5 grid_lon1 = -46.5 grid_lat2 = 38.9 grid_lon2 = 46.9 grid_mode = 136 grid_dx = 9.99990041077413E24 grid_dy = 9.99990041077413E24 grid_scan = 64 isUVEastNorth = false grid_latsp = -25.0 grid_lonsp = 0.0 grid_rotang = 0.0 2) Then, from this Grib GDS section, the projection was instantiated with: new RotLatLon(grid_lonsp, grid_latsp, grid_rotang) 3) And here is the resulting instantiated projection code: public RotLatLon(double aLonsp, double aLatsp, double aRotsp) { lonsp = aLonsp; latsp = aLatsp; rotsp = aRotsp; double dlat_rad = (latsp - (-90)) * DEG2RAD; //delta latitude sinDlat = Math.sin(dlat_rad); cosDlat = Math.cos(dlat_rad); } static final double DEG2RAD = Math.PI / 180; static final double RAD2DEG = 180 / Math.PI; public double[] fwd(double[] lonlat) { return transform(lonlat, lonsp, rotsp, sinDlat); } public double[] inv(double[] lonlat) { return transform(lonlat, -rotsp, -lonsp, -sinDlat); } private double[] transform(double [] lonlat, double rot1, double rot2, double s) { double e = DEG2RAD * (lonlat[0] - rot1); //east double n = DEG2RAD * lonlat[1]; //north double cn = Math.cos(n); double x = cn * Math.cos(e); double y = cn * Math.sin(e); double z = Math.sin(n); double x2 = cosDlat * x + s * z; double z2 = -s * x + cosDlat * z; double R = Math.sqrt(x2 * x2 + y * y); double e2 = Math.atan2(y, x2); double n2 = Math.atan2(z2, R); double rlon = RAD2DEG * e2 - rot2; double rlat = RAD2DEG * n2; return new double[]{rlon, rlat}; } Thanks, Tor C-----Original Message----- From: Robb Kambic [mailto:address@hidden] Sent: 24. november 2008 02:51 To: Tor Christian Bekkvik Subject: RE: [Fwd: [netcdf-java] FW: Adding a new projection to GRIB netcdf decoder] (fwd) Tor, thanks for the example file, it helped. I'm still having problems with the projection, it seems that it's -90 off. If I enter 0, 0 into the projection transform I would expect that the results would be the La1, Lo1 and that's not the case. Is this the formula that you used in creating the grid? Any help would be appreciated. I get all the parameters correct. THnaks, Robb... On Thu, 20 Nov 2008, Tor Christian Bekkvik wrote:Hi, I have been a little busy; Sorry for the bad example file I sent first, the attached file should have more meaningful values. Some file info, logged from a java-application using JGRIB (JGRIB supports rotated latlon for GRIB-1) : Gribfile: rotatedlatlon.grb Records: 1 Grid: GDS section: Rotated LatLon Grid (468x378) lon: -46.5 to 46.9 (dx 9.99990041077413E24) lat: -36.5 to 38.9 (dy 9.99990041077413E24) south pole: lon 0.0 lat -25.0 rot angle: 0.0 mode(17):0x88 type(6) :0xa scan(28):0x40 PDS cent tab proc: 96 1 1 -> 96 -1 1 ncep_opn.tab Z type(105): fixed height above ground Z #levels: 1 Dates(1): min Mon Dec 15 19:00:00 CET 2003 max Mon Dec 15 19:00:00 CET 2003 Param: 7 HGT [gpm] Geopotential height [-275.23438 4048.7656] PS. Hoping that rotated latlon projection could be supported for GRIB-2 Thanks, Tor C-----Original Message----- From: Robb Kambic [mailto:address@hidden] Sent: 20. november 2008 00:44 To: Tor Christian Bekkvik Subject: Re: [Fwd: [netcdf-java] FW: Adding a newprojection to GRIBnetcdf decoder] (fwd) Tor, An image of the result would be appreciated. Thanks, Robb... On Wed, 19 Nov 2008, Robb Kambic wrote:Tor, I've been working on this problem but the sample file yousent has allthe data values equal to 0. could you make availableanother samplefile with a description so I can tell if I am displayingcorrectly.Robb... -------- Original Message -------- Subject: [netcdf-java] FW: Adding a new projection to GRIB netcdf decoder Date: Wed, 29 Oct 2008 09:35:55 +0100 From: Tor Christian Bekkvik <address@hidden> To: netcdf-java <address@hidden> -----Original Message----- From: Tor Christian Bekkvik Sent: 27. oktober 2008 14:48 To: 'address@hidden' Cc: Kjell Røang Subject: Adding a new projection to GRIB netcdf decoder Hi, We would like to add support for Rotated LatLon projection when reading GRIB files. (using netcdf decoders) To open the gribfile, we have used class NetcdfDataset: NetcdfDataset.openDataset(gribFile,...) But this fails for with rotated latlong projection, as itis not supported:ucar.grib.NoValidGribException: GDS: Unknown Grid Type :10) is notsupported. atucar.grib.grib1.Grib1GridDefinitionSection.<init>(Grib1GridDef initionSection.java:401)at ucar.grib.grib1.Grib1Input.scan(Grib1Input.java:155) atucar.grib.grib1.Grib1Indexer.writeFileIndex(Grib1Indexer.java:90)atucar.nc2.iosp.grib.GribServiceProvider.writeIndex(GribServiceP rovider.java:163)atucar.nc2.iosp.grib.GribServiceProvider.getIndex(GribServicePro vider.java:139)atucar.nc2.iosp.grib.GribServiceProvider.open(GribServiceProvide r.java:75)at ucar.nc2.NetcdfFile.<init>(NetcdfFile.java:1092) at ucar.nc2.NetcdfFile.open(NetcdfFile.java:485) at ucar.nc2.NetcdfFile.open(NetcdfFile.java:278) atucar.nc2.dataset.NetcdfDataset.openFile(NetcdfDataset.java:341)atucar.nc2.dataset.NetcdfDataset.openDataset(NetcdfDataset.java:175)atucar.nc2.dataset.NetcdfDataset.openDataset(NetcdfDataset.java:161)To add this new projection, it seems that at least thefollowing codeneeds update: 1) ProjectionImpl (I can implement RotLatLonProjection extends ProjectionImpl) 2) Grib1GridDefinitionSection (to fixNoValidGribException, listed above)3) GribHorizCoordSys.makeProjection (read necessaryparameters from GribGDS section and return RotLatLonProjection): private void makeProjection(NetcdfFile ncfile) { switch (lookup.getProjectionType(gdsIndex)) { case TableLookup.PolarStereographic: makePS(); break; case TableLookup.LambertConformal: makeLC(); break; case TableLookup.Mercator: makeMercator(); break; case TableLookup.Orthographic: makeSpaceViewOrOthographic(); break; default: throw new UnsupportedOperationException("unknownprojection = "+ gdsIndex.grid_type); } Currently, I don't have sufficient GRIB/NetCDF knowledge(or sourcecode) to do this. But I have already added Rotated LatLon projection for the JGRIB project (http://jgrib.sourceforge.net/), and we would liketo do thesame with NetCDF/GRIB. - Any comments how to best do this ? Thanks, Tor Chr Bekkvik _______________________________________________ netcdf-java mailing list address@hidden For list information or to unsubscribe, visit: http://www.unidata.ucar.edu/mailing_lists/============================================================== ================= Robb Kambic Unidata Program Center Software Engineer III Univ. Corp for Atmospheric Research address@hidden WWW: http://www.unidata.ucar.edu/ ============================================================== =============================================================================== ================= Robb Kambic Unidata Program Center Software Engineer III Univ. Corp for Atmospheric Research address@hidden WWW: http://www.unidata.ucar.edu/ ============================================================== ================================================================================================ Robb Kambic Unidata Program Center Software Engineer III Univ. Corp for Atmospheric Research address@hidden WWW: http://www.unidata.ucar.edu/ ============================================================================================================================================================== Robb Kambic Unidata Program Center Software Engineer III Univ. Corp for Atmospheric Research address@hidden WWW: http://www.unidata.ucar.edu/ ===============================================================================
=============================================================================== Robb Kambic Unidata Program Center Software Engineer III Univ. Corp for Atmospheric Research address@hidden WWW: http://www.unidata.ucar.edu/ ===============================================================================