[gdal-dev] ERROR 1: No inverse operation

Even Rouault even.rouault at spatialys.com
Wed Mar 8 09:30:58 PST 2023


The file is indeed fine and was recognized in earlier GDAL versions (<= 
3.4). Bug fixed per https://github.com/OSGeo/gdal/pull/7382


Le 08/03/2023 à 14:22, Micha Silver a écrit :
>
> Inline...
>
>
> On 08/03/2023 14:33, Even Rouault wrote:
>>
>> Ah, the issue here is that the GeoTIFF keys in the .tif partially 
>> define the projection (or possibly the .tif.aux.xml contain this 
>> invalid WKT). But METHOD["unnamed"] means that the projection method 
>> is unknown, so it is logical that GDAL/PROJ can't compute geographic 
>> coordinates from projected coordinates.
>>
>>
>> The output of "listgeo HiBar_21022022_A1.tif" might be useful to 
>> further diagnose.
>>
>>
> I can't see any problem here:
>
>
> listgeo HiBar_21022022_A1.tif
> Geotiff_Information:
>   Version: 1
>   Key_Revision: 1.0
>   Tagged_Information:
>      ModelTiepointTag (2,3):
>         0                 0                 0
>         202245.240786848  421915.720555977  0
>      ModelPixelScaleTag (1,3):
>         0.125             0.125             1
>      End_Of_Tags.
>   Keyed_Information:
>      GTModelTypeGeoKey (Short,1): ModelTypeProjected
>      GTRasterTypeGeoKey (Short,1): RasterPixelIsArea
>      GTCitationGeoKey (Ascii,26): "PCS Name = Israel_TM_Grid"
>      GeographicTypeGeoKey (Short,1): Code-4141 (Israel 1993)
>      GeogAngularUnitsGeoKey (Short,1): Angular_Degree
>      GeogEllipsoidGeoKey (Short,1): Ellipse_GRS_1980
>      GeogSemiMajorAxisGeoKey (Double,1): 6378137
>      GeogSemiMinorAxisGeoKey (Double,1): 6356752.31414036
>      GeogInvFlatteningGeoKey (Double,1): 298.257222101
>      GeogTOWGS84GeoKey (Double,7): -48              55               52
> 0                0                0
> 0
>      ProjectedCSTypeGeoKey (Short,1): Code-2039 (Israel 1993 / Israeli 
> TM Grid)
>      PCSCitationGeoKey (Ascii,15): "Israel_TM_Grid"
>      ProjLinearUnitsGeoKey (Short,1): Linear_Meter
>      End_Of_Keys.
>   End_Of_Geotiff.
>
> PCS = 2039 (Israel 1993 / Israeli TM Grid)
> Projection = 18204 (Israeli TM)
> Projection Method: CT_TransverseMercator
>   ProjNatOriginLatGeoKey: 31.734394 ( 31d44' 3.82"N)
>   ProjNatOriginLongGeoKey: 35.204517 ( 35d12'16.26"E)
>   ProjScaleAtNatOriginGeoKey: 1.000007
>   ProjFalseEastingGeoKey: 219529.584000 m
>   ProjFalseNorthingGeoKey: 626907.390000 m
> GCS: 4141/Israel 1993
> Datum: 6141/Israel 1993
> Ellipsoid: 7019/GRS 1980 (6378137.00,6356752.31)
> Prime Meridian: 8901/Greenwich (0.000000/  0d 0' 0.00"E)
> TOWGS84: -48,55,52,0,0,0,0
> Projection Linear Units: 9001/metre (1.000000m)
>
> Corner Coordinates:
> Upper Left    (  202245.241,  421915.721)  ( 35d 1'32.11"E, 29d53' 
> 7.03"N)
> Lower Left    (  202245.241,  418790.721)  ( 35d 1'32.29"E, 
> 29d51'25.54"N)
> Upper Right   (  205370.241,  421915.721)  ( 35d 3'28.57"E, 29d53' 
> 7.17"N)
> Lower Right   (  205370.241,  418790.721)  ( 35d 3'28.72"E, 
> 29d51'25.68"N)
> Center        (  203807.741,  420353.221)  ( 35d 2'30.42"E, 29d52'16.36"N)
>
> Then, after:
>
>
> gdal_translate -a_srs EPSG:2039 HiBar_21022022_A1.tif 
> HiBar_21022022_A1b.tif
>
>
> The ERROR goes away and I get:
>
>
> gdalinfo HiBar_21022022_A1b.tif
> Driver: GTiff/GeoTIFF
> Files: HiBar_21022022_A1b.tif
> Size is 25000, 25000
> Coordinate System is:
> PROJCRS["Israel 1993 / Israeli TM Grid",
>    BASEGEOGCRS["Israel 1993",
>        DATUM["Israel 1993",
>            ELLIPSOID["GRS 1980",6378137,298.257222101,
>                LENGTHUNIT["metre",1]]],
>        PRIMEM["Greenwich",0,
>            ANGLEUNIT["degree",0.0174532925199433]],
>        ID["EPSG",4141]],
>    CONVERSION["Israeli TM",
>        METHOD["Transverse Mercator",
>            ID["EPSG",9807]],
>        PARAMETER["Latitude of natural origin",31.7343936111111,
>            ANGLEUNIT["degree",0.0174532925199433],
>            ID["EPSG",8801]],
>        PARAMETER["Longitude of natural origin",35.2045169444444,
>            ANGLEUNIT["degree",0.0174532925199433],
>            ID["EPSG",8802]],
>        PARAMETER["Scale factor at natural origin",1.0000067,
>            SCALEUNIT["unity",1],
>            ID["EPSG",8805]],
>        PARAMETER["False easting",219529.584,
>            LENGTHUNIT["metre",1],
>            ID["EPSG",8806]],
>        PARAMETER["False northing",626907.39,
>            LENGTHUNIT["metre",1],
>            ID["EPSG",8807]]],
>    CS[Cartesian,2],
>        AXIS["(E)",east,
>            ORDER[1],
>            LENGTHUNIT["metre",1]],
>        AXIS["(N)",north,
>            ORDER[2],
>            LENGTHUNIT["metre",1]],
>    USAGE[
>        SCOPE["Cadastre, engineering survey, topographic mapping (large 
> and medium scale)."],
>        AREA["Israel - onshore; Palestine Territory - onshore."],
>        BBOX[29.45,34.17,33.28,35.69]],
>    ID["EPSG",2039]]
> Data axis to CRS axis mapping: 1,2
> Origin = (202245.240786847687559,421915.720555976789910)
> Pixel Size = (0.125000000000000,-0.125000000000000)
> Metadata:
>  AREA_OR_POINT=Area
> Image Structure Metadata:
>  INTERLEAVE=PIXEL
> Corner Coordinates:
> Upper Left  (  202245.241,  421915.721) ( 35d 1'32.11"E, 29d53' 7.03"N)
> Lower Left  (  202245.241,  418790.721) ( 35d 1'32.29"E, 29d51'25.54"N)
> Upper Right (  205370.241,  421915.721) ( 35d 3'28.57"E, 29d53' 7.17"N)
> Lower Right (  205370.241,  418790.721) ( 35d 3'28.72"E, 29d51'25.68"N)
> Center      (  203807.741,  420353.221) ( 35d 2'30.42"E, 29d52'16.36"N)
> Band 1 Block=25000x1 Type=Byte, ColorInterp=Red
> Band 2 Block=25000x1 Type=Byte, ColorInterp=Green
> Band 3 Block=25000x1 Type=Byte, ColorInterp=Blue
> Band 4 Block=25000x1 Type=Byte, ColorInterp=Undefined
>
> and the listgeo output is also different:
>
>
> listgeo HiBar_21022022_A1b.tif
> Geotiff_Information:
>   Version: 1
>   Key_Revision: 1.0
>   Tagged_Information:
>      ModelTiepointTag (2,3):
>         0                 0                 0
>         202245.240786848  421915.720555977  0
>      ModelPixelScaleTag (1,3):
>         0.125             0.125             0
>      End_Of_Tags.
>   Keyed_Information:
>      GTModelTypeGeoKey (Short,1): ModelTypeProjected
>      GTRasterTypeGeoKey (Short,1): RasterPixelIsArea
>      GTCitationGeoKey (Ascii,30): "Israel 1993 / Israeli TM Grid"
>      GeogCitationGeoKey (Ascii,12): "Israel 1993"
>      GeogAngularUnitsGeoKey (Short,1): Angular_Degree
>      ProjectedCSTypeGeoKey (Short,1): Code-2039 (Israel 1993 / Israeli 
> TM Grid)
>      ProjLinearUnitsGeoKey (Short,1): Linear_Meter
>      End_Of_Keys.
>   End_Of_Geotiff.
>
> PCS = 2039 (Israel 1993 / Israeli TM Grid)
> Projection = 18204 (Israeli TM)
> Projection Method: CT_TransverseMercator
>   ProjNatOriginLatGeoKey: 31.734394 ( 31d44' 3.82"N)
>   ProjNatOriginLongGeoKey: 35.204517 ( 35d12'16.26"E)
>   ProjScaleAtNatOriginGeoKey: 1.000007
>   ProjFalseEastingGeoKey: 219529.584000 m
>   ProjFalseNorthingGeoKey: 626907.390000 m
> GCS: 4141/Israel 1993
> Datum: 6141/Israel 1993
> 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    (  202245.241,  421915.721)  ( 35d 1'32.11"E, 29d53' 
> 7.03"N)
> Lower Left    (  202245.241,  418790.721)  ( 35d 1'32.29"E, 
> 29d51'25.54"N)
> Upper Right   (  205370.241,  421915.721)  ( 35d 3'28.57"E, 29d53' 
> 7.17"N)
> Lower Right   (  205370.241,  418790.721)  ( 35d 3'28.72"E, 
> 29d51'25.68"N)
> Center        (  203807.741,  420353.221)  ( 35d 2'30.42"E, 29d52'16.36"N)
>
>
> So I guess that solves the problem, allowing for reprojection, and 
> mapping together with other datasets.
>
> Thanks for you help (spot on, as always!)
>
>
>
>> You likely need to gdal_translate -a_srs EPSG:2039 to do something 
>> useful with that file, and report to the data producer that the CRS 
>> encoding is broken
>>
>>
>> Le 08/03/2023 à 13:18, Micha Silver a écrit :
>>>
>>> Hi Even:
>>>
>>>
>>> Here's the full gdalinfo, with the reported CRS:
>>>
>>>
>>> gdalinfo HiBar_21022022_A1.tif
>>> Driver: GTiff/GeoTIFF
>>> Files: HiBar_21022022_A1.tif
>>>       HiBar_21022022_A1.tif.aux.xml
>>> Size is 25000, 25000
>>> Coordinate System is:
>>> PROJCRS["Israel 1993 / Israeli TM Grid",
>>>    BASEGEOGCRS["WGS 84",
>>>        DATUM["World Geodetic System 1984",
>>>            ELLIPSOID["WGS 84",6378137,298.257223563,
>>>                LENGTHUNIT["metre",1,
>>>                    ID["EPSG",9001]]]],
>>>        PRIMEM["Greenwich",0,
>>>            ANGLEUNIT["degree",0.0174532925199433,
>>>                ID["EPSG",9122]]]],
>>>    CONVERSION["unnamed",
>>>        METHOD["unnamed"]],
>>>    CS[Cartesian,2],
>>>        AXIS["(E)",east,
>>>            ORDER[1],
>>>            LENGTHUNIT["metre",1,
>>>                ID["EPSG",9001]]],
>>>        AXIS["(N)",north,
>>>            ORDER[2],
>>>            LENGTHUNIT["metre",1,
>>>                ID["EPSG",9001]]]]
>>> Data axis to CRS axis mapping: 1,2
>>> Origin = (202245.240786847687559,421915.720555976789910)
>>> Pixel Size = (0.125000000000000,-0.125000000000000)
>>> Metadata:
>>>  AREA_OR_POINT=Area
>>> Image Structure Metadata:
>>>  INTERLEAVE=PIXEL
>>> Corner Coordinates:
>>> ERROR 1: No inverse operation
>>> Upper Left  (  202245.241,  421915.721)
>>> ERROR 1: No inverse operation
>>> Lower Left  (  202245.241,  418790.721)
>>> ERROR 1: No inverse operation
>>> Upper Right (  205370.241,  421915.721)
>>> ERROR 1: No inverse operation
>>> Lower Right (  205370.241,  418790.721)
>>> ERROR 1: No inverse operation
>>> Center      (  203807.741,  420353.221)
>>>
>>> and proj version:
>>>
>>>
>>> proj
>>> Rel. 9.1.1, December 1st, 2022
>>>
>>>
>>> (a newly installed debian 12)
>>>
>>>
>>> Thanks,
>>>
>>> Micha
>>>
>>>
>>> On 08/03/2023 13:47, Even Rouault wrote:
>>>>
>>>> Micha,
>>>>
>>>>
>>>> which CRS is reported and which PROJ version do you use ? It is 
>>>> presumably one projection method for which PROJ only implements the 
>>>> forward conversion, that is from geographic to projected 
>>>> coordinates, but not the reverse way
>>>>
>>>>
>>>> Even
>>>>
>>>>
>>>> Le 08/03/2023 à 09:37, Micha Silver a écrit :
>>>>>
>>>>> I received some Geotiff files, and gdalinfo shows:
>>>>>
>>>>>
>>>>> Corner Coordinates:
>>>>> ERROR 1: No inverse operation
>>>>> Upper Left  (  205370.241,  418790.721)
>>>>> ERROR 1: No inverse operation
>>>>> Lower Left  (  205370.241,  415665.721)
>>>>> ERROR 1: No inverse operation
>>>>> Upper Right (  206198.866,  418790.721)
>>>>> ERROR 1: No inverse operation
>>>>> Lower Right (  206198.866,  415665.721)
>>>>> ERROR 1: No inverse operation
>>>>> Center      (  205784.553,  417228.221)
>>>>>
>>>>> What does that mean?
>>>>>
>>>>>
>>>>> gdalinfo --version
>>>>> GDAL 3.6.2, released 2023/01/02
>>>>>
>>>>> Thanks
>>>>>
>>>>>
>>>>> -- 
>>>>> Micha Silver
>>>>> Remote Sensing Lab, Sde Boker
>>>>> Ben Gurion University
>>>>> +972-523-665918
>>>>>
>>>>> _______________________________________________
>>>>> gdal-dev mailing list
>>>>> gdal-dev at lists.osgeo.org
>>>>> https://lists.osgeo.org/mailman/listinfo/gdal-dev
>>>> -- 
>>>> http://www.spatialys.com
>>>> My software is free, but my time generally not.
>>> -- 
>>> Micha Silver
>>> Remote Sensing Lab, Sde Boker
>>> Ben Gurion University
>>> +972-523-665918
>> -- 
>> http://www.spatialys.com
>> My software is free, but my time generally not.
> -- 
> Micha Silver
> Remote Sensing Lab, Sde Boker
> Ben Gurion University
> +972-523-665918

-- 
http://www.spatialys.com
My software is free, but my time generally not.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/gdal-dev/attachments/20230308/9471d84d/attachment-0001.htm>


More information about the gdal-dev mailing list