[gdal-dev] [PATCH v2] Support Mercator_2SP in GeoTIFF
Antti Castrén
antti.castren at iki.fi
Fri Nov 22 14:18:27 PST 2013
Hi Trent et al.
Thanks for your work on the issue.
I tested the files, and results were not satisfactory. But at least we
got more information on the issue.
FME opened both files.
The georeferencing of 18440_lowres_nol.tif was recognized by FME as:
DESC_NM: Global Mercator
DT_NAME: NAD83
PARM2: 47.667
PROJ: MRCAT
UNIT: METER
This should place the file at correct location.
The georeferencing of 18440_lowres_sp1.tif was recognized by FME as:
DESC_NM: Global Mercator
DT_NAME: NAD83
PROJ: MRCAT
UNIT: METER
To compare to the first file which you provided last week, its FME
parameters were:
DESC_NM: Global Mercator
DT_NAME: NAD83
PROJ: MRCAT
UNIT: METER
It seems FME doesn't read StdParallel1 or at least it doesn't use it.
Geomedia did NOT recognize either of the files as valid GeoTIFF.
To get projection (and of course location) right GeoMedia reads following keys:
ProjCoordTransGeoKey (Short,1): CT_Mercator
ProjLinearUnitsGeoKey (Short,1): Linear_Meter
ProjNatOriginLongGeoKey (Double,1): 0
ProjNatOriginLatGeoKey (Double,1): 0
ProjFalseEastingGeoKey (Double,1): 0
ProjFalseNorthingGeoKey (Double,1): 0
ProjCenterLatGeoKey (Double,1): 60.7666667
ProjScaleAtNatOriginGeoKey (Double,1): 0.48961699
These keys are (in my humble opinion) logical and the most correct
ones. Of course there is some redundancy. Unfortunately FME and ArcGIS
don't read last two at all, and therefore they place the file as if it
was projected from the Equator with scale 0.
ArcGIS opened all your three files at correct location.
By the way, I have seen that many sofware which are able to
interchange latitude of true scale and scalefactor at the Equator fail
to do the calculation with enough decimal places. My experience is
that you need at least 8 decimal places in scale factor to align your
raster data correctly. (As seen in the GeoMedia example above)
Br,
Antti
2013/11/20 Trent Piepho <tpiepho at gmail.com>:
> Antii,
>
> Thank you for the testing. It looks like Arc/GIS works while FME and
> Intergraph don't.
>
> That file had both StdParallel1 and ScaleAtNatOrigin defined.
>
> I've uploaded two new files. This one has only StdParallel1 defined with no
> ScaleAtNatOrigin.
> https://drive.google.com/file/d/0BybuTedE9CLxRnNlNUl2Zm5JRzg/edit?usp=sharing
>
> This file has the lat_ts stored in NatOriginLat and no ScaleAtNatOrigin
> defined.
> https://drive.google.com/file/d/0BybuTedE9CLxMzltYlBFWU51RzA/edit?usp=sharing
>
> I think the latter might have better compatibility. Yet it seems to clearly
> be the incorrect way to do it.
>
>
> On Mon, Nov 18, 2013 at 4:54 AM, Antti Castrén <antti.castren at iki.fi> wrote:
>>
>> Hi Trent,
>>
>> The chart opens in ArcMap (ArcGIS 10.0) well, and it is in right
>> location (Seattle).
>>
>> Relevant properties of the file as seen by ArcGIS:
>> Cell Size: 240.003787, 240.003787
>> Extent
>> Top:4143530.20898
>> Left: -9278434.53415
>> Right: -9165872.75807
>> Bottom: 3984167.69444
>> Spatial Reference: Global Mercator
>> Linear Unit: Meter
>> false_easting: 0
>> false_northing: 0
>> central_meridian: 0
>> standard_parallel_1: 47.667
>> Datum: D_North_American_1983
>>
>>
>> I opened the file in FME Data Inspector (FME Desktop 2013 SP1) also.
>> It did not read the latitude of true scale from the
>> ProjStdParallel1GeoKey as you can see from the following parameters.
>> Therefore it places the chart to NE Georgia/NW South Carolina as if
>> the latitude of true scale was the Equator.
>>
>> <FME Parameters>
>> Coordinate System Parameters
>> CS_NAME: _FME_0
>> DESC_NM: Global Mercator
>> DT_NAME: NAD83
>> PROJ: MRCAT
>> UNIT: METER
>>
>> Datum Parameters
>> DESC_NM: NAD 1983, Alaska, Canada,Continental US,Mexico,Central America
>> ELLIPSOID: GRS1980
>> SOURCE: US Defense Mapping Agency, TR-8350.2-B, August 1993
>> USE: NAD83
>>
>> Ellipsoid Parameters
>> DESC_NM: Geodetic Reference System of 1980
>> E_RAD: 6378137
>> P_RAD: 6356752.3141403478
>> SOURCE: Stem, L.E., Jan 1989, State Plane Coordinate System of 1983
>>
>> OGC WKT Description
>> PROJCS["Global Mercator",
>> GEOGCS["NAD83",
>> DATUM["North_American_Datum_1983",
>> SPHEROID["Geodetic Reference System of
>> 1980",6378137,298.2572221008916,
>> AUTHORITY["EPSG","7019"]],
>> AUTHORITY["EPSG","6269"]],
>> PRIMEM["Greenwich",0],
>> UNIT["degree",0.0174532925199433],
>> AUTHORITY["EPSG","4269"]],
>> PROJECTION["Mercator_1SP"],
>> PARAMETER["central_meridian",0],
>> PARAMETER["scale_factor",1],
>> PARAMETER["false_easting",0],
>> PARAMETER["false_northing",0],
>> UNIT["METER",1]]
>>
>> Esri WKT Description
>>
>> PROJCS["Global_Mercator",GEOGCS["GCS_North_American_1983",DATUM["D_North_American_1983",SPHEROID["Geodetic_Reference_System_of_1980",6378137,298.2572221008916]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Mercator"],PARAMETER["central_meridian",0],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["Meter",1],PARAMETER["standard_parallel_1",0.0]]
>> </FME Parameters>
>>
>> Intergraph's Geomedia Professional recognises the file as GeoTIFF, but
>> the chart appears to be in the same false place as with FME.
>>
>> If you want, I could update the file's georeferencing by FME and
>> Geomedia so that location is correct. Then I could send listgeo
>> information or even put the actual files available if desired.
>>
>> The goal is, of course, to get various software vendors to use the
>> GeoTIFF-information the same way. Hopefully the common way would also
>> be logically correct, but I guess that's too much to ask.
>>
>> Best regards,
>>
>> Antti
>>
>>
>> Here is the listgeo info from the original file for those who are
>> interested:
>>
>> Geotiff_Information:
>> Version: 1
>> Key_Revision: 1.0
>> Tagged_Information:
>> ModelTiepointTag (2,3):
>> 0 0 0
>> -9278434.53 4143530.21 0
>> ModelPixelScaleTag (1,3):
>> 240.003787 240.003787 0
>> End_Of_Tags.
>> Keyed_Information:
>> GTModelTypeGeoKey (Short,1): ModelTypeProjected
>> GTRasterTypeGeoKey (Short,1): RasterPixelIsArea
>> GTCitationGeoKey (Ascii,16): "Global Mercator"
>> GeographicTypeGeoKey (Short,1): GCS_NAD83
>> GeogCitationGeoKey (Ascii,6): "NAD83"
>> GeogAngularUnitsGeoKey (Short,1): Angular_Degree
>> GeogSemiMajorAxisGeoKey (Double,1): 6378137
>> GeogInvFlatteningGeoKey (Double,1): 298.257222
>> Unknown-2062 (Double,3): 0 0 0
>> ProjectedCSTypeGeoKey (Short,1): User-Defined
>> ProjectionGeoKey (Short,1): User-Defined
>> ProjCoordTransGeoKey (Short,1): CT_Mercator
>> ProjLinearUnitsGeoKey (Short,1): Linear_Meter
>> ProjStdParallel1GeoKey (Double,1): 47.667
>> ProjNatOriginLongGeoKey (Double,1): 0
>> ProjNatOriginLatGeoKey (Double,1): 0
>> ProjFalseEastingGeoKey (Double,1): 0
>> ProjFalseNorthingGeoKey (Double,1): 0
>> ProjScaleAtNatOriginGeoKey (Double,1): 1
>> End_Of_Keys.
>> End_Of_Geotiff.
>> Projection Method: CT_Mercator
>> ProjNatOriginLatGeoKey: 0.000000 ( 0d 0' 0.00"N)
>> ProjNatOriginLongGeoKey: 0.000000 ( 0d 0' 0.00"E)
>> ProjScaleAtNatOriginGeoKey: 1.000000
>> ProjFalseEastingGeoKey: 0.000000 m
>> ProjFalseNorthingGeoKey: 0.000000 m
>> GCS: 4269/NAD83
>> Datum: 6269/North American Datum 1983
>> Ellipsoid: 7019/GRS 1980 (6378137.00,6356752.31)
>> Prime Meridian: 8901/Greenwich (0.000000/ 0d 0' 0.00"E)
>> Projection Linear Units: 9001/metre (1.000000m)
>> Corner Coordinates:
>> Upper Left (-9278434.534, 4143530.209)
>> Lower Left (-9278434.534, 3984167.694)
>> Upper Right (-9165872.758, 4143530.209)
>> Lower Right (-9165872.758, 3984167.694)
>> Center (-9222153.646, 4063848.952)
>>
>>
>> 2013/11/14 Trent Piepho <tpiepho at gmail.com>:
>> > It looks like arcgis does not support Mercator_1SP in geotiff files.
>> > It just ignores the scale factor and uses 1.0. Which means in GDAL
>> > and libgeotiff as they are now, there is no way to export Mercator
>> > projections that don't have a scale of 1.0 to arcgis.
>> >
>> > ArcGIS does support Mercator_2SP, but instead of placing the latitude
>> > of true scale in ProjStdParallel1GeoKey, it uses
>> > ProjNatOriginLatGeoKey. Which seems to be plainly wrong, yet that's
>> > what it does.
>> >
>> > I wonder if ESRI used libgeotiff? With unpatched libgeotiff, the only
>> > way to get a latitude from a Mercator geotiff it to use one of the
>> > origin latitude GeoKeys. It's not correct, but it does allow one to
>> > stick the number in the file and get it back out without having to
>> > modify or understand any libgeotiff code.
>> >
>> > It's possible that ArcGIS supports getting the latitude of true scale
>> > from ProjStdParallel1GeoKey, even it it doesn't store it there.
>> > Antii, could you test this file with ArcGIS?
>> >
>> >
>> > https://drive.google.com/file/d/0BybuTedE9CLxak9iakZRTC00cHM/edit?usp=sharing
>> >
>> > Pick file->download to get the geotiff file.
>> >
>> > It should have projection information like this:
>> > PROJECTION["Mercator_2SP"],
>> > PARAMETER["standard_parallel_1",47.667],
>> > PARAMETER["central_meridian",0],
>> > PARAMETER["false_easting",0],
>> > PARAMETER["false_northing",0],
>> > UNIT["metre",1,
>> > AUTHORITY["EPSG","9001"]]]
>> > Origin = (-9278434.534151820465922,4143530.208984711207449)
>> > Pixel Size = (240.003786961970206,-240.003786961970206)
>> > Corner Coordinates:
>> > Upper Left (-9278434.534, 4143530.209) (123d32'26.22"W, 48d23'56.59"N)
>> > Lower Left (-9278434.534, 3984167.694) (123d32'26.22"W, 46d57'58.64"N)
>> > Upper Right (-9165872.758, 4143530.209) (122d 2'30.76"W, 48d23'56.59"N)
>> > Lower Right (-9165872.758, 3984167.694) (122d 2'30.76"W, 46d57'58.64"N)
>> > Center (-9222153.646, 4063848.952) (122d47'28.49"W, 47d41'15.43"N)
>
>
--
Antti Castrén
+358505997538
antti.castren at iki.fi
More information about the gdal-dev
mailing list