xref: /plan9/sys/src/cmd/map/libmap/lambert.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 static struct coord stdp0, stdp1;
6219b2ee8SDavid du Colombier static double k;
73e12c5d1SDavid du Colombier 
83e12c5d1SDavid du Colombier static int
Xlambert(struct place * place,double * x,double * y)9219b2ee8SDavid du Colombier Xlambert(struct place *place, double *x, double *y)
103e12c5d1SDavid du Colombier {
11219b2ee8SDavid du Colombier 	double r;
123e12c5d1SDavid du Colombier 	if(place->nlat.l < -80.*RAD)
133e12c5d1SDavid du Colombier 		return(-1);
143e12c5d1SDavid du Colombier 	if(place->nlat.l > 89.*RAD)
153e12c5d1SDavid du Colombier 		r = 0;	/* slovenly */
163e12c5d1SDavid du Colombier 	else
173e12c5d1SDavid du Colombier 		r = stdp0.c*exp(0.5*k*log(
183e12c5d1SDavid du Colombier 		   (1+stdp0.s)*(1-place->nlat.s)/((1-stdp0.s)*(1+place->nlat.s))));
193e12c5d1SDavid du Colombier 	if(stdp1.l<0.)
203e12c5d1SDavid du Colombier 		r = -r;
213e12c5d1SDavid du Colombier 	*x = - r*sin(k * place->wlon.l);
223e12c5d1SDavid du Colombier 	*y = - r*cos(k * place->wlon.l);
233e12c5d1SDavid du Colombier 	return(1);
243e12c5d1SDavid du Colombier }
253e12c5d1SDavid du Colombier 
263e12c5d1SDavid du Colombier proj
lambert(double par0,double par1)27219b2ee8SDavid du Colombier lambert(double par0, double par1)
283e12c5d1SDavid du Colombier {
29219b2ee8SDavid du Colombier 	double temp;
303e12c5d1SDavid du Colombier 	if(fabs(par0)>fabs(par1)){
313e12c5d1SDavid du Colombier 		temp = par0;
323e12c5d1SDavid du Colombier 		par0 = par1;
333e12c5d1SDavid du Colombier 		par1 = temp;
343e12c5d1SDavid du Colombier 	}
353e12c5d1SDavid du Colombier 	deg2rad(par0, &stdp0);
363e12c5d1SDavid du Colombier 	deg2rad(par1, &stdp1);
373e12c5d1SDavid du Colombier 	if(fabs(par1+par0)<.1)
383e12c5d1SDavid du Colombier 		return(mercator());
393e12c5d1SDavid du Colombier 	if(fabs(par1-par0)<.1)
403e12c5d1SDavid du Colombier 		return(perspective(-1.));
413e12c5d1SDavid du Colombier 	if(fabs(par0)>89.5||fabs(par1)>89.5)
423e12c5d1SDavid du Colombier 		return(0);
433e12c5d1SDavid du Colombier 	k = 2*log(stdp1.c/stdp0.c)/log(
443e12c5d1SDavid du Colombier 		(1+stdp0.s)*(1-stdp1.s)/((1-stdp0.s)*(1+stdp1.s)));
453e12c5d1SDavid du Colombier 	return(Xlambert);
463e12c5d1SDavid du Colombier }
47