From b08b9580ab0aca70c8762b3f8f0039484ddaca60 Mon Sep 17 00:00:00 2001 From: Even Rouault Date: Wed, 20 Feb 2019 11:02:53 +0100 Subject: vgridshift: propagate multiplier to avoid false-positive detection of nodata values in the grids with US VERTCON grids that are in millimeters --- src/apply_vgridshift.cpp | 20 ++++++++++---------- 1 file changed, 10 insertions(+), 10 deletions(-) (limited to 'src/apply_vgridshift.cpp') diff --git a/src/apply_vgridshift.cpp b/src/apply_vgridshift.cpp index 61e0c528..951bdf6d 100644 --- a/src/apply_vgridshift.cpp +++ b/src/apply_vgridshift.cpp @@ -35,16 +35,16 @@ #include "proj_internal.h" #include "proj_internal.h" -static int is_nodata(float value) +static int is_nodata(float value, double vmultiplier) { /* nodata? */ /* GTX official nodata value if -88.88880f, but some grids also */ /* use other big values for nodata (e.g naptrans2008.gtx has */ /* nodata values like -2147479936), so test them too */ - return value > 1000 || value < -1000 || value == -88.88880f; + return value * vmultiplier > 1000 || value * vmultiplier < -1000 || value == -88.88880f; } -static double read_vgrid_value( PJ *defn, PJ_LP input, int *gridlist_count_p, PJ_GRIDINFO **tables, struct CTABLE *ct) { +static double read_vgrid_value( PJ *defn, PJ_LP input, double vmultiplier, int *gridlist_count_p, PJ_GRIDINFO **tables, struct CTABLE *ct) { int itable = 0; double value = HUGE_VAL; double grid_x, grid_y; @@ -129,28 +129,28 @@ static double read_vgrid_value( PJ *defn, PJ_LP input, int *gridlist_count_p, PJ double total_weight = 0.0; int n_weights = 0; value = 0.0f; - if( !is_nodata(value_a) ) + if( !is_nodata(value_a, vmultiplier) ) { double weight = (1.0-grid_x) * (1.0-grid_y); value += value_a * weight; total_weight += weight; n_weights ++; } - if( !is_nodata(value_b) ) + if( !is_nodata(value_b, vmultiplier) ) { double weight = (grid_x) * (1.0-grid_y); value += value_b * weight; total_weight += weight; n_weights ++; } - if( !is_nodata(value_c) ) + if( !is_nodata(value_c, vmultiplier) ) { double weight = (1.0-grid_x) * (grid_y); value += value_c * weight; total_weight += weight; n_weights ++; } - if( !is_nodata(value_d) ) + if( !is_nodata(value_d, vmultiplier) ) { double weight = (grid_x) * (grid_y); value += value_d * weight; @@ -218,7 +218,7 @@ int pj_apply_vgridshift( PJ *defn, const char *listname, input.phi = y[io]; input.lam = x[io]; - value = read_vgrid_value(defn, input, gridlist_count_p, tables, &ct); + value = read_vgrid_value(defn, input, 1.0, gridlist_count_p, tables, &ct); if( inverse ) z[io] -= value; @@ -310,7 +310,7 @@ int proj_vgrid_init(PJ* P, const char *grids) { } /***********************************************/ -double proj_vgrid_value(PJ *P, PJ_LP lp){ +double proj_vgrid_value(PJ *P, PJ_LP lp, double vmultiplier){ /*********************************************** Read grid value at position lp in grids loaded @@ -324,7 +324,7 @@ double proj_vgrid_value(PJ *P, PJ_LP lp){ double value; memset(&used_grid, 0, sizeof(struct CTABLE)); - value = read_vgrid_value(P, lp, &(P->vgridlist_geoid_count), P->vgridlist_geoid, &used_grid); + value = read_vgrid_value(P, lp, vmultiplier, &(P->vgridlist_geoid_count), P->vgridlist_geoid, &used_grid); proj_log_trace(P, "proj_vgrid_value: (%f, %f) = %f", lp.lam*RAD_TO_DEG, lp.phi*RAD_TO_DEG, value); return value; -- cgit v1.2.3 From 04844ac495f65e824a0bd9f9e49ea3360f2c063f Mon Sep 17 00:00:00 2001 From: Even Rouault Date: Wed, 20 Feb 2019 13:21:15 +0100 Subject: Apply multiplier in proj_vgrid_value() --- src/apply_vgridshift.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) (limited to 'src/apply_vgridshift.cpp') diff --git a/src/apply_vgridshift.cpp b/src/apply_vgridshift.cpp index 951bdf6d..ae23e39a 100644 --- a/src/apply_vgridshift.cpp +++ b/src/apply_vgridshift.cpp @@ -165,7 +165,7 @@ static double read_vgrid_value( PJ *defn, PJ_LP input, double vmultiplier, int * } - return value; + return value * vmultiplier; } /************************************************************************/ -- cgit v1.2.3 From 63224ed0e81b0e5795610c10e44578c932bbc33b Mon Sep 17 00:00:00 2001 From: Even Rouault Date: Tue, 26 Mar 2019 14:30:33 +0100 Subject: read_vgrid_value(): make it clear that nullptr deref of ct->cvs cannot happen. Coverity CID 193525 --- src/apply_vgridshift.cpp | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) (limited to 'src/apply_vgridshift.cpp') diff --git a/src/apply_vgridshift.cpp b/src/apply_vgridshift.cpp index ae23e39a..9be65f08 100644 --- a/src/apply_vgridshift.cpp +++ b/src/apply_vgridshift.cpp @@ -98,10 +98,13 @@ static double read_vgrid_value( PJ *defn, PJ_LP input, double vmultiplier, int * } /* load the grid shift info if we don't have it. */ - if( ct->cvs == nullptr && !pj_gridinfo_load( pj_get_ctx(defn), gi ) ) + if( ct->cvs == nullptr ) { - pj_ctx_set_errno( defn->ctx, PJD_ERR_FAILED_TO_LOAD_GRID ); - return PJD_ERR_FAILED_TO_LOAD_GRID; + if( !pj_gridinfo_load( pj_get_ctx(defn), gi ) || ct->cvs == nullptr ) + { + pj_ctx_set_errno( defn->ctx, PJD_ERR_FAILED_TO_LOAD_GRID ); + return PJD_ERR_FAILED_TO_LOAD_GRID; + } } -- cgit v1.2.3 From 095d2204f8bb05d172936aebbb1e9e44852c049f Mon Sep 17 00:00:00 2001 From: Chris Mayo Date: Fri, 29 Mar 2019 19:17:37 +0000 Subject: Remove duplicate instances of #include "proj_internal.h" Introduced by "Merge projects.h into proj_internal.h" 8ab6f683. --- src/apply_vgridshift.cpp | 1 - 1 file changed, 1 deletion(-) (limited to 'src/apply_vgridshift.cpp') diff --git a/src/apply_vgridshift.cpp b/src/apply_vgridshift.cpp index 9be65f08..e73c6da7 100644 --- a/src/apply_vgridshift.cpp +++ b/src/apply_vgridshift.cpp @@ -33,7 +33,6 @@ #include "proj_math.h" #include "proj_internal.h" -#include "proj_internal.h" static int is_nodata(float value, double vmultiplier) { -- cgit v1.2.3 From 7c543c30ba77b4a722b5b193f99af569b4ade41b Mon Sep 17 00:00:00 2001 From: Even Rouault Date: Tue, 16 Apr 2019 22:42:25 +0200 Subject: vgridshift: handle longitude wrap-around for grids with 360deg longitude extent Like egm96_15.gtx Fixes #1415 Technically, a similar fix could be done for horizontal grids, but world extent is less common for them. --- src/apply_vgridshift.cpp | 48 +++++++++++++++++++++++++++++++++++++++++------- 1 file changed, 41 insertions(+), 7 deletions(-) (limited to 'src/apply_vgridshift.cpp') diff --git a/src/apply_vgridshift.cpp b/src/apply_vgridshift.cpp index e73c6da7..377e36e0 100644 --- a/src/apply_vgridshift.cpp +++ b/src/apply_vgridshift.cpp @@ -28,6 +28,7 @@ #define PJ_LIB__ +#include #include #include @@ -62,9 +63,19 @@ static double read_vgrid_value( PJ *defn, PJ_LP input, double vmultiplier, int * ct = gi->ct; - /* skip tables that don't match our point at all. */ - if( ct->ll.phi > input.phi || ct->ll.lam > input.lam - || ct->ll.phi + (ct->lim.phi-1) * ct->del.phi < input.phi + /* skip tables that don't match our point at all (latitude check). */ + if( ct->ll.phi > input.phi + || ct->ll.phi + (ct->lim.phi-1) * ct->del.phi < input.phi ) + continue; + + bool fullWorldLongExtent = false; + if( fabs(ct->lim.lam * ct->del.lam - 2 * M_PI) < 1e-10 ) + { + fullWorldLongExtent = true; + } + + /* skip tables that don't match our point at all (longitude check). */ + else if( ct->ll.lam > input.lam || ct->ll.lam + (ct->lim.lam-1) * ct->del.lam < input.lam ) continue; @@ -77,8 +88,17 @@ static double read_vgrid_value( PJ *defn, PJ_LP input, double vmultiplier, int * { struct CTABLE *ct1 = child->ct; - if( ct1->ll.phi > input.phi || ct1->ll.lam > input.lam - || ct1->ll.phi+(ct1->lim.phi-1)*ct1->del.phi < input.phi + fullWorldLongExtent = false; + + if( ct1->ll.phi > input.phi + || ct1->ll.phi+(ct1->lim.phi-1)*ct1->del.phi < input.phi) + continue; + + if( fabs(ct1->lim.lam * ct1->del.lam - 2 * M_PI) < 1e-10 ) + { + fullWorldLongExtent = true; + } + else if( ct1->ll.lam > input.lam || ct1->ll.lam+(ct1->lim.lam-1)*ct1->del.lam < input.lam) continue; @@ -109,15 +129,29 @@ static double read_vgrid_value( PJ *defn, PJ_LP input, double vmultiplier, int * /* Interpolation a location within the grid */ grid_x = (input.lam - ct->ll.lam) / ct->del.lam; + if( fullWorldLongExtent ) { + // The first fmod goes to ]-lim, lim[ range + // So we add lim again to be in ]0, 2*lim[ and fmod again + grid_x = fmod( + fmod(grid_x + ct->lim.lam, ct->lim.lam) + ct->lim.lam, + ct->lim.lam); + } grid_y = (input.phi - ct->ll.phi) / ct->del.phi; grid_ix = lround(floor(grid_x)); + assert(grid_ix >= 0 && grid_ix < ct->lim.lam); grid_iy = lround(floor(grid_y)); + assert(grid_iy >= 0 && grid_iy < ct->lim.phi); grid_x -= grid_ix; grid_y -= grid_iy; grid_ix2 = grid_ix + 1; - if( grid_ix2 >= ct->lim.lam ) - grid_ix2 = ct->lim.lam - 1; + if( grid_ix2 >= ct->lim.lam ) { + if( fullWorldLongExtent ) { + grid_ix2 = 0; + } else { + grid_ix2 = ct->lim.lam - 1; + } + } grid_iy2 = grid_iy + 1; if( grid_iy2 >= ct->lim.phi ) grid_iy2 = ct->lim.phi - 1; -- cgit v1.2.3