/*
 *  geoloc.c
 *
 *	This file contains code that manages values of type geoloc.
 *	Values of type geoloc are georeferenced coordinates in one
 *	of several different projections.  We use a table to store
 *	function entries for the various types, and execute those
 *	functions in order to do work.
 */

#include <stdio.h>
#include <string.h>

#define DBL_DIG	52
#define FLT_DIG	25

#include "tmp/c.h"
#include "tmp/postgres.h"

#include "utils/rel.h"
#include "utils/log.h"

#include "access/htup.h"
#include "access/skey.h"
#include "access/heapam.h"

#include "catalog/syscache.h"
#include "catalog/pg_proc.h"
#include "catalog/pg_type.h"

#include "geoloc.h"
#include "geounits.h"

/* routines declared here */
extern ObjectId		findfunc ARGS((char *s, int nargs));
extern ObjectId		get_units_func ARGS((char *projname, int attnum));
extern HeapTuple	findtuple ARGS(( Relation reln, int attnum, ObjectId func, Datum value));

/*
 *  geoloc_in() -- Read a geoloc string and return a pointer to the internal
 *		   representation.
 *
 *	The input representation is
 *
 *		"<projection>,<precision>,<x>,<y>,<state>"
 *
 *	and we use the projection name to look up the appropriate GEO_UNITS
 *	tuple.
 *
 *	As a side effect, this routine folds case on the projection name
 *	so that all projections have lower-case names.
 */

GEOLOC *
geoloc_in(s)
	char *s;
{
	char *comma;
	ObjectId fid;

	if (s == (char *) NULL)
		elog(WARN, "Cannot interpret null geoloc string");

	/* extract the projection name */
	if ((comma = index(s, ',')) == (char *) NULL)
		elog(WARN, "'%s' is not a valid geoloc string", s);
	*comma = '\0';

	/* find the in function for this procection */
	fid = get_units_func(s, Anum_geo_units_inproc);

	/* restore the string to its original state */
	*comma++ = ',';

	/* call the input function on the string, excluding the name */
	return ((GEOLOC *) fmgr(fid, comma));
}

/*
 *  geoloc_out(), geoloc_fullout() -- Turn a GEOLOC value into a printable
 *				      string.
 *
 *	To get a simplified version of the original input string, use
 *	geoloc_out().  To get the full input string, use geoloc_fullout().
 *
 *	As with geoloc_in, we use the projection name (stored in the
 *	geoloc struct) to search the GEO_UNITS table for a matching
 *	tuple, and call the desired function there.
 */

char *
geoloc_out(v)
	GEOLOC *v;
{
	ObjectId fid;

	/* find the in function for this procection */
	fid = get_units_func(v->geo_name, Anum_geo_units_outproc);

	/* call the output function on the supplied value */
	return ((char *) fmgr(fid, v));
}

text *
geoloc_fullout(v)
	GEOLOC *v;
{
	ObjectId fid;

	/* find the in function for this procection */
	fid = get_units_func(v->geo_name, Anum_geo_units_fulloutproc);

	/* call the output function on the supplied value */
	return ((text *) fmgr(fid, v));
}

/*
 *  geoloc_tostd(), geoloc_fromstd() -- Convert a geoloc value to or from
 *					"standard" geoloc units.
 *
 *	By convention, standard units are latlon.  When converting to
 *	standard units, if the supplied geoloc value is not currently
 *	a latlon value, we look up the conversion procedure and call it.
 *
 *	When converting from standard units, we expect the user to supply
 *	the name of a projection.  We find that projection and call the
 *	fromstd proc stored in the associated tuple.  If necessary, we
 *	first convert to standard units.
 */

GEOLOC *
geoloc_tostd(v)
	GEOLOC *v;
{
	ObjectId fid;

	if (v == (GEOLOC *) NULL)
		return (v);

	/* try to avoid doing any work */
	if (strcmp(v->geo_name, "latlon") == 0)
		return (v);

	/* find the tostd function for this procection */
	fid = get_units_func(v->geo_name, Anum_geo_units_tostdproc);

	/* call the output function on the supplied value */
	return ((GEOLOC *) fmgr(fid, v));
}

GEOLOC *
geoloc_fromstd(v, newproj)
	GEOLOC *v;
	char *newproj;
{
	ObjectId fid;

	if (v == (GEOLOC *) NULL)
		return (v);

	/* try to avoid doing any work */
	if (strcmp(v->geo_name, newproj) == 0)
		return (v);

	/* if we're not yet in standard units, get there */
	if (strcmp(v->geo_name, "latlon") != 0)
		v = geoloc_tostd(v);

	/* find the tostd function for this procection */
	fid = get_units_func(newproj, Anum_geo_units_fromstdproc);

	/* call the output function on the supplied value */
	return ((GEOLOC *) fmgr(fid, v));
}

/*
 *  get_units_func -- Find the desired regproc id in the desired GEO_UNITS
 *		      tuple.
 *
 *	The caller supplies a projection name and the attribute number of
 *	the desired function, and we return the regproc id.
 */

ObjectId
get_units_func(projname, attnum)
	char *projname;
	int attnum;
{
	Relation unitsreln;
	HeapTuple units_tup;
	Form_geo_units units;
	ObjectId char8eq;
	ObjectId fid;

	/* find the function id for char8eq */
	char8eq = findfunc("char8eq", 2);

	/* find the geo_units tuple of interest */
	if ((unitsreln = heap_openr(Name_geo_units)) == (Relation) NULL)
		elog(WARN, "cannot open %s to find %s",
			    Name_geo_units, projname);

	units_tup = findtuple(unitsreln, Anum_geo_units_name, char8eq,
			      PointerGetDatum(projname));

	if (!HeapTupleIsValid(units_tup))
		elog(WARN, "Cannot find tuple for projection %s in %s",
			   projname, Name_geo_units);

	/* skip the tuple header and reach right into the tuple struct */
	units = (Form_geo_units) GETSTRUCT(units_tup);

	/* extract the desired fid */
	switch (attnum) {
		case Anum_geo_units_inproc:
			fid = (ObjectId) units->gu_inproc;
			break;

		case Anum_geo_units_outproc:
			fid = (ObjectId) units->gu_outproc;
			break;

		case Anum_geo_units_fulloutproc:
			fid = (ObjectId) units->gu_fulloutproc;
			break;

		case Anum_geo_units_fromstdproc:
			fid = (ObjectId) units->gu_fromstdproc;
			break;

		case Anum_geo_units_tostdproc:
			fid = (ObjectId) units->gu_tostdproc;
			break;

		default:
			elog(WARN, "Attribute %d of %s is not a regproc",
				   attnum, Name_geo_units);
	}

	/* be tidy */
	heap_close(unitsreln);

	return (fid);
}


/*
 *  findfunc() -- Find a function by function name.
 *
 *	This routine looks up the named function in the pg_proc system
 *	cache and returns its function id, for use by the function manager.
 *	We could supply the actual argument types to the function, but
 *	that's a fair amount of trouble.  Instead, we allow the system
 *	to find the function without any type information.  This only
 *	works if there's exactly one such named function in the catalogs.
 */

ObjectId
findfunc(funcname, nargs)
	char *funcname;
	int nargs;
{
	int i;
	ObjectId fid;
	bool retset;
	ObjectId rettype;
	ObjectId *actualtypes;
	ObjectId argtypes[8];

	for (i = 0; i < 8; i++)
		argtypes[i] = UNKNOWNOID;

	if (!func_get_detail(funcname, nargs, argtypes, &fid, &rettype,
			     &retset, &actualtypes))
		elog(WARN, "cannot find function %s in %s", funcname,
			   Name_pg_proc);

	return (fid);
}

/*
 *  findtuple() -- A simplified scan interface.
 *
 *	This routine takes a relation, an attribute, a function, and
 *	a value, and searches the relation for the first tuple such that
 *	the func(attribute, value) returns true.  This just abstracts the
 *	scan management code to a single place.
 */

HeapTuple
findtuple(reln, attnum, func, value)
	Relation reln;
	int attnum;
	ObjectId func;
	Datum value;
{
	HeapScanDesc sdesc;
	HeapTuple htup;
	ScanKeyEntryData skey;

	/* set up a scan key */
	ScanKeyEntryInitialize(&skey, (bits16) 0x0, (AttributeNumber) attnum,
			       (RegProcedure) func, (Datum) value);

	/* find the first qualifying tuple */
	sdesc = heap_beginscan(reln, 0, NowTimeQual, 1, (ScanKey) &skey);
	htup = heap_getnext(sdesc, 0, (Buffer *) NULL);

	/* be tidy */
	heap_endscan(sdesc);
	heap_close(reln);

	/* that's it */
	return (htup);
}

/*
 *  latlon_in() -- Read a Latitude/Longitude spec.
 *
 *	By here, we know that the input string begins with "latlon,", so
 *	we skip over everything up to and including the first comma.
 */

GEOLOC *
latlon_in(s)
	char *s;
{
	int i;
	int n;
	GEOLOC *v;

	v = (GEOLOC *) palloc(sizeof(GEOLOC));

	v->geo_varlen = sizeof(GEOLOC);
	strcpy(&(v->geo_name[0]), "latlon");

	n = sscanf(s, "%f, %lg, %lg", &v->geo_prec, &v->geo_x, &v->geo_y);

	if (n != 3)
		elog(WARN, "bad latlon spec %s", s);

	return (v);
}

char *
latlon_out(v)
	GEOLOC *v;
{
	char *p;

	if (v == (GEOLOC *) NULL)
		return ((char *) NULL);

	if (strcmp(v->geo_name, "latlon") != 0)
		v = geoloc_tostd(v);

	/* shouldn't have to be this stupid */
	p = (char *) palloc(128);

	sprintf(p, "%lg, %lg latlon", v->geo_x, v->geo_y);

	return (p);
}

text *
latlon_fullout(v)
	GEOLOC *v;
{
	text *p;

	if (v == (GEOLOC *) NULL)
		return ((text *) NULL);

	if (strcmp(v->geo_name, "latlon") != 0)
		v = geoloc_tostd(v);

	/* shouldn't have to be this stupid */
	p = (text *) palloc(256);

	sprintf(VARDATA(p), "latlon, %g, %lg, %lg",
		   v->geo_prec, v->geo_x, v->geo_y);

	VARSIZE(p) = strlen(VARDATA(p)) + sizeof(VARSIZE(p));
	return (p);
}

/*
 *  lazea_in() -- Read a Lambert Azimuthal Equal-Area projection spec.
 *
 *	By here, we know that the input string begins with "lazea,", so
 *	we skip over everything up to and including the first comma.
 */

LAZEA *
lazea_in(s)
	char *s;
{
	int i;
	int n;
	int sphcode;
	LAZEA *v;

	v = (LAZEA *) palloc(sizeof(LAZEA));

	v->la_g.geo_varlen = sizeof(LAZEA);
	strcpy(&(v->la_g.geo_name[0]), "lazea");

	n = sscanf(s, "%g, %lg, %lg, %d, %lg, %lg, %lg, %lg",
			&v->la_g.geo_prec, &v->la_g.geo_x,
			&v->la_g.geo_y, &sphcode, &v->la_lon0,
			&v->la_lat0, &v->la_falseE, &v->la_falseN);

	if (n != 8)
		elog(WARN, "bad lazea spec %s", s);

	if (sphcode > 256 ||  sphcode < 0)
		elog(WARN, "spheroid code %d out of range in lazea spec '%s'",
			sphcode, s);
	v->la_sphcode = (uint8) sphcode;

	return (v);
}

char *
lazea_out(v)
	LAZEA *v;
{
	char *p;

	/* shouldn't have to be this stupid */
	p = (char *) palloc(128);

	if (v == (LAZEA *) NULL)
		return ((char *) NULL);

	if (strcmp(v->la_g.geo_name, "lazea") != 0)
		v = (LAZEA *) geoloc_fromstd(geoloc_tostd(v), "lazea");

	sprintf(p, "%lg, %lg lazea", v->la_g.geo_x, v->la_g.geo_y);

	return (p);
}

text *
lazea_fullout(v)
	LAZEA *v;
{
	text *p;

	if (v == (LAZEA *) NULL)
		return ((text *) NULL);

	if (strcmp(v->la_g.geo_name, "lazea") != 0)
		v = (LAZEA *) geoloc_fromstd(geoloc_tostd(v), "lazea");

	/* shouldn't have to be this stupid */
	p = (text *) palloc(256);

	sprintf(VARDATA(p),
		"lazea, %g, %lg, %lg, %d, %lg, %lg, %lg, %lg",
		v->la_g.geo_prec,
		v->la_g.geo_x, v->la_g.geo_y,
		v->la_sphcode,
		v->la_lat0, v->la_lon0, v->la_falseE, v->la_falseN);

	VARSIZE(p) = strlen(VARDATA(p)) + sizeof(VARSIZE(p));
	return (p);
}

/*
 *  utm_in() -- Read a UTM projection spec.
 *
 *	By here, we know that the input string begins with "utm,", so
 *	we skip over everything up to and including the first comma.
 */

UTM *
utm_in(s)
	char *s;
{
	int n;
	int i;
	int sphcode, zone;
	UTM *v;

	v = (UTM *) palloc(sizeof(UTM));

	v->utm_g.geo_varlen = sizeof(UTM);
	strcpy(&(v->utm_g.geo_name[0]), "utm");

	n = sscanf(s, "%f, %lg, %lg, %d, %d",
			&v->utm_g.geo_prec, &v->utm_g.geo_x, &v->utm_g.geo_y,
			&sphcode, &zone);

	if (n != 5)
		elog(WARN, "bad utm spec %s", s);

	if (sphcode > 256 ||  sphcode < 0)
		elog(WARN, "spheroid code %d out of range in utm spec '%s'",
			sphcode, s);
	if (zone > 256 ||  zone < 0)
		elog(WARN, "zone %d out of range in utm spec '%s'",
			sphcode, s);
	v->utm_sphcode = (uint8) sphcode;
	v->utm_zone = (uint8) zone;
	return (v);
}

char *
utm_out(v)
	UTM *v;
{
	char *p;

	if (v == (UTM *) NULL)
		return ((char *) NULL);

	if (strcmp(v->utm_g.geo_name, "utm") != 0)
		v = (UTM *) geoloc_fromstd(geoloc_tostd(v), "utm");

	/* shouldn't have to be this stupid */
	p = (char *) palloc(128);

	sprintf(p, "%lg, %lg utm", v->utm_g.geo_x, v->utm_g.geo_y);

	return (p);
}

text *
utm_fullout(v)
	UTM *v;
{
	text *p;

	if (v == (UTM *) NULL)
		return ((text *) NULL);

	if (strcmp(v->utm_g.geo_name, "utm") != 0)
		v = (UTM *) geoloc_fromstd(geoloc_tostd(v), "utm");

	/* shouldn't have to be this stupid */
	p = (text *) palloc(256);

	sprintf(VARDATA(p), "utm, %g, %lg, %lg, %d, %d",
		v->utm_g.geo_prec,
		v->utm_g.geo_x, v->utm_g.geo_y,
		v->utm_sphcode, v->utm_zone);

	VARSIZE(p) = strlen(VARDATA(p)) + sizeof(VARSIZE(p));
	return (p);
}
