diff options
| -rw-r--r-- | src/grids.cpp | 20 | ||||
| -rw-r--r-- | test/cli/ntv2_out.dist | 8 | ||||
| -rw-r--r-- | test/cli/td_out.dist | 2 | ||||
| -rwxr-xr-x | test/cli/testntv2 | 11 |
4 files changed, 30 insertions, 11 deletions
diff --git a/src/grids.cpp b/src/grids.cpp index 3007fedc..b4c3ad24 100644 --- a/src/grids.cpp +++ b/src/grids.cpp @@ -2360,12 +2360,14 @@ bool HorizontalShiftGridSet::reopen(PJ_CONTEXT *ctx) { // --------------------------------------------------------------------------- +#define REL_TOLERANCE_HGRIDSHIFT 1e-5 + const HorizontalShiftGrid *HorizontalShiftGrid::gridAt(double lon, double lat) const { for (const auto &child : m_children) { const auto &extentChild = child->extentAndRes(); - const double epsilon = - (extentChild.resLon + extentChild.resLat) / 10000.0; + const double epsilon = (extentChild.resLon + extentChild.resLat) * + REL_TOLERANCE_HGRIDSHIFT; if ((extentChild.fullWorldLongitude() || (lon + epsilon >= extentChild.westLon && lon - epsilon <= extentChild.eastLon)) && @@ -2385,7 +2387,8 @@ const HorizontalShiftGrid *HorizontalShiftGridSet::gridAt(double lon, return grid.get(); } const auto &extent = grid->extentAndRes(); - const double epsilon = (extent.resLon + extent.resLat) / 10000.0; + const double epsilon = + (extent.resLon + extent.resLat) * REL_TOLERANCE_HGRIDSHIFT; if ((extent.fullWorldLongitude() || (lon + epsilon >= extent.westLon && lon - epsilon <= extent.eastLon)) && @@ -2857,26 +2860,26 @@ static PJ_LP pj_hgrid_interpolate(PJ_LP t, const HorizontalShiftGrid *grid, frct.phi = t.phi - indx.phi; val.lam = val.phi = HUGE_VAL; if (indx.lam < 0) { - if (indx.lam == -1 && frct.lam > 0.99999999999) { + if (indx.lam == -1 && frct.lam > 1 - 10 * REL_TOLERANCE_HGRIDSHIFT) { ++indx.lam; frct.lam = 0.; } else return val; } else if ((in = indx.lam + 1) >= grid->width()) { - if (in == grid->width() && frct.lam < 1e-11) { + if (in == grid->width() && frct.lam < 10 * REL_TOLERANCE_HGRIDSHIFT) { --indx.lam; frct.lam = 1.; } else return val; } if (indx.phi < 0) { - if (indx.phi == -1 && frct.phi > 0.99999999999) { + if (indx.phi == -1 && frct.phi > 1 - 10 * REL_TOLERANCE_HGRIDSHIFT) { ++indx.phi; frct.phi = 0.; } else return val; } else if ((in = indx.phi + 1) >= grid->height()) { - if (in == grid->height() && frct.phi < 1e-11) { + if (in == grid->height() && frct.phi < 10 * REL_TOLERANCE_HGRIDSHIFT) { --indx.phi; frct.phi = 1.; } else @@ -2937,8 +2940,6 @@ static PJ_LP pj_hgrid_apply_internal(projCtx ctx, PJ_LP in, tb.lam -= extent->westLon; tb.phi -= extent->southLat; - tb.lam = adjlon(tb.lam - M_PI) + M_PI; - t = pj_hgrid_interpolate(tb, grid, true); if (grid->hasChanged()) { shouldRetry = gridset->reopen(ctx); @@ -2981,7 +2982,6 @@ static PJ_LP pj_hgrid_apply_internal(projCtx ctx, PJ_LP in, tb = in; tb.lam -= extent->westLon; tb.phi -= extent->southLat; - tb.lam = adjlon(tb.lam - M_PI) + M_PI; dif.lam = std::numeric_limits<double>::max(); dif.phi = std::numeric_limits<double>::max(); continue; diff --git a/test/cli/ntv2_out.dist b/test/cli/ntv2_out.dist index 650a69d8..ce866f37 100644 --- a/test/cli/ntv2_out.dist +++ b/test/cli/ntv2_out.dist @@ -12,3 +12,11 @@ Try with NTv2 and NTv1 together ... falls back to NTv1 ############################################################## Switching between NTv2 subgrids -112.5839956 49.4914451 0 -112.58307487 49.49145197 0.00000000 +############################################################## +Interpolating very close (and sometimes a bit outside) to the edges a NTv2 subgrid (#209) +-115.58333333 51.25000000 0 -115.58228512 51.24997866 0.00000000 +-115.58333333 51.25000010 0 -115.58228512 51.24997876 0.00000000 +-115.58333334 51.25000000 0 -115.58228513 51.24997866 0.00000000 +-115.49166667 51.07500000 0 -115.49062909 51.07497666 0.00000000 +-115.49166668 51.07500000 0 -115.49062910 51.07497666 0.00000000 +-115.49166667 51.07499990 0 -115.49062909 51.07497656 0.00000000 diff --git a/test/cli/td_out.dist b/test/cli/td_out.dist index 82b64321..76c6e6ec 100644 --- a/test/cli/td_out.dist +++ b/test/cli/td_out.dist @@ -21,7 +21,7 @@ edge or even a wee bit outside (#141). -5.5001 52.0 -5.500100000000 52.000000000000 0.000000000000 -5.5 52.0 -5.498893534472 52.000109529716 0.000000000000 -5.5000000000001 52.0000000000001 -5.498893534472 52.000109529717 0.000000000000 --5.4999 51.9999 -5.498793541695 52.000009529743 0.000000000000 +-5.4999 51.9999 -5.498793593803 52.000009531513 0.000000000000 -5.5001 52.0 -5.500100000000 52.000000000000 0.000000000000 ############################################################## NAD27 -> NAD83: 1st through ntv1, 2nd through conus diff --git a/test/cli/testntv2 b/test/cli/testntv2 index 2a31304e..d15ad7aa 100755 --- a/test/cli/testntv2 +++ b/test/cli/testntv2 @@ -60,6 +60,17 @@ $EXE +proj=latlong +datum=NAD83 +to +proj=latlong +ellps=clrk66 +nadgrids=ntv2_0 -112.5839956 49.4914451 0 EOF +echo "##############################################################" >> ${OUT} +echo "Interpolating very close (and sometimes a bit outside) to the edges a NTv2 subgrid (#209)" >> ${OUT} +$EXE +proj=latlong +datum=NAD83 +to +proj=latlong +ellps=clrk66 +nadgrids=ntv2_0.gsb -E -d 8 >>${OUT} <<EOF +-115.58333333 51.25000000 0 +-115.58333333 51.25000010 0 +-115.58333334 51.25000000 0 +-115.49166667 51.07500000 0 +-115.49166668 51.07500000 0 +-115.49166667 51.07499990 0 +EOF + # ############################################################################## # Done! |