Blame


1 28994509 2004-04-21 devnull #include <u.h>
2 28994509 2004-04-21 devnull #include <libc.h>
3 28994509 2004-04-21 devnull #include "map.h"
4 28994509 2004-04-21 devnull
5 28994509 2004-04-21 devnull static struct place gywhem, gyehem;
6 28994509 2004-04-21 devnull static struct coord gytwist;
7 28994509 2004-04-21 devnull static double gyconst, gykc, gyside;
8 28994509 2004-04-21 devnull
9 28994509 2004-04-21 devnull
10 28994509 2004-04-21 devnull static void
11 28994509 2004-04-21 devnull dosquare(double z1, double z2, double *x, double *y)
12 28994509 2004-04-21 devnull {
13 28994509 2004-04-21 devnull double w1,w2;
14 28994509 2004-04-21 devnull w1 = z1 -1;
15 28994509 2004-04-21 devnull if(fabs(w1*w1+z2*z2)>.000001) {
16 28994509 2004-04-21 devnull cdiv(z1+1,z2,w1,z2,&w1,&w2);
17 28994509 2004-04-21 devnull w1 *= gyconst;
18 28994509 2004-04-21 devnull w2 *= gyconst;
19 28994509 2004-04-21 devnull if(w1<0)
20 28994509 2004-04-21 devnull w1 = 0;
21 28994509 2004-04-21 devnull elco2(w1,w2,gykc,1.,1.,x,y);
22 28994509 2004-04-21 devnull } else {
23 28994509 2004-04-21 devnull *x = gyside;
24 28994509 2004-04-21 devnull *y = 0;
25 28994509 2004-04-21 devnull }
26 28994509 2004-04-21 devnull }
27 28994509 2004-04-21 devnull
28 28994509 2004-04-21 devnull int
29 28994509 2004-04-21 devnull Xguyou(struct place *place, double *x, double *y)
30 28994509 2004-04-21 devnull {
31 28994509 2004-04-21 devnull int ew; /*which hemisphere*/
32 28994509 2004-04-21 devnull double z1,z2;
33 28994509 2004-04-21 devnull struct place pl;
34 28994509 2004-04-21 devnull ew = place->wlon.l<0;
35 28994509 2004-04-21 devnull copyplace(place,&pl);
36 28994509 2004-04-21 devnull norm(&pl,ew?&gyehem:&gywhem,&gytwist);
37 28994509 2004-04-21 devnull Xstereographic(&pl,&z1,&z2);
38 28994509 2004-04-21 devnull dosquare(z1/2,z2/2,x,y);
39 28994509 2004-04-21 devnull if(!ew)
40 28994509 2004-04-21 devnull *x -= gyside;
41 28994509 2004-04-21 devnull return(1);
42 28994509 2004-04-21 devnull }
43 28994509 2004-04-21 devnull
44 28994509 2004-04-21 devnull proj
45 28994509 2004-04-21 devnull guyou(void)
46 28994509 2004-04-21 devnull {
47 28994509 2004-04-21 devnull double junk;
48 28994509 2004-04-21 devnull gykc = 1/(3+2*sqrt(2.));
49 28994509 2004-04-21 devnull gyconst = -(1+sqrt(2.));
50 28994509 2004-04-21 devnull elco2(-gyconst,0.,gykc,1.,1.,&gyside,&junk);
51 28994509 2004-04-21 devnull gyside *= 2;
52 28994509 2004-04-21 devnull latlon(0.,90.,&gywhem);
53 28994509 2004-04-21 devnull latlon(0.,-90.,&gyehem);
54 28994509 2004-04-21 devnull deg2rad(0.,&gytwist);
55 28994509 2004-04-21 devnull return(Xguyou);
56 28994509 2004-04-21 devnull }
57 28994509 2004-04-21 devnull
58 28994509 2004-04-21 devnull int
59 28994509 2004-04-21 devnull guycut(struct place *g, struct place *og, double *cutlon)
60 28994509 2004-04-21 devnull {
61 28994509 2004-04-21 devnull int c;
62 28994509 2004-04-21 devnull c = picut(g,og,cutlon);
63 28994509 2004-04-21 devnull if(c!=1)
64 28994509 2004-04-21 devnull return(c);
65 28994509 2004-04-21 devnull *cutlon = 0.;
66 28994509 2004-04-21 devnull if(g->nlat.c<.7071||og->nlat.c<.7071)
67 28994509 2004-04-21 devnull return(ckcut(g,og,0.));
68 28994509 2004-04-21 devnull return(1);
69 28994509 2004-04-21 devnull }
70 28994509 2004-04-21 devnull
71 28994509 2004-04-21 devnull static int
72 28994509 2004-04-21 devnull Xsquare(struct place *place, double *x, double *y)
73 28994509 2004-04-21 devnull {
74 28994509 2004-04-21 devnull double z1,z2;
75 28994509 2004-04-21 devnull double r, theta;
76 28994509 2004-04-21 devnull struct place p;
77 28994509 2004-04-21 devnull copyplace(place,&p);
78 28994509 2004-04-21 devnull if(place->nlat.l<0) {
79 28994509 2004-04-21 devnull p.nlat.l = -p.nlat.l;
80 28994509 2004-04-21 devnull p.nlat.s = -p.nlat.s;
81 28994509 2004-04-21 devnull }
82 28994509 2004-04-21 devnull if(p.nlat.l<FUZZ && fabs(p.wlon.l)>PI-FUZZ){
83 28994509 2004-04-21 devnull *y = -gyside/2;
84 28994509 2004-04-21 devnull *x = p.wlon.l>0?0:gyside;
85 28994509 2004-04-21 devnull return(1);
86 28994509 2004-04-21 devnull }
87 28994509 2004-04-21 devnull Xstereographic(&p,&z1,&z2);
88 28994509 2004-04-21 devnull r = sqrt(sqrt(hypot(z1,z2)/2));
89 28994509 2004-04-21 devnull theta = atan2(z1,-z2)/4;
90 28994509 2004-04-21 devnull dosquare(r*sin(theta),-r*cos(theta),x,y);
91 28994509 2004-04-21 devnull if(place->nlat.l<0)
92 28994509 2004-04-21 devnull *y = -gyside - *y;
93 28994509 2004-04-21 devnull return(1);
94 28994509 2004-04-21 devnull }
95 28994509 2004-04-21 devnull
96 28994509 2004-04-21 devnull proj
97 28994509 2004-04-21 devnull square(void)
98 28994509 2004-04-21 devnull {
99 28994509 2004-04-21 devnull guyou();
100 28994509 2004-04-21 devnull return(Xsquare);
101 28994509 2004-04-21 devnull }