mirror
/
OpenDroneMap-ODM
kopia lustrzana https://github.com/OpenDroneMap/ODM
1
0
Forkuj 0
Kod Zgłoszenia Projekty Wydania Wiki Aktywność

z-occlusion filtering, removed outlier removal, still need testing

pull/542/head
Piero Toffanin 2017-04-17 17:39:47 -04:00
rodzic 0a50f5d35a
commit 9f4da6271d
5 zmienionych plików z 121 dodań i 142 usunięć
Pokaż statystyki Ściągnij plik aktualizacji Ściągnij plik porównania Expand all files Collapse all files

240
modules/odm_25dmeshing/src/Odm25dMeshing.cpp
Wyświetl plik

@ -1,25 +1,14 @@
#include "Odm25dMeshing.hpp"
//We define a vertex_base with info. The "info" (size_t) allow us to keep track of the original point index.
//typedef CGAL::Triangulation_vertex_base_with_info_2<size_t, Kernel> Vb;
//typedef CGAL::Triangulation_data_structure_2<Vb> Tds;
//typedef CGAL::Delaunay_triangulation_2<Kernel, Tds> DT;
//typedef DT::Point cgalPoint;
//typedef DT::Vertex_circulator Vertex_circulator;
typedef CGAL::Triangulation_vertex_base_with_info_2<size_t, Kernel> Vb;
typedef CGAL::Triangulation_data_structure_2<Vb> Tds;
typedef CGAL::Delaunay_triangulation_2<Kernel, Tds> DT;
typedef DT::Point cgalPoint;
typedef DT::Vertex_circulator Vertex_circulator;
typedef CGAL::Projection_traits_xy_3<Kernel> Gt;
typedef CGAL::Triangulation_vertex_base_2<Gt> Vb;
typedef CGAL::Delaunay_mesh_face_base_2<Gt> Fb;
typedef CGAL::Triangulation_data_structure_2<Vb, Fb> Tds;
typedef CGAL::Constrained_Delaunay_triangulation_2<Gt, Tds> CDT;
typedef CGAL::Delaunay_mesh_size_criteria_2<CDT> Criteria;
typedef CGAL::Delaunay_mesher_2<CDT, Criteria> Mesher;
typedef CDT::Vertex_handle Vertex_handle;
typedef CDT::Vertex_circulator Vertex_circulator;
typedef CDT::Point cgalPoint;
typedef CGAL::First_of_pair_property_map<Pwn> Point_map;
typedef CGAL::Second_of_pair_property_map<Pwn> Normal_map;
// Concurrency
#ifdef CGAL_LINKED_WITH_TBB
@ -41,7 +30,7 @@ int Odm25dMeshing::run(int argc, char **argv) {
parseArguments(argc, argv);
// loadPointCloud();
loadPointCloud();
buildMesh();
@ -83,15 +72,6 @@ void Odm25dMeshing::parseArguments(int argc, char **argv) {
if (ss.bad()) throw Odm25dMeshingException("Argument '" + argument + "' has a bad value (wrong type).");
maxVertexCount = std::max<unsigned int>(maxVertexCount, 0);
log << "Vertex count was manually set to: " << maxVertexCount << "\n";
} else if (argument == "-outliersRemovalPercentage" && argIndex < argc) {
++argIndex;
if (argIndex >= argc) throw Odm25dMeshingException("Argument '" + argument + "' expects 1 more input following it, but no more inputs were provided.");
std::stringstream ss(argv[argIndex]);
ss >> outliersRemovalPercentage;
if (ss.bad()) throw Odm25dMeshingException("Argument '" + argument + "' has a bad value (wrong type).");
outliersRemovalPercentage = std::min<double>(99.99, std::max<double>(outliersRemovalPercentage, 0));
log << "Outliers removal was manually set to: " << outliersRemovalPercentage << "\n";
} else if (argument == "-wlopIterations" && argIndex < argc) {
++argIndex;
if (argIndex >= argc) throw Odm25dMeshingException("Argument '" + argument + "' expects 1 more input following it, but no more inputs were provided.");
@ -167,90 +147,53 @@ void Odm25dMeshing::loadPointCloud(){
}
void Odm25dMeshing::buildMesh(){
// const unsigned int NEIGHBORS = 24;
//
// if (outliersRemovalPercentage > 0){
// size_t pointCountBeforeOutRemoval = points.size();
// log << "Removing outliers\n";
//
// points.erase(CGAL::remove_outliers(points.begin(), points.end(),
// CGAL::First_of_pair_property_map<Pwn>(), NEIGHBORS, outliersRemovalPercentage),
// points.end());
// std::vector<Pwn>(points).swap(points);
//
// log << "Removed " << (pointCountBeforeOutRemoval - points.size()) << " points\n";
// }
//
// // Filter points by normal orientation
// log << "Applying normal filtering\n";
//
// double COSINE_LIMIT = std::cos(M_PI / 4); // ~ 45 degrees
//
// Vector3 up(0, 0, 1);
// if (flipFaces) up = Vector3(0, 0, -1);
//
// size_t pointCountBeforeNormalFilter = points.size();
//
// points.erase(std::remove_if(points.begin(), points.end(),[&](Pwn pointWithNormal){
// double cosine = pointWithNormal.second * up / CGAL::sqrt(pointWithNormal.second*pointWithNormal.second) / CGAL::sqrt(up * up);
// return cosine <= COSINE_LIMIT;
// }), points.end());
//
const unsigned int NEIGHBORS = 16;
size_t pointCount = points.size();
//
// log << "Removed " << (pointCountBeforeNormalFilter - pointCount) << " points\n";
//
// const double RETAIN_PERCENTAGE = std::min<double>(((100. * (double)maxVertexCount) / (double)pointCount), 80.); // percentage of points to retain.
std::vector<Point3> simplifiedPoints;
//
// log << "Computing average spacing... ";
//
// FT avgSpacing = CGAL::compute_average_spacing<Concurrency_tag>(
// points.begin(),
// points.end(),
// CGAL::First_of_pair_property_map<Pwn>(),
// NEIGHBORS);
//
// log << avgSpacing << "\n";
//
// log << "Performing weighted locally optimal projection simplification and regularization (retain: " << RETAIN_PERCENTAGE << "%, iterate: " << wlopIterations << ")" << "\n";
//
// CGAL::wlop_simplify_and_regularize_point_set<Concurrency_tag>(
// points.begin(),
// points.end(),
// std::back_inserter(simplifiedPoints),
// CGAL::First_of_pair_property_map<Pwn>(),
// RETAIN_PERCENTAGE,
// 8 * avgSpacing,
// wlopIterations,
// false); // require_uniform_sampling
// for (size_t c = 0; c < pointCount; c++){
// simplifiedPoints.push_back(points[c].first);
// }
// std::ofstream f("out.bin");
// f << simplifiedPoints.size() * 3 << " ";
//
// for (size_t i = 0; i < simplifiedPoints.size(); i++){
// f << simplifiedPoints[i].x() << " ";
// f << simplifiedPoints[i].y() << " ";
// f << simplifiedPoints[i].z() << " ";
// }
// f.close();
// exit(1);
std::ifstream fin("out.bin");
pointCount = 0;
fin >> pointCount;
for (size_t i = 0; i < pointCount / 3; i++){
float x, y, z;
fin >> x;
fin >> y;
fin >> z;
simplifiedPoints.push_back(Point3(x, y, z));
log << "Performing bilateral smoothing... ";
for (int i = 0; i < 3; i++){
double err = CGAL::bilateral_smooth_point_set<Concurrency_tag>(
points.begin(),
points.end(),
Point_map(),
Normal_map(),
NEIGHBORS,
45.0,
Kernel());
log << err << " ";
if (err < 0.001){
log << "stopping early\n";
}else if (i >= 2){
log << "iteration limit reached\n";
}
}
const double RETAIN_PERCENTAGE = std::min<double>(((100. * (double)maxVertexCount) / (double)pointCount), 80.); // percentage of points to retain.
std::vector<Point3> simplifiedPoints;
log << "Computing average spacing... ";
FT avgSpacing = CGAL::compute_average_spacing<Concurrency_tag>(
points.begin(),
points.end(),
Point_map(),
NEIGHBORS);
log << avgSpacing << "\n";
log << "Performing weighted locally optimal projection simplification and regularization (retain: " << RETAIN_PERCENTAGE << "%, iterate: " << wlopIterations << ")" << "\n";
CGAL::wlop_simplify_and_regularize_point_set<Concurrency_tag>(
points.begin(),
points.end(),
std::back_inserter(simplifiedPoints),
Point_map(),
RETAIN_PERCENTAGE,
8 * avgSpacing,
wlopIterations,
true);
pointCount = simplifiedPoints.size();
if (pointCount < 3){
@ -259,21 +202,74 @@ void Odm25dMeshing::buildMesh(){
log << "Vertex count is " << pointCount << "\n";
// std::vector< std::pair<cgalPoint, size_t > > pts;
// try{
// pts.reserve(pointCount);
// } catch (const std::bad_alloc&){
// throw Odm25dMeshingException("Not enough memory");
// }
log << "Z-occlusion filtering\n";
size_t pointCountBeforeOcclusionFiltering = pointCount;
const double gridStep = 0.05;
Kernel::Iso_cuboid_3 bbox = CGAL::bounding_box(simplifiedPoints.begin(), simplifiedPoints.end());
Vector3 boxDiag = bbox.max() - bbox.min();
int gridWidth = 1 + static_cast<unsigned>(boxDiag.x() / gridStep + 0.5);
int gridHeight = 1 + static_cast<unsigned>(boxDiag.y() / gridStep + 0.5);
std::unordered_map<int, Point3> grid;
for (size_t c = 0; c < pointCount; c++){
const Point3 &p = simplifiedPoints[c];
Vector3 relativePos = p - bbox.min();
int i = static_cast<int>((relativePos.x() / gridStep + 0.5));
int j = static_cast<int>((relativePos.y() / gridStep + 0.5));
if ((i >= 0 && i < gridWidth) && (j >= 0 && j < gridHeight)){
int key = i * gridWidth + j;
if (grid.find(key) == grid.end()){
grid[key] = p;
}else if (p.z() > grid[key].z()){
grid[key] = p;
}
//
// for (size_t i = 0; i < pointCount; ++i){
// pts.push_back(std::make_pair(cgalPoint(simplifiedPoints[i].x(), simplifiedPoints[i].y()), i));
// if (grid[i][j] == NULL){
// grid[i][j] = static_cast<Point3*>(&p);
// }else if (p.z() > grid[i][j]->z()){
// grid[i][j] = static_cast<Point3*>(&p);
// }
}
}
std::vector<Point3> gridPoints;
for ( auto it = grid.begin(); it != grid.end(); ++it ){
gridPoints.push_back(it->second);
}
// for (int i = 0; i < gridWidth; i++){
// for (int j = 0; j < gridHeight; j++){
// if (grid[i][j] != NULL){
// gridPoints.push_back(*grid[i][j]);
// }
// }
// }
pointCount = gridPoints.size();
log << "Removed " << (pointCountBeforeOcclusionFiltering - pointCount) << " points\n";
std::vector< std::pair<cgalPoint, size_t > > pts;
try{
pts.reserve(pointCount);
} catch (const std::bad_alloc&){
throw Odm25dMeshingException("Not enough memory");
}
for (size_t i = 0; i < pointCount; ++i){
pts.push_back(std::make_pair(cgalPoint(gridPoints[i].x(), gridPoints[i].y()), i));
}
log << "Computing delaunay triangulation\n";
//The delaunay triangulation is built according to the 2D point cloud
// DT dt(pts.begin(), pts.end());
DT dt(pts.begin(), pts.end());
unsigned int numberOfTriangles = static_cast<unsigned >(dt.number_of_faces());
unsigned int triIndexes = dt.number_of_faces()*3;
@ -295,9 +291,9 @@ void Odm25dMeshing::buildMesh(){
for (size_t i = 0; i < pointCount; ++i){
vertices.push_back(simplifiedPoints[i].x());
vertices.push_back(simplifiedPoints[i].y());
vertices.push_back(simplifiedPoints[i].z());
vertices.push_back(gridPoints[i].x());
vertices.push_back(gridPoints[i].y());
vertices.push_back(gridPoints[i].z());
}
for (DT::Face_iterator face = dt.faces_begin(); face != dt.faces_end(); ++face) {
@ -314,12 +310,11 @@ void Odm25dMeshing::buildMesh(){
log << "Removing spikes\n";
const float THRESHOLD = 0.1;
const float THRESHOLD = 0.05;
std::vector<float> heights;
unsigned int spikesRemoved = 0;
for (DT::Vertex_iterator vertex = dt.vertices_begin(); vertex != dt.vertices_end(); ++vertex){
// Check if the height between this vertex and its
// incident vertices is greater than THRESHOLD
@ -387,7 +382,6 @@ void Odm25dMeshing::printHelp() {
<< " -outputFile <path> where the output PLY 2.5D mesh should be saved (default: " << outputFile << ")\n"
<< " -logFile <path> log file path (default: " << logFilePath << ")\n"
<< " -verbose whether to print verbose output (default: " << (printInCoutPop ? "true" : "false") << ")\n"
<< " -outliersRemovalPercentage <0 - 99.99> percentage of outliers to remove. Set to 0 to disable. (default: " << outliersRemovalPercentage << ")\n"
<< " -maxVertexCount <0 - N> Maximum number of vertices in the output mesh. The mesh might have fewer vertices, but will not exceed this limit. (default: " << maxVertexCount << ")\n"
<< " -wlopIterations <1 - 1000> Iterations of the Weighted Locally Optimal Projection (WLOP) simplification algorithm. Higher values take longer but produce a smoother mesh. (default: " << wlopIterations << ")\n"

8
modules/odm_25dmeshing/src/Odm25dMeshing.hpp
Wyświetl plik

@ -3,17 +3,16 @@
// STL
#include <string>
#include <iostream>
#include <unordered_map>
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Constrained_Delaunay_triangulation_2.h>
#include <CGAL/Triangulation_vertex_base_with_info_2.h>
#include <CGAL/Delaunay_triangulation_2.h>
#include <CGAL/Triangulation_2.h>
#include <CGAL/Projection_traits_xy_3.h>
#include <CGAL/remove_outliers.h>
#include <CGAL/wlop_simplify_and_regularize_point_set.h>
#include <CGAL/nearest_neighbor_delaunay_2.h>
#include <CGAL/compute_average_spacing.h>
#include <CGAL/bilateral_smooth_point_set.h>
#include <CGAL/bounding_box.h>
#include "Logger.hpp"
#include "PlyInterpreter.hpp"
@ -65,7 +64,6 @@ private:
std::string outputFile = "odm_25dmesh.ply";
std::string logFilePath = "odm_25dmeshing_log.txt";
unsigned int maxVertexCount = 100000;
double outliersRemovalPercentage = 2;
unsigned int wlopIterations = 35;
std::vector<Pwn> points;
bool flipFaces = false;

8
opendm/config.py
Wyświetl plik

@ -262,14 +262,6 @@ def config():
'Applies to 2.5D mesh only. '
'Default: %(default)s'))
parser.add_argument('--mesh-remove-outliers',
metavar='<percent>',
default=2,
type=float,
help=('Percentage of outliers to remove from the point set. Set to 0 to disable. '
'Applies to 2.5D mesh only. '
'Default: %(default)s'))
parser.add_argument('--mesh-wlop-iterations',
metavar='<positive integer>',
default=70,

1
scripts/odm_app.py
Wyświetl plik

@ -59,7 +59,6 @@ class ODMApp(ecto.BlackBox):
samples=p.args.mesh_samples,
solver=p.args.mesh_solver_divide,
max_vertex_25d=getattr(p.args, '25d_mesh_size'),
remove_outliers=p.args.mesh_remove_outliers,
wlop_iterations=p.args.mesh_wlop_iterations,
verbose=p.args.verbose),
'texturing': ODMMvsTexCell(data_term=p.args.texturing_data_term,

6
scripts/odm_meshing.py
Wyświetl plik

@ -21,9 +21,6 @@ class ODMeshingCell(ecto.Cell):
'times slightly but helps reduce memory usage.', 9)
params.declare("max_vertex_25d", 'The maximum vertex count of the 2.5D output mesh.', 50000)
params.declare("remove_outliers", 'Percentage of outliers to remove from the point set. Set to 0 to disable. '
'Applies to 2.5D mesh only.', 2)
params.declare("wlop_iterations", 'Iterations of the Weighted Locally Optimal Projection (WLOP) simplification algorithm. '
'Higher values take longer but produce a smoother mesh. '
'Applies to 2.5D mesh only. ', 70)
@ -98,14 +95,13 @@ class ODMeshingCell(ecto.Cell):
'log': tree.odm_25dmeshing_log,
'verbose': verbose,
'max_vertex': self.params.max_vertex_25d,
'remove_outliers': self.params.remove_outliers,
'wlop_iterations': self.params.wlop_iterations
}
# run 2.5D meshing binary
system.run('{bin}/odm_25dmeshing -inputFile {infile} '
'-outputFile {outfile} -logFile {log} '
'-maxVertexCount {max_vertex} -outliersRemovalPercentage {remove_outliers} '
'-maxVertexCount {max_vertex} '
'-wlopIterations {wlop_iterations} {verbose}'.format(**kwargs))
else:
log.ODM_WARNING('Found a valid ODM 2.5D Mesh file in: %s' %