sys/src/cmd/map/libmap/lambert.c

   1 #include <u.h>
   2 #include <libc.h>
   3 #include "map.h"
   4
   5 static struct coord stdp0, stdp1;
   6 static double k;
   7
   8 static int
   9 Xlambert(struct place *place, double *x, double *y)
  10 {
  11         double r;
  12         if(place->nlat.l < -80.*RAD)
  13                 return(-1);
  14         if(place->nlat.l > 89.*RAD)
  15                 r = 0;  /* slovenly */
  16         else
  17                 r = stdp0.c*exp(0.5*k*log(
  18                    (1+stdp0.s)*(1-place->nlat.s)/((1-stdp0.s)*(1+place->nlat.s))));
  19         if(stdp1.l<0.)
  20                 r = -r;
  21         *x = - r*sin(k * place->wlon.l);
  22         *y = - r*cos(k * place->wlon.l);
  23         return(1);
  24 }
  25
  26 proj
  27 lambert(double par0, double par1)
  28 {
  29         double temp;
  30         if(fabs(par0)>fabs(par1)){
  31                 temp = par0;
  32                 par0 = par1;
  33                 par1 = temp;
  34         }
  35         deg2rad(par0, &stdp0);
  36         deg2rad(par1, &stdp1);
  37         if(fabs(par1+par0)<.1)
  38                 return(mercator());
  39         if(fabs(par1-par0)<.1)
  40                 return(perspective(-1.));
  41         if(fabs(par0)>89.5||fabs(par1)>89.5)
  42                 return(0);
  43         k = 2*log(stdp1.c/stdp0.c)/log(
  44                 (1+stdp0.s)*(1-stdp1.s)/((1-stdp0.s)*(1+stdp1.s)));
  45         return(Xlambert);
  46 }