plan9port

fork of plan9port with libvec, libstr and libsdb
Log | Files | Refs | README | LICENSE

perspective.c (1445B)


      1 #include <u.h>
      2 #include <libc.h>
      3 #include "map.h"
      4 
      5 #define ORTHRAD 1000
      6 static double viewpt;
      7 
      8 static int
      9 Xperspective(struct place *place, double *x, double *y)
     10 {
     11 	double r;
     12 	if(viewpt<=1+FUZZ && fabs(place->nlat.s)<=viewpt+.01)
     13 		return(-1);
     14 	r = place->nlat.c*(viewpt - 1.)/(viewpt - place->nlat.s);
     15 	*x = - r*place->wlon.s;
     16 	*y = - r*place->wlon.c;
     17 	if(r>4.)
     18 		return(-1);
     19 	if(fabs(viewpt)>1 && place->nlat.s<1/viewpt ||
     20 	   fabs(viewpt)<=1 && place->nlat.s<viewpt)
     21 			return 0;
     22 	return(1);
     23 }
     24 
     25 proj
     26 perspective(double radius)
     27 {
     28 	viewpt = radius;
     29 	if(viewpt >= ORTHRAD)
     30 		return(Xorthographic);
     31 	if(fabs(viewpt-1.)<.0001)
     32 		return(0);
     33 	return(Xperspective);
     34 }
     35 
     36 	/* called from various conformal projections,
     37            but not from stereographic itself */
     38 int
     39 Xstereographic(struct place *place, double *x, double *y)
     40 {
     41 	double v = viewpt;
     42 	int retval;
     43 	viewpt = -1;
     44 	retval = Xperspective(place, x, y);
     45 	viewpt = v;
     46 	return retval;
     47 }
     48 
     49 proj
     50 stereographic(void)
     51 {
     52 	viewpt = -1.;
     53 	return(Xperspective);
     54 }
     55 
     56 proj
     57 gnomonic(void)
     58 {
     59 	viewpt = 0.;
     60 	return(Xperspective);
     61 }
     62 
     63 int
     64 plimb(double *lat, double *lon, double res)
     65 {
     66 	static int first = 1;
     67 	if(viewpt >= ORTHRAD)
     68 		return olimb(lat, lon, res);
     69 	if(first) {
     70 		first = 0;
     71 		*lon = -180;
     72 		if(fabs(viewpt) < .01)
     73 			*lat = 0;
     74 		else if(fabs(viewpt)<=1)
     75 			*lat = asin(viewpt)/RAD;
     76 		else
     77 			*lat = asin(1/viewpt)/RAD;
     78 	} else
     79 		*lon += res;
     80 	if(*lon <= 180)
     81 		return 1;
     82 	first = 1;
     83 	return -1;
     84 }