[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