Orthorectified WorldView rasters offset after 'grid based image resampling'

OTB Usage
1

Using WorldView panchromatics rasters

<PRODUCTLEVEL>Stereo OR2A</PRODUCTLEVEL>

assumed to be ‘ortho ready’ from above product level in the xml file, so I anticipated having to orthorectify before carrying out Stereoscopic reconstruction creation.

First I carry out Orthorectification of source rasters:

otbcli_OrthoRectification
-elev.dem F:\GIS\opentopo
-elev.geoid F:\OTB-6.6.1-Win64\geoid\egm96.grd
-interpolator “bco”
-interpolator.bco.radius 2
-io.in “F:/GIS/056078906030_01_P001_PAN/16JAN23133047-P2AS-056078906030_01_P001.TIF?&skipcarto=true”
-map “epsg”
-map.epsg.code 32723
-opt.gridspacing 4.0
-opt.rpc 10 -outputs.isotropic True
-outputs.mode “auto”
-io.out “F:/GIS/orthoR.tif”

which seems to work fine when briefly matching against google sat in QGIS -

image
image1064×778 412 KB

I’ve then carried out Stereo rectification grid generator

otbcli_StereoRectificationGridGenerator
-io.inleft “F:/GIS/orthoL.tif?&skipcarto=true”
-io.inright “F:/GIS/orthoR.tif?&skipcarto=true”
-epi.step 2
-epi.elevation.dem F:/GIS/opentopo
-epi.elevation.default 45.35 -epi.elevation.geoid F:\OTB-6.6.1-Win64\geoid\egm96.grd
-io.outleft F:/GIS/grid_stereo_rio_small_L.tif
-io.outright F:/GIS/grid_stereo_rio_small_R.tif

and then Grid based image resampling for left and right (only showing right below):

otbcli_GridBasedImageResampling
-io.in “F:/GIS/orthoR.tif?&skipcarto=true”
-io.out F:/GIS/epi_stereo_rio_small_R.tif
-grid.in F:/GIS/grid_stereo_rio_small_R.tif
-out.sizex 14002
-out.sizey 16053

But my issue is I’m getting an offset in the 2 rasters?

image
image664×735 194 KB

image
image932×766 238 KB

Are ‘ortho ready’ products sometimes already orthorectified?

Below is result of ReadImageInfo on source WorldView raster:

Input parameters:
{ ‘in’ : ‘F:/GIS/OTB_Help/Rio_de_Janeiro/056078906030_01_P001_PAN/16JAN23133047-P2AS-056078906030_01_P001.TIF’, ‘keywordlist’ : False, ‘outkwl’ : ‘TEMPORARY_OUTPUT’ }

2020-01-23 01:10:18 (INFO): Default RAM limit for OTB is 128 MB
2020-01-23 01:10:18 (INFO): GDAL maximum cache size is 812 MB
2020-01-23 01:10:18 (INFO): OTB will use at most 8 threads
2020-01-23 01:10:18 (INFO): Loading kwl metadata from official product in file F:/GIS/OTB_Help/Rio_de_Janeiro/056078906030_01_P001_PAN/16JAN23133047-P2AS-056078906030_01_P001.TIF
2020-01-23 01:10:18 (INFO):
Image general information:
Number of bands : 1
No data flags : Not found
Start index : [0,0]
Size : [13632,11244]
Origin : [684717,7.46199e+06]
Spacing : [0.3,-0.3]
Estimated ground spacing (in meters): [0.299544,0.301234]

Image acquisition information:
Sensor : WV03
Image identification number: 16JAN23133047-P2AS-056078906030_01_P001
Image projection : PROJCS[“WGS 84 / UTM zone 23S”,
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”]],
PROJECTION[“Transverse_Mercator”],
PARAMETER[“latitude_of_origin”,0],
PARAMETER[“central_meridian”,-45],
PARAMETER[“scale_factor”,0.9996],
PARAMETER[“false_easting”,500000],
PARAMETER[“false_northing”,10000000],
UNIT[“metre”,1,
AUTHORITY[“EPSG”,“9001”]],
AUTHORITY[“EPSG”,“32723”]]

Image footprint coordinates:
Upper left corner (latitude, longitude) = [-22.94,-43.1986]
Upper right corner (latitude, longitude) = [-22.9395,-43.1587]
Lower left corner (latitude, longitude) = [-22.9704,-43.1982]
Lower right corner (latitude, longitude) = [-22.97,-43.1583]

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: 13632
sizey: 11244
spacingx: 0.3000000119
spacingy: -0.3000000119
originx: 684716.5625
originy: 7461990
estimatedgroundspacingx: 0.299544096
estimatedgroundspacingy: 0.3012338281
numberbands: 1
sensor: WV03
id: 16JAN23133047-P2AS-056078906030_01_P001
time:
ullat: -22.9399929
ullon: -43.19857788
urlat: -22.93953514
urlon: -43.15871429
lrlat: -22.96998787
lrlon: -43.15830231
lllat: -22.97044563
lllon: -43.19817734
town: Not available
country: Not available
rgb.r: 0
rgb.g: 1
rgb.b: 2
projectionref: PROJCS[“WGS 84 / UTM zone 23S”,
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”]],
PROJECTION[“Transverse_Mercator”],
PARAMETER[“latitude_of_origin”,0],
PARAMETER[“central_meridian”,-45],
PARAMETER[“scale_factor”,0.9996],
PARAMETER[“false_easting”,500000],
PARAMETER[“false_northing”,10000000],
UNIT[“metre”,1,
AUTHORITY[“EPSG”,“9001”]],
AUTHORITY[“EPSG”,“32723”]]
keyword:
gcp.count: 0
gcp.proj:
gcp.ids:
gcp.info:
gcp.imcoord:
gcp.geocoord:

And then ReadImageInfo on raster generated after the Grid based image resampling step:

Input parameters:
{ ‘in’ : ‘F:/GIS/OTB_Help/Rio_de_Janeiro/22_01_20b/epi_stereo_rio_small_R.tif’, ‘keywordlist’ : False, ‘outkwl’ : ‘TEMPORARY_OUTPUT’ }

2020-01-23 01:13:59 (INFO): Default RAM limit for OTB is 128 MB
2020-01-23 01:13:59 (INFO): GDAL maximum cache size is 812 MB
2020-01-23 01:13:59 (INFO): OTB will use at most 8 threads
2020-01-23 01:13:59 (INFO): Loading kwl metadata from attached geom file F:/GIS/OTB_Help/Rio_de_Janeiro/22_01_20b/epi_stereo_rio_small_R.geom
2020-01-23 01:13:59 (INFO):
Image general information:
Number of bands : 1
No data flags :
Band 1: 0
Start index : [0,0]
Size : [14002,16053]
Origin : [0,0]
Spacing : [1,1]
Estimated ground spacing (in meters): [0.299487,0.301214]

Image acquisition information:
Sensor : WV03
Image identification number: 16JAN23133047-P2AS-056078906030_01_P001

Image footprint coordinates:
Upper left corner (latitude, longitude) = [-22.94,-43.1986]
Upper right corner (latitude, longitude) = [-22.9395,-43.1587]
Lower left corner (latitude, longitude) = [-22.9704,-43.1982]
Lower right corner (latitude, longitude) = [-22.97,-43.1583]

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: 14002
sizey: 16053
spacingx: 1
spacingy: 1
originx: 0
originy: 0
estimatedgroundspacingx: 0.2994868755
estimatedgroundspacingy: 0.3012136221
numberbands: 1
sensor: WV03
id: 16JAN23133047-P2AS-056078906030_01_P001
time:
ullat: -22.9399929
ullon: -43.19857788
urlat: -22.93953514
urlon: -43.15871429
lrlat: -22.96998787
lrlon: -43.15830231
lllat: -22.97044563
lllon: -43.19817734
town: Not available
country: Not available
rgb.r: 0
rgb.g: 1
rgb.b: 2
projectionref:
keyword:
gcp.count: 0
gcp.proj:
gcp.ids:
gcp.info:
gcp.imcoord:
gcp.geocoord:

2

Do you discovered if the OrthoReady image was already orthorectified or what was the problem?

3

Hi Andrea,

Looking in the associated XML file 056078906030_01_P001_PAN\16JAN23133047-P2AS-056078906030_01_P001.XML

I noted the type of product:

<PRODUCTLEVEL>Stereo OR2A</PRODUCTLEVEL>

And reading here - https://www.c-agg.org/wp-content/uploads/DigitalGlobe-Base-Product-FAQ.pdf

I see they mention the following re OR2A products:

Ortho Ready Standard Imagery can be orthorectified using COTS packages like ERDAS IMAGINE, PCI Geomatics, ENVI, and SOCET SET/GXP

So I have to assume it’s not been orthorectified.

However when reading notes in orfeo cookbook documentation:
https://www.orfeo-toolbox.org/CookBook/recipes/optpreproc.html?highlight=ortho%20ready

They warn of ‘ortho ready’ products:

There are some image products, called “ortho-ready”, that should be processed carefully. They are actual products in raw geometry, but their metadata also contains projection data:

Having checked the rasters with OTB’s ReadImageInfo they do indeed hold projection information, so I tried to get around this by using the recommended ?&skipcarto=true suffix within the commands.

But I still get the offset issue just after Grid based image resampling step, as shown above :frowning:

Hopefully someone has seen something similar and has some idea how to get around this?

I have tested the whole process without attempting orthorectification and I don’t get the offset issue, however the results are not good when attempting the full Stereoscopic reconstruction process.

Results are very noisy -

GetImage (2)
GetImage (2)997×761 370 KB

Hence why I’m attempting to first orthorectify source rasters, with the hope that it will improve resulting DEM…

4

Hi all,

I confirm that:

  • StereoRectificationGridGenerator is not meant to be used with ortho-rectified images,
  • skipcarto option should be used whenever dealing with ortho-ready products (I do not know who invented that ortho-ready thing, but it is definitely confusing for users and developers).

If the results are noisy with the full pipeline, you may want to tweak downstream parameters such as output.fusionmethod (try the mean option). postproc.med (activate it), postproc.metrict (increase it), mask.variancet (activate it).

That said, remember that OTB is not a tool dedicated to stereo-restitution, and as such it will not compete with the quality you will get from other open-source products such as MicMac, S2P or AMES. In particular, it is using traditional block matching when all the other makes use of the SGM algorithm.

1 Like
5

Thank you,

Very helpful tips.

I have not re-visited attempts using the single application ‘StereoFramework’ as I had trouble with it so have since resorted to breaking down the task using the individual processes.

With that in mind, I did notice that -mask.variancet float is available within otbcli_BlockMatching.

However re the other 3 you mentioned:
-output.fusionmethod [max|min|mean|acc] Default value: max
-postproc.med bool Default value: false
-postproc.metrict float Default value: 0.6

I can’t seem to find those parameters within the individual processes?

Are they only available within the ‘StereoFramework’?

Regards,

Ross