#include "map.h" #define NTRACK 10	/* max number of -t and -u files */#define NFILE 30	/* max number of map files */#define NVERT 20	/* max number of vertices in a -v polygon */#define HALFWIDTH 8192	/* output scaled to fit in -HALFWIDTH,HALFWIDTH */#define LONGLINES (HALFWIDTH*4)	/* permissible segment lengths */#define SHORTLINES (HALFWIDTH/8)#define SCALERATIO 10	/* of abs to rel data (see map(5)) */#define RESOL 2.	/* coarsest resolution for tracing grid (degrees) */#define TWO_THRD 0.66666666666666667 int normproj(float, float, float *, float *);int posproj(float, float, float *, float *);int picut(struct place *, struct place *, float *);float reduce(float);short getshort(FILE *);char *mapindex(char *);proj projection;  static char *mapdir = "/lib/map";	/* default map directory */static char *file[NFILE+1] = {	/* list of map files */	"world",		/* default map */	0};  extern struct index index[];int halfwidth = HALFWIDTH; static int (*cut)(struct place *, struct place *, float *);static int poles;static float orientation[3] = { 90., 0., 0. };	/* -o option */static oriented;	/* nonzero if -o option occurred */static upright;		/* 1 if orientation[0]==90, -1 if -90, else 0*/static int delta = 1;	/* -d setting */static float limits[4] = {	/* -l parameters */	-90., 90., -180., 180.};static float klimits[4] = {	/* -k parameters */	-90., 90., -180., 180.};static int limcase;static float rlimits[4];	/* limits expressed in radians */static float lolat, hilat, lolon, hilon;static float window[4] = {	/* option -w */	-90., 90., -180., 180.};static windowed;	/* nozero if option -w */static struct vert { float x, y; } v[NVERT+2];	/*clipping polygon*/static struct edge { float a, b, c; } e[NVERT]; /* coeffs for linear inequality */static int nvert;	/* number of vertices in clipping polygon */ static float rwindow[4];	/* window, expressed in radians */static float params[2];		/* projection params *//* bounds on output values before scaling; found by coarse survey */static float xmin = 100.;static float xmax = -100.;static float ymin = 100.;static float ymax = -100.;static float xcent, ycent;static float xoff, yoff;float xrange, yrange;static int left = -HALFWIDTH;static int right = HALFWIDTH;static int bottom = -HALFWIDTH;static int top = HALFWIDTH;static int longlines = SHORTLINES; /* drop longer segments */static int shortlines = SHORTLINES;static int bflag = 1;	/* 0 for option -b */static int sflag = 0;	/* 1 for option -s */static int rflag = 0;	/* 1 for option -r */static int mflag = 0;	/* 1 if option -m occurred */static int kflag = 0;	/* 1 if option -k occurred */static float position[3];	/* option -p */static float center[3] = {0., 0., 0.};	/* option -c */static struct coord crot;		/* option -c */static float grid[3] = { 10., 10., RESOL };	/* option -g */static float dlat, dlon;	/* resolution for tracing grid in lat and lon */static float scaling;	/* to compute final integer output */static struct track {	/* options -t and -u */	int tracktyp;	/* 't' or 'u' */	char *tracknam;	/* name of input file */} track[NTRACK];static int ntrack;	/* number of tracks present */static char *symbolfile;	/* option -y */ void	clamp(float *px, float v);void	clipinit(void);float	diddle(struct place *, float, float);float	diddle(struct place *, float, float);void	dobounds(float, float, float, float, int);void	dogrid(float, float, float, float);int	duple(struct place *, float);double	fmax(double, double);double	fmin(double, double);void	getdata(char *);int	gridpt(float, float, int);int	inpoly(float, float);int	inwindow(struct place *);void	pathnames(void);int	pnorm(float);void	radbds(float *w, float *rw);void	revlon(struct place *, float);void	satellite(struct track *);int	seeable(float, float);void	windlim(void);void	realcut(void); intoption(char *s) { 	if(s[0]=='-' && (s[1]<'0'||s[1]>'9'))		return(s[1]!='.'&&s[1]!=0);	else		return(0);} voidconv(int k, struct coord *g){	g->l = (0.0001/SCALERATIO)*k;	sincos(g);} intmain(int argc, char *argv[]){	int i,k;	char *s, *t;	float x, y;	double lat, lon;	float *wlim;	double dd;	if(sizeof(short)!=2)		abort();	/* getshort() won't work */	s = getenv("MAP");	if(s)		file[0] = s;	s = getenv("MAPDIR");	if(s)		mapdir = s;	if(argc<=1) 		error("usage: map projection params options");	for(k=0;index[k].name;k++) {		s = index[k].name;		t = argv[1];		while(*s == *t){			if(*s==0) goto found;			s++;			t++;		}	}	fprintf(stderr,"projections:\n");	for(i=0;index[i].name;i++)  {		fprintf(stderr,"%s",index[i].name);		for(k=0; k<index[i].npar; k++)			fprintf(stderr," p%d", k);		fprintf(stderr,"\n");	}	exits("error");found:	argv += 2;	argc -= 2;	cut = index[k].cut;	poles = index[k].poles;	for(i=0;i<index[k].npar;i++) {		if(i>=argc||option(argv[i])) {			fprintf(stderr,"%s needs %d params\n",index[k].name,index[k].npar);			exits("error");		}		params[i] = atof(argv[i]);	}	argv += i;	argc -= i;	while(argc>0&&option(argv[0])) {		argc--;		argv++;		switch(argv[-1][1]) {		case 'm':			i = 0;			if(!mflag) while(file[i]!=0)					i++;			for(i=0;i<NFILE&&argc>i&&!option(argv[i]);i++)				file[i] = argv[i];			file[i] = 0;			mflag++;			argc -= i;			argv += i;			break;		case 'b':			bflag = 0;			for(nvert=0;nvert<NVERT&&argc>=2;nvert++) {				if(option(*argv))					break;				v[nvert].x = atof(*argv++);				argc--;				if(option(*argv))					break;				v[nvert].y = atof(*argv++);				argc--;			}			if(nvert>=NVERT)				error("too many clipping vertices");			break;		case 'g':			for(i=0;i<3&&argc>i&&!option(argv[i]);i++)				grid[i] = atof(argv[i]);			switch(i) {			case 0:				grid[0] = grid[1] = 0.;				break;			case 1:				grid[1] = grid[0];			}			argc -= i;			argv += i;			break;		case 't':		case 'u':			for(i=0;ntrack<NTRACK&&argc>i&&!option(argv[i]);i++) {				track[ntrack].tracktyp = argv[-1][1];				track[ntrack++].tracknam = argv[i];			}			argc -= i;			argv +=i;			break;		case 'r':			rflag++;			break;		case 's':			sflag++;			break;		case 'o':			for(i=0;i<3&&i<argc&&!option(argv[i]);i++)				orientation[i] = atof(argv[i]);			oriented++;			argv += i;			argc -= i;			break;		case 'l':			for(i=0;i<argc&&i<4&&!option(argv[i]);i++)				limits[i] = atof(argv[i]);			argv += i;			argc -= i;			break;		case 'k':			kflag++;			for(i=0;i<argc&&i<4&&!option(argv[i]);i++)				klimits[i] = atof(argv[i]);			argv += i;			argc -= i;			break;		case 'd':			if(argc>0&&!option(argv[0])) {				delta = atoi(argv[0]);				argv++;				argc--;			}			break;		case 'w':			windowed++;			for(i=0;i<argc&&i<4&&!option(argv[i]);i++)				window[i] = atof(argv[i]);			argv += i;			argc -= i;			break;		case 'c':			for(i=0;i<3&&argc>i&&!option(argv[i]);i++) 				center[i] = atof(argv[i]);			argc -= i;			argv += i;			break;		case 'p':			for(i=0;i<3&&argc>i&&!option(argv[i]);i++)				position[i] = atof(argv[i]);			argc -= i;			argv += i;			if(i!=3||position[2]<=0) 				error("incomplete positioning");			break;		case 'y':			if(argc>0&&!option(argv[0]))				symbolfile = argv[0];			argc--;			argv++;			break;		}	}	if(argc>0)		error("error in arguments");	pathnames();	clamp(&limits[0],-90.);	clamp(&limits[1],90.);	clamp(&klimits[0],-90.);	clamp(&klimits[1],90.);	clamp(&window[0],-90.);	clamp(&window[1],90.);	radbds(limits,rlimits);	limcase = limits[2]<-180.?0:		  limits[3]>180.?2:		  1;	if(		window[0]>=window[1]||		window[2]>=window[3]||		window[0]>90.||		window[1]<-90.||		window[2]>180.||		window[3]<-180.)		error("unreasonable window");	windlim();	radbds(window,rwindow);	upright = orientation[0]==90? 1: orientation[0]==-90? -1: 0;	if(index[k].spheroid && !upright)		error("can't tilt the spheroid");	if(limits[2]>limits[3])		limits[3] += 360;	if(!oriented)		orientation[2] = (limits[2]+limits[3])/2;	orient(orientation[0],orientation[1],orientation[2]);	projection = (*index[k].prog)(params[0],params[1]);	if(projection == 0)		error("unreasonable projection parameters");	clipinit();	grid[0] = fabs(grid[0]);	grid[1] = fabs(grid[1]);	if(!kflag)		for(i=0;i<4;i++)			klimits[i] = limits[i];	if(klimits[2]>klimits[3])		klimits[3] += 360;	lolat = limits[0];	hilat = limits[1];	lolon = limits[2];	hilon = limits[3];	if(lolon>=hilon||lolat>=hilat||lolat<-90.||hilat>90.)		error("unreasonable limits");	wlim = kflag? klimits: window;	dlat = fmin(hilat-lolat,wlim[1]-wlim[0])/16;	dlon = fmin(hilon-lolon,wlim[3]-wlim[2])/32;	dd = fmax(dlat,dlon);	while(grid[2]>fmin(dlat,dlon)/2)		grid[2] /= 2;	realcut();	if(nvert<=0) {		for(lat=klimits[0];lat<klimits[1]+dd-FUZZ;lat+=dd) {			if(lat>klimits[1])				lat = klimits[1];			for(lon=klimits[2];lon<klimits[3]+dd-FUZZ;lon+=dd) {				if((kflag?posproj:normproj)					(lat,lon+(lon<klimits[3]?FUZZ:-FUZZ),				   &x,&y)<=0)					continue;				if(x<xmin) xmin = x;				if(x>xmax) xmax = x;				if(y<ymin) ymin = y;				if(y>ymax) ymax = y;			}		}	} else {		for(i=0; i<nvert; i++) {			x = v[i].x;			y = v[i].y;			if(x<xmin) xmin = x;			if(x>xmax) xmax = x;			if(y<ymin) ymin = y;			if(y>ymax) ymax = y;		}	}	xrange = xmax - xmin;	yrange = ymax - ymin;	if(xrange<=0||yrange<=0)		error("map seems to be empty");	scaling = 2;	/*plotting area from -1 to 1*/	if(position[2]!=0) {		if(posproj(position[0]-.5,position[1],&xcent,&ycent)<=0||		   posproj(position[0]+.5,position[1],&x,&y)<=0)			error("unreasonable position");		scaling /= (position[2]*hypot(x-xcent,y-ycent));		if(posproj(position[0],position[1],&xcent,&ycent)<=0)			error("unreasonable position");	} else {		scaling /= (xrange>yrange?xrange:yrange);		xcent = (xmin+xmax)/2;		ycent = (ymin+ymax)/2;	}	xoff = center[0]/scaling;	yoff = center[1]/scaling;	crot.l = center[2]*RAD;	sincos(&crot);	scaling *= HALFWIDTH*0.9;	if(symbolfile) 		getsyms(symbolfile);	if(!sflag) {		openpl();		range(left,bottom,right,top);		erase();	}	pen("dotted");	if(grid[0]>0.)		for(lat=ceil(lolat/grid[0])*grid[0];		    lat<=hilat;lat+=grid[0]) 			dogrid(lat,lat,lolon,hilon);	if(grid[1]>0.)		for(lon=ceil(lolon/grid[1])*grid[1];		    lon<=hilon;lon+=grid[1]) 			dogrid(lolat,hilat,lon,lon);	pen("solid");	if(bflag) {		dobounds(lolat,hilat,lolon,hilon,0);		dobounds(window[0],window[1],window[2],window[3],1);	}	lolat = floor(limits[0]/10)*10;	hilat = ceil(limits[1]/10)*10;	lolon = floor(limits[2]/10)*10;	hilon = ceil(limits[3]/10)*10;	if(lolon>hilon)		hilon += 360.;	/*do tracks first so as not to lose the standard input*/	for(i=0;i<ntrack;i++) {		longlines = LONGLINES;		satellite(&track[i]);		longlines = shortlines;	}	pen("solid");	for(i=0;file[i];i++)		getdata(file[i]); 	move(right,bottom);	return 0;}  intnormproj(float lat, float lon, float *x, float *y){	int i;	struct place geog;	latlon(lat,lon,&geog);/*	printp(&geog);*/	normalize(&geog);	if(!inwindow(&geog))		return(-1);	i = (*projection)(&geog,x,y);	if(rflag) 		*x = -*x;/*	printp(&geog);	fprintf(stderr,"%d %.3f %.3f\n",i,*x,*y);*/	return(i);} intposproj(float lat, float lon, float *x, float *y){	int i;	struct place geog;	latlon(lat,lon,&geog);	normalize(&geog);	i = (*projection)(&geog,x,y);	if(rflag) 		*x = -*x;	return(i);} intinwindow(struct place *geog){	if(geog->nlat.l<rwindow[0]||	   geog->nlat.l>rwindow[1]||	   geog->wlon.l<rwindow[2]||	   geog->wlon.l>rwindow[3])		return(0);	else return(1);} intinlimits(struct place *g){	if(rlimits[0]>g->nlat.l||	   rlimits[1]<g->nlat.l)		return(0);	switch(limcase) {	case 0:		if(rlimits[2]+TWOPI>g->wlon.l&&		   rlimits[3]<g->wlon.l)			return(0);		break;	case 1:		if(rlimits[2]>g->wlon.l||		   rlimits[3]<g->wlon.l)			return(0);		break;	case 2:		if(rlimits[2]>g->wlon.l&&

Lines 1–500 of 1121. Load more