Clean up various accumulated cruft.
This commit is contained in:
curt
parent
1d71909929
commit
f048b73a9a
2 changed files with 4 additions and 439 deletions
|
|
@ -61,21 +61,7 @@ SG_USING_STD(cout);
|
|||
SG_USING_STD(endl);
|
||||
|
||||
|
||||
#define MAX_EX_NODES 10000
|
||||
|
||||
#if 0
|
||||
#ifdef WIN32
|
||||
# ifdef __BORLANDC__
|
||||
# include <dir.h>
|
||||
# define MKDIR(a) mkdir(a)
|
||||
# else
|
||||
# define MKDIR(a) mkdir(a,S_IRWXU) // I am just guessing at this flag (NHV)
|
||||
# endif // __BORLANDC__
|
||||
#endif // WIN32
|
||||
#endif //0
|
||||
|
||||
|
||||
FGDem::FGDem( void ) :
|
||||
FGDem::FGDem() :
|
||||
z_units(2) // meters
|
||||
{
|
||||
// cout << "class FGDem CONstructor called." << endl;
|
||||
|
|
@ -476,406 +462,7 @@ FGDem::write_area( const string& root, SGBucket& b, bool compress ) {
|
|||
}
|
||||
|
||||
|
||||
#if 0
|
||||
|
||||
// return the current altitude based on grid data. We should rewrite
|
||||
// this to interpolate exact values, but for now this is good enough
|
||||
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;
|
||||
int x1, x2, x3, y1, y2, y3;
|
||||
float z1, z2, z3;
|
||||
int xindex, yindex;
|
||||
|
||||
/* determine if we are in the lower triangle or the upper triangle
|
||||
______
|
||||
| /|
|
||||
| / |
|
||||
| / |
|
||||
|/ |
|
||||
------
|
||||
|
||||
then calculate our end points
|
||||
*/
|
||||
|
||||
xlocal = (lon - originx) / col_step;
|
||||
ylocal = (lat - originy) / row_step;
|
||||
|
||||
xindex = (int)(xlocal);
|
||||
yindex = (int)(ylocal);
|
||||
|
||||
// printf("xindex = %d yindex = %d\n", xindex, yindex);
|
||||
|
||||
if ( xindex + 1 == cols ) {
|
||||
xindex--;
|
||||
}
|
||||
|
||||
if ( yindex + 1 == rows ) {
|
||||
yindex--;
|
||||
}
|
||||
|
||||
if ( (xindex < 0) || (xindex + 1 >= cols) ||
|
||||
(yindex < 0) || (yindex + 1 >= rows) ) {
|
||||
return(-9999);
|
||||
}
|
||||
|
||||
dx = xlocal - xindex;
|
||||
dy = ylocal - yindex;
|
||||
|
||||
if ( dx > dy ) {
|
||||
// lower triangle
|
||||
// printf(" Lower triangle\n");
|
||||
|
||||
x1 = xindex;
|
||||
y1 = yindex;
|
||||
z1 = dem_data[x1][y1];
|
||||
|
||||
x2 = xindex + 1;
|
||||
y2 = yindex;
|
||||
z2 = dem_data[x2][y2];
|
||||
|
||||
x3 = xindex + 1;
|
||||
y3 = yindex + 1;
|
||||
z3 = dem_data[x3][y3];
|
||||
|
||||
// printf(" dx = %.2f dy = %.2f\n", dx, dy);
|
||||
// printf(" (x1,y1,z1) = (%d,%d,%d)\n", x1, y1, z1);
|
||||
// printf(" (x2,y2,z2) = (%d,%d,%d)\n", x2, y2, z2);
|
||||
// printf(" (x3,y3,z3) = (%d,%d,%d)\n", x3, y3, z3);
|
||||
|
||||
zA = dx * (z2 - z1) + z1;
|
||||
zB = dx * (z3 - z1) + z1;
|
||||
|
||||
// printf(" zA = %.2f zB = %.2f\n", zA, zB);
|
||||
|
||||
if ( dx > SG_EPSILON ) {
|
||||
elev = dy * (zB - zA) / dx + zA;
|
||||
} else {
|
||||
elev = zA;
|
||||
}
|
||||
} else {
|
||||
// upper triangle
|
||||
// printf(" Upper triangle\n");
|
||||
|
||||
x1 = xindex;
|
||||
y1 = yindex;
|
||||
z1 = dem_data[x1][y1];
|
||||
|
||||
x2 = xindex;
|
||||
y2 = yindex + 1;
|
||||
z2 = dem_data[x2][y2];
|
||||
|
||||
x3 = xindex + 1;
|
||||
y3 = yindex + 1;
|
||||
z3 = dem_data[x3][y3];
|
||||
|
||||
// printf(" dx = %.2f dy = %.2f\n", dx, dy);
|
||||
// printf(" (x1,y1,z1) = (%d,%d,%d)\n", x1, y1, z1);
|
||||
// printf(" (x2,y2,z2) = (%d,%d,%d)\n", x2, y2, z2);
|
||||
// printf(" (x3,y3,z3) = (%d,%d,%d)\n", x3, y3, z3);
|
||||
|
||||
zA = dy * (z2 - z1) + z1;
|
||||
zB = dy * (z3 - z1) + z1;
|
||||
|
||||
// printf(" zA = %.2f zB = %.2f\n", zA, zB );
|
||||
// printf(" xB - xA = %.2f\n", col_step * dy / row_step);
|
||||
|
||||
if ( dy > SG_EPSILON ) {
|
||||
elev = dx * (zB - zA) / dy + zA;
|
||||
} else {
|
||||
elev = zA;
|
||||
}
|
||||
}
|
||||
|
||||
return(elev);
|
||||
}
|
||||
|
||||
|
||||
// Use least squares to fit a simpler data set to dem data
|
||||
void FGDem::fit( double error, SGBucket& p ) {
|
||||
double x[DEM_SIZE_1], y[DEM_SIZE_1];
|
||||
double m, b, ave_error, max_error;
|
||||
double cury, lasty;
|
||||
int n, row, start, end;
|
||||
int colmin, colmax, rowmin, rowmax;
|
||||
bool good_fit;
|
||||
// FILE *dem, *fit, *fit1;
|
||||
|
||||
printf("Initializing output mesh structure\n");
|
||||
outputmesh_init();
|
||||
|
||||
// determine dimensions
|
||||
colmin = p.get_x() * ( (cols - 1) / 8);
|
||||
colmax = colmin + ( (cols - 1) / 8);
|
||||
rowmin = p.get_y() * ( (rows - 1) / 8);
|
||||
rowmax = rowmin + ( (rows - 1) / 8);
|
||||
printf("Fitting region = %d,%d to %d,%d\n", colmin, rowmin, colmax, rowmax);
|
||||
|
||||
// include the corners explicitly
|
||||
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");
|
||||
|
||||
for ( row = rowmin; row <= rowmax; row++ ) {
|
||||
// fit = fopen("fit.dat", "w");
|
||||
// fit1 = fopen("fit1.dat", "w");
|
||||
|
||||
start = colmin;
|
||||
|
||||
// printf(" fitting row = %d\n", row);
|
||||
|
||||
while ( start < colmax ) {
|
||||
end = start + 1;
|
||||
good_fit = true;
|
||||
|
||||
x[(end - start) - 1] = 0.0 + ( start * col_step );
|
||||
y[(end - start) - 1] = dem_data[start][row];
|
||||
|
||||
while ( (end <= colmax) && good_fit ) {
|
||||
n = (end - start) + 1;
|
||||
// printf("Least square of first %d points\n", n);
|
||||
x[end - start] = 0.0 + ( end * col_step );
|
||||
y[end - start] = dem_data[end][row];
|
||||
least_squares(x, y, n, &m, &b);
|
||||
ave_error = least_squares_error(x, y, n, m, b);
|
||||
max_error = least_squares_max_error(x, y, n, m, b);
|
||||
|
||||
/*
|
||||
printf("%d - %d ave error = %.2f max error = %.2f y = %.2f*x + %.2f\n",
|
||||
start, end, ave_error, max_error, m, b);
|
||||
|
||||
f = fopen("gnuplot.dat", "w");
|
||||
for ( j = 0; j <= end; j++) {
|
||||
fprintf(f, "%.2f %.2f\n", 0.0 + ( j * col_step ),
|
||||
dem_data[row][j]);
|
||||
}
|
||||
for ( j = start; j <= end; j++) {
|
||||
fprintf(f, "%.2f %.2f\n", 0.0 + ( j * col_step ),
|
||||
dem_data[row][j]);
|
||||
}
|
||||
fclose(f);
|
||||
|
||||
printf("Please hit return: "); gets(junk);
|
||||
*/
|
||||
|
||||
if ( max_error > error ) {
|
||||
good_fit = false;
|
||||
}
|
||||
|
||||
end++;
|
||||
}
|
||||
|
||||
if ( !good_fit ) {
|
||||
// error exceeded the threshold, back up
|
||||
end -= 2; // back "end" up to the last good enough fit
|
||||
n--; // back "n" up appropriately too
|
||||
} else {
|
||||
// we popped out of the above loop while still within
|
||||
// the error threshold, so we must be at the end of
|
||||
// the data set
|
||||
end--;
|
||||
}
|
||||
|
||||
least_squares(x, y, n, &m, &b);
|
||||
ave_error = least_squares_error(x, y, n, m, b);
|
||||
max_error = least_squares_max_error(x, y, n, m, b);
|
||||
|
||||
/*
|
||||
printf("\n");
|
||||
printf("%d - %d ave error = %.2f max error = %.2f y = %.2f*x + %.2f\n",
|
||||
start, end, ave_error, max_error, m, b);
|
||||
printf("\n");
|
||||
|
||||
fprintf(fit1, "%.2f %.2f\n", x[0], m * x[0] + b);
|
||||
fprintf(fit1, "%.2f %.2f\n", x[end-start], m * x[end-start] + b);
|
||||
*/
|
||||
|
||||
if ( start > colmin ) {
|
||||
// skip this for the first line segment
|
||||
cury = m * x[0] + b;
|
||||
outputmesh_set_pt(start, row, (lasty + cury) / 2);
|
||||
// fprintf(fit, "%.2f %.2f\n", x[0], (lasty + cury) / 2);
|
||||
}
|
||||
|
||||
lasty = m * x[end-start] + b;
|
||||
start = end;
|
||||
}
|
||||
|
||||
/*
|
||||
fclose(fit);
|
||||
fclose(fit1);
|
||||
|
||||
dem = fopen("gnuplot.dat", "w");
|
||||
for ( j = 0; j < DEM_SIZE_1; j++) {
|
||||
fprintf(dem, "%.2f %.2f\n", 0.0 + ( j * col_step ),
|
||||
dem_data[j][row]);
|
||||
}
|
||||
fclose(dem);
|
||||
*/
|
||||
|
||||
// NOTICE, this is for testing only. This instance of
|
||||
// output_nodes should be removed. It should be called only
|
||||
// once at the end once all the nodes have been generated.
|
||||
// newmesh_output_nodes(&nm, "mesh.node");
|
||||
// printf("Please hit return: "); gets(junk);
|
||||
}
|
||||
|
||||
// outputmesh_output_nodes(fg_root, p);
|
||||
}
|
||||
|
||||
|
||||
// Initialize output mesh structure
|
||||
void FGDem::outputmesh_init( void ) {
|
||||
int i, j;
|
||||
|
||||
for ( j = 0; j < DEM_SIZE_1; j++ ) {
|
||||
for ( i = 0; i < DEM_SIZE_1; i++ ) {
|
||||
output_data[i][j] = -9999.0;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
// Get the value of a mesh node
|
||||
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( 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.
|
||||
// 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( const string& fg_root, SGBucket& p )
|
||||
{
|
||||
double exnodes[MAX_EX_NODES][3];
|
||||
struct stat stat_buf;
|
||||
string dir;
|
||||
char file[256], exfile[256];
|
||||
#ifdef WIN32
|
||||
char tmp_path[256];
|
||||
#endif
|
||||
string command;
|
||||
FILE *fd;
|
||||
long int index;
|
||||
int colmin, colmax, rowmin, rowmax;
|
||||
int i, j, count, excount, result;
|
||||
|
||||
// determine dimensions
|
||||
colmin = p.get_x() * ( (cols - 1) / 8);
|
||||
colmax = colmin + ( (cols - 1) / 8);
|
||||
rowmin = p.get_y() * ( (rows - 1) / 8);
|
||||
rowmax = rowmin + ( (rows - 1) / 8);
|
||||
cout << " dumping region = " << colmin << "," << rowmin << " to " <<
|
||||
colmax << "," << rowmax << "\n";
|
||||
|
||||
// generate the base directory
|
||||
string base_path = p.gen_base_path();
|
||||
cout << "fg_root = " << fg_root << " Base Path = " << base_path << endl;
|
||||
dir = fg_root + "/" + base_path;
|
||||
cout << "Dir = " << dir << endl;
|
||||
|
||||
// stat() directory and create if needed
|
||||
errno = 0;
|
||||
result = stat(dir.c_str(), &stat_buf);
|
||||
if ( result != 0 && errno == ENOENT ) {
|
||||
cout << "Creating directory\n";
|
||||
|
||||
#ifdef _MSC_VER
|
||||
fg_mkdir( dir.c_str() );
|
||||
#else
|
||||
command = "mkdir -p " + dir + "\n";
|
||||
system( command.c_str() );
|
||||
#endif
|
||||
|
||||
} else {
|
||||
// assume directory exists
|
||||
}
|
||||
|
||||
// get index and generate output file name
|
||||
index = p.gen_index();
|
||||
sprintf(file, "%s/%ld.node", dir.c_str(), 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 ) {
|
||||
int junki;
|
||||
fscanf(fd, "%d %d %d %d", &excount, &junki, &junki, &junki);
|
||||
|
||||
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", &junki,
|
||||
&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 regular nodes to generate header
|
||||
count = 0;
|
||||
for ( j = rowmin; j <= rowmax; j++ ) {
|
||||
for ( i = colmin; i <= colmax; i++ ) {
|
||||
if ( output_data[i][j] > -9000.0 ) {
|
||||
count++;
|
||||
}
|
||||
}
|
||||
// printf(" count = %d\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]);
|
||||
}
|
||||
|
||||
// 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 ) {
|
||||
fprintf(fd, "%d %.2f %.2f %.2f\n",
|
||||
count++,
|
||||
originx + (double)i * col_step,
|
||||
originy + (double)j * row_step,
|
||||
output_data[i][j]);
|
||||
}
|
||||
}
|
||||
// printf(" count = %d\n", count);
|
||||
}
|
||||
|
||||
fclose(fd);
|
||||
}
|
||||
#endif
|
||||
|
||||
|
||||
FGDem::~FGDem( void ) {
|
||||
FGDem::~FGDem() {
|
||||
// printf("class FGDem DEstructor called.\n");
|
||||
delete [] dem_data;
|
||||
delete [] output_data;
|
||||
|
|
|
|||
|
|
@ -96,11 +96,11 @@ private:
|
|||
public:
|
||||
|
||||
// Constructor
|
||||
FGDem( void );
|
||||
FGDem();
|
||||
FGDem( const string& file );
|
||||
|
||||
// Destructor
|
||||
~FGDem( void );
|
||||
~FGDem();
|
||||
|
||||
// open a DEM file (use "-" if input is coming from stdin)
|
||||
int open ( const string& file );
|
||||
|
|
@ -122,28 +122,6 @@ public:
|
|||
// hand corner.
|
||||
int write_area( const string& root, SGBucket& b, bool compress );
|
||||
|
||||
#if 0
|
||||
// return the current altitude based on grid data. We should
|
||||
// rewrite this to interpolate exact values, but for now this is
|
||||
// good enough
|
||||
double interpolate_altitude( double lon, double lat );
|
||||
|
||||
// Use least squares to fit a simpler data set to dem data
|
||||
void fit( double error, SGBucket& p );
|
||||
|
||||
// Initialize output mesh structure
|
||||
void outputmesh_init( void );
|
||||
|
||||
// Get the value of a mesh node
|
||||
double outputmesh_get_pt( int i, int j );
|
||||
|
||||
// Set the value of a mesh node
|
||||
void outputmesh_set_pt( int i, int j, double value );
|
||||
|
||||
// Write out a node file that can be used by the "triangle" program
|
||||
void outputmesh_output_nodes( const string& fg_root, SGBucket& p );
|
||||
#endif
|
||||
|
||||
// Informational methods
|
||||
inline double get_originx() const { return originx; }
|
||||
inline double get_originy() const { return originy; }
|
||||
|
|
|
|||
Loading…
Add table
Reference in a new issue