otbcli_ReadImageInfo does not write info to external geom file

OTB Usage
1

Configuration setup

My system: Ubuntu 18.04
Version of the OTB: 7.2

If relevant, you may also provide:
QGIS version: 3.16
Python version: 3.7

This issue had been previously raised here and here.

Simply put, extracting image metadata as mentioned in the cookbook here does not work.

Here are the logs

$ otbcli_ReadImageInfo -in 2014_test.tif -outkwl TheGeom.geom

2021-05-17 15:30:13 (INFO) ReadImageInfo: Default RAM limit for OTB is 256 MB
2021-05-17 15:30:13 (INFO) ReadImageInfo: GDAL maximum cache size is 3203 MB
2021-05-17 15:30:13 (INFO) ReadImageInfo: OTB will use at most 56 threads
2021-05-17 15:30:13 (INFO) ReadImageInfo: 
Image general information:
	Number of bands : 3
	Data type : float
	No data flags :
		Band 1: 0
		Band 2: 0
		Band 3: 0
	Start index :  [0,0]
	Size :  [4097,8193]
	Origin :  [780671,2.58396e+06]
	Spacing :  [2,-2]
	Estimated ground spacing (in meters): [1.99546,2.00688]

Image acquisition information:
	Sensor : 
	Image identification number: 
	Image projection : PROJCS["WGS 84 / UTM zone 39N",
    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",51],
    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","32639"]]

Image default RGB composition:
	[R, G, B] = [0,1,2]

Ground control points information:
	Number of GCPs = 0
	GCPs projection = 

Output parameters value:
indexx: 0
indexy: 0
sizex: 4097
sizey: 8193
spacingx: 2
spacingy: -2
originx: 780671
originy: 2583955
estimatedgroundspacingx: 1.995457888
estimatedgroundspacingy: 2.00687623
numberbands: 3
datatype: float
sensor: 
id: 
time: 
town: 
country: 
rgb.r: 0
rgb.g: 1
rgb.b: 2
projectionref: PROJCS["WGS 84 / UTM zone 39N",
    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",51],
    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","32639"]]
keyword: 
gcp.count: 0
gcp.proj: 
gcp.ids: 
gcp.info: 
gcp.imcoord: 
gcp.geocoord:
2

Dear @ashnair1,
Thank you for the feedback. If you need a bypass, you can try to use the ExtractROI application, it should generate the geom file.

By the way, we are currently working on the version 8.0 of OTB that won’t write the geom files anymore (but it will always be able to read them).

3

Dear @julienosman

My goal is to use Residual Registration to register two images. For the Residual Registration recipe, we require the geom file as per the documentation. This seems to be the case in the version 8.0 as well. Is the Residual Registration documentation out of date?

4

Indeed, the documentation haven’t been updated yet. I will make sure it is updated before the release. The new metadata framework allows to read the geometry metadata directly from the image without having to use a geom file.

5

That would be great. So how should residual registration be done according to the current setup? I’m still using OTB v7.2.0 and would appreciate any input.

6

Instead of using otbcli_ReadImageInfo (which appears to be broken), try otbcli_ExtractROI -in 2014_test.tif -out out.tif. It should generate the geom file.

7

I tried running the command and unfortunately it did not generate the geom file. Only the out.tif was created.

Logs

2021-05-18 16:57:55 (INFO) ExtractROI: Default RAM limit for OTB is 256 MB
2021-05-18 16:57:55 (INFO) ExtractROI: GDAL maximum cache size is 3203 MB
2021-05-18 16:57:55 (INFO) ExtractROI: OTB will use at most 56 threads
2021-05-18 16:57:55 (INFO): Estimated memory for full processing: 1152.31MB (avail.: 256 MB), optimal image partitioning: 5 blocks
2021-05-18 16:57:55 (INFO): File out.tif will be written in 6 blocks of 4097x1366 pixels
Writing out.tif...: 100% [**************************************************] (1s)
8

And with the extended filename ?

otbcli_ExtractROI -in 2014_test.tif?&writegeom=true

9

Adding ?&writegeom=true throws an error

otbcli_ExtractROI -in 2014_test.tif?&writegeom=true -out out.tif
[1] 2253

Command '-out' not found, did you mean:

  command 'gout' from deb scotch
  command 'lout' from deb lout

Try: sudo apt install <deb name>

ERROR: Missing mandatory parameter -out.


This is the ExtractROI application, version 7.2.0

Extract a ROI defined by the user.
Complete documentation: https://www.orfeo-toolbox.org/CookBook/Applications/app_ExtractROI.html or -help

Parameters: 
(base) ashwin@ashwin-laptop:~/Desktop/Projects/changedet/data$         -in                <string>         Input Image  (mandatory, default value is 2014_test.tif?)
MISSING -out               <string> [pixel] Output Image  [pixel=uint8/uint16/int16/uint32/int32/float/double/cint16/cint32/cfloat/cdouble] (default value is float) (mandatory)
        -mode              <string>         Extraction mode [standard/fit/extent/radius] (mandatory, default value is standard)
        -mode.fit.im       <string>         Reference image  (mandatory)
        -mode.fit.vect     <string>         Reference vector  (mandatory)
        -mode.extent.ulx   <float>          X coordinate of the Upper left corner  (mandatory, default value is 0)
        -mode.extent.uly   <float>          Y coordinate of Upper Left corner point  (mandatory, default value is 0)
        -mode.extent.lrx   <float>          X coordinate of Lower Right corner point  (mandatory, default value is 1000)
        -mode.extent.lry   <float>          Y coordinate of Lower Right corner point  (mandatory, default value is 1000)
        -mode.extent.unit  <string>         Unit [pxl/phy/lonlat] (mandatory, default value is pxl)
        -mode.radius.r     <float>          Radius  (mandatory, default value is 499)
        -mode.radius.unitr <string>         Radius unit [pxl/phy] (mandatory, default value is pxl)
        -mode.radius.cx    <float>          X coordinate of the center  (mandatory, default value is 499)
        -mode.radius.cy    <float>          Y coordinate of the center  (mandatory, default value is 499)
        -mode.radius.unitc <string>         Center unit [pxl/phy/lonlat] (mandatory, default value is pxl)
        -startx            <int32>          Start X  (mandatory, default value is 0)
        -starty            <int32>          Start Y  (mandatory, default value is 0)
        -sizex             <int32>          Size X  (mandatory, default value is 1000)
        -sizey             <int32>          Size Y  (mandatory, default value is 1000)
        -cl                <string list>    Output Image channels  (mandatory, default value is channel1)
        -elev              <group>          Elevation management 
        -elev.dem          <string>         DEM directory  (optional, off by default)
        -elev.geoid        <string>         Geoid File  (optional, off by default)
        -elev.default      <float>          Default elevation  (mandatory, default value is 0)
        -ram               <int32>          Available RAM (MB)  (optional, off by default, default value is 256)
        -progress          <boolean>        Report progress 
        -help              <string list>    Display long help (empty list), or help for given parameters keys

Use -help param1 [... paramN] to see detailed documentation of those parameters.

Examples: 
otbcli_ExtractROI -in VegetationIndex.hd -mode extent -mode.extent.ulx 40 -mode.extent.uly 40 -mode.extent.lrx 150 -mode.extent.lry 150 -out ExtractROI.tif
10

My bad. The extended filename goes to the -out parameter…

otbcli_ExtractROI -in 2014_test.tif -out out.tif?&writegeom=true

11

Still no luck I’m afraid. Running
otbcli_ExtractROI -in 2014_test.tif -out out.tif?&writegeom=true

results in out.tif file being written and then hanging for some reason (I’m assuming resolving the ?&writegeom=true bit)

2021-05-18 21:09:14 (INFO) ExtractROI: GDAL maximum cache size is 1602 MB
2021-05-18 21:09:14 (INFO) ExtractROI: OTB will use at most 12 threads
2021-05-18 21:09:14 (INFO): Estimated memory for full processing: 68.4357MB (avail.: 256 MB), optimal image partitioning: 1 blocks
2021-05-18 21:09:14 (INFO): File out.tif will be written in 1 blocks of 1000x1000 pixels
Writing out.tif?...: 100% [**************************************************] (0s)
^C
[1]+  Done                    /home/ashwin/Downloads/OTB-7.2.0-Linux64/bin/otbcli_ExtractROI -in 2014_test.tif -out out.tif?
12

It seems that Bash doesn’t like the & (it is interpreted as a “run it in background”, and what follows is ignored). To avoid this, write the parameter in quotes:

otbcli_ExtractROI -in 2014_test.tif -out "out.tif?&writegeom=true"

13

Now it’s gone back to just writing out.tif. Seems the writegeom is ignored.

$ otbcli_ExtractROI -in 2014_test.tif -out "out.tif?&writegeom=true"

2021-05-19 12:39:58 (INFO) ExtractROI: Default RAM limit for OTB is 256 MB
2021-05-19 12:39:58 (INFO) ExtractROI: GDAL maximum cache size is 3203 MB
2021-05-19 12:39:58 (INFO) ExtractROI: OTB will use at most 56 threads
2021-05-19 12:39:58 (INFO): Estimated memory for full processing: 1152.31MB (avail.: 256 MB), optimal image partitioning: 5 blocks
2021-05-19 12:39:58 (INFO): File out.tif will be written in 6 blocks of 4097x1366 pixels
Writing out.tif?&writegeom=true...: 100% [**************************************************] (1s)
14

Are you sure your image contains geometric metadata? What is the output of gdalinfo 2014_test.tif?

15

Output of gdalinfo 2014_test.tif

Driver: GTiff/GeoTIFF
Files: 2014_test.tif
       2014_test.tif.ovr
       2014_test.tif.ovr.aux.xml
       2014_test.tif.aux.xml
Size is 4097, 8193
Coordinate System is:
PROJCRS["WGS 84 / UTM zone 39N",
    BASEGEOGCRS["WGS 84",
        DATUM["World Geodetic System 1984",
            ELLIPSOID["WGS 84",6378137,298.257223563,
                LENGTHUNIT["metre",1]]],
        PRIMEM["Greenwich",0,
            ANGLEUNIT["degree",0.0174532925199433]],
        ID["EPSG",4326]],
    CONVERSION["UTM zone 39N",
        METHOD["Transverse Mercator",
            ID["EPSG",9807]],
        PARAMETER["Latitude of natural origin",0,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8801]],
        PARAMETER["Longitude of natural origin",51,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8802]],
        PARAMETER["Scale factor at natural origin",0.9996,
            SCALEUNIT["unity",1],
            ID["EPSG",8805]],
        PARAMETER["False easting",500000,
            LENGTHUNIT["metre",1],
            ID["EPSG",8806]],
        PARAMETER["False northing",0,
            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["unknown"],
        AREA["World - N hemisphere - 48°E to 54°E - by country"],
        BBOX[0,48,84,54]],
    ID["EPSG",32639]]
Data axis to CRS axis mapping: 1,2
Origin = (780670.000000000000000,2583956.000000000000000)
Pixel Size = (2.000000000000000,-2.000000000000000)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left  (  780670.000, 2583956.000) ( 53d44'41.86"E, 23d20'28.30"N)
Lower Left  (  780670.000, 2567570.000) ( 53d44'30.96"E, 23d11'36.06"N)
Upper Right (  788864.000, 2583956.000) ( 53d49'30.09"E, 23d20'23.17"N)
Lower Right (  788864.000, 2567570.000) ( 53d49'18.87"E, 23d11'30.96"N)
Center      (  784767.000, 2575763.000) ( 53d47' 0.43"E, 23d15'59.64"N)
Band 1 Block=4097x1 Type=Float32, ColorInterp=Gray
  NoData Value=0
  Overviews: 2049x4097
Band 2 Block=4097x1 Type=Float32, ColorInterp=Undefined
  NoData Value=0
  Overviews: 2049x4097
Band 3 Block=4097x1 Type=Float32, ColorInterp=Undefined
  NoData Value=0
  Overviews: 2049x4097
16

The image doesn’t seems to contain geometric metadata, so there is nothing to write in the geom file.
If the image doesn’t contain a Sensor Model, you won’t be able to run the RefineSensorModel application.

17

Hi @julienosman, in this context what do you mean by geometric metadata? I was assuming you were talking about the projection information but it seems that was incorrect. Could you give me an example? Maybe the output of gdalinfo for a tiff that has the metadata?

I ask because I have a lot of tiffs like this and want to know what is required to do Residual Registration on them?

18

Here is an example of metadata for a QUICKBIRD image (the one used in the recipe you are pointing at):

gdalinfo ~/QUICKBIRD/TOULOUSE/000000128955_01_P001_MUL/02APR01105228-M1BS-000000128955_01_P001.TIF

Driver: GTiff/GeoTIFF
Files: QUICKBIRD/TOULOUSE/000000128955_01_P001_MUL/02APR01105228-M1BS-000000128955_01_P001.TIF
       QUICKBIRD/TOULOUSE/000000128955_01_P001_MUL/02APR01105228-M1BS-000000128955_01_P001.TIF.ovr
       QUICKBIRD/TOULOUSE/000000128955_01_P001_MUL/02APR01105228-M1BS-000000128955_01_P001.TIF.ovr.aux.xml
       QUICKBIRD/TOULOUSE/000000128955_01_P001_MUL/02APR01105228-M1BS-000000128955_01_P001.IMD
       QUICKBIRD/TOULOUSE/000000128955_01_P001_MUL/02APR01105228-M1BS-000000128955_01_P001.RPB
Size is 6876, 7265
Coordinate System is `'
GCP Projection = 
GEOGCS["WGS 84",
    DATUM["WGS_1984",
        SPHEROID["WGS 84",6378137,298.257223563,
            AUTHORITY["EPSG","7030"]],
        AUTHORITY["EPSG","6326"]],
    PRIMEM["Greenwich",0],
    UNIT["degree",0.0174532925199433],
    AUTHORITY["EPSG","4326"]]
GCP[  0]: Id=1, Info=
          (0.5,0.5) -> (1.3544523884966,43.6633899110959,0)
GCP[  1]: Id=2, Info=
          (6875.5,0.5) -> (1.57254683562334,43.663858603131,0)
GCP[  2]: Id=3, Info=
          (6875.5,7264.5) -> (1.57224677372287,43.4879451348962,0)
GCP[  3]: Id=4, Info=
          (0.5,7264.5) -> (1.35617289802566,43.4876035537,0)
Metadata:
  AREA_OR_POINT=Area
  METADATATYPE=DG
  TIFFTAG_COPYRIGHT=(C) COPYRIGHT 2002 EarthWatch Incorporated, Longmont CO USA 80501-6700
  TIFFTAG_DATETIME=2004:05:25 19:46:40
  TIFFTAG_IMAGEDESCRIPTION={
  bandList = 
  [
    1;
    2;
    3;
    4;
  ]
}
  TIFFTAG_MAXSAMPLEVALUE=2047
  TIFFTAG_MINSAMPLEVALUE=1
Image Structure Metadata:
  INTERLEAVE=BAND
RPC Metadata:
  HEIGHT_OFF=241
  HEIGHT_SCALE=500
  LAT_OFF=43.5757
  LAT_SCALE=0.0881
  LINE_DEN_COEFF=+1.000000E+00 +3.073151E-05 +6.767593E-04 +2.759762E-04 +1.707011E-05 -4.215218E-07 -2.031982E-06 +6.281466E-06 -2.098686E-05 +1.091061E-05 +3.716818E-08 +5.435656E-08 -4.356190E-06 +1.090220E-08 +1.785190E-08 +4.342067E-04 +6.371344E-08 +0.000000E+00 -1.416592E-05 +1.653665E-08
  LINE_NUM_COEFF=-1.183432E-03 +1.540212E-03 -9.858373E-01 +1.594806E-02 +2.730892E-05 -9.785795E-07 +2.528461E-04 -7.762763E-04 +1.998298E-03 -5.855501E-06 -8.293341E-07 +4.496584E-08 +2.925877E-05 +0.000000E+00 -8.818227E-06 -4.465644E-05 -1.061420E-05 +4.172296E-07 -4.632173E-06 +1.915274E-07
  LINE_OFF=3630
  LINE_SCALE=3691
  LONG_OFF=1.4635
  LONG_SCALE=0.1091
  SAMP_DEN_COEFF=+1.000000E+00 +1.117681E-03 +2.334335E-03 -6.207135E-04 -5.962135E-05 +1.552717E-06 -1.840090E-06 -4.876876E-06 +3.158963E-05 -1.201072E-05 -1.579758E-07 -9.944747E-08 +1.942372E-06 +0.000000E+00 -8.945496E-08 -1.480163E-06 -6.008504E-08 +3.055564E-08 +1.110320E-07 +1.314046E-08
  SAMP_NUM_COEFF=-2.142200E-03 +1.003674E+00 +3.177080E-03 +2.667587E-04 -2.315384E-03 +6.203277E-04 -4.161892E-04 +1.035653E-03 -1.198854E-03 +7.437719E-06 +2.503507E-06 -2.185972E-05 -4.893912E-05 -1.149441E-05 +8.351017E-05 -7.642065E-05 -2.684169E-07 +2.340974E-08 +7.001246E-06 -3.846245E-08
  SAMP_OFF=3438
  SAMP_SCALE=3443
Corner Coordinates:
Upper Left  (    0.0,    0.0)
Lower Left  (    0.0, 7265.0)
Upper Right ( 6876.0,    0.0)
Lower Right ( 6876.0, 7265.0)
Center      ( 3438.0, 3632.5)
Band 1 Block=6876x76 Type=UInt16, ColorInterp=Red
  Overviews: 3438x3633, 1719x1817, 860x909, 430x455, 215x228, 108x114, 54x57, 27x29, 14x15, 7x8, 4x4, 2x2, 1x1
Band 2 Block=6876x76 Type=UInt16, ColorInterp=Green
  Overviews: 3438x3633, 1719x1817, 860x909, 430x455, 215x228, 108x114, 54x57, 27x29, 14x15, 7x8, 4x4, 2x2, 1x1
Band 3 Block=6876x76 Type=UInt16, ColorInterp=Blue
  Overviews: 3438x3633, 1719x1817, 860x909, 430x455, 215x228, 108x114, 54x57, 27x29, 14x15, 7x8, 4x4, 2x2, 1x1
Band 4 Block=6876x76 Type=UInt16, ColorInterp=Undefined
  Overviews: 3438x3633, 1719x1817, 860x909, 430x455, 215x228, 108x114, 54x57, 27x29, 14x15, 7x8, 4x4, 2x2, 1x1

This image contains RPC Metadata which are used to model the sensor.
Usually, once the geometric corrections are done, this model is not needed anymore, so it is removed from the image. If your images are already orthorectified, it is already registered, so you don’t really need to apply the Residual Registration on them.
In the recipe, the process is applied to a level1 Quickbird image. Most of the time, end users have access to level2 or level3 images, and don’t need to use this process.

1 Like
19

I see. That makes sense. Thanks for the detailed explanation.

Just to provide some context, my objective was to actually align two sets of images. One from 2014 and one from 2016 so that I can do Change Detection on them. However when I tried to do this via OTB, it said the images were not in the same physical space.

2021-05-20 13:28:54 (FATAL) MultivariateAlterationDetector: itk::ERROR: MultivariateAlterationDetectorImageFilter(0x7f0f9801fdb0): Inputs do not occupy the same physical space! 
InputImage Origin: [7.8067100e+05, 2.5839550e+06], InputImage_1 Origin: [7.7809500e+05, 2.5808890e+06]
	Tolerance: 2.0000000e-06

So I thought I had to co-register my images to make it work. We actually had a discussion previously here regarding a similar issue.

So if this is not the way to go, how can I prepare my two tiffs so that it will work with OTB’s change detection?

20

In the case of change detection, I would use SuperImpose to project one image into the other. Then, I could use FineRegistration to be sure the images are well registred (that’s probably what you wanted do do with RefineSensorModel).