/*--------------------------------------------------------------------
 *	$Id: grdclip.c,v 1.6 2004/01/02 22:45:13 pwessel Exp $
 *
 *	Copyright (c) 1991-2004 by P. Wessel and W. H. F. Smith
 *	See COPYING file for copying and redistribution conditions.
 *
 *	This program is free software; you can redistribute it and/or modify
 *	it under the terms of the GNU General Public License as published by
 *	the Free Software Foundation; version 2 of the License.
 *
 *	This program is distributed in the hope that it will be useful,
 *	but WITHOUT ANY WARRANTY; without even the implied warranty of
 *	MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *	GNU General Public License for more details.
 *
 *	Contact info: gmt.soest.hawaii.edu
 *--------------------------------------------------------------------*/
/*
   grdclip read a grdfile and sets all values < the user-supplied
   lower limit to the value <below>, and all values > the user-
   supplied upper limit to the value <above>.  above/below can
   be any number including NaN.

   WHF Smith, 22 Aug 88
   Modified for v2.0 4/26/91 P.Wessel
   Modified for v3.1 4/29/98 P.Wessel
   Modified for v3.3.5 06/29/2000 P.Wessel
   Version:	4
 
 */


#include "gmt.h"

float *data;

main (int argc, char **argv)
{

	BOOLEAN	error, set_above, set_below;

	int i, nxy, n_above, n_below, n;

	float	above, below, low, high;

	char *infile, *outfile, format[BUFSIZ], txt[50];

	struct GRD_HEADER header;

	argc = GMT_begin (argc, argv);

	error = set_above = set_below = FALSE;
	n_above = n_below = 0;
	low = high = 0.0;

	for (i = 1; i < argc; i++) {
		if (argv[i][0] == '-') {
			switch (argv[i][1]) {
				/* Common parameters */
			
				case 'V':
				case '\0':
					error += GMT_get_common_args (argv[i], 0, 0, 0, 0);
					break;
				
				/* Supplemental parameters */
			
				case 'S':
					if (argv[i][2] == 'a') {
						set_above = TRUE;
						n = sscanf (&argv[i][3], "%f/%s", &high, txt);
						if (n != 2) {
							fprintf (stderr, "%s: GMT SYNTAX ERROR -Sa option.  Correct syntax:\n", GMT_program);
							fprintf (stderr, "\t-Sa<high>/<above>, <above> may be set to NaN\n");
							error = TRUE;
						}
						else 
							above = (txt[0] == 'N' || txt[0] == 'n') ? GMT_f_NaN : (float)atof (txt);
					}
					else if (argv[i][2] == 'b') {
						set_below = TRUE;
						n = sscanf (&argv[i][3], "%f/%s", &low, txt);
						if (n != 2) {
							fprintf (stderr, "%s: GMT SYNTAX ERROR -Sb option.  Correct syntax:\n", GMT_program);
							fprintf (stderr, "\t-Sb<low>/<below>, <below> may be set to NaN\n");
							error = TRUE;
						}
						else
							below = (txt[0] == 'N' || txt[0] == 'n') ? GMT_f_NaN : (float)atof (txt);
					}
					else {
						fprintf (stderr, "%s: GMT SYNTAX ERROR -S option.  Correct syntax:\n", GMT_program);
						fprintf (stderr, "\t-Sa<high>/<above> or -Sb<low>/<below>\n");
						error = TRUE;
					}
					break;
				case 'G':
					outfile = &argv[i][2];
					break;

				default:
					error = TRUE;
					GMT_default_error (argv[i][1]);
					break;
			}
		}
		else
			infile = argv[i];
	}
	
	if (GMT_quick || argc == 1) {
		fprintf (stderr, "grdclip %s - Clipping of range in grdfiles\n\n", GMT_VERSION);
		fprintf (stderr, "usage: grdclip <input_grdfile> -G<output_grdfile>\n");
		fprintf (stderr, "\t[-Sa<high/above>] [-Sb<low/below>] [-V]\n");
		
		if (GMT_quick) exit (EXIT_FAILURE);
		
		fprintf (stderr, "\tYou must choose at least one -S option.\n");
		fprintf (stderr, "\t-G name of output grid.\n");
		fprintf (stderr, "\n\tOPTIONS:\n");
		fprintf (stderr, "\t-Sa will set all data > high to the <above> value.\n");
		fprintf (stderr, "\t-Sb will set all data < low to the <below> value.\n");
		fprintf (stderr, "\t    above/below can be any number including NaN\n");
		GMT_explain_option ('V');
		GMT_explain_option ('.');
		exit (EXIT_FAILURE);
	}

	if (!outfile) {
		fprintf (stderr, "%s: GMT SYNTAX ERROR -G option:.  Must specify output file\n", GMT_program);
		error++;
	}
	if (!infile) {
		fprintf (stderr, "%s: GMT SYNTAX ERROR:.  Must specify input file\n", GMT_program);
		error++;
	}
	
	if (error) exit (EXIT_FAILURE);

	GMT_put_history (argc, argv);	/* Update .gmtcommands4 */

		
	if (GMT_read_grd_info (infile, &header)) {
		fprintf (stderr, "%s: Error opening file %s\n", GMT_program, infile);
		exit (EXIT_FAILURE);
	}
	
	GMT_grd_init (&header, argc, argv, TRUE);
	
	nxy = header.nx * header.ny;
	data = (float *) GMT_memory (VNULL, (size_t)nxy, sizeof (float), GMT_program);
	
	if (GMT_read_grd (infile, &header, data, 0.0, 0.0, 0.0, 0.0, GMT_pad, FALSE)) {
		fprintf (stderr, "%s: Error reading file %s\n", GMT_program, infile);
		exit (EXIT_FAILURE);
	}
	
	for (i = 0; i < nxy; i++) {
		if (set_above && (!GMT_is_fnan (data[i])) && data[i] > high) {
			data[i] = above;
			n_above++;
		}
		else if (set_below && (!GMT_is_fnan (data[i])) && data[i] < low) {
			data[i] = below;
			n_below++;
		}
	}
	
	if (GMT_write_grd (outfile, &header, data, 0.0, 0.0, 0.0, 0.0, GMT_pad, FALSE)) {
		fprintf (stderr, "%s: Error writing file %s\n", GMT_program, infile);
		exit (EXIT_FAILURE);
	}
	
	GMT_free ((void *) data);

	if (gmtdefs.verbose) {
		sprintf (format, "%s set to %s\n", gmtdefs.d_format, gmtdefs.d_format);
		if (set_below) {
			fprintf (stderr, "%s: %d values < ", GMT_program, n_below);
			fprintf (stderr, format, (double)low, (double)below);
		}
		if (set_above) {
			fprintf (stderr, "%s: %d values > ", GMT_program, n_above);
			fprintf (stderr, format, (double)high, (double)above);
		}
	}
	
	GMT_end (argc, argv);
}
