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

   1 #include <u.h>
   2 #include <libc.h>
   3 #include <stdio.h>
   4 #include "map.h"
   5
   6 static double cirmod(double);
   7
   8 static struct place pole;       /* map pole is tilted to here */
   9 static struct coord twist;      /* then twisted this much */
  10 static struct place ipole;      /* inverse transfrom */
  11 static struct coord itwist;
  12
  13 void
  14 orient(double lat, double lon, double theta)
  15 {
  16         lat = cirmod(lat);
  17         if(lat>90.) {
  18                 lat = 180. - lat;
  19                 lon -= 180.;
  20                 theta -= 180.;
  21         } else if(lat < -90.) {
  22                 lat = -180. - lat;
  23                 lon -= 180.;
  24                 theta -= 180;
  25         }
  26         latlon(lat,lon,&pole);
  27         deg2rad(theta, &twist);
  28         latlon(lat,180.-theta,&ipole);
  29         deg2rad(180.-lon, &itwist);
  30 }
  31
  32 void
  33 latlon(double lat, double lon, struct place *p)
  34 {
  35         lat = cirmod(lat);
  36         if(lat>90.) {
  37                 lat = 180. - lat;
  38                 lon -= 180.;
  39         } else if(lat < -90.) {
  40                 lat = -180. - lat;
  41                 lon -= 180.;
  42         }
  43         deg2rad(lat,&p->nlat);
  44         deg2rad(lon,&p->wlon);
  45 }
  46
  47 void
  48 deg2rad(double theta, struct coord *coord)
  49 {
  50         theta = cirmod(theta);
  51         coord->l = theta*RAD;
  52         if(theta==90) {
  53                 coord->s = 1;
  54                 coord->c = 0;
  55         } else if(theta== -90) {
  56                 coord->s = -1;
  57                 coord->c = 0;
  58         } else
  59                 sincos(coord);
  60 }
  61
  62 static double
  63 cirmod(double theta)
  64 {
  65         while(theta >= 180.)
  66                 theta -= 360;
  67         while(theta<-180.)
  68                 theta += 360.;
  69         return(theta);
  70 }
  71
  72 void
  73 sincos(struct coord *coord)
  74 {
  75         coord->s = sin(coord->l);
  76         coord->c = cos(coord->l);
  77 }
  78
  79 void
  80 normalize(struct place *gg)
  81 {
  82         norm(gg,&pole,&twist);
  83 }
  84
  85 void
  86 invert(struct place *g)
  87 {
  88         norm(g,&ipole,&itwist);
  89 }
  90
  91 void
  92 norm(struct place *gg, struct place *pp, struct coord *tw)
  93 {
  94         register struct place *g;       /*geographic coords */
  95         register struct place *p;       /* new pole in old coords*/
  96         struct place m;                 /* standard map coords*/
  97         g = gg;
  98         p = pp;
  99         if(p->nlat.s == 1.) {
 100                 if(p->wlon.l+tw->l == 0.)
 101                         return;
 102                 g->wlon.l -= p->wlon.l+tw->l;
 103         } else {
 104                 if(p->wlon.l != 0) {
 105                         g->wlon.l -= p->wlon.l;
 106                         sincos(&g->wlon);
 107                 }
 108                 m.nlat.s = p->nlat.s * g->nlat.s
 109                         + p->nlat.c * g->nlat.c * g->wlon.c;
 110                 m.nlat.c = sqrt(1. - m.nlat.s * m.nlat.s);
 111                 m.nlat.l = atan2(m.nlat.s, m.nlat.c);
 112                 m.wlon.s = g->nlat.c * g->wlon.s;
 113                 m.wlon.c = p->nlat.c * g->nlat.s
 114                         - p->nlat.s * g->nlat.c * g->wlon.c;
 115                 m.wlon.l = atan2(m.wlon.s, - m.wlon.c)
 116                         - tw->l;
 117                 *g = m;
 118         }
 119         sincos(&g->wlon);
 120         if(g->wlon.l>PI)
 121                 g->wlon.l -= 2*PI;
 122         else if(g->wlon.l<-PI)
 123                 g->wlon.l += 2*PI;
 124 }
 125
 126 void
 127 printp(struct place *g)
 128 {
 129 printf("%.3f %.3f %.3f %.3f %.3f %.3f\n",
 130 g->nlat.l,g->nlat.s,g->nlat.c,g->wlon.l,g->wlon.s,g->wlon.c);
 131 }
 132
 133 void
 134 copyplace(struct place *g1, struct place *g2)
 135 {
 136         *g2 = *g1;
 137 }