TransWikia.com

Warp vector WGS84 polygons to an orthorectified sub image

Geographic Information Systems Asked by Erotemic on July 14, 2021

I’m using GDAL Python bindings (and gdal command line tools via subprocessing) and I’m having trouble warping WGS84 polygons to an orthorectified subimage. I’m fairly new to GIS, so please excuse any incorrect terminology. I’ll try to explain in as much detail as possible.

First, here is what I have:

  • A large geotiff with RPC metadata.
  • A WGS84 region of interest.
  • Multiple WGS84 polygons inside that region of interest.

What I have been able to do

  • Use gdalwarp to extract an the orthorectified ROI from the large geotiff. I did this via gdalwap:
gdalwarp 
-te {xmin} {ymin} {xmax} {ymax} 
-te_srs epsg:4326 
-s_srs epsg:4326 
-rpc -et 0 
-to RPC_DEM={dem_fpath} 
-overwrite 
{orig_geotiff} {roi_geotiff}

Where xmin/ymin/xmax/ymax are the lat/lon bounds of the ROI, dem_fpath is the path to the corresponding DEM tif, orig_geotiff is the large source geotiff and roi_geotiff is where the new much smaller orthorectified image is written.

  • I can warp the WGS84 polygons into pixelspace of the original orig_geotiff large geotiff.
import gdal
src_ref = gdal.Open(orig_geotiff)
tr = gdal.Transformer(src_ref, None, ['METHOD=RPC', 'RPC_PIXEL_ERROR_THRESHOLD=0.00' f'RPC_DEM={dem_fpath}'])
orig_pixel_polys = []
for poly in geo_poly_list:
    pts_in = poly.data['exterior'].data
    pts_out = []
    lats, lons = pts_in.T
    # I can use the DEM to find a better
    # initialization for this but I'm keeping it
    # simple in this example
    elevs = [0] * len(pts_in)
    for x_in, y_in, z_in in zip(lons, lats, elevs):
        bDstToSrc = 1
        success, pnt = tr.TransformPoint(bDstToSrc, x_in, y_in, z_in)
        x_out, y_out, _ = pnt
        pts_out.append((x_out, y_out))
    pts_out = np.array(pts_out)
    orig_pixel_polys.append(pts_out)

What I need help with

But what I’m having an incredibly hard time doing is trying to either (1) warp those pixel-polygon points in the original image to the orthorectified image or (2) warp the WGS84-polygon points directly into the orthorectified image.

In the first case, I’m not sure how to find whatever the transform gdalwarp is using and apply it to a vector dataset.

In the second case, the file written by gdalwarp does not recompute new RPCs for the orthorectified subimage; it only brings over the affine geotransform. Unfortunately, this affine transform does not do a good job of mapping the world coordinate system to the pixel coordinates. If I look at the alignment between the data and the pixel coordinates in the original image it looks fine (because it used the RPCs to do map WGS84 to original pixel coordinates) but when I try warping using the affine geotransform in the orthorectified image (or in the original image for that matter) they do not align well.

Summary

In summary I’m looking for an answer for one of the following:

  • How do I take pixel coordinates in an original image and map them to the corresponding pixel coordinates in a a cropped subimage that was orthorectified with RPC metadata?

  • How do I take WGS84 coordinates and map them into a pixel coordinates in a cropped subimage that was orthorectified with RPC metadata?

A method for either of these would answer my one question of: "How do I warp vector WGS84 polygons to an orthorectified sub image".

Add your own answers!

Ask a Question

Get help from others!

© 2024 TransWikia.com. All rights reserved. Sites we Love: PCI Database, UKBizDB, Menu Kuliner, Sharing RPP