|
RE: [vtp] converting lat/lon to terrain coords (fwd): msg#00099gis.gdal.devel
Hi GDALers, I'm new on the list - the Virtual Terrain Project folks suggested that I try asking this question here. I'm trying to convert from some lat/lon coords coming out of our gazetteer into the projection that the terrain I'm visualizing in VTP is using. My code is below, as is the text description of one of the (several) source projections I've tried, and the target projection (below the code). The problem is that I seem to be getting weird results out of the tranformation. (the method returns zero) This is a typical one... lat 34.400830 lon -119.738330 -> -2.089828, 0.600408 Inside VTP it looks like the numbers I should be getting are in the millions - maybe 6.00M for x and 1.97M for y. Someone suggested that it looks like I'm just getting the result in radians instead of feet, so I tried setting the linear units manually (HERE in the code) even though it looks ok in the text description. I've tried that and a bunch of different source projections (both WellKnown and EPSG codes) and all seem to yield the same results. Thanks in advance for your help. I got this far impressively fast - C++ is not my first language and when I get something working as quickly as I got into this, it means the documentation is good and the architecture is clean. Certainly seems true for GDAL. Dan void AClient::projectPoints(const vtProjection &targetProj){ OGRSpatialReference sourceProj; // will be geog OGRSpatialReference targetOGR; OGRCoordinateTransformation *trans; double x,y; DPoint3 p; char *ss = NULL; char *st = NULL; int res = 0; targetOGR = (OGRSpatialReference)targetProj; res = sourceProj.SetWellKnownGeogCS( "EPSG:4326" ); //res = sourceProj.importFromEPSG(4269); // HERE trying to set linear units manually... res = targetOGR.SetLinearUnits(SRS_UL_US_FOOT,1.0); sourceProj.exportToWkt(&ss); targetOGR.exportToWkt(&st); VTLOG("from proj: %s\nto proj:%s\n",ss,st); trans = OGRCreateCoordinateTransformation( &sourceProj, &targetOGR); if(!trans) VTLOG("argh."); for (int i=0; i < m_gazList.n; ++i) { x = m_gazList.gazItems[i].lon; y = m_gazList.gazItems[i].lat; int result = trans->Transform(1, &x, &y); if (result < 1) VTLOG("Transformation failed.\n"); else { p.x = x; p.y = y; p.z = 0; m_gazList.gazItems[i].pProj = p; VTLOG("lat %f lon %f -> %f, %f\n", m_gazList.gazItems[i].lat,m_gazList.gazItems[i].lon ,x,y); } } } source projection ----------------- GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563, AUTHORITY["EPSG","7030"]], AUTHORITY["EPSG","6326"]], PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]], UNIT["degree",0.0174532925199433],AUTHORITY["EPSG","4326"]] target projection ----------------- PROJCS["NAD83 / California zone 5", GEOGCS["NAD83", DATUM["North_American_Datum_1983", SPHEROID["GRS 1980",6378137,298.2572221010002, AUTHORITY["EPSG","7019"]], AUTHORITY["EPSG","6269"]], PRIMEM["Greenwich",0], UNIT["degree",0.0174532925199433],AUTHORITY["EPSG","4269"]], UNIT["US survey foot",0.3048006096012192, AUTHORITY["EPSG","9003"]],AUTHORITY["EPSG","26945"]] To unsubscribe from this group, send an email to: vtp-unsubscribe@xxxxxxxxxxx Your use of Yahoo! Groups is subject to http://docs.yahoo.com/info/terms/ |
|
| <Prev in Thread] | Current Thread | [Next in Thread> |
|---|---|---|
| Previous by Date: | Re: gdal_translate can now add GCPs to image, gdalwarptest renamed: 00099, Frank Warmerdam |
|---|---|
| Next by Date: | Re: Running gdal_merge.py: 00099, Sean Gillies |
| Previous by Thread: | Running gdal_merge.pyi: 00099, Reinaldo Escada Chohfi |
| Next by Thread: | RE: [vtp] converting lat/lon to terrain coords (fwd): 00099, dan ancona |
| Indexes: | [Date] [Thread] [Top] [All Lists] |
| News | FAQ | advertise |