Annotation of researchv10no/cmd/view2d/vdata.c, revision 1.1.1.1

1.1       root        1: #include       <stdio.h>
                      2: #include       <math.h>
                      3: #include       "view2d.h"
                      4: #define MAXPTS 4096
                      5: 
                      6: int NX = 16, NY = 16;
                      7: char *progname;
                      8: 
                      9: typedef struct Vector
                     10: {
                     11:        float x, y, z;
                     12: } Vector;
                     13: 
                     14: typedef struct Frame
                     15: {
                     16:        Vector data[MAXPTS];
                     17:        double time;
                     18:         int npts;
                     19: } Frame;
                     20: Frame frame[40];
                     21: int interp = 1;
                     22: 
                     23: float grid[200][200];
                     24: float xmin, xmax, ymin, ymax;
                     25: float xstep, ystep;
                     26: short *data;
                     27: int verbose;
                     28: 
                     29: main(argc, argv)
                     30:        char **argv;
                     31: {
                     32:        register i, j;
                     33:        int nframes;
                     34:        float x;
                     35: 
                     36:        verbose = 0;
                     37:        progname = *argv;
                     38:        for(argc--, argv++; **argv == '-'; argv++)
                     39:                switch(argv[0][1])
                     40:                {
                     41:                case 'n':
                     42:                        NY = -1;
                     43:                        sscanf(&argv[0][2], "%d,%d", &NX, &NY);
                     44:                        if(NY < 0) NY = NX;
                     45:                        break;
                     46:                case 'c':   /* piecewise constant */
                     47:                        interp = 2;
                     48:                        break;
                     49:                case 's':   /* scatterplot on black background */
                     50:                        interp = 3;
                     51:                        break;
                     52:                case 'i':
                     53:                        interp = 0;
                     54:                        break;
                     55:                case 'v':
                     56:                        verbose = 1;
                     57:                        break;
                     58:                default:
                     59:                        error("usage: vdata -c -s -i -v -nNX,NY file");
                     60:                }
                     61:        if(*argv)
                     62:                if(freopen(*argv, "r", stdin)==NULL)
                     63:                        error("freopen(%s) failed",*argv);
                     64:        if((NX>200)||(NY>200)) error("NX, NY must be <= 200");
                     65:        data = (short *)Malloc(NX*NY*sizeof(short));
                     66: 
                     67:        if(interp)
                     68:        {
                     69:                xmax = ymax = -HUGE;
                     70:                xmin = ymin = HUGE;
                     71:        }
                     72:        for(nframes = 0; readf(&frame[nframes]); nframes++){
                     73:          if(nframes>38) error("too many frames");
                     74:         }
                     75:        if(nframes<=0) error("empty input");
                     76:        if(verbose)fprintf(stderr, "nx=%d ny=%d nframes=%d\n", NX, NY, nframes);
                     77: 
                     78:        if(interp)
                     79:        {
                     80:                xstep = (xmax-xmin)/(NX-1);
                     81:                ystep = (ymax-ymin)/(NY-1);
                     82:                if(verbose)fprintf(stderr, "%.5g<=x<=%.5g  %.5g<=y<=%.5g\n",
                     83:                        xmin, xmax, ymin, ymax);
                     84:        }
                     85:        else
                     86:        {
                     87:                xmin = 0; xstep = 1; xmax = xmin + (NX-1);
                     88:                ymin = 0; ystep = 1; ymax = ymin + (NY-1);
                     89:        }
                     90:        for(i = 0; i < nframes; i++){
                     91:                fit(&frame[i]);
                     92:                out(frame[i].time,frame[i].npts);
                     93:        }
                     94:        exit(0);
                     95: }
                     96: 
                     97: readf(s)
                     98:        Frame *s;
                     99: {
                    100:        register Vector *v;
                    101:        register i;
                    102:        float f;
                    103: 
                    104:        if(scanf("%d %E", &s->npts, &s->time) != 2) return(0);
                    105:        if(verbose) fprintf(stderr,"%d %g\n",s->npts,s->time);
                    106:        if(s->npts>MAXPTS) error("too many data points");
                    107:        if(interp)
                    108:                for(v = s->data, i = 0; i < s->npts; i++, v++)
                    109:                {
                    110:                        scanf("%f %f %f", &v->x, &v->y, &v->z);
                    111:                        if(v->x < xmin) xmin = v->x;
                    112:                        if(v->x > xmax) xmax = v->x;
                    113:                        if(v->y < ymin) ymin = v->y;
                    114:                        if(v->y > ymax) ymax = v->y;
                    115:                }
                    116:        else
                    117:                for(v = s->data, i = 0; i < s->npts; i++, v++)
                    118:                {
                    119:                        scanf("%f", &v->z);
                    120:                }
                    121:        return(1);
                    122: }
                    123: 
                    124: fit(s)
                    125:        Frame *s;
                    126: {
                    127:        register Vector *v;
                    128:        Vector h;
                    129:        float g, x, y;
                    130:        register i, j, k;
                    131:        Vector *c1, *c2, *c3;
                    132:        float d, d1, d2, d3;
                    133: 
                    134:        if( (interp==1) || (interp==2) ) {
                    135:                for(j = 0; j < NY; j++)
                    136:                        for(i = 0; i < NX; i++) {
                    137:                                x = xmin + i*xstep;
                    138:                                y = ymin + j*ystep;
                    139:                                g = grid[i][j];
                    140:                                d1 = d2 = d3 = HUGE;
                    141:                                for(v = s->data, k = 0; k < s->npts; k++, v++) {
                    142:                                        d = (v->x - x)*(v->x - x) +
                    143:                                                (v->y - y)*(v->y - y);
                    144:                                        if(d < d1) {
                    145:                                                c3 = c2; d3 = d2;
                    146:                                                c2 = c1; d2 = d1;
                    147:                                                c1 = v; d1 = d;
                    148:                                        } else if(d < d2) {
                    149:                                                c3 = c2; d3 = d2;
                    150:                                                c2 = v; d2 = d;
                    151:                                        } else if(d < d3) {
                    152:                                                c3 = v; d3 = d;
                    153:                                        }
                    154:                                }
                    155:                                /* c1 c2 c3 three closest points */
                    156:                                if(interp==1) {
                    157:                                        h.x = x; h.y = y; h.z = g;
                    158:                                        interp_(&h, c1, c2, c3);
                    159:                                        grid[i][j] = h.z;
                    160:                                } else if(interp==2)
                    161:                                        grid[i][j] = c1->z;
                    162:                        }
                    163:        } else if( interp==3 ) {
                    164:                for(j = 0; j < NY; j++)
                    165:                        for(i = 0; i < NX; i++)
                    166:                                grid[i][j] = -1e26;
                    167:                for(v = s->data, k = 0; k < s->npts; k++, v++) {
                    168:                        i = (NX-1)*(v->x-xmin)/(xmax-xmin);
                    169:                        j = (NY-1)*(v->y-ymin)/(ymax-ymin);
                    170:                        if((i<0)||(i>=NX)||(j<0)||(j>=NY))
                    171:                                error("can't: i,j=%d,%d\n",i,j);
                    172:                        grid[i][j] = v->z;
                    173:                }
                    174:        } else {  /* option -i */
                    175:                for(v = s->data, j = 0; j < NY; j++)
                    176:                        for(i = 0; i < NX; i++, v++)
                    177:                                grid[i][j] = v->z;
                    178:        }
                    179: }
                    180: 
                    181: out(tim,npt)
                    182:   double tim;
                    183:   int npt;
                    184: {
                    185:        register short *d, i, j;
                    186:        register float zmax, zmin, zrange, f, v2;
                    187:        short u, v;
                    188:        extern double pow2();
                    189: 
                    190:        zmin = HUGE;
                    191:        zmax = -HUGE;
                    192:        for(j = 0; j < NY; j++)
                    193:                for(i = 0; i < NX; i++){
                    194:                        if(grid[i][j] < -1e25 ) continue;
                    195:                        if(grid[i][j] < zmin) zmin = grid[i][j];
                    196:                        if(grid[i][j] > zmax) zmax = grid[i][j];
                    197:                }
                    198: 
                    199:        if(verbose) fprintf(stderr, "time=%g npts=%d %g<=z<=%g\n",
                    200:                 tim, npt, zmin, zmax);
                    201:        zmax = (zmax<0)?0:zmax;
                    202:        zmin = (zmin>0)?0:zmin;
                    203:        zrange = zmax - zmin;
                    204:        if(zrange<=0)
                    205:                v=0, u=0, v2=1, f=0;
                    206:        else{
                    207:                f = 2*BIG-2;
                    208:                v = log(zrange)/log(2.)-14;
                    209:                v2 = pow2(v);
                    210:                while( 2*zrange < f*v2 ){ v--; v2 /= 2; }
                    211:                f = -(zmin+zmax)/(2*v2);
                    212:                u = f;
                    213:                if( (f<-BIG)||(f>BIG) ){
                    214:                        if(verbose)fprintf(stderr,"zmin=%g zmax=%g u=%d v=%d\n",
                    215:                                zmin,zmax,u,v);
                    216:                        error("u out of bounds f=%g",f);
                    217:                }
                    218:        }
                    219:        d = data;
                    220:        for(j = 0; j < NY; j++)
                    221:                for(i = 0; i < NX; i++)
                    222:                        if(grid[i][j] < -1e25){ *d++ = -BIG-1; }
                    223:                         else{ *d++ = f + grid[i][j]/v2; }
                    224:        view2d(1, NX, NY, tim, u, v, 0, 0, 0, data);
                    225: }

unix.superglobalmegacorp.com

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