From b85239dfacb343220a932f476f295038ed1616ff Mon Sep 17 00:00:00 2001 From: Paul Wessel Date: Thu, 2 Dec 2021 15:35:31 -1000 Subject: [PATCH 1/4] Update grdcontour.c --- src/grdcontour.c | 23 +++++++++++++++++++++-- 1 file changed, 21 insertions(+), 2 deletions(-) diff --git a/src/grdcontour.c b/src/grdcontour.c index ed329030dad..e3c8e2cb256 100644 --- a/src/grdcontour.c +++ b/src/grdcontour.c @@ -521,6 +521,19 @@ static int parse (struct GMT_CTRL *GMT, struct GRDCONTOUR_CTRL *Ctrl, struct GMT /* Three sub functions used by GMT_grdcontour */ +#ifdef DEBUG +/* If we need to dump out a polygon then */ +GMT_LOCAL void dumppol (uint64_t n, double *x, double *y, int id) { + uint64_t i; + FILE *fp = NULL; + char line[GMT_LEN64]; + snprintf (line, GMT_LEN64, "dump_%d.txt", id); + fp = fopen (line, "w"); + for (i = 0; i < n; i++) fprintf (fp, "%g\t%g\n", x[i], y[i]); + fclose (fp); +} +#endif + GMT_LOCAL void grdcontour_sort_and_plot_ticks (struct GMT_CTRL *GMT, struct PSL_CTRL *PSL, struct GRDCONTOUR_SAVE *save, size_t n, struct GMT_GRID *G, struct CONTOUR_CLOSED *I, unsigned int mode, struct GMT_DATASET *T) { /* Labeling and ticking of inner-most contours cannot happen until all contours are found and we can determine which are the innermost ones. Here, all the n candidate contours are passed via the save array. @@ -532,7 +545,7 @@ GMT_LOCAL void grdcontour_sort_and_plot_ticks (struct GMT_CTRL *GMT, struct PSL_ Note on mode bitflags: mode = 1 (plot only), 2 (save labels only), 3 (both) */ - int np, j, k, inside, col, row, stop, n_ticks, way, form; + int np, j, k, inside, col, row, stop, n_ticks, way, form, dump = 1; uint64_t ij; size_t pol, pol2; bool done, match, found; @@ -547,7 +560,7 @@ GMT_LOCAL void grdcontour_sort_and_plot_ticks (struct GMT_CTRL *GMT, struct PSL_ lbl[0] = (I->txt[0]) ? I->txt[0] : def[0]; lbl[1] = (I->txt[1]) ? I->txt[1] : def[1]; - /* The x/y coordinates in SAVE in original cooordinates */ + /* The x/y coordinates in SAVE in original cooordinates and hence are closed polygons in lon/lat (if geographic) */ for (pol = 0; pol < n; pol++) { /* Set y min/max for polar caps */ if (abs (save[pol].kind) < 3) continue; /* Skip all but polar caps */ @@ -666,6 +679,11 @@ GMT_LOCAL void grdcontour_sort_and_plot_ticks (struct GMT_CTRL *GMT, struct PSL_ np = (int)gmt_clip_to_map (GMT, save[pol].x, save[pol].y, np, &xp, &yp); /* Convert to inches */ if (np == 0) continue; + if (gmt_polygon_is_open (GMT, xp, yp, np)) { /* Closed polygon got projected to an open one */ + save[pol].do_it = false; + gmt_M_free (GMT, s); gmt_M_free (GMT, xp); gmt_M_free (GMT, yp); + continue; + } s = gmt_M_memory (GMT, NULL, np, double); /* Compute distance along the contour */ for (j = 1, s[0] = 0.0; j < np; j++) @@ -677,6 +695,7 @@ GMT_LOCAL void grdcontour_sort_and_plot_ticks (struct GMT_CTRL *GMT, struct PSL_ continue; } + if (dump) dumppol (np, xp, yp, pol); way = gmt_polygon_centroid (GMT, xp, yp, np, &save[pol].xlabel, &save[pol].ylabel); /* -1 is CCW, +1 is CW */ /* Compute mean location of closed contour ~hopefully a good point inside to place label. */ From 07b5c1093217b27bdee3e9f2b770e545905ac9b8 Mon Sep 17 00:00:00 2001 From: Paul Wessel Date: Thu, 2 Dec 2021 16:12:24 -1000 Subject: [PATCH 2/4] Update constants --- src/gmt_constants.h | 4 ++-- src/gmt_plot.c | 4 ++-- src/grdcontour.c | 29 +++++++---------------------- src/pscontour.c | 7 ++++++- 4 files changed, 17 insertions(+), 27 deletions(-) diff --git a/src/gmt_constants.h b/src/gmt_constants.h index 0b963b4be0a..1d1c25fdf2b 100644 --- a/src/gmt_constants.h +++ b/src/gmt_constants.h @@ -374,8 +374,8 @@ enum GMT_enum_inside { /*! Return codes from gmt_polygon_orientation */ enum GMT_enum_polorient { - GMT_POL_IS_CCW = 0, - GMT_POL_IS_CW = 1}; + GMT_POL_IS_CCW = -1, + GMT_POL_IS_CW = +1}; /*! Codes for -q selections */ enum GMT_enum_skiprows { diff --git a/src/gmt_plot.c b/src/gmt_plot.c index baa9cff4fa5..e60c5340ce9 100644 --- a/src/gmt_plot.c +++ b/src/gmt_plot.c @@ -9421,7 +9421,7 @@ void gmt_geo_polygons (struct GMT_CTRL *GMT, struct GMT_DATASEGMENT *S) { struct GMT_DATASEGMENT *S2 = NULL; struct GMT_DATASEGMENT_HIDDEN *SH = gmt_get_DS_hidden (S); bool add_pole, separate; - int outline = 0, P_handedness = GMT_NOTSET, H_handedness; + int outline = 0, P_handedness = 0, H_handedness; int geo = gmt_M_is_geographic (GMT, GMT_IN); uint64_t used = 0; char *type[2] = {"Perimeter", "Polar cap perimeter"}; @@ -9458,7 +9458,7 @@ void gmt_geo_polygons (struct GMT_CTRL *GMT, struct GMT_DATASEGMENT *S) { if (PSL->internal.comments) snprintf (comment, GMT_LEN64, "%s polygon for %s\n", type[add_pole], use[PSL->current.outline]); used = gmtplot_geo_polygon_segment (GMT, S, add_pole, true, comment); /* First lay down perimeter */ for (S2 = gmt_get_next_S (S); S2; S2 = gmt_get_next_S (S2)) { /* Process all holes [none processed if there aren't any holes] */ - if (P_handedness == GMT_NOTSET) + if (P_handedness == 0) P_handedness = gmt_polygon_orientation (GMT, S->data[GMT_X], S->data[GMT_Y], S->n_rows, geo); H_handedness = gmt_polygon_orientation (GMT, S2->data[GMT_X], S2->data[GMT_Y], S2->n_rows, geo); diff --git a/src/grdcontour.c b/src/grdcontour.c index c051e54ff0e..7bbb2fd0c08 100644 --- a/src/grdcontour.c +++ b/src/grdcontour.c @@ -521,19 +521,6 @@ static int parse (struct GMT_CTRL *GMT, struct GRDCONTOUR_CTRL *Ctrl, struct GMT /* Three sub functions used by GMT_grdcontour */ -#ifdef DEBUG -/* If we need to dump out a polygon then */ -GMT_LOCAL void dumppol (uint64_t n, double *x, double *y, int id) { - uint64_t i; - FILE *fp = NULL; - char line[GMT_LEN64]; - snprintf (line, GMT_LEN64, "dump_%d.txt", id); - fp = fopen (line, "w"); - for (i = 0; i < n; i++) fprintf (fp, "%g\t%g\n", x[i], y[i]); - fclose (fp); -} -#endif - GMT_LOCAL void grdcontour_sort_and_plot_ticks (struct GMT_CTRL *GMT, struct PSL_CTRL *PSL, struct GRDCONTOUR_SAVE *save, size_t n, struct GMT_GRID *G, struct CONTOUR_CLOSED *I, unsigned int mode, struct GMT_DATASET *T) { /* Labeling and ticking of inner-most contours cannot happen until all contours are found and we can determine which are the innermost ones. Here, all the n candidate contours are passed via the save array. @@ -545,7 +532,7 @@ GMT_LOCAL void grdcontour_sort_and_plot_ticks (struct GMT_CTRL *GMT, struct PSL_ Note on mode bitflags: mode = 1 (plot only), 2 (save labels only), 3 (both) */ - int np, j, k, inside, col, row, stop, n_ticks, way, form, dump = 1; + int np, j, k, inside, col, row, stop, n_ticks, way, form; uint64_t ij; size_t pol, pol2; bool done, match, found; @@ -560,7 +547,7 @@ GMT_LOCAL void grdcontour_sort_and_plot_ticks (struct GMT_CTRL *GMT, struct PSL_ lbl[0] = (I->txt[0]) ? I->txt[0] : def[0]; lbl[1] = (I->txt[1]) ? I->txt[1] : def[1]; - /* The x/y coordinates in SAVE in original cooordinates and hence are closed polygons in lon/lat (if geographic) */ + /* The x/y coordinates in SAVE in original coordinates and hence are closed polygons in lon/lat (if geographic) */ for (pol = 0; pol < n; pol++) { /* Set y min/max for polar caps */ if (abs (save[pol].kind) < 3) continue; /* Skip all but polar caps */ @@ -679,11 +666,6 @@ GMT_LOCAL void grdcontour_sort_and_plot_ticks (struct GMT_CTRL *GMT, struct PSL_ np = (int)gmt_clip_to_map (GMT, save[pol].x, save[pol].y, np, &xp, &yp); /* Convert to inches */ if (np == 0) continue; - if (gmt_polygon_is_open (GMT, xp, yp, np)) { /* Closed polygon got projected to an open one */ - save[pol].do_it = false; - gmt_M_free (GMT, s); gmt_M_free (GMT, xp); gmt_M_free (GMT, yp); - continue; - } s = gmt_M_memory (GMT, NULL, np, double); /* Compute distance along the contour */ for (j = 1, s[0] = 0.0; j < np; j++) @@ -695,9 +677,12 @@ GMT_LOCAL void grdcontour_sort_and_plot_ticks (struct GMT_CTRL *GMT, struct PSL_ continue; } - if (dump) dumppol (np, xp, yp, pol); - way = gmt_polygon_centroid (GMT, xp, yp, np, &save[pol].xlabel, &save[pol].ylabel); /* -1 is CCW, +1 is CW */ /* Compute mean location of closed contour ~hopefully a good point inside to place label. */ + way = gmt_polygon_centroid (GMT, xp, yp, np, &save[pol].xlabel, &save[pol].ylabel); /* -1 is CCW, +1 is CW */ + if (way == GMT_POL_IS_CW) { /* So far this has been found to be the wrong way so we switch */ + GMT_Report (GMT->parent, GMT_MSG_WARNING, "gmt_polygon_centroid found CW polygon (by mistake?), switch to CCW\n"); + way = GMT_POL_IS_CCW; + } if (mode & 1) { /* Tick the innermost contour */ x_back = xp[np-1]; /* Last point along contour */ diff --git a/src/pscontour.c b/src/pscontour.c index c9f670463c8..79dc2d4a901 100644 --- a/src/pscontour.c +++ b/src/pscontour.c @@ -345,8 +345,13 @@ GMT_LOCAL void pscontour_sort_and_plot_ticks (struct GMT_CTRL *GMT, struct PSL_C if (n_ticks == 0) continue; /* Too short to be ticked or labeled */ gmt_setpen (GMT, &save[pol].pen); + /* Compute mean location of closed contour ~hopefully a good point inside to place label. */ way = gmt_polygon_centroid (GMT, save[pol].x, save[pol].y, np, &x_mean, &y_mean); /* -1 is CCW, +1 is CW */ - if (I->label) { /* Compute mean location of closed contour ~hopefully a good point inside to place label. */ + if (way == GMT_POL_IS_CW) { /* So far this has been found to be the wrong way so we switch */ + GMT_Report (GMT->parent, GMT_MSG_WARNING, "gmt_polygon_centroid found CW polygon (by mistake?), switch to CCW\n"); + way = GMT_POL_IS_CCW; + } + if (I->label) { if (mode & 1) { form = gmt_setfont (GMT, &save[pol].font); PSL_plottext (PSL, x_mean, y_mean, GMT->current.setting.font_annot[GMT_PRIMARY].size, lbl[save[pol].high], 0.0, PSL_MC, form); From 03054c3f37366baa2bd1db62da866fc9894fadea Mon Sep 17 00:00:00 2001 From: Paul Wessel Date: Thu, 2 Dec 2021 16:13:34 -1000 Subject: [PATCH 3/4] Update grdcontour.c --- src/grdcontour.c | 1 + 1 file changed, 1 insertion(+) diff --git a/src/grdcontour.c b/src/grdcontour.c index 7bbb2fd0c08..82d3cac6134 100644 --- a/src/grdcontour.c +++ b/src/grdcontour.c @@ -680,6 +680,7 @@ GMT_LOCAL void grdcontour_sort_and_plot_ticks (struct GMT_CTRL *GMT, struct PSL_ /* Compute mean location of closed contour ~hopefully a good point inside to place label. */ way = gmt_polygon_centroid (GMT, xp, yp, np, &save[pol].xlabel, &save[pol].ylabel); /* -1 is CCW, +1 is CW */ if (way == GMT_POL_IS_CW) { /* So far this has been found to be the wrong way so we switch */ + /* See https://github.com/GenericMappingTools/gmt/issues/6080 which used way as is for 13000 km contour */ GMT_Report (GMT->parent, GMT_MSG_WARNING, "gmt_polygon_centroid found CW polygon (by mistake?), switch to CCW\n"); way = GMT_POL_IS_CCW; } From 5d91ebf762d8e3fe31e4e8abb21bd6bafa7dad74 Mon Sep 17 00:00:00 2001 From: Paul Wessel Date: Thu, 2 Dec 2021 16:21:19 -1000 Subject: [PATCH 4/4] Update closed.sh --- test/grdcontour/closed.sh | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/test/grdcontour/closed.sh b/test/grdcontour/closed.sh index d00cc1fd759..087b16adf51 100755 --- a/test/grdcontour/closed.sh +++ b/test/grdcontour/closed.sh @@ -1,6 +1,7 @@ #!/usr/bin/env bash # - +# GRAPHICSMAGICK_RMS = 0.02 +# See https://github.com/GenericMappingTools/gmt/issues/6080 for why we up the RMS here (to avoid updating PS with minor changes) ps=closed.ps # Make a grid with closed contours at N pole, one crossing the periodic boundary, and one safely in middle