Skip to content
Snippets Groups Projects
Unverified Commit 88eb5896 authored by Austin Sanders's avatar Austin Sanders Committed by GitHub
Browse files

Cahvor fix (#490)

* Removed tolerance check from optical shift addition

* Updated changelog
parent b861a956
No related branches found
No related tags found
No related merge requests found
......@@ -35,6 +35,9 @@ release.
## [Unreleased]
### Fixed
- Fixed CAHVOR model optical shifts by removing tolerance check [#488](https://github.com/DOI-USGS/usgscsm/issues/488)
## [2.0.1] - 2024-01-23
### Changed
......@@ -51,4 +54,4 @@ release.
- Updated ALE submodule [#470](https://github.com/DOI-USGS/usgscsm/pull/470)
### Fixed
- Fixed issue with radial distortion computation [#464](https://github.com/DOI-USGS/usgscsm/pull/464)
\ No newline at end of file
- Fixed issue with radial distortion computation [#464](https://github.com/DOI-USGS/usgscsm/pull/464)
......@@ -328,16 +328,13 @@ void removeDistortion(double dx, double dy, double &ux, double &uy,
double shiftedDy = dy - opticalDistCoeffs[4];
double rr = shiftedDx * shiftedDx + shiftedDy * shiftedDy;
if (rr > tolerance)
{
double dr = opticalDistCoeffs[0] +
(rr * (opticalDistCoeffs[1] + rr * opticalDistCoeffs[2]));
double dr = opticalDistCoeffs[0] +
(rr * (opticalDistCoeffs[1] + rr * opticalDistCoeffs[2]));
ux = shiftedDx * (1.0 - dr);
uy = shiftedDy * (1.0 - dr);
ux += opticalDistCoeffs[3];
uy += opticalDistCoeffs[4];
}
ux = shiftedDx * (1.0 - dr);
uy = shiftedDy * (1.0 - dr);
ux += opticalDistCoeffs[3];
uy += opticalDistCoeffs[4];
}
break;
......@@ -587,51 +584,48 @@ void applyDistortion(double ux, double uy, double &dx, double &dy,
double shiftedUx = ux - opticalDistCoeffs[3];
double shiftedUy = uy - opticalDistCoeffs[4];
double rp2 = (ux * ux) + (uy * uy);
double rp = sqrt(rp2);
// Compute first fractional distortion using rp
double drOverR =
opticalDistCoeffs[0] +
(rp2 * (opticalDistCoeffs[1] + (rp2 * opticalDistCoeffs[2])));
if (rp2 > tolerance)
// Compute first distorted point estimate, r
double r = rp + (drOverR * rp);
double r_prev, r2_prev;
int iteration = 0;
do
{
double rp = sqrt(rp2);
// Compute first fractional distortion using rp
double drOverR =
opticalDistCoeffs[0] +
(rp2 * (opticalDistCoeffs[1] + (rp2 * opticalDistCoeffs[2])));
// Compute first distorted point estimate, r
double r = rp + (drOverR * rp);
double r_prev, r2_prev;
int iteration = 0;
do
// Don't get in an end-less loop. This algorithm should
// converge quickly. If not then we are probably way outside
// of the focal plane. Just set the distorted position to the
// undistorted position. Also, make sure the focal plane is less
// than 1km, it is unreasonable for it to grow larger than that.
if (iteration >= 20 || r > 1E9)
{
// Don't get in an end-less loop. This algorithm should
// converge quickly. If not then we are probably way outside
// of the focal plane. Just set the distorted position to the
// undistorted position. Also, make sure the focal plane is less
// than 1km, it is unreasonable for it to grow larger than that.
if (iteration >= 20 || r > 1E9)
{
drOverR = 0.0;
break;
}
r_prev = r;
r2_prev = r * r;
// Compute new fractional distortion:
drOverR = opticalDistCoeffs[0] +
(r2_prev *
(opticalDistCoeffs[1] + (r2_prev * opticalDistCoeffs[2])));
// Compute new estimate of r
r = rp + (drOverR * r_prev);
iteration++;
} while (fabs(r - r_prev) > desiredPrecision);
dx = shiftedUx / (1.0 - drOverR);
dy = shiftedUy / (1.0 - drOverR);
dx += opticalDistCoeffs[3];
dy += opticalDistCoeffs[4];
}
drOverR = 0.0;
break;
}
r_prev = r;
r2_prev = r * r;
// Compute new fractional distortion:
drOverR = opticalDistCoeffs[0] +
(r2_prev *
(opticalDistCoeffs[1] + (r2_prev * opticalDistCoeffs[2])));
// Compute new estimate of r
r = rp + (drOverR * r_prev);
iteration++;
} while (fabs(r - r_prev) > desiredPrecision);
dx = shiftedUx / (1.0 - drOverR);
dy = shiftedUy / (1.0 - drOverR);
dx += opticalDistCoeffs[3];
dy += opticalDistCoeffs[4];
}
break;
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment