Skip to content
Permalink
Browse files

Fix transformCoordinates() for rasters with RPC georeferencing (#38350)

Unfortunately for raster with RPC georeferencing, the transform
used in GDAL provider's transformCoordinates() was using RPC_HEIGHT
at zero elevation, while rendering pipeline was using RPC_HEIGHT
at the correct elevation from HEIGHT_OFF metadata. This was causing
shift in what was rendered in map canvas vs what transform function
returned. This is now fixed by also the transform function using
the correct elevation.
  • Loading branch information
wonder-sk committed Aug 21, 2020
1 parent b98cfe1 commit 418bcf3d57e50f65a626f79e0a2cd67f4a340704
Showing with 30 additions and 1 deletion.
  1. +5 −1 src/core/providers/gdal/qgsgdalprovider.cpp
  2. +15 −0 src/core/qgsgdalutils.cpp
  3. +10 −0 src/core/qgsgdalutils.h
@@ -2920,6 +2920,7 @@ void QgsGdalProvider::initBaseDataset()
mGdalDataset =
QgsGdalUtils::rpcAwareAutoCreateWarpedVrt( mGdalBaseDataset, nullptr, nullptr,
GRA_NearestNeighbour, 0.2, nullptr );
mGdalTransformerArg = QgsGdalUtils::rpcAwareCreateTransformer( mGdalBaseDataset );
}
else
{
@@ -2943,7 +2944,10 @@ void QgsGdalProvider::initBaseDataset()
mGdalDataset = mGdalBaseDataset;
}

mGdalTransformerArg = GDALCreateGenImgProjTransformer( mGdalBaseDataset, nullptr, nullptr, nullptr, TRUE, 1.0, 0 );
if ( !mGdalTransformerArg )
{
mGdalTransformerArg = GDALCreateGenImgProjTransformer( mGdalBaseDataset, nullptr, nullptr, nullptr, TRUE, 1.0, 0 );
}

if ( !hasGeoTransform )
{
@@ -490,6 +490,21 @@ GDALDatasetH QgsGdalUtils::rpcAwareAutoCreateWarpedVrt(
return GDALAutoCreateWarpedVRTEx( hSrcDS, pszSrcWKT, pszDstWKT, eResampleAlg, dfMaxError, psOptionsIn, opts );
}

void *QgsGdalUtils::rpcAwareCreateTransformer( GDALDatasetH hSrcDS, GDALDatasetH hDstDS, char **papszOptions )
{
char **opts = CSLDuplicate( papszOptions );
if ( GDALGetMetadata( hSrcDS, "RPC" ) )
{
// well-behaved RPC should have height offset a good value for RPC_HEIGHT
const char *heightOffStr = GDALGetMetadataItem( hSrcDS, "HEIGHT_OFF", "RPC" );
if ( heightOffStr )
opts = CSLAddNameValue( opts, "RPC_HEIGHT", heightOffStr );
}
void *transformer = GDALCreateGenImgProjTransformer2( hSrcDS, hDstDS, opts );
CSLDestroy( opts );
return transformer;
}

#ifndef QT_NO_NETWORKPROXY
void QgsGdalUtils::setupProxy()
{
@@ -116,6 +116,16 @@ class CORE_EXPORT QgsGdalUtils
double dfMaxError,
const GDALWarpOptions *psOptionsIn );

/**
* This is a wrapper around GDALCreateGenImgProjTransformer2() that takes into account RPC
* georeferencing (it sets RPC_HEIGHT in GDALCreateGenImgProjTransformer2 based on HEIGHT_OFF).
* By default GDAL would assume that the imagery has zero elevation - if that is not the case,
* the image would not be shown in the correct location.
*
* \since QGIS 3.16
*/
static void *rpcAwareCreateTransformer( GDALDatasetH hSrcDS, GDALDatasetH hDstDS = nullptr, char **papszOptions = nullptr );

#ifndef QT_NO_NETWORKPROXY
//! Sets the gdal proxy variables
static void setupProxy();

0 comments on commit 418bcf3

Please sign in to comment.
You can’t perform that action at this time.