|
|
1.1 ! root 1: #include "map.h" ! 2: ! 3: #define NVERT 20 /* max number of vertices in a -v polygon */ ! 4: #define HALFWIDTH 8192 /* output scaled to fit in -HALFWIDTH,HALFWIDTH */ ! 5: #define LONGLINES (HALFWIDTH*4) /* permissible segment lengths */ ! 6: #define SHORTLINES (HALFWIDTH/8) ! 7: #define SCALERATIO 10 /* of abs to rel data (see map(5)) */ ! 8: #define RESOL 2. /* coarsest resolution for tracing grid (degrees) */ ! 9: #define TWO_THRD 0.66666666666666667 ! 10: ! 11: int normproj(double, double, double *, double *); ! 12: int posproj(double, double, double *, double *); ! 13: int picut(struct place *, struct place *, double *); ! 14: double reduce(double); ! 15: short getshort(FILE *); ! 16: char *mapindex(char *); ! 17: proj projection; ! 18: ! 19: ! 20: static char *mapdir = "/usr/dict"; /* default map directory */ ! 21: struct file { ! 22: char *name; ! 23: char *color; ! 24: char *style; ! 25: }; ! 26: static struct file dfltfile = { ! 27: "world", BLACK, SOLID /* default map */ ! 28: }; ! 29: static struct file *file = &dfltfile; /* list of map files */ ! 30: static int nfile = 1; /* length of list */ ! 31: static char *currcolor = BLACK; /* current color */ ! 32: static char *gridcolor = BLACK; ! 33: static char *bordcolor = BLACK; ! 34: ! 35: extern struct index index[]; ! 36: int halfwidth = HALFWIDTH; ! 37: ! 38: static int (*cut)(struct place *, struct place *, double *); ! 39: static int (*limb)(double*, double*, double); ! 40: static void dolimb(void); ! 41: static int onlimb; ! 42: static int poles; ! 43: static double orientation[3] = { 90., 0., 0. }; /* -o option */ ! 44: static oriented; /* nonzero if -o option occurred */ ! 45: static upright; /* 1 if orientation[0]==90, -1 if -90, else 0*/ ! 46: static int delta = 1; /* -d setting */ ! 47: static double limits[4] = { /* -l parameters */ ! 48: -90., 90., -180., 180. ! 49: }; ! 50: static double klimits[4] = { /* -k parameters */ ! 51: -90., 90., -180., 180. ! 52: }; ! 53: static int limcase; ! 54: static double rlimits[4]; /* limits expressed in radians */ ! 55: static double lolat, hilat, lolon, hilon; ! 56: static double window[4] = { /* option -w */ ! 57: -90., 90., -180., 180. ! 58: }; ! 59: static windowed; /* nozero if option -w */ ! 60: static struct vert { double x, y; } v[NVERT+2]; /*clipping polygon*/ ! 61: static struct edge { double a, b, c; } e[NVERT]; /* coeffs for linear inequality */ ! 62: static int nvert; /* number of vertices in clipping polygon */ ! 63: ! 64: static double rwindow[4]; /* window, expressed in radians */ ! 65: static double params[2]; /* projection params */ ! 66: /* bounds on output values before scaling; found by coarse survey */ ! 67: static double xmin = 100.; ! 68: static double xmax = -100.; ! 69: static double ymin = 100.; ! 70: static double ymax = -100.; ! 71: static double xcent, ycent; ! 72: static double xoff, yoff; ! 73: double xrange, yrange; ! 74: static int left = -HALFWIDTH; ! 75: static int right = HALFWIDTH; ! 76: static int bottom = -HALFWIDTH; ! 77: static int top = HALFWIDTH; ! 78: static int longlines = SHORTLINES; /* drop longer segments */ ! 79: static int shortlines = SHORTLINES; ! 80: static int bflag = 1; /* 0 for option -b */ ! 81: static int s1flag = 0; /* 1 for option -s1 */ ! 82: static int s2flag = 0; /* 1 for option -s2 */ ! 83: static int rflag = 0; /* 1 for option -r */ ! 84: static int kflag = 0; /* 1 if option -k occurred */ ! 85: static int xflag = 0; /* 1 if option -x occurred */ ! 86: int vflag = 1; /* -1 if option -v occurred */ ! 87: static double position[3]; /* option -p */ ! 88: static double center[3] = {0., 0., 0.}; /* option -c */ ! 89: static struct coord crot; /* option -c */ ! 90: static double grid[3] = { 10., 10., RESOL }; /* option -g */ ! 91: static double dlat, dlon; /* resolution for tracing grid in lat and lon */ ! 92: static double scaling; /* to compute final integer output */ ! 93: static struct file *track; /* options -t and -u */ ! 94: static int ntrack; /* number of tracks present */ ! 95: static char *symbolfile; /* option -y */ ! 96: ! 97: void clamp(double *px, double v); ! 98: void clipinit(void); ! 99: double diddle(struct place *, double, double); ! 100: double diddle(struct place *, double, double); ! 101: void dobounds(double, double, double, double, int); ! 102: void dogrid(double, double, double, double); ! 103: int duple(struct place *, double); ! 104: double fmax(double, double); ! 105: double fmin(double, double); ! 106: void getdata(char *); ! 107: int gridpt(double, double, int); ! 108: int inpoly(double, double); ! 109: int inwindow(struct place *); ! 110: void pathnames(void); ! 111: int pnorm(double); ! 112: void radbds(double *w, double *rw); ! 113: void revlon(struct place *, double); ! 114: void satellite(struct file *); ! 115: int seeable(double, double); ! 116: void windlim(void); ! 117: void realcut(void); ! 118: ! 119: int ! 120: option(char *s) ! 121: { ! 122: ! 123: if(s[0]=='-' && (s[1]<'0'||s[1]>'9')) ! 124: return(s[1]!='.'&&s[1]!=0); ! 125: else ! 126: return(0); ! 127: } ! 128: ! 129: void ! 130: conv(int k, struct coord *g) ! 131: { ! 132: g->l = (0.0001/SCALERATIO)*k; ! 133: sincos(g); ! 134: } ! 135: ! 136: int ! 137: main(int argc, char *argv[]) ! 138: { ! 139: int i,k; ! 140: char *s, *t, *style; ! 141: double x, y; ! 142: double lat, lon; ! 143: double *wlim; ! 144: double dd; ! 145: if(sizeof(short)!=2) ! 146: abort(); /* getshort() won't work */ ! 147: s = getenv("MAP"); ! 148: if(s) ! 149: file[0].name = s; ! 150: s = getenv("MAPDIR"); ! 151: if(s) ! 152: mapdir = s; ! 153: if(argc<=1) ! 154: error("usage: map projection params options"); ! 155: for(k=0;index[k].name;k++) { ! 156: s = index[k].name; ! 157: t = argv[1]; ! 158: while(*s == *t){ ! 159: if(*s==0) goto found; ! 160: s++; ! 161: t++; ! 162: } ! 163: } ! 164: fprintf(stderr,"projections:\n"); ! 165: for(i=0;index[i].name;i++) { ! 166: fprintf(stderr,"%s",index[i].name); ! 167: for(k=0; k<index[i].npar; k++) ! 168: fprintf(stderr," p%d", k); ! 169: fprintf(stderr,"\n"); ! 170: } ! 171: exit(1); ! 172: found: ! 173: argv += 2; ! 174: argc -= 2; ! 175: cut = index[k].cut; ! 176: limb = index[k].limb; ! 177: poles = index[k].poles; ! 178: for(i=0;i<index[k].npar;i++) { ! 179: if(i>=argc||option(argv[i])) { ! 180: fprintf(stderr,"%s needs %d params\n",index[k].name,index[k].npar); ! 181: exit(1); ! 182: } ! 183: params[i] = atof(argv[i]); ! 184: } ! 185: argv += i; ! 186: argc -= i; ! 187: while(argc>0&&option(argv[0])) { ! 188: argc--; ! 189: argv++; ! 190: switch(argv[-1][1]) { ! 191: case 'm': ! 192: if(file == &dfltfile) { ! 193: file = 0; ! 194: nfile = 0; ! 195: } ! 196: while(argc && !option(*argv)) { ! 197: file = realloc(file,(nfile+1)*sizeof(*file)); ! 198: file[nfile].name = *argv; ! 199: file[nfile].color = currcolor; ! 200: file[nfile].style = SOLID; ! 201: nfile++; ! 202: argv++; ! 203: argc--; ! 204: } ! 205: break; ! 206: case 'b': ! 207: bflag = 0; ! 208: for(nvert=0;nvert<NVERT&&argc>=2;nvert++) { ! 209: if(option(*argv)) ! 210: break; ! 211: v[nvert].x = atof(*argv++); ! 212: argc--; ! 213: if(option(*argv)) ! 214: break; ! 215: v[nvert].y = atof(*argv++); ! 216: argc--; ! 217: } ! 218: if(nvert>=NVERT) ! 219: error("too many clipping vertices"); ! 220: break; ! 221: case 'g': ! 222: gridcolor = currcolor; ! 223: for(i=0;i<3&&argc>i&&!option(argv[i]);i++) ! 224: grid[i] = atof(argv[i]); ! 225: switch(i) { ! 226: case 0: ! 227: grid[0] = grid[1] = 0.; ! 228: break; ! 229: case 1: ! 230: grid[1] = grid[0]; ! 231: } ! 232: argc -= i; ! 233: argv += i; ! 234: break; ! 235: case 't': ! 236: style = SOLID; ! 237: goto casetu; ! 238: case 'u': ! 239: style = DOTDASH; ! 240: casetu: ! 241: while(argc && !option(*argv)) { ! 242: track = realloc(track,(ntrack+1)*sizeof(*track)); ! 243: track[ntrack].name = *argv; ! 244: track[ntrack].color = currcolor; ! 245: track[ntrack].style = style; ! 246: ntrack++; ! 247: argv++; ! 248: argc--; ! 249: } ! 250: break; ! 251: case 'r': ! 252: rflag++; ! 253: break; ! 254: case 's': ! 255: switch(argv[-1][2]) { ! 256: case '1': ! 257: s1flag++; ! 258: break; ! 259: case 0: /* compatibility */ ! 260: case '2': ! 261: s2flag++; ! 262: } ! 263: break; ! 264: case 'o': ! 265: for(i=0;i<3&&i<argc&&!option(argv[i]);i++) ! 266: orientation[i] = atof(argv[i]); ! 267: oriented++; ! 268: argv += i; ! 269: argc -= i; ! 270: break; ! 271: case 'l': ! 272: bordcolor = currcolor; ! 273: for(i=0;i<argc&&i<4&&!option(argv[i]);i++) ! 274: limits[i] = atof(argv[i]); ! 275: argv += i; ! 276: argc -= i; ! 277: break; ! 278: case 'k': ! 279: kflag++; ! 280: for(i=0;i<argc&&i<4&&!option(argv[i]);i++) ! 281: klimits[i] = atof(argv[i]); ! 282: argv += i; ! 283: argc -= i; ! 284: break; ! 285: case 'd': ! 286: if(argc>0&&!option(argv[0])) { ! 287: delta = atoi(argv[0]); ! 288: argv++; ! 289: argc--; ! 290: } ! 291: break; ! 292: case 'w': ! 293: bordcolor = currcolor; ! 294: windowed++; ! 295: for(i=0;i<argc&&i<4&&!option(argv[i]);i++) ! 296: window[i] = atof(argv[i]); ! 297: argv += i; ! 298: argc -= i; ! 299: break; ! 300: case 'c': ! 301: for(i=0;i<3&&argc>i&&!option(argv[i]);i++) ! 302: center[i] = atof(argv[i]); ! 303: argc -= i; ! 304: argv += i; ! 305: break; ! 306: case 'p': ! 307: for(i=0;i<3&&argc>i&&!option(argv[i]);i++) ! 308: position[i] = atof(argv[i]); ! 309: argc -= i; ! 310: argv += i; ! 311: if(i!=3||position[2]<=0) ! 312: error("incomplete positioning"); ! 313: break; ! 314: case 'y': ! 315: if(argc>0&&!option(argv[0])) { ! 316: symbolfile = argv[0]; ! 317: argc--; ! 318: argv++; ! 319: } ! 320: break; ! 321: case 'v': ! 322: if(index[k].limb == 0) ! 323: error("-v does not apply here"); ! 324: vflag = -1; ! 325: break; ! 326: case 'x': ! 327: xflag = 1; ! 328: break; ! 329: case 'C': ! 330: if(argc && !option(*argv)) { ! 331: currcolor = *argv; ! 332: argc--; ! 333: argv++; ! 334: } ! 335: break; ! 336: } ! 337: } ! 338: if(argc>0) ! 339: error("error in arguments"); ! 340: pathnames(); ! 341: clamp(&limits[0],-90.); ! 342: clamp(&limits[1],90.); ! 343: clamp(&klimits[0],-90.); ! 344: clamp(&klimits[1],90.); ! 345: clamp(&window[0],-90.); ! 346: clamp(&window[1],90.); ! 347: radbds(limits,rlimits); ! 348: limcase = limits[2]<-180.?0: ! 349: limits[3]>180.?2: ! 350: 1; ! 351: if( ! 352: window[0]>=window[1]|| ! 353: window[2]>=window[3]|| ! 354: window[0]>90.|| ! 355: window[1]<-90.|| ! 356: window[2]>180.|| ! 357: window[3]<-180.) ! 358: error("unreasonable window"); ! 359: windlim(); ! 360: radbds(window,rwindow); ! 361: upright = orientation[0]==90? 1: orientation[0]==-90? -1: 0; ! 362: if(index[k].spheroid && !upright) ! 363: error("can't tilt the spheroid"); ! 364: if(limits[2]>limits[3]) ! 365: limits[3] += 360; ! 366: if(!oriented) ! 367: orientation[2] = (limits[2]+limits[3])/2; ! 368: orient(orientation[0],orientation[1],orientation[2]); ! 369: projection = (*index[k].prog)(params[0],params[1]); ! 370: if(projection == 0) ! 371: error("unreasonable projection parameters"); ! 372: clipinit(); ! 373: grid[0] = fabs(grid[0]); ! 374: grid[1] = fabs(grid[1]); ! 375: if(!kflag) ! 376: for(i=0;i<4;i++) ! 377: klimits[i] = limits[i]; ! 378: if(klimits[2]>klimits[3]) ! 379: klimits[3] += 360; ! 380: lolat = limits[0]; ! 381: hilat = limits[1]; ! 382: lolon = limits[2]; ! 383: hilon = limits[3]; ! 384: if(lolon>=hilon||lolat>=hilat||lolat<-90.||hilat>90.) ! 385: error("unreasonable limits"); ! 386: wlim = kflag? klimits: window; ! 387: dlat = fmin(hilat-lolat,wlim[1]-wlim[0])/16; ! 388: dlon = fmin(hilon-lolon,wlim[3]-wlim[2])/32; ! 389: dd = fmax(dlat,dlon); ! 390: while(grid[2]>fmin(dlat,dlon)/2) ! 391: grid[2] /= 2; ! 392: realcut(); ! 393: if(nvert<=0) { ! 394: for(lat=klimits[0];lat<klimits[1]+dd-FUZZ;lat+=dd) { ! 395: if(lat>klimits[1]) ! 396: lat = klimits[1]; ! 397: for(lon=klimits[2];lon<klimits[3]+dd-FUZZ;lon+=dd) { ! 398: i = (kflag?posproj:normproj) ! 399: (lat,lon+(lon<klimits[3]?FUZZ:-FUZZ), ! 400: &x,&y); ! 401: if(i*vflag <= 0) ! 402: continue; ! 403: if(x<xmin) xmin = x; ! 404: if(x>xmax) xmax = x; ! 405: if(y<ymin) ymin = y; ! 406: if(y>ymax) ymax = y; ! 407: } ! 408: } ! 409: } else { ! 410: for(i=0; i<nvert; i++) { ! 411: x = v[i].x; ! 412: y = v[i].y; ! 413: if(x<xmin) xmin = x; ! 414: if(x>xmax) xmax = x; ! 415: if(y<ymin) ymin = y; ! 416: if(y>ymax) ymax = y; ! 417: } ! 418: } ! 419: xrange = xmax - xmin; ! 420: yrange = ymax - ymin; ! 421: if(xrange<=0||yrange<=0) ! 422: error("map seems to be empty"); ! 423: scaling = 2; /*plotting area from -1 to 1*/ ! 424: if(position[2]!=0) { ! 425: if(posproj(position[0]-.5,position[1],&xcent,&ycent)==0|| ! 426: posproj(position[0]+.5,position[1],&x,&y)==0) ! 427: error("unreasonable position"); ! 428: scaling /= (position[2]*hypot(x-xcent,y-ycent)); ! 429: if(posproj(position[0],position[1],&xcent,&ycent)==0) ! 430: error("unreasonable position"); ! 431: } else { ! 432: scaling /= (xrange>yrange?xrange:yrange); ! 433: xcent = (xmin+xmax)/2; ! 434: ycent = (ymin+ymax)/2; ! 435: } ! 436: xoff = center[0]/scaling; ! 437: yoff = center[1]/scaling; ! 438: crot.l = center[2]*RAD; ! 439: sincos(&crot); ! 440: scaling *= HALFWIDTH*0.9; ! 441: if(symbolfile) ! 442: getsyms(symbolfile); ! 443: if(!s2flag) { ! 444: openpl(); ! 445: erase(); ! 446: } ! 447: range(left,bottom,right,top); ! 448: color(gridcolor); ! 449: pen(DOTTED); ! 450: if(grid[0]>0.) ! 451: for(lat=ceil(lolat/grid[0])*grid[0]; ! 452: lat<=hilat;lat+=grid[0]) ! 453: dogrid(lat,lat,lolon,hilon); ! 454: if(grid[1]>0.) ! 455: for(lon=ceil(lolon/grid[1])*grid[1]; ! 456: lon<=hilon;lon+=grid[1]) ! 457: dogrid(lolat,hilat,lon,lon); ! 458: color(bordcolor); ! 459: pen(SOLID); ! 460: if(bflag) { ! 461: dolimb(); ! 462: dobounds(lolat,hilat,lolon,hilon,0); ! 463: dobounds(window[0],window[1],window[2],window[3],1); ! 464: } ! 465: lolat = floor(limits[0]/10)*10; ! 466: hilat = ceil(limits[1]/10)*10; ! 467: lolon = floor(limits[2]/10)*10; ! 468: hilon = ceil(limits[3]/10)*10; ! 469: if(lolon>hilon) ! 470: hilon += 360.; ! 471: /*do tracks first so as not to lose the standard input*/ ! 472: for(i=0;i<ntrack;i++) { ! 473: longlines = LONGLINES; ! 474: satellite(&track[i]); ! 475: longlines = shortlines; ! 476: } ! 477: for(i=0;i<nfile;i++) { ! 478: pen(file[i].style); ! 479: color(file[i].color); ! 480: getdata(file[i].name); ! 481: } ! 482: move(right,bottom); ! 483: if(!s1flag) ! 484: closepl(); ! 485: return 0; ! 486: } ! 487: ! 488: /* Out of perverseness (really to recover from a dubious, ! 489: but documented, convention) the returns from projection ! 490: functions (-1 unplottable, 0 wrong sheet, 1 good) are ! 491: recoded into -1 wrong sheet, 0 unplottable, 1 good. */ ! 492: ! 493: int ! 494: fixproj(struct place *g, double *x, double *y) ! 495: { ! 496: int i = (*projection)(g,x,y); ! 497: return i<0? 0: i==0? -1: 1; ! 498: } ! 499: ! 500: int ! 501: normproj(double lat, double lon, double *x, double *y) ! 502: { ! 503: int i; ! 504: struct place geog; ! 505: latlon(lat,lon,&geog); ! 506: /* ! 507: printp(&geog); ! 508: */ ! 509: normalize(&geog); ! 510: if(!inwindow(&geog)) ! 511: return(-1); ! 512: i = fixproj(&geog,x,y); ! 513: if(rflag) ! 514: *x = -*x; ! 515: /* ! 516: printp(&geog); ! 517: fprintf(stderr,"%d %.3f %.3f\n",i,*x,*y); ! 518: */ ! 519: return(i); ! 520: } ! 521: ! 522: int ! 523: posproj(double lat, double lon, double *x, double *y) ! 524: { ! 525: int i; ! 526: struct place geog; ! 527: latlon(lat,lon,&geog); ! 528: normalize(&geog); ! 529: i = fixproj(&geog,x,y); ! 530: if(rflag) ! 531: *x = -*x; ! 532: return(i); ! 533: } ! 534: ! 535: int ! 536: inwindow(struct place *geog) ! 537: { ! 538: if(geog->nlat.l<rwindow[0]-FUZZ|| ! 539: geog->nlat.l>rwindow[1]+FUZZ|| ! 540: geog->wlon.l<rwindow[2]-FUZZ|| ! 541: geog->wlon.l>rwindow[3]+FUZZ) ! 542: return(0); ! 543: else return(1); ! 544: } ! 545: ! 546: int ! 547: inlimits(struct place *g) ! 548: { ! 549: if(rlimits[0]-FUZZ>g->nlat.l|| ! 550: rlimits[1]+FUZZ<g->nlat.l) ! 551: return(0); ! 552: switch(limcase) { ! 553: case 0: ! 554: if(rlimits[2]+TWOPI-FUZZ>g->wlon.l&& ! 555: rlimits[3]+FUZZ<g->wlon.l) ! 556: return(0); ! 557: break; ! 558: case 1: ! 559: if(rlimits[2]-FUZZ>g->wlon.l|| ! 560: rlimits[3]+FUZZ<g->wlon.l) ! 561: return(0); ! 562: break; ! 563: case 2: ! 564: if(rlimits[2]>g->wlon.l&& ! 565: rlimits[3]-TWOPI+FUZZ<g->wlon.l) ! 566: return(0); ! 567: break; ! 568: } ! 569: return(1); ! 570: } ! 571: ! 572: ! 573: long patch[18][36]; ! 574: ! 575: void ! 576: getdata(char *mapfile) ! 577: { ! 578: char *indexfile; ! 579: int kx,ky,c; ! 580: int k; ! 581: long b; ! 582: long *p; ! 583: int ip, jp; ! 584: int n; ! 585: struct place g; ! 586: int i, j; ! 587: double lat, lon; ! 588: int conn; ! 589: FILE *ifile, *xfile; ! 590: ! 591: indexfile = mapindex(mapfile); ! 592: xfile = fopen(indexfile,"r"); ! 593: if(xfile==NULL) ! 594: filerror("can't find map index", indexfile); ! 595: free(indexfile); ! 596: for(i=0,p=patch[0];i<18*36;i++,p++) ! 597: *p = 1; ! 598: while(!feof(xfile) && fscanf(xfile,"%d%d%ld",&i,&j,&b)==3) ! 599: patch[i+9][j+18] = b; ! 600: fclose(xfile); ! 601: ifile = fopen(mapfile,"r"); ! 602: if(ifile==NULL) ! 603: filerror("can't find map data", mapfile); ! 604: for(lat=lolat;lat<hilat;lat+=10.) ! 605: for(lon=lolon;lon<hilon;lon+=10.) { ! 606: if(!seeable(lat,lon)) ! 607: continue; ! 608: i = pnorm(lat); ! 609: j = pnorm(lon); ! 610: if((b=patch[i+9][j+18])&1) ! 611: continue; ! 612: fseek(ifile,b,0); ! 613: while((ip=getc(ifile))>=0&&(jp=getc(ifile))>=0){ ! 614: if(ip!=(i&0377)||jp!=(j&0377)) ! 615: break; ! 616: n = getshort(ifile); ! 617: conn = 0; ! 618: if(n > 0) { /* absolute coordinates */ ! 619: for(k=0;k<n;k++){ ! 620: kx = SCALERATIO*getshort(ifile); ! 621: ky = SCALERATIO*getshort(ifile); ! 622: if (((k%delta) != 0) && (k != (n-1))) ! 623: continue; ! 624: conv(kx,&g.nlat); ! 625: conv(ky,&g.wlon); ! 626: conn = plotpt(&g,conn); ! 627: } ! 628: } else { /* differential, scaled by SCALERATI0 */ ! 629: n = -n; ! 630: kx = SCALERATIO*getshort(ifile); ! 631: ky = SCALERATIO*getshort(ifile); ! 632: for(k=0; k<n; k++) { ! 633: c = getc(ifile); ! 634: if(c&0200) c|= ~0177; ! 635: kx += c; ! 636: c = getc(ifile); ! 637: if(c&0200) c|= ~0177; ! 638: ky += c; ! 639: if(k%delta!=0&&k!=n-1) ! 640: continue; ! 641: conv(kx,&g.nlat); ! 642: conv(ky,&g.wlon); ! 643: conn = plotpt(&g,conn); ! 644: } ! 645: } ! 646: if(k==1) { ! 647: conv(kx,&g.nlat); ! 648: conv(ky,&g.wlon); ! 649: conn = plotpt(&g,conn); ! 650: } ! 651: } ! 652: } ! 653: fclose(ifile); ! 654: } ! 655: ! 656: int ! 657: seeable(double lat0, double lon0) ! 658: { ! 659: double x, y; ! 660: double lat, lon; ! 661: for(lat=lat0;lat<=lat0+10;lat+=2*grid[2]) ! 662: for(lon=lon0;lon<=lon0+10;lon+=2*grid[2]) ! 663: if(normproj(lat,lon,&x,&y)*vflag>0) ! 664: return(1); ! 665: return(0); ! 666: } ! 667: ! 668: void ! 669: satellite(struct file *t) ! 670: { ! 671: char sym[50]; ! 672: char lbl[50]; ! 673: double scale; ! 674: register conn; ! 675: double lat,lon; ! 676: struct place place; ! 677: static FILE *ifile = stdin; ! 678: if(t->name[0]!='-'||t->name[1]!=0) { ! 679: fclose(ifile); ! 680: if((ifile=fopen(t->name,"r"))==NULL) ! 681: filerror("can't find track", t->name); ! 682: } ! 683: color(t->color); ! 684: pen(t->style); ! 685: for(;;) { ! 686: conn = 0; ! 687: while(!feof(ifile) && fscanf(ifile,"%lf%lf",&lat,&lon)==2){ ! 688: latlon(lat,lon,&place); ! 689: if(fscanf(ifile,"%1s",lbl) == 1) { ! 690: if(strchr("+-.0123456789",*lbl)==0) ! 691: break; ! 692: ungetc(*lbl,ifile); ! 693: } ! 694: conn = plotpt(&place,conn); ! 695: } ! 696: if(feof(ifile)) ! 697: return; ! 698: fscanf(ifile,"%[^\n]",lbl+1); ! 699: switch(*lbl) { ! 700: case '"': ! 701: if(plotpt(&place,conn)) ! 702: text(lbl+1); ! 703: break; ! 704: case ':': ! 705: case '!': ! 706: if(sscanf(lbl+1,"%s %lf",sym,&scale) <= 1) ! 707: scale = 1; ! 708: if(plotpt(&place,conn?conn:-1)) { ! 709: int r = *lbl=='!'?0:rflag?-1:1; ! 710: if(putsym(&place,sym,scale,r) == 0) ! 711: text(lbl); ! 712: } ! 713: break; ! 714: default: ! 715: if(plotpt(&place,conn)) ! 716: text(lbl); ! 717: break; ! 718: } ! 719: } ! 720: } ! 721: ! 722: int ! 723: pnorm(double x) ! 724: { ! 725: int i; ! 726: i = x/10.; ! 727: i %= 36; ! 728: if(i>=18) return(i-36); ! 729: if(i<-18) return(i+36); ! 730: return(i); ! 731: } ! 732: ! 733: void ! 734: error(char *s) ! 735: { ! 736: fprintf(stderr,"map: \r\n%s\n",s); ! 737: exit(1); ! 738: } ! 739: ! 740: void ! 741: filerror(char *s, char *f) ! 742: { ! 743: fprintf(stderr,"\r\n%s %s\n",s,f); ! 744: exit(1); ! 745: } ! 746: ! 747: char * ! 748: mapindex(char *s) ! 749: { ! 750: char *t = malloc(strlen(s)+3); ! 751: strcpy(t,s); ! 752: strcat(t,".x"); ! 753: return t; ! 754: } ! 755: ! 756: #define NOPT 32767 ! 757: static ox = NOPT, oy = NOPT; ! 758: ! 759: int ! 760: cpoint(int xi, int yi, int conn) ! 761: { ! 762: int dx = abs(ox-xi); ! 763: int dy = abs(oy-yi); ! 764: if(!xflag && (xi<left||xi>=right || yi<bottom||yi>=top)) { ! 765: ox = oy = NOPT; ! 766: return 0; ! 767: } ! 768: if(conn == -1) /* isolated plotting symbol */ ! 769: ; ! 770: else if(!conn) ! 771: move(xi,yi); ! 772: else { ! 773: if(dx+dy>longlines) { ! 774: ox = oy = NOPT; /* don't leap across cuts */ ! 775: return 0; ! 776: } ! 777: if(dx || dy) ! 778: vec(xi,yi); ! 779: } ! 780: ox = xi, oy = yi; ! 781: return dx+dy<=2? 2: 1; /* 2=very near; see dogrid */ ! 782: } ! 783: ! 784: ! 785: struct place oldg; ! 786: ! 787: int ! 788: plotpt(struct place *g, int conn) ! 789: { ! 790: int kx,ky; ! 791: int ret; ! 792: double cutlon; ! 793: if(!inlimits(g)) { ! 794: return(0); ! 795: } ! 796: normalize(g); ! 797: if(!inwindow(g)) { ! 798: return(0); ! 799: } ! 800: switch((*cut)(g,&oldg,&cutlon)) { ! 801: case 2: ! 802: if(conn) { ! 803: ret = duple(g,cutlon)|duple(g,cutlon); ! 804: oldg = *g; ! 805: return(ret); ! 806: } ! 807: case 0: ! 808: conn = 0; ! 809: default: /* prevent diags about bad return value */ ! 810: case 1: ! 811: oldg = *g; ! 812: ret = doproj(g,&kx,&ky); ! 813: if(ret==0 || !onlimb && ret*vflag<=0) ! 814: return(0); ! 815: ret = cpoint(kx,ky,conn); ! 816: return ret; ! 817: } ! 818: } ! 819: ! 820: int ! 821: doproj(struct place *g, int *kx, int *ky) ! 822: { ! 823: int i; ! 824: double x,y,x1,y1; ! 825: /*fprintf(stderr,"dopr1 %f %f \n",g->nlat.l,g->wlon.l);*/ ! 826: i = fixproj(g,&x,&y); ! 827: if(i == 0) ! 828: return(0); ! 829: if(rflag) ! 830: x = -x; ! 831: /*fprintf(stderr,"dopr2 %f %f\n",x,y);*/ ! 832: if(!inpoly(x,y)) { ! 833: return 0; ! 834: } ! 835: x1 = x - xcent; ! 836: y1 = y - ycent; ! 837: x = (x1*crot.c - y1*crot.s + xoff)*scaling; ! 838: y = (x1*crot.s + y1*crot.c + yoff)*scaling; ! 839: *kx = x + (x>0?.5:-.5); ! 840: *ky = y + (y>0?.5:-.5); ! 841: return(i); ! 842: } ! 843: ! 844: int ! 845: duple(struct place *g, double cutlon) ! 846: { ! 847: int kx,ky; ! 848: int okx,oky; ! 849: struct place ig; ! 850: revlon(g,cutlon); ! 851: revlon(&oldg,cutlon); ! 852: ig = *g; ! 853: invert(&ig); ! 854: if(!inlimits(&ig)) ! 855: return(0); ! 856: if(doproj(g,&kx,&ky)*vflag<=0 || ! 857: doproj(&oldg,&okx,&oky)*vflag<=0) ! 858: return(0); ! 859: cpoint(okx,oky,0); ! 860: cpoint(kx,ky,1); ! 861: return(1); ! 862: } ! 863: ! 864: void ! 865: revlon(struct place *g, double cutlon) ! 866: { ! 867: g->wlon.l = reduce(cutlon-reduce(g->wlon.l-cutlon)); ! 868: sincos(&g->wlon); ! 869: } ! 870: ! 871: ! 872: /* recognize problems of cuts ! 873: * move a point across cut to side of its predecessor ! 874: * if its very close to the cut ! 875: * return(0) if cut interrupts the line ! 876: * return(1) if line is to be drawn normally ! 877: * return(2) if line is so close to cut as to ! 878: * be properly drawn on both sheets ! 879: */ ! 880: ! 881: int ! 882: picut(struct place *g, struct place *og, double *cutlon) ! 883: { ! 884: *cutlon = PI; ! 885: return(ckcut(g,og,PI)); ! 886: } ! 887: ! 888: int ! 889: nocut(struct place *g, struct place *og, double *cutlon) ! 890: { ! 891: #pragma ref g ! 892: #pragma ref og ! 893: #pragma ref cutlon ! 894: return(1); ! 895: } ! 896: ! 897: int ! 898: ckcut(struct place *g1, struct place *g2, double lon) ! 899: { ! 900: double d1, d2; ! 901: double f1, f2; ! 902: int kx,ky; ! 903: d1 = reduce(g1->wlon.l -lon); ! 904: d2 = reduce(g2->wlon.l -lon); ! 905: if((f1=fabs(d1))<FUZZ) ! 906: d1 = diddle(g1,lon,d2); ! 907: if((f2=fabs(d2))<FUZZ) { ! 908: d2 = diddle(g2,lon,d1); ! 909: if(doproj(g2,&kx,&ky)*vflag>0) ! 910: cpoint(kx,ky,0); ! 911: } ! 912: if(f1<FUZZ&&f2<FUZZ) ! 913: return(2); ! 914: if(f1>PI*TWO_THRD||f2>PI*TWO_THRD) ! 915: return(1); ! 916: return(d1*d2>=0); ! 917: } ! 918: ! 919: double ! 920: diddle(struct place *g, double lon, double d) ! 921: { ! 922: double d1; ! 923: d1 = FUZZ/2; ! 924: if(d<0) ! 925: d1 = -d1; ! 926: g->wlon.l = reduce(lon+d1); ! 927: sincos(&g->wlon); ! 928: return(d1); ! 929: } ! 930: ! 931: double ! 932: reduce(double lon) ! 933: { ! 934: if(lon>PI) ! 935: lon -= 2*PI; ! 936: else if(lon<-PI) ! 937: lon += 2*PI; ! 938: return(lon); ! 939: } ! 940: ! 941: ! 942: double tetrapt = 35.26438968; /* atan(1/sqrt(2)) */ ! 943: ! 944: void ! 945: dogrid(double lat0, double lat1, double lon0, double lon1) ! 946: { ! 947: double slat,slon,tlat,tlon; ! 948: register int conn, oconn; ! 949: slat = tlat = slon = tlon = 0; ! 950: if(lat1>lat0) ! 951: slat = tlat = fmin(grid[2],dlat); ! 952: else ! 953: slon = tlon = fmin(grid[2],dlon);; ! 954: conn = oconn = 0; ! 955: while(lat0<=lat1&&lon0<=lon1) { ! 956: conn = gridpt(lat0,lon0,conn); ! 957: if(projection==Xguyou&&slat>0) { ! 958: if(lat0<-45&&lat0+slat>-45) ! 959: conn = gridpt(-45.,lon0,conn); ! 960: else if(lat0<45&&lat0+slat>45) ! 961: conn = gridpt(45.,lon0,conn); ! 962: } else if(projection==Xtetra&&slat>0) { ! 963: if(lat0<-tetrapt&&lat0+slat>-tetrapt) { ! 964: gridpt(-tetrapt-.001,lon0,conn); ! 965: conn = gridpt(-tetrapt+.001,lon0,0); ! 966: } ! 967: else if(lat0<tetrapt&&lat0+slat>tetrapt) { ! 968: gridpt(tetrapt-.001,lon0,conn); ! 969: conn = gridpt(tetrapt+.001,lon0,0); ! 970: } ! 971: } ! 972: if(conn==0 && oconn!=0) { ! 973: if(slat+slon>.05) { ! 974: lat0 -= slat; /* steps too big */ ! 975: lon0 -= slon; /* or near bdry */ ! 976: slat /= 2; ! 977: slon /= 2; ! 978: conn = oconn = gridpt(lat0,lon0,conn); ! 979: } else ! 980: oconn = 0; ! 981: } else { ! 982: if(conn==2) { ! 983: slat = tlat; ! 984: slon = tlon; ! 985: conn = 1; ! 986: } ! 987: oconn = conn; ! 988: } ! 989: lat0 += slat; ! 990: lon0 += slon; ! 991: } ! 992: gridpt(lat1,lon1,conn); ! 993: } ! 994: ! 995: static gridinv; /* nonzero when doing window bounds */ ! 996: ! 997: int ! 998: gridpt(double lat, double lon, int conn) ! 999: { ! 1000: struct place g; ! 1001: /*fprintf(stderr,"%f %f\n",lat,lon);*/ ! 1002: latlon(lat,lon,&g); ! 1003: if(gridinv) ! 1004: invert(&g); ! 1005: return(plotpt(&g,conn)); ! 1006: } ! 1007: ! 1008: /* win=0 ordinary grid lines, win=1 window lines */ ! 1009: ! 1010: void ! 1011: dobounds(double lolat, double hilat, double lolon, double hilon, int win) ! 1012: { ! 1013: gridinv = win; ! 1014: if(lolat>-90 || win && (poles&1)!=0) ! 1015: dogrid(lolat+FUZZ,lolat+FUZZ,lolon,hilon); ! 1016: if(hilat<90 || win && (poles&2)!=0) ! 1017: dogrid(hilat-FUZZ,hilat-FUZZ,lolon,hilon); ! 1018: if(hilon-lolon<360 || win && cut==picut) { ! 1019: dogrid(lolat,hilat,lolon+FUZZ,lolon+FUZZ); ! 1020: dogrid(lolat,hilat,hilon-FUZZ,hilon-FUZZ); ! 1021: } ! 1022: gridinv = 0; ! 1023: } ! 1024: ! 1025: static void ! 1026: dolimb(void) ! 1027: { ! 1028: double lat, lon; ! 1029: double res = fmin(dlat, dlon)/4; ! 1030: int conn = 0; ! 1031: int newconn; ! 1032: if(limb == 0) ! 1033: return; ! 1034: onlimb = gridinv = 1; ! 1035: for(;;) { ! 1036: newconn = (*limb)(&lat, &lon, res); ! 1037: if(newconn == -1) ! 1038: break; ! 1039: conn = gridpt(lat, lon, conn*newconn); ! 1040: } ! 1041: onlimb = gridinv = 0; ! 1042: } ! 1043: ! 1044: ! 1045: void ! 1046: radbds(double *w, double *rw) ! 1047: { ! 1048: int i; ! 1049: for(i=0;i<4;i++) ! 1050: rw[i] = w[i]*RAD; ! 1051: rw[0] -= FUZZ; ! 1052: rw[1] += FUZZ; ! 1053: rw[2] -= FUZZ; ! 1054: rw[3] += FUZZ; ! 1055: } ! 1056: ! 1057: void ! 1058: windlim(void) ! 1059: { ! 1060: double center = orientation[0]; ! 1061: double colat; ! 1062: if(center>90) ! 1063: center = 180 - center; ! 1064: if(center<-90) ! 1065: center = -180 - center; ! 1066: if(fabs(center)>90) ! 1067: error("unreasonable orientation"); ! 1068: colat = 90 - window[0]; ! 1069: if(center-colat>limits[0]) ! 1070: limits[0] = center - colat; ! 1071: if(center+colat<limits[1]) ! 1072: limits[1] = center + colat; ! 1073: } ! 1074: ! 1075: ! 1076: short ! 1077: getshort(FILE *f) ! 1078: { ! 1079: int c, r; ! 1080: c = getc(f); ! 1081: r = (c | getc(f)<<8); ! 1082: if (r&0x8000) ! 1083: r |= ~0xFFFF; /* in case short > 16 bits */ ! 1084: return r; ! 1085: } ! 1086: ! 1087: double ! 1088: fmin(double x, double y) ! 1089: { ! 1090: return(x<y?x:y); ! 1091: } ! 1092: ! 1093: double ! 1094: fmax(double x, double y) ! 1095: { ! 1096: return(x>y?x:y); ! 1097: } ! 1098: ! 1099: void ! 1100: clamp(double *px, double v) ! 1101: { ! 1102: *px = (v<0?fmax:fmin)(*px,v); ! 1103: } ! 1104: ! 1105: void ! 1106: pathnames(void) ! 1107: { ! 1108: int i; ! 1109: char *t, *indexfile, *name; ! 1110: FILE *f, *fx; ! 1111: for(i=0; i<nfile; i++) { ! 1112: name = file[i].name; ! 1113: if(*name=='/') ! 1114: continue; ! 1115: indexfile = mapindex(name); ! 1116: /* ansi equiv of unix access() call */ ! 1117: f = fopen(name, "r"); ! 1118: fx = fopen(indexfile, "r"); ! 1119: if(f) fclose(f); ! 1120: if(fx) fclose(fx); ! 1121: free(indexfile); ! 1122: if(f && fx) ! 1123: continue; ! 1124: t = malloc(strlen(name)+strlen(mapdir)+2); ! 1125: strcpy(t,mapdir); ! 1126: strcat(t,"/"); ! 1127: strcat(t,name); ! 1128: file[i].name = t; ! 1129: } ! 1130: } ! 1131: ! 1132: void ! 1133: clipinit(void) ! 1134: { ! 1135: register i; ! 1136: double s,t; ! 1137: if(nvert<=0) ! 1138: return; ! 1139: for(i=0; i<nvert; i++) { /*convert latlon to xy*/ ! 1140: if(normproj(v[i].x,v[i].y,&v[i].x,&v[i].y)==0) ! 1141: error("invisible clipping vertex"); ! 1142: } ! 1143: if(nvert==2) { /*rectangle with diag specified*/ ! 1144: nvert = 4; ! 1145: v[2] = v[1]; ! 1146: v[1].x=v[0].x, v[1].y=v[2].y, v[3].x=v[2].x, v[3].y=v[0].y; ! 1147: } ! 1148: v[nvert] = v[0]; ! 1149: v[nvert+1] = v[1]; ! 1150: s = 0; ! 1151: for(i=1; i<=nvert; i++) { /*test for convexity*/ ! 1152: t = (v[i-1].x-v[i].x)*(v[i+1].y-v[i].y) - ! 1153: (v[i-1].y-v[i].y)*(v[i+1].x-v[i].x); ! 1154: if(t<-FUZZ && s>=0) s = 1; ! 1155: if(t>FUZZ && s<=0) s = -1; ! 1156: if(-FUZZ<=t&&t<=FUZZ || t*s>0) { ! 1157: s = 0; ! 1158: break; ! 1159: } ! 1160: } ! 1161: if(s==0) ! 1162: error("improper clipping polygon"); ! 1163: for(i=0; i<nvert; i++) { /*edge equation ax+by=c*/ ! 1164: e[i].a = s*(v[i+1].y - v[i].y); ! 1165: e[i].b = s*(v[i].x - v[i+1].x); ! 1166: e[i].c = s*(v[i].x*v[i+1].y - v[i].y*v[i+1].x); ! 1167: } ! 1168: } ! 1169: ! 1170: int ! 1171: inpoly(double x, double y) ! 1172: { ! 1173: register i; ! 1174: for(i=0; i<nvert; i++) { ! 1175: register struct edge *ei = &e[i]; ! 1176: double val = x*ei->a + y*ei->b - ei->c; ! 1177: if(val>10*FUZZ) ! 1178: return(0); ! 1179: } ! 1180: return 1; ! 1181: } ! 1182: ! 1183: void ! 1184: realcut() ! 1185: { ! 1186: struct place g; ! 1187: double lat; ! 1188: ! 1189: if(cut != picut) /* punt on unusual cuts */ ! 1190: return; ! 1191: for(lat=window[0]; lat<=window[1]; lat+=grid[2]) { ! 1192: g.wlon.l = PI; ! 1193: sincos(&g.wlon); ! 1194: g.nlat.l = lat*RAD; ! 1195: sincos(&g.nlat); ! 1196: if(!inwindow(&g)) { ! 1197: break; ! 1198: } ! 1199: invert(&g); ! 1200: if(inlimits(&g)) { ! 1201: return; ! 1202: } ! 1203: } ! 1204: longlines = shortlines = LONGLINES; ! 1205: cut = nocut; /* not necessary; small eff. gain */ ! 1206: }
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.