<html xmlns:v="urn:schemas-microsoft-com:vml" xmlns:o="urn:schemas-microsoft-com:office:office" xmlns:w="urn:schemas-microsoft-com:office:word" xmlns:m="http://schemas.microsoft.com/office/2004/12/omml" xmlns="http://www.w3.org/TR/REC-html40"><head><meta http-equiv=Content-Type content="text/html; charset=us-ascii"><meta name=Generator content="Microsoft Word 12 (filtered medium)"><style><!--
/* Font Definitions */
@font-face
{font-family:"Arial Unicode MS";
panose-1:2 11 6 4 2 2 2 2 2 4;}
@font-face
{font-family:Calibri;
panose-1:2 15 5 2 2 2 4 3 2 4;}
@font-face
{font-family:"\@Arial Unicode MS";
panose-1:2 11 6 4 2 2 2 2 2 4;}
/* Style Definitions */
p.MsoNormal, li.MsoNormal, div.MsoNormal
{margin:0in;
margin-bottom:.0001pt;
font-size:11.0pt;
font-family:"Calibri","sans-serif";}
a:link, span.MsoHyperlink
{mso-style-priority:99;
color:blue;
text-decoration:underline;}
a:visited, span.MsoHyperlinkFollowed
{mso-style-priority:99;
color:purple;
text-decoration:underline;}
span.EmailStyle17
{mso-style-type:personal-compose;
font-family:"Calibri","sans-serif";
color:windowtext;}
.MsoChpDefault
{mso-style-type:export-only;}
@page WordSection1
{size:8.5in 11.0in;
margin:1.0in 1.0in 1.0in 1.0in;}
div.WordSection1
{page:WordSection1;}
--></style><!--[if gte mso 9]><xml>
<o:shapedefaults v:ext="edit" spidmax="1026" />
</xml><![endif]--><!--[if gte mso 9]><xml>
<o:shapelayout v:ext="edit">
<o:idmap v:ext="edit" data="1" />
</o:shapelayout></xml><![endif]--></head><body lang=EN-US link=blue vlink=purple><div class=WordSection1><p class=MsoNormal>Hi All,<o:p></o:p></p><p class=MsoNormal><o:p> </o:p></p><p class=MsoNormal>I am experiencing a bit of confusion regarding the use of GDALCreateGenImgProjTransformer.<o:p></o:p></p><p class=MsoNormal><o:p> </o:p></p><p class=MsoNormal>I am attempting to create a ProjTransformer in 2 different ways to prove to myself that I would get the same result…The problem is I don’t but don’t understand exactly why.<o:p></o:p></p><p class=MsoNormal><o:p> </o:p></p><p class=MsoNormal>The first method I am trying is to open the GRIB2 dataset using GDALOpen to get a GDAL dataset handle. This handle is then passed to GDALCreateGenImgProjTransformer as the first argument with the second argument being passed as NULL. The destination dataset is one created using the well-known “WGS84” Tag…For example<o:p></o:p></p><p class=MsoNormal><o:p> </o:p></p><p class=MsoNormal>OGRSpatialReferenceH hDEST;<o:p></o:p></p><p class=MsoNormal> hDEST = OSRNewSpatialReference(NULL);<o:p></o:p></p><p class=MsoNormal> OSRSetWellKnownGeogCS(hDEST, "WGS84");<o:p></o:p></p><p class=MsoNormal> OSRExportToWkt(hDEST, &pszDstWKT);<o:p></o:p></p><p class=MsoNormal> OSRDestroySpatialReference(hDEST);<o:p></o:p></p><p class=MsoNormal><o:p> </o:p></p><p class=MsoNormal>void *my_transformer = GDALCreateGenImgProjTransformer(hSrcDS,<o:p></o:p></p><p class=MsoNormal> NULL,<o:p></o:p></p><p class=MsoNormal> NULL,<o:p></o:p></p><p class=MsoNormal> pszDstWKT,<o:p></o:p></p><p class=MsoNormal> FALSE, 0.0, 0);<o:p></o:p></p><p class=MsoNormal><o:p> </o:p></p><p class=MsoNormal>The second method passes NULL for the first argument and instead passed the WKT instead. The destination data set is the same as created above.<o:p></o:p></p><p class=MsoNormal><o:p> </o:p></p><p class=MsoNormal> void *my_second_transformer = GDALCreateGenImgProjTransformer( NULL,<o:p></o:p></p><p class=MsoNormal> GDALGetProjectionRef( hSrcDS ),<o:p></o:p></p><p class=MsoNormal> NULL,<o:p></o:p></p><p class=MsoNormal> pszDstWKT,<o:p></o:p></p><p class=MsoNormal> FALSE,<o:p></o:p></p><p class=MsoNormal> 0.0, 0 );<o:p></o:p></p><p class=MsoNormal><o:p> </o:p></p><p class=MsoNormal>The result of GDALGetProjectionRef( hSrcDS ) for reference :<o:p></o:p></p><p class=MsoNormal> <o:p></o:p></p><p class=MsoNormal>PROJCS["unnamed",<o:p></o:p></p><p class=MsoNormal> GEOGCS["Coordinate System imported from GRIB file",<o:p></o:p></p><p class=MsoNormal> DATUM["unknown",<o:p></o:p></p><p class=MsoNormal> SPHEROID["Sphere",6371229,0]],<o:p></o:p></p><p class=MsoNormal> PRIMEM["Greenwich",0],<o:p></o:p></p><p class=MsoNormal> UNIT["degree",0.0174532925199433]],<o:p></o:p></p><p class=MsoNormal> PROJECTION["Lambert_Conformal_Conic_2SP"],<o:p></o:p></p><p class=MsoNormal> PARAMETER["standard_parallel_1",25],<o:p></o:p></p><p class=MsoNormal> PARAMETER["standard_parallel_2",25],<o:p></o:p></p><p class=MsoNormal> PARAMETER["latitude_of_origin",0],<o:p></o:p></p><p class=MsoNormal> PARAMETER["central_meridian",265],<o:p></o:p></p><p class=MsoNormal> PARAMETER["false_easting",0],<o:p></o:p></p><p class=MsoNormal> PARAMETER["false_northing",0]]<o:p></o:p></p><p class=MsoNormal><o:p> </o:p></p><p class=MsoNormal>Lambert Conformal Grid Information as extracted from gdalinfo:<o:p></o:p></p><p class=MsoNormal>Origin = (-4226106.996915472671390,7240925.337282469496131)<o:p></o:p></p><p class=MsoNormal>Pixel Size = (12191.000000000000000,-12191.000000000000000)<o:p></o:p></p><p class=MsoNormal>Corner Coordinates:<o:p></o:p></p><p class=MsoNormal>Upper Left (-4226106.997, 7240925.337) (152d52'43.04"W, 54d33'55.24"N)<o:p></o:p></p><p class=MsoNormal>Lower Left (-4226106.997, 2023177.337) (133d25'42.19"W, 12d 5'14.28"N)<o:p></o:p></p><p class=MsoNormal>Upper Right ( 3259167.003, 7240925.337) ( 49d15'27.05"W, 57d17'54.49"N)<o:p></o:p></p><p class=MsoNormal>Lower Right ( 3259167.003, 2023177.337) ( 65d 2'28.98"W, 14d12'48.02"N)<o:p></o:p></p><p class=MsoNormal>Center ( -483469.997, 4632051.337) (100d30'19.44"W, 40d34'21.84"N)<o:p></o:p></p><p class=MsoNormal><o:p> </o:p></p><p class=MsoNormal><o:p> </o:p></p><p class=MsoNormal>To test both of these I pass in the known coordinates up the upper left hand corner of the grid ( -152.877, 54.567 ) using code:<o:p></o:p></p><p class=MsoNormal><o:p> </o:p></p><p class=MsoNormal>double x_coord = atof( argv[2] );<o:p></o:p></p><p class=MsoNormal>double y_coord = atof( argv[3] );<o:p></o:p></p><p class=MsoNormal><o:p> </o:p></p><p class=MsoNormal> std::cout << std::fixed << "original x_coord: " << x_coord << " y_coord: " << y_coord << std::endl;<o:p></o:p></p><p class=MsoNormal><o:p> </o:p></p><p class=MsoNormal> result = GDALGenImgProjTransform( my_transformer, TRUE, 1, &x_coord, &y_coord,<o:p></o:p></p><p class=MsoNormal> NULL, &panSuccess);<o:p></o:p></p><p class=MsoNormal><o:p> </o:p></p><p class=MsoNormal> std::cout << std::fixed << "calculated x_coord: " << x_coord << " y_coord: " << y_coord << std::endl;<o:p></o:p></p><p class=MsoNormal><o:p> </o:p></p><p class=MsoNormal>// Reset x-coord and y-coord from command line.<o:p></o:p></p><p class=MsoNormal> x_coord = atof( argv[2] );<o:p></o:p></p><p class=MsoNormal> y_coord = atof( argv[3] );<o:p></o:p></p><p class=MsoNormal><o:p> </o:p></p><p class=MsoNormal> std::cout << std::fixed << "second original x_coord: " << x_coord << " y_coord: " << y_coord << std::endl;<o:p></o:p></p><p class=MsoNormal><o:p> </o:p></p><p class=MsoNormal> esult =<o:p></o:p></p><p class=MsoNormal> GDALGenImgProjTransform( my_second_transformer, TRUE, 1, &x_coord, &y_coord,<o:p></o:p></p><p class=MsoNormal> NULL, &panSuccess);<o:p></o:p></p><p class=MsoNormal><o:p> </o:p></p><p class=MsoNormal> std::cout << std::fixed << "second calculated x_coord: " << x_coord << " y_coord: " << y_coord << std::endl;<o:p></o:p></p><p class=MsoNormal><o:p> </o:p></p><p class=MsoNormal>The results:<o:p></o:p></p><p class=MsoNormal><o:p> </o:p></p><p class=MsoNormal><b>// This result is right on the money and exactly what I expect<o:p></o:p></b></p><p class=MsoNormal>original x_coord: -152.877000 y_coord: 54.567000<o:p></o:p></p><p class=MsoNormal>calculated x_coord: 0.016427 y_coord: -0.011915<o:p></o:p></p><p class=MsoNormal><o:p> </o:p></p><p class=MsoNormal><b>// This result is not what I expect. I would have assumed since I created the transformer in a different way but using the WKT information it should have been the same.<o:p></o:p></b></p><p class=MsoNormal><b>// instead I am getting the Lambert Easting and Northing in Meters.<o:p></o:p></b></p><p class=MsoNormal>second original x_coord: -152.877000 y_coord: 54.567000<o:p></o:p></p><p class=MsoNormal>second calculated x_coord: -4225906.734228 y_coord: 7241070.590061<o:p></o:p></p><p class=MsoNormal><o:p> </o:p></p><p class=MsoNormal>Am I missing a step here?<o:p></o:p></p><p class=MsoNormal><o:p> </o:p></p><p class=MsoNormal>Thanks,<o:p></o:p></p><p class=MsoNormal>Bill<o:p></o:p></p><p class=MsoNormal><o:p> </o:p></p><p class=MsoNormal><span style='font-size:9.0pt'>++++<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:9.0pt'><o:p> </o:p></span></p><p class=MsoNormal><span style='font-size:9.0pt'>William Cassanova | Senior GFS Developer | The Weather Channel | 770.226.2368 | bcassanova@weather.com<o:p></o:p></span></p><p class=MsoNormal><o:p> </o:p></p></div></body></html>