xref: /plan9/sys/src/cmd/map/libmap/elliptic.c (revision 59cc4ca53493a3c6d2349fe2b7f7c40f7dce7294)
1*59cc4ca5SDavid du Colombier #include <u.h>
2*59cc4ca5SDavid du Colombier #include <libc.h>
33e12c5d1SDavid du Colombier #include "map.h"
43e12c5d1SDavid du Colombier 
53e12c5d1SDavid du Colombier struct coord center;
63e12c5d1SDavid du Colombier 
73e12c5d1SDavid du Colombier static int
Xelliptic(struct place * place,double * x,double * y)8219b2ee8SDavid du Colombier Xelliptic(struct place *place, double *x, double *y)
93e12c5d1SDavid du Colombier {
10219b2ee8SDavid du Colombier 	double r1,r2;
113e12c5d1SDavid du Colombier 	r1 = acos(place->nlat.c*(place->wlon.c*center.c
123e12c5d1SDavid du Colombier 		- place->wlon.s*center.s));
133e12c5d1SDavid du Colombier 	r2 = acos(place->nlat.c*(place->wlon.c*center.c
143e12c5d1SDavid du Colombier 		+ place->wlon.s*center.s));
153e12c5d1SDavid du Colombier 	*x = -(r1*r1 - r2*r2)/(4*center.l);
163e12c5d1SDavid du Colombier 	*y = (r1*r1+r2*r2)/2 - (center.l*center.l+*x**x);
173e12c5d1SDavid du Colombier 	if(*y < 0)
183e12c5d1SDavid du Colombier 		*y = 0;
193e12c5d1SDavid du Colombier 	*y = sqrt(*y);
203e12c5d1SDavid du Colombier 	if(place->nlat.l<0)
213e12c5d1SDavid du Colombier 		*y = -*y;
223e12c5d1SDavid du Colombier 	return(1);
233e12c5d1SDavid du Colombier }
243e12c5d1SDavid du Colombier 
253e12c5d1SDavid du Colombier proj
elliptic(double l)26219b2ee8SDavid du Colombier elliptic(double l)
273e12c5d1SDavid du Colombier {
283e12c5d1SDavid du Colombier 	l = fabs(l);
293e12c5d1SDavid du Colombier 	if(l>89)
303e12c5d1SDavid du Colombier 		return(0);
313e12c5d1SDavid du Colombier 	if(l<1)
323e12c5d1SDavid du Colombier 		return(Xazequidistant);
333e12c5d1SDavid du Colombier 	deg2rad(l,&center);
343e12c5d1SDavid du Colombier 	return(Xelliptic);
353e12c5d1SDavid du Colombier }
36