libtin/src/tool/grd2tin.cpp

274 lines
9.3 KiB
C++

#include "../lib/tin.h"
#include "display_help.h"
#include "iostream"
#include "fstream"
#include "sstream"
#include "string"
#include "cmath"
#include "getopt.h"
void display_help(std::string exe_name)
{
dispHelp dh;
dh.changeLayerOut(0, 10);
dh.addHeadInfo(exe_name, "1.0", "Generating the Triangular Irregular Networks (TIN) from regular DEM grids.",
"Yi Zhang (zhangyiss@icloud.com)");
dh.addUsage(exe_name+" -f<grid-file> -r<xmin>/<xmax>/<ymin>/<ymax> -i<dx>/<dy> -m<mesh-file> [-p<poly-file>] [-t<threshold>] [-l<log-file>] [-n<neighbor-file>] [-z]");
dh.addOption("Input DEM grid file.", "-f", "--grid-file");
dh.addOption("Input DEM grid's range.", "-r", "--range");
dh.addOption("Input DEM grid's interval.", "-i", "--interval");
dh.addOption("Output Gmsh file (.msh) of the generated TIN.", "-m", "--mesh-file");
dh.addOption("Input text file of a polygon to control the outline shape of the generated TIN.", "-p", "--polygon-file");
dh.addOption("Threshold of the maximal error of the generated TIN with respect to the input DEM grid. (Default is 1.0)", "-t", "--threshold");
dh.addOption("Output text file of a log file of the maximal error of the generated TIN.", "-l", "--log-file");
dh.addOption("Output text file of neighborships of the generated TIN.", "-n", "--neighbor-file");
dh.addOption("The input DEM grid is a 1-column z-table.", "-z", "--z-table");
dh.addOption("Display help information.", "-h", "--help");
dh.show();
return;
}
int main(int argc, char* argv[])
{
try
{
static struct option long_opts[] =
{
{"grid-file", required_argument, NULL, 'f'},
{"range", required_argument, NULL, 'r'},
{"interval", required_argument, NULL, 'i'},
{"mesh-file", required_argument, NULL, 'm'},
{"polygon-file", required_argument, NULL, 'p'},
{"threshold", required_argument, NULL, 't'},
{"log-file", required_argument, NULL, 'l'},
{"neighbor-file", required_argument, NULL, 'n'},
{"z-table", no_argument, NULL, 'z'},
{"help", no_argument, NULL, 'h'},
{0, 0, 0, 0},
};
if (argc == 1)
{
display_help(argv[0]);
return 0;
}
double range[4] = {0, 1000, 0, 1000};
double interval[2] = {10, 10};
double threshold = 1.0;
std::string grid_file, mesh_file, poly_file, log_file, neigh_file;
grid_file = mesh_file = "NULL";
poly_file = log_file = neigh_file = "NULL";
bool z_table = false;
int curr;
while (1)
{
int optIndex = 0;
curr = getopt_long(argc, argv, "zhf:r:i:m:p:t:l:n:", long_opts, &optIndex);
if (curr == -1) break;
switch (curr)
{
case 'h':
display_help(argv[0]);
return 0;
case 'z':
z_table = true;
break;
case 'f':
grid_file = optarg; break;
case 'r':
if (4 != sscanf(optarg, "%lf/%lf/%lf/%lf", &range[0], &range[1], &range[2], &range[3]))
{
throw "Invalid range parameters";
}
break;
case 'i':
if (2 != sscanf(optarg, "%lf/%lf", &interval[0], &interval[1])) //格式化读入参数
{
throw "Invalid interval parameters";
}
break;
case 'm':
mesh_file = optarg; break;
case 'p':
poly_file = optarg; break;
case 't':
if (1 != sscanf(optarg,"%lf", &threshold)) //格式化读入参数
{
throw "Invalid threshold";
}
break;
case 'l':
log_file = optarg; break;
case 'n':
neigh_file = optarg; break;
case '?':
display_help(argv[0]);
return 0;
default:
abort();
}
}
if (grid_file == "NULL" || mesh_file == "NULL")
{
throw "No input grid text file or Gmsh mesh file";
}
// Prepare DEM parameters
int xnum = round((range[1] - range[0])/interval[0]) + 1;
int ynum = round((range[3] - range[2])/interval[1]) + 1;
std::vector<double> topo(xnum*ynum);
// Read DEM grid
double tmp_d;
std::ifstream infile;
infile.open(grid_file);
if (z_table)
{
for (int i = 0; i < xnum*ynum; i++)
{
infile >> topo[i];
}
}
else
{
for (int i = 0; i < xnum*ynum; i++)
{
infile >> tmp_d >> tmp_d >> topo[i];
}
}
infile.close();
// Read Polygon file
int tmp_count;
double tmp_x, tmp_y;
std::vector<vertex2dc> poly_vert;
if (poly_file != "NULL")
{
infile.open(poly_file);
infile >> tmp_count;
poly_vert.resize(tmp_count);
for (int i = 0; i < tmp_count; i++)
{
infile >> tmp_x >> tmp_y;
poly_vert[i].set(tmp_x, tmp_y, 0.0, i);
}
infile.close();
}
std::vector<double> err_records;
std::vector<vertex2dc*> tin_vert;
std::vector<triangle*> tin_ele;
if (poly_file == "NULL" && log_file == "NULL")
{
grd2tin(topo, range[0], range[1], range[2], range[3], interval[0], interval[1],
tin_vert, tin_ele, threshold, nullptr, nullptr);
}
else if (poly_file == "NULL")
{
grd2tin(topo, range[0], range[1], range[2], range[3], interval[0], interval[1],
tin_vert, tin_ele, threshold, nullptr, &err_records);
}
else if (log_file == "NULL")
{
grd2tin(topo, range[0], range[1], range[2], range[3], interval[0], interval[1],
tin_vert, tin_ele, threshold, &poly_vert, nullptr);
}
else
{
grd2tin(topo, range[0], range[1], range[2], range[3], interval[0], interval[1],
tin_vert, tin_ele, threshold, &poly_vert, &err_records);
}
// Write a log file
if (log_file != "NULL")
{
std::ofstream logfile(log_file);
logfile << "# Insertion Maxi-Error\n";
for (int i = 0; i < err_records.size(); ++i)
{
logfile << i+1 << " " << err_records[i] << std::endl;
}
logfile.close();
}
// Write a Gmsh's .msh file
std::ofstream outfile(mesh_file);
outfile << "$MeshFormat" << std::endl << "2.2 0 8" << std::endl << "$EndMeshFormat "<<std::endl;
outfile << "$Nodes" << std::endl << tin_vert.size() << std::endl;
for (int i = 0; i < tin_vert.size(); i++)
{
outfile << tin_vert[i]->id + 1 << " " << std::setprecision(16)
<< tin_vert[i]->x << " " << tin_vert[i]->y << " " << tin_vert[i]->elev << std::endl;
}
outfile<<"$EndNodes"<<std::endl;
outfile << "$Elements" << std::endl << tin_ele.size() <<std::endl;
for (int i = 0; i < tin_ele.size(); i++)
{
outfile << i + 1 << " 2 0";
for (int j = 0; j < 3; j++)
{
outfile << " " << tin_ele[i]->vert[j]->id + 1;
}
outfile << std::endl;
}
outfile << "$EndElements"<< std::endl;
outfile<<"$NodeData"<<std::endl;
outfile<<1<<std::endl
<<"\"Topography (m)\"" <<std::endl
<< 1 <<std::endl<< 0.0 <<std::endl
<< 3 <<std::endl<< 0<<std::endl
<< 1 <<std::endl<< tin_vert.size() <<std::endl;
for (int i = 0; i < tin_vert.size(); i++)
{
outfile << tin_vert[i]->id + 1 << " " << std::setprecision(16) << tin_vert[i]->elev << std::endl;
}
outfile << "$EndNodeData" << std::endl;
outfile.close();
// write a neighbor file
if (neigh_file != "NULL")
{
outfile.open(neigh_file);
outfile << tin_ele.size() << std::endl;
for (int i = 0; i < tin_ele.size(); i++)
{
outfile << i + 1;
for (int j = 0; j < 3; j++)
{
if (tin_ele[i]->neigh[j] != nullptr)
{
outfile << " " << tin_ele[i]->neigh[j]->id + 1;
}
else outfile << " -1";
}
outfile << std::endl;
}
outfile.close();
}
// Destroy memories allocated by the dem2tin function
for (int i = 0; i < tin_vert.size(); ++i)
{
delete tin_vert[i];
}
for (int i = 0; i < tin_ele.size(); ++i)
{
delete tin_ele[i];
}
}
catch(const char* err_char)
{
std::cerr << err_char << '\n';
}
return 0;
}