/*--------------------------------------------------------------------
 *	$Id: psclip.c,v 1.19 2004/09/10 23:16:15 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
 *--------------------------------------------------------------------*/
/*
 * psclip reads one or many xy-files and draws polygons just like psxy
 * with the exception that these polygons are set up to act as clip
 * paths for subsequent plotting.  psclip uses the even-odd winding
 * rule to decide what's inside and outside the path.
 *
 * Author:	Paul Wessel
 * Date:	6-FEB-1991-2000
 * Version:	4
 *
 */

#include "gmt.h"

main (int argc, char **argv)
{
	int i, n, fno, n_alloc, n_args, n_files = 0, n_fields, n_expected_fields = 0;
	
	BOOLEAN error = FALSE, nofile = TRUE, end_of_clip = FALSE;
	BOOLEAN invert = FALSE, done, first, set_z = FALSE;
	
	char line[BUFSIZ];
	
	double west, east, north, south, x0, y0, new_z_level = 0.0, *in;
	double *xx, *yy;	/* Arrays holding the contour xy values */

	FILE *fp = NULL;
	
	west = east = north = south = 0.0;

	argc = GMT_begin (argc, argv);

	for (i = 1; i < argc; i++) {
		if (argv[i][0] == '-') {
			switch (argv[i][1]) {
		
				/* Common parameters */
			
				case 'B':
				case 'H':
				case 'J':
				case 'K':
				case 'O':
				case 'P':
				case 'R':
				case 'U':
				case 'V':
				case 'X':
				case 'x':
				case 'Y':
				case 'y':
				case 'b':
				case 'c':
				case 'f':
				case ':':
				case '\0':
					error += GMT_get_common_args (argv[i], &west, &east, &south, &north);
					break;
				
				/* Supplemental parameters */
			
				case 'S':	/* Undocumented backwards compatible option */
				case 'C':
					end_of_clip = TRUE;
					break;
				case 'E':
					sscanf (&argv[i][2], "%lf/%lf",&z_project.view_azimuth, &z_project.view_elevation);
					break;
				case 'N':
					invert = TRUE;
					break;
				case 'M':               /* Multiple line segments */
					GMT_multisegment (&argv[i][2]);
					break;
				case 'Z':
					if (argv[i][2]) {
						new_z_level = atof (&argv[i][2]);
						set_z = TRUE;
					}
					break;
				default:
					error = TRUE;
					GMT_default_error (argv[i][1]);
					break;
			}
		}
		else
			n_files++;
	}
	
	if (argc == 1 || GMT_quick) {
		fprintf (stderr,"psclip %s - To set up polygonal clip paths\n\n", GMT_VERSION);
		fprintf (stderr, "usage: psclip -C [-K] [-O]  OR\n");
		fprintf (stderr, "\tpsclip <xy-files> -J<params> -R<west/east/south/north> [-B<tickinfo>]\n");
		fprintf (stderr, "\t[-Eaz/el] [-H[<nrec>]] [-K] [-M<flag>] [-N] [-O] [-P] [-U[label] [-V]\n");
		fprintf (stderr, "\t[-X<x_shift>] [-Y<y_shift>] [-c<ncopies>] [-:] [-bi[s][<n>]] [-f[i|o]<colinfo>]\n\n");
		
		if (GMT_quick) exit (EXIT_FAILURE);
		
		fprintf (stderr, "\t-C means stop existing clip-path.  -R, -J are not required\n");
		fprintf (stderr, "\t<xy-files> is one or more polygon files.  If none, standard input is read\n");
		GMT_explain_option ('j');
		GMT_explain_option ('R');
		fprintf (stderr, "\n\tOPTIONS:\n");
		GMT_explain_option ('b');
		fprintf (stderr, "\t-E set azimuth and elevation of viewpoint for 3-D perspective [180/90]\n");
		GMT_explain_option ('H');
		GMT_explain_option ('K');
		GMT_explain_option ('M');
		fprintf (stderr, "\t-N uses the outside of the polygons as clip paths\n");
		GMT_explain_option ('O');
		GMT_explain_option ('P');
		GMT_explain_option ('U');
		GMT_explain_option ('V');
		GMT_explain_option ('X');
		GMT_explain_option ('c');
		GMT_explain_option (':');
		GMT_explain_option ('i');
		GMT_explain_option ('n');
		fprintf (stderr, "\t   Default is 2 input columns\n");
		GMT_explain_option ('f');
		GMT_explain_option ('.');
		exit (EXIT_FAILURE);
	}
	
	if (!end_of_clip) {
		if (!project_info.region_supplied) {
			fprintf (stderr, "%s: GMT SYNTAX ERROR:  Must specify -R option\n", GMT_program);
			error++;
		}
		if (z_project.view_azimuth > 360.0 || z_project.view_elevation <= 0.0 || z_project.view_elevation > 90.0) {
			fprintf (stderr, "%s: GMT SYNTAX ERROR -E option:  Enter azimuth in 0-360 range, elevation in 0-90 range\n", GMT_program);
			error++;
		}
	}
	
	if (GMT_io.binary[GMT_IN] && gmtdefs.io_header[GMT_IN]) {
		fprintf (stderr, "%s: GMT SYNTAX ERROR.  Binary input data cannot have header -H\n", GMT_program);
		error++;
	}
        if (GMT_io.binary[GMT_IN] && GMT_io.ncol[GMT_IN] == 0) GMT_io.ncol[GMT_IN] = 2;
        if (GMT_io.binary[GMT_IN] && GMT_io.ncol[GMT_IN] < 2) {
                fprintf (stderr, "%s: GMT SYNTAX ERROR.  Binary input data (-bi) must have at least 2 columns\n", GMT_program);
		error++;
	}

	if (error) exit (EXIT_FAILURE);

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

	if (GMT_io.binary[GMT_IN] && gmtdefs.verbose) {
		char *type[2] = {"double", "single"};
		fprintf (stderr, "%s: Expects %d-column %s-precision binary data\n", GMT_program, GMT_io.ncol[GMT_IN], type[GMT_io.single_precision[GMT_IN]]);
	}

	if (!project_info.x_off_supplied && gmtdefs.overlay) gmtdefs.x_origin = 0.0;
	if (!project_info.y_off_supplied && gmtdefs.overlay) gmtdefs.y_origin = 0.0;

	if (end_of_clip)
		gmtdefs.page_orientation |= GMT_CLIP_OFF;	/* Signal that this program terminates clipping that initiated prior to this process */
	else
		gmtdefs.page_orientation |= GMT_CLIP_ON;	/* Signal that this program initiates clipping that wil outlive this process */
		
	ps_plotinit (CNULL, gmtdefs.overlay, gmtdefs.page_orientation, gmtdefs.x_origin, gmtdefs.y_origin,
		gmtdefs.global_x_scale, gmtdefs.global_y_scale, gmtdefs.n_copies,
		gmtdefs.dpi, GMT_INCH, gmtdefs.paper_width, gmtdefs.page_rgb, gmtdefs.encoding.name, GMT_epsinfo (argv[0]));
	GMT_echo_command (argc, argv);
	if (gmtdefs.unix_time) GMT_timestamp (argc, argv);

	if (end_of_clip && !frame_info.plot) {	/* Just undo previous clip-path, no basemap needed */
		ps_clipoff ();
		ps_plotend (gmtdefs.last_page);
	
		if (gmtdefs.verbose) fprintf (stderr, "%s: Done!\n", GMT_program);
	
		GMT_end (argc, argv);
	}

	GMT_map_setup (west, east, south, north);
		
	if (project_info.three_D) ps_transrotate (-z_project.xmin, -z_project.ymin, 0.0);

	if (set_z) project_info.z_level = new_z_level;

	n_expected_fields = (GMT_io.ncol[GMT_IN]) ? GMT_io.ncol[GMT_IN] : 2;

	if (!end_of_clip) {	/* Start new clip_path */
	
		GMT_map_basemap ();
		
		xx = (double *) GMT_memory (VNULL, (size_t)GMT_CHUNK, sizeof (double), GMT_program);
		yy = (double *) GMT_memory (VNULL, (size_t)GMT_CHUNK, sizeof (double), GMT_program);
		n_alloc = GMT_CHUNK;
		
		if (invert) GMT_map_clip_on (GMT_no_rgb, 1);	/* Must clip map */
		first = (invert) ? FALSE : TRUE;
		done = FALSE;
		n_args = (argc > 1) ? argc : 2;
		if (n_files > 0) nofile = FALSE;
		
		for (fno = 1; !done && fno < n_args; fno++) {	/* Loop over all input files */
			if (!nofile && argv[fno][0] == '-') continue;
			if (nofile) {
				fp = GMT_stdin;
				done = TRUE;
#ifdef SET_IO_MODE
				GMT_setmode (GMT_IN);
#endif
			}
			else if ((fp = GMT_fopen (argv[fno], GMT_io.r_mode)) == NULL) {
				fprintf (stderr, "%s: Cannot open file %s\n", GMT_program, argv[fno]);
				continue;
			}

			if (!nofile && gmtdefs.verbose) {
				fprintf (stderr, "%s: Working on file %s\n", GMT_program, argv[fno]);
				sprintf (line, "File: %s", argv[fno]);
				ps_comment (line);
			}
			if (gmtdefs.io_header[GMT_IN]) for (i = 0; i < gmtdefs.n_header_recs; i++) GMT_fgets (line, BUFSIZ, fp);
			if (GMT_io.multi_segments) {
				GMT_fgets (line, BUFSIZ, fp);
				if (gmtdefs.verbose) ps_comment (line);
			}
			n_fields = GMT_input (fp, &n_expected_fields, &in);

			while (! (GMT_io.status & GMT_IO_EOF)) {	/* Not yet EOF */

				while (GMT_io.status & GMT_IO_SEGMENT_HEADER) {
					if (gmtdefs.verbose) ps_comment (GMT_io.segment_header);
					n_fields = GMT_input (fp, &n_expected_fields, &in);
				}
				if ((GMT_io.status & GMT_IO_EOF)) continue;	/* At EOF */

				/* Start a new segment path */

				n = 0;

				while (! (GMT_io.status & (GMT_IO_SEGMENT_HEADER | GMT_IO_EOF))) {	/* Keep going until FALSE */
					if (GMT_io.status & GMT_IO_MISMATCH) {
						fprintf (stderr, "%s: Mismatch between actual (%d) and expected (%d) fields near line %d\n", GMT_program, n_fields, 2, n);
						exit (EXIT_FAILURE);
					}

					xx[n] = in[0];	yy[n] = in[1];
					n++;
					if (n == n_alloc) {
						n_alloc += GMT_CHUNK;
						xx = (double *) GMT_memory ((void *)xx, (size_t)n_alloc, sizeof (double), GMT_program);
						yy = (double *) GMT_memory ((void *)yy, (size_t)n_alloc, sizeof (double), GMT_program);
					}
					n_fields = GMT_input (fp, &n_expected_fields, &in);
				}
		
				xx = (double *) GMT_memory ((void *)xx, (size_t)n, sizeof (double), GMT_program);
				yy = (double *) GMT_memory ((void *)yy, (size_t)n, sizeof (double), GMT_program);
				n_alloc = n;
				
				for (i = 0; i < n; i++) {
					GMT_geo_to_xy (xx[i], yy[i], &x0, &y0);
					xx[i] = x0; yy[i] = y0;
				}
				
				if (project_info.three_D) GMT_2D_to_3D (xx, yy, project_info.z_level, n);
				
				ps_clipon (xx, yy, n, GMT_no_rgb, first);
				first = FALSE;
			}
			if (fp != GMT_stdin) GMT_fclose (fp);
		}

		ps_clipon (xx, yy, 0, GMT_no_rgb, 2);	/* Finish up the composite clip path */

		GMT_free ((void *)xx);
		GMT_free ((void *)yy);
	}

	else {	/* Undo previous clip-path and draw basemap */
		ps_clipoff ();
		GMT_map_basemap ();
	}
	
	if (project_info.three_D) ps_rotatetrans (z_project.xmin, z_project.ymin, 0.0);

	ps_plotend (gmtdefs.last_page);
	
	if (gmtdefs.verbose) fprintf (stderr, "%s: Done!\n", GMT_program);
	
	GMT_end (argc, argv);
}
