diff --git a/doc/rst/source/grdinfo.rst b/doc/rst/source/grdinfo.rst index 642349f9e5c..ec8ef2db727 100644 --- a/doc/rst/source/grdinfo.rst +++ b/doc/rst/source/grdinfo.rst @@ -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| ] @@ -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), @@ -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 diff --git a/src/grdinfo.c b/src/grdinfo.c index 183478962f4..d9184a4049b 100644 --- a/src/grdinfo.c +++ b/src/grdinfo.c @@ -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; @@ -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; @@ -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[[/]][+i]] [-E[x|y][+l|L|u|U]] [-F] [-G] [-I[[/]|b|i|r]] [-L[a|0|1|2|p]] " - "[-M] [%s] [-T[][+a[]][+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[][+a[]][+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); @@ -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[][+a[]][+s]"); GMT_Usage (API, -2, "Print global -Tvmin/vmax[/dv] (in rounded multiples of , if given). Optional modifiers:"); @@ -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; @@ -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 */ @@ -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]; @@ -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); @@ -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; @@ -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 */ @@ -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; } @@ -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) { @@ -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); } @@ -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); @@ -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); @@ -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);