[gdal-dev] gdalwarp with RPCs causing bogus output

Even Rouault even.rouault at spatialys.com
Wed Apr 17 15:20:10 PDT 2024


Hi,


I haven't done myself that exercice but I know that computing RPC values 
that are stable enough might be challenging, so perhaps the issue is 
related to that


A few other remarks:

- your gdalwarp command line refers to  RPC_DEM_SRS and 
RPC_DEM_MISSING_VALUE but doesn't include a RPC_DEM itself, hence those 
are likely to be non effective

- you generally don't want to use both RPC_HEIGHT and RPC_DEM. 
RPC_HEIGHT is essentially useful when you don't have a DEM available, 
and thus fallback taking an average elevation

- there's a subtlety regarding DEM. RPC reference for altitudes is WGS84 
ellipsoidal height, not orthometric/MSL altitude. But DEM values use 
orthometric/MSL altitude, hence a geoid correction must be applied. So 
you'd rather want to use RPC_DEM_SRS=EPSG:4326+5773 for example to use 
the EGM96 geoid (cf https://github.com/OSGeo/gdal/issues/3298). That 
said, misusing orthometric altitude vs ellipsoidal height generally 
accounts for small shifts, not big ones

- you could use gdaltransform with all your -to parameter and 
input_with_RPCs.tiff to check at least that the forward RPC 
transformation path works correctly (the one from longitude, latitude -> 
column, line). And with -i to check the inverse RPC transformer. The 
inverse RPC transformer can have a hard time converging in montainous 
areas or for images off-nadir


Even


Le 17/04/2024 à 23:30, Joseph McGlinchy via gdal-dev a écrit :
> Hello,
>
> I am attempting to implement georegistration through RPC. I have the 
> following information I've used to calibrate the RPC coefficients, 
> using all terms for numerator and denominator for both sample and line 
> equations.
>
>  *
>     image grid, stored in tiff format with no geo-information
>     associated with it, so it reflects the imaging orientation
>  *
>     a 'ground' grid, which corresponds to the longitude/latitude
>     coordinates for each pixel determined from a line-of-sight vector
>     and imaging system coordinates where the pseudo-rays intersect the
>     WGS84 geoid
>  *
>     a random sample of up to 4000 image coordinate / object space
>     coordinate pairs (I arrived at this number through trial and
>     error; using all pixels explodes RAM)
>  *
>     elevation extracted at the longitude/latitude coordinates from a
>     DEM, in this case, SRTM, but have also tried using NASADEM
>
>
> I am able to write out an image to EPSG:4326 by populating the RPC 
> metadata and using the non-georeferenced image data. However, I am 
> struggling to use |gdalwarp|​ in a reliable way to orthorectify that 
> data, let alone writing it out in a way that 'bakes in' the 
> georegistration with the RPCs whether that is in EPSG:4326 or the 
> local UTM zone. I see strips of the image "ripped out" with odd curves 
> in various places throughout the image.
>
> The only way I've been able to use |gdalwarp|​ to write the image at 
> all is with the following parameters (any DEM reference is to the SRTM 
> DEM):
>
> gdalwarp --config CPL_DEBUG ON -t_srs EPSG:32610 -rpc -to 
> "RPC_DEM_SRS=+proj=longlat +datum=WGS84 +no_def" -to "RPC_HEIGHT=350" 
> -to "RPC_DEM_MISSING_VALUE=0.001" -to "RPC_FOOTPRINT='POLYGON ((list 
> of polygon coordinates comprising the long/lat grid))'" -to 
> "RPC_MAX_ITERATIONS=101" input_with_RPCs.tiff output.tiff
>
> This is the only configuration i can use to run |gdalwarp|​ 
> successfully.  removing any single RPC_X tranformer option gives me 
> bogus output. The RPC_HEIGHT value i specify above is not close to the 
> mean or median elevation of the extent of my data; mean is ~195m and 
> median is ~150m.
>
> With the debug turned on, any other set of parameters gives me failed 
> RPC convergence on several points. I am able to reproduce this 
> regularly by specifying RPC_DEM=dem.tif, where dem.tif is the same 
> data I used to extract elevation values when calibrating the RPCs. I 
> am seeing normalized latitude and longitude values with magnitude > 1 
> (I checked every location in the image, based on the metadata, the 
> range is not outside of [-1,1]), as well as normalized altitude values 
> with magnitude > 1 (there are some, not many, that have magnitude of 
> 1.75).
>
> My workflow can be summarized as:
>
> 1.
>     load grids (image data, longitude, latitude)
> 2.
>     randomly sample up to 4000 points in image coordinates, object
>     coordinates
>     1.
>         assign z-value from SRTM DEM
>     2.
>         evaluate if any of the points are in NODATA areas of SRTM
>         (image is coastal, so there are NODATA areas for SRTM here),
>         if so, remove those and generate more points
> 3.
>     normalize coordinates of grids to be in [-1,1], recording offsets
>     and scale
> 4.
>      calibrate RPC coefficients using all terms
> 5.
>     write out GeoTIFF with image grid for pixels, along with RPC
>     required metadata fields and CRS EPSG:4326
>
>
> System information (please let me know if more is needed)
> OS: Ubuntu 20.04 LTS (GNU/Linux 5.10.16.3-microsoft-standard-WSL2 
> x86_64) (Windows Subsystem for Linux)
> GDAL 3.6.0 (python)
>
> Thank you in advance for any insight into this process! I am happy to 
> package up any of the data I am using, as well. I have placed the 
> initial data from the end of Step 5 described above, along with some 
> additional files, a sample gdalwarp call, and a file-list.txt, at 
> https://drive.google.com/drive/folders/1BfevhKQa4ZHi_OQfiX_rUk2sqeoNVhyM?usp=sharing 
> <https://drive.google.com/drive/folders/1BfevhKQa4ZHi_OQfiX_rUk2sqeoNVhyM?usp=sharing>
>
> If more is needed for anyone interested in having a look, please let 
> me know and I'll upload.
>
>
> Cheers,
> Joe
>
>
> _______________________________________________
> 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.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/gdal-dev/attachments/20240418/70863c28/attachment.htm>


More information about the gdal-dev mailing list