Skip to content

Introduce directives -Mc and -Mf in grdinfo #6767

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
PaulWessel merged 4 commits into from
Jun 6, 2022
Merged

Introduce directives -Mc and -Mf in grdinfo #6767

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
4 commits
Select commit Hold shift + click to select a range
a3907c6
Introduce directives -Mc and -Mf in grdinfo
PaulWessel Jun 5, 2022
598985c
Use enums
PaulWessel Jun 5, 2022
ba294ea
polish
PaulWessel Jun 5, 2022
1b3052d
Update grdinfo.c
PaulWessel Jun 5, 2022
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
13 changes: 8 additions & 5 deletions doc/rst/source/grdinfo.rst
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -19,7 +19,7 @@ Synopsis
[ |-G| ]
[ |-I|\ [*dx*\ [/*dy*]\|\ **b**\|\ **i**\|\ **r**] ]
[ |-L|\ [**0**\|\ **1**\|\ **2**\|\ **p**\|\ **a**] ]
[ |-M| ]
[ |-M|\ [**c**\|\ **f**] ]
[ |SYN_OPT-R| ]
[ |-T|\ [*dv*]\ [**+a**\ [*alpha*]]\ [**+s**] ]
[ |SYN_OPT-V| ]
Expand Down Expand Up @@ -65,8 +65,8 @@ Optional Arguments
*name w e s n {b t} v0 v1 dx dy {dz} nx ny {nz}* [*x0 y0 {z0} x1 y1 {z1}*] [*med scale*] [*mean std rms*] [*n\_nan*] *registration gtype*

The data in brackets are
output only if the corresponding options |-M|, **-L1**, **-L2**,
and |-M| are used, respectively, while the data in braces only apply if
output only if the corresponding options |-M| (with no directive), **-L1**, and **-L2**
are used, respectively, while the data in braces only apply if
used with 3-D data cubes. Use **-Ct** to place file *name*
at the end of the output record or **-Cn** to only output numerical
columns. The *registration* is either 0 (gridline) or 1 (pixel),
Expand Down Expand Up @@ -145,9 +145,12 @@ Optional Arguments

.. _-M:

**-M**
**-M**\ [**c**\|\ **f**]
Find and report the location of min/max *v*-values, and count and
report the number of nodes set to NaN, if any.
report the number of nodes set to NaN, if any [Default]. Use directive
**f** to instead force an update of the *v*-value min/max by reading the
matrix, or use **c** for conditionally doing so if the header information
does not contain a valid *v* range.

.. |Add_-R| replace:: Using the |-R| option will select a subsection of the input grid(s). If this subsection
exceeds the boundaries of the grid, only the common region will be extracted. For cubes you must also
Expand Down
59 changes: 44 additions & 15 deletions src/grdinfo.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -48,6 +48,11 @@ enum Opt_C_modes {
GRDINFO_NUMERICAL = 1,
GRDINFO_TRAILING = 2};

enum Opt_M_modes {
GRDINFO_FORCE_REPORT = 1,
GRDINFO_FORCE = 2,
GRDINFO_CONDITIONAL = 3};

struct GRDINFO_CTRL {
unsigned int n_files; /* How many grids given */
bool is_cube;
Expand Down Expand Up @@ -77,8 +82,9 @@ struct GRDINFO_CTRL {
unsigned int status;
double inc[2];
} I;
struct GRDINFO_M { /* -M */
struct GRDINFO_M { /* -M[c|f] */
bool active;
unsigned int mode; /* 1 is force, 2 is conditional */
} M;
struct GRDINFO_L { /* -L[0|1|2|p] */
bool active;
Expand Down Expand Up @@ -112,7 +118,7 @@ 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[n|t]] [-D[<offx>[/<offy>]][+i]] [-E[x|y][+l|L|u|U]] [-F] [-G] [-I[<dx>[/<dy>]|b|i|r]] [-L[a|0|1|2|p]] "
"[-M] [%s] [-T[<dv>][+a[<alpha>]][+s]] [%s] [%s] [%s] [%s] [%s]\n", name, GMT_INGRID, GMT_Rgeo_OPT, GMT_V_OPT, GMT_f_OPT, GMT_ho_OPT, GMT_o_OPT, GMT_PAR_OPT);
"[-M[c|f]] [%s] [-T[<dv>][+a[<alpha>]][+s]] [%s] [%s] [%s] [%s] [%s]\n", name, GMT_INGRID, GMT_Rgeo_OPT, GMT_V_OPT, GMT_f_OPT, GMT_ho_OPT, GMT_o_OPT, GMT_PAR_OPT);

if (level == GMT_SYNOPSIS) return (GMT_MODULE_SYNOPSIS);

Expand Down Expand Up @@ -158,7 +164,9 @@ static int usage (struct GMTAPI_CTRL *API, int level) {
GMT_Usage (API, 3, "p: Report mode (lms) and LMS-scale (MAD w.r.t. mode) of data set.");
GMT_Usage (API, 3, "a: All of the above.");
GMT_Usage (API, -2, "Note: If grid is geographic then we report area-weighted statistics.");
GMT_Usage (API, 1, "\n-M Search for the global min and max locations (x0,y0{,z0}) and (x1,y1{,z1}).");
GMT_Usage (API, 1, "\n-M[c|f]");
GMT_Usage (API, -2, "\nSearch for the global data min and max locations (x0,y0{,z0}) and (x1,y1{,z1}) [Default]."
" Append c to only determine data min/max range if missing from the header, or f to force that search to override the header range.");
GMT_Option (API, "R");
GMT_Usage (API, 1, "\n-T[<dv>][+a[<alpha>]][+s]");
GMT_Usage (API, -2, "Print global -Tvmin/vmax[/dv] (in rounded multiples of <dv>, if given). Optional modifiers:");
Expand Down Expand Up @@ -310,9 +318,17 @@ static int parse (struct GMT_CTRL *GMT, struct GRDINFO_CTRL *Ctrl, struct GMT_OP
Ctrl->L.norm |= (1+2+4); break;
}
break;
case 'M': /* Global extrema */
case 'M': /* Global extrema and|or update missing header data range */
n_errors += gmt_M_repeated_module_option (API, Ctrl->M.active);
n_errors += gmt_get_no_argument (GMT, opt->arg, opt->option, 0);
switch (opt->arg[0]) {
case 'c': Ctrl->M.mode = GRDINFO_CONDITIONAL; break;
case 'f': Ctrl->M.mode = GRDINFO_FORCE; break;
case '\0': Ctrl->M.mode = GRDINFO_FORCE_REPORT; break;
default:
GMT_Report (API, GMT_MSG_ERROR, "Option -M: Unrecognized modifier %s\n", opt->arg);
n_errors++;
break;
}
break;
case 'Q': /* Expect cubes is no longer an option as we detect it automatically */
break;
Expand Down Expand Up @@ -623,7 +639,7 @@ EXTERN_MSC int GMT_grdinfo (void *V_API, int mode, void *args) {
n_cols = (is_cube) ? 8 : 6; /* w e s n {b t} w0 w1 */
if (!Ctrl->I.active) {
n_cols += (is_cube) ? 6 : 4; /* Add dx dy {dz} n_columns n_rows {n_layers} */
if (Ctrl->M.active) n_cols += (is_cube) ? 7 : 5; /* Add x0 y0 {z0} x1 y1 {z1} nnan */
if (Ctrl->M.mode == GRDINFO_FORCE_REPORT) n_cols += (is_cube) ? 7 : 5; /* Add x0 y0 {z0} x1 y1 {z1} nnan */
if (Ctrl->L.norm & 1) n_cols += 2; /* Add median scale */
if (Ctrl->L.norm & 2) n_cols += 3; /* Add mean stdev rms */
if (Ctrl->L.norm & 4) n_cols += 2; /* Add mode lmsscale */
Expand Down Expand Up @@ -705,7 +721,7 @@ EXTERN_MSC int GMT_grdinfo (void *V_API, int mode, void *args) {
GMT->current.io.col_type[GMT_OUT][n] = gmt_M_type (GMT, GMT_IN, n); /* Since grids may differ in types */

if (num_report && n_grds == 0) { /* Need to prepare col_type for output */
if (Ctrl->M.active) { /* Will write x,y[,z] min and max locations */
if (Ctrl->M.mode == GRDINFO_FORCE_REPORT) { /* Will write x,y[,z] min and max locations */
GMT->current.io.col_type[GMT_OUT][x_col_min] = GMT->current.io.col_type[GMT_OUT][x_col_max] = GMT->current.io.col_type[GMT_OUT][GMT_X];
GMT->current.io.col_type[GMT_OUT][y_col_min] = GMT->current.io.col_type[GMT_OUT][y_col_max] = GMT->current.io.col_type[GMT_OUT][GMT_Y];
if (is_cube) GMT->current.io.col_type[GMT_OUT][z_col_min] = GMT->current.io.col_type[GMT_OUT][z_col_max] = GMT->current.io.col_type[GMT_OUT][GMT_Z];
Expand All @@ -716,7 +732,16 @@ EXTERN_MSC int GMT_grdinfo (void *V_API, int mode, void *args) {
}
n_grds++;

if (Ctrl->E.active || Ctrl->M.active || Ctrl->L.active || subset || Ctrl->D.mode || (Ctrl->T.mode & 2)) { /* Need to read the data (all or subset) */
if (Ctrl->M.mode == GRDINFO_CONDITIONAL) { /* Conditionally read and find min/max if not in header */
if (doubleAlmostEqualZero (header->z_min, header->z_max)) { /* Force search */
GMT_Report (API, GMT_MSG_WARNING, "No actual range for data in header - must read matrix and determine range\n");
Ctrl->M.mode = GRDINFO_FORCE;
}
else /* else the header has the data range; no need to read file */
GMT_Report (API, GMT_MSG_DEBUG, "Found valid actual data range in header - no need to read matrix\n");
}

if (Ctrl->E.active || Ctrl->M.mode || Ctrl->L.active || subset || Ctrl->D.mode || (Ctrl->T.mode & 2)) { /* Need to read the data (all or subset) */
if (is_cube) {
if (GMT_Read_Data (API, GMT_IS_CUBE, GMT_IS_FILE, GMT_IS_VOLUME, GMT_DATA_ONLY, wesn, opt->arg, U) == NULL) {
Return (API->error);
Expand Down Expand Up @@ -779,7 +804,7 @@ EXTERN_MSC int GMT_grdinfo (void *V_API, int mode, void *args) {

if (Ctrl->T.mode & 2) strncpy (grdfile, opt->arg, PATH_MAX-1);

if (Ctrl->M.active || Ctrl->L.active) { /* Must determine the location of global min and max values */
if (Ctrl->M.mode == GRDINFO_FORCE_REPORT || Ctrl->M.mode == GRDINFO_FORCE || Ctrl->L.active) { /* Must determine the [location of] global min and max values */
uint64_t ij_min, ij_max, here = 0, node;
unsigned int col, row;
gmt_grdfloat *data = (is_cube) ? U->data : G->data;
Expand Down Expand Up @@ -816,6 +841,10 @@ EXTERN_MSC int GMT_grdinfo (void *V_API, int mode, void *args) {
}
else /* Not a single valid node */
x_min = x_max = y_min = y_max = z_min = z_max = GMT->session.d_NaN;
if (Ctrl->M.mode == GRDINFO_FORCE) { /* Update header */
header->z_min = v_min;
header->z_max = v_max;
}
}

if (Ctrl->L.norm && gmt_M_is_geographic (GMT, GMT_IN)) { /* Must use spherical weights */
Expand Down Expand Up @@ -913,7 +942,7 @@ EXTERN_MSC int GMT_grdinfo (void *V_API, int mode, void *args) {
out[col++] = header->z_min; out[col++] = header->z_max;
out[col++] = header->inc[GMT_X]; out[col++] = header->inc[GMT_Y]; if (is_cube) out[col++] = U->z_inc;
out[col++] = header->n_columns; out[col++] = header->n_rows; if (is_cube) out[col++] = header->n_bands;
if (Ctrl->M.active) {
if (Ctrl->M.mode == GRDINFO_FORCE_REPORT) {
out[col++] = x_min; out[col++] = y_min; if (is_cube) out[col++] = z_min;
out[col++] = x_max; out[col++] = y_max; if (is_cube) out[col++] = z_max;
}
Expand All @@ -923,7 +952,7 @@ EXTERN_MSC int GMT_grdinfo (void *V_API, int mode, void *args) {
if (Ctrl->L.norm & 2) {
out[col++] = z_mean; out[col++] = z_stdev; out[col++] = z_rms;
}
if (Ctrl->M.active) {
if (Ctrl->M.mode == GRDINFO_FORCE_REPORT) {
out[col++] = (double)n_nan;
}
if (Ctrl->L.norm & 4) {
Expand Down Expand Up @@ -975,7 +1004,7 @@ EXTERN_MSC int GMT_grdinfo (void *V_API, int mode, void *args) {
gmt_ascii_format_col (GMT, text, (double)header->n_bands, GMT_OUT, GMT_W);
strcat (record, text); strcat (record, sep);
}
if (Ctrl->M.active) {
if (Ctrl->M.mode == GRDINFO_FORCE_REPORT) {
strcat (record, sep); gmt_ascii_format_col (GMT, text, x_min, GMT_OUT, GMT_X); strcat (record, text);
strcat (record, sep); gmt_ascii_format_col (GMT, text, y_min, GMT_OUT, GMT_Y); strcat (record, text);
if (is_cube) { strcat (record, sep); gmt_ascii_format_col (GMT, text, z_min, GMT_OUT, GMT_Z); strcat (record, text); }
Expand All @@ -992,7 +1021,7 @@ EXTERN_MSC int GMT_grdinfo (void *V_API, int mode, void *args) {
strcat (record, sep); gmt_ascii_format_col (GMT, text, z_stdev, GMT_OUT, GMT_Z); strcat (record, text);
strcat (record, sep); gmt_ascii_format_col (GMT, text, z_rms, GMT_OUT, GMT_Z); strcat (record, text);
}
if (Ctrl->M.active) { sprintf (text, "%s%" PRIu64, sep, n_nan); strcat (record, text); }
if (Ctrl->M.mode == GRDINFO_FORCE_REPORT) { sprintf (text, "%s%" PRIu64, sep, n_nan); strcat (record, text); }
if (Ctrl->L.norm & 4) {
strcat (record, sep); gmt_ascii_format_col (GMT, text, z_mode, GMT_OUT, GMT_Z); strcat (record, text);
strcat (record, sep); gmt_ascii_format_col (GMT, text, z_lmsscl, GMT_OUT, GMT_Z); strcat (record, text);
Expand Down Expand Up @@ -1117,7 +1146,7 @@ EXTERN_MSC int GMT_grdinfo (void *V_API, int mode, void *args) {
}
}

if (Ctrl->M.active) {
if (Ctrl->M.mode == GRDINFO_FORCE_REPORT) {
if (v_min == DBL_MAX) v_min = GMT->session.d_NaN;
if (v_max == -DBL_MAX) v_max = GMT->session.d_NaN;
sprintf (record, "%s: v_min: ", HH->name);
Expand Down Expand Up @@ -1159,7 +1188,7 @@ EXTERN_MSC int GMT_grdinfo (void *V_API, int mode, void *args) {
(header->z_max - header->z_add_offset) / header->z_scale_factor);
}
GMT_Put_Record (API, GMT_WRITE_DATA, Out);
if (n_nan || Ctrl->M.active) {
if (n_nan || Ctrl->M.mode == GRDINFO_FORCE_REPORT) {
double percent = 100.0 * n_nan / header->nm;
sprintf (record, "%s: %" PRIu64 " nodes (%.1f%%) set to NaN", HH->name, n_nan, percent);
GMT_Put_Record (API, GMT_WRITE_DATA, Out);
Expand Down