Add -Ccpt option to sample colors from sample1d and grdxyz

doc/rst/source/dump_rgb.rst_

@@ -0,0 +1,33 @@
+**-C**\ [*section*/]\ *master*\|\ *cpt*\|\ *color*\ :math:`_1`,\ *color*\ :math:`_2`\ [,\ *color*\ :math:`_3`\ ,...]\ [**+h**\ [*hinge*]][**+i**\ *dz*][**+u**\|\ **U**\ *unit*][**+s**\ *fname*]
+    Determine the color components based on the *z* values.
+    Append name of a master CPT, an input CPT file or a comma-separated list of colors from which to build a CPT.
+    If no argument is given then under modern mode we select the current CPT, if it is available.
+    Generally, the input can be many things:
+
+    #. A standard GMT *master* CPT file, e.g., *earth* 
+       (see :ref:`Of Colors and Color Legends`) and can be either addressed
+       by *master* or *section*/*master* (without any **.cpt** extension).
+    #. File name of already custom-made *cpt* file (e.g, *my_colors.cpt*).
+    #. Build a linear continuous CPT from *color*\ :math:`_1`,\ *color*\ :math:`_2`\ [,\ *color*\ :math:`_3`\ ,...]
+       automatically, where *z* starts at 0 and is incremented by one for each color. In this case,
+       *color*\ :math:`_i` can be a *r*/*g*/*b* (e.g., 255/100/75), *h*-*s*-*v* triplet (e.g., 180-0.8-1),
+       a *c*/*m*/*y*/*k* quadruple (e.g., 80/50/40/30), an HTML hexadecimal color (e.g. *#aabbcc*),
+       or a :ref:`color name <RGBchart>`. No spaces between commas are allowed.
+
+    A few modifiers pertains to hinges and units:
+
+    - **+h** - If given a master CPT with soft hinges then you can enable the hinge at
+      data value *hinge* [0], whereas for hard-hinge CPTs you can adjust the location
+      of the hinge [0] but not disable it.
+    - **+i** - Append increment *dz* to quantize the grid range [Default uses the exact grid range].
+    - **+s** - Append *fname* to save the finalized CPT in a disk file.
+      This is useful when the CPT is generated automatically, but if used, **must** be
+      at the end of the **-C** option.
+    - **+u** - For any other *master* CPT, you may convert their *z*-values from 
+      other distance `Units`_ to meter by appending the original unit code.
+    - **+U** - Likewise, you may convert their *z*-values from meter to other `Units`_
+      by appending the desired unit code.
+
+    **Note**: We append four extra output columns at the end of the output records.  These
+    are the red, green, and blue components color components (all in 0-255 range) and the
+    transparency value (0-100% range, with 0 meaning opaque).

doc/rst/source/grd2xyz.rst

@@ -13,7 +13,8 @@ Synopsis
 .. include:: common_SYN_OPTs.rst_
 
 **gmt grd2xyz** *ingrid*
-[ |-C|\ [**f**\|\ **i**] ]
+[ |-C|\ [*section*/]\ *master*\|\ *cpt*\|\ *color*\ :math:`_1`,\ *color*\ :math:`_2`\ [,\ *color*\ :math:`_3`\ ,...]\ [**+h**\ [*hinge*]][**+i**\ *dz*][**+u**\|\ **U**\ *unit*][**+s**\ *fname*] ]
+[ |-F|\ [**f**\|\ **i**] ]
 [ |-L|\ [**c**\|\ **r**\|\ **x**\|\ **y**]\ *value* ]
 [ |SYN_OPT-R| ]
 [ |SYN_OPT-V| ]
@@ -40,7 +41,8 @@ precision of the ASCII output format by editing the
 **--FORMAT_FLOAT_OUT**\ =\ *format* on the command line, or choose binary
 output using single or double precision storage. As an option you may
 output z-values without the (*x, y*) coordinates (see |-Z| below) or you can
-save the grid in the STL format for 3-D printers.
+save the grid in the STL format for 3-D printers. Also, by giving a CPT via
+|-C| we will add *r*, *g*, *b*, *a* columns to the output based on *z* values.
 
 Required Arguments
 ------------------
@@ -55,7 +57,11 @@ Optional Arguments
 
 .. _-C:
 
-**-C**\ [**f**\|\ **i**]
+.. include:: dump_rgb.rst_
+
+.. _-F:
+
+**-F**\ [**f**\|\ **i**]
     Replace the x- and y-coordinates on output with the corresponding
     column and row numbers. These start at 0 (C-style counting); append
     **f** to start at 1 (FORTRAN-style counting). Alternatively, append

doc/rst/source/sample1d.rst

@@ -14,6 +14,7 @@ Synopsis
 
 **gmt sample1d** [ *table* ]
 [ |-A|\ [**f**\|\ **p**\|\ **m**\|\ **r**\|\ **R**][**+d**][**+l**] ]
+[ |-C|\ [*section*/]\ *master*\|\ *cpt*\|\ *color*\ :math:`_1`,\ *color*\ :math:`_2`\ [,\ *color*\ :math:`_3`\ ,...]\ [**+h**\ [*hinge*]][**+i**\ *dz*][**+u**\|\ **U**\ *unit*][**+s**\ *fname*] ]
 [ |-E| ]
 [ |-F|\ **a**\|\ **c**\|\ **e**\|\ **l**\|\ **n**\|\ **s**\ *p*\ [**+d1**\|\ **2**] ]
 [ |-N|\ *col* ]
@@ -49,7 +50,8 @@ Equidistant or arbitrary sampling can be selected. All columns
 are resampled based on the new sampling interval. Several interpolation
 schemes are available, in addition to a *smoothing* spline which trades off misfit
 for curvature. Extrapolation outside the range of the input data
-is not supported.
+is not supported. By giving a CPT via |-C| we will add *r*, *g*, *b*, *a* columns to the output
+based on the last input data column.
 
 Required Arguments
 ------------------
@@ -87,6 +89,10 @@ Optional Arguments
     **Note**: Calculation mode for loxodromes is spherical, hence **-je**
     cannot be used in combination with **+l**.
 
+.. _-C:
+
+.. include:: dump_rgb.rst_
+
 .. _-E:
 
 **-E**

src/grd2xyz.c

@@ -42,15 +42,21 @@ enum grd2xyz_mode {
 };
 
 struct GRD2XYZ_CTRL {
-	struct GRD2XYZ_C {	/* -C[f|i] */
+	struct GRD2XYZ_C {	/* -C<cpt> or -C<color1>,<color2>[,<color3>,...][+i<dz>] */
 		bool active;
-		unsigned int mode;
+		double dz;	/* Rounding for min/max determined from data */
+		char *file;
+		char *savecpt;	/* For when we want to save the automatically generated CPT */
 	} C;
 	struct GRD2XYZ_E {	/* -E[f][<nodata>] */
 		bool active;
 		bool floating;
 		double nodata;
 	} E;
+	struct GRD2XYZ_F {	/* -F[f|i] */
+		bool active;
+		unsigned int mode;
+	} F;
 	struct GRD2XYZ_L {	/* -Lc|r|x|y<value> */
 		bool active;
 		unsigned int mode;
@@ -90,23 +96,26 @@ static void *New_Ctrl (struct GMT_CTRL *GMT) {	/* Allocate and initialize a new
 
 static void Free_Ctrl (struct GMT_CTRL *GMT, struct GRD2XYZ_CTRL *C) {	/* Deallocate control structure */
 	if (!C) return;
+	gmt_M_str_free (C->C.file);
+	gmt_M_str_free (C->C.savecpt);
 	gmt_M_free (GMT, C);
 }
 
 static int usage (struct GMTAPI_CTRL *API, int level) {
 	const char *name = gmt_show_name_and_purpose (API, THIS_MODULE_LIB, THIS_MODULE_CLASSIC_NAME, THIS_MODULE_PURPOSE);
 	if (level == GMT_MODULE_PURPOSE) return (GMT_NOERROR);
-	GMT_Usage (API, 0, "usage: %s %s [-C[f|i]] [-Lc|r|x|y<value>] [%s] [-T[a|b][<base>]] [%s] [-W[a[+u<unit>]|<weight>]] "
+	GMT_Usage (API, 0, "usage: %s %s  [-C%s] [-F[f|i]] [-Lc|r|x|y<value>] [%s] [-T[a|b][<base>]] [%s] [-W[a[+u<unit>]|<weight>]] "
 		"[-Z[<flags>]] [%s] [%s] [%s] [%s] [%s] [%s] [%s] [%s] [%s]\n",
-		name, GMT_INGRID, GMT_Rgeo_OPT, GMT_V_OPT, GMT_bo_OPT, GMT_d_OPT, GMT_f_OPT, GMT_ho_OPT,
+		name, GMT_INGRID, CPT_OPT_ARGS, GMT_Rgeo_OPT, GMT_V_OPT, GMT_bo_OPT, GMT_d_OPT, GMT_f_OPT, GMT_ho_OPT,
 		GMT_o_OPT, GMT_qo_OPT, GMT_s_OPT, GMT_colon_OPT, GMT_PAR_OPT);
 
 	if (level == GMT_SYNOPSIS) return (GMT_MODULE_SYNOPSIS);
 
 	GMT_Message (API, GMT_TIME_NONE, "  REQUIRED ARGUMENTS:\n");
 	gmt_ingrid_syntax (API, 0, "Name of one or more grid files");
 	GMT_Message (API, GMT_TIME_NONE, "\n  OPTIONAL ARGUMENTS:\n");
-	GMT_Usage (API, 1, "\n-C[f|i]");
+	gmt_explain_cpt_input (API, 'C');
+	GMT_Usage (API, 1, "\n-F[f|i]");
 	GMT_Usage (API, -2, "Write row,col instead of x,y. Append f to start at 1, else 0 [Default]. "
 		"Use -Ci to write grid index instead of (x,y).");
 	GMT_Usage (API, 1, "\n-Lc|r|x|y<value>");
@@ -161,6 +170,7 @@ static int parse (struct GMT_CTRL *GMT, struct GRD2XYZ_CTRL *Ctrl, struct GMT_Z_
 	 */
 
 	unsigned int n_errors = 0, n_files = 0, k = 0;
+	char *f = NULL;
 	struct GMT_OPTION *opt = NULL;
 	struct GMTAPI_CTRL *API = GMT->parent;
 
@@ -175,11 +185,22 @@ static int parse (struct GMT_CTRL *GMT, struct GRD2XYZ_CTRL *Ctrl, struct GMT_Z_
 
 			/* Processes program-specific parameters */
 
-			case 'C':	/* Write row,col or index instead of x,y */
+			case 'C':	/* CPT */
+				if (opt->arg[0] && strchr ("if", opt->arg[0])) goto the_F;	/* Backwars compatibility for -Ci|f */
 				n_errors += gmt_M_repeated_module_option (API, Ctrl->C.active);
-				if (opt->arg[0] == 'c') Ctrl->C.mode = 0;
-				else if (opt->arg[0] == 'f') Ctrl->C.mode = 1;
-				else if (opt->arg[0] == 'i') Ctrl->C.mode = 2;
+				gmt_M_str_free (Ctrl->C.file);
+				if (opt->arg[0]) Ctrl->C.file = strdup (opt->arg);
+				if (opt->arg[0] && (f = gmt_strrstr (Ctrl->C.file, "+s")) != NULL) {	/* Filename has a +s<outname>, extract that part */
+					Ctrl->C.savecpt = &f[2];
+					f[0] = '\0';		/* Remove the +s<outname> from Ctrl->C.file */
+				}
+				gmt_cpt_interval_modifier (GMT, &(Ctrl->C.file), &(Ctrl->C.dz));
+				break;
+the_F:		case 'F':	/* Write row,col or index instead of x,y */
+				n_errors += gmt_M_repeated_module_option (API, Ctrl->F.active);
+				if (opt->arg[0] == 'c') Ctrl->F.mode = 0;
+				else if (opt->arg[0] == 'f') Ctrl->F.mode = 1;
+				else if (opt->arg[0] == 'i') Ctrl->F.mode = 2;
 				break;
 			case 'E':	/* Old ESRI option */
 				if (gmt_M_compat_check (GMT, 4)) {
@@ -412,7 +433,7 @@ void place_base_triangles (struct GMT_CTRL *GMT, FILE *fp, struct GMT_GRID *G, b
 #define Return(code) {Free_Ctrl (GMT, Ctrl); gmt_end_module (GMT, GMT_cpy); bailout (code);}
 
 EXTERN_MSC int GMT_grd2xyz (void *V_API, int mode, void *args) {
-	bool first = true, first_geo = true;
+	bool first = true, first_geo = true, first_grd = true;
 	openmp_int row, col;
 	unsigned int n_output, w_col = 3, orig_mode;
 	int error = 0, write_error = 0;
@@ -421,13 +442,15 @@ EXTERN_MSC int GMT_grd2xyz (void *V_API, int mode, void *args) {
 
 	char header[GMT_BUFSIZ];
 
-	double wesn[4], d_value, A_scale = 1.0, A_geo = 0.0, out[4], *x = NULL, *y = NULL;
+	double wesn[4], d_value, A_scale = 1.0, A_geo = 0.0, out[7], *x = NULL, *y = NULL;
+	double rgb[4] = {0.0, 0.0, 0.0, 0.0};
 
 	struct GMT_GRID *G = NULL, *W = NULL;
 	struct GMT_GRID_HEADER_HIDDEN *H = NULL;
 	struct GMT_RECORD *Out = NULL;
 	struct GMT_Z_IO io;
 	struct GMT_OPTION *opt = NULL;
+	struct GMT_PALETTE *P = NULL;
 	struct GRD2XYZ_CTRL *Ctrl = NULL;
 	struct GMT_CTRL *GMT = NULL, *GMT_cpy = NULL;
 	struct GMT_OPTION *options = NULL;
@@ -466,13 +489,15 @@ EXTERN_MSC int GMT_grd2xyz (void *V_API, int mode, void *args) {
 	}
 	else if (io.binary) GMT->common.b.active[GMT_OUT] = true;
 
-	n_output = (Ctrl->Z.active) ? 1 : ((Ctrl->W.active) ? 4 : ((Ctrl->C.mode == 2) ? 2 : 3));
+	n_output = (Ctrl->Z.active) ? 1 : ((Ctrl->W.active) ? 4 : ((Ctrl->F.mode == 2) ? 2 : 3));
 	if (Ctrl->Z.active)
 		n_output = 1;
+	else if (Ctrl->C.active)	/* x, y, z, r, g, b, a */
+		n_output = 7;
 	else {
-		n_output = (Ctrl->C.mode == 2) ? 2 : 3;
+		n_output = (Ctrl->F.mode == 2) ? 2 : 3;
 		if (Ctrl->W.active) n_output++;
-		if (Ctrl->C.mode == 2) w_col = 2;
+		if (Ctrl->F.mode == 2) w_col = 2;
 	}
 	if ((error = GMT_Set_Columns (API, GMT_OUT, n_output, GMT_COL_FIX_NO_TEXT)) != GMT_NOERROR) Return (error);
 
@@ -488,6 +513,14 @@ EXTERN_MSC int GMT_grd2xyz (void *V_API, int mode, void *args) {
 
 	out[w_col] = Ctrl->W.weight;
 
+	if (Ctrl->C.active) {	/* Just sample the CPT for the grid z values and append three columns r g b a to the x y z */
+		if ((P = gmt_get_palette (GMT, Ctrl->C.file, GMT_CPT_OPTIONAL, 0.0, 255.0, Ctrl->C.dz)) == NULL) {
+			GMT_Report (API, GMT_MSG_ERROR, "Failed to read CPT %s.\n", Ctrl->C.file);
+			Return (API->error);	/* Well, that did not go so well... */
+		}
+		Out = gmt_new_record (GMT, out, NULL);	/* Since we only need to worry about numerics in this module */
+	}
+
 	for (opt = options; opt; opt = opt->next) {	/* Loop over arguments, skip options */
 
 		if (opt->option != '<') continue;	/* We are only processing input files here */
@@ -505,6 +538,25 @@ EXTERN_MSC int GMT_grd2xyz (void *V_API, int mode, void *args) {
 			Return (API->error);	/* Get subset */
 		}
 
+		if (Ctrl->C.active) {	/* Just sample the CPT for the z value and append four columns r g b a */
+			unsigned int k, geo = gmt_M_is_geographic (GMT, GMT_IN);
+			static char *coord[2] = {"x\ty\tz\tred\tgreen\tblue\ttransparency", "lon\tlat\tz\tred\tgreen\tblue\ttransparency"};
+			if (first_grd) GMT_Put_Record (API, GMT_WRITE_TABLE_HEADER, coord[geo]);
+			first_grd = false;
+			x = (G->x) ? G->x : gmt_grd_coord (GMT, G->header, GMT_X);
+			y = (G->y) ? G->y : gmt_grd_coord (GMT, G->header, GMT_Y);
+			gmt_M_grd_loop (GMT, G, row, col, ij) {
+				out[GMT_X] = x[col];	out[GMT_Y] = y[row];	out[GMT_Z] = G->data[ij];
+				gmt_get_rgb_from_z (GMT, P, G->data[ij], rgb);
+				for (k = 0; k < 4; k++) out[GMT_Z + 1 + k] = rgb[k];
+				write_error = GMT_Put_Record (API, GMT_WRITE_DATA, Out);
+				if (write_error == GMT_NOTSET) n_suppressed++;	/* Bad value caught by -s[r] */
+			}
+			if (GMT_Destroy_Data (API, &G) != GMT_NOERROR) {	/* Free the grid since there may be more of them */
+				Return (API->error);
+			}
+			continue;	/* Back to next grid */
+		}
 		if (Ctrl->L.active) {	/* Check that limits are within range and warn otherwise */
 			orig_mode = Ctrl->L.mode;	/* In case we have many grids */
 			if (Ctrl->L.mode == GRD2XYZ_Y) {	/* Convert y coordinate to nearest row */
@@ -735,10 +787,10 @@ EXTERN_MSC int GMT_grd2xyz (void *V_API, int mode, void *args) {
 			x = (G->x) ? G->x : gmt_grd_coord (GMT, G->header, GMT_X);
 			y = (G->y) ? G->y : gmt_grd_coord (GMT, G->header, GMT_Y);
 			Out = gmt_new_record (GMT, out, NULL);	/* Since we only need to worry about numerics in this module */
-			if (Ctrl->C.active) {	/* Replace x,y with col,row */
-				if (Ctrl->C.mode < 2) {
-					gmt_M_row_loop  (GMT, G, row) y[row] = row + Ctrl->C.mode;
-					gmt_M_col_loop2 (GMT, G, col) x[col] = col + Ctrl->C.mode;
+			if (Ctrl->F.active) {	/* Replace x,y with col,row */
+				if (Ctrl->F.mode < 2) {
+					gmt_M_row_loop  (GMT, G, row) y[row] = row + Ctrl->F.mode;
+					gmt_M_col_loop2 (GMT, G, col) x[col] = col + Ctrl->F.mode;
 				}
 				else
 					GMT->current.io.io_nan_col[0] = GMT_Y;	/* Since that is where z will go now */
@@ -748,7 +800,7 @@ EXTERN_MSC int GMT_grd2xyz (void *V_API, int mode, void *args) {
 			if (GMT->current.setting.io_header[GMT_OUT] && first) {
 				if (!G->header->x_units[0]) strcpy (G->header->x_units, "x");
 				if (!G->header->y_units[0]) strcpy (G->header->y_units, "y");
-				if (Ctrl->C.active) {
+				if (Ctrl->F.active) {
 					strcpy (G->header->x_units, "col");
 					strcpy (G->header->y_units, "row");
 				}
@@ -772,7 +824,7 @@ EXTERN_MSC int GMT_grd2xyz (void *V_API, int mode, void *args) {
 					if (Ctrl->L.mode == GRD2XYZ_Y && !doubleAlmostEqualZero (y[row], Ctrl->L.value)) continue;	/* Skip these rows */
 					if (Ctrl->L.mode == GRD2XYZ_X && !doubleAlmostEqualZero (x[col], Ctrl->L.value)) continue;	/* Skip these columns */
 				}
-				if (Ctrl->C.mode == 2) {	/* Write index, z */
+				if (Ctrl->F.mode == 2) {	/* Write index, z */
 					out[GMT_X] = (double)gmt_M_ij0 (G->header, row, col);
 					out[GMT_Y] = G->data[ij];
 					if (GMT->common.d.active[GMT_OUT] && gmt_M_is_dnan (out[GMT_Y]))	/* Input matched no-data setting, so change to NaN */
@@ -807,6 +859,8 @@ EXTERN_MSC int GMT_grd2xyz (void *V_API, int mode, void *args) {
 		Return (API->error);
 	}
 
+	if (Ctrl->C.active) gmt_M_free (GMT, Out);
+
 	if (!Ctrl->T.active) GMT_Report (API, GMT_MSG_INFORMATION, "%" PRIu64 " values extracted\n", n_total - n_suppressed);
 	if (n_suppressed) {
 		if (GMT->current.setting.io_nan_mode & GMT_IO_NAN_KEEP)

src/longopt/grd2xyz_inc.h

@@ -25,7 +25,8 @@ static struct GMT_KEYWORD_DICTIONARY module_kw[] = {
 		  short_directives,    long_directives,
 		  short_modifiers,     long_modifiers,
 		  transproc_mask */
-	{ 0, 'C', "",
+	GMT_C_CPT_KW,
+	{ 0, 'F', "",
 	          "",                  "",
 	          "",                  "",
 		  GMT_TP_STANDARD },

src/longopt/sample1d_inc.h

@@ -29,6 +29,7 @@ static struct GMT_KEYWORD_DICTIONARY module_kw[] = {
 	          "",                  "",
 	          "",                  "",
 		  GMT_TP_STANDARD },
+	GMT_C_CPT_KW,
 	{ 0, 'E', "",
 	          "",                  "",
 	          "",                  "",