Return to map.c CVS log

Up to [Research Unix] / researchv10dc / cmd / map

Annotation of researchv10dc/cmd/map/map.c, revision 1.1

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

This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.