There is now an additional wrinkle to the project. We need to now extract the 4 closest values for a given lat-lon point to complete a bilinear interpolation for added accuracy.
Would the correct method be: (1) Use the suggestion that Frank made below about a hybrid approach to my original (2). From this value I can easily determine the 4 closest neighbors in raster space. Since I need to know the distance between my original query point and these 4 corner points I could create a reverse transformation and push the 4 corner points back through the reverse transformation and extract the lat-lon to calculate distance. It would seem as rich as GDAL is that there may be a better way to do this but perhaps not. In a second related question. Has a book ever been produced by anyone that describes "GDAL cookbook" type of solutions? It would seem that the type of question I am asking come up often and a book written to an audience that includes both GIS professionals as well as general programmers with limited GIS knowledge would be helpful. Thanks, Bill -----Original Message----- From: Frank Warmerdam [mailto:[email protected]] Sent: Monday, July 20, 2009 3:52 PM To: Cassanova, Bill Cc: [email protected] Subject: Re: [gdal-dev] General Questions about gdal usage with Grib1 Data Set Cassanova, Bill wrote: > Hi, > > I am new to GDal and just beginning to dig into the intricacies of the > library and would like some guidance on the best way to proceed. > > I have a data set that is in a grib1 format but with a Lambert Conformal > Projection. If need be this file can be converted to either grib2 or > shape file using an application called degrib. > > My task is given this data-set I need to be able to extract values from > the file for a given lat-lon and don't necessarily want to code the > intricacies of the Lambert-Conformal-SP2 projection. So in other words > given a lat-lon of 38.59, - 97.61 I want to be able to interrogate the > raster and extract a value for that given coordinate. > > There seem to be a few options and just need advice from the season GDAL > guru's of the best way to proceed: > > (1) Convert the lambert-conformal grib1 to a cylindrical equidistant > grib1 or tiff using gdalwarp and do the simpler ratio based math on the > rectangular grid as opposed to the projected grid. > > (2) Use OGRCreateCoordinateTransformation specifying the source as > Lambert Conformal and the destination as cylindrical equidistant. My > presumption at this point is that the > OGRCoordinateTransformation->Transform( 1, &x, &y ) method call could > take a longitude as the x-argument, and latitude as the y-argument. A > successful transform call would replace the original lat-lon with the > new lat-lon and from this point the data grid ( I, J ) offsets could be > calculated using simply ratios. > > (3) Are they any more direct simpler ways of accomplishing the same > thing? Bill, I would suggest a variation on (2). You can call GDALCreateGenImgProjTransformer to create a transformer from image pixel/line to long/lat. The call is: void CPL_DLL * GDALCreateGenImgProjTransformer( GDALDatasetH hSrcDS, const char *pszSrcWKT, GDALDatasetH hDstDS, const char *pszDstWKT, int bGCPUseOK, double dfGCPErrorThreshold, int nOrder ); But you can leave the hDstDS as NULL, and provide the WGS84 pszDstWKT yourself. If the source file coordinate system is already known by GDAL (check gdalinfo report) then you could leave pszSrcWKT NULL and it will be picked up from the dataset. The bGCPUseOK can be FALSE, and the GCPErrorThreshold 0.0, and nOrder 0. The pointer you get back can be used with GDALGenImgProjTransform(): int CPL_DLL GDALGenImgProjTransform( void *pTransformArg, int bDstToSrc, int nPointCount, double *x, double *y, double *z, int *panSuccess ); In your case you want to pass in long/lat, and get back image pixel/line so you would pass bDstToSrc=TRUE. This approach is somewhat easier than OGRCreateCoordinateTransformation because the GenImgProj stuff knows also how to apply the geotransform step. Best regards, -- ---------------------------------------+-------------------------------- ------ I set the clouds in motion - turn up | Frank Warmerdam, [email protected] light and sound - activate the windows | http://pobox.com/~warmerdam and watch the world go round - Rush | Geospatial Programmer for Rent _______________________________________________ gdal-dev mailing list [email protected] http://lists.osgeo.org/mailman/listinfo/gdal-dev
