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