// (Log is kept at end of this file)
+#ifdef HAVE_CONFIG_H
+# include <config.h>
+#endif
+
#include <ctype.h> // isspace()
+#include <stdlib.h> // atoi()
#include <math.h> // rint()
#include <stdio.h>
-#include <stdlib.h> // atoi()
#include <string.h>
#include <sys/stat.h> // stat()
#include <unistd.h> // stat()
+#include <zlib/zlib.h>
+
#include "dem.hxx"
#include "leastsqs.hxx"
#include <Include/fg_constants.h>
+#define MAX_EX_NODES 10000
+
#ifdef WIN32
# define MKDIR(a) mkdir(a,S_IRWXU) // I am just guessing at this flag (NHV)
#endif // WIN32
fgDEM::fgDEM( void ) {
// printf("class fgDEM CONstructor called.\n");
+ dem_data = new float[DEM_SIZE_1][DEM_SIZE_1];
+ output_data = new float[DEM_SIZE_1][DEM_SIZE_1];
}
// open input file (or read from stdin)
if ( strcmp(file, "-") == 0 ) {
printf("Loading DEM data file: stdin\n");
- fd = stdin;
+ // fd = stdin;
+ fd = gzdopen(STDIN_FILENO, "r");
} else {
- if ( (fd = fopen(file, "r")) == NULL ) {
- printf("Cannot open %s\n", file);
+ if ( (fd = gzopen(file, "rb")) == NULL ) {
+ printf("Cannot gzopen %s\n", file);
return(0);
}
printf("Loading DEM data file: %s\n", file);
// close a DEM file
int fgDEM::close ( void ) {
- fclose(fd);
+ gzclose(fd);
return(1);
}
// return next token from input stream
-static void next_token(FILE *fd, char *token) {
- int result;
-
- result = fscanf(fd, "%s", token);
+static void next_token(gzFile fd, char *token) {
+ int i, result;
+ char c;
+
+ i = 0;
+ c = gzgetc(fd);
+ // skip past spaces
+ while ( (c != -1) && (c == ' ') ) {
+ c = gzgetc(fd);
+ }
+ while ( (c != -1) && (c != ' ') && (c != '\n') ){
+ token[i] = c;
+ i++;
+ c = gzgetc(fd);
+ }
+ token[i] = '\0';
- if ( result == EOF ) {
+ if ( c == -1 ) {
strcpy(token, "__END_OF_FILE__");
printf(" Warning: Reached end of file!\n");
}
// return next integer from input stream
-static int next_int(FILE *fd) {
+static int next_int(gzFile fd) {
char token[80];
next_token(fd, token);
// return next double from input stream
-static double next_double(FILE *fd) {
+static double next_double(gzFile fd) {
char token[80];
next_token(fd, token);
// return next exponential num from input stream
-static int next_exp(FILE *fd) {
+static int next_exp(gzFile fd) {
+ char token[80];
double mantissa;
int exp, acc;
int i;
- fscanf(fd, "%lfD%d", &mantissa, &exp);
+ next_token(fd, token);
+
+ sscanf(token, "%lfD%d", &mantissa, &exp);
// printf(" Mantissa = %.4f Exp = %d\n", mantissa, exp);
// read and parse DEM "A" record
-void fgDEM::read_a_record( void ) {
+int fgDEM::read_a_record( void ) {
int i, inum;
double dnum;
- char name[144];
+ char name[256];
char token[80];
char *ptr;
// get the name field (144 characters)
for ( i = 0; i < 144; i++ ) {
- name[i] = fgetc(fd);
+ name[i] = gzgetc(fd);
}
name[i+1] = '\0';
inum = next_int(fd);
printf(" DEM level code = %d\n", inum);
+ if ( inum > 3 ) {
+ return(0);
+ }
+
// Pattern code, 1 indicates a regular elevation pattern
inum = next_int(fd);
printf(" Pattern code = %d\n", inum);
// number of profiles
dem_num_profiles = cols = next_int(fd);
printf(" Expecting %d profiles\n", dem_num_profiles);
+
+ return(1);
}
// read and parse DEM "B" record
-void fgDEM::read_b_record(float dem_data[DEM_SIZE_1][DEM_SIZE_1])
-{
+void fgDEM::read_b_record( void ) {
char token[80];
int i;
// parse dem file
-int fgDEM::parse( float dem_data[DEM_SIZE_1][DEM_SIZE_1] ) {
+int fgDEM::parse( void ) {
int i;
- cur_row = 0;
+ cur_col = 0;
- read_a_record();
+ if ( !read_a_record() ) {
+ return(0);
+ }
for ( i = 0; i < dem_num_profiles; i++ ) {
- read_b_record( dem_data );
+ // printf("Ready to read next b record\n");
+ read_b_record();
cur_col++;
if ( cur_col % 100 == 0 ) {
printf(" Done parsing\n");
- return(0);
+ return(1);
}
// return the current altitude based on mesh data. We should rewrite
// this to interpolate exact values, but for now this is good enough
-double fgDEM::interpolate_altitude( float dem_data[DEM_SIZE_1][DEM_SIZE_1],
- double lon, double lat)
-{
+double fgDEM::interpolate_altitude( double lon, double lat ) {
// we expect incoming (lon,lat) to be in arcsec for now
double xlocal, ylocal, dx, dy, zA, zB, elev;
// Use least squares to fit a simpler data set to dem data
-void fgDEM::fit( float dem_data[DEM_SIZE_1][DEM_SIZE_1],
- float output_data[DEM_SIZE_1][DEM_SIZE_1],
- char *fg_root, double error, struct fgBUCKET *p )
-{
+void fgDEM::fit( double error, fgBUCKET *p ) {
double x[DEM_SIZE_1], y[DEM_SIZE_1];
double m, b, ave_error, max_error;
double cury, lasty;
// FILE *dem, *fit, *fit1;
printf("Initializing output mesh structure\n");
- outputmesh_init( output_data );
+ outputmesh_init();
// determine dimensions
colmin = p->x * ( (cols - 1) / 8);
printf("Fitting region = %d,%d to %d,%d\n", colmin, rowmin, colmax, rowmax);
// include the corners explicitly
- outputmesh_set_pt(output_data, colmin, rowmin, dem_data[colmin][rowmin]);
- outputmesh_set_pt(output_data, colmin, rowmax, dem_data[colmin][rowmax]);
- outputmesh_set_pt(output_data, colmax, rowmax, dem_data[colmax][rowmax]);
- outputmesh_set_pt(output_data, colmax, rowmin, dem_data[colmax][rowmin]);
+ outputmesh_set_pt(colmin, rowmin, dem_data[colmin][rowmin]);
+ outputmesh_set_pt(colmin, rowmax, dem_data[colmin][rowmax]);
+ outputmesh_set_pt(colmax, rowmax, dem_data[colmax][rowmax]);
+ outputmesh_set_pt(colmax, rowmin, dem_data[colmax][rowmin]);
printf("Beginning best fit procedure\n");
if ( start > colmin ) {
// skip this for the first line segment
cury = m * x[0] + b;
- outputmesh_set_pt(output_data, start, row, (lasty + cury) / 2);
+ outputmesh_set_pt(start, row, (lasty + cury) / 2);
// fprintf(fit, "%.2f %.2f\n", x[0], (lasty + cury) / 2);
}
// printf("Please hit return: "); gets(junk);
}
- outputmesh_output_nodes(output_data, fg_root, p);
+ // outputmesh_output_nodes(fg_root, p);
}
// Initialize output mesh structure
-void fgDEM::outputmesh_init( float output_data[DEM_SIZE_1][DEM_SIZE_1] ) {
+void fgDEM::outputmesh_init( void ) {
int i, j;
for ( j = 0; j < DEM_SIZE_1; j++ ) {
// Get the value of a mesh node
-double fgDEM::outputmesh_get_pt( float output_data[DEM_SIZE_1][DEM_SIZE_1],
- int i, int j )
-{
+double fgDEM::outputmesh_get_pt( int i, int j ) {
return ( output_data[i][j] );
}
// Set the value of a mesh node
-void fgDEM::outputmesh_set_pt( float output_data[DEM_SIZE_1][DEM_SIZE_1],
- int i, int j, double value )
-{
+void fgDEM::outputmesh_set_pt( int i, int j, double value ) {
// printf("Setting data[%d][%d] = %.2f\n", i, j, value);
output_data[i][j] = value;
}
-// Write out a node file that can be used by the "triangle" program
-void fgDEM::outputmesh_output_nodes( float output_data[DEM_SIZE_1][DEM_SIZE_1],
- char *fg_root, struct fgBUCKET *p )
-{
+// Write out a node file that can be used by the "triangle" program.
+// Check for an optional "index.node.ex" file in case there is a .poly
+// file to go along with this node file. Include these nodes first
+// since they are referenced by position from the .poly file.
+void fgDEM::outputmesh_output_nodes( char *fg_root, fgBUCKET *p ) {
+ double exnodes[MAX_EX_NODES][3];
+ double junk1, junk2, junk3;
struct stat stat_buf;
- char base_path[256], dir[256], file[256];
+ char base_path[256], dir[256], file[256], exfile[256];
#ifdef WIN32
char tmp_path[256];
#endif
FILE *fd;
long int index;
int colmin, colmax, rowmin, rowmax;
- int i, j, count, result;
+ int i, j, count, excount, result;
// determine dimensions
colmin = p->x * ( (cols - 1) / 8);
index = fgBucketGenIndex(p);
sprintf(file, "%s/%ld.node", dir, index);
+ // get (optional) extra node file name (in case there is matching
+ // .poly file.
+ strcpy(exfile, file);
+ strcat(exfile, ".ex");
+
+ // load extra nodes if they exist
+ excount = 0;
+ if ( (fd = fopen(exfile, "r")) != NULL ) {
+ fscanf(fd, "%d %d %d %d", &excount, &junk1, &junk2, &junk3);
+
+ if ( excount > MAX_EX_NODES - 1 ) {
+ printf("Error, too many 'extra' nodes, increase array size\n");
+ exit(-1);
+ } else {
+ printf(" Expecting %d 'extra' nodes\n", excount);
+ }
+
+ for ( i = 1; i <= excount; i++ ) {
+ fscanf(fd, "%d %lf %lf %lf\n", &junk1,
+ &exnodes[i][0], &exnodes[i][1], &exnodes[i][2]);
+ printf("(extra) %d %.2f %.2f %.2f\n",
+ i, exnodes[i][0], exnodes[i][1], exnodes[i][2]);
+ }
+ fclose(fd);
+ }
+
printf("Creating node file: %s\n", file);
fd = fopen(file, "w");
- // first count nodes to generate header
+ // first count regular nodes to generate header
count = 0;
for ( j = rowmin; j <= rowmax; j++ ) {
for ( i = colmin; i <= colmax; i++ ) {
}
// printf(" count = %d\n", count);
}
- fprintf(fd, "%d 2 1 0\n", count);
+ fprintf(fd, "%d 2 1 0\n", count + excount);
+
+ // now write out extra node data
+ for ( i = 1; i <= excount; i++ ) {
+ fprintf(fd, "%d %.2f %.2f %.2f\n",
+ i, exnodes[i][0], exnodes[i][1], exnodes[i][2]);
+ }
- // now write out actual node data
- count = 1;
+ // write out actual node data
+ count = excount + 1;
for ( j = rowmin; j <= rowmax; j++ ) {
for ( i = colmin; i <= colmax; i++ ) {
if ( output_data[i][j] > -9000.0 ) {
fgDEM::~fgDEM( void ) {
// printf("class fgDEM DEstructor called.\n");
+ delete(dem_data);
+ delete(output_data);
}
// $Log$
+// Revision 1.13 1998/09/09 16:24:04 curt
+// Fixed a bug in the handling of exclude files which was causing
+// a crash by calling fclose() on an invalid file handle.
+//
+// Revision 1.12 1998/08/24 20:03:31 curt
+// Eliminated a possible memory overrun error.
+// Use the proper free() rather than the incorrect delete().
+//
+// Revision 1.11 1998/07/20 12:46:11 curt
+// When outputing to a .node file, first check for an optional
+// "index.node.ex" file in case there is a .poly file to go along with this
+// node file. Include these nodes first since they are referenced by position
+// from the .poly file. This is my first pass at adding an area "cutout"
+// feature to the terrain generation pipeline.
+//
+// Revision 1.10 1998/07/13 20:58:02 curt
+// .
+//
+// Revision 1.9 1998/07/13 15:29:49 curt
+// Added #ifdef HAVE_CONFIG_H
+//
+// Revision 1.8 1998/07/04 00:47:18 curt
+// typedef'd struct fgBUCKET.
+//
+// Revision 1.7 1998/06/05 18:14:39 curt
+// Abort out early when reading the "A" record if it doesn't look like
+// a proper DEM file.
+//
+// Revision 1.6 1998/05/02 01:49:21 curt
+// Fixed a bug where the wrong variable was being initialized.
+//
+// Revision 1.5 1998/04/25 15:00:32 curt
+// Changed "r" to "rb" in gzopen() options. This fixes bad behavior in win32.
+//
+// Revision 1.4 1998/04/22 13:14:46 curt
+// Fixed a bug in zlib usage.
+//
+// Revision 1.3 1998/04/18 03:53:05 curt
+// Added zlib support.
+//
+// Revision 1.2 1998/04/14 02:43:27 curt
+// Used "new" to auto-allocate large DEM parsing arrays in class constructor.
+//
// Revision 1.1 1998/04/08 22:57:22 curt
// Adopted Gnu automake/autoconf system.
//