[gdal-dev] Unable to read PDS3 (img + lbl) Moon DEM
I've tried to access a lunar elevation model, but gdalinfo and gdal_translate give me a strange error message. Error1: Sample_bits of 0 is not supported in this gdal PDS reader The elevation model is a PDS3 pair (.img with seperate .lbl file), and pointing gdalinfo to the lbl file gives me the above error message (pointing to the .img directly gives not a recognized file format): IMG (binary) http://imbrium.mit.edu/DATA/LOLA_GDR/CYLINDRICAL/IMG/LDEM_256_00N_90N_000_180.IMG LBL (plaintext label) http://imbrium.mit.edu/DATA/LOLA_GDR/CYLINDRICAL/IMG/LDEM_256_00N_90N_000_180.LBL the strange thing is - the label file explicitely specifies says SAMPLE_BITS = 16 in line 52. Is this a bug in the reader, an invalid label file or am i doing something wrong? -- View this message in context: http://osgeo-org.1560.x6.nabble.com/Unable-to-read-PDS3-img-lbl-Moon-DEM-tp5150432.html Sent from the GDAL - Dev mailing list archive at Nabble.com. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] python scripts - create layer problem
Hi Martin, I would say that your Zsj layer is a non-spatial layer. Is it the case ? If yes, then the output of the log and the content of geometry_columns are consistant. You may not see it directly when running ogrinfo on the PostGIS DB since by default only spatial layers are listed. You can however query it directly with GetLayerByName(), or define PG_LIST_ALL_TABLES=YES as a configuration option (see http://gdal.org/drv_pg_advanced.html) If your layer is supposed to have a geometry column, then the problem is likely on the read side. Even ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Unable to get features from OpenFileGDB data source
Nik,
By using GetFeature() the way you do, you assume that the passed parameter is an
index number whereas it is a feature identifier. In most drivers that don't
really support the concept of feature identifier, both concepts are identical.
But in FileGDB, Feature Identifiers start at 1, and not 0. And you can have
holes in the numbering due to deleted records (this is the case with layer
PointBonnetHill that has a single feature whose FID is 2)
For efficiency and robustness with all formats, I'd encourage you to use
GetNextFeature() instead. And you don't need to call OGR_L_GetFeatureCount()
which can be a slow operation on some drivers. Just iterate over features with
GetNextFeature() until it returns a NULL pointer.
Best regards,
Even
Hi GDAL devs,
I'm trying to support ESRI File Geodatabase support (read-only) in my iOS app
and have recently got GDAL 1.11.0 compiled and running in the app for its new
OpenFileGDB driver. However, I've been unable to extract features from GDBs
using the exact same code that works fine for other formats/drivers.
An extract of the code in question is:
ds = OGROpen(path, FALSE, NULL);
layerCount = OGR_DS_GetLayerCount(ds);
for ( int i = 0; i layerCount; i++ )
{
OGRLayerH sourceLayer = OGR_DS_GetLayer(ds, i);
int sourceFeatureCount = OGR_L_GetFeatureCount(sourceLayer,
YES);
for ( long j = 0; j sourceFeatureCount; j++ )
{
OGRFeatureH sourceFeature =
OGR_L_GetFeature(sourceLayer, j);
if ( ! sourceFeature )
NSLog(@Failed to get feature at index '%ld'
from layer at index '%d'!,
j, i);
The NSLog line is always run for my GDB data sources but never for data
sources of other formats. So for GDB data sources 'sourceFeature' is always
NULL and of course this causes all subsequent OGR functions to which
'sourceFeature' is passed to fail with errors like:
ERROR 10: Pointer 'hFeat' is NULL in 'OGR_...'.
Using the CLI utility 'ogrinfo' (different build of GDAL) on the same
datasource works fine showing all the features are there as expected and it
works fine in ESRI ArcGIS. And the fact that it gets into the second 'for'
loop at all indicates that there are features there, and yet it cannot get
the feature at index '0'.
The File Geodatabase datasource can be downloaded for testing from:
https://dl.dropboxusercontent.com/u/12436846/Bandicoots.gdb.zip
Can anybody advise me on how I should trouble shoot this further, or how the
problem can be resolved?
Cheers,
Nik.
___
gdal-dev mailing list
gdal-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/gdal-dev
___
gdal-dev mailing list
gdal-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Unable to read PDS3 (img + lbl) Moon DEM
Selon GrJo gregor.jochm...@rt.rif-ev.de: I've tried to access a lunar elevation model, but gdalinfo and gdal_translate give me a strange error message. Error1: Sample_bits of 0 is not supported in this gdal PDS reader The elevation model is a PDS3 pair (.img with seperate .lbl file), and pointing gdalinfo to the lbl file gives me the above error message (pointing to the .img directly gives not a recognized file format): IMG (binary) http://imbrium.mit.edu/DATA/LOLA_GDR/CYLINDRICAL/IMG/LDEM_256_00N_90N_000_180.IMG LBL (plaintext label) http://imbrium.mit.edu/DATA/LOLA_GDR/CYLINDRICAL/IMG/LDEM_256_00N_90N_000_180.LBL the strange thing is - the label file explicitely specifies says SAMPLE_BITS = 16 in line 52. Is this a bug in the reader, an invalid label file or am i doing something wrong? Nothing wrong. You are just using an old version of GDAL. Your issue is http://trac.osgeo.org/gdal/ticket/3943 that was fixed in GDAL 1.9.0. -- View this message in context: http://osgeo-org.1560.x6.nabble.com/Unable-to-read-PDS3-img-lbl-Moon-DEM-tp5150432.html Sent from the GDAL - Dev mailing list archive at Nabble.com. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
[gdal-dev] warp taking time
Dear All, I am trying to generate mosaic of 28 input images using *gdalwarp*. All are in same projection. It has taken 5 days to generate the final output. Please find the gdalinfo of input and output attached. The input files are passed to warp in the order in which the gdalinfo is present in the attachement. Any input on improving the time would be really helpful. Thanks input_files_gdalinfo.txt http://osgeo-org.1560.x6.nabble.com/file/n5150453/input_files_gdalinfo.txt final_output_gdalinfo.txt http://osgeo-org.1560.x6.nabble.com/file/n5150453/final_output_gdalinfo.txt -- View this message in context: http://osgeo-org.1560.x6.nabble.com/warp-taking-time-tp5150453.html Sent from the GDAL - Dev mailing list archive at Nabble.com. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
[gdal-dev] warp taking time
Dear All, I am trying to generate mosaic of 28 input images using gdalwarp. All are in same projection. It has taken 5 days to generate the final output. Please find the gdalinfo of input and output attached. The input files are passed to warp in the order in which the gdalinfo is present in the attachement. Any input on improving the time would be really helpful. Thanks input_files_gdalinfo.txt http://osgeo-org.1560.x6.nabble.com/file/n5150454/input_files_gdalinfo.txt final_output_gdalinfo.txt http://osgeo-org.1560.x6.nabble.com/file/n5150454/final_output_gdalinfo.txt -- View this message in context: http://osgeo-org.1560.x6.nabble.com/gdal-dev-warp-taking-time-tp5150454.html Sent from the GDAL - Dev mailing list archive at Nabble.com. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] warp taking time
garfy, Please provide the gdalwarp command you used. Did you modify GDAL_CACHEMAX? On Thu, Jul 10, 2014 at 3:26 PM, garfy adiba.niz...@gmail.com wrote: Dear All, I am trying to generate mosaic of 28 input images using *gdalwarp*. All are in same projection. It has taken 5 days to generate the final output. Please find the gdalinfo of input and output attached. The input files are passed to warp in the order in which the gdalinfo is present in the attachement. Any input on improving the time would be really helpful. Thanks input_files_gdalinfo.txt http://osgeo-org.1560.x6.nabble.com/file/n5150453/input_files_gdalinfo.txt final_output_gdalinfo.txt http://osgeo-org.1560.x6.nabble.com/file/n5150453/final_output_gdalinfo.txt -- View this message in context: http://osgeo-org.1560.x6.nabble.com/warp-taking-time-tp5150453.html Sent from the GDAL - Dev mailing list archive at Nabble.com. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev -- Best regards, Chaitanya kumar CH. +91-9494447584 17.2416N 80.1426E ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Unable to read PDS3 (img + lbl) Moon DEM
Thanks a lot, i was using a version that came with the outdated FWTools. Works great with the new binaries (1.11) from gisinternals.com -- View this message in context: http://osgeo-org.1560.x6.nabble.com/Unable-to-read-PDS3-img-lbl-Moon-DEM-tp5150432p5150457.html Sent from the GDAL - Dev mailing list archive at Nabble.com. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] Can't compile 1.11.0 for iOS with 'expected statement' error for gdalgrid.cpp
Selon Nik Sands nix...@nixanz.com: PS. (This always happens just AFTER I send to the list, no matter how long I mull it over first...) Actually, I think I see the issue now. The closing bracket is the end of an if () stanza in which there is no content at all if neither 'HAVE_AVX_AT_COMPILE_TIME' or 'HAVE_SSE_AT_COMPILE_TIME' are defined, which must be the case in my situation. I guess I could either delete the entire stanza, or add in a statement that does nothing useful. Anyhow, I thought you guys should know this situation exists so that you may be able to deal with it yourselves for some future release. :-) Fixed per http://trac.osgeo.org/gdal/ticket/5566 In fact the error was when HAVE_AVX_AT_COMPILE_TIME is available but not HAVE_SSE_AT_COMPILE_TIME. Which shouldn't happen normally. So I guess there's an issue in the SSE detection case with that enivronment. Even -- Geospatial professional services http://even.rouault.free.fr/services.html ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] warp taking time
I am wondering why you are using gdalwarp for mosaicing images having same projection instead of gdal_merge.py? Bas Smit On Thu, Jul 10, 2014 at 7:14 PM, Chaitanya kumar CH chaitanya...@gmail.com wrote: garfy, Please provide the gdalwarp command you used. Did you modify GDAL_CACHEMAX? On Thu, Jul 10, 2014 at 3:26 PM, garfy adiba.niz...@gmail.com wrote: Dear All, I am trying to generate mosaic of 28 input images using *gdalwarp*. All are in same projection. It has taken 5 days to generate the final output. Please find the gdalinfo of input and output attached. The input files are passed to warp in the order in which the gdalinfo is present in the attachement. Any input on improving the time would be really helpful. Thanks input_files_gdalinfo.txt http://osgeo-org.1560.x6.nabble.com/file/n5150453/input_files_gdalinfo.txt final_output_gdalinfo.txt http://osgeo-org.1560.x6.nabble.com/file/n5150453/final_output_gdalinfo.txt -- View this message in context: http://osgeo-org.1560.x6.nabble.com/warp-taking-time-tp5150453.html Sent from the GDAL - Dev mailing list archive at Nabble.com. ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev -- Best regards, Chaitanya kumar CH. +91-9494447584 17.2416N 80.1426E ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] warp taking time
Hello, I’m interested in this thread. One comment I can add about gdal_merge.py is that it does not properly handle nodata values from the full (padded) extent of the input files. gdalwarp, which can be VERY slow, does handle the nodata values correctly. Michele Thornton From: gdal-dev-boun...@lists.osgeo.org [mailto:gdal-dev-boun...@lists.osgeo.org] On Behalf Of bas smit Sent: Thursday, July 10, 2014 7:52 AM To: Chaitanya kumar CH Cc: gdal dev Subject: Re: [gdal-dev] warp taking time I am wondering why you are using gdalwarp for mosaicing images having same projection instead of gdal_merge.py? Bas Smit On Thu, Jul 10, 2014 at 7:14 PM, Chaitanya kumar CH chaitanya...@gmail.commailto:chaitanya...@gmail.com wrote: garfy, Please provide the gdalwarp command you used. Did you modify GDAL_CACHEMAX? On Thu, Jul 10, 2014 at 3:26 PM, garfy adiba.niz...@gmail.commailto:adiba.niz...@gmail.com wrote: Dear All, I am trying to generate mosaic of 28 input images using *gdalwarp*. All are in same projection. It has taken 5 days to generate the final output. Please find the gdalinfo of input and output attached. The input files are passed to warp in the order in which the gdalinfo is present in the attachement. Any input on improving the time would be really helpful. Thanks input_files_gdalinfo.txt http://osgeo-org.1560.x6.nabble.com/file/n5150453/input_files_gdalinfo.txt final_output_gdalinfo.txt http://osgeo-org.1560.x6.nabble.com/file/n5150453/final_output_gdalinfo.txt -- View this message in context: http://osgeo-org.1560.x6.nabble.com/warp-taking-time-tp5150453.html Sent from the GDAL - Dev mailing list archive at Nabble.com. ___ gdal-dev mailing list gdal-dev@lists.osgeo.orgmailto:gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev -- Best regards, Chaitanya kumar CH. +91-9494447584tel:%2B91-9494447584 17.2416N 80.1426E ___ gdal-dev mailing list gdal-dev@lists.osgeo.orgmailto:gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
[gdal-dev] State Plane to WGS84 Transformation problem with false easting and northing values
Frank, Even or whoever, I have a source geotiff that has the following projection: PROJCS[State Plane,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],AUTHORIT Y[EPSG,4269]],PROJECTION[Lambert_Conformal_Conic_2SP],PARAMETER[sta ndard_parallel_1,35.47],PARAMETER[standard_parallel_2,34.03 ],PARAMETER[latitude_of_origin,33.5],PARAMETER[central_meri dian,-118],PARAMETER[false_easting,21527734.7222],PARAMETER[false_ northing,5381933.68055],UNIT[feet,0.3048006096012192],AUTHORITY[EP SG,26945]] I want to convert the bounds to WGS84. When I setup my coordinate transform object the OGRProj4CT::InitializeNoLock() function in ogrct.cpp calls poSRSSource-exportToProj4() and returns the following: +proj=lcc +lat_1=35.47 +lat_2=34.03 +lat_0=33.5 +lon_0=-118 +x_0=6561666.6 +y_0=1640416.7 +datum=NAD83 +units=us-ft +no_defs Then, when calling TransformEx() on the bounds the resulting coordinates are wrong. It appears that somewhere in the transformation process the UNIT RADIANS value is not being factored correctly with the false eating and northing values. If I translate the source image into a new geotiff with the exact same parameters as the source geotiff then the new projection is as follows: 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],AUTHORIT Y[EPSG,4269]],PROJECTION[Lambert_Conformal_Conic_2SP],PARAMETER[sta ndard_parallel_1,35.47],PARAMETER[standard_parallel_2,34.03 ],PARAMETER[latitude_of_origin,33.5],PARAMETER[central_meri dian,-118],PARAMETER[false_easting,6561666.6],PARAMETER[false_ northing,1640416.7],UNIT[US survey foot,0.3048006096012192,AUTHORITY[EPSG,9003]],AUTHORITY[EPSG,26945 ]] And the exportToProj4() returns this: +proj=lcc +lat_1=35.47 +lat_2=34.03 +lat_0=33.5 +lon_0=-118 +x_0=200 +y_0=50.01 +datum=NAD83 +units=us-ft +no_defs Which produces a correct result when using TransformEx(). When I get the SRS_PP_FALSE_EASTING and SRS_PP_FALSE_NORTHING from both projections using the OGRSpatialReference::GetProjParm() function I notice the difference between the two is a factor of the GEOGCS unit radians: SOURCE IMAGE: false_easting: 21527734.7 false_northing: 5381933.680556 TRANSLATED IMAGE: false_easting: 6561666.67 false_northing: 1640416.67 The translated offsets divided by the source offsets seem to be a factor of the unit radians: 0.3048 I think the reason that the translate to geotiff corrects the problem is that the geotiff keys for the projection are calculated in a different manner than when using the OGRCoordinateTransform object. Any help on where to apply the fix in the GDAL code if you think one is required would be much appreciated. If you think I need to accommodate this in my code any help would be greatly appreciated. Best regards, Martin ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
[gdal-dev] State Plane to WGS84 Transformation problem with false easting and northing values
Frank, Even or whoever, I have a source geotiff that has the following projection: PROJCS[State Plane,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],AUTHORIT Y[EPSG,4269]],PROJECTION[Lambert_Conformal_Conic_2SP],PARAMETER[sta ndard_parallel_1,35.47],PARAMETER[standard_parallel_2,34.03 ],PARAMETER[latitude_of_origin,33.5],PARAMETER[central_meri dian,-118],PARAMETER[false_easting,21527734.7222],PARAMETER[false_ northing,5381933.68055],UNIT[feet,0.3048006096012192],AUTHORITY[EP SG,26945]] I want to convert the bounds to WGS84. When I setup my coordinate transform object the OGRProj4CT::InitializeNoLock() function in ogrct.cpp calls poSRSSource-exportToProj4() and returns the following: +proj=lcc +lat_1=35.47 +lat_2=34.03 +lat_0=33.5 +lon_0=-118 +x_0=6561666.6 +y_0=1640416.7 +datum=NAD83 +units=us-ft +no_defs Then, when calling TransformEx() on the bounds the resulting coordinates are wrong. It appears that somewhere in the transformation process the UNIT RADIANS value is not being factored correctly with the false eating and northing values. If I translate the source image into a new geotiff with the exact same parameters as the source geotiff then the new projection is as follows: 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],AUTHORIT Y[EPSG,4269]],PROJECTION[Lambert_Conformal_Conic_2SP],PARAMETER[sta ndard_parallel_1,35.47],PARAMETER[standard_parallel_2,34.03 ],PARAMETER[latitude_of_origin,33.5],PARAMETER[central_meri dian,-118],PARAMETER[false_easting,6561666.6],PARAMETER[false_ northing,1640416.7],UNIT[US survey foot,0.3048006096012192,AUTHORITY[EPSG,9003]],AUTHORITY[EPSG,26945 ]] And the exportToProj4() returns this: +proj=lcc +lat_1=35.47 +lat_2=34.03 +lat_0=33.5 +lon_0=-118 +x_0=200 +y_0=50.01 +datum=NAD83 +units=us-ft +no_defs Which produces a correct result when using TransformEx(). When I get the SRS_PP_FALSE_EASTING and SRS_PP_FALSE_NORTHING from both projections using the OGRSpatialReference::GetProjParm() function I notice the difference between the two is a factor of the GEOGCS unit radians: SOURCE IMAGE: false_easting: 21527734.7 false_northing: 5381933.680556 TRANSLATED IMAGE: false_easting: 6561666.67 false_northing: 1640416.67 The translated offsets divided by the source offsets seem to be a factor of the unit radians: 0.3048 I think the reason that the translate to geotiff corrects the problem is that the geotiff keys for the projection are calculated in a different manner than when using the OGRCoordinateTransform object. Any help on where to apply the fix in the GDAL code if you think one is required would be much appreciated. If you think I need to accommodate this in my code any help would be greatly appreciated. Best regards, Martin ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] State Plane to WGS84 Transformation problem with false easting and northing values
Martin, Yes this issue has been reported several times (last time was http://lists.osgeo.org/pipermail/gdal-dev/2014-May/039251.html) and is a GDAL bug that the user code shouldn't compensate for. I think it is https://trac.osgeo.org/gdal/ticket/3901 / https://trac.osgeo.org/gdal/ticket/4954 The 0.3048 ratio you see is the length of a foot in meters (not sure why you mention unit radians here). Basically the conversion from foot to meters is wrongly applied twice in the geotiff SRS decoding. 21527734.7222 * 0.3048 * 0.3048 = 200 Something I've looked at a bit, but did not pursue because of the lack of incentive to do so. Even Le jeudi 10 juillet 2014 22:12:04, Martin Chapman a écrit : Frank, Even or whoever, I have a source geotiff that has the following projection: PROJCS[State Plane,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],AUTHORIT Y[EPSG,4269]],PROJECTION[Lambert_Conformal_Conic_2SP],PARAMETER[sta ndard_parallel_1,35.47],PARAMETER[standard_parallel_2,34.03 ],PARAMETER[latitude_of_origin,33.5],PARAMETER[central_meri dian,-118],PARAMETER[false_easting,21527734.7222],PARAMETER[false_ northing,5381933.68055],UNIT[feet,0.3048006096012192],AUTHORITY[EP SG,26945]] I want to convert the bounds to WGS84. When I setup my coordinate transform object the OGRProj4CT::InitializeNoLock() function in ogrct.cpp calls poSRSSource-exportToProj4() and returns the following: +proj=lcc +lat_1=35.47 +lat_2=34.03 +lat_0=33.5 +lon_0=-118 +x_0=6561666.6 +y_0=1640416.7 +datum=NAD83 +units=us-ft +no_defs Then, when calling TransformEx() on the bounds the resulting coordinates are wrong. It appears that somewhere in the transformation process the UNIT RADIANS value is not being factored correctly with the false eating and northing values. If I translate the source image into a new geotiff with the exact same parameters as the source geotiff then the new projection is as follows: 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],AUTHORIT Y[EPSG,4269]],PROJECTION[Lambert_Conformal_Conic_2SP],PARAMETER[sta ndard_parallel_1,35.47],PARAMETER[standard_parallel_2,34.03 ],PARAMETER[latitude_of_origin,33.5],PARAMETER[central_meri dian,-118],PARAMETER[false_easting,6561666.6],PARAMETER[false_ northing,1640416.7],UNIT[US survey foot,0.3048006096012192,AUTHORITY[EPSG,9003]],AUTHORITY[EPSG,26945 ]] And the exportToProj4() returns this: +proj=lcc +lat_1=35.47 +lat_2=34.03 +lat_0=33.5 +lon_0=-118 +x_0=200 +y_0=50.01 +datum=NAD83 +units=us-ft +no_defs Which produces a correct result when using TransformEx(). When I get the SRS_PP_FALSE_EASTING and SRS_PP_FALSE_NORTHING from both projections using the OGRSpatialReference::GetProjParm() function I notice the difference between the two is a factor of the GEOGCS unit radians: SOURCE IMAGE: false_easting: 21527734.7 false_northing: 5381933.680556 TRANSLATED IMAGE: false_easting: 6561666.67 false_northing: 1640416.67 The translated offsets divided by the source offsets seem to be a factor of the unit radians: 0.3048 I think the reason that the translate to geotiff corrects the problem is that the geotiff keys for the projection are calculated in a different manner than when using the OGRCoordinateTransform object. Any help on where to apply the fix in the GDAL code if you think one is required would be much appreciated. If you think I need to accommodate this in my code any help would be greatly appreciated. Best regards, Martin -- Geospatial professional services http://even.rouault.free.fr/services.html ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
Re: [gdal-dev] State Plane to WGS84 Transformation problem with false easting and northing values
Even, Thanks for the quick response. Yeah, I goofed, I meant to say units, guess I am just getting old. Anyway, thanks for pointing me in the right direction! Best regards, Martin -Original Message- From: Even Rouault [mailto:even.roua...@mines-paris.org] Sent: Thursday, July 10, 2014 2:35 PM To: gdal-dev@lists.osgeo.org Cc: Martin Chapman Subject: Re: [gdal-dev] State Plane to WGS84 Transformation problem with false easting and northing values Martin, Yes this issue has been reported several times (last time was http://lists.osgeo.org/pipermail/gdal-dev/2014-May/039251.html) and is a GDAL bug that the user code shouldn't compensate for. I think it is https://trac.osgeo.org/gdal/ticket/3901 / https://trac.osgeo.org/gdal/ticket/4954 The 0.3048 ratio you see is the length of a foot in meters (not sure why you mention unit radians here). Basically the conversion from foot to meters is wrongly applied twice in the geotiff SRS decoding. 21527734.7222 * 0.3048 * 0.3048 = 200 Something I've looked at a bit, but did not pursue because of the lack of incentive to do so. Even Le jeudi 10 juillet 2014 22:12:04, Martin Chapman a écrit : Frank, Even or whoever, I have a source geotiff that has the following projection: PROJCS[State Plane,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],AUTH ORIT Y[EPSG,4269]],PROJECTION[Lambert_Conformal_Conic_2SP],PARAMETER[ sta ndard_parallel_1,35.47],PARAMETER[standard_parallel_2,3 4.03 ],PARAMETER[latitude_of_origin,33.5],PARAMETER[central_ meri dian,-118],PARAMETER[false_easting,21527734.7222],PARAMETER[fa lse_ northing,5381933.68055],UNIT[feet,0.3048006096012192],AUTHORITY [EP SG,26945]] I want to convert the bounds to WGS84. When I setup my coordinate transform object the OGRProj4CT::InitializeNoLock() function in ogrct.cpp calls poSRSSource-exportToProj4() and returns the following: +proj=lcc +lat_1=35.47 +lat_2=34.03 ++lat_0=33.5 +lon_0=-118 +x_0=6561666.6 +y_0=1640416.7 +datum=NAD83 +units=us-ft +no_defs Then, when calling TransformEx() on the bounds the resulting coordinates are wrong. It appears that somewhere in the transformation process the UNIT RADIANS value is not being factored correctly with the false eating and northing values. If I translate the source image into a new geotiff with the exact same parameters as the source geotiff then the new projection is as follows: 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],AUTH ORIT Y[EPSG,4269]],PROJECTION[Lambert_Conformal_Conic_2SP],PARAMETER[ sta ndard_parallel_1,35.47],PARAMETER[standard_parallel_2,3 4.03 ],PARAMETER[latitude_of_origin,33.5],PARAMETER[central_ meri dian,-118],PARAMETER[false_easting,6561666.6],PARAMETER[fa lse_ northing,1640416.7],UNIT[US survey foot,0.3048006096012192,AUTHORITY[EPSG,9003]],AUTHORITY[EPSG,2 6945 ]] And the exportToProj4() returns this: +proj=lcc +lat_1=35.47 +lat_2=34.03 ++lat_0=33.5 +lon_0=-118 +x_0=200 +y_0=50.01 +datum=NAD83 ++units=us-ft no_defs Which produces a correct result when using TransformEx(). When I get the SRS_PP_FALSE_EASTING and SRS_PP_FALSE_NORTHING from both projections using the OGRSpatialReference::GetProjParm() function I notice the difference between the two is a factor of the GEOGCS unit radians: SOURCE IMAGE: false_easting: 21527734.7 false_northing: 5381933.680556 TRANSLATED IMAGE: false_easting: 6561666.67 false_northing: 1640416.67 The translated offsets divided by the source offsets seem to be a factor of the unit radians: 0.3048 I think the reason that the translate to geotiff corrects the problem is that the geotiff keys for the projection are calculated in a different manner than when using the OGRCoordinateTransform object. Any help on where to apply the fix in the GDAL code if you think one is required would be much appreciated. If you think I need to accommodate this in my code any help would be greatly appreciated. Best regards, Martin -- Geospatial professional services http://even.rouault.free.fr/services.html ___ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev
