[gdal-dev] conversion of imagery which has an empty coordinate
system to geotif format via gdalwarp
Vincent Schut
schut at sarvision.nl
Wed Feb 27 05:00:52 EST 2008
John Mitchell wrote:
> Hi,
>
> I would like to convert imagery which has an empty coordinate system to
> geotif format via gdalwarp:
John,
this is a file with gcp's, and the gcp's do have a projection defined.
You should just run gdalwarp on this file without any additional
parameters regarding coordinates, then gdalwarp should pick up the gcp's
and use those to warp the file to the crs of the gcp's, which indeed
looks like epsg 4326, but I'm no expert on that.
So just run gdalwarp like this: 'gdalwarp -rc <sourcefile> <dstfile>',
which should do the trick.
Cheers,
Vincent.
>
> Below is part of a gdalinfo dump of the original imagery. Notice that
> the Coordinate System is empty
> ----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
> Size is 12074, 13232
> Coordinate System is `'
> GCP Projection = GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS
> 84",6378137,298.257223563,AUTHORITY["EPSG",7030]],TOWGS84[0,0,0,0,0,0,0],AUTHORITY["EPSG",6326]],PRIMEM["Greenwich",0,AUTHORITY["EPSG",8901]],UNIT["DMSH",0.0174532925199433,AUTHORITY["EPSG",9108]],AXIS["Lat",NORTH],AXIS["Long",EAST],AUTHORITY["EPSG",4326]]
> GCP[ 0]: Id=1, Info=
> (0.5,0.5) -> (36.77378,-17.270641,0)
> GCP[ 1]: Id=2, Info=
> (6037,0.5) -> (35.398528,-16.938865,0)
> GCP[ 2]: Id=3, Info=
> (12073.5,0.5) -> (34.027899,-16.597882,0)
> GCP[ 3]: Id=4, Info=
> (0.5,3307.5) -> (36.58987,-17.99724,0)
> GCP[ 4]: Id=5, Info=
> (6037,3307.5) -> (35.209644,-17.663643,0)
> GCP[ 5]: Id=6, Info=
> (12073.5,3307.5) -> (33.834311,-17.320427,0)
> GCP[ 6]: Id=7, Info=
> (0.5,6614.5) -> (36.404318,-18.723611,0)
> GCP[ 7]: Id=8, Info=
> (6037,6614.5) -> (35.01887,-18.388159,0)
> GCP[ 8]: Id=9, Info=
> (12073.5,6614.5) -> (33.638591,-18.042675,0)
> GCP[ 9]: Id=10, Info=
> (0.5,9921.5) -> (36.21706,-19.449749,0)
> GCP[ 10]: Id=11, Info=
> (6037,9921.5) -> (34.826137,-19.112408,0)
> GCP[ 11]: Id=12, Info=
> (12073.5,9921.5) -> (33.440671,-18.764619,0)
> GCP[ 12]: Id=13, Info=
> (0.5,13228.5) -> (36.028031,-20.175626,0)
> GCP[ 13]: Id=14, Info=
> (6037,13228.5) -> (34.631381,-19.836361,0)
> GCP[ 14]: Id=15, Info=
> (12073.5,13228.5) -> (33.24048,-19.486228,0)
> Metadata:
> ELLIPSOID=WGS84
> Corner Coordinates:
> Upper Left ( 0.0, 0.0)
> Lower Left ( 0.0,13232.0)
> Upper Right (12074.0, 0.0)
> Lower Right (12074.0,13232.0)
> Center ( 6037.0, 6616.0)
> Band 1 Block=12074x1 Type=Byte, ColorInterp=Undefined
> -----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
>
>
> Based on the above gdalinfo data I assumed that the projection should be
> wgs84/EPSG:4326
>
> So I run the following command against the imagery to convert it to
> geotiff with an EPSG code of 4326:
> gdalwarp -rc -s_srs epsg:4326 -t_srs epsg:4326 source_imagery
> destination_geotiff
>
> Next I run gdalinfo against the converted geotiff file which is listed
> below:
> ------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
> Driver: GTiff/GeoTIFF
> Files: C:\dat_01.tif
> Size is 15174, 15369
> Coordinate System is:
> GEOGCS["WGS 84",
> DATUM["WGS_1984",
> SPHEROID["WGS 84",6378137,298.2572235630016,
> AUTHORITY["EPSG","7030"]],
> AUTHORITY["EPSG","6326"]],
> PRIMEM["Greenwich",0],
> UNIT["degree",0.0174532925199433],
> AUTHORITY["EPSG","4326"]]
> Origin = (33.240318507099211,-16.597638311748319)
> Pixel Size = (0.000232863689106,-0.000232863689106)
> Metadata:
> AREA_OR_POINT=Area
> Image Structure Metadata:
> INTERLEAVE=BAND
> Corner Coordinates:
> Upper Left ( 33.2403185, -16.5976383) ( 33d14'25.15"E, 16d35'51.50"S)
> Lower Left ( 33.2403185, -20.1765203) ( 33d14'25.15"E, 20d10'35.47"S)
> Upper Right ( 36.7737921, -16.5976383) ( 36d46'25.65"E, 16d35'51.50"S)
> Lower Right ( 36.7737921, -20.1765203) ( 36d46'25.65"E, 20d10'35.47"S)
> Center ( 35.0070553, -18.3870793) ( 35d 0'25.40"E, 18d23'13.49"S)
> Band 1 Block=15174x1 Type=Byte, ColorInterp=Gray
> -----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
>
> The readme file has the CORNER COORDINATES of:
> -----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
>
> 17°16'14.30" S 16°35'52.37" S
> 36°46'25.60" E 34°01'40.43" E
>
> 20°10'34.62" S 19°29'12.77" S
> 36°01'40.29" E 33°14'25.06" E
>
> CENTER 18°22'S 035°01'E
> ----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
>
> In comparing the readme file to converted geotiff the resolution
> increased slightly and the coordinates were slightly different.
> Am I correct in assuming that the EPSG code should be 4326 when
> converting to geotiff format via gdalwarp?
>
> Thanks,
>
> John J. Mitchell
>
>
> ------------------------------------------------------------------------
>
> _______________________________________________
> gdal-dev mailing list
> gdal-dev at lists.osgeo.org
> http://lists.osgeo.org/mailman/listinfo/gdal-dev
More information about the gdal-dev
mailing list