Menu
▾
▴
[Meshdb-cvs] src/geo/alt main.c,1.3.2.1,1.3.2.2
|
From: <le...@pr...> - 2005-01-21 14:12:45
|
Update of /cvsroot/meshdb/src/geo/alt In directory sc8-pr-cvs1.sourceforge.net:/tmp/cvs-serv22710 Modified Files: Tag: leonard-dev main.c Log Message: add "-l <height>" flag to do line-of-sight. make texturing greyscale. Index: main.c =================================================================== RCS file: /cvsroot/meshdb/src/geo/alt/main.c,v retrieving revision 1.3.2.1 retrieving revision 1.3.2.2 diff -u -d -r1.3.2.1 -r1.3.2.2 --- main.c 7 Nov 2004 09:49:42 -0000 1.3.2.1 +++ main.c 21 Jan 2005 14:12:37 -0000 1.3.2.2 @@ -25,6 +25,12 @@ #include "compat/compat.h" #include "compat/err.h" +/* +#define image_begin(w, h) printf("P6\n%u %u 255\n",w,h) +#define image_write_rgb(r) fwrite(r,3,1,stdout) +#define image_end() fflush(stdout) +*/ + #define CROSS 2 #define CROSS2 3 @@ -37,6 +43,8 @@ struct cross *cross = NULL; int ncross = 0; +static int debug = 0; + static const float NaN = 0.0f/0.0f; static const float Inf = 1.0f/0.0f; @@ -93,26 +101,52 @@ int tflag = 0; int i; int vflag = 0; + int lflag = 0; + double lvalue = 0; double e0, n0, scale, e, n; int w, h, w2, h2; char *Shiftfile = NULL; char *Altfile = NULL; struct cross *cr; int zone = 56; + float shade; + float *hm, *azm, *azmaxm; + char *visiblem; - while ((ch = getopt(argc, argv, "S:A:to:x:X:F:v")) != -1) + fprintf(stderr, "ALT:"); + for (i = 0; i < argc; i++) fprintf(stderr, " \"%s\"", argv[i]); + fprintf(stderr, "\n"); + + while ((ch = getopt(argc, argv, "A:F:S:dl:tvz:X:o:x:")) != -1) switch (ch) { + case 'A': + Altfile = optarg; + break; + case 'F': + loadfile(optarg); + break; case 'S': Shiftfile = optarg; break; - case 'A': - Altfile = optarg; + case 'd': + debug++; + break; + case 'l': + /* Line-of-sight initial altitude */ + lflag = 1; + sscanf(optarg, "%lf", &lvalue); break; case 't': tflag = 1; break; - case 'o': + case 'v': + vflag = 1; + break; + case 'z': + zone = atoi(optarg); + break; case 'X': + case 'o': case 'x': if (argv[optind] == NULL) { error = 1; @@ -124,15 +158,6 @@ sscanf(argv[optind], "%lf", &cr->n); optind++; break; - case 'F': - loadfile(optarg); - break; - case 'z': - zone = atoi(optarg); - break; - case 'v': - vflag = 1; - break; default: error = 1; } @@ -143,13 +168,17 @@ if (error) { fprintf(stderr, "usage: %s" - " [-S shift.gsb] [-A alt.desc] [-t] [-v]" + " [-dtvl]" + " [-S shift.gsb] [-A alt.desc]" " [-F filename]*" " [-[x|X|o] e n]*" " e n w h [scale]\n", argv[0]); exit(1); } + if (lflag && (vflag || tflag)) + errx(1, "-l incompatible with -v and -t"); + /* Read the rest of the arguments */ sscanf(argv[optind], "%lf", &e0); sscanf(argv[optind+1], "%lf", &n0); @@ -190,8 +219,9 @@ cross[i].kind = 0; /* make invisible */ } - /* generate PPM output by reading individual pixels */ + if (vflag) { + /* Generate VRML output */ printf("#VRML V2.0 utf8\n"); printf("# $Id$\n"); printf("# Brisbane Mesh www.brismesh.org\n" @@ -214,92 +244,181 @@ " height [\n", w2 * 2, h2 * 2, 1.0/scale, 1.0/scale); - } else + } else { + /* generate PNG output */ image_begin(w, h); + } + + if (debug > 1) + { extern int ntv2_verbose, qtdebug, altdebug; + fprintf(stderr, "alt: "); + ntv2_verbose++; + alt_print_at(&alt, e0, n0); + ntv2_verbose--; + qtdebug++; altdebug++; + fprintf(stderr, " = %f", alt_at(&alt, e0, n0)); + qtdebug--; altdebug--; + fprintf(stderr, "\n"); + } + + /* Allocate the height map as (w+2)*(h+2) floats. */ + hm = (float *)malloc((w+2) * (h+2) * sizeof (float)); + + /* Shift hm so that H[-h2-1,-w2-1] == h[0] (step of w+1) */ +#define MAP(m,x,y) (m)[(x) + (w+2)*(y) - (w+2)*(h+2)/2] +#define H(x,y) MAP(hm,x,y) + + /* Compute the height map */ + for (y = -h2 -1; y <= h2; y++) { + n = n0 - y / scale; + for (x = -w2 -1; x <= w2; x++) { + e = e0 + x / scale; + H(x,y) = alt_fastat(&alt, e, n); + } + } + + if (lflag) { + int n, nmax = w2 > h2 ? w2 : h2; + + azm = (float *)malloc((w+2) * (h+2) * sizeof (float)); + azmaxm = (float *)malloc((w+2) * (h+2) * sizeof (float)); + visiblem = (char *)malloc((w+2) * (h+2) * sizeof (char)); +#define AZ(x,y) MAP(azm,x,y) +#define AZMAX(x,y) MAP(azmaxm,x,y) +#define VISIBLE(x,y) MAP(visiblem,x,y) +#define ZENITH (-Inf) + + AZ(0,0) = AZMAX(0,0) = ZENITH; + VISIBLE(0,0) = 1; + + for (n = 1; n < nmax; n++) { + int p; + float d; + +#define LOS(x,y,xp,yp,xq,yq) \ + if (-w2<(x) && (x)<w2 && -h2<(y) && (y)<h2) { \ + float az, azprev, azmaxprev; \ + az = (H(x,y) - lvalue) / d; \ + if (isinf(az)) az = ZENITH; \ + azprev = AZ(xp,yp)*(p/n) + \ + AZ(xq,yq)*(1-p/n); \ + azmaxprev = AZMAX(xp,yp)*(p/n) + \ + AZMAX(xq,yq)*(1-p/n); \ + if (az < azmaxprev) { \ + AZMAX(x,y) = azmaxprev; \ + VISIBLE(x,y) = 0; \ + } else { \ + AZMAX(x,y) = az; \ + VISIBLE(x,y) = 1; \ + } \ + AZ(x,y) = az; \ + } + + for (p = 1; p <= n; p++) { + LOS(n,p,n-1,p-1,n-1,p) /* ESE */ + LOS(-p,n,-(p-1),n-1,-(p-1),n) /* SSW */ + LOS(-n,-p,-(n-1),-(p-1),-(n-1),-p) /* WNW */ + LOS(p,-n,p-1,-(n-1),p-1,-n) /* NNE */ + } + for (p = 0; p < n; p++) { + LOS(n,-p,n-1,-(p-1),n-1,-p) /* ENE */ + LOS(p,n,p-1,n-1,p-1,n) /* SSE */ + LOS(-n,p,-(n-1),p-1,-(n-1),p) /* WSW */ + LOS(-p,-n,-(p-1),-(n-1),-(p-1),-n) /* NNW */ + } + } + } + for (y = -h2; y < h2; y++) { - n = n0 - y / scale; for (x = -w2; x < w2; x++) { int iscross = 0; - float f; - e = e0 + x / scale; - f = alt_fastat(&alt, e, n); + float f = H(x,y); + /* VRML is easy */ if (vflag) { if (isnanf(f)) printf(" 0"); else printf(" %.2f", f); - } else { + continue; + } - for (i = 0; i < ncross && !iscross; i++) { - int dx, dy; - switch (cross[i].kind) { - case 0: - break; + /* Is this pixel on a cross/circle? */ + for (i = 0; i < ncross && !iscross; i++) { + int dx, dy; + switch (cross[i].kind) { + case 0: + break; - case 'X': - dx = x - cross[i].x; - dy = y - cross[i].y; - if ((dx == dy || dx == -dy) && - (-CROSS2 < dx && dx < CROSS2)) - iscross = 1; - break; + case 'X': + dx = x - cross[i].x; + dy = y - cross[i].y; + if ((dx == dy || dx == -dy) && + (-CROSS2 < dx && dx < CROSS2)) + iscross = 1; + break; - case 'x': - dx = x - cross[i].x; - dy = y - cross[i].y; - if ((dx == dy || dx == -dy) && - (-CROSS < dx && dx < CROSS)) - iscross = 1; - break; + case 'x': + dx = x - cross[i].x; + dy = y - cross[i].y; + if ((dx == dy || dx == -dy) && + (-CROSS < dx && dx < CROSS)) + iscross = 1; + break; - case 'o': - dx = x - cross[i].x; - dy = y - cross[i].y; - if (dx*dx+dy*dy < (CROSS2+1)*(CROSS2+1) && - dx*dx+dy*dy >= (CROSS2-0)*(CROSS2-0)) - iscross = 1; - break; - } + case 'o': + dx = x - cross[i].x; + dy = y - cross[i].y; + if (dx*dx+dy*dy < (CROSS2+1)*(CROSS2+1) && + dx*dx+dy*dy >= (CROSS2-0)*(CROSS2-0)) + iscross = 1; + break; } + } + + /* Crosses (and circles in red) */ + if (iscross) { + rgb[0] = 255; + rgb[1] = rgb[2] = 0; + goto plot; + } + + if (isnanf(f)) { + /* Unknown data (sea) is black */ + rgb[0] = rgb[1] = rgb[2] = 0; + goto plot; + } + + /* + * Shade the pixel after determining the + * cosine of the normal to the top-left sun + */ + shade = (f - H(x-1,y-1)) / ((fmax-fmin) * 0.3) + 0.5; + if (shade < 0) shade = 0; + if (shade > 1) shade = 1; + + rgb[0] = rgb[1] = rgb[2] = 255 * shade; + + /* Tinge with blue if visible */ + if (lflag && VISIBLE(x,y)) + rgb[2] = 255 - (255-rgb[2])/4; - if (iscross) - rgb[0] = rgb[1] = rgb[2] = 255; - else { - if (isnanf(f)) - /* Unknown data (sea) is black */ - rgb[0] = rgb[1] = rgb[2] = 0; - else if (f < fmin) { - /* sub-minima is blue */ - rgb[0] = 0; - rgb[1] = 0; - rgb[2] = 128; - } else if (f > fmax) { - /* super-maxima is green */ - rgb[0] = 0; - rgb[1] = 255; - rgb[2] = 0; - } else { - /* everything else is in between */ - float i = (f - fmin)/(fmax - fmin); - rgb[0] = 0; - rgb[1] = i * 255.0; - rgb[2] = (1.0 - i) * (i*0.5+0.5) * 255.0; - } /* fprintf(stderr, "[%d,%d] = %lf,%lf -> %f -> <%u,%u,%u>\n", - x,y,lon,lat,f,rgb[0],rgb[1],rgb[2]); +x,y,lon,lat,f,rgb[0],rgb[1],rgb[2]); */ - } - image_write_rgb(&rgb); - } + + plot: + image_write_rgb(&rgb); } } if (vflag) printf("\n ]\n" " }\n" "}\n"); - else + else { image_end(); + fflush(stdout); + } exit(0); } |
