Problem reading metadata with WV-1, WV-2 and Formosat-2

Matthijs-Oosterhuis · April 13, 2021, 2:09pm

Dear all,
Im trying to optically calibrate a WV-1, WV-2 and a Formosat-2 image, but it seems there is an issue reading the metadata. My directory and metadata is all present in the folder of these 3 missions, as provided by the mission provider. The .IMD, .TIL and .XML data is all present. An example of how I call it is this:

otbcli_OpticalCalibration -in /home/working_folder/012403626010_01/012403626010_01_P001_PSH/20MAY02102340-S3DS_R2C1-012403626010_01_P001.TIF -level toa -out /home/processed_data/OptiCali_EW01_PAN_20171225 float -ram 16000

With all 3 of these missions I get this error:

    Sensor unknown!
[2021-04-13 15:48:47,712] {bash.py:173} INFO - Additional parameters are necessary, please provide them (cf. documentation)!
[2021-04-13 15:48:47,712] {bash.py:173} INFO - 
[2021-04-13 15:48:47,712] {bash.py:173} INFO - 2021-04-13 15:48:47 (DEBUG): file ../Modules/Wrappers/ApplicationEngine/src/otbWrapperApplication.cxx, line 1545, OpticalCalibration (0x1678cb0): MISSING : acqui.gainbias (Level 2)
[2021-04-13 15:48:47,712] {bash.py:173} INFO - ERROR: Missing mandatory parameter -acqui.gainbias.
[2021-04-13 15:48:47,712] {bash.py:173} INFO -

Any thoughts or suggestions of what I’m doing wrong?

Thanks in advance,
Thijs

Cedric · April 14, 2021, 8:33am

Hello,

I don’t think WV-1 products are supported in the OpticalCalibration application. I am not sure though, maybe OTB is able to read WV1 data as WV2 metadatas. But I don’t have a WV1 product on hand to test it.

Regarding WV2 and Formasat, are you able to read the products with gdalinfo ? In particular it should return METADATATYPE=DIMAP for Formosat and METADATATYPE=DG for WV2.

Can you run the ReadImageInformation application on these products

otbcli_ReadImageInformation -in /home/working_folder/012403626010_01/012403626010_01_P001_PSH/20MAY02102340-S3DS_R2C1-012403626010_01_P001.TIF -outkwl out.geom

the produced geom file will contain the parsed metadata, it can be helpful to identify what is going wrong here.

Cédric

Matthijs-Oosterhuis · April 20, 2021, 2:43pm

Dear Cedric,
Strangely enough METADATATYPE is not given in the result of gdalinfo.

This is my result of gdalinfo:

Driver: GTiff/GeoTIFF
Files: 20MAY02102340-S3DS_R1C1-012403626010_01_P001.TIF
Size is 14808, 16384
Coordinate System is:
PROJCRS["WGS 84 / UTM zone 32N",
    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["Transverse Mercator",
        METHOD["Transverse Mercator",
            ID["EPSG",9807]],
        PARAMETER["Latitude of natural origin",0,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8801]],
        PARAMETER["Longitude of natural origin",9,
            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["easting",east,
            ORDER[1],
            LENGTHUNIT["metre",1]],
        AXIS["northing",north,
            ORDER[2],
            LENGTHUNIT["metre",1]],
    ID["EPSG",32632]]
Data axis to CRS axis mapping: 1,2
Origin = (548162.500000000000000,5064503.500000000000000)
Pixel Size = (0.500000000000000,-0.500000000000000)
Metadata:
  AREA_OR_POINT=Area
  TIFFTAG_COPYRIGHT=(C) COPYRIGHT 2019 DigitalGlobe, Inc., Longmont CO USA 80503
  TIFFTAG_DATETIME=2020:05:02 11:33:06
  TIFFTAG_IMAGEDESCRIPTION={
  bandList =
  [
    4;
    3;
    2;
  ]
}
  TIFFTAG_MAXSAMPLEVALUE=2047
  TIFFTAG_MINSAMPLEVALUE=0
Image Structure Metadata:
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  (  548162.500, 5064503.500) (  9d37' 8.46"E, 45d43'56.68"N)
Lower Left  (  548162.500, 5056311.500) (  9d37' 5.54"E, 45d39'31.25"N)
Upper Right (  555566.500, 5064503.500) (  9d42'51.02"E, 45d43'54.68"N)
Lower Right (  555566.500, 5056311.500) (  9d42'47.64"E, 45d39'29.26"N)
Center      (  551864.500, 5060407.500) (  9d39'58.17"E, 45d41'43.00"N)
Band 1 Block=14808x35 Type=UInt16, ColorInterp=Red
Band 2 Block=14808x35 Type=UInt16, ColorInterp=Green
Band 3 Block=14808x35 Type=UInt16, ColorInterp=Blue

After performing the otbcli_ReadImageInfo command i get this result:
2021-04-20 16:26:04 (INFO) ReadImageInfo: Default RAM limit for OTB is 256 MB
2021-04-20 16:26:04 (INFO) ReadImageInfo: GDAL maximum cache size is 782 MB
2021-04-20 16:26:04 (INFO) ReadImageInfo: OTB will use at most 8 threads
2021-04-20 16:26:04 (INFO) ReadImageInfo:
Image general information:
Number of bands : 3
Data type : unsigned_short
No data flags : Not found
Start index : [0,0]
Size : [6946,628]
Origin : [548163,5.05631e+06]
Spacing : [0.5,-0.5]
Estimated ground spacing (in meters): [0.498746,0.500754]

Image acquisition information:
        Sensor :
        Image identification number:
        Image projection : PROJCS["WGS 84 / UTM zone 32N",
    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","9122"]],
        AUTHORITY["EPSG","4326"]],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",0],
    PARAMETER["central_meridian",9],
    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","32632"]]

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: 6946
sizey: 628
spacingx: 0.5
spacingy: -0.5
originx: 548162.75
originy: 5056311
estimatedgroundspacingx: 0.4987464845
estimatedgroundspacingy: 0.5007536411
numberbands: 3
datatype: unsigned_short
sensor:
id:
time:
town:
country:
rgb.r: 0
rgb.g: 1
rgb.b: 2
projectionref: PROJCS["WGS 84 / UTM zone 32N",
    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","9122"]],
        AUTHORITY["EPSG","4326"]],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",0],
    PARAMETER["central_meridian",9],
    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","32632"]]
keyword:
gcp.count: 0
gcp.proj:
gcp.ids:
gcp.info:
gcp.imcoord:
gcp.geocoord:

The folder structure of my Worldview-2 image is this, might be of help:

20MAY02102340-S3DS-012403626010_01_P001-BROWSE.JPG
2020 20MAY02102340-S3DS-012403626010_01_P001.IMD
2020 20MAY02102340-S3DS-012403626010_01_P001_README.TXT
20MAY02102340-S3DS-012403626010_01_P001.TIL
20MAY02102340-S3DS-012403626010_01_P001.XML
20MAY02102340-S3DS_R1C1-012403626010_01_P001.TIF
20MAY02102340-S3DS_R2C1-012403626010_01_P001.TIF
BASE.TXT

Thanks for your thoughts!
Thijs

Cedric · April 21, 2021, 7:45am

Hello,

I don’t why, but it seems that the RiCi tile naming convention is not supported by OTB nor by GDAL. For an input product “product.TIF”, both readers are looking for:

product.IMD
product.RPB (this one appears to be missing in your case, this is the file containing RPC coefficients)
product XML (I’m not sure this one is actually needed)

you can try to rename the files to follow this convention (i think the "2020 " prefix before the .IMD file should also be removed), but I don’t think you will be able to open the dataset as a WV2 product without the .RPB file.

Cédric

Matthijs-Oosterhuis · May 4, 2021, 2:33pm

Thanks for the help Cedric. I found that when using a different product type (OrthoReady), it did include .RPB and ran succesfully. Also surprisingly enough for Worldview-1.

I wonder if you have any advise how to orthorectify the product with .TIL to a product with .RPB included.

Thanks for the support!
Thijs