WIP Improve detection of closed but split contours

src/gmt_constants.h

@@ -433,8 +433,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 {

src/gmt_plot.c

@@ -9755,7 +9755,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"};
@@ -9792,7 +9792,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);
 

src/grdcontour.c

@@ -535,7 +535,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 coordinates */
+	/* 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 */
@@ -670,8 +670,13 @@ GMT_LOCAL void grdcontour_sort_and_plot_ticks (struct GMT_CTRL *GMT, struct PSL_
 			continue;
 		}
 
-		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 */
+			/* 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;
+		}
 
 		if (mode & 1) {	/* Tick the innermost contour */
 			x_back = xp[np-1];	/* Last point along contour */

src/pscontour.c

@@ -351,8 +351,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);

test/grdcontour/closed.sh

@@ -1,5 +1,4 @@
 #!/usr/bin/env bash
-
 ps=closed.ps
 
 # Make a grid with closed contours at N pole, one crossing the periodic boundary, and one safely in middle
@@ -8,4 +7,3 @@ gmt grdmath -Rg -I1 0 0 SDIST KM2DEG 35 DIV 2 POW NEG EXP 0 90 SDIST KM2DEG 50 D
 contour="gmt grdcontour -A2 -C1 -L8.5/10.5 -Gd4 tmp.nc -Bxa60g30 -By30g30 -BWS -T+lLH -Wa1p,red -Wc1p,blue"
 $contour -JN180/7i -P -K > $ps
 $contour -JG30/35/5i -O -Y4i -X1i >> $ps
-