[gdal-dev] map info and projections from GeoTIFF file to create ENVI header

nguler at jlab.org nguler at jlab.org
Mon Jan 18 09:27:50 PST 2016


It is solved. Thank you very much Even and Dmitry for your helps. I have
missed the GDAL_DATA variable and it caused me all that trouble.

Now I need to get the map info string:
map info = {UTM, 1, 1, 145185, 2189715, 30, 30, 19, North,WGS-84}

Using hDS->GetGeoTransform(adfGeoTransform) does not give me all the
information in the map info string. It only gives me elements [4-7], how to
get the others?

Best regards, - Nevzat

> Yes, Thank you Even, I was missing the GDAL_DATA environmental variable. Now I
> get the same gdalinfo output. Do I need to set this variable while using the
> C++ classes you suggested before: GetGeoTransform and GetProjectionRef()..?
> Because, there also I had the same problem.
>
> Best regards, - Nevzat
>
>
>> Le lundi 18 janvier 2016 17:13:14, nguler at jlab.org a écrit :
>>> Dmitry and All,
>>> I just downloaded the set of tiff files from earthexplorer and run gdalinfo
>>> on it. Please see below. I am missing some fields in the output
>>> information. It is different from what you got. I don't know what can
>>> cause this? My gdalinfo.exe is not good? Any idea?
>>
>> Have a look at the FAQ entries about GDAL_DATA :
>> https://trac.osgeo.org/gdal/wiki/FAQInstallationAndBuilding#WhatisGDAL_DATAenvironmentvariable
>>
>>>
>>> C:\Cubes\LC80080472013174LGN00.tar>gdalinfo LC80080472013174LGN00_B1.TIF
>>> Driver: GTiff/GeoTIFF
>>> Files: LC80080472013174LGN00_B1.TIF
>>>        .\LC80080472013174LGN00_MTL.txt
>>> Size is 7571, 7411
>>> Coordinate System is:
>>> PROJCS["unnamed",
>>>     PROJECTION["Transverse_Mercator"],
>>>     PARAMETER["latitude_of_origin",0],
>>>     PARAMETER["central_meridian",-69],
>>>     PARAMETER["scale_factor",0.9996],
>>>     PARAMETER["false_easting",500000],
>>>     PARAMETER["false_northing",0],
>>>     UNIT["metre",1,
>>>         AUTHORITY["EPSG","9001"]]]
>>> Origin = (145185.000000000000000,2189715.000000000000000)
>>> Pixel Size = (30.000000000000000,-30.000000000000000)
>>> Metadata:
>>>   AREA_OR_POINT=Point
>>>   METADATATYPE=ODL
>>> Image Structure Metadata:
>>>   INTERLEAVE=BAND
>>> Corner Coordinates:
>>> Upper Left  (  145185.000, 2189715.000)
>>> Lower Left  (  145185.000, 1967385.000)
>>> Upper Right (  372315.000, 2189715.000)
>>> Lower Right (  372315.000, 1967385.000)
>>> Center      (  258750.000, 2078550.000)
>>> Band 1 Block=7571x1 Type=UInt16, ColorInterp=Gray
>>>
>>> C:\Cubes\LC80080472013174LGN00.tar>gdalinfo --version
>>> GDAL 2.0.1, released 2015/09/15
>>>
>>> C:\Cubes\LC80080472013174LGN00.tar>
>>>
>>> > I have the 2.0.1 version:
>>> > C:\warmerda\bld\bin>gdalinfo --version
>>> > GDAL 2.0.1, released 2015/09/15
>>> >
>>> > I will check the history of the file I have. Maybe it is modified. I will
>>> > also download the same file from earthexplorer and try my routines on
>>> > it. Thanks for the suggestions. I appreciate.
>>> > - Nevzat
>>> >
>>> > -----Original Message-----
>>> > From: Dmitry Baryshnikov [mailto:bishop.dev at gmail.com]
>>> > Sent: Sunday, January 17, 2016 12:57 PM
>>> > To: Nevzat Guler <nguler at jlab.org>; gdal-dev at lists.osgeo.org
>>> > Subject: Re: [gdal-dev] map info and projections from GeoTIFF file to
>>> > create ENVI header
>>> >
>>> > I downloaded this scene from earthexplorer, and everything is fine:
>>> >
>>> > $ gdalinfo LC80080472013174LGN00_B1.TIF
>>> > Driver: GTiff/GeoTIFF
>>> > Files: LC80080472013174LGN00_B1.TIF
>>> >
>>> >         ./LC80080472013174LGN00_MTL.txt
>>> >
>>> > Size is 7571, 7411
>>> > Coordinate System is:
>>> > PROJCS["WGS 84 / UTM zone 19N",
>>> >
>>> >      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","9122"]],
>>> >
>>> >          AUTHORITY["EPSG","4326"]],
>>> >
>>> >      PROJECTION["Transverse_Mercator"],
>>> >      PARAMETER["latitude_of_origin",0],
>>> >      PARAMETER["central_meridian",-69],
>>> >      PARAMETER["scale_factor",0.9996],
>>> >      PARAMETER["false_easting",500000],
>>> >      PARAMETER["false_northing",0],
>>> >      UNIT["metre",1,
>>> >
>>> >          AUTHORITY["EPSG","9001"]],
>>> >
>>> >      AXIS["Easting",EAST],
>>> >      AXIS["Northing",NORTH],
>>> >      AUTHORITY["EPSG","32619"]]
>>> >
>>> > Origin = (145185.000000000000000,2189715.000000000000000)
>>> > Pixel Size = (30.000000000000000,-30.000000000000000)
>>> >
>>> > Metadata:
>>> >    AREA_OR_POINT=Point
>>> >    METADATATYPE=ODL
>>> >
>>> > Image Structure Metadata:
>>> >    INTERLEAVE=BAND
>>> >
>>> > Corner Coordinates:
>>> > Upper Left  (  145185.000, 2189715.000) ( 72d23' 7.97"W, 19d46'16.43"N)
>>> > Lower Left  (  145185.000, 1967385.000) ( 72d20'44.51"W, 17d45'55.15"N)
>>> > Upper Right (  372315.000, 2189715.000) ( 70d13' 8.50"W, 19d47'56.92"N)
>>> > Lower Right ( 372315.000, 1967385.000) ( 70d12'16.74"W, 17d47'24.74"N)
>>> > Center      (  258750.000, 2078550.000) ( 71d17'19.50"W, 18d47' 4.68"N)
>>> > Band 1 Block=7571x1 Type=UInt16, ColorInterp=Gray
>>> >
>>> > $ gdalinfo --version
>>> > GDAL 2.1.0dev, released 2015/99/99
>>> >
>>> > What is your gdal version? Maybe the Landsat geotiff was modified from
>>> > some software?
>>> >
>>> > Best regards,
>>> >
>>> >      Dmitry
>>> >
>>> > 17.01.2016 20:15, Nevzat Guler пишет:
>>> >> Thank you Dmitry,
>>> >> What is the best way to extract the information we got using
>>> >> gdalinfo.exe by using GDAL classes, within a C++ code? As I mentioned,
>>> >> I am using GetGeoTransform and GetProjectionRef() but the results miss
>>> >> many of the information (specifically the map info elements as I
>>> >> explained in the previous email).
>>> >>
>>> >> On another note, I actually get a different output from gdalinfo.exe
>>> >> using my file (please see the printed result below). I got the geotiff
>>> >> file from my colleague. When I open the file with ENVI, I get all the
>>> >> information, including GEOGCS information.
>>> >>
>>> >> Can it be that my file is just corrupted and that's why I don't get all
>>> >> the information from GetGeoTransform and GetProjectionRef()? But then
>>> >> how would ENVI is able to extract this information if it were not in
>>> >> the file?
>>> >>
>>> >> Thanks again, - Nevzat
>>> >>
>>> >> C:\Cubes\Landsat\LC80080472013174LGN00>gdalinfo
>>> >> LC80080472013174LGN00_B1.TIF
>>> >> Driver: GTiff/GeoTIFF
>>> >> Files: LC80080472013174LGN00_B1.TIF
>>> >>
>>> >>         .\LC80080472013174LGN00_MTL.txt Size is 7571, 7411 Coordinate
>>> >>
>>> >> System is:
>>> >> PROJCS["unnamed",
>>> >>
>>> >>      PROJECTION["Transverse_Mercator"],
>>> >>      PARAMETER["latitude_of_origin",0],
>>> >>      PARAMETER["central_meridian",-69],
>>> >>      PARAMETER["scale_factor",0.9996],
>>> >>      PARAMETER["false_easting",500000],
>>> >>      PARAMETER["false_northing",0],
>>> >>      UNIT["metre",1,
>>> >>
>>> >>          AUTHORITY["EPSG","9001"]]]
>>> >>
>>> >> Origin = (145185.000000000000000,2189715.000000000000000)
>>> >> Pixel Size = (30.000000000000000,-30.000000000000000)
>>> >>
>>> >> Metadata:
>>> >>    AREA_OR_POINT=Point
>>> >>    METADATATYPE=ODL
>>> >>
>>> >> Image Structure Metadata:
>>> >>    INTERLEAVE=BAND
>>> >>
>>> >> Corner Coordinates:
>>> >> Upper Left  (  145185.000, 2189715.000) Lower Left  (  145185.000,
>>> >> 1967385.000) Upper Right (  372315.000, 2189715.000) Lower Right (
>>> >> 372315.000, 1967385.000)
>>> >> Center      (  258750.000, 2078550.000)
>>> >> Band 1 Block=7571x1 Type=UInt16, ColorInterp=Gray
>>> >>
>>> >> -----Original Message-----
>>> >> From: gdal-dev [mailto:gdal-dev-bounces at lists.osgeo.org] On Behalf Of
>>> >> Dmitry Baryshnikov
>>> >> Sent: Sunday, January 17, 2016 6:57 AM
>>> >> To: gdal-dev at lists.osgeo.org
>>> >> Subject: Re: [gdal-dev] map info and projections from GeoTIFF file to
>>> >> create ENVI header
>>> >>
>>> >> Hi Nevzat,
>>> >>
>>> >> The landsat geotiff already have CRS embeded in file
>>> >>
>>> >> $ gdalinfo LC80150362013097LGN03_B1.TIF
>>> >> Driver: GTiff/GeoTIFF
>>> >> Files: LC80150362013097LGN03_B1.TIF
>>> >>
>>> >>          ./LC80150362013097LGN03_MTL.txt Size is 7281, 7261 Coordinate
>>> >>
>>> >> System is:
>>> >> PROJCS["WGS 84 / UTM zone 17N",
>>> >>
>>> >>       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","9122"]],
>>> >>
>>> >>           AUTHORITY["EPSG","4326"]],
>>> >>
>>> >>       PROJECTION["Transverse_Mercator"],
>>> >>       PARAMETER["latitude_of_origin",0],
>>> >>       PARAMETER["central_meridian",-81],
>>> >>       PARAMETER["scale_factor",0.9996],
>>> >>       PARAMETER["false_easting",500000],
>>> >>       PARAMETER["false_northing",0],
>>> >>       UNIT["metre",1,
>>> >>
>>> >>           AUTHORITY["EPSG","9001"]],
>>> >>
>>> >>       AXIS["Easting",EAST],
>>> >>       AXIS["Northing",NORTH],
>>> >>       AUTHORITY["EPSG","32617"]]
>>> >>
>>> >> Origin = (706185.000000000000000,3943515.000000000000000)
>>> >> Pixel Size = (30.000000000000000,-30.000000000000000)
>>> >>
>>> >> Metadata:
>>> >>     AREA_OR_POINT=Point
>>> >>     METADATATYPE=ODL
>>> >>
>>> >> Image Structure Metadata:
>>> >>     INTERLEAVE=BAND
>>> >>
>>> >> Corner Coordinates:
>>> >> Upper Left  (  706185.000, 3943515.000) ( 78d43'24.74"W, 35d36'50.15"N)
>>> >> Lower Left  (  706185.000, 3725685.000) ( 78d46'35.65"W, 33d39' 3.82"N)
>>> >> Upper Right (  924615.000, 3943515.000) ( 76d19' 2.05"W, 35d32'39.74"N)
>>> >> Lower Right (  924615.000, 3725685.000) ( 76d25'33.46"W, 33d35'10.99"N)
>>> >> Center      (  815400.000, 3834600.000) ( 77d33'38.95"W, 34d36'17.35"N)
>>> >> Band 1 Block=7281x1 Type=UInt16, ColorInterp=Gray
>>> >>
>>> >> And this is right CRS.
>>> >> Where do you get you geotiff?
>>> >>
>>> >> Best regards,
>>> >>
>>> >>       Dmitry
>>> >>
>>> >> 17.01.2016 03:41, Nevzat Guler пишет:
>>> >>> Dear experts,
>>> >>>
>>> >>> I am reading LANDSAT GeoTIFF files and combining bands to create and
>>> >>> ENVI cube.
>>> >>> I am using GDAL to read the GeoTIFF file, and my own custom cube
>>> >>> writer to write the ENVI style cube. Transferring the map information
>>> >>> proved to be a bit convoluted. I seek your help to figure out this.
>>> >>>
>>> >>> My first problem is to access map info from GeoTIFF file. The map
>>> >>> info
>>> >>> includes:
>>> >>>
>>> >>> map info = {UTM, 1.000, 1.000, 145185.000, 2189715.000,
>>> >>> 3.0000000000e+001, 3.0000000000e+001, 19, North, WGS-84,
>>> >>> units=Meters}
>>> >>>
>>> >>> I am able to get the location (x,y), and pixel size (x,y) from :
>>> >>> double adfGeoTransform[6];
>>> >>> srcDS->GetGeoTransform(adfGeoTransform);
>>> >>>
>>> >>> But, I could not find any way to get the other variables in the
>>> >>> string: UTM, pixel tie points, Unit and DATUM. Is there any way to get
>>> >>> them?
>>> >>>
>>> >>> Finally, for the Coordinate System String, GDAL gives me:
>>> >>>
>>> >>> PROJCS["unnamed",PROJECTION["Transverse_Mercator"],PARAMETER["latitud
>>> >>> e
>>> >>> _of_or
>>> >>> igin
>>> >>> ",0],PARAMETER["central_meridian",-69],PARAMETER["scale_factor",0.999
>>> >>> 6
>>> >>> ],PARA
>>> >>> METE
>>> >>> R["false_easting",500000],PARAMETER["false_northing",0],UNIT["Meter",
>>> >>> 1
>>> >>> ]]
>>> >>>
>>> >>> But I should actually get
>>> >>>
>>> >>> PROJCS["UTM_Zone_19N",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS
>>> >>> 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],TOWGS84[0,0,0,0,0
>>> >>> ,
>>> >>> 0,0],A
>>> >>> UTHO
>>> >>> RITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],U
>>> >>> N
>>> >>> IT["de
>>> >>> gree
>>> >>> ",0.0174532925199433,AUTHORITY["EPSG","9108"]],AUTHORITY["EPSG","4326"
>>> >>> ]],PRO
>>> >>> JECT
>>> >>> ION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETE
>>> >>> R
>>> >>> ["cent
>>> >>> ral_
>>> >>> meridian",-69],PARAMETER["scale_factor",0.9996],PARAMETER["false_east
>>> >>> i
>>> >>> ng",50
>>> >>> 0000
>>> >>> ],PARAMETER["false_northing",0],UNIT["Meter",1]]
>>> >>>
>>> >>> GDAL gives me PROJCS, "unnamed". And no GEOGCS information is given.
>>> >>> I am able to get the coordinate system properties by using the
>>> >>> following commands but I first need to set few things to use
>>> >>> exportToWTK.
>>> >>>
>>> >>>>>         OGRSpatialReference oSRS;
>>> >>>>>         char *pszWKT = NULL;
>>> >>>>>         oSRS.SetProjCS("UTM_Zone_19N");
>>> >>>>>         oSRS.SetWellKnownGeogCS("WGS84");
>>> >>>>>         oSRS.SetUTM(19, TRUE);
>>> >>>>>         oSRS.SetWellKnownGeogCS("WGS84");
>>> >>>>>         oSRS.exportToWkt(&pszWKT);
>>> >>>>>         printf("%s\n", pszWKT);
>>> >>>
>>> >>> However, this requires me to set SetProjCS, SetWellKnownGeogCS and
>>> >>> SetUTM first. I don't know these before I read the file. Is there any
>>> >>> way to get them from the file so that I don't need prior knowledge to
>>> >>> create the header. When I pull the GeoTIFF file into ENVI, all these
>>> >>> information is already there. If this list is only or development
>>> >>> purposes, please let me know if there is any other email list that
>>> >>> can give me guidance on this quest.
>>> >>>
>>> >>> I appreciate your help and suggestions.
>>> >>> Best regards,
>>> >>> - Nevzat
>>> >>>
>>> >>> _______________________________________________
>>> >>> gdal-dev mailing list
>>> >>> gdal-dev at lists.osgeo.org
>>> >>> http://lists.osgeo.org/mailman/listinfo/gdal-dev
>>> >>
>>> >> _______________________________________________
>>> >> gdal-dev mailing list
>>> >> gdal-dev at lists.osgeo.org
>>> >> http://lists.osgeo.org/mailman/listinfo/gdal-dev
>>> >
>>> > _______________________________________________
>>> > gdal-dev mailing list
>>> > gdal-dev at lists.osgeo.org
>>> > http://lists.osgeo.org/mailman/listinfo/gdal-dev
>>>
>>> Nevzat Guler
>>> Senior Scientist
>>> Spectral Sciences
>>> Phone: 781-273-4770
>>> http://www.spectral.com/index.shtml
>>> http://www.linkedin.com/in/nevzatguler
>>>
>>>
>>> _______________________________________________
>>> gdal-dev mailing list
>>> gdal-dev at lists.osgeo.org
>>> http://lists.osgeo.org/mailman/listinfo/gdal-dev
>>
>> --
>> Spatialys - Geospatial professional services
>> http://www.spatialys.com
>>
>
>
> Nevzat Guler
> Senior Scientist
> Spectral Sciences
> Phone: 781-273-4770
> http://www.spectral.com/index.shtml
> http://www.linkedin.com/in/nevzatguler
>
>


Nevzat Guler
Senior Scientist
Spectral Sciences
Phone: 781-273-4770
http://www.spectral.com/index.shtml
http://www.linkedin.com/in/nevzatguler




More information about the gdal-dev mailing list