Skip to content

Commit

Permalink
Merge pull request #10602 from rouault/fix_10600
Browse files Browse the repository at this point in the history
gdal_translate/GDALOverviewDataset: fix half-pixel shift issue when rescaling RPC
  • Loading branch information
rouault authored Aug 25, 2024
2 parents db4b439 + 64ab5eb commit 0266bc0
Show file tree
Hide file tree
Showing 5 changed files with 74 additions and 33 deletions.
22 changes: 18 additions & 4 deletions apps/gdal_translate_lib.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -1826,10 +1826,24 @@ GDALDatasetH GDALTranslate(const char *pszDest, GDALDatasetH hSrcDataset,

const double df2 = adfSrcWinOri[2];
const double df3 = adfSrcWinOri[3];
dfSAMP_OFF *= nOXSize / df2;
dfLINE_OFF *= nOYSize / df3;
dfSAMP_SCALE *= nOXSize / df2;
dfLINE_SCALE *= nOYSize / df3;
const double dfXRatio = nOXSize / df2;
const double dfYRatio = nOYSize / df3;

// For line offset and pixel offset, we need to convert from RPC
// pixel center registration convention to GDAL pixel top-left corner
// registration convention by adding an initial 0.5 shift, and un-apply
// it after scaling.

dfSAMP_OFF += 0.5;
dfSAMP_OFF *= dfXRatio;
dfSAMP_OFF -= 0.5;

dfLINE_OFF += 0.5;
dfLINE_OFF *= dfYRatio;
dfLINE_OFF -= 0.5;

dfSAMP_SCALE *= dfXRatio;
dfLINE_SCALE *= dfYRatio;

CPLString osField;
osField.Printf("%.15g", dfLINE_OFF);
Expand Down
13 changes: 6 additions & 7 deletions autotest/gcore/overviewds.py
Original file line number Diff line number Diff line change
Expand Up @@ -225,13 +225,12 @@ def test_overviewds_4(tmp_path):

for key in rpc_md:
assert ds.GetMetadataItem(key, "RPC") == got_md[key]
if (
key == "LINE_SCALE"
or key == "SAMP_SCALE"
or key == "LINE_OFF"
or key == "SAMP_OFF"
):
assert float(got_md[key]) == myfloat(rpc_md[key]) / 2
if key == "LINE_SCALE" or key == "SAMP_SCALE":
assert float(got_md[key]) == pytest.approx(myfloat(rpc_md[key]) / 2)
elif key == "LINE_OFF" or key == "SAMP_OFF":
assert float(got_md[key]) == pytest.approx(
(myfloat(rpc_md[key]) + 0.5) / 2 - 0.5
)
elif got_md[key] != rpc_md[key]:
print(got_md[key])
print(rpc_md[key])
Expand Down
28 changes: 22 additions & 6 deletions autotest/utilities/test_gdal_translate.py
Original file line number Diff line number Diff line change
Expand Up @@ -872,16 +872,32 @@ def test_gdal_translate_32(gdal_translate_path, tmp_path):

dst_tif = str(tmp_path / "out.tif")

src_ds = gdal.Open("../gcore/data/byte_rpc.tif")
src_md = src_ds.GetMetadata("RPC")
srcxoff = 1
srcyoff = 2
srcwidth = 13
srcheight = 14
widthratio = 200
heightratio = 300
gdaltest.runexternal(
f"{gdal_translate_path} ../gcore/data/byte_rpc.tif {dst_tif} -srcwin 1 2 13 14 -outsize 150% 300%"
f"{gdal_translate_path} ../gcore/data/byte_rpc.tif {dst_tif} -srcwin {srcxoff} {srcyoff} {srcwidth} {srcheight} -outsize {widthratio}% {heightratio}%"
)
widthratio /= 100.0
heightratio /= 100.0
ds = gdal.Open(dst_tif)
md = ds.GetMetadata("RPC")
assert (
float(md["LINE_OFF"]) == pytest.approx(47496, abs=1e-5)
and float(md["LINE_SCALE"]) == pytest.approx(47502, abs=1e-5)
and float(md["SAMP_OFF"]) == pytest.approx(19676.6923076923, abs=1e-5)
and float(md["SAMP_SCALE"]) == pytest.approx(19678.1538461538, abs=1e-5)
assert float(md["LINE_OFF"]) == pytest.approx(
(float(src_md["LINE_OFF"]) - srcyoff + 0.5) * heightratio - 0.5, abs=1e-5
)
assert float(md["LINE_SCALE"]) == pytest.approx(
float(src_md["LINE_SCALE"]) * heightratio, abs=1e-5
)
assert float(md["SAMP_OFF"]) == pytest.approx(
(float(src_md["SAMP_OFF"]) - srcxoff + 0.5) * widthratio - 0.5, abs=1e-5
)
assert float(md["SAMP_SCALE"]) == pytest.approx(
float(src_md["SAMP_SCALE"]) * widthratio, abs=1e-5
)


Expand Down
8 changes: 6 additions & 2 deletions autotest/utilities/test_gdal_translate_lib.py
Original file line number Diff line number Diff line change
Expand Up @@ -1281,11 +1281,15 @@ def test_gdal_translate_ovr_rpc():
src_rpc = src_ds.GetMetadata("RPC")
ovr_rpc = ds.GetMetadata("RPC")
assert ovr_rpc
assert float(ovr_rpc["LINE_OFF"]) == pytest.approx(0.5 * float(src_rpc["LINE_OFF"]))
assert float(ovr_rpc["LINE_OFF"]) == pytest.approx(
0.5 * (float(src_rpc["LINE_OFF"]) + 0.5) - 0.5
)
assert float(ovr_rpc["LINE_SCALE"]) == pytest.approx(
0.5 * float(src_rpc["LINE_SCALE"])
)
assert float(ovr_rpc["SAMP_OFF"]) == pytest.approx(0.5 * float(src_rpc["SAMP_OFF"]))
assert float(ovr_rpc["SAMP_OFF"]) == pytest.approx(
0.5 * (float(src_rpc["SAMP_OFF"]) + 0.5) - 0.5
)
assert float(ovr_rpc["SAMP_SCALE"]) == pytest.approx(
0.5 * float(src_rpc["SAMP_SCALE"])
)
Expand Down
36 changes: 22 additions & 14 deletions gcore/gdaloverviewdataset.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -71,7 +71,8 @@ class GDALOverviewDataset final : public GDALDataset
GDALOverviewBand *m_poMaskBand = nullptr;

static void Rescale(char **&papszMD, const char *pszItem, double dfRatio,
double dfDefaultVal);
double dfDefaultVal, double dfPreShift = 0,
double dfPostShift = 0);

protected:
CPLErr IRasterIO(GDALRWFlag, int, int, int, int, void *, int, int,
Expand Down Expand Up @@ -455,12 +456,17 @@ const GDAL_GCP *GDALOverviewDataset::GetGCPs()
/* Rescale() */
/************************************************************************/

/* static */
void GDALOverviewDataset::Rescale(char **&papszMD, const char *pszItem,
double dfRatio, double dfDefaultVal)
double dfRatio, double dfDefaultVal,
double dfPreShift /*= 0*/,
double dfPostShift /*= 0*/)
{
double dfVal = CPLAtofM(CSLFetchNameValueDef(
papszMD, pszItem, CPLSPrintf("%.18g", dfDefaultVal)));
dfVal += dfPreShift;
dfVal *= dfRatio;
dfVal += dfPostShift;
papszMD = CSLSetNameValue(papszMD, pszItem, CPLSPrintf("%.18g", dfVal));
}

Expand All @@ -487,18 +493,20 @@ char **GDALOverviewDataset::GetMetadata(const char *pszDomain)
return papszMD_RPC;
papszMD_RPC = CSLDuplicate(papszMD);

Rescale(papszMD_RPC, RPC_LINE_OFF,
static_cast<double>(nRasterYSize) / poMainDS->GetRasterYSize(),
0.0);
Rescale(papszMD_RPC, RPC_LINE_SCALE,
static_cast<double>(nRasterYSize) / poMainDS->GetRasterYSize(),
1.0);
Rescale(papszMD_RPC, RPC_SAMP_OFF,
static_cast<double>(nRasterXSize) / poMainDS->GetRasterXSize(),
0.0);
Rescale(papszMD_RPC, RPC_SAMP_SCALE,
static_cast<double>(nRasterXSize) / poMainDS->GetRasterXSize(),
1.0);
const double dfXRatio =
static_cast<double>(nRasterXSize) / poMainDS->GetRasterXSize();
const double dfYRatio =
static_cast<double>(nRasterYSize) / poMainDS->GetRasterYSize();

// For line offset and pixel offset, we need to convert from RPC
// pixel center registration convention to GDAL pixel top-left corner
// registration convention by adding an initial 0.5 shift, and un-apply
// it after scaling.

Rescale(papszMD_RPC, RPC_LINE_OFF, dfYRatio, 0.0, 0.5, -0.5);
Rescale(papszMD_RPC, RPC_LINE_SCALE, dfYRatio, 1.0);
Rescale(papszMD_RPC, RPC_SAMP_OFF, dfXRatio, 0.0, 0.5, -0.5);
Rescale(papszMD_RPC, RPC_SAMP_SCALE, dfXRatio, 1.0);

papszMD = papszMD_RPC;
}
Expand Down

0 comments on commit 0266bc0

Please sign in to comment.