diff --git a/.gitattributes b/.gitattributes new file mode 100644 index 0000000..68bc67b --- /dev/null +++ b/.gitattributes @@ -0,0 +1,20 @@ +# Set the default behavior, in case people don't have core.autocrlf set. +* text=auto + +# Explicitly declare text files you want to always be normalized and converted +# to native line endings on checkout. +*.tex text +*.bib text +*.svg text +*.py text +*.vbs text +*.cpp text +*.hpp text +Makefile text + +# Declare files that will always have CRLF line endings on checkout. +*.sln text eol=crlf + +# Denote all files that are truly binary and should not be modified. +*.png binary +*.jpg binary diff --git a/.gitignore b/.gitignore index 6e81a50..ba969e9 100644 --- a/.gitignore +++ b/.gitignore @@ -1,3 +1,8 @@ +# Extern +/extern/.cache +/extern/spdlog + +# Build src/ispc/energy.h *build*/* .idea/* diff --git a/CMakeLists.txt b/CMakeLists.txt index e41648b..e36c9cb 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -8,7 +8,6 @@ project(TetWild) # Options ################################################################################ # tetwild -option(TETWILD_QUIET "Mute log and unnecassary checks" OFF) option(TETWILD_WITH_ISPC "Use ISPC" OFF) # libigl library option(LIBIGL_USE_STATIC_LIBRARY "Use libigl as static library" OFF) @@ -44,6 +43,9 @@ endif() set(TETWILD_EXTERNAL "${CMAKE_CURRENT_SOURCE_DIR}/extern") list(APPEND CMAKE_MODULE_PATH ${CMAKE_CURRENT_SOURCE_DIR}/cmake) +# Color output +include(UseColors) + # Use folder in Visual Studio set_property(GLOBAL PROPERTY USE_FOLDERS ON) @@ -56,6 +58,16 @@ set(CMAKE_POSITION_INDEPENDENT_CODE ON) ################################################################################ # 3rd party libraries ################################################################################ + +# download external dependencies +include(TetWildDownloadExternal) + +# spdlog +if(NOT TARGET spdlog::spdlog) + tetwild_download_spdlog() + add_subdirectory(${TETWILD_EXTERNAL}/spdlog) +endif() + # cgal find_package(CGAL REQUIRED) @@ -71,44 +83,45 @@ add_subdirectory(${TETWILD_EXTERNAL}/pymesh) # CL11 add_subdirectory(${TETWILD_EXTERNAL}/CLI) -# cout -if(TETWILD_QUIET) - add_definitions(-DMUTE_COUT) -endif() - ################################################################################ # TetWild ################################################################################ # Build static library for executable add_library(libTetWild STATIC - src/tetwild/Preprocess.h - src/tetwild/Preprocess.cpp - src/tetwild/DelaunayTetrahedralization.h - src/tetwild/DelaunayTetrahedralization.cpp - src/tetwild/MeshConformer.h - src/tetwild/MeshConformer.cpp - src/tetwild/BSPSubdivision.h src/tetwild/BSPSubdivision.cpp - src/tetwild/SimpleTetrahedralization.h - src/tetwild/SimpleTetrahedralization.cpp - src/tetwild/MeshRefinement.h - src/tetwild/MeshRefinement.cpp - src/tetwild/LocalOperations.h - src/tetwild/LocalOperations.cpp - src/tetwild/EdgeSplitter.h - src/tetwild/EdgeSplitter.cpp - src/tetwild/EdgeCollapser.h + src/tetwild/BSPSubdivision.h + src/tetwild/CGALCommon.cpp + src/tetwild/CGALCommon.h + src/tetwild/Common.cpp + src/tetwild/Common.h + src/tetwild/DelaunayTetrahedralization.cpp + src/tetwild/DelaunayTetrahedralization.h src/tetwild/EdgeCollapser.cpp - src/tetwild/EdgeRemover.h + src/tetwild/EdgeCollapser.h src/tetwild/EdgeRemover.cpp - src/tetwild/VertexSmoother.h - src/tetwild/VertexSmoother.cpp - src/tetwild/InoutFiltering.h + src/tetwild/EdgeRemover.h + src/tetwild/EdgeSplitter.cpp + src/tetwild/EdgeSplitter.h src/tetwild/InoutFiltering.cpp - src/tetwild/heads.h - src/tetwild/heads.cpp - src/tetwild/tetwild.cpp + src/tetwild/InoutFiltering.h + src/tetwild/LocalOperations.cpp + src/tetwild/LocalOperations.h + src/tetwild/Logger.cpp + src/tetwild/Logger.h + src/tetwild/MeshConformer.cpp + src/tetwild/MeshConformer.h + src/tetwild/MeshRefinement.cpp + src/tetwild/MeshRefinement.h + src/tetwild/Preprocess.cpp + src/tetwild/Preprocess.h + src/tetwild/SimpleTetrahedralization.cpp + src/tetwild/SimpleTetrahedralization.h + src/tetwild/State.cpp + src/tetwild/State.h src/tetwild/TetmeshElements.h + src/tetwild/tetwild.cpp + src/tetwild/VertexSmoother.cpp + src/tetwild/VertexSmoother.h include/tetwild/tetwild.h ) target_include_directories(libTetWild @@ -123,6 +136,7 @@ target_link_libraries(libTetWild igl::core igl::cgal pymesh::pymesh + spdlog::spdlog ) set_target_properties(libTetWild PROPERTIES OUTPUT_NAME "tetwild") diff --git a/cmake/DownloadProject.CMakeLists.cmake.in b/cmake/DownloadProject.CMakeLists.cmake.in new file mode 100644 index 0000000..89be4fd --- /dev/null +++ b/cmake/DownloadProject.CMakeLists.cmake.in @@ -0,0 +1,17 @@ +# Distributed under the OSI-approved MIT License. See accompanying +# file LICENSE or https://github.com/Crascit/DownloadProject for details. + +cmake_minimum_required(VERSION 2.8.2) + +project(${DL_ARGS_PROJ}-download NONE) + +include(ExternalProject) +ExternalProject_Add(${DL_ARGS_PROJ}-download + ${DL_ARGS_UNPARSED_ARGUMENTS} + SOURCE_DIR "${DL_ARGS_SOURCE_DIR}" + BINARY_DIR "${DL_ARGS_BINARY_DIR}" + CONFIGURE_COMMAND "" + BUILD_COMMAND "" + INSTALL_COMMAND "" + TEST_COMMAND "" +) diff --git a/cmake/DownloadProject.cmake b/cmake/DownloadProject.cmake new file mode 100644 index 0000000..e300f42 --- /dev/null +++ b/cmake/DownloadProject.cmake @@ -0,0 +1,182 @@ +# Distributed under the OSI-approved MIT License. See accompanying +# file LICENSE or https://github.com/Crascit/DownloadProject for details. +# +# MODULE: DownloadProject +# +# PROVIDES: +# download_project( PROJ projectName +# [PREFIX prefixDir] +# [DOWNLOAD_DIR downloadDir] +# [SOURCE_DIR srcDir] +# [BINARY_DIR binDir] +# [QUIET] +# ... +# ) +# +# Provides the ability to download and unpack a tarball, zip file, git repository, +# etc. at configure time (i.e. when the cmake command is run). How the downloaded +# and unpacked contents are used is up to the caller, but the motivating case is +# to download source code which can then be included directly in the build with +# add_subdirectory() after the call to download_project(). Source and build +# directories are set up with this in mind. +# +# The PROJ argument is required. The projectName value will be used to construct +# the following variables upon exit (obviously replace projectName with its actual +# value): +# +# projectName_SOURCE_DIR +# projectName_BINARY_DIR +# +# The SOURCE_DIR and BINARY_DIR arguments are optional and would not typically +# need to be provided. They can be specified if you want the downloaded source +# and build directories to be located in a specific place. The contents of +# projectName_SOURCE_DIR and projectName_BINARY_DIR will be populated with the +# locations used whether you provide SOURCE_DIR/BINARY_DIR or not. +# +# The DOWNLOAD_DIR argument does not normally need to be set. It controls the +# location of the temporary CMake build used to perform the download. +# +# The PREFIX argument can be provided to change the base location of the default +# values of DOWNLOAD_DIR, SOURCE_DIR and BINARY_DIR. If all of those three arguments +# are provided, then PREFIX will have no effect. The default value for PREFIX is +# CMAKE_BINARY_DIR. +# +# The QUIET option can be given if you do not want to show the output associated +# with downloading the specified project. +# +# In addition to the above, any other options are passed through unmodified to +# ExternalProject_Add() to perform the actual download, patch and update steps. +# The following ExternalProject_Add() options are explicitly prohibited (they +# are reserved for use by the download_project() command): +# +# CONFIGURE_COMMAND +# BUILD_COMMAND +# INSTALL_COMMAND +# TEST_COMMAND +# +# Only those ExternalProject_Add() arguments which relate to downloading, patching +# and updating of the project sources are intended to be used. Also note that at +# least one set of download-related arguments are required. +# +# If using CMake 3.2 or later, the UPDATE_DISCONNECTED option can be used to +# prevent a check at the remote end for changes every time CMake is run +# after the first successful download. See the documentation of the ExternalProject +# module for more information. It is likely you will want to use this option if it +# is available to you. Note, however, that the ExternalProject implementation contains +# bugs which result in incorrect handling of the UPDATE_DISCONNECTED option when +# using the URL download method or when specifying a SOURCE_DIR with no download +# method. Fixes for these have been created, the last of which is scheduled for +# inclusion in CMake 3.8.0. Details can be found here: +# +# https://gitlab.kitware.com/cmake/cmake/commit/bdca68388bd57f8302d3c1d83d691034b7ffa70c +# https://gitlab.kitware.com/cmake/cmake/issues/16428 +# +# If you experience build errors related to the update step, consider avoiding +# the use of UPDATE_DISCONNECTED. +# +# EXAMPLE USAGE: +# +# include(DownloadProject) +# download_project(PROJ googletest +# GIT_REPOSITORY https://github.com/google/googletest.git +# GIT_TAG master +# UPDATE_DISCONNECTED 1 +# QUIET +# ) +# +# add_subdirectory(${googletest_SOURCE_DIR} ${googletest_BINARY_DIR}) +# +#======================================================================================== + + +set(_DownloadProjectDir "${CMAKE_CURRENT_LIST_DIR}") + +include(CMakeParseArguments) + +function(download_project) + + set(options QUIET) + set(oneValueArgs + PROJ + PREFIX + DOWNLOAD_DIR + SOURCE_DIR + BINARY_DIR + # Prevent the following from being passed through + CONFIGURE_COMMAND + BUILD_COMMAND + INSTALL_COMMAND + TEST_COMMAND + ) + set(multiValueArgs "") + + cmake_parse_arguments(DL_ARGS "${options}" "${oneValueArgs}" "${multiValueArgs}" ${ARGN}) + + # Hide output if requested + if (DL_ARGS_QUIET) + set(OUTPUT_QUIET "OUTPUT_QUIET") + else() + unset(OUTPUT_QUIET) + message(STATUS "Downloading/updating ${DL_ARGS_PROJ}") + endif() + + # Set up where we will put our temporary CMakeLists.txt file and also + # the base point below which the default source and binary dirs will be. + # The prefix must always be an absolute path. + if (NOT DL_ARGS_PREFIX) + set(DL_ARGS_PREFIX "${CMAKE_BINARY_DIR}") + else() + get_filename_component(DL_ARGS_PREFIX "${DL_ARGS_PREFIX}" ABSOLUTE + BASE_DIR "${CMAKE_CURRENT_BINARY_DIR}") + endif() + if (NOT DL_ARGS_DOWNLOAD_DIR) + set(DL_ARGS_DOWNLOAD_DIR "${DL_ARGS_PREFIX}/${DL_ARGS_PROJ}-download") + endif() + + # Ensure the caller can know where to find the source and build directories + if (NOT DL_ARGS_SOURCE_DIR) + set(DL_ARGS_SOURCE_DIR "${DL_ARGS_PREFIX}/${DL_ARGS_PROJ}-src") + endif() + if (NOT DL_ARGS_BINARY_DIR) + set(DL_ARGS_BINARY_DIR "${DL_ARGS_PREFIX}/${DL_ARGS_PROJ}-build") + endif() + set(${DL_ARGS_PROJ}_SOURCE_DIR "${DL_ARGS_SOURCE_DIR}" PARENT_SCOPE) + set(${DL_ARGS_PROJ}_BINARY_DIR "${DL_ARGS_BINARY_DIR}" PARENT_SCOPE) + + # The way that CLion manages multiple configurations, it causes a copy of + # the CMakeCache.txt to be copied across due to it not expecting there to + # be a project within a project. This causes the hard-coded paths in the + # cache to be copied and builds to fail. To mitigate this, we simply + # remove the cache if it exists before we configure the new project. It + # is safe to do so because it will be re-generated. Since this is only + # executed at the configure step, it should not cause additional builds or + # downloads. + file(REMOVE "${DL_ARGS_DOWNLOAD_DIR}/CMakeCache.txt") + + # Create and build a separate CMake project to carry out the download. + # If we've already previously done these steps, they will not cause + # anything to be updated, so extra rebuilds of the project won't occur. + # Make sure to pass through CMAKE_MAKE_PROGRAM in case the main project + # has this set to something not findable on the PATH. + configure_file("${_DownloadProjectDir}/DownloadProject.CMakeLists.cmake.in" + "${DL_ARGS_DOWNLOAD_DIR}/CMakeLists.txt") + execute_process(COMMAND ${CMAKE_COMMAND} -G "${CMAKE_GENERATOR}" + -D "CMAKE_MAKE_PROGRAM:FILE=${CMAKE_MAKE_PROGRAM}" + . + RESULT_VARIABLE result + ${OUTPUT_QUIET} + WORKING_DIRECTORY "${DL_ARGS_DOWNLOAD_DIR}" + ) + if(result) + message(FATAL_ERROR "CMake step for ${DL_ARGS_PROJ} failed: ${result}") + endif() + execute_process(COMMAND ${CMAKE_COMMAND} --build . + RESULT_VARIABLE result + ${OUTPUT_QUIET} + WORKING_DIRECTORY "${DL_ARGS_DOWNLOAD_DIR}" + ) + if(result) + message(FATAL_ERROR "Build step for ${DL_ARGS_PROJ} failed: ${result}") + endif() + +endfunction() diff --git a/cmake/TetWildDownloadExternal.cmake b/cmake/TetWildDownloadExternal.cmake new file mode 100644 index 0000000..6e8cdf3 --- /dev/null +++ b/cmake/TetWildDownloadExternal.cmake @@ -0,0 +1,22 @@ +################################################################################ +include(DownloadProject) + +# Shortcut function +function(tetwild_download_project name) + download_project( + PROJ ${name} + SOURCE_DIR ${TETWILD_EXTERNAL}/${name} + DOWNLOAD_DIR ${TETWILD_EXTERNAL}/.cache/${name} + ${ARGN} + ) +endfunction() + +################################################################################ + +## spdlog +function(tetwild_download_spdlog) + tetwild_download_project(spdlog + GIT_REPOSITORY https://github.com/gabime/spdlog.git + GIT_TAG v1.1.0 + ) +endfunction() diff --git a/cmake/UseColors.cmake b/cmake/UseColors.cmake new file mode 100644 index 0000000..7801756 --- /dev/null +++ b/cmake/UseColors.cmake @@ -0,0 +1,47 @@ +################################################################################ +# When using Clang, there is nothing to do: colors are enabled by default +# When using GCC >= 4.9, colored diagnostics can be enabled natively +# When using an older version, one can use gccfilter (a perl script) +# +# I do not recommend using gccfilter as of now (May 2014), because it seems to +# be bugged. But if you still want to try, here is how to install it on Ubuntu: +# +# +# 1) Download the perl script and add it to you $PATH +# mkdir -p ~/.local/bin +# wget -P ~/.local/bin http://www.mixtion.org/gccfilter/gccfilter +# chmod +x ~/local/bin/gccfilter +# echo 'PATH="$HOME/.local/bin:$PATH"' >> ~/.bashrc +# +# 2) Install the dependencies +# * Term::ANSIColor +# sudo cpan +# cpan> install Term::ANSIColor +# * The module "Getopt::Long" is included in "perl-base" +# * For Getopt::ArgvFile and Regexp::Common ... +# sudo apt-get install libgetopt-argvfile-perl libregexp-common-perl +# +################################################################################ + +if(CMAKE_COMPILER_IS_GNUCXX) + # If GCC >= 4.9, just activate the right option + # We enable colorized diagnostics always instead of using "auto" so that + # they're still colored when using Ninja. + if(NOT CMAKE_CXX_COMPILER_VERSION VERSION_LESS 4.9) + message(STATUS "GCC >= 4.9 detected, enabling colored diagnostics") + add_definitions(-fdiagnostics-color=always) + set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -fdiagnostics-color=always") + set(CMAKE_C_FLAGS "${CMAKE_C_FLAGS} -fdiagnostics-color=always") + return() + endif() + # If GCC < 4.9, maybe we can use gccfilter + find_program(GCC_FILTER gccfilter) + if(GCC_FILTER) + option(COLOR_GCC "Use GCCFilter to color compiler output messages" OFF) + set(COLOR_GCC_OPTIONS "-c -r -w" CACHE STRING "Arguments that are passed to gccfilter when output coloring is switchend on. Defaults to -c -r -w.") + if(COLOR_GCC) + set_property(GLOBAL PROPERTY RULE_LAUNCH_COMPILE "${GCC_FILTER} ${COLOR_GCC_OPTIONS}") + message(STATUS "Using gccfilter for colored diagnostics") + endif() + endif() +endif() diff --git a/extern/pymesh/MshLoader.h b/extern/pymesh/MshLoader.h index 99fc67c..3b08311 100644 --- a/extern/pymesh/MshLoader.h +++ b/extern/pymesh/MshLoader.h @@ -8,12 +8,13 @@ #include #include + +namespace PyMesh { + typedef double Float; typedef Eigen::VectorXd VectorF; typedef Eigen::VectorXi VectorI; -namespace PyMesh { - class MshLoader { public: typedef std::map FieldMap; @@ -82,4 +83,4 @@ class MshLoader { FieldMap m_element_fields; }; -} \ No newline at end of file +} diff --git a/extern/pymesh/MshSaver.cpp b/extern/pymesh/MshSaver.cpp index ff7c63f..ed92f01 100644 --- a/extern/pymesh/MshSaver.cpp +++ b/extern/pymesh/MshSaver.cpp @@ -370,4 +370,4 @@ void MshSaver::save_elem_tensor_field(const std::string& fieldname, const Vector fout << "$EndElementData" << std::endl; fout.flush(); -} \ No newline at end of file +} diff --git a/extern/pymesh/MshSaver.h b/extern/pymesh/MshSaver.h index 81419c2..105727c 100644 --- a/extern/pymesh/MshSaver.h +++ b/extern/pymesh/MshSaver.h @@ -5,12 +5,13 @@ #include #include #include + +namespace PyMesh { + typedef double Float; typedef Eigen::VectorXd VectorF; typedef Eigen::VectorXi VectorI; -namespace PyMesh { - class MshSaver { public: MshSaver(const std::string& filename, bool binary=true); @@ -52,4 +53,4 @@ class MshSaver { std::ofstream fout; }; -} \ No newline at end of file +} diff --git a/include/tetwild/tetwild.h b/include/tetwild/tetwild.h index 9d04a7a..721a2e6 100644 --- a/include/tetwild/tetwild.h +++ b/include/tetwild/tetwild.h @@ -16,7 +16,8 @@ #include namespace tetwild { - struct Args{ + + struct Args { double i_ideal_edge_length = 20; double i_epsilon = 1000; int stage = 1; @@ -28,12 +29,12 @@ namespace tetwild { std::string bg_mesh = ""; bool is_quiet = false; }; - extern Args parameters; void tetrahedralization(const std::vector>& V_in, const std::vector>& F_in, std::vector>& V_out, - std::vector>& T_out); + std::vector>& T_out, + Args params = Args()); } #endif //TETWILD_TETWILD_H diff --git a/src/main.cpp b/src/main.cpp index 3985dd9..91ed327 100644 --- a/src/main.cpp +++ b/src/main.cpp @@ -5,21 +5,25 @@ // This Source Code Form is subject to the terms of the Mozilla Public License // v. 2.0. If a copy of the MPL was not distributed with this file, You can // obtain one at http://mozilla.org/MPL/2.0/. -#include +#include #include #include #include #include #include #include +#include +#include #include +using namespace tetwild; + MeshRefinement MR; void outputFinalQuality(double time, const std::vector& tet_vertices, const std::vector>& tets, const std::vector &t_is_removed, const std::vector& tet_qualities, const std::vector& v_ids) { - cout << "final quality:" << endl; + logger().debug("final quality:"); double min = 10, max = 0; double min_avg = 0, max_avg = 0; // double max_asp_ratio = 0, avg_asp_ratio = 0; @@ -54,20 +58,10 @@ void outputFinalQuality(double time, const std::vector& tet_vertices, cmp_cnt[j + 3]++; } } - cout << "min_d_angle = " << min - << ", max_d_angle = " << max -// << ", max_aspect_ratio = " << max_asp_ratio - << ", max_slim_energy = " << max_slim_energy - << endl; - cout << "avg_min_d_angle = " << min_avg / cnt - << ", avg_max_d_angle = " << max_avg / cnt -// << ", avg_aspect_ratio = " << avg_asp_ratio / cnt - << ", avg_slim_energy = " << avg_slim_energy / cnt - << endl; - cout << "min_d_angle: <6 " << cmp_cnt[0] / cnt << "; <12 " << cmp_cnt[1] / cnt << "; <18 " << cmp_cnt[2] / cnt - << endl; - cout << "max_d_angle: >174 " << cmp_cnt[5] / cnt << "; >168 " << cmp_cnt[4] / cnt << "; >162 " << cmp_cnt[3] / cnt - << endl; + logger().debug("min_d_angle = {}, max_d_angle = {}, max_slim_energy = {}", min, max, max_slim_energy); + logger().debug("avg_min_d_angle = {}, avg_max_d_angle = {}, avg_slim_energy = {}", min_avg / cnt, max_avg / cnt, avg_slim_energy / cnt); + logger().debug("min_d_angle: <6 {}; <12 {}; <18 {}", cmp_cnt[0] / cnt, cmp_cnt[1] / cnt, cmp_cnt[2] / cnt); + logger().debug("max_d_angle: >174 {}; >168 {}; >162 {}", cmp_cnt[5] / cnt, cmp_cnt[4] / cnt, cmp_cnt[3] / cnt); addRecord(MeshRecord(MeshRecord::OpType::OP_WN, time, v_ids.size(), cnt, min, min_avg / cnt, max, max_avg / cnt, max_slim_energy, avg_slim_energy / cnt)); @@ -78,7 +72,7 @@ void outputFinalQuality(double time, const std::vector& tet_vertices, if (!tet_vertices[v_id].is_rounded) cnt++; } - cout << cnt << "/" << v_ids.size() << " vertices are unrounded!!!" << endl; + logger().debug("{}/{} vertices are unrounded!!!", cnt, v_ids.size()); addRecord(MeshRecord(MeshRecord::OpType::OP_UNROUNDED, -1, cnt, -1)); } @@ -140,7 +134,7 @@ void outputFinalSurface(MeshRefinement& MR){ break; } } - igl::writeOBJ(g_working_dir+g_postfix+"_sf.obj", V_sf, F_sf); + igl::writeOBJ(State::state().g_working_dir+State::state().g_postfix+"_sf.obj", V_sf, F_sf); } void outputFinalTetmesh(MeshRefinement& MR) { @@ -151,23 +145,19 @@ void outputFinalTetmesh(MeshRefinement& MR) { std::vector &tet_qualities = MR.tet_qualities; int t_cnt = std::count(t_is_removed.begin(), t_is_removed.end(), false); double tmp_time = 0; - if (!args.is_laplacian) { + if (!GArgs::args().is_laplacian) { InoutFiltering IOF(tet_vertices, tets, MR.is_surface_fs, v_is_removed, t_is_removed, tet_qualities); -#ifndef MUTE_COUT igl::Timer igl_timer; igl_timer.start(); -#endif IOF.filter(); t_cnt = std::count(t_is_removed.begin(), t_is_removed.end(), false); -#ifndef MUTE_COUT tmp_time = igl_timer.getElapsedTime(); - cout << "time = " << tmp_time << "s" << endl; - cout << t_cnt << " tets inside!" << endl; -#endif + logger().info("time = {}s", tmp_time); + logger().debug("{} tets inside!", t_cnt); } //output result - cout<<"Writing mesh to "< v_ids; for (int i = 0; i < tets.size(); i++) { if (t_is_removed[i]) @@ -195,12 +185,12 @@ void outputFinalTetmesh(MeshRefinement& MR) { oT(cnt * 4 + j) = map_ids[tets[i][j]]; cnt++; } - cout << "#v = " << oV.rows() / 3 << endl; - cout << "#t = " << oT.rows() / 4 << endl; + logger().debug("#v = {}", oV.rows() / 3); + logger().debug("#t = {}", oT.rows() / 4); - std::string output_format = g_output_file.substr(g_output_file.size() - 4, 4); + std::string output_format = State::state().g_output_file.substr(State::state().g_output_file.size() - 4, 4); if (output_format == "mesh") { - std::fstream f(g_output_file, std::ios::out); + std::fstream f(State::state().g_output_file, std::ios::out); f.precision(std::numeric_limits::digits10 + 1); f << "MeshVersionFormatted 1" << std::endl; f << "Dimension 3" << std::endl; @@ -208,8 +198,8 @@ void outputFinalTetmesh(MeshRefinement& MR) { f << "Vertices" << std::endl << oV.rows() / 3 << std::endl; for (int i = 0; i < oV.rows() / 3; i++) f << oV(i * 3) << " " << oV(i * 3 + 1) << " " << oV(i * 3 + 2) << " " << 0 << std::endl; - f << "Triangles" << endl << 0 < m_vertices; std::vector> m_faces; pp.process(MR.geo_sf_mesh, m_vertices, m_faces); -#ifndef MUTE_COUT tmp_time = igl_timer.getElapsedTime(); addRecord(MeshRecord(MeshRecord::OpType::OP_PREPROCESSING, tmp_time, m_vertices.size(), m_faces.size())); sum_time += tmp_time; - cout << "time = " << tmp_time << "s" << endl; - cout << endl; -#endif + logger().info("time = {}s", tmp_time); //delaunay tetrahedralization -#ifndef MUTE_COUT igl_timer.start(); - cout << "Delaunay tetrahedralizing..." << endl; -#endif + logger().info("Delaunay tetrahedralizing..."); DelaunayTetrahedralization DT; std::vector raw_e_tags; std::vector> raw_conn_e4v; @@ -310,90 +285,67 @@ void gtet_new() { std::vector bsp_faces; std::vector bsp_nodes; DT.tetra(m_vertices, MR.geo_sf_mesh, bsp_vertices, bsp_edges, bsp_faces, bsp_nodes); -#ifndef MUTE_COUT - cout << "# bsp_vertices = " << bsp_vertices.size() << endl; - cout << "# bsp_edges = " << bsp_edges.size() << endl; - cout << "# bsp_faces = " << bsp_faces.size() << endl; - cout << "# bsp_nodes = " << bsp_nodes.size() << endl; - cout << "Delaunay tetrahedralization done!" << endl; + logger().debug("# bsp_vertices = {}", bsp_vertices.size()); + logger().debug("# bsp_edges = {}", bsp_edges.size()); + logger().debug("# bsp_faces = {}", bsp_faces.size()); + logger().debug("# bsp_nodes = {}", bsp_nodes.size()); + logger().info("Delaunay tetrahedralization done!"); tmp_time = igl_timer.getElapsedTime(); addRecord(MeshRecord(MeshRecord::OpType::OP_DELAUNEY_TETRA, tmp_time, bsp_vertices.size(), bsp_nodes.size())); sum_time += tmp_time; - cout << "time = " << tmp_time << "s" << endl; - cout << endl; -#endif + logger().info("time = {}s", tmp_time); //mesh conforming -#ifndef MUTE_COUT igl_timer.start(); - cout << "Divfaces matching..." << endl; -#endif + logger().info("Divfaces matching..."); MeshConformer MC(m_vertices, m_faces, bsp_vertices, bsp_edges, bsp_faces, bsp_nodes); MC.match(); -#ifndef MUTE_COUT - cout << "Divfaces matching done!" << endl; + logger().info("Divfaces matching done!"); tmp_time = igl_timer.getElapsedTime(); addRecord(MeshRecord(MeshRecord::OpType::OP_DIVFACE_MATCH, tmp_time, bsp_vertices.size(), bsp_nodes.size())); - cout << "time = " << tmp_time << "s" << endl; - cout << endl; -#endif + logger().info("time = {}s", tmp_time); //bsp subdivision -#ifndef MUTE_COUT igl_timer.start(); - cout << "BSP subdivision ..." << endl; -#endif + logger().info("BSP subdivision ..."); BSPSubdivision BS(MC); BS.init(); BS.subdivideBSPNodes(); -#ifndef MUTE_COUT - cout << "Output: " << endl; - cout << "# node = " << MC.bsp_nodes.size() << endl; - cout << "# face = " << MC.bsp_faces.size() << endl; - cout << "# edge = " << MC.bsp_edges.size() << endl; - cout << "# vertex = " << MC.bsp_vertices.size() << endl; - cout << "BSP subdivision done!" << endl; + logger().debug("Output: "); + logger().debug("# node = {}", MC.bsp_nodes.size()); + logger().debug("# face = {}", MC.bsp_faces.size()); + logger().debug("# edge = {}", MC.bsp_edges.size()); + logger().debug("# vertex = {}", MC.bsp_vertices.size()); + logger().info("BSP subdivision done!"); tmp_time = igl_timer.getElapsedTime(); addRecord(MeshRecord(MeshRecord::OpType::OP_BSP, tmp_time, bsp_vertices.size(), bsp_nodes.size())); sum_time += tmp_time; - cout << "time = " << tmp_time << "s" << endl; - cout << endl; -#endif + logger().info("time = {}s", tmp_time); //simple tetrahedralization -#ifndef MUTE_COUT igl_timer.start(); - cout << "Tetrehedralizing ..." << endl; -#endif + logger().info("Tetrehedralizing ..."); SimpleTetrahedralization ST(MC); ST.tetra(MR.tet_vertices, MR.tets); ST.labelSurface(m_f_tags, raw_e_tags, raw_conn_e4v, MR.tet_vertices, MR.tets, MR.is_surface_fs); ST.labelBbox(MR.tet_vertices, MR.tets); - if (!g_is_close)//if input is an open mesh + if (!State::state().g_is_close)//if input is an open mesh ST.labelBoundary(MR.tet_vertices, MR.tets, MR.is_surface_fs); -#ifndef MUTE_COUT - cout << "# tet_vertices = " << MR.tet_vertices.size() << endl; - cout << "# tets = " << MR.tets.size() << endl; - cout << "Tetrahedralization done!" << endl; + logger().debug("# tet_vertices = {}", MR.tet_vertices.size()); + logger().debug("# tets = {}", MR.tets.size()); + logger().info("Tetrahedralization done!"); tmp_time = igl_timer.getElapsedTime(); addRecord(MeshRecord(MeshRecord::OpType::OP_SIMPLE_TETRA, tmp_time, MR.tet_vertices.size(), MR.tets.size())); sum_time += tmp_time; - cout << "time = " << tmp_time << "s" << endl; - cout << endl; - cout << "Total time for the first stage = " << sum_time << endl; -#endif + logger().info("time = {}s", tmp_time); + logger().info("Total time for the first stage = {}", sum_time); } /// STAGE 2 //init -#ifndef MUTE_COUT - cout << "Refinement initializing..." << endl; -#endif + logger().info("Refinement initializing..."); MR.prepareData(); -#ifndef MUTE_COUT - cout << "Refinement intialization done!" << endl; - cout << endl; -#endif + logger().info("Refinement intialization done!"); //improvement MR.refine(energy_type); @@ -408,29 +360,34 @@ void gtet_new_slz(const std::string& sf_file, const std::string& slz_file, int m MR.deserialization(sf_file, slz_file); // MR.is_dealing_unrounded = true; - MR.refine(ENERGY_AMIPS, ops, false, true); + MR.refine(State::state().ENERGY_AMIPS, ops, false, true); outputFinalTetmesh(MR); } int main(int argc, char *argv[]) { #ifdef MUTE_COUT - cout<<"Unnecessary checks are muted."<required(); - app.add_option("--postfix", args.postfix, "--postfix P. Postfix P for output files. (string, optinal, default: '_')"); - app.add_option("--output", args.output, "--output OUTPUT. Output tetmesh OUTPUT in .msh format. (string, optional, default: input_file+postfix+'.msh')"); - app.add_option("--ideal-edge-length", args.i_ideal_edge_length, "--ideal-edge-length L. ideal_edge_length = diag_of_bbox / L. (double, optional, default: 20)"); - app.add_option("--epsilon", args.i_epsilon, "--epsilon EPS. epsilon = diag_of_bbox / EPS. (double, optional, default: 1000)"); - app.add_option("--stage", args.stage, "--stage STAGE. Run pipeline in stage STAGE. (integer, optional, default: 1)"); - app.add_option("--filter-energy", args.filter_energy, "--filter-energy ENERGY. Stop mesh improvement when the maximum energy is smaller than ENERGY. (double, optional, default: 10)"); - app.add_option("--max-pass", args.max_pass, "--max-pass PASS. Do PASS mesh improvement passes in maximum. (integer, optional, default: 80)"); - - app.add_option("--is-laplacian", args.is_laplacian, "--is-laplacian ISLAP. Do Laplacian smoothing for the surface of output on the holes of input, if ISLAP = 1. Otherwise, ISLAP = 0. (integer, optinal, default: 0)"); - app.add_option("--targeted-num-v", args.targeted_num_v, "--targeted-num-v TV. Output tetmesh that contains TV vertices. (integer, optinal, tolerance: 5%)"); - app.add_option("--bg-mesh", args.bg_mesh, "--bg-mesh BGMESH. Background tetmesh BGMESH in .msh format for applying sizing field. (string, optional)"); - app.add_option("--is-quiet", args.is_quiet, "--is-quiet Q. Mute log info and only output tetmesh if Q = 1. (integer, optional, default: 0)"); + app.add_option("input,--input", GArgs::args().input, "Input surface mesh INPUT in .off/.obj/.stl/.ply format. (string, required)")->required(); + app.add_option("output,--output", GArgs::args().output, "Output tetmesh OUTPUT in .msh format. (string, optional, default: input_file+postfix+'.msh')"); + app.add_option("--postfix", GArgs::args().postfix, "Postfix P for output files. (string, optional, default: '_')"); + app.add_option("-l,--ideal-edge-length", GArgs::args().i_ideal_edge_length, "ideal_edge_length = diag_of_bbox / L. (double, optional, default: 20)"); + app.add_option("-e,--epsilon", GArgs::args().i_epsilon, "epsilon = diag_of_bbox / EPS. (double, optional, default: 1000)"); + app.add_option("--stage", GArgs::args().stage, "Run pipeline in stage STAGE. (integer, optional, default: 1)"); + app.add_option("--filter-energy", GArgs::args().filter_energy, "Stop mesh improvement when the maximum energy is smaller than ENERGY. (double, optional, default: 10)"); + app.add_option("--max-pass", GArgs::args().max_pass, "Do PASS mesh improvement passes in maximum. (integer, optional, default: 80)"); + + app.add_option("--is-laplacian", GArgs::args().is_laplacian, "Do Laplacian smoothing for the surface of output on the holes of input, if ISLAP = 1. Otherwise, ISLAP = 0. (integer, optinal, default: 0)"); + app.add_option("--targeted-num-v", GArgs::args().targeted_num_v, "Output tetmesh that contains TV vertices. (integer, optinal, tolerance: 5%)"); + app.add_option("--bg-mesh", GArgs::args().bg_mesh, "Background tetmesh BGMESH in .msh format for applying sizing field. (string, optional)"); + app.add_option("-q,--is-quiet", GArgs::args().is_quiet, "Mute console output. (integer, optional, default: 0)"); + app.add_option("--log", log_filename, "Log info to given file."); + app.add_option("--level", log_level, "Log level (0 = most verbose, 6 = off)."); try { app.parse(argc, argv); @@ -438,36 +395,46 @@ int main(int argc, char *argv[]) { return app.exit(e); } + Logger::init(!GArgs::args().is_quiet, log_filename); + log_level = std::max(0, std::min(6, log_level)); + spdlog::set_level(static_cast(log_level)); + spdlog::flush_every(std::chrono::seconds(3)); + + // logger().info("this is a test"); + // logger().debug("debug stuff"); + //initalization GEO::initialize(); - g_postfix = args.postfix; - if(args.slz_file != "") - g_working_dir = args.input.substr(0, args.slz_file.size() - 4); + State::state().g_postfix = GArgs::args().postfix; + if(GArgs::args().slz_file != "") + State::state().g_working_dir = GArgs::args().input.substr(0, GArgs::args().slz_file.size() - 4); else - g_working_dir = args.input.substr(0, args.input.size() - 4); + State::state().g_working_dir = GArgs::args().input.substr(0, GArgs::args().input.size() - 4); - if(args.csv_file == "") - g_stat_file = g_working_dir + g_postfix + ".csv"; + if(GArgs::args().csv_file == "") + State::state().g_stat_file = State::state().g_working_dir + State::state().g_postfix + ".csv"; else - g_stat_file = args.csv_file; + State::state().g_stat_file = GArgs::args().csv_file; - if(args.output == "") - g_output_file = g_working_dir + g_postfix + ".msh"; + if(GArgs::args().output == "") + State::state().g_output_file = State::state().g_working_dir + State::state().g_postfix + ".msh"; else - g_output_file = args.output; + State::state().g_output_file = GArgs::args().output; - if(args.is_quiet) { - args.is_output_csv = false; - std::streambuf* orig_buf = cout.rdbuf(); - cout.rdbuf(NULL); -// cout.setstate(std::ios_base::failbit);//use std::cout.clear() to get it back + if(GArgs::args().is_quiet) { + GArgs::args().is_output_csv = false; + std::streambuf* orig_buf = std::cout.rdbuf(); + std::cout.rdbuf(NULL); +// std::cout.setstate(std::ios_base::failbit);//use std::std::cout.clear() to get it back } //do tetrahedralization - if(args.slz_file != "") - gtet_new_slz(args.input, args.slz_file, args.max_pass, std::array({true, false, true, true})); + if(GArgs::args().slz_file != "") + gtet_new_slz(GArgs::args().input, GArgs::args().slz_file, GArgs::args().max_pass, std::array({true, false, true, true})); else gtet_new(); + spdlog::shutdown(); + return 0; } diff --git a/src/tetwild/BSPElements.h b/src/tetwild/BSPElements.h index 59c6f71..1bd3545 100644 --- a/src/tetwild/BSPElements.h +++ b/src/tetwild/BSPElements.h @@ -1,9 +1,9 @@ // This file is part of TetWild, a software for generating tetrahedral meshes. -// +// // Copyright (C) 2018 Yixin Hu -// -// This Source Code Form is subject to the terms of the Mozilla Public License -// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can // obtain one at http://mozilla.org/MPL/2.0/. // // Created by yihu on 8/22/17. @@ -14,6 +14,8 @@ #include #include +namespace tetwild { + class BSPEdge{ public: std::vector vertices; @@ -42,4 +44,6 @@ class BSPtreeNode{ std::unordered_set div_faces; }; +} // namespace tetwild + #endif //NEW_GTET_BSPELEMENTS_H diff --git a/src/tetwild/BSPSubdivision.cpp b/src/tetwild/BSPSubdivision.cpp index 24b6df0..9ea62f6 100644 --- a/src/tetwild/BSPSubdivision.cpp +++ b/src/tetwild/BSPSubdivision.cpp @@ -1,9 +1,9 @@ // This file is part of TetWild, a software for generating tetrahedral meshes. -// +// // Copyright (C) 2018 Yixin Hu -// -// This Source Code Form is subject to the terms of the Mozilla Public License -// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can // obtain one at http://mozilla.org/MPL/2.0/. // // Created by Yixin Hu on 4/3/17. @@ -11,6 +11,8 @@ #include +namespace tetwild { + void BSPSubdivision::init() { for (int old_n_id = 0; old_n_id < MC.bsp_nodes.size(); old_n_id++) { if (MC.bsp_nodes[old_n_id].div_faces.size() == 0) { @@ -26,10 +28,7 @@ void BSPSubdivision::init() { } nf+=1; } -#ifndef MUTE_COUT - cout << "# nodes need subdivision = " << nf<<"/"<< processing_n_ids.size() << "/" - << MC.bsp_nodes.size() << endl; -#endif + logger().debug("# nodes need subdivision = {}/{}/{}", nf, processing_n_ids.size(), MC.bsp_nodes.size()); } void BSPSubdivision::subdivideBSPNodes() { @@ -205,7 +204,7 @@ void BSPSubdivision::subdivideBSPNodes() { v_sides[new_v_id] = V_ON;//fixed } else { - cout << "error cal p!" << endl; + logger().debug("error cal p!"); exit(250); } @@ -379,4 +378,6 @@ void BSPSubdivision::getVertices(BSPFace& face){ } for(auto it=vs.begin();it!=vs.end();it++) face.vertices.push_back(*it); -} \ No newline at end of file +} + +} // namespace tetwild diff --git a/src/tetwild/BSPSubdivision.h b/src/tetwild/BSPSubdivision.h index b74e713..6f903c1 100644 --- a/src/tetwild/BSPSubdivision.h +++ b/src/tetwild/BSPSubdivision.h @@ -1,9 +1,9 @@ // This file is part of TetWild, a software for generating tetrahedral meshes. -// +// // Copyright (C) 2018 Yixin Hu -// -// This Source Code Form is subject to the terms of the Mozilla Public License -// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can // obtain one at http://mozilla.org/MPL/2.0/. // // Created by Yixin Hu on 4/3/17. @@ -13,6 +13,9 @@ #define NEW_GTET_BSPSUBDIVISION_H #include +#include + +namespace tetwild { class BSPSubdivision { public: @@ -38,5 +41,6 @@ class BSPSubdivision { }; +} // namespace tetwild #endif //NEW_GTET_BSPSUBDIVISION_H diff --git a/src/tetwild/CGALCommon.cpp b/src/tetwild/CGALCommon.cpp new file mode 100644 index 0000000..e69de29 diff --git a/src/tetwild/CGALCommon.h b/src/tetwild/CGALCommon.h new file mode 100644 index 0000000..655ca16 --- /dev/null +++ b/src/tetwild/CGALCommon.h @@ -0,0 +1,121 @@ +// This file is part of TetWild, a software for generating tetrahedral meshes. +// +// Copyright (C) 2018 Jeremie Dumas +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// obtain one at http://mozilla.org/MPL/2.0/. +// +// Created by Jeremie Dumas on 09/04/18. +// + +#pragma once + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +namespace tetwild { + +typedef CGAL::Exact_predicates_exact_constructions_kernel K; +typedef K::Point_2 Point_2; +typedef K::Segment_2 Segment_2; +typedef K::Line_2 Line_2; +typedef K::Iso_rectangle_2 Iso_rectangle_2; +typedef K::Triangle_2 Triangle_2; +typedef K::Intersect_2 Intersect_2; +//typedef CGAL::Polygon_2 Polygon_2; + +typedef K::Point_3 Point_3; +typedef K::Vector_3 Vector_3; +typedef K::Segment_3 Segment_3; +typedef K::Line_3 Line_3; +typedef K::Plane_3 Plane_3; +typedef K::Triangle_3 Triangle_3; +typedef K::Intersect_3 Intersect_3; +typedef K::Tetrahedron_3 Tetrahedron_3; +typedef K::Direction_3 Direction_3; + +typedef CGAL::Exact_predicates_inexact_constructions_kernel Kf; +typedef Kf::Point_3 Point_3f; +typedef Kf::Vector_3 Vector_3f; +typedef Kf::Plane_3 Plane_3f; +typedef Kf::Triangle_3 Triangle_3f; +typedef Kf::Segment_3 Segment_3f; +typedef Kf::Line_3 Line_3f; + +typedef CGAL::Epeck::FT CGAL_FT; +//#include +//typedef CGAL::Simple_cartesian::FT CGAL_FT; + +typedef K::Iso_cuboid_3 Bbox_3; + +} // namespace tetwild + +//for serialization +namespace igl { + namespace serialization { + template<> + inline void serialize(const tetwild::Point_3 &p, std::vector &buffer) { + ::igl::serialize(STR(CGAL::exact(p[0])), std::string("x"), buffer); + ::igl::serialize(STR(CGAL::exact(p[1])), std::string("y"), buffer); + ::igl::serialize(STR(CGAL::exact(p[2])), std::string("z"), buffer); + } + + template<> + inline void deserialize(tetwild::Point_3 &p, const std::vector &buffer) { + std::string s1, s2, s3; + ::igl::deserialize(s1, std::string("x"), buffer); + ::igl::deserialize(s2, std::string("y"), buffer); + ::igl::deserialize(s3, std::string("z"), buffer); +// p=Point_3(CGAL_FT(s1), CGAL_FT(s2), CGAL_FT(s3)); + } + + template<> + inline void serialize(const tetwild::Point_3f &p, std::vector &buffer) { + ::igl::serialize(p[0], std::string("x"), buffer); + ::igl::serialize(p[1], std::string("y"), buffer); + ::igl::serialize(p[2], std::string("z"), buffer); + } + + template<> + inline void deserialize(tetwild::Point_3f &p, const std::vector &buffer) { + double x, y, z; + ::igl::deserialize(x, std::string("x"), buffer); + ::igl::deserialize(y, std::string("y"), buffer); + ::igl::deserialize(z, std::string("z"), buffer); + p=tetwild::Point_3f(x, y, z); + } + + template<> + inline void serialize(const std::array &arr, std::vector &buffer) { + for(int i=0;i<3;i++) + ::igl::serialize(arr[i], std::to_string(i), buffer); + } + + template<> + inline void deserialize(std::array &arr, const std::vector &buffer) { + for(int i=0;i<3;i++) + ::igl::deserialize(arr[i], std::to_string(i), buffer); + } + + template<> + inline void serialize(const std::array &arr, std::vector &buffer) { + for(int i=0;i<4;i++) + ::igl::serialize(arr[i], std::to_string(i), buffer); + } + + template<> + inline void deserialize(std::array &arr, const std::vector &buffer) { + for(int i=0;i<4;i++) + ::igl::deserialize(arr[i], std::to_string(i), buffer); + } + } +} diff --git a/src/tetwild/heads.cpp b/src/tetwild/Common.cpp similarity index 76% rename from src/tetwild/heads.cpp rename to src/tetwild/Common.cpp index ba699f8..a8b124b 100644 --- a/src/tetwild/heads.cpp +++ b/src/tetwild/Common.cpp @@ -1,54 +1,30 @@ // This file is part of TetWild, a software for generating tetrahedral meshes. -// -// Copyright (C) 2018 Yixin Hu -// -// This Source Code Form is subject to the terms of the Mozilla Public License -// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// +// Copyright (C) 2018 Jeremie Dumas +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can // obtain one at http://mozilla.org/MPL/2.0/. // -// Created by Yixin Hu on 3/29/17. +// Created by Jeremie Dumas on 09/04/18. // -#include - -//Eigen::MatrixXd G; -std::string g_working_dir = ""; -std::string g_stat_file = ""; -std::string g_postfix = ""; -std::string g_output_file = ""; - -bool is_using_energy_max = true; - -bool is_using_sampling = true; -bool is_use_project = false; +#include +#include +#include -int NOT_SURFACE = 0; +namespace tetwild { -double g_eps = 0; -double g_eps_2 = 0; -double g_dd = 0; -double g_ideal_l = 0; -double g_diag_l = 0; -bool g_is_close = true; - -double g_eps_input = 0; -double g_eps_delta = 0; -int g_cur_stage = 1; - -bool is_print_tmp = false; - -Args args; -//std::vector mesh_records = std::vector(); -bool is_app_csv=false; void addRecord(const MeshRecord& record) { - if (!args.is_output_csv) + if (!GArgs::args().is_output_csv) return; + static bool first_time = true; std::ofstream f; - if (is_app_csv) - f.open(g_stat_file, std::ios::app); - else { - f.open(g_stat_file); - is_app_csv = true; + if (first_time) { + f.open(State::state().g_stat_file); + first_time = false; + } else { + f.open(State::state().g_stat_file, std::ios::app); } f << record.op << "," << record.timing << "," << record.n_v << "," << record.n_t << "," << record.min_min_d_angle << "," << record.avg_min_d_angle << "," @@ -58,9 +34,9 @@ void addRecord(const MeshRecord& record) { } void pausee(){ - cout<<"Is pausing... (Enter '0' to exit and other characters to continue.)"<>c; + std::cin>>c; if(c=='0') exit(0); } @@ -134,7 +110,7 @@ void sampleTriangle(const std::array& vs, std::vector& auto min_max = std::minmax_element(ls.begin(), ls.end()); int min_i = min_max.first - ls.begin(); int max_i = min_max.second - ls.begin(); - double N = sqrt(ls[max_i]) / g_dd; + double N = sqrt(ls[max_i]) / State::state().g_dd; if (N <= 1) { for (int i = 0; i < 3; i++) ps.push_back(vs[i]); @@ -149,12 +125,12 @@ void sampleTriangle(const std::array& vs, std::vector& GEO::vec3 n_v0v1 = GEO::normalize(v1 - v0); for (int n = 0; n <= N; n++) { - ps.push_back(v0 + n_v0v1 * g_dd * n); + ps.push_back(v0 + n_v0v1 * State::state().g_dd * n); } ps.push_back(v1); double h = GEO::distance(GEO::dot((v2 - v0), (v1 - v0)) * (v1 - v0) / ls[max_i] + v0, v2); - int M = h / (sqrt3_2 * g_dd); + int M = h / (sqrt3_2 * State::state().g_dd); if (M < 1) { ps.push_back(v2); return; @@ -174,45 +150,45 @@ void sampleTriangle(const std::array& vs, std::vector& if (m % 2 == 0 && n == n1) { n += 1; } - GEO::vec3 v0_m = v0 + m * sqrt3_2 * g_dd / sin_v0 * n_v0v2; - GEO::vec3 v1_m = v1 + m * sqrt3_2 * g_dd / sin_v1 * n_v1v2; - if (GEO::distance(v0_m, v1_m) <= g_dd) + GEO::vec3 v0_m = v0 + m * sqrt3_2 * State::state().g_dd / sin_v0 * n_v0v2; + GEO::vec3 v1_m = v1 + m * sqrt3_2 * State::state().g_dd / sin_v1 * n_v1v2; + if (GEO::distance(v0_m, v1_m) <= State::state().g_dd) break; - double delta_d = ((n + (m % 2) / 2.0) - m * sqrt3_2 / tan_v0) * g_dd; + double delta_d = ((n + (m % 2) / 2.0) - m * sqrt3_2 / tan_v0) * State::state().g_dd; GEO::vec3 v = v0_m + delta_d * n_v0v1; - int N1 = GEO::distance(v, v1_m) / g_dd; + int N1 = GEO::distance(v, v1_m) / State::state().g_dd; // ps.push_back(v0_m); for (int i = 0; i <= N1; i++) { - ps.push_back(v + i * n_v0v1 * g_dd); + ps.push_back(v + i * n_v0v1 * State::state().g_dd); } // ps.push_back(v1_m); } ps.push_back(v2); //sample edges - N = sqrt(ls[(max_i + 1) % 3]) / g_dd; + N = sqrt(ls[(max_i + 1) % 3]) / State::state().g_dd; if (N > 1) { if (N == int(N)) N -= 1; GEO::vec3 n_v1v2 = GEO::normalize(v2 - v1); for (int n = 1; n <= N; n++) { - ps.push_back(v1 + n_v1v2 * g_dd * n); + ps.push_back(v1 + n_v1v2 * State::state().g_dd * n); } } - N = sqrt(ls[(max_i + 2) % 3]) / g_dd; + N = sqrt(ls[(max_i + 2) % 3]) / State::state().g_dd; if (N > 1) { if (N == int(N)) N -= 1; GEO::vec3 n_v2v0 = GEO::normalize(v0 - v2); for (int n = 1; n <= N; n++) { - ps.push_back(v2 + n_v2v0 * g_dd * n); + ps.push_back(v2 + n_v2v0 * State::state().g_dd * n); } } -// cout << "ps.size = " << ps.size() << endl; -// cout << "is output samples?" << endl; +// logger().debug("ps.size = {}", ps.size()); +// logger().debug("is output samples?"); // int anw = 0; // cin >> anw; // if (anw != 0) { @@ -232,6 +208,8 @@ void sampleTriangle(const std::array& vs, std::vector& // F_tmp(1 + i, k) = (1 + i) * 3 + k; // } // } -// igl::writeSTL(g_working_dir + "_sample.stl", V_tmp, F_tmp); +// igl::writeSTL(State::state().g_working_dir + "_sample.stl", V_tmp, F_tmp); // } -} \ No newline at end of file +} + +} // namespace tetwild diff --git a/src/tetwild/Common.h b/src/tetwild/Common.h new file mode 100644 index 0000000..6a3ab61 --- /dev/null +++ b/src/tetwild/Common.h @@ -0,0 +1,35 @@ +// This file is part of TetWild, a software for generating tetrahedral meshes. +// +// Copyright (C) 2018 Jeremie Dumas +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// obtain one at http://mozilla.org/MPL/2.0/. +// +// Created by Jeremie Dumas on 09/04/18. +// + +#pragma once + +#include +#include +#include +#include +#include +#include +#include + +#define TIMING_BREAKDOWN true + +namespace tetwild { + +void pausee(); + +bool isHaveCommonEle(const std::unordered_set& v1, const std::unordered_set& v2); +void setIntersection(const std::unordered_set& s1, const std::unordered_set& s2, std::unordered_set& s); +void setIntersection(const std::unordered_set& s1, const std::unordered_set& s2, std::vector& s); +void sampleTriangle(const std::array& vs, std::vector& ps); + +void addRecord(const MeshRecord& record); + +} // namespace tetwild diff --git a/src/tetwild/DelaunayTetrahedralization.cpp b/src/tetwild/DelaunayTetrahedralization.cpp index 0d26bb4..0d40c56 100644 --- a/src/tetwild/DelaunayTetrahedralization.cpp +++ b/src/tetwild/DelaunayTetrahedralization.cpp @@ -1,9 +1,9 @@ // This file is part of TetWild, a software for generating tetrahedral meshes. -// +// // Copyright (C) 2018 Yixin Hu -// -// This Source Code Form is subject to the terms of the Mozilla Public License -// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can // obtain one at http://mozilla.org/MPL/2.0/. // // Created by Yixin Hu on 3/29/17. @@ -15,14 +15,14 @@ #include #include #include - #include #include #include +#include #include -//#include +#include -#include +namespace tetwild { void DelaunayTetrahedralization::init(const std::vector& m_vertices, const std::vector>& m_faces, std::vector& m_f_tags, std::vector& raw_e_tags, std::vector>& raw_conn_e4v) { @@ -56,10 +56,10 @@ void DelaunayTetrahedralization::getVoxelPoints(const Point_3& p_min, const Poin GEO::MeshFacetsAABB geo_face_tree(geo_surface_mesh); double voxel_resolution; - if(args.i_ideal_edge_length > 20) - voxel_resolution = g_diag_l / 20.0; + if(GArgs::args().i_ideal_edge_length > 20) + voxel_resolution = State::state().g_diag_l / 20.0; else - voxel_resolution = g_diag_l / args.i_ideal_edge_length; + voxel_resolution = State::state().g_diag_l / GArgs::args().i_ideal_edge_length; std::array d; std::array N; for (int i = 0; i < 3; i++) { @@ -77,7 +77,7 @@ void DelaunayTetrahedralization::getVoxelPoints(const Point_3& p_min, const Poin } double min_dis = voxel_resolution * voxel_resolution / 4; -// double min_dis = g_ideal_l * g_ideal_l;//epsilon*2 +// double min_dis = State::state().g_ideal_l * State::state().g_ideal_l;//epsilon*2 for (int i = 0; i < ds[0].size(); i++) { for (int j = 0; j < ds[1].size(); j++) { for (int k = 0; k < ds[2].size(); k++) { @@ -108,16 +108,16 @@ void DelaunayTetrahedralization::tetra(const std::vector& m_vertices, G Point_3 p_min = bbox.min(); Point_3 p_max = bbox.max(); - double dis = g_eps * 2;//todo: use epsilon to determine the size of bbx - if (dis < g_diag_l / 20) - dis = g_diag_l / 20; + double dis = State::state().g_eps * 2;//todo: use epsilon to determine the size of bbx + if (dis < State::state().g_diag_l / 20) + dis = State::state().g_diag_l / 20; else - dis = g_eps * 1.1; + dis = State::state().g_eps * 1.1; p_min = Point_3(p_min[0] - dis, p_min[1] - dis, p_min[2] - dis); p_max = Point_3(p_max[0] + dis, p_max[1] + dis, p_max[2] + dis); -// cout<<"p_min: "< p; @@ -132,18 +132,16 @@ void DelaunayTetrahedralization::tetra(const std::vector& m_vertices, G } ///add voxel points std::vector voxel_points; - if(args.is_using_voxel) + if(GArgs::args().is_using_voxel) getVoxelPoints(p_min, p_max, geo_surface_mesh, voxel_points); for(int i=0;i& m_vertices, G void DelaunayTetrahedralization::outputTetmesh(const std::vector& m_vertices, std::vector>& cells, const std::string& output_file){ - std::streambuf* coutBuf = std::cout.rdbuf(); std::ofstream of(output_file); - std::cout.rdbuf(of.rdbuf()); - cout< -// -// This Source Code Form is subject to the terms of the Mozilla Public License -// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can // obtain one at http://mozilla.org/MPL/2.0/. // // Created by Yixin Hu on 3/29/17. @@ -12,11 +12,14 @@ #ifndef GTET_DELAUNAYTETRAHEDRALIZATION_H #define GTET_DELAUNAYTETRAHEDRALIZATION_H -#include - +#include #include +#include #include #include + +namespace tetwild { + typedef CGAL::Triangulation_vertex_base_with_info_3 Vb; typedef CGAL::Triangulation_data_structure_3 Tds; typedef CGAL::Delaunay_triangulation_3 Delaunay; @@ -39,5 +42,6 @@ class DelaunayTetrahedralization { const std::string& output_file); }; +} // namespace tetwild #endif //GTET_DELAUNAYTETRAHEDRALIZATION_H diff --git a/src/tetwild/EdgeCollapser.cpp b/src/tetwild/EdgeCollapser.cpp index bf8f285..1efad27 100644 --- a/src/tetwild/EdgeCollapser.cpp +++ b/src/tetwild/EdgeCollapser.cpp @@ -1,16 +1,19 @@ // This file is part of TetWild, a software for generating tetrahedral meshes. -// +// // Copyright (C) 2018 Yixin Hu -// -// This Source Code Form is subject to the terms of the Mozilla Public License -// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can // obtain one at http://mozilla.org/MPL/2.0/. // // Created by Yixin Hu on 4/11/17. // #include - +#include + +namespace tetwild { + void EdgeCollapser::init() { energy_time = 0; @@ -64,9 +67,7 @@ void EdgeCollapser::init() { void EdgeCollapser::collapse() { tet_tss.assign(tets.size(), 0); int cnt = 0; -#ifndef MUTE_COUT - cout << "edge queue size = " << ec_queue.size() << endl; -#endif + logger().debug("edge queue size = {}", ec_queue.size()); while (!ec_queue.empty()) { std::array v_ids = ec_queue.top().v_ids; double old_weight = ec_queue.top().weight; @@ -140,15 +141,13 @@ void EdgeCollapser::collapse() { counter++; } -#ifndef MUTE_COUT - cout << suc_counter << " " << counter << " " << inf_es.size() << endl; - cout << "envelop accept = " << envelop_accept_cnt << endl; -#endif + logger().debug("{} {} {}", suc_counter, counter, inf_es.size()); + logger().debug("envelop accept = {}", envelop_accept_cnt); if (suc_counter == 0 || inf_es.size() == 0) { -// cout<<"checking......."< v_ids = ec_queue.top().v_ids; @@ -185,27 +184,25 @@ void EdgeCollapser::collapse() { // else // cnt_suc++; // } -// cout< tet_qs; @@ -383,14 +378,14 @@ int EdgeCollapser::collapseAnEdge(int v1_id, int v2_id) { calTetQualities(new_tets, tet_qs); energy_time+=tmp_timer.getElapsedTime(); - if (energy_type != ENERGY_NA && is_check_quality) { + if (energy_type != State::state().ENERGY_NA && is_check_quality) { TetQuality old_tq, new_tq; getCheckQuality(old_t_ids, old_tq); getCheckQuality(tet_qs, new_tq); // if (is_edge_too_short) -// cout << "old " << old_tq.slim_energy << " new " << new_tq.slim_energy << endl; +// logger().debug("old {} new {}", old_tq.slim_energy, new_tq.slim_energy); // if (is_soft && old_tq.slim_energy < soft_energy) { -// old_tq.slim_energy = args.filter_energy; +// old_tq.slim_energy = GArgs::args().filter_energy; // } if(is_soft) old_tq.slim_energy = soft_energy; @@ -398,7 +393,7 @@ int EdgeCollapser::collapseAnEdge(int v1_id, int v2_id) { new_tq.slim_energy = 0; if (!is_edge_degenerate && !new_tq.isBetterOrEqualThan(old_tq, energy_type)) { // if (is_edge_too_short) -// cout << "quality" << endl; +// logger().debug("quality"); return QUALITY; } } @@ -413,7 +408,7 @@ int EdgeCollapser::collapseAnEdge(int v1_id, int v2_id) { tet_vertices[v1_id].posf = old_pf; tet_vertices[v1_id].pos = old_p; // if (is_edge_too_short) -// cout << "boundary" << endl; +// logger().debug("boundary"); return ENVELOP; } tet_vertices[v1_id].posf = old_pf; @@ -422,24 +417,22 @@ int EdgeCollapser::collapseAnEdge(int v1_id, int v2_id) { //check 3 bool is_envelop_suc = false; - if (g_eps != EPSILON_NA && g_eps != EPSILON_INFINITE && tet_vertices[v1_id].is_on_surface) { + if (State::state().g_eps != State::state().EPSILON_NA && State::state().g_eps != State::state().EPSILON_INFINITE && tet_vertices[v1_id].is_on_surface) { if (!is_edge_degenerate && !isCollapsable_epsilon(v1_id, v2_id)) { // if (is_edge_too_short) -// cout << "envelop" << endl; +// logger().debug("envelop"); return ENVELOP; } is_envelop_suc = true; envelop_accept_cnt++; -#ifndef MUTE_COUT if (envelop_accept_cnt % 1000 == 0) - cout << "1000 accepted!" << endl; -#endif + logger().debug("1000 accepted!"); } //real update // if(is_edge_too_short) -// cout<<"success"<"; +// logger().debug("{}{}jt==tri.end()", n1_v_ids.size(), "->"; std::vector n1_v_ids_vec, n12_v_ids_vec; n1_v_ids_vec.reserve(n1_v_ids.size()); n12_v_ids_vec.reserve(n12_v_ids.size()); @@ -664,7 +657,7 @@ bool EdgeCollapser::isCollapsable_cd1(int v1_id, int v2_id) { } //check the surface tags //if the vertex is on the surface -// if (g_eps != EPSILON_NA) { +// if (State::state().g_eps != State::state().EPSILON_NA) { // return true; // } return true; @@ -719,11 +712,11 @@ bool EdgeCollapser::isCollapsable_epsilon(int v1_id, int v2_id) { // std::vector tris; // for (auto it = tet_vertices[v1_id].conn_tets.begin(); it != tet_vertices[v1_id].conn_tets.end(); it++) { // for (int j = 0; j < 4; j++) { -// if (tets[*it][j] == v2_id && is_surface_fs[*it][j] != NOT_SURFACE) { +// if (tets[*it][j] == v2_id && is_surface_fs[*it][j] != State::state().NOT_SURFACE) { // std::array tri = {tets[*it][(j + 1) % 4], tets[*it][(j + 2) % 4], tets[*it][(j + 3) % 4]}; // auto jt = std::find(tri.begin(), tri.end(), v1_id); // if(jt==tri.end()){ -// cout<<"jt==tri.end()"<> tri_ids; for (auto it = tet_vertices[v1_id].conn_tets.begin(); it != tet_vertices[v1_id].conn_tets.end(); it++) { for (int j = 0; j < 4; j++) { - if (tets[*it][j] != v1_id && is_surface_fs[*it][j] != NOT_SURFACE) { + if (tets[*it][j] != v1_id && is_surface_fs[*it][j] != State::state().NOT_SURFACE) { std::array tri = {tets[*it][(j + 1) % 4], tets[*it][(j + 2) % 4], tets[*it][(j + 3) % 4]}; std::sort(tri.begin(), tri.end()); tri_ids.push_back(tri); @@ -789,3 +782,5 @@ bool EdgeCollapser::isEdgeValid(const std::array& e){ // ec_queue.push(ele); // } //} + +} // namespace tetwild diff --git a/src/tetwild/EdgeCollapser.h b/src/tetwild/EdgeCollapser.h index c4c96f6..2477334 100644 --- a/src/tetwild/EdgeCollapser.h +++ b/src/tetwild/EdgeCollapser.h @@ -14,6 +14,8 @@ #include +namespace tetwild { + class ElementInQueue_ec{ public: std::array v_ids; @@ -94,5 +96,6 @@ class EdgeCollapser: public LocalOperations { igl::Timer igl_timer; }; +} // namespace tetwild #endif //NEW_GTET_EDGECOLLAPSER_H diff --git a/src/tetwild/EdgeRemover.cpp b/src/tetwild/EdgeRemover.cpp index 3c0949f..de9cb6c 100644 --- a/src/tetwild/EdgeRemover.cpp +++ b/src/tetwild/EdgeRemover.cpp @@ -1,15 +1,18 @@ // This file is part of TetWild, a software for generating tetrahedral meshes. -// +// // Copyright (C) 2018 Yixin Hu -// -// This Source Code Form is subject to the terms of the Mozilla Public License -// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can // obtain one at http://mozilla.org/MPL/2.0/. // // Created by Yixin Hu on 4/17/17. // #include +#include + +namespace tetwild { void EdgeRemover::init() { energy_time = 0; @@ -79,7 +82,7 @@ void EdgeRemover::swap(){ std::array v_ids=ele.v_ids; er_queue.pop(); -// cout< tmp_v_ids = er_queue.top().v_ids; @@ -102,21 +105,19 @@ void EdgeRemover::swap(){ } // if(is_fail){ -// cout<<"f"<& old_t_ids) { @@ -199,7 +200,7 @@ bool EdgeRemover::removeAnEdge_32(int v1_id, int v2_id, const std::vector& auto it = std::find(fs.begin(), fs.end(), tmp); is_surface_fs[t_ids[0]][i] = is_sf_fs[it - fs.begin()]; } else - is_surface_fs[t_ids[0]][i] = NOT_SURFACE; + is_surface_fs[t_ids[0]][i] = State::state().NOT_SURFACE; if (tets[t_ids[1]][i] != v1_id) { std::array tmp = {tets[t_ids[1]][(i + 1) % 4], tets[t_ids[1]][(i + 2) % 4], @@ -208,7 +209,7 @@ bool EdgeRemover::removeAnEdge_32(int v1_id, int v2_id, const std::vector& auto it = std::find(fs.begin(), fs.end(), tmp); is_surface_fs[t_ids[1]][i] = is_sf_fs[it - fs.begin()]; } else - is_surface_fs[t_ids[1]][i] = NOT_SURFACE; + is_surface_fs[t_ids[1]][i] = State::state().NOT_SURFACE; } tet_vertices[v_ids[0]].conn_tets.erase(std::find(tet_vertices[v_ids[0]].conn_tets.begin(), @@ -401,7 +402,7 @@ bool EdgeRemover::removeAnEdge_44(int v1_id, int v2_id, const std::vector& for (int i = 0; i < old_t_ids.size(); i++) {//old_t_ids contains new tets for (int j = 0; j < 4; j++) { - is_surface_fs[old_t_ids[i]][j] = NOT_SURFACE; + is_surface_fs[old_t_ids[i]][j] = State::state().NOT_SURFACE; if (tets[old_t_ids[i]][j] == v_ids[0] || tets[old_t_ids[i]][j] == v_ids[1]) { std::array tmp = {tets[old_t_ids[i]][(j + 1) % 4], tets[old_t_ids[i]][(j + 2) % 4], tets[old_t_ids[i]][(j + 3) % 4]}; @@ -572,7 +573,7 @@ bool EdgeRemover::removeAnEdge_56(int v1_id, int v2_id, const std::vector& qs.push_back(tet_qs[(i - 1 + 5) % 5][j]); } if(qs.size() != 6){ - cout<<"ERROR: qs.size() != 6"<& //update on_surface -- 2 for (int i = 0; i < new_t_ids.size(); i++) { for (int j = 0; j < 4; j++) { - is_surface_fs[new_t_ids[i]][j] = NOT_SURFACE; + is_surface_fs[new_t_ids[i]][j] = State::state().NOT_SURFACE; if (tets[new_t_ids[i]][j] != v1_id && tets[new_t_ids[i]][j] != v2_id && tets[new_t_ids[i]][j] != n12_v_ids[(selected_id + 1) % 5] && tets[new_t_ids[i]][j] != n12_v_ids[(selected_id - 1 + 5) % 5]) { @@ -767,4 +768,6 @@ void EdgeRemover::addNewEdge(const std::array& e){ } } } -} \ No newline at end of file +} + +} // namespace tetwild diff --git a/src/tetwild/EdgeRemover.h b/src/tetwild/EdgeRemover.h index d11cd21..d05bed5 100644 --- a/src/tetwild/EdgeRemover.h +++ b/src/tetwild/EdgeRemover.h @@ -1,9 +1,9 @@ // This file is part of TetWild, a software for generating tetrahedral meshes. -// +// // Copyright (C) 2018 Yixin Hu -// -// This Source Code Form is subject to the terms of the Mozilla Public License -// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can // obtain one at http://mozilla.org/MPL/2.0/. // // Created by Yixin Hu on 4/17/17. @@ -14,6 +14,8 @@ #include +namespace tetwild { + class ElementInQueue_er{ public: std::array v_ids; @@ -70,5 +72,6 @@ class EdgeRemover:public LocalOperations { double energy_time = 0; }; +} // namespace tetwild #endif //NEW_GTET_EDGEREMOVER_H diff --git a/src/tetwild/EdgeSplitter.cpp b/src/tetwild/EdgeSplitter.cpp index 357f97e..d423a2e 100644 --- a/src/tetwild/EdgeSplitter.cpp +++ b/src/tetwild/EdgeSplitter.cpp @@ -1,9 +1,9 @@ // This file is part of TetWild, a software for generating tetrahedral meshes. -// +// // Copyright (C) 2018 Yixin Hu -// -// This Source Code Form is subject to the terms of the Mozilla Public License -// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can // obtain one at http://mozilla.org/MPL/2.0/. // // Created by Yixin Hu on 4/11/17. @@ -11,6 +11,8 @@ #include +namespace tetwild { + void EdgeSplitter::getMesh_ui(const std::vector>& tets, Eigen::MatrixXd& V, Eigen::MatrixXi& F) { ///get V, F, C V.resize(tets.size() * 4, 3); @@ -92,23 +94,21 @@ void EdgeSplitter::split() { if (t_slot_size < es_queue.size() * 6 * 2) tets.reserve(es_queue.size() * 6 * 2 - t_slot_size + 1); } -#ifndef MUTE_COUT - cout << es_queue.size() << endl; - cout << "ideal_weight = " << ideal_weight << endl; -#endif + logger().debug("{}", es_queue.size()); + logger().debug("ideal_weight = {}", ideal_weight); while (!es_queue.empty()) { const ElementInQueue_es &ele = es_queue.top(); std::array v_ids = ele.v_ids; -// if (is_print_tmp) -// cout << v_ids[0] << ' ' << v_ids[1] -// << " " << std::sqrt(calEdgeLength(v_ids)) -// << " " << std::sqrt(ideal_weight) * +// if (State::state().is_print_tmp) +// logger().debug("{}{}{} {} {} {} {}", v_ids[0], ' ', v_ids[1] +//, std::sqrt(calEdgeLength(v_ids)) +//, std::sqrt(ideal_weight) * // (tet_vertices[v_ids[0]].adaptive_scale + tet_vertices[v_ids[1]].adaptive_scale) / 2.0 -// << " " << tet_vertices[v_ids[0]].adaptive_scale -// << " " << tet_vertices[v_ids[1]].adaptive_scale -// << endl; +//, tet_vertices[v_ids[0]].adaptive_scale +//, tet_vertices[v_ids[1]].adaptive_scale +//); es_queue.pop(); if (splitAnEdge(v_ids)) suc_counter++; @@ -146,12 +146,6 @@ void EdgeSplitter::split() { } -void printConnTets(const TetVertex& v){ - for(auto it=v.conn_tets.begin();it!=v.conn_tets.end();it++) - cout<<*it<<' '; - cout<& edge) { int v1_id = edge[0]; int v2_id = edge[1]; @@ -245,10 +239,10 @@ bool EdgeSplitter::splitAnEdge(const std::array& edge) { } //update surface tags - if (g_eps != EPSILON_INFINITE) { + if (State::state().g_eps != State::state().EPSILON_INFINITE) { if (isEdgeOnSurface(v1_id, v2_id)) { tet_vertices[v_id].is_on_surface = true; - if (g_eps == EPSILON_NA) { + if (State::state().g_eps == State::state().EPSILON_NA) { setIntersection(tet_vertices[v1_id].on_edge, tet_vertices[v2_id].on_edge, tet_vertices[v_id].on_edge); setIntersection(tet_vertices[v1_id].on_face, tet_vertices[v2_id].on_face, tet_vertices[v_id].on_face); } @@ -273,7 +267,7 @@ bool EdgeSplitter::splitAnEdge(const std::array& edge) { for (int i = 0; i < new_t_ids.size(); i++) { for (int j = 0; j < 4; j++) {//v1->v if (tets[new_t_ids[i]][j] == v2_id) - is_surface_fs[new_t_ids[i]][j] = NOT_SURFACE; + is_surface_fs[new_t_ids[i]][j] = State::state().NOT_SURFACE; // else if(tets[new_t_ids[i]][j]==v_id)//no need to change // is_surface_fs[new_t_ids[i]][j]=is_surface_fs[old_t_ids[i]][j]; } @@ -281,7 +275,7 @@ bool EdgeSplitter::splitAnEdge(const std::array& edge) { for (int i = 0; i < old_t_ids.size(); i++) { for (int j = 0; j < 4; j++) {//v2->v if (tets[old_t_ids[i]][j] == v1_id) - is_surface_fs[old_t_ids[i]][j] = NOT_SURFACE; + is_surface_fs[old_t_ids[i]][j] = State::state().NOT_SURFACE; } } @@ -349,7 +343,7 @@ int EdgeSplitter::getOverRefineScale(int v1_id, int v2_id){ std::vector n12_t_ids; setIntersection(tet_vertices[v1_id].conn_tets, tet_vertices[v2_id].conn_tets, n12_t_ids); for(int i=0;i 500) {//todo: add || for other types of energy int scale = 1; scale = (tet_qualities[n12_t_ids[i]].slim_energy - 500) / 500.0; @@ -376,7 +370,7 @@ bool EdgeSplitter::isSplittable_cd1(double weight) { bool EdgeSplitter::isSplittable_cd1(int v1_id, int v2_id, double weight) { double adaptive_scale = (tet_vertices[v1_id].adaptive_scale + tet_vertices[v2_id].adaptive_scale) / 2.0; // if(adaptive_scale==0){ -// cout<<"adaptive_scale==0!!!"< ideal_weight * adaptive_scale * adaptive_scale) return true; @@ -414,3 +408,4 @@ void EdgeSplitter::getNewTetSlots(int n, std::vector& new_conn_tets) { } } +} // namespace tetwild diff --git a/src/tetwild/EdgeSplitter.h b/src/tetwild/EdgeSplitter.h index 5cca7f9..7cc4f29 100644 --- a/src/tetwild/EdgeSplitter.h +++ b/src/tetwild/EdgeSplitter.h @@ -1,9 +1,9 @@ // This file is part of TetWild, a software for generating tetrahedral meshes. -// +// // Copyright (C) 2018 Yixin Hu -// -// This Source Code Form is subject to the terms of the Mozilla Public License -// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can // obtain one at http://mozilla.org/MPL/2.0/. // // Created by Yixin Hu on 4/11/17. @@ -15,6 +15,8 @@ #include //#include +namespace tetwild { + class ElementInQueue_es{ public: std::array v_ids; @@ -62,4 +64,6 @@ class EdgeSplitter:public LocalOperations { unsigned int budget = 0; }; +} // namespace tetwild + #endif //NEW_GTET_EDGESPLITTER_H diff --git a/src/tetwild/InoutFiltering.cpp b/src/tetwild/InoutFiltering.cpp index 938f69a..5a9a5f2 100644 --- a/src/tetwild/InoutFiltering.cpp +++ b/src/tetwild/InoutFiltering.cpp @@ -1,23 +1,23 @@ // This file is part of TetWild, a software for generating tetrahedral meshes. -// +// // Copyright (C) 2018 Yixin Hu -// -// This Source Code Form is subject to the terms of the Mozilla Public License -// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can // obtain one at http://mozilla.org/MPL/2.0/. // // Created by yihu on 6/13/17. // #include - +#include #include #include +namespace tetwild { + void InoutFiltering::filter() { -#ifndef MUTE_COUT - cout << "In/out filtering..." << endl; -#endif + logger().debug("In/out filtering..."); Eigen::MatrixXd C(std::count(t_is_removed.begin(), t_is_removed.end(), false), 3); int cnt = 0; @@ -52,13 +52,13 @@ void InoutFiltering::filter() { //if the surface is totally reversed //TODO: test the correctness // if(std::count(tmp_t_is_removed.begin(), tmp_t_is_removed.end(), false)==0) { -// cout<<"Winding number gives a empty mesh! trying again"< 0) {//outside + if (is_surface_fs[i][j] != State::state().NOT_SURFACE && is_surface_fs[i][j] > 0) {//outside std::array v_ids = {tets[i][(j + 1) % 4], tets[i][(j + 2) % 4], tets[i][(j + 3) % 4]}; if (CGAL::orientation(tet_vertices[v_ids[0]].pos, tet_vertices[v_ids[1]].pos, tet_vertices[v_ids[2]].pos, tet_vertices[tets[i][j]].pos) != CGAL::POSITIVE) { @@ -117,7 +115,7 @@ void InoutFiltering::getSurface(Eigen::MatrixXd& V, Eigen::MatrixXi& F){ F(i, j)=map_ids[fs[i][j]]; } -// igl::writeSTL(g_working_dir+g_postfix+"_debug.stl", V, F); +// igl::writeSTL(State::state().g_working_dir+State::state().g_postfix+"_debug.stl", V, F); } void InoutFiltering::outputWindingNumberField(const Eigen::VectorXd& W){ @@ -136,7 +134,7 @@ void InoutFiltering::outputWindingNumberField(const Eigen::VectorXd& W){ for (int i = 0; i < v_ids.size(); i++) map_ids[v_ids[i]] = i; - PyMesh::MshSaver mSaver(g_working_dir+g_postfix+"_wn.msh", true); + PyMesh::MshSaver mSaver(State::state().g_working_dir+State::state().g_postfix+"_wn.msh", true); Eigen::VectorXd oV(v_ids.size() * 3); Eigen::VectorXi oT(t_cnt * 4); for (int i = 0; i < v_ids.size(); i++) { @@ -152,8 +150,10 @@ void InoutFiltering::outputWindingNumberField(const Eigen::VectorXd& W){ cnt++; } mSaver.save_mesh(oV, oT, 3, mSaver.TET); - cout << "#v = " << oV.rows() / 3 << endl; - cout << "#t = " << oT.rows() / 4 << endl; + logger().debug("#v = {}", oV.rows() / 3); + logger().debug("#t = {}", oT.rows() / 4); mSaver.save_elem_scalar_field("winding number", W); -} \ No newline at end of file +} + +} // namespace tetwild diff --git a/src/tetwild/InoutFiltering.h b/src/tetwild/InoutFiltering.h index dc7ed34..cca5647 100644 --- a/src/tetwild/InoutFiltering.h +++ b/src/tetwild/InoutFiltering.h @@ -1,9 +1,9 @@ // This file is part of TetWild, a software for generating tetrahedral meshes. -// +// // Copyright (C) 2018 Yixin Hu -// -// This Source Code Form is subject to the terms of the Mozilla Public License -// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can // obtain one at http://mozilla.org/MPL/2.0/. // // Created by yihu on 6/13/17. @@ -12,9 +12,11 @@ #ifndef NEW_GTET_INOUTFILTERING_H #define NEW_GTET_INOUTFILTERING_H -#include +#include #include +namespace tetwild { + class InoutFiltering { public: std::vector& tet_vertices; @@ -37,5 +39,6 @@ class InoutFiltering { void outputWindingNumberField(const Eigen::VectorXd& W); }; +} // namespace tetwild #endif //NEW_GTET_INOUTFILTERING_H diff --git a/src/tetwild/LocalOperations.cpp b/src/tetwild/LocalOperations.cpp index 22bee74..fb9e5d2 100644 --- a/src/tetwild/LocalOperations.cpp +++ b/src/tetwild/LocalOperations.cpp @@ -12,13 +12,17 @@ //#include #include +#include #include +#include //#include //#include #define FULL_LOG false #define CHECK_ENVELOP true +namespace tetwild { + double LocalOperations::comformalAMIPSEnergy_new(const std::vector& T) { double helper_0[12]; helper_0[0] = T[0]; @@ -294,15 +298,15 @@ void LocalOperations::check() { for (auto it = tet_vertices[i].conn_tets.begin(); it != tet_vertices[i].conn_tets.end(); it++) { if (t_is_removed[*it]) - cout << "t " << *it << " is removed!" << endl; + logger().debug("t {} is removed!", *it); auto jt=std::find(tets[*it].begin(), tets[*it].end(), i); if(jt==tets[*it].end()){ - cout<<"t "<<*it<<" is not a conn_tet for v "<0); } @@ -327,11 +331,11 @@ void LocalOperations::check() { is_found = true; } if (t_is_removed[*it]) - cout << "tet " << *it << " is removed!" << endl; + logger().debug("tet {} is removed!", *it); } if (!is_found) { - cout << tets[i][0] << " " << tets[i][1] << " " << tets[i][2] << " " << tets[i][3] << endl; - cout << "tet " << i << " should be conn to v " << tets[i][j] << endl; + logger().debug("{} {} {} {}", tets[i][0], tets[i][1], tets[i][2], tets[i][3]); + logger().debug("tet {} should be conn to v {}", i, tets[i][j]); } } } @@ -340,14 +344,14 @@ void LocalOperations::check() { } void LocalOperations::outputInfo(int op_type, double time, bool is_log) { - cout<<"outputing info"<> fs; // for (int i = 0; i < tets.size(); i++) { // if (t_is_removed[i]) // continue; // for (int j = 0; j < 4; j++) { -// if (is_surface_fs[i][j] != NOT_SURFACE) { +// if (is_surface_fs[i][j] != State::state().NOT_SURFACE) { // std::array f = {tets[i][(j + 1) % 4], tets[i][(j + 2) % 4], tets[i][(j + 3) % 4]}; // std::sort(f.begin(), f.end()); // fs.push_back(std::array({f[0], f[1], f[2], is_surface_fs[i][j]})); @@ -408,21 +412,19 @@ void LocalOperations::outputInfo(int op_type, double time, bool is_log) { // } // } // if (fs.size() % 2 != 0) { -// cout << "fs.size()%2!=0" << endl; +// logger().debug("fs.size()%2!=0"); // } // std::sort(fs.begin(), fs.end()); // for (int i = 0; i < fs.size() - 1; i += 2) { // if (fs[i][0] == fs[i + 1][0] && fs[i][1] == fs[i + 1][1] && fs[i][2] == fs[i + 1][2] && // fs[i][3] + fs[i + 1][3] == 0);//good // else { -// cout << i << endl; -// cout << "hehehehe" << endl; +// logger().debug("{}", i); +// logger().debug("hehehehe"); // for (int j = 0; j < 4; j++) -// cout << fs[i][j] << " "; -// cout << endl; +// logger().debug("{}{}{}{}#vertices outside of envelop = {}", fs[i][j], " "; // for (int j = 0; j < 4; j++) -// cout << fs[i + 1][j] << " "; -// cout << endl; +// std::cout, fs[i + 1][j], " "; // pausee(); // } // } @@ -434,11 +436,11 @@ void LocalOperations::outputInfo(int op_type, double time, bool is_log) { // if (!v_is_removed[i] && tet_vertices[i].is_on_surface) { // double dis = geo_sf_tree.squared_distance( // GEO::vec3(tet_vertices[i].posf[0], tet_vertices[i].posf[1], tet_vertices[i].posf[2])); -// if (dis > g_eps_2) +// if (dis > State::state().g_eps_2) // cnt++; // } // } -// cout << "#vertices outside of envelop = " << cnt << endl; +// std::cout, cnt); // } int cnt = 0; @@ -450,10 +452,10 @@ void LocalOperations::outputInfo(int op_type, double time, bool is_log) { // if (tet_vertices[i].pos[0] != tet_vertices[i].posf[0] // || tet_vertices[i].pos[1] != tet_vertices[i].posf[1] // || tet_vertices[i].pos[2] != tet_vertices[i].posf[2]) { -// cout << "tet_vertices[i].pos!=tet_vertices[i].posf" << endl; -// cout << tet_vertices[i].pos[0] - tet_vertices[i].posf[0] << " " -// << tet_vertices[i].pos[1] - tet_vertices[i].posf[1] << " " -// << tet_vertices[i].pos[2] - tet_vertices[i].posf[2] << endl; +// logger().debug("tet_vertices[i].pos!=tet_vertices[i].posf"); +// logger().debug("{}{}{}{}{}", tet_vertices[i].pos[0] - tet_vertices[i].posf[0], " " +//, tet_vertices[i].pos[1] - tet_vertices[i].posf[1], " " +//, tet_vertices[i].pos[2] - tet_vertices[i].posf[2]); // // } r_cnt++; @@ -462,26 +464,26 @@ void LocalOperations::outputInfo(int op_type, double time, bool is_log) { // if (CGAL::to_double(tet_vertices[i].pos[0]) != tet_vertices[i].posf[0] // || CGAL::to_double(tet_vertices[i].pos[1]) != tet_vertices[i].posf[1] // || CGAL::to_double(tet_vertices[i].pos[2]) != tet_vertices[i].posf[2]) { -// cout << "CGAL::to_double(tet_vertices[i].pos)!=tet_vertices[i].posf" << endl; -// cout << CGAL::to_double(tet_vertices[i].pos[0]) - tet_vertices[i].posf[0] << " " -// << CGAL::to_double(tet_vertices[i].pos[1]) - tet_vertices[i].posf[1] << " " -// << CGAL::to_double(tet_vertices[i].pos[2]) - tet_vertices[i].posf[2] << endl; +// logger().debug("CGAL::to_double(tet_vertices[i].pos)!=tet_vertices[i].posf"); +// logger().debug("{}{}{}{}{}", CGAL::to_double(tet_vertices[i].pos[0]) - tet_vertices[i].posf[0], " " +//, CGAL::to_double(tet_vertices[i].pos[1]) - tet_vertices[i].posf[1], " " +//, CGAL::to_double(tet_vertices[i].pos[2]) - tet_vertices[i].posf[2]); // } // } } } - cout << "# vertices = " << cnt << "(" << tet_vertices.size() << ") " << r_cnt << "(r)" << endl; + logger().debug("# vertices = {}({}) {}(r)", cnt, tet_vertices.size(), r_cnt); cnt = 0; for (int i = 0; i < tets.size(); i++) { if (!t_is_removed[i]) cnt++; } - cout << "# tets = " << cnt << "(" << tets.size() << ")" << endl; - cout << "# total operations = " << counter << endl; - cout << "# accepted operations = " << suc_counter << endl; + logger().debug("# tets = {}({})", cnt, tets.size()); + logger().debug("# total operations = {}", counter); + logger().debug("# accepted operations = {}", suc_counter); double min = 10, max = 0; @@ -517,18 +519,10 @@ void LocalOperations::outputInfo(int op_type, double time, bool is_log) { } } - cout << "min_d_angle = " << min - << ", max_d_angle = " << max - << ", max_slim_energy = " << max_slim_energy - << endl; - cout << "avg_min_d_angle = " << min_avg / cnt - << ", avg_max_d_angle = " << max_avg / cnt - << ", avg_slim_energy = " << avg_slim_energy / cnt - << endl; - cout << "min_d_angle: <6 " << cmp_cnt[0] / cnt << "; <12 " << cmp_cnt[1] / cnt << "; <18 " << cmp_cnt[2] / cnt - << endl; - cout << "max_d_angle: >174 " << cmp_cnt[5] / cnt << "; >168 " << cmp_cnt[4] / cnt << "; >162 " << cmp_cnt[3] / cnt - << endl; + logger().debug("min_d_angle = {}, max_d_angle = {}, max_slim_energy = {}", min, max, max_slim_energy); + logger().debug("avg_min_d_angle = {}, avg_max_d_angle = {}, avg_slim_energy = {}", min_avg / cnt, max_avg / cnt, avg_slim_energy / cnt); + logger().debug("min_d_angle: <6 {}; <12 {}; <18 {}", cmp_cnt[0] / cnt, cmp_cnt[1] / cnt, cmp_cnt[2] / cnt); + logger().debug("max_d_angle: >174 {}; >168 {}; >162 {}", cmp_cnt[5] / cnt, cmp_cnt[4] / cnt, cmp_cnt[3] / cnt); if(is_log) addRecord(MeshRecord(op_type, time, std::count(v_is_removed.begin(), v_is_removed.end(), false), cnt, @@ -578,13 +572,13 @@ bool LocalOperations::isFlip(const std::vector>& new_tets) { void LocalOperations::getCheckQuality(const std::vector& tet_qs, TetQuality& tq) { double slim_sum = 0, slim_max = 0; for (int i = 0; i < tet_qs.size(); i++) { - if (is_using_energy_max) { + if (State::state().is_using_energy_max) { if (tet_qs[i].slim_energy > slim_max) slim_max = tet_qs[i].slim_energy; } else slim_sum += tet_qs[i].slim_energy * tet_qs[i].volume; } - if (is_using_energy_max) + if (State::state().is_using_energy_max) tq.slim_energy = slim_max; else tq.slim_energy = slim_sum; @@ -593,13 +587,13 @@ void LocalOperations::getCheckQuality(const std::vector& tet_qs, Tet void LocalOperations::getCheckQuality(const std::vector& t_ids, TetQuality& tq){ double slim_sum = 0, slim_max = 0; for (int i = 0; i < t_ids.size(); i++) { - if (is_using_energy_max) { + if (State::state().is_using_energy_max) { if (tet_qualities[t_ids[i]].slim_energy > slim_max) slim_max = tet_qualities[t_ids[i]].slim_energy; } else slim_sum += tet_qualities[t_ids[i]].slim_energy * tet_qualities[t_ids[i]].volume; } - if (is_using_energy_max) + if (State::state().is_using_energy_max) tq.slim_energy = slim_max; else tq.slim_energy = slim_sum; @@ -621,7 +615,7 @@ void LocalOperations::getAvgMaxEnergy(double& avg_tq, double& max_tq) { } avg_tq /= cnt; if(std::isinf(avg_tq)) - avg_tq = MAX_ENERGY; + avg_tq = State::state().MAX_ENERGY; } double LocalOperations::getMaxEnergy(){ @@ -642,7 +636,7 @@ double LocalOperations::getSecondMaxEnergy(double max_energy){ for (unsigned int i = 0; i < tet_qualities.size(); i++) { if (t_is_removed[i]) continue; - if(tet_qualities[i].slim_energy == MAX_ENERGY) + if(tet_qualities[i].slim_energy == State::state().MAX_ENERGY) continue; if(isTetLocked_ui(i)) continue; @@ -659,12 +653,12 @@ double LocalOperations::getFilterEnergy(bool& is_clean_up) { for (unsigned int i = 0; i < tet_qualities.size(); i++) { if (t_is_removed[i]) continue; - if (tet_qualities[i].slim_energy > args.filter_energy - 1 + 1e10) + if (tet_qualities[i].slim_energy > GArgs::args().filter_energy - 1 + 1e10) buckets[10]++; else { for (int j = 0; j < 10; j++) { - if (tet_qualities[i].slim_energy > args.filter_energy - 1 + pow(10, j) - && tet_qualities[i].slim_energy <= args.filter_energy - 1 + pow(10, j + 1)) { + if (tet_qualities[i].slim_energy > GArgs::args().filter_energy - 1 + pow(10, j) + && tet_qualities[i].slim_energy <= GArgs::args().filter_energy - 1 + pow(10, j + 1)) { buckets[j]++; break; } @@ -688,7 +682,7 @@ double LocalOperations::getFilterEnergy(bool& is_clean_up) { for (int i = 0; i < 8; i++) { if (tmps1[i] < tmps2[i] && tmps1[i + 1] > tmps2[i + 1]){ - return args.filter_energy - 1 + 5 * pow(10, i+1); + return GArgs::args().filter_energy - 1 + 5 * pow(10, i+1); } } @@ -718,7 +712,7 @@ void LocalOperations::calTetQualities(const std::vector>& new static double* energy = 0; if (T0 == 0) { - std::cerr << "Initial ISPC allocation: n = " << n << endl; + logger().warn("Initial ISPC allocation: n = {}", n); current_max_size = n; T0 = new double[n]; T1 = new double[n]; @@ -736,7 +730,7 @@ void LocalOperations::calTetQualities(const std::vector>& new } if (current_max_size < n) { - std::cerr << "ISPC reallocation: n = " << n << endl; + logger.warn("ISPC reallocation: n = {}", n); free(T0); free(T1); free(T2); @@ -790,13 +784,13 @@ void LocalOperations::calTetQualities(const std::vector>& new tet_vertices[new_tets[i][2]].posf, tet_vertices[new_tets[i][3]].posf); if (ori != CGAL::POSITIVE) { - tet_qs[i].slim_energy = MAX_ENERGY; + tet_qs[i].slim_energy = State::state().MAX_ENERGY; continue; } else tet_qs[i].slim_energy = energy[i]; if (std::isinf(energy[i]) || std::isnan(energy[i])) - tet_qs[i].slim_energy = MAX_ENERGY; + tet_qs[i].slim_energy = State::state().MAX_ENERGY; } #else for (int i = 0; i < new_tets.size(); i++) { @@ -836,15 +830,15 @@ void LocalOperations::calTetQuality_AD(const std::array& tet, TetQuality heights[i] = CGAL::squared_distance(tet_vertices[tet[i]].posf, tmp_p); // if(std::isnan(heights[i])){//because pln is degenerate -// cout<& tet, TetQuality if(*tmp == 0 || heights[i] == 0){ t_quality.min_d_angle = 0; t_quality.max_d_angle = M_PI; -// t_quality.asp_ratio_2 = MAX_ENERGY; +// t_quality.asp_ratio_2 = State::state().MAX_ENERGY; return; } else if (*tmp < 1e-5) { nv[i] = Vector_3f(nv[i][0] / *tmp, nv[i][1] / *tmp, nv[i][2] / *tmp); @@ -893,13 +887,13 @@ void LocalOperations::calTetQuality_AD(const std::array& tet, TetQuality } void LocalOperations::calTetQuality_AMIPS(const std::array& tet, TetQuality& t_quality) { - if (energy_type == ENERGY_AMIPS) { + if (energy_type == State::state().ENERGY_AMIPS) { CGAL::Orientation ori = CGAL::orientation(tet_vertices[tet[0]].posf, tet_vertices[tet[1]].posf, tet_vertices[tet[2]].posf, tet_vertices[tet[3]].posf); if (ori != CGAL::POSITIVE) {//degenerate in floats - t_quality.slim_energy = MAX_ENERGY; + t_quality.slim_energy = State::state().MAX_ENERGY; } else { std::vector T; T.reserve(12); @@ -909,11 +903,11 @@ void LocalOperations::calTetQuality_AMIPS(const std::array& tet, TetQual } t_quality.slim_energy = comformalAMIPSEnergy_new(T); if (std::isinf(t_quality.slim_energy) || std::isnan(t_quality.slim_energy)) - t_quality.slim_energy = MAX_ENERGY; + t_quality.slim_energy = State::state().MAX_ENERGY; } } if(std::isinf(t_quality.slim_energy) || std::isnan(t_quality.slim_energy) || t_quality.slim_energy <= 0) - t_quality.slim_energy = MAX_ENERGY; + t_quality.slim_energy = State::state().MAX_ENERGY; } bool LocalOperations::isEdgeOnSurface(int v1_id, int v2_id) { @@ -939,7 +933,7 @@ bool LocalOperations::isEdgeOnSurface(int v1_id, int v2_id, const std::vector 2) return false; @@ -998,7 +992,7 @@ bool LocalOperations::isEdgeOnBoundary(int v1_id, int v2_id) { bool LocalOperations::isFaceOutEnvelop(const Triangle_3f& tri) { #if CHECK_ENVELOP - if(is_using_sampling){ + if(State::state().is_using_sampling){ return isFaceOutEnvelop_sampling(tri); } return true; @@ -1018,9 +1012,9 @@ bool LocalOperations::isFaceOutEnvelop(const Triangle_3f& tri) { // tris_queue.push(tri); // cnt++; // while (!tris_queue.empty()) { -//// cout << "depth = " << depth << ", " -//// << "cnt = " << cnt << ", " -//// << "sub_cnt = " << sub_cnt << endl; +//// logger().debug("depth = {}{}cnt = {}{}sub_cnt = {}", depth, ", " +////, cnt, ", " +////, sub_cnt); // if (depth == 6) // return true; //// return false; @@ -1037,7 +1031,7 @@ bool LocalOperations::isFaceOutEnvelop(const Triangle_3f& tri) { // tris[3] = Triangle_3f(mps[0], mps[1], mps[2]); // // for (int j = 0; j < 4; j++) { -//// cout << j << ": "; +//// logger().debug("{}{}depth = {}{}cnt = {}{}sub_cnt = {}", j, ": "; // side = getUpperLowerBounds(tris[j]); // if (side == EnvelopSide::OUTSIDE) // return true; @@ -1049,9 +1043,9 @@ bool LocalOperations::isFaceOutEnvelop(const Triangle_3f& tri) { // // tris_queue.pop(); // cnt--; -//// cout << "depth = " << depth << ", " -//// << "cnt = " << cnt << ", " -//// << "sub_cnt = " << sub_cnt << endl; +//// std::cout, depth, ", " +////, cnt, ", " +////, sub_cnt); // if (cnt == 0) { // cnt = sub_cnt; // sub_cnt = 0; @@ -1066,7 +1060,7 @@ bool LocalOperations::isFaceOutEnvelop(const Triangle_3f& tri) { bool LocalOperations::isPointOutEnvelop(const Point_3f& p) { #if CHECK_ENVELOP GEO::vec3 geo_p(p[0], p[1], p[2]); - if (geo_sf_tree.squared_distance(geo_p) > g_eps_2) + if (geo_sf_tree.squared_distance(geo_p) > State::state().g_eps_2) return true; return false; @@ -1108,7 +1102,7 @@ bool LocalOperations::isFaceOutEnvelop_sampling(const Triangle_3f& tri) { sq_dist = current_point.distance2(nearest_point); geo_sf_tree.nearest_facet_with_hint(current_point, prev_facet, nearest_point, sq_dist); double dis = current_point.distance2(nearest_point); - if (dis > g_eps_2) { + if (dis > State::state().g_eps_2) { #if TIMING_BREAKDOWN breakdown_timing0[id_aabb] += igl_timer0.getElapsedTime(); #endif @@ -1132,7 +1126,7 @@ bool LocalOperations::isFaceOutEnvelop_sampling(const Triangle_3f& tri) { bool LocalOperations::isPointOutBoundaryEnvelop(const Point_3f& p) { #if CHECK_ENVELOP GEO::vec3 geo_p(p[0], p[1], p[2]); - if (geo_b_tree.squared_distance(geo_p) > g_eps_2) { + if (geo_b_tree.squared_distance(geo_p) > State::state().g_eps_2) { return true; } return false; @@ -1145,7 +1139,7 @@ bool LocalOperations::isBoundarySlide(int v1_id, int v2_id, Point_3f& old_pf){ return false; #if CHECK_ENVELOP - if(g_is_close) + if(State::state().g_is_close) return false; std::unordered_set n_v_ids; @@ -1169,7 +1163,7 @@ bool LocalOperations::isBoundarySlide(int v1_id, int v2_id, Point_3f& old_pf){ GEO::vec3 p2(tet_vertices[v_id].posf[0], tet_vertices[v_id].posf[1], tet_vertices[v_id].posf[2]); b_points.push_back(p1); b_points.push_back(p2); - int n = GEO::distance(p1, p2) / g_dd + 1; + int n = GEO::distance(p1, p2) / State::state().g_dd + 1; if (n == 1) continue; b_points.reserve(b_points.size() + n + 1); @@ -1197,7 +1191,7 @@ bool LocalOperations::isBoundarySlide(int v1_id, int v2_id, Point_3f& old_pf){ if (!is_12_on_boundary) { GEO::vec3 p1(tet_vertices[v1_id].posf[0], tet_vertices[v1_id].posf[1], tet_vertices[v1_id].posf[2]); GEO::vec3 p2(old_pf[0], old_pf[1], old_pf[2]); - int n = GEO::distance(p1, p2) / g_dd + 1; + int n = GEO::distance(p1, p2) / State::state().g_dd + 1; b_points.reserve(b_points.size() + n + 1); b_points.push_back(p1); for (int k = 1; k <= n - 1; k++) @@ -1238,7 +1232,7 @@ bool LocalOperations::isBoundarySlide(int v1_id, int v2_id, Point_3f& old_pf){ sq_dist = current_point.distance2(nearest_point); geo_b_tree.nearest_facet_with_hint(current_point, prev_facet, nearest_point, sq_dist); double dis = current_point.distance2(nearest_point); - if (dis > g_eps_2) { + if (dis > State::state().g_eps_2) { #if TIMING_BREAKDOWN breakdown_timing0[id_aabb] += igl_timer0.getElapsedTime(); #endif @@ -1260,7 +1254,7 @@ bool LocalOperations::isBoundarySlide(int v1_id, int v2_id, Point_3f& old_pf){ bool LocalOperations::isTetOnSurface(int t_id){ for(int i=0;i<4;i++){ - if(is_surface_fs[t_id][i]!=NOT_SURFACE) + if(is_surface_fs[t_id][i]!=State::state().NOT_SURFACE) return false; } return true; @@ -1297,7 +1291,7 @@ void LocalOperations::getFaceConnTets(int v1_id, int v2_id, int v3_id, std::vect bool LocalOperations::isIsolated(int v_id) { for (auto it = tet_vertices[v_id].conn_tets.begin(); it != tet_vertices[v_id].conn_tets.end(); it++) { for (int j = 0; j < 4; j++) { - if (tets[*it][j] != v_id && is_surface_fs[*it][j] != NOT_SURFACE) + if (tets[*it][j] != v_id && is_surface_fs[*it][j] != State::state().NOT_SURFACE) return false; } } @@ -1332,10 +1326,8 @@ void LocalOperations::checkUnrounded() { if(!is_output) return; - std::streambuf *coutBuf = std::cout.rdbuf(); std::ofstream of; - of.open(g_working_dir + "unrounded_check.txt"); - std::cout.rdbuf(of.rdbuf()); + of.open(State::state().g_working_dir + "unrounded_check.txt"); int cnt_sf = 0; int cnt_b = 0; int cnt_all = 0; @@ -1360,7 +1352,7 @@ void LocalOperations::checkUnrounded() { for (int t_id:tet_vertices[i].conn_tets) { for (int j = 0; j < 4; j++) { if (tets[t_id][j] == i) { - if (is_surface_fs[t_id][j] != NOT_SURFACE) { + if (is_surface_fs[t_id][j] != State::state().NOT_SURFACE) { cnt_sf1++; is_found = true; } @@ -1378,19 +1370,15 @@ void LocalOperations::checkUnrounded() { diss.push_back(dis); } - cout << "# all = " << cnt_all << endl; - cout << "# surface = " << cnt_sf << endl; - cout << "# boundary = " << cnt_b << endl; -// cout<<"Is closed? "< -#include #include +#include +#include #ifdef TETWILD_WITH_ISPC #include #endif +namespace tetwild { + enum class EnvelopSide{ OUTSIDE=0, INSIDE=1, @@ -112,4 +115,6 @@ class LocalOperations { bool isTetLocked_ui(int tid); }; +} // namespace tetwild + #endif //NEW_GTET_LOCALOPERATIONS_H diff --git a/src/tetwild/Logger.cpp b/src/tetwild/Logger.cpp new file mode 100644 index 0000000..0708346 --- /dev/null +++ b/src/tetwild/Logger.cpp @@ -0,0 +1,49 @@ +// This file is part of TetWild, a software for generating tetrahedral meshes. +// +// Copyright (C) 2018 Jeremie Dumas +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// obtain one at http://mozilla.org/MPL/2.0/. +// +// Created by Jeremie Dumas on 09/04/18. +// + +#include +#include +#include +#include +#include +#include +#include +#include + +namespace tetwild { + +std::shared_ptr Logger::logger_; + +// Some code was copied over from +void Logger::init(bool use_cout, const std::string &filename, bool truncate) { + std::vector sinks; + if (use_cout) { + sinks.emplace_back(std::make_shared()); + } + if (!filename.empty()) { + sinks.emplace_back(std::make_shared(filename, truncate)); + } + + auto ®istry_inst = spdlog::details::registry::instance(); + + // create global thread pool if not already exists.. + std::lock_guard tp_lock(registry_inst.tp_mutex()); + auto tp = registry_inst.get_tp(); + if (tp == nullptr) { + tp = std::make_shared(spdlog::details::default_async_q_size, 1); + registry_inst.set_tp(tp); + } + + logger_ = std::make_shared("tetwild", sinks.begin(), sinks.end(), std::move(tp), spdlog::async_overflow_policy::block); + registry_inst.register_and_init(logger_); +} + +} // namespace tetwild diff --git a/src/tetwild/Logger.h b/src/tetwild/Logger.h new file mode 100644 index 0000000..d02da43 --- /dev/null +++ b/src/tetwild/Logger.h @@ -0,0 +1,34 @@ +// This file is part of TetWild, a software for generating tetrahedral meshes. +// +// Copyright (C) 2018 Jeremie Dumas +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// obtain one at http://mozilla.org/MPL/2.0/. +// +// Created by Jeremie Dumas on 09/04/18. +// + +#pragma once + +#include +#include + +namespace tetwild { + +struct Logger { + static std::shared_ptr logger_; + + // By default, write to stdout, but don't write to any file + static void init(bool use_cout = true, const std::string &filename = "", bool truncate = true); +}; + +// Retrieve current logger, or create one if not available +inline spdlog::async_logger & logger() { + if (!Logger::logger_) { + Logger::init(); + } + return *Logger::logger_; +} + +} // namespace tetwild diff --git a/src/tetwild/MeshConformer.cpp b/src/tetwild/MeshConformer.cpp index 1b7e2c4..c718f3e 100644 --- a/src/tetwild/MeshConformer.cpp +++ b/src/tetwild/MeshConformer.cpp @@ -1,9 +1,9 @@ // This file is part of TetWild, a software for generating tetrahedral meshes. -// +// // Copyright (C) 2018 Yixin Hu -// -// This Source Code Form is subject to the terms of the Mozilla Public License -// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can // obtain one at http://mozilla.org/MPL/2.0/. // // Created by Yixin Hu on 3/29/17. @@ -11,6 +11,8 @@ #include +namespace tetwild { + void MeshConformer::match() { is_matched.assign(m_faces.size(), false); matchDivFaces(); @@ -173,9 +175,7 @@ void MeshConformer::matchDivFaces() { seed_nids.insert(new_nids.begin(), new_nids.end());//c++11 } } -#ifndef MUTE_COUT - cout << std::count(is_matched.begin(), is_matched.end(), true) << " faces matched!" << endl; -#endif + logger().debug("{} faces matched!", std::count(is_matched.begin(), is_matched.end(), true)); } void MeshConformer::getOrientedVertices(int bsp_f_id){ @@ -208,9 +208,8 @@ void MeshConformer::getOrientedVertices(int bsp_f_id){ bsp_faces[bsp_f_id].vertices=vertices; if(vertices.size()!=bsp_faces[bsp_f_id].vertices.size()){ - cout<<"error!"<(&*result); return *p; } else { - cout << "error to3d!" << endl; + logger().debug("error to3d!"); exit(250); } } - - +} // namespace tetwild diff --git a/src/tetwild/MeshConformer.h b/src/tetwild/MeshConformer.h index c21f8e6..64f4cbe 100644 --- a/src/tetwild/MeshConformer.h +++ b/src/tetwild/MeshConformer.h @@ -1,9 +1,9 @@ // This file is part of TetWild, a software for generating tetrahedral meshes. -// +// // Copyright (C) 2018 Yixin Hu -// -// This Source Code Form is subject to the terms of the Mozilla Public License -// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can // obtain one at http://mozilla.org/MPL/2.0/. // // Created by Yixin Hu on 3/29/17. @@ -12,9 +12,11 @@ #ifndef GTET_MESHCONFORMER_H #define GTET_MESHCONFORMER_H -//#include -#include +#include #include +//#include + +namespace tetwild { class MeshConformer { public: @@ -49,4 +51,6 @@ class MeshConformer { Point_3 to3d(const Point_2 &p, const Plane_3 &pln); }; +} // namespace tetwild + #endif //GTET_MESHCONFORMER_H diff --git a/src/tetwild/MeshRefinement.cpp b/src/tetwild/MeshRefinement.cpp index e46e8d0..ba5fe32 100644 --- a/src/tetwild/MeshRefinement.cpp +++ b/src/tetwild/MeshRefinement.cpp @@ -1,28 +1,29 @@ // This file is part of TetWild, a software for generating tetrahedral meshes. -// +// // Copyright (C) 2018 Yixin Hu -// -// This Source Code Form is subject to the terms of the Mozilla Public License -// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can // obtain one at http://mozilla.org/MPL/2.0/. // // Created by Yixin Hu on 4/11/17. // #include -#include +#include +#include +#include +#include #include #include #include #include -#include #include +namespace tetwild { void MeshRefinement::prepareData(bool is_init) { -#ifndef MUTE_COUT igl_timer.start(); -#endif if (is_init) { t_is_removed = std::vector(tets.size(), false);//have to v_is_removed = std::vector(tet_vertices.size(), false); @@ -38,13 +39,11 @@ void MeshRefinement::prepareData(bool is_init) { getSimpleMesh(simple_mesh); GEO::MeshFacetsAABB simple_tree(simple_mesh); LocalOperations localOperation(tet_vertices, tets, is_surface_fs, v_is_removed, t_is_removed, tet_qualities, - ENERGY_AMIPS, simple_tree, simple_tree); + State::state().ENERGY_AMIPS, simple_tree, simple_tree); localOperation.calTetQualities(tets, tet_qualities, true);//cal all measure -#ifndef MUTE_COUT double tmp_time = igl_timer.getElapsedTime(); - cout << tmp_time << "s" << endl; + logger().debug("{}s", tmp_time); localOperation.outputInfo(MeshRecord::OpType::OP_OPT_INIT, tmp_time); -#endif } void MeshRefinement::round() { @@ -86,9 +85,7 @@ void MeshRefinement::round() { sub_cnt++; } } -#ifndef MUTE_COUT - cout << "round: " << cnt << "(" << tet_vertices.size() << ")" << endl; -#endif + logger().debug("round: {}({})", cnt, tet_vertices.size()); //for check // for (int i = 0; i < tets.size(); i++) { @@ -97,7 +94,7 @@ void MeshRefinement::round() { // CGAL::Orientation ori = CGAL::orientation(tet_vertices[tets[i][0]].pos, tet_vertices[tets[i][1]].pos, // tet_vertices[tets[i][2]].pos, tet_vertices[tets[i][3]].pos); // if (ori != CGAL::POSITIVE) { -// cout << "round hehe" << endl; +// logger().debug("round hehe"); // } // } } @@ -124,66 +121,46 @@ int MeshRefinement::doOperations(EdgeSplitter& splitter, EdgeCollapser& collapse double tmp_time; if (ops[0]) { -#ifndef MUTE_COUT igl_timer.start(); - cout << "edge splitting..." << endl; -#endif + logger().info("edge splitting..."); splitter.init(); splitter.split(); -#ifndef MUTE_COUT tmp_time = igl_timer.getElapsedTime(); splitter.outputInfo(MeshRecord::OpType::OP_SPLIT, tmp_time, is_log); - cout << "edge splitting done!" << endl; - cout << "time = " << tmp_time << "s" << endl; - cout << endl; -#endif + logger().info("edge splitting done!"); + logger().info("time = {}s", tmp_time); } if (ops[1]) { -#ifndef MUTE_COUT igl_timer.start(); - cout << "edge collasing..." << endl; -#endif + logger().info("edge collapsing..."); collapser.init(); collapser.collapse(); -#ifndef MUTE_COUT tmp_time = igl_timer.getElapsedTime(); collapser.outputInfo(MeshRecord::OpType::OP_COLLAPSE, tmp_time, is_log); - cout << "edge collasing done!" << endl; - cout << "time = " << tmp_time << "s" << endl; - cout << endl; -#endif + logger().info("edge collasing done!"); + logger().info("time = {}s", tmp_time); } if (ops[2]) { -#ifndef MUTE_COUT igl_timer.start(); - cout << "edge removing..." << endl; -#endif + logger().info("edge removing..."); edge_remover.init(); edge_remover.swap(); -#ifndef MUTE_COUT tmp_time = igl_timer.getElapsedTime(); edge_remover.outputInfo(MeshRecord::OpType::OP_SWAP, tmp_time, is_log); - cout << "edge removal done!" << endl; - cout << "time = " << tmp_time << "s" << endl; - cout << endl; -#endif + logger().info("edge removal done!"); + logger().info("time = {}s", tmp_time); } if (ops[3]) { -#ifndef MUTE_COUT igl_timer.start(); - cout << "vertex smoothing..." << endl; -#endif + logger().info("vertex smoothing..."); smoother.smooth(); -#ifndef MUTE_COUT tmp_time = igl_timer.getElapsedTime(); smoother.outputInfo(MeshRecord::OpType::OP_SMOOTH, tmp_time, is_log); - cout << "vertex smooth done!" << endl; - cout << "time = " << tmp_time << "s" << endl; - cout << endl; -#endif + logger().info("vertex smooth done!"); + logger().info("time = {}s", tmp_time); } round(); @@ -209,8 +186,8 @@ int MeshRefinement::doOperationLoops(EdgeSplitter& splitter, EdgeCollapser& coll double tmp_avg_energy, tmp_max_energy; splitter.getAvgMaxEnergy(tmp_avg_energy, tmp_max_energy); - if (std::abs(tmp_avg_energy - avg_energy) < args.delta_energy - && std::abs(tmp_max_energy - max_energy) < args.delta_energy) + if (std::abs(tmp_avg_energy - avg_energy) < GArgs::args().delta_energy + && std::abs(tmp_max_energy - max_energy) < GArgs::args().delta_energy) break; avg_energy = tmp_avg_energy; max_energy = tmp_max_energy; @@ -227,26 +204,26 @@ void MeshRefinement::refine(int energy_type, const std::array& ops, boo GEO::MeshFacetsAABB geo_b_tree(geo_b_mesh); if (is_dealing_unrounded) - min_adaptive_scale = g_eps / g_ideal_l * 0.5; //min to eps/2 + min_adaptive_scale = State::state().g_eps / State::state().g_ideal_l * 0.5; //min to eps/2 else -// min_adaptive_scale = g_eps_input / g_ideal_l; // g_eps_input / g_ideal_l * 0.5 is too small - min_adaptive_scale = (g_diag_l / 1000) / g_ideal_l; // set min_edge_length to diag / 1000 would be better +// min_adaptive_scale = State::state().g_eps_input / State::state().g_ideal_l; // State::state().g_eps_input / State::state().g_ideal_l * 0.5 is too small + min_adaptive_scale = (State::state().g_diag_l / 1000) / State::state().g_ideal_l; // set min_edge_length to diag / 1000 would be better LocalOperations localOperation(tet_vertices, tets, is_surface_fs, v_is_removed, t_is_removed, tet_qualities, energy_type, geo_sf_tree, geo_b_tree); - EdgeSplitter splitter(localOperation, g_ideal_l * (4.0 / 3.0) * g_ideal_l * (4.0 / 3.0)); - EdgeCollapser collapser(localOperation, g_ideal_l * (4.0 / 5.0) * g_ideal_l * (4.0 / 5.0)); - EdgeRemover edge_remover(localOperation, g_ideal_l * (4.0 / 3.0) * g_ideal_l * (4.0 / 3.0)); + EdgeSplitter splitter(localOperation, State::state().g_ideal_l * (4.0 / 3.0) * State::state().g_ideal_l * (4.0 / 3.0)); + EdgeCollapser collapser(localOperation, State::state().g_ideal_l * (4.0 / 5.0) * State::state().g_ideal_l * (4.0 / 5.0)); + EdgeRemover edge_remover(localOperation, State::state().g_ideal_l * (4.0 / 3.0) * State::state().g_ideal_l * (4.0 / 3.0)); VertexSmoother smoother(localOperation); collapser.is_check_quality = true; - if (args.mid_result == 1) + if (GArgs::args().mid_result == 1) outputMidResult(false, 1); -// double old_g_eps = g_eps; -// g_eps = 0.5 * old_g_eps; -// g_eps_2 = g_eps * g_eps; +// double old_State::state().g_eps = State::state().g_eps; +// State::state().g_eps = 0.5 * old_State::state().g_eps; +// State::state().g_eps_2 = State::state().g_eps * State::state().g_eps; if (is_pre) refine_pre(splitter, collapser, edge_remover, smoother); @@ -267,17 +244,15 @@ void MeshRefinement::refine(int energy_type, const std::array& ops, boo int update_cnt = 0; int is_output = true; // const double eps_s = 0.8; -// g_eps *= eps_s; -// g_eps_2 *= eps_s*eps_s; +// State::state().g_eps *= eps_s; +// State::state().g_eps_2 *= eps_s*eps_s; bool is_split = true; - for (int pass = old_pass; pass < old_pass + args.max_pass; pass++) { + for (int pass = old_pass; pass < old_pass + GArgs::args().max_pass; pass++) { if (is_dealing_unrounded && pass == old_pass) { - updateScalarField(false, false, args.filter_energy); + updateScalarField(false, false, GArgs::args().filter_energy); } -#ifndef MUTE_COUT - cout << "////////////////Pass " << pass << " ////////////////" << endl; -#endif + logger().info("//////////////// Pass {} ////////////////", pass); if (is_dealing_unrounded) collapser.is_limit_length = false; doOperations(splitter, collapser, edge_remover, smoother, @@ -293,25 +268,25 @@ void MeshRefinement::refine(int energy_type, const std::array& ops, boo is_finished = false; } if (is_finished) { - cout << "all vertices rounded!!" << endl; + logger().debug("all vertices rounded!!"); // break; } } - if (localOperation.getMaxEnergy() < args.filter_energy) + if (localOperation.getMaxEnergy() < GArgs::args().filter_energy) break; //check and mark is_bad_element double avg_energy, max_energy; localOperation.getAvgMaxEnergy(avg_energy, max_energy); - if (pass > 0 && pass < old_pass + args.max_pass - 1 - && avg_energy0 - avg_energy < args.delta_energy && max_energy0 - max_energy < args.delta_energy) { - -// if (args.targeted_num_v > 0 && getInsideVertexSize() > 1.05 * args.targeted_num_v && isRegionFullyRounded()) { -// if (g_cur_stage < args.stage) { -// g_eps += g_eps_delta; -// g_eps_2 = g_eps * g_eps; -// g_cur_stage++; + if (pass > 0 && pass < old_pass + GArgs::args().max_pass - 1 + && avg_energy0 - avg_energy < GArgs::args().delta_energy && max_energy0 - max_energy < GArgs::args().delta_energy) { + +// if (GArgs::args().targeted_num_v > 0 && getInsideVertexSize() > 1.05 * GArgs::args().targeted_num_v && isRegionFullyRounded()) { +// if (State::state().g_cur_stage < GArgs::args().stage) { +// State::state().g_eps += State::state().g_eps_delta; +// State::state().g_eps_2 = State::state().g_eps * State::state().g_eps; +// State::state().g_cur_stage++; // avg_energy0 = avg_energy; // max_energy0 = max_energy; // continue; @@ -330,12 +305,12 @@ void MeshRefinement::refine(int energy_type, const std::array& ops, boo continue; } if (update_buget == 0) { - if (g_cur_stage > 1 && g_cur_stage < args.stage) { - g_eps += g_eps_delta; - g_eps_2 = g_eps * g_eps; - g_cur_stage++; + if (State::state().g_cur_stage > 1 && State::state().g_cur_stage < GArgs::args().stage) { + State::state().g_eps += State::state().g_eps_delta; + State::state().g_eps_2 = State::state().g_eps * State::state().g_eps; + State::state().g_cur_stage++; update_buget = 2; -// cout<<"[[[[[[[[[[[[[[UPDATE EPSILON "<& ops, boo //get target energy double target_energy = localOperation.getMaxEnergy() / 100; target_energy = std::min(target_energy, target_energy0 / 10); - target_energy = std::max(target_energy, args.filter_energy * 0.8); + target_energy = std::max(target_energy, GArgs::args().filter_energy * 0.8); target_energy0 = target_energy; updateScalarField(false, false, target_energy); - if (g_cur_stage == 1 && g_cur_stage < args.stage - && target_energy < args.filter_energy) { - g_eps += g_eps_delta; - g_eps_2 = g_eps * g_eps; - g_cur_stage++; -// cout<<"[[[[[[[[[[[[[[UPDATE EPSILON "<& ops, boo // Eigen::MatrixXd V_tmp; // Eigen::MatrixXi F_tmp; // getSurface(V_tmp, F_tmp); -// localOperation.outputSurfaceColormap(V_tmp, F_tmp, old_g_eps); +// localOperation.outputSurfaceColormap(V_tmp, F_tmp, old_State::state().g_eps); // } } avg_energy0 = avg_energy; max_energy0 = max_energy; } - old_pass = old_pass + args.max_pass; + old_pass = old_pass + GArgs::args().max_pass; // if (!isRegionFullyRounded()) { // refine_unrounded(splitter, collapser, edge_remover, smoother); // } // if (max_energy0 > 1e3) { -// refine_local(splitter, collapser, edge_remover, smoother, args.filter_energy); +// refine_local(splitter, collapser, edge_remover, smoother, GArgs::args().filter_energy); // } // if (!isRegionFullyRounded() || max_energy0 > 1e3) -// serialization(g_working_dir + g_postfix + ".slz"); +// serialization(State::state().g_working_dir + State::state().g_postfix + ".slz"); - if (!args.is_quiet) { + if (!GArgs::args().is_quiet) { double max_e = localOperation.getMaxEnergy(); if (max_e > 100) { bool is_print = false; std::ofstream f; - f.open(g_working_dir + args.postfix + ".tmp"); + f.open(State::state().g_working_dir + GArgs::args().postfix + ".tmp"); for (int i = 0; i < tet_qualities.size(); i++) { if (t_is_removed[i]) continue; @@ -429,24 +404,24 @@ void MeshRefinement::refine(int energy_type, const std::array& ops, boo << "v2 " << tets[i][v2_id] << " " << tet_vertices[tets[i][v2_id]].is_on_surface << " " << tet_vertices[tets[i][v2_id]].is_on_boundary << " " << localOperation.isPointOutEnvelop(tet_vertices[tets[i][v2_id]].posf) << " " - << localOperation.isPointOutBoundaryEnvelop(tet_vertices[tets[i][v2_id]].posf) << endl; + << localOperation.isPointOutBoundaryEnvelop(tet_vertices[tets[i][v2_id]].posf) << std::endl; } } if (is_print) - f << g_eps << endl; + f << State::state().g_eps << std::endl; f.close(); } } if (is_post) { - if (args.targeted_num_v > 0) { + if (GArgs::args().targeted_num_v > 0) { double n = getInsideVertexSize(); - if (n > args.targeted_num_v) { + if (n > GArgs::args().targeted_num_v) { collapser.is_limit_length = false; collapser.is_soft = true; collapser.soft_energy = localOperation.getMaxEnergy(); collapser.budget = - (n - args.targeted_num_v) * std::count(v_is_removed.begin(), v_is_removed.end(), false) / n * + (n - GArgs::args().targeted_num_v) * std::count(v_is_removed.begin(), v_is_removed.end(), false) / n * 1.5; } } @@ -454,37 +429,35 @@ void MeshRefinement::refine(int energy_type, const std::array& ops, boo } - if (args.targeted_num_v > 0) + if (GArgs::args().targeted_num_v > 0) applyTargetedVertexNum(splitter, collapser, edge_remover, smoother); - if (args.bg_mesh != "") { + if (GArgs::args().bg_mesh != "") { applySizingField(splitter, collapser, edge_remover, smoother); } - if (args.mid_result == 2) + if (GArgs::args().mid_result == 2) outputMidResult(true, 2);//mark in/out -// if (!args.is_quiet) { +// if (!GArgs::args().is_quiet) { //// Eigen::MatrixXd V_tmp; //// Eigen::MatrixXi F_tmp; //// getTrackedSurface(V_tmp, F_tmp); -//// igl::writeOBJ(g_working_dir + g_postfix + "_tracked_sf1.obj", V_tmp, F_tmp); +//// igl::writeOBJ(State::state().g_working_dir + State::state().g_postfix + "_tracked_sf1.obj", V_tmp, F_tmp); //// getSurface(V_tmp, F_tmp); -//// igl::writeOBJ(g_working_dir + g_postfix + "_tracked_sf2.obj", V_tmp, F_tmp); -////// localOperation.outputSurfaceColormap(V_tmp, F_tmp, g_eps_input);//compared with user input eps +//// igl::writeOBJ(State::state().g_working_dir + State::state().g_postfix + "_tracked_sf2.obj", V_tmp, F_tmp); +////// localOperation.outputSurfaceColormap(V_tmp, F_tmp, State::state().g_eps_input);//compared with user input eps //// localOperation.checkUnrounded(); // } - if (args.is_laplacian) + if (GArgs::args().is_laplacian) postProcess(smoother); } void MeshRefinement::refine_pre(EdgeSplitter& splitter, EdgeCollapser& collapser, EdgeRemover& edge_remover, VertexSmoother& smoother){ -#ifndef MUTE_COUT - cout<<"////////////////// Pre-processing //////////////////"<{false, true, false, false}); collapser.is_limit_length = true; @@ -492,9 +465,7 @@ void MeshRefinement::refine_pre(EdgeSplitter& splitter, EdgeCollapser& collapser void MeshRefinement::refine_post(EdgeSplitter& splitter, EdgeCollapser& collapser, EdgeRemover& edge_remover, VertexSmoother& smoother){ -#ifndef MUTE_COUT - cout<<"////////////////// Post-processing //////////////////"< 0 && pass < args.max_pass - 1 - && avg_energy0 - avg_energy < args.delta_energy && max_energy - max_energy0 < args.delta_energy) { + if (pass > 0 && pass < GArgs::args().max_pass - 1 + && avg_energy0 - avg_energy < GArgs::args().delta_energy && max_energy - max_energy0 < GArgs::args().delta_energy) { updateScalarField(false, true, target_energy); } } @@ -544,18 +515,18 @@ bool MeshRefinement::refine_unrounded(EdgeSplitter& splitter, EdgeCollapser& col EdgeSplitter &localOperation = splitter; int scalar_update = 3; double old_min_adaptive_scale = min_adaptive_scale; - min_adaptive_scale = g_eps / g_ideal_l * 0.5; + min_adaptive_scale = State::state().g_eps / State::state().g_ideal_l * 0.5; collapser.is_limit_length = false; updateScalarField(true, false, -1, true); for (int pass = 0; pass < 5 * scalar_update; pass++) { - cout << "////////////////// Local Pass " << pass << " //////////////////" << endl; + logger().info("////////////////// Local Pass {} //////////////////", pass); doOperations(splitter, collapser, edge_remover, smoother); if (isRegionFullyRounded()) break; - if (scalar_update > 0 && pass % scalar_update == scalar_update - 1 && pass < args.max_pass * scalar_update - 1) { + if (scalar_update > 0 && pass % scalar_update == scalar_update - 1 && pass < GArgs::args().max_pass * scalar_update - 1) { updateScalarField(true, false, -1); } } @@ -582,7 +553,7 @@ void MeshRefinement::refine_revert(EdgeSplitter& splitter, EdgeCollapser& collap int n_v0 = std::count(v_is_removed.begin(), v_is_removed.end(), false); for (int pass = 0; pass < 10; pass++) { - cout << "////////////////// Local (revert) Pass " << pass << " //////////////////" << endl; + logger().info("////////////////// Local (revert) Pass {} //////////////////", pass); doOperations(splitter, collapser, edge_remover, smoother, std::array({false, true, true, true})); // doOperations(splitter, collapser, edge_remover, smoother); @@ -630,9 +601,9 @@ void MeshRefinement::markInOut(std::vector& tmp_t_is_removed){ Eigen::MatrixXi F; getSurface(V, F); Eigen::VectorXd W; - cout << "winding number..." << endl; + logger().debug("winding number..."); igl::winding_number(V, F, C, W); - cout << "winding number done" << endl; + logger().debug("winding number done"); cnt = 0; for (int i = 0; i < tets.size(); i++) { @@ -645,14 +616,14 @@ void MeshRefinement::markInOut(std::vector& tmp_t_is_removed){ void MeshRefinement::applySizingField(EdgeSplitter& splitter, EdgeCollapser& collapser, EdgeRemover& edge_remover, VertexSmoother& smoother) { - PyMesh::MshLoader mshLoader(args.bg_mesh); + PyMesh::MshLoader mshLoader(GArgs::args().bg_mesh); Eigen::VectorXd V_in = mshLoader.get_nodes(); Eigen::VectorXi T_in = mshLoader.get_elements(); Eigen::VectorXd values = mshLoader.get_node_field("values"); if (V_in.rows() == 0 || T_in.rows() == 0 || values.rows() == 0) return; - cout << "Applying sizing field..." << endl; + logger().debug("Applying sizing field..."); GEO::Mesh bg_mesh; bg_mesh.vertices.clear(); @@ -697,7 +668,7 @@ void MeshRefinement::applySizingField(EdgeSplitter& splitter, EdgeCollapser& col value += weight * values(T_in(bg_t_id * 4 + (j + 3) % 4)); } - tet_vertices[i].adaptive_scale = value / g_ideal_l; //we allow .adaptive_scale > 1 + tet_vertices[i].adaptive_scale = value / State::state().g_ideal_l; //we allow .adaptive_scale > 1 } // for debugging @@ -707,19 +678,19 @@ void MeshRefinement::applySizingField(EdgeSplitter& splitter, EdgeCollapser& col collapser.is_limit_length = true; collapser.is_soft = true; collapser.soft_energy = splitter.getMaxEnergy(); -// is_print_tmp = true;//debugging splitting +// State::state().is_print_tmp = true;//debugging splitting doOperationLoops(splitter, collapser, edge_remover, smoother, 20); } void MeshRefinement::applyTargetedVertexNum(EdgeSplitter& splitter, EdgeCollapser& collapser, EdgeRemover& edge_remover, VertexSmoother& smoother) { - if (args.targeted_num_v < 0) + if (GArgs::args().targeted_num_v < 0) return; - if (args.targeted_num_v == 0) + if (GArgs::args().targeted_num_v == 0) for (int i = 0; i < t_is_removed.size(); i++) t_is_removed[i] = true; - double N = args.targeted_num_v; //targeted #v + double N = GArgs::args().targeted_num_v; //targeted #v //marking in/out std::vector tmp_t_is_removed; @@ -744,7 +715,7 @@ void MeshRefinement::applyTargetedVertexNum(EdgeSplitter& splitter, EdgeCollapse if (std::abs(cnt - N) / N < size_threshold) return; - cout< target "< target {}", cnt, N); if (cnt > N) {//reduce vertices double max_energy = splitter.getMaxEnergy(); @@ -753,8 +724,8 @@ void MeshRefinement::applyTargetedVertexNum(EdgeSplitter& splitter, EdgeCollapse collapser.is_soft = true; collapser.soft_energy = max_energy; -// g_eps *= 1.5; -// g_eps_2 *= 1.5 * 1.5; +// State::state().g_eps *= 1.5; +// State::state().g_eps_2 *= 1.5 * 1.5; // for (int i = 0; i < tet_vertices.size(); i++) // tet_vertices[i].is_locked = false; @@ -790,25 +761,21 @@ bool MeshRefinement::isRegionFullyRounded(){ } void MeshRefinement::updateScalarField(bool is_clean_up_unrounded, bool is_clean_up_local, double filter_energy, bool is_lock) { -#ifndef MUTE_COUT igl_timer.start(); - cout << "marking adaptive scales..." << endl; -#endif + logger().debug("marking adaptive scales..."); double tmp_time = 0; - double radius0 = g_ideal_l * 1.8;//increasing the radius would increase the #v in output + double radius0 = State::state().g_ideal_l * 1.8;//increasing the radius would increase the #v in output if(is_hit_min) radius0 *= 2; if (is_clean_up_local) - radius0 = g_ideal_l; + radius0 = State::state().g_ideal_l; if (is_clean_up_unrounded) radius0 *= 2; -#ifndef MUTE_COUT - cout << "filter_energy = " << filter_energy << endl; -#endif + logger().debug("filter_energy = {}", filter_energy); std::vector adap_tmp(tet_vertices.size(), 1.5); - double dynamic_adaptive_scale = args.adaptive_scalar; + double dynamic_adaptive_scale = GArgs::args().adaptive_scalar; const int N = -int(std::log2(min_adaptive_scale) - 1); std::vector> v_ids(N, std::vector()); @@ -907,14 +874,12 @@ void MeshRefinement::updateScalarField(bool is_clean_up_unrounded, bool is_clean tet_vertices[i].adaptive_scale = new_scale; } if (is_clean_up_unrounded && is_lock) - cout << cnt << " vertices locked" << endl; + logger().debug("{} vertices locked", cnt); -#ifndef MUTE_COUT - cout << "marked!" << endl; + logger().debug("marked!"); tmp_time = igl_timer.getElapsedTime(); - cout << "time = " << tmp_time << "s" << endl << endl; + logger().debug("time = {}s", tmp_time); addRecord(MeshRecord(MeshRecord::OpType::OP_ADAP_UPDATE, tmp_time, -1, -1)); -#endif // outputMidResult(true); } @@ -968,8 +933,8 @@ void MeshRefinement::postProcess(VertexSmoother& smoother) { if (!tet_vertices[i].is_on_surface) b_v_ids.push_back(i); } - cout << "tmp_is_on_surface.size = " << tmp_is_on_surface.size() << endl; - cout << "b_v_ids.size = " << b_v_ids.size() << endl; + logger().debug("tmp_is_on_surface.size = {}", tmp_is_on_surface.size()); + logger().debug("b_v_ids.size = {}", b_v_ids.size()); for (int i = 0; i < 20; i++) { if (smoother.laplacianBoundary(b_v_ids, tmp_is_on_surface, tmp_t_is_removed) == 0) { break; @@ -988,9 +953,9 @@ void MeshRefinement::check() { continue; for (auto it = tet_vertices[i].conn_tets.begin(); it != tet_vertices[i].conn_tets.end(); it++) if (t_is_removed[*it]) - cout << "v " << i << " should have conn_tet t" << *it << endl; + logger().debug("v {} should have conn_tet t{}", i, *it); if (tet_vertices[i].conn_tets.size() == 0) { - cout << "empty conn_tets: v " << i << endl; + logger().debug("empty conn_tets: v {}", i); } } @@ -1005,9 +970,9 @@ void MeshRefinement::check() { tet_vertices[tets[i][3]].pos); if (ori == CGAL::COPLANAR) { - cout << "tet " << i << " is degenerate!" << endl; + logger().debug("tet {} is degenerate!", i); } else if (ori == CGAL::NEGATIVE) { - cout << "tet " << i << " is flipped!" << endl; + logger().debug("tet {} is flipped!", i); } for (int j = 0; j < 4; j++) { @@ -1018,11 +983,11 @@ void MeshRefinement::check() { is_found = true; } if (t_is_removed[*it]) - cout << "tet " << *it << " is removed!" << endl; + logger().debug("tet {} is removed!", *it); } if (!is_found) { - cout << tets[i][0] << " " << tets[i][1] << " " << tets[i][2] << " " << tets[i][3] << endl; - cout << "tet " << i << " should be conn to v " << tets[i][j] << endl; + logger().debug("{} {} {} {}", tets[i][0], tets[i][1], tets[i][2], tets[i][3]); + logger().debug("tet {} should be conn to v {}", i, tets[i][j]); } std::array f = {tets[i][j], tets[i][(j + 1) % 4], tets[i][(j + 2) % 4]}; @@ -1039,7 +1004,7 @@ void MeshRefinement::check() { setIntersection(tet_vertices[tet_faces[i][2]].conn_tets, tmp, tmp); if (tmp.size() != 1 && tmp.size() != 2) - cout << tmp.size() << endl; + logger().debug("{}", tmp.size()); } } @@ -1066,9 +1031,9 @@ void MeshRefinement::outputMidResult(bool is_with_bbox, double id) { Eigen::MatrixXi F; getSurface(V, F); Eigen::VectorXd W; - cout<<"winding number..."< 0) {//outside + if (is_surface_fs[i][j] != State::state().NOT_SURFACE && is_surface_fs[i][j] > 0) {//outside std::array v_ids = {tets[i][(j + 1) % 4], tets[i][(j + 2) % 4], tets[i][(j + 3) % 4]}; if (CGAL::orientation(tet_vertices[v_ids[0]].pos, tet_vertices[v_ids[1]].pos, tet_vertices[v_ids[2]].pos, tet_vertices[tets[i][j]].pos) != CGAL::POSITIVE) { @@ -1189,7 +1154,7 @@ void MeshRefinement::getTrackedSurface(Eigen::MatrixXd& V, Eigen::MatrixXi& F) { if (t_is_removed[i]) continue; for (int j = 0; j < 4; j++) { - if (is_surface_fs[i][j] != NOT_SURFACE && is_surface_fs[i][j] >= 0) {//outside + if (is_surface_fs[i][j] != State::state().NOT_SURFACE && is_surface_fs[i][j] >= 0) {//outside std::array v_ids = {tets[i][(j + 1) % 4], tets[i][(j + 2) % 4], tets[i][(j + 3) % 4]}; if (CGAL::orientation(tet_vertices[v_ids[0]].pos, tet_vertices[v_ids[1]].pos, tet_vertices[v_ids[2]].pos, tet_vertices[tets[i][j]].pos) != CGAL::POSITIVE) { @@ -1342,7 +1307,7 @@ void getBoudnaryMesh(const Eigen::MatrixXd& V_sf, const Eigen::MatrixXi& F_sf, G } bool MeshRefinement::deserialization(const std::string& sf_file, const std::string& slz_file) { - cout<<"deserializing ..."< -// -// This Source Code Form is subject to the terms of the Mozilla Public License -// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can // obtain one at http://mozilla.org/MPL/2.0/. // // Created by Yixin Hu on 4/11/17. @@ -16,10 +16,11 @@ #include #include #include - #include #include +namespace tetwild { + class MeshRefinement { public: //init @@ -89,5 +90,6 @@ class MeshRefinement { void serialization(const std::string& slz_file); }; +} // namespace tetwild #endif //NEW_GTET_MESHREFINEMENT_H diff --git a/src/tetwild/Preprocess.cpp b/src/tetwild/Preprocess.cpp index c8a40fb..4b5569f 100644 --- a/src/tetwild/Preprocess.cpp +++ b/src/tetwild/Preprocess.cpp @@ -1,28 +1,31 @@ // This file is part of TetWild, a software for generating tetrahedral meshes. -// +// // Copyright (C) 2018 Yixin Hu -// -// This Source Code Form is subject to the terms of the Mozilla Public License -// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can // obtain one at http://mozilla.org/MPL/2.0/. // // Created by Yixin Hu on 10/12/17. // #include -#include +#include +#include #include #include - #include #include #include #include #include #include +#include + +namespace tetwild { void checkBoundary(const Eigen::MatrixXd& V, const Eigen::MatrixXi& F) { - PyMesh::MshSaver mSaver(g_working_dir+args.postfix+"_boundary.msh", true); + PyMesh::MshSaver mSaver(State::state().g_working_dir+GArgs::args().postfix+"_boundary.msh", true); Eigen::VectorXd oV; Eigen::VectorXi oF; oV.resize(V.rows() * 3); @@ -60,12 +63,12 @@ void checkBoundary(const Eigen::MatrixXd& V, const Eigen::MatrixXi& F) { std::back_inserter(tmp)); if (tmp.size() == 1) { bF(tmp[0]) = 1; -// cout << "boundary tri! " << tmp[0] << endl; +// logger().debug("boundary tri! {}", tmp[0]); // Triangle_3f tri(Point_3f(V(F(i, 0), 0), V(F(i, 0), 1), V(F(i, 0), 2)), // Point_3f(V(F(i, 1), 0), V(F(i, 1), 1), V(F(i, 1), 2)), // Point_3f(V(F(i, 2), 0), V(F(i, 2), 1), V(F(i, 2), 2))); // if(tri.is_degenerate()) -// cout<<"degenerate"< " << V_in.rows() << endl; - cout << "#f = " << F_tmp.rows() << " -> " << F_in.rows() << endl; -#endif + logger().debug("#v = {} -> {}", V_tmp.rows(), V_in.rows()); + logger().debug("#f = {} -> {}", F_tmp.rows(), F_in.rows()); // checkBoundary(V_in, F_in); ////set global parameters - g_diag_l = igl::bounding_box_diagonal(V_in); - g_eps_input = g_diag_l / args.i_epsilon; - - if (args.i_dd > 0) {//for testing only - g_dd = g_diag_l / args.i_dd; - g_eps = g_diag_l / args.i_epsilon; - g_eps_2 = g_eps * g_eps; - args.stage = 1; + State::state().g_diag_l = igl::bounding_box_diagonal(V_in); + State::state().g_eps_input = State::state().g_diag_l / GArgs::args().i_epsilon; + + if (GArgs::args().i_dd > 0) {//for testing only + State::state().g_dd = State::state().g_diag_l / GArgs::args().i_dd; + State::state().g_eps = State::state().g_diag_l / GArgs::args().i_epsilon; + State::state().g_eps_2 = State::state().g_eps * State::state().g_eps; + GArgs::args().stage = 1; } else { -// g_dd = g_eps_input / args.stage; -// g_cur_stage = 1; -// g_eps = g_eps_input - g_dd / std::sqrt(3) * (args.stage + 1 - g_cur_stage); -// g_eps_delta = g_dd / std::sqrt(3); -// g_eps_2 = g_eps * g_eps; - - g_dd = g_eps_input / args.stage; - g_cur_stage = 1; - g_eps = g_eps_input - g_dd / std::sqrt(3) * (args.stage + 1 - g_cur_stage); - g_eps_delta = g_dd / std::sqrt(3); - g_eps_2 = g_eps * g_eps; +// State::state().g_dd = State::state().g_eps_input / GArgs::args().stage; +// State::state().g_cur_stage = 1; +// State::state().g_eps = State::state().g_eps_input - State::state().g_dd / std::sqrt(3) * (GArgs::args().stage + 1 - State::state().g_cur_stage); +// State::state().g_eps_delta = State::state().g_dd / std::sqrt(3); +// State::state().g_eps_2 = State::state().g_eps * State::state().g_eps; + + State::state().g_dd = State::state().g_eps_input / GArgs::args().stage; + State::state().g_cur_stage = 1; + State::state().g_eps = State::state().g_eps_input - State::state().g_dd / std::sqrt(3) * (GArgs::args().stage + 1 - State::state().g_cur_stage); + State::state().g_eps_delta = State::state().g_dd / std::sqrt(3); + State::state().g_eps_2 = State::state().g_eps * State::state().g_eps; } - g_ideal_l = g_diag_l / args.i_ideal_edge_length; + State::state().g_ideal_l = State::state().g_diag_l / GArgs::args().i_ideal_edge_length; -// cout << "eps = " << g_eps << endl; -// cout << "ideal_l = " << g_ideal_l << endl; +// logger().debug("eps = {}", State::state().g_eps); +// logger().debug("ideal_l = {}", State::state().g_ideal_l); ////get GEO meshes geo_sf_mesh.vertices.clear(); @@ -191,7 +192,7 @@ bool Preprocess::init(const Eigen::MatrixXd& V_tmp, const Eigen::MatrixXi& F_tmp geo_sf_mesh.facets.compute_borders(); getBoudnaryMesh(geo_b_mesh); - g_is_close = (geo_b_mesh.vertices.nb() == 0); + State::state().g_is_close = (geo_b_mesh.vertices.nb() == 0); return true; } @@ -208,7 +209,7 @@ void Preprocess::getBoudnaryMesh(GEO::Mesh& b_mesh) { //check isolated vertices // for(int i=0;i> b_edges; @@ -260,9 +261,9 @@ void Preprocess::process(GEO::Mesh& geo_sf_mesh, std::vector& m_vertice double eps_scalar = 0.8; double eps_scalar_2 = eps_scalar*eps_scalar; - g_eps *= eps_scalar; - g_eps_2 *= eps_scalar_2; -// g_dd *= eps_scalar*2; + State::state().g_eps *= eps_scalar; + State::state().g_eps_2 *= eps_scalar_2; +// State::state().g_dd *= eps_scalar*2; conn_fs.resize(V_in.size()); for (int i = 0; i < F_in.rows(); i++) { @@ -321,11 +322,9 @@ void Preprocess::process(GEO::Mesh& geo_sf_mesh, std::vector& m_vertice F_out(cnt, j) = new_v_ids[F_in(i, j)]; cnt++; } -// igl::writeSTL(g_working_dir+args.postfix+"_simplified.stl", V_out, F_out); -#ifndef MUTE_COUT - cout<<"#v = "<& m_vertice conn_fs[F_in(i, j)].insert(i); } swap(geo_face_tree); - if(args.mid_result == 0) - igl::writeSTL(g_working_dir+args.postfix+"_swapped.stl", V_in, F_in); + if(GArgs::args().mid_result == 0) + igl::writeSTL(State::state().g_working_dir+GArgs::args().postfix+"_swapped.stl", V_in, F_in); // checkBoundary(V_in, F_in); @@ -352,14 +351,12 @@ void Preprocess::process(GEO::Mesh& geo_sf_mesh, std::vector& m_vertice if (!tr.is_degenerate())//delete all degenerate triangles m_faces.push_back(f); } -#ifndef MUTE_COUT - cout << "#v = " << m_vertices.size() << endl; - cout << "#f = " << F_in.rows()<<"->"<{}", F_in.rows(), m_faces.size()); - g_eps /= eps_scalar; - g_eps_2 /= eps_scalar_2; -// g_dd /= eps_scalar*2; + State::state().g_eps /= eps_scalar; + State::state().g_eps_2 /= eps_scalar_2; +// State::state().g_dd /= eps_scalar*2; //output colormap // outputSurfaceColormap(geo_face_tree, geo_sf_mesh); @@ -459,9 +456,7 @@ void Preprocess::swap(GEO::MeshFacetsAABB& face_aabb_tree) { if (is_swapped) cnt++; } -#ifndef MUTE_COUT - cout << cnt << " faces are swapped!!" << endl; -#endif + logger().debug("{} faces are swapped!!", cnt); } double Preprocess::getCosAngle(int v_id, int v1_id, int v2_id) { @@ -473,7 +468,7 @@ double Preprocess::getCosAngle(int v_id, int v1_id, int v2_id) { void Preprocess::simplify(GEO::MeshFacetsAABB& face_aabb_tree) { int cnt = 0; -// cout << "queue.size() = " << sm_queue.size() << endl; +// logger().debug("queue.size() = {}", sm_queue.size()); while (!sm_queue.empty()) { std::array v_ids = sm_queue.top().v_ids; double old_weight = sm_queue.top().weight; @@ -487,25 +482,19 @@ void Preprocess::simplify(GEO::MeshFacetsAABB& face_aabb_tree) { inf_e_tss.push_back(ts); } else { cnt++; -#ifndef MUTE_COUT if (cnt % 1000 == 0) - cout << "1000 vertices removed" << endl; -#endif + logger().debug("1000 vertices removed"); } } -#ifndef MUTE_COUT - cout << cnt << endl; - cout << c << endl; -#endif + logger().debug("{}", cnt); + logger().debug("{}", c); if (cnt > 0) postProcess(face_aabb_tree); } void Preprocess::postProcess(GEO::MeshFacetsAABB& face_aabb_tree){ -#ifndef MUTE_COUT - cout << "postProcess!" << endl; -#endif + logger().debug("postProcess!"); std::vector> tmp_inf_es; const int inf_es_size = inf_es.size(); @@ -541,7 +530,7 @@ bool Preprocess::removeAnEdge(int v1_id, int v2_id, GEO::MeshFacetsAABB& face_aa std::vector n12_f_ids; setIntersection(conn_fs[v1_id], conn_fs[v2_id], n12_f_ids); if (n12_f_ids.size() != 2) {//!!! -// cout << "error: n12_f_ids.size()!=2" << endl; +// logger().debug("error: n12_f_ids.size()!=2"); return false; } @@ -715,8 +704,8 @@ bool Preprocess::isOutEnvelop(const std::unordered_set& new_f_ids, GEO::Mes sampleTriangle(vs, ps); -// cout << "ps.size = " << ps.size() << endl; -// cout << "is output samples?" << endl; +// logger().debug("ps.size = {}", ps.size()); +// logger().debug("is output samples?"); // int anw = 0; // cin >> anw; // if (anw != 0) { @@ -736,7 +725,7 @@ bool Preprocess::isOutEnvelop(const std::unordered_set& new_f_ids, GEO::Mes // F_tmp(1 + i, k) = (1 + i) * 3 + k; // } // } -// igl::writeSTL(g_working_dir + "_sample.stl", V_tmp, F_tmp); +// igl::writeSTL(State::state().g_working_dir + "_sample.stl", V_tmp, F_tmp); // } @@ -748,7 +737,7 @@ bool Preprocess::isOutEnvelop(const std::unordered_set& new_f_ids, GEO::Mes // int min_i = min_max.first - ls.begin(); // int max_i = min_max.second - ls.begin(); // -// double n = ls[max_i] / g_dd; +// double n = ls[max_i] / State::state().g_dd; // if (n <= 1) { // for (int i = 0; i < 3; i++) // ps.push_back(vs[i]); @@ -761,10 +750,10 @@ bool Preprocess::isOutEnvelop(const std::unordered_set& new_f_ids, GEO::Mes // break; // ps.push_back(j / n * vs[(min_i + 2) % 3] + (n - j) / n * vs[(min_i + 1) % 3]); // } -// if (ls[min_i] > g_dd) { +// if (ls[min_i] > State::state().g_dd) { // const int ps_size = ps.size(); // for (int i = 0; i < ps_size - 1; i += 2) { -// double m = GEO::length(ps[i] - ps[i + 1]) / g_dd; +// double m = GEO::length(ps[i] - ps[i + 1]) / State::state().g_dd; // if (m < 1) // break; // m = int(m) + 1; @@ -784,7 +773,7 @@ bool Preprocess::isOutEnvelop(const std::unordered_set& new_f_ids, GEO::Mes sq_dist = current_point.distance2(nearest_point); geo_face_tree.nearest_facet_with_hint(current_point, prev_facet, nearest_point, sq_dist); double dis = current_point.distance2(nearest_point); - if (dis > g_eps_2) + if (dis > State::state().g_eps_2) return true; } } @@ -793,7 +782,7 @@ bool Preprocess::isOutEnvelop(const std::unordered_set& new_f_ids, GEO::Mes } bool Preprocess::isPointOutEnvelop(int v_id, GEO::MeshFacetsAABB& geo_face_tree){ - if (geo_face_tree.squared_distance(GEO::vec3(V_in(v_id, 0), V_in(v_id, 1), V_in(v_id, 2))) > g_eps_2) + if (geo_face_tree.squared_distance(GEO::vec3(V_in(v_id, 0), V_in(v_id, 1), V_in(v_id, 2))) > State::state().g_eps_2) return true; return false; } @@ -818,13 +807,13 @@ int calEuclidean(const std::vector>& fs){ } bool Preprocess::isEuclideanValid(int v1_id, int v2_id){ -// cout<<"v1:"<> fs; @@ -840,9 +829,9 @@ bool Preprocess::isEuclideanValid(int v1_id, int v2_id){ } std::sort(fs.begin(), fs.end()); fs.erase(std::unique(fs.begin(), fs.end()), fs.end()); -// cout<<"fs.size() = "<> fs1; for(int i=0;i f = {fs[i][0], fs[i][1], fs[i][2]}; std::sort(f.begin(), f.end()); @@ -861,9 +850,9 @@ bool Preprocess::isEuclideanValid(int v1_id, int v2_id){ } std::sort(fs1.begin(), fs1.end()); fs1.erase(std::unique(fs1.begin(), fs1.end()), fs1.end()); -// cout<<"fs1.size() = "< g_dd) { +// if(ls[min_i] > State::state().g_dd) { // int ps_size = ps.size(); // for (int i = 0; i < ps_size; i += 2) { -// double m = int(GEO::length(ps[i] - ps[i + 1]) / g_dd + 1); +// double m = int(GEO::length(ps[i] - ps[i + 1]) / State::state().g_dd + 1); // if(m==0) // break; // for (int j = 1; j < m; j++) @@ -942,7 +931,7 @@ void Preprocess::outputSurfaceColormap(GEO::MeshFacetsAABB& geo_face_tree, GEO:: // geo_face_tree.nearest_facet_with_hint(current_point, prev_facet, nearest_point, sq_dist); // double dis = current_point.distance2(nearest_point); // if(f_id==2514) -// cout< max_dis) { max_dis = dis; pp=current_point; @@ -953,46 +942,42 @@ void Preprocess::outputSurfaceColormap(GEO::MeshFacetsAABB& geo_face_tree, GEO:: cnt = 0; if(f_id==1681) { for (const GEO::vec3 &p:ps) { - cout << cnt << ": " << p[0] << ", " << p[1] << ", " << p[2] << "; " << fs[cnt] << "; " - << geo_face_tree.squared_distance(p) << endl; + logger().debug("{}: {}, {}, {}; {}; {}", cnt, p[0], p[1], p[2], fs[cnt], geo_face_tree.squared_distance(p)); cnt++; } } - eps_dis(f_id) = sqrt(max_dis / g_eps_2); + eps_dis(f_id) = sqrt(max_dis / State::state().g_eps_2); if(eps_dis(f_id)>1) { - cout << "ERROR: simplified input goes outside of the envelop" << endl; - cout< vf={1681, 1675, 1671, 1666}; for(int j=0;j v_ids = {v1_id ,v2_id, v3_id}; for (int k = 0; k < 3; k++) { - cout << v_ids[k] << ": " << geo_sf_mesh.vertices.point(v_ids[k])[0] << " " - << geo_sf_mesh.vertices.point(v_ids[k])[1] - << " " << geo_sf_mesh.vertices.point(v_ids[k])[2] << endl; + logger().debug("{}: {} {} {}", v_ids[k], geo_sf_mesh.vertices.point(v_ids[k])[0], geo_sf_mesh.vertices.point(v_ids[k])[1], geo_sf_mesh.vertices.point(v_ids[k])[2]); } } @@ -1020,20 +1005,20 @@ void Preprocess::outputSurfaceColormap(GEO::MeshFacetsAABB& geo_face_tree, GEO:: } for(int i=0;i -// -// This Source Code Form is subject to the terms of the Mozilla Public License -// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can // obtain one at http://mozilla.org/MPL/2.0/. // // Created by Yixin Hu on 10/12/17. @@ -12,17 +12,20 @@ #ifndef NEW_GTET_PREPROCESS_H #define NEW_GTET_PREPROCESS_H -#include +#include +#include +#include #include #include #include #include - #include #include #include #include +namespace tetwild { + class ElementInQueue_sm{ public: std::array v_ids; @@ -82,5 +85,6 @@ class Preprocess { void outputSurfaceColormap(GEO::MeshFacetsAABB& geo_face_tree, GEO::Mesh& geo_sf_mesh); }; +} // namespace tetwild #endif //NEW_GTET_PREPROCESS_H diff --git a/src/tetwild/SimpleTetrahedralization.cpp b/src/tetwild/SimpleTetrahedralization.cpp index 68fb468..1cf2efd 100644 --- a/src/tetwild/SimpleTetrahedralization.cpp +++ b/src/tetwild/SimpleTetrahedralization.cpp @@ -1,34 +1,36 @@ // This file is part of TetWild, a software for generating tetrahedral meshes. -// +// // Copyright (C) 2018 Yixin Hu -// -// This Source Code Form is subject to the terms of the Mozilla Public License -// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can // obtain one at http://mozilla.org/MPL/2.0/. // // Created by Yixin Hu on 3/31/17. // #include +#include +#include +#include //triangulation #include #include #include -typedef CGAL::Constrained_Delaunay_triangulation_2 CDT; +typedef CGAL::Constrained_Delaunay_triangulation_2 CDT; typedef CDT::Point Point_cdt_2; -typedef CGAL::Polygon_2 Polygon_2; +typedef CGAL::Polygon_2 Polygon_2; //arrangement #include #include -typedef CGAL::Arr_segment_traits_2 Traits_2; +typedef CGAL::Arr_segment_traits_2 Traits_2; typedef Traits_2::Point_2 Point_arr_2; typedef Traits_2::X_monotone_curve_2 Segment_arr_2; typedef CGAL::Arrangement_2 Arrangement_2; -#include - +namespace tetwild { void SimpleTetrahedralization::tetra(std::vector& tet_vertices, std::vector>& tets) { std::vector &faces = MC.bsp_faces; @@ -36,9 +38,7 @@ void SimpleTetrahedralization::tetra(std::vector& tet_vertices, std:: ///cal arrangement & tetrahedralization triangulation(tet_vertices, tets); -#ifndef MUTE_COUT - cout<<"#v = "<& tet_vertice bsp_faces[i].edges.push_back(bsp_edges.size() - 1); } } -#ifndef MUTE_COUT - cout<<"2D arr "<& tet_vertice } } -#ifndef MUTE_COUT - cout<<"improvement "<>> cdt_faces(bsp_faces.size(), std::vector>()); @@ -281,7 +277,7 @@ void SimpleTetrahedralization::triangulation(std::vector& tet_vertice MC.to2d(bsp_vertices[bsp_edges[bsp_faces[i].edges[j]].vertices[1]])); } if(cdt.number_of_vertices() != bsp_faces[i].vertices.size()){ - cout<<"error: cdt.number_of_vertices() != bsp_faces[i].vertices.size()"< vs_cdt2bsp; for (int j = 0; j < bsp_faces[i].vertices.size(); j++) { @@ -364,12 +360,10 @@ void SimpleTetrahedralization::triangulation(std::vector& tet_vertice //todo: calculate a new position } } -#ifndef MUTE_COUT - cout<<"all_cnt = "<& m_f_tags, const std::vector& m_e_tags, @@ -472,9 +466,9 @@ void SimpleTetrahedralization::labelSurface(const std::vector& m_f_tags, co } ////is face on surface//// - NOT_SURFACE = m_faces.size()+1; + State::state().NOT_SURFACE = m_faces.size()+1; is_surface_fs=std::vector>(tets.size(), - std::array({NOT_SURFACE, NOT_SURFACE, NOT_SURFACE, NOT_SURFACE})); + std::array({State::state().NOT_SURFACE, State::state().NOT_SURFACE, State::state().NOT_SURFACE, State::state().NOT_SURFACE})); // std::vector> is_visited(tets.size(), std::array({false, false, false, false})); for(unsigned int i = 0; i < tets.size(); i++) { @@ -502,26 +496,26 @@ void SimpleTetrahedralization::labelSurface(const std::vector& m_f_tags, co if (!tet_vertices[tets[i][(j + 1) % 4]].is_on_surface || !tet_vertices[tets[i][(j + 2) % 4]].is_on_surface || !tet_vertices[tets[i][(j + 3) % 4]].is_on_surface) { - is_surface_fs[i][j] = NOT_SURFACE; + is_surface_fs[i][j] = State::state().NOT_SURFACE; // if (opp_i >= 0) -// is_visited[opp_i][opp_j] = NOT_SURFACE; +// is_visited[opp_i][opp_j] = State::state().NOT_SURFACE; continue; } std::unordered_set sf_faces_tmp; setIntersection(tet_vertices[tets[i][(j + 1) % 4]].on_face, tet_vertices[tets[i][(j + 2) % 4]].on_face, sf_faces_tmp); if (sf_faces_tmp.size() == 0) { - is_surface_fs[i][j] = NOT_SURFACE; + is_surface_fs[i][j] = State::state().NOT_SURFACE; // if (opp_i >= 0) -// is_visited[opp_i][opp_j] = NOT_SURFACE; +// is_visited[opp_i][opp_j] = State::state().NOT_SURFACE; continue; } std::vector sf_faces; setIntersection(sf_faces_tmp, tet_vertices[tets[i][(j + 3) % 4]].on_face, sf_faces); if (sf_faces.size() == 0) { - is_surface_fs[i][j] = NOT_SURFACE; + is_surface_fs[i][j] = State::state().NOT_SURFACE; // if (opp_i >= 0) -// is_visited[opp_i][opp_j] = NOT_SURFACE; +// is_visited[opp_i][opp_j] = State::state().NOT_SURFACE; continue; } @@ -539,7 +533,7 @@ void SimpleTetrahedralization::labelSurface(const std::vector& m_f_tags, co CGAL::Oriented_side side = pln.oriented_side(tet_vertices[tets[i][j]].pos); if (side == CGAL::ON_ORIENTED_BOUNDARY) { - cout << "ERROR: side == CGAL::ON_ORIENTED_BOUNDARY!!" << endl; + logger().debug("ERROR: side == CGAL::ON_ORIENTED_BOUNDARY!!"); exit(250); } if (side == CGAL::ON_POSITIVE_SIDE)//outside @@ -567,7 +561,7 @@ void SimpleTetrahedralization::labelSurface(const std::vector& m_f_tags, co else is_surface_fs[i][j] -= delta; // else { -// cout << "wrong direction!!" << endl; +// logger().debug("wrong direction!!"); // pausee(); // } } @@ -579,7 +573,7 @@ void SimpleTetrahedralization::labelSurface(const std::vector& m_f_tags, co // CGAL::Oriented_side side = pln.oriented_side(tet_vertices[tets[i][j]].pos); // // if (side == CGAL::ON_ORIENTED_BOUNDARY) { -// cout << "ERROR: side == CGAL::ON_ORIENTED_BOUNDARY!!" << endl; +// logger().debug("ERROR: side == CGAL::ON_ORIENTED_BOUNDARY!!"); // exit(250); // } // if (side == CGAL::ON_POSITIVE_SIDE)//outside @@ -657,9 +651,7 @@ void SimpleTetrahedralization::labelBbox(std::vector& tet_vertices, s i7 = I; i++; } -#ifndef MUTE_COUT - cout<<"#v on bbox = "<& tet_vertices, std::vector>& tets, @@ -695,9 +687,9 @@ void SimpleTetrahedralization::labelBoundary(std::vector& tet_vertice opp_js.push_back(j); } if (opp_js.size() == 2) { - if (is_surface_fs[t_id][opp_js[0]] != NOT_SURFACE) + if (is_surface_fs[t_id][opp_js[0]] != State::state().NOT_SURFACE) cnt++; - if (is_surface_fs[t_id][opp_js[1]] != NOT_SURFACE) + if (is_surface_fs[t_id][opp_js[1]] != State::state().NOT_SURFACE) cnt++; if (cnt > 2) break; @@ -716,10 +708,8 @@ void SimpleTetrahedralization::labelBoundary(std::vector& tet_vertice if (tet_vertices[i].is_on_surface) cnt_surface++; } -#ifndef MUTE_COUT - cout << cnt_boundary << " vertices on boundary" << endl; - cout << cnt_surface << " vertices on surface" << endl; -#endif + logger().debug("{} vertices on boundary", cnt_boundary); + logger().debug("{} vertices on surface", cnt_surface); } void SimpleTetrahedralization::constructPlane(int bsp_f_id, Plane_3& pln) { @@ -735,4 +725,4 @@ void SimpleTetrahedralization::constructPlane(int bsp_f_id, Plane_3& pln) { assert(!(pln.is_degenerate())); } - +} // namespace tetwild diff --git a/src/tetwild/SimpleTetrahedralization.h b/src/tetwild/SimpleTetrahedralization.h index 38d7225..3563cf7 100644 --- a/src/tetwild/SimpleTetrahedralization.h +++ b/src/tetwild/SimpleTetrahedralization.h @@ -1,9 +1,9 @@ // This file is part of TetWild, a software for generating tetrahedral meshes. -// +// // Copyright (C) 2018 Yixin Hu -// -// This Source Code Form is subject to the terms of the Mozilla Public License -// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can // obtain one at http://mozilla.org/MPL/2.0/. // // Created by Yixin Hu on 3/31/17. @@ -15,6 +15,8 @@ #include #include +namespace tetwild { + class SimpleTetrahedralization { public: MeshConformer& MC; @@ -40,5 +42,6 @@ class SimpleTetrahedralization { void constructPlane(int bsp_f_id, Plane_3& pln); }; +} // namespace tetwild #endif //NEW_GTET_SIMPLETETRAHEDRALIZATION_H diff --git a/src/tetwild/State.cpp b/src/tetwild/State.cpp new file mode 100644 index 0000000..854f6a2 --- /dev/null +++ b/src/tetwild/State.cpp @@ -0,0 +1,10 @@ +// This file is part of TetWild, a software for generating tetrahedral meshes. +// +// Copyright (C) 2018 Jeremie Dumas +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// obtain one at http://mozilla.org/MPL/2.0/. +// +// Created by Jeremie Dumas on 09/04/18. +// diff --git a/src/tetwild/State.h b/src/tetwild/State.h new file mode 100644 index 0000000..378f118 --- /dev/null +++ b/src/tetwild/State.h @@ -0,0 +1,148 @@ +// This file is part of TetWild, a software for generating tetrahedral meshes. +// +// Copyright (C) 2018 Jeremie Dumas +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// obtain one at http://mozilla.org/MPL/2.0/. +// +// Created by Jeremie Dumas on 09/04/18. +// + +#pragma once + +#include + +namespace tetwild { + +struct State { + const int EPSILON_INFINITE=-2; + const int EPSILON_NA=-1; + const int ENERGY_NA=0; + const int ENERGY_AD=1; + const int ENERGY_AMIPS=2; + const int ENERGY_DIRICHLET=3; + const double MAX_ENERGY = 1e50; + int NOT_SURFACE = 0; + + //global parameters + std::string g_working_dir; + std::string g_stat_file; + std::string g_postfix; + std::string g_output_file; + + double g_eps = 0; + double g_eps_2 = 0; + double g_dd = 0; + double g_ideal_l = 0; + double g_diag_l = 0; + bool g_is_close = 0; + + double g_eps_input = 0; + double g_eps_delta = 0; + int g_cur_stage = 1; + + //for test + bool is_using_energy_max = true; + bool is_using_sampling = true; + bool is_use_project = false; + bool is_print_tmp = false; + + static State & state() { + static State st; + return st; + } + +private: + State() = default; +}; + + +struct GArgs { + std::string input; + std::string output = ""; + std::string postfix = "_"; + double i_ideal_edge_length = 20; + double i_epsilon = 1000; + int i_dd = -1; + int stage = 1; + double adaptive_scalar = 0.6; + double filter_energy = 10; + double delta_energy = 0.1; + int max_pass = 80; + int is_output_csv = true; + std::string csv_file = ""; + std::string slz_file = ""; + + int mid_result = -1; + bool is_using_voxel = true; + bool is_laplacian = false; + + int targeted_num_v = -1; + std::string bg_mesh = ""; + + bool is_quiet = false; + + static GArgs & args() { + static GArgs ag; + return ag; + } + +private: + GArgs() = default; +}; + + +struct MeshRecord { + enum OpType { + OP_INIT = 0, + OP_PREPROCESSING, + OP_DELAUNEY_TETRA, + OP_DIVFACE_MATCH, + OP_BSP, + OP_SIMPLE_TETRA, + + OP_OPT_INIT, + OP_SPLIT, + OP_COLLAPSE, + OP_SWAP, + OP_SMOOTH, + OP_ADAP_UPDATE, + OP_WN, + OP_UNROUNDED + }; + + int op; + double timing; + int n_v; + int n_t; + double min_min_d_angle = -1; + double avg_min_d_angle = -1; + double max_max_d_angle = -1; + double avg_max_d_angle = -1; + double max_energy = -1; + double avg_energy = -1; + + MeshRecord(int op, double timing, int n_v, int n_t, double min_min_d_angle, double avg_min_d_angle, + double max_max_d_angle, double avg_max_d_angle, double max_energy, double avg_energy) { + this->op = op; + this->timing = timing; + this->n_v = n_v; + this->n_t = n_t; + this->min_min_d_angle = min_min_d_angle; + this->avg_min_d_angle = avg_min_d_angle; + this->max_max_d_angle = max_max_d_angle; + this->avg_max_d_angle = avg_max_d_angle; + this->max_energy = max_energy; + this->avg_energy = avg_energy; + } + + MeshRecord(int op, double timing, int n_v, int n_t) { + this->op = op; + this->timing = timing; + this->n_v = n_v; + this->n_t = n_t; + } +}; + +} // namespace tetwild diff --git a/src/tetwild/TetmeshElements.h b/src/tetwild/TetmeshElements.h index 2389e26..b890839 100644 --- a/src/tetwild/TetmeshElements.h +++ b/src/tetwild/TetmeshElements.h @@ -1,9 +1,9 @@ // This file is part of TetWild, a software for generating tetrahedral meshes. -// +// // Copyright (C) 2018 Yixin Hu -// -// This Source Code Form is subject to the terms of the Mozilla Public License -// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can // obtain one at http://mozilla.org/MPL/2.0/. // // Created by yihu on 8/22/17. @@ -12,7 +12,9 @@ #ifndef NEW_GTET_TETMESHELEMENTS_H #define NEW_GTET_TETMESHELEMENTS_H -#include +#include + +namespace tetwild { //const int ON_SURFACE_FALSE = 0;//delete //const int ON_SURFACE_TRUE_INSIDE = 1;//delete @@ -53,13 +55,12 @@ class TetVertex { } void printInfo(){ - cout<<"is_on_surface = "< tq.min_d_angle && max_d_angle < tq.max_d_angle; } else return false; } bool isBetterOrEqualThan(const TetQuality& tq, int energy_type) { - if (energy_type == ENERGY_AMIPS || energy_type == ENERGY_DIRICHLET) { + if (energy_type == State::state().ENERGY_AMIPS || energy_type == State::state().ENERGY_DIRICHLET) { return slim_energy <= tq.slim_energy; } - else if (energy_type == ENERGY_AD) { + else if (energy_type == State::state().ENERGY_AD) { return min_d_angle >= tq.min_d_angle && max_d_angle <= tq.max_d_angle; } else @@ -156,10 +157,12 @@ class Stage{ } }; +} // namespace tetwild + namespace igl { namespace serialization { template<> - inline void serialize(const TetVertex &v, std::vector &buffer) { + inline void serialize(const tetwild::TetVertex &v, std::vector &buffer) { ::igl::serialize(v.pos, std::string("pos"), buffer); ::igl::serialize(v.posf, std::string("posf"), buffer); @@ -186,7 +189,7 @@ namespace igl { } template<> - inline void deserialize(TetVertex &v, const std::vector &buffer) { + inline void deserialize(tetwild::TetVertex &v, const std::vector &buffer) { ::igl::deserialize(v.pos, std::string("pos"), buffer); ::igl::deserialize(v.posf, std::string("posf"), buffer); diff --git a/src/tetwild/VertexSmoother.cpp b/src/tetwild/VertexSmoother.cpp index 03336a7..f178f2a 100644 --- a/src/tetwild/VertexSmoother.cpp +++ b/src/tetwild/VertexSmoother.cpp @@ -10,6 +10,9 @@ // #include +#include + +namespace tetwild { void VertexSmoother::smooth() { tets_tss = std::vector(tets.size(), 1); @@ -24,24 +27,20 @@ void VertexSmoother::smooth() { double suc_surface = 0; smoothSingle(); suc_in = suc_counter; - if (g_eps >= 0) { + if (State::state().g_eps >= 0) { smoothSurface(); suc_surface = suc_counter; } -#ifndef MUTE_COUT - cout << (suc_in + suc_surface) / v_cnt << endl; + logger().debug("{}", (suc_in + suc_surface) / v_cnt); if (suc_in + suc_surface < v_cnt * 0.1) { - cout << i << endl; + logger().debug("{}", i); break; } -#endif } -#ifndef MUTE_COUT for (int i = 0; i < breakdown_timing.size(); i++) { - cout << breakdown_name[i] << ": " << breakdown_timing[i] << "s" << endl; + logger().debug("{}: {}s", breakdown_name[i], breakdown_timing[i]); breakdown_timing[i] = 0;//reset } -#endif } bool VertexSmoother::smoothSingleVertex(int v_id, bool is_cal_energy){ @@ -77,7 +76,7 @@ bool VertexSmoother::smoothSingleVertex(int v_id, bool is_cal_energy){ return false; } else { Point_3f pf; - if(energy_type == ENERGY_AMIPS) { + if(energy_type == State::state().ENERGY_AMIPS) { if (!NewtonsMethod(t_ids, new_tets, v_id, pf)) return false; } @@ -91,7 +90,7 @@ bool VertexSmoother::smoothSingleVertex(int v_id, bool is_cal_energy){ tet_vertices[v_id].posf = pf; tet_vertices[v_id].is_rounded = true; if (isFlip(new_tets)) {//TODO: why it happens? - cout << "flip in the end" << endl; + logger().debug("flip in the end"); tet_vertices[v_id].pos = old_p; tet_vertices[v_id].posf = old_pf; tet_vertices[v_id].is_rounded = old_is_rounded; @@ -119,7 +118,7 @@ void VertexSmoother::smoothSingle() { continue; if (tet_vertices[v_id].is_on_bbox) continue; - if (g_eps != EPSILON_INFINITE && tet_vertices[v_id].is_on_surface) + if (State::state().g_eps != State::state().EPSILON_INFINITE && tet_vertices[v_id].is_on_surface) continue; if (tet_vertices[v_id].is_locked) @@ -177,7 +176,7 @@ void VertexSmoother::smoothSingle() { continue; } else { Point_3f pf; - if (energy_type == ENERGY_AMIPS) { + if (energy_type == State::state().ENERGY_AMIPS) { if (!NewtonsMethod(t_ids, new_tets, v_id, pf)) continue; } @@ -193,7 +192,7 @@ void VertexSmoother::smoothSingle() { tet_vertices[v_id].posf = pf; tet_vertices[v_id].is_rounded = true; if (isFlip(new_tets)) {//TODO: why it happens? - cout << "flip in the end" << endl; + logger().debug("flip in the end"); tet_vertices[v_id].pos = old_p; tet_vertices[v_id].posf = old_pf; tet_vertices[v_id].is_rounded = old_is_rounded; @@ -301,7 +300,7 @@ void VertexSmoother::smoothSurface() {//smoothing surface using two methods if (!is_valid) { continue; } else { - if (energy_type == ENERGY_AMIPS) { + if (energy_type == State::state().ENERGY_AMIPS) { if (!NewtonsMethod(old_t_ids, new_tets, v_id, pf_out)) continue; } @@ -315,7 +314,7 @@ void VertexSmoother::smoothSurface() {//smoothing surface using two methods std::vector> tri_ids; for (auto it = tet_vertices[v_id].conn_tets.begin(); it != tet_vertices[v_id].conn_tets.end(); it++) { for (int j = 0; j < 4; j++) { - if (tets[*it][j] != v_id && is_surface_fs[*it][j] != NOT_SURFACE) { + if (tets[*it][j] != v_id && is_surface_fs[*it][j] != State::state().NOT_SURFACE) { std::array tri = {tets[*it][(j + 1) % 4], tets[*it][(j + 2) % 4], tets[*it][(j + 3) % 4]}; std::sort(tri.begin(), tri.end()); tri_ids.push_back(tri); @@ -327,7 +326,7 @@ void VertexSmoother::smoothSurface() {//smoothing surface using two methods Point_3f pf; Point_3 p; - if (is_use_project) {//we have to use exact construction here. Or the projecting points may be not exactly on the plane. + if (State::state().is_use_project) {//we have to use exact construction here. Or the projecting points may be not exactly on the plane. std::vector tris; for (int i = 0; i < tri_ids.size(); i++) { tris.push_back(Triangle_3(tet_vertices[tri_ids[i][0]].pos, tet_vertices[tri_ids[i][1]].pos, @@ -455,14 +454,10 @@ void VertexSmoother::smoothSurface() {//smoothing surface using two methods suc_counter++; sf_suc_counter++; -#ifndef MUTE_COUT if (sf_suc_counter % 1000 == 0) - cout << "1000 accepted!" << endl; -#endif + logger().debug("1000 accepted!"); } -#ifndef MUTE_COUT - cout << "Totally " << sf_suc_counter << "(" << sf_counter << ")" << " vertices on surface are smoothed." << endl; -#endif + logger().debug("Totally {}({}) vertices on surface are smoothed.", sf_suc_counter, sf_counter); } bool VertexSmoother::NewtonsMethod(const std::vector& t_ids, const std::vector>& new_tets, @@ -567,7 +562,7 @@ double VertexSmoother::getNewEnergy(const std::vector& t_ids) { static double* energy = 0; if (T0 == 0) { - std::cerr << "Initial ISPC allocation: n = " << n << endl; + logger().warn("Initial ISPC allocation: n = {}", n); current_max_size = n; T0 = new double[n]; T1 = new double[n]; @@ -585,7 +580,7 @@ double VertexSmoother::getNewEnergy(const std::vector& t_ids) { } if (current_max_size < n) { - std::cerr << "ISPC reallocation: n = " << n << endl; + logger().warn("ISPC reallocation: n = {}", n); free(T0); free(T1); free(T2); @@ -644,16 +639,14 @@ double VertexSmoother::getNewEnergy(const std::vector& t_ids) { for (int k = 0; k < 3; k++) t.push_back(tet_vertices[tets[t_ids[i]][j]].posf[k]); } - if (energy_type == ENERGY_AMIPS) { + if (energy_type == State::state().ENERGY_AMIPS) { s_energy += comformalAMIPSEnergy_new(t); } } #endif - if (std::isinf(s_energy) || std::isnan(s_energy) || s_energy <= 0 || s_energy > MAX_ENERGY) { -#ifndef MUTE_COUT - cout << "new E inf" << endl; -#endif - s_energy = MAX_ENERGY; + if (std::isinf(s_energy) || std::isnan(s_energy) || s_energy <= 0 || s_energy > State::state().MAX_ENERGY) { + logger().debug("new E inf"); + s_energy = State::state().MAX_ENERGY; } return s_energy; @@ -712,33 +705,23 @@ bool VertexSmoother::NewtonsUpdate(const std::vector& t_ids, int v_id, #endif if (std::isinf(energy)) { -#ifndef MUTE_COUT - cout << v_id << " E inf" << endl; -#endif - energy = MAX_ENERGY; + logger().debug("{} E inf", v_id); + energy = State::state().MAX_ENERGY; } if (std::isnan(energy)) { -#ifndef MUTE_COUT - cout << v_id << " E nan" << endl; -#endif + logger().debug("{} E nan", v_id); return false; } if (energy <= 0) { -#ifndef MUTE_COUT - cout << v_id << " E < 0" << endl; -#endif + logger().debug("{} E < 0", v_id); return false; } if (!J.allFinite()) { -#ifndef MUTE_COUT - cout << v_id << " J inf/nan" << endl; -#endif + logger().debug("{} J inf/nan", v_id); return false; } if (!H.allFinite()) { -#ifndef MUTE_COUT - cout << v_id << " H inf/nan" << endl; -#endif + logger().debug("{} H inf/nan", v_id); return false; } @@ -813,7 +796,7 @@ int VertexSmoother::laplacianBoundary(const std::vector& b_v_ids, const std //give stop condition bool is_stop = true; for (int j = 0; j < 3; j++) - if (vec[j] * a > g_eps) + if (vec[j] * a > State::state().g_eps) is_stop = false; if (is_stop) break; @@ -859,8 +842,8 @@ int VertexSmoother::laplacianBoundary(const std::vector& b_v_ids, const std } // do normal smoothing on neighbor vertices -// cout<<"n_v_ids.size = "<& b_v_ids, const std // } } -#ifndef MUTE_COUT - cout<<"suc.size = "< v_ids; std::vector new_ids(tet_vertices.size(), -1); for(int t_id: tet_vertices[v_id].conn_tets){ @@ -925,3 +906,4 @@ void VertexSmoother::outputOneRing(int v_id, std::string s){ mSaver.save_elem_scalar_field("quality", q); } +} // namespace tetwild diff --git a/src/tetwild/VertexSmoother.h b/src/tetwild/VertexSmoother.h index afd406a..69b3ac8 100644 --- a/src/tetwild/VertexSmoother.h +++ b/src/tetwild/VertexSmoother.h @@ -1,9 +1,9 @@ // This file is part of TetWild, a software for generating tetrahedral meshes. -// +// // Copyright (C) 2018 Yixin Hu -// -// This Source Code Form is subject to the terms of the Mozilla Public License -// v. 2.0. If a copy of the MPL was not distributed with this file, You can +// +// This Source Code Form is subject to the terms of the Mozilla Public License +// v. 2.0. If a copy of the MPL was not distributed with this file, You can // obtain one at http://mozilla.org/MPL/2.0/. // // Created by Yixin Hu on 4/11/17. @@ -14,6 +14,8 @@ #include +namespace tetwild { + class VertexSmoother:public LocalOperations { public: VertexSmoother(LocalOperations& lo): LocalOperations(lo){} @@ -48,5 +50,6 @@ class VertexSmoother:public LocalOperations { igl::Timer igl_timer; }; +} // namespace tetwild #endif //NEW_GTET_VERTEXSMOOTHER_H diff --git a/src/tetwild/heads.h b/src/tetwild/heads.h deleted file mode 100644 index 68e668c..0000000 --- a/src/tetwild/heads.h +++ /dev/null @@ -1,277 +0,0 @@ -// This file is part of TetWild, a software for generating tetrahedral meshes. -// -// Copyright (C) 2018 Yixin Hu -// -// This Source Code Form is subject to the terms of the Mozilla Public License -// v. 2.0. If a copy of the MPL was not distributed with this file, You can -// obtain one at http://mozilla.org/MPL/2.0/. -// -// Created by Yixin Hu on 3/29/17. -// - -#ifndef GTET_HEADS_H -#define GTET_HEADS_H - -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include - -#include -#include -#include -#include -#include -#include - -// PyMesh -#include -#include - -using std::cerr; -using std::cout; -using std::cin; -using std::endl; - -#include -#include -#include -#include - -typedef CGAL::Exact_predicates_exact_constructions_kernel K; -typedef K::Point_2 Point_2; -typedef K::Segment_2 Segment_2; -typedef K::Line_2 Line_2; -typedef K::Iso_rectangle_2 Iso_rectangle_2; -typedef K::Triangle_2 Triangle_2; -typedef K::Intersect_2 Intersect_2; -//typedef CGAL::Polygon_2 Polygon_2; - -typedef K::Point_3 Point_3; -typedef K::Vector_3 Vector_3; -typedef K::Segment_3 Segment_3; -typedef K::Line_3 Line_3; -typedef K::Plane_3 Plane_3; -typedef K::Triangle_3 Triangle_3; -typedef K::Intersect_3 Intersect_3; -typedef K::Tetrahedron_3 Tetrahedron_3; -typedef K::Direction_3 Direction_3; - -typedef CGAL::Exact_predicates_inexact_constructions_kernel Kf; -typedef Kf::Point_3 Point_3f; -typedef Kf::Vector_3 Vector_3f; -typedef Kf::Plane_3 Plane_3f; -typedef Kf::Triangle_3 Triangle_3f; -typedef Kf::Segment_3 Segment_3f; -typedef Kf::Line_3 Line_3f; - -typedef CGAL::Epeck::FT CGAL_FT; -//#include -//typedef CGAL::Simple_cartesian::FT CGAL_FT; - -#include -#include -#include - -#include -typedef K::Iso_cuboid_3 Bbox_3; - -#include - -void pausee(); -bool isHaveCommonEle(const std::unordered_set& v1, const std::unordered_set& v2); -void setIntersection(const std::unordered_set& s1, const std::unordered_set& s2, std::unordered_set& s); -void setIntersection(const std::unordered_set& s1, const std::unordered_set& s2, std::vector& s); -void sampleTriangle(const std::array& vs, std::vector& ps); - -#define TIMING_BREAKDOWN true - -const int EPSILON_INFINITE=-2; -const int EPSILON_NA=-1; - -const int ENERGY_NA=0; -const int ENERGY_AD=1; -const int ENERGY_AMIPS=2; -const int ENERGY_DIRICHLET=3; - -const double MAX_ENERGY = 1e50; - -extern int NOT_SURFACE; - -//global parameters -extern std::string g_working_dir; -extern std::string g_stat_file; -extern std::string g_postfix; -extern std::string g_output_file; - -struct Args{ - std::string input; - std::string output = ""; - std::string postfix = "_"; - double i_ideal_edge_length = 20; - double i_epsilon = 1000; - int i_dd = -1; - int stage = 1; - double adaptive_scalar = 0.6; - double filter_energy = 10; - double delta_energy = 0.1; - int max_pass = 80; - int is_output_csv = true; - std::string csv_file = ""; - std::string slz_file = ""; - - int mid_result = -1; - bool is_using_voxel = true; - bool is_laplacian = false; - - int targeted_num_v = -1; - std::string bg_mesh = ""; - - bool is_quiet = false; -}; - -extern Args args; - -extern double g_eps; -extern double g_eps_2; -extern double g_dd; -extern double g_ideal_l; -extern double g_diag_l; -extern bool g_is_close; - -extern double g_eps_input; -extern double g_eps_delta; -extern int g_cur_stage; - -//for test -extern bool is_using_energy_max; -extern bool is_using_sampling; -extern bool is_use_project; - -extern bool is_print_tmp; - -struct MeshRecord { - enum OpType { - OP_INIT = 0, - OP_PREPROCESSING, - OP_DELAUNEY_TETRA, - OP_DIVFACE_MATCH, - OP_BSP, - OP_SIMPLE_TETRA, - - OP_OPT_INIT, - OP_SPLIT, - OP_COLLAPSE, - OP_SWAP, - OP_SMOOTH, - OP_ADAP_UPDATE, - OP_WN, - OP_UNROUNDED - }; - int op; - double timing; - int n_v; - int n_t; - double min_min_d_angle = -1; - double avg_min_d_angle = -1; - double max_max_d_angle = -1; - double avg_max_d_angle = -1; - double max_energy = -1; - double avg_energy = -1; - MeshRecord(int op, double timing, int n_v, int n_t, double min_min_d_angle, double avg_min_d_angle, - double max_max_d_angle, double avg_max_d_angle, double max_energy, double avg_energy) { - this->op = op; - this->timing = timing; - this->n_v = n_v; - this->n_t = n_t; - this->min_min_d_angle = min_min_d_angle; - this->avg_min_d_angle = avg_min_d_angle; - this->max_max_d_angle = max_max_d_angle; - this->avg_max_d_angle = avg_max_d_angle; - this->max_energy = max_energy; - this->avg_energy = avg_energy; - } - MeshRecord(int op, double timing, int n_v, int n_t) { - this->op = op; - this->timing = timing; - this->n_v = n_v; - this->n_t = n_t; - } -}; -//extern std::vector mesh_records; -extern bool is_app_csv; -void addRecord(const MeshRecord& record); - -//for serialization -namespace igl { - namespace serialization { - template<> - inline void serialize(const Point_3 &p, std::vector &buffer) { - ::igl::serialize(STR(CGAL::exact(p[0])), std::string("x"), buffer); - ::igl::serialize(STR(CGAL::exact(p[1])), std::string("y"), buffer); - ::igl::serialize(STR(CGAL::exact(p[2])), std::string("z"), buffer); - } - - template<> - inline void deserialize(Point_3 &p, const std::vector &buffer) { - std::string s1, s2, s3; - ::igl::deserialize(s1, std::string("x"), buffer); - ::igl::deserialize(s2, std::string("y"), buffer); - ::igl::deserialize(s3, std::string("z"), buffer); -// p=Point_3(CGAL_FT(s1), CGAL_FT(s2), CGAL_FT(s3)); - } - - template<> - inline void serialize(const Point_3f &p, std::vector &buffer) { - ::igl::serialize(p[0], std::string("x"), buffer); - ::igl::serialize(p[1], std::string("y"), buffer); - ::igl::serialize(p[2], std::string("z"), buffer); - } - - template<> - inline void deserialize(Point_3f &p, const std::vector &buffer) { - double x, y, z; - ::igl::deserialize(x, std::string("x"), buffer); - ::igl::deserialize(y, std::string("y"), buffer); - ::igl::deserialize(z, std::string("z"), buffer); - p=Point_3f(x, y, z); - } - - template<> - inline void serialize(const std::array &arr, std::vector &buffer) { - for(int i=0;i<3;i++) - ::igl::serialize(arr[i], std::to_string(i), buffer); - } - - template<> - inline void deserialize(std::array &arr, const std::vector &buffer) { - for(int i=0;i<3;i++) - ::igl::deserialize(arr[i], std::to_string(i), buffer); - } - - template<> - inline void serialize(const std::array &arr, std::vector &buffer) { - for(int i=0;i<4;i++) - ::igl::serialize(arr[i], std::to_string(i), buffer); - } - - template<> - inline void deserialize(std::array &arr, const std::vector &buffer) { - for(int i=0;i<4;i++) - ::igl::deserialize(arr[i], std::to_string(i), buffer); - } - } -} -#endif //GTET_HEADS_H diff --git a/src/tetwild/tetwild.cpp b/src/tetwild/tetwild.cpp index 67bdd37..9f83a2d 100644 --- a/src/tetwild/tetwild.cpp +++ b/src/tetwild/tetwild.cpp @@ -8,6 +8,7 @@ // // Created by Yixin Hu on 5/31/18. // + #include #include #include @@ -15,351 +16,306 @@ #include #include #include +#include namespace tetwild { + // MeshRefinement MR; - Args parameters; - - void outputFinalQuality(double time, const std::vector& tet_vertices, const std::vector>& tets, - const std::vector &t_is_removed, const std::vector& tet_qualities, - const std::vector& v_ids) { - cout << "final quality:" << endl; - double min = 10, max = 0; - double min_avg = 0, max_avg = 0; + +void outputFinalQuality(double time, const std::vector& tet_vertices, const std::vector>& tets, + const std::vector &t_is_removed, const std::vector& tet_qualities, + const std::vector& v_ids) { + logger().debug("final quality:"); + double min = 10, max = 0; + double min_avg = 0, max_avg = 0; // double max_asp_ratio = 0, avg_asp_ratio = 0; - double max_slim_energy = 0, avg_slim_energy = 0; - std::array cmp_cnt = {0, 0, 0, 0, 0, 0}; - std::array cmp_d_angles = {6 / 180.0 * M_PI, 12 / 180.0 * M_PI, 18 / 180.0 * M_PI, - 162 / 180.0 * M_PI, 168 / 180.0 * M_PI, 174 / 180.0 * M_PI}; - int cnt = 0; - for (int i = 0; i < tet_qualities.size(); i++) { - if (t_is_removed[i]) - continue; - cnt++; - if (tet_qualities[i].min_d_angle < min) - min = tet_qualities[i].min_d_angle; - if (tet_qualities[i].max_d_angle > max) - max = tet_qualities[i].max_d_angle; + double max_slim_energy = 0, avg_slim_energy = 0; + std::array cmp_cnt = {0, 0, 0, 0, 0, 0}; + std::array cmp_d_angles = {6 / 180.0 * M_PI, 12 / 180.0 * M_PI, 18 / 180.0 * M_PI, + 162 / 180.0 * M_PI, 168 / 180.0 * M_PI, 174 / 180.0 * M_PI}; + int cnt = 0; + for (int i = 0; i < tet_qualities.size(); i++) { + if (t_is_removed[i]) + continue; + cnt++; + if (tet_qualities[i].min_d_angle < min) + min = tet_qualities[i].min_d_angle; + if (tet_qualities[i].max_d_angle > max) + max = tet_qualities[i].max_d_angle; // if (tet_qualities[i].asp_ratio_2 > max_asp_ratio) // max_asp_ratio = tet_qualities[i].asp_ratio_2; - if (tet_qualities[i].slim_energy > max_slim_energy) - max_slim_energy = tet_qualities[i].slim_energy; - min_avg += tet_qualities[i].min_d_angle; - max_avg += tet_qualities[i].max_d_angle; + if (tet_qualities[i].slim_energy > max_slim_energy) + max_slim_energy = tet_qualities[i].slim_energy; + min_avg += tet_qualities[i].min_d_angle; + max_avg += tet_qualities[i].max_d_angle; // avg_asp_ratio += tet_qualities[i].asp_ratio_2; - avg_slim_energy += tet_qualities[i].slim_energy; - - for (int j = 0; j < 3; j++) { - if (tet_qualities[i].min_d_angle < cmp_d_angles[j]) - cmp_cnt[j]++; - } - for (int j = 0; j < 3; j++) { - if (tet_qualities[i].max_d_angle > cmp_d_angles[j + 3]) - cmp_cnt[j + 3]++; - } + avg_slim_energy += tet_qualities[i].slim_energy; + + for (int j = 0; j < 3; j++) { + if (tet_qualities[i].min_d_angle < cmp_d_angles[j]) + cmp_cnt[j]++; } - cout << "min_d_angle = " << min - << ", max_d_angle = " << max - // << ", max_aspect_ratio = " << max_asp_ratio - << ", max_slim_energy = " << max_slim_energy - << endl; - cout << "avg_min_d_angle = " << min_avg / cnt - << ", avg_max_d_angle = " << max_avg / cnt - // << ", avg_aspect_ratio = " << avg_asp_ratio / cnt - << ", avg_slim_energy = " << avg_slim_energy / cnt - << endl; - cout << "min_d_angle: <6 " << cmp_cnt[0] / cnt << "; <12 " << cmp_cnt[1] / cnt << "; <18 " << cmp_cnt[2] / cnt - << endl; - cout << "max_d_angle: >174 " << cmp_cnt[5] / cnt << "; >168 " << cmp_cnt[4] / cnt << "; >162 " << cmp_cnt[3] / cnt - << endl; - - addRecord(MeshRecord(MeshRecord::OpType::OP_WN, time, v_ids.size(), cnt, - min, min_avg / cnt, max, max_avg / cnt, max_slim_energy, avg_slim_energy / cnt)); - - ////output unrounded vertices: - cnt = 0; - for (int v_id: v_ids) { - if (!tet_vertices[v_id].is_rounded) - cnt++; + for (int j = 0; j < 3; j++) { + if (tet_qualities[i].max_d_angle > cmp_d_angles[j + 3]) + cmp_cnt[j + 3]++; } - cout << cnt << "/" << v_ids.size() << " vertices are unrounded!!!" << endl; - addRecord(MeshRecord(MeshRecord::OpType::OP_UNROUNDED, -1, cnt, -1)); } + logger().debug("min_d_angle = {}, max_d_angle = {}, max_slim_energy = {}", min, max, max_slim_energy); + logger().debug("avg_min_d_angle = {}, avg_max_d_angle = {}, avg_slim_energy = {}", min_avg / cnt, max_avg / cnt, avg_slim_energy / cnt); + logger().debug("min_d_angle: <6 {}; <12 {}; <18 {}", cmp_cnt[0] / cnt, cmp_cnt[1] / cnt, cmp_cnt[2] / cnt); + logger().debug("max_d_angle: >174 {}; >168 {}; >162 {}", cmp_cnt[5] / cnt, cmp_cnt[4] / cnt, cmp_cnt[3] / cnt); - void outputFinalTetmesh(MeshRefinement& MR, - std::vector>& out_vertices, - std::vector>& out_tets) { - std::vector &tet_vertices = MR.tet_vertices; - std::vector> &tets = MR.tets; - std::vector &v_is_removed = MR.v_is_removed; - std::vector &t_is_removed = MR.t_is_removed; - std::vector &tet_qualities = MR.tet_qualities; - int t_cnt = std::count(t_is_removed.begin(), t_is_removed.end(), false); - double tmp_time = 0; - if (!args.is_laplacian) { - InoutFiltering IOF(tet_vertices, tets, MR.is_surface_fs, v_is_removed, t_is_removed, tet_qualities); -#ifndef MUTE_COUT - igl::Timer igl_timer; - igl_timer.start(); -#endif - IOF.filter(); -#ifndef MUTE_COUT - tmp_time = igl_timer.getElapsedTime(); - cout << "time = " << tmp_time << "s" << endl; - t_cnt = std::count(t_is_removed.begin(), t_is_removed.end(), false); - cout << t_cnt << " tets inside!" << endl; -#endif - } + addRecord(MeshRecord(MeshRecord::OpType::OP_WN, time, v_ids.size(), cnt, + min, min_avg / cnt, max, max_avg / cnt, max_slim_energy, avg_slim_energy / cnt)); - //output result - std::vector v_ids; - for (int i = 0; i < tets.size(); i++) { - if (t_is_removed[i]) - continue; - for (int j = 0; j < 4; j++) - v_ids.push_back(tets[i][j]); - } - std::sort(v_ids.begin(), v_ids.end()); - v_ids.erase(std::unique(v_ids.begin(), v_ids.end()), v_ids.end()); - std::unordered_map map_ids; - for (int i = 0; i < v_ids.size(); i++) - map_ids[v_ids[i]] = i; - - out_vertices.reserve(v_ids.size()); - for (int i = 0; i < v_ids.size(); i++) { - out_vertices.push_back(std::array({tet_vertices[v_ids[i]].posf[0], - tet_vertices[v_ids[i]].posf[1], - tet_vertices[v_ids[i]].posf[2]})); - } - out_tets.reserve(std::count(t_is_removed.begin(), t_is_removed.end(), false)); - for (int i = 0; i < tets.size(); i++) { - if (t_is_removed[i]) - continue; - out_tets.push_back(std::array({map_ids[tets[i][0]], map_ids[tets[i][1]], map_ids[tets[i][2]], - map_ids[tets[i][3]]})); - } + ////output unrounded vertices: + cnt = 0; + for (int v_id: v_ids) { + if (!tet_vertices[v_id].is_rounded) + cnt++; + } + logger().debug("{}/{} vertices are unrounded!!!", cnt, v_ids.size()); + addRecord(MeshRecord(MeshRecord::OpType::OP_UNROUNDED, -1, cnt, -1)); +} - if(args.is_quiet) - return; +void outputFinalTetmesh(MeshRefinement& MR, + std::vector>& out_vertices, + std::vector>& out_tets) { + std::vector &tet_vertices = MR.tet_vertices; + std::vector> &tets = MR.tets; + std::vector &v_is_removed = MR.v_is_removed; + std::vector &t_is_removed = MR.t_is_removed; + std::vector &tet_qualities = MR.tet_qualities; + int t_cnt = std::count(t_is_removed.begin(), t_is_removed.end(), false); + double tmp_time = 0; + if (!GArgs::args().is_laplacian) { + InoutFiltering IOF(tet_vertices, tets, MR.is_surface_fs, v_is_removed, t_is_removed, tet_qualities); + igl::Timer igl_timer; + igl_timer.start(); + IOF.filter(); + tmp_time = igl_timer.getElapsedTime(); + logger().debug("time = {}s", tmp_time); + t_cnt = std::count(t_is_removed.begin(), t_is_removed.end(), false); + logger().debug("{} tets inside!", t_cnt); + } + + //output result + std::vector v_ids; + for (int i = 0; i < tets.size(); i++) { + if (t_is_removed[i]) + continue; + for (int j = 0; j < 4; j++) + v_ids.push_back(tets[i][j]); + } + std::sort(v_ids.begin(), v_ids.end()); + v_ids.erase(std::unique(v_ids.begin(), v_ids.end()), v_ids.end()); + std::unordered_map map_ids; + for (int i = 0; i < v_ids.size(); i++) + map_ids[v_ids[i]] = i; -#ifndef MUTE_COUT - outputFinalQuality(tmp_time, tet_vertices, tets, t_is_removed, tet_qualities, v_ids); -#endif + out_vertices.reserve(v_ids.size()); + for (int i = 0; i < v_ids.size(); i++) { + out_vertices.push_back(std::array({tet_vertices[v_ids[i]].posf[0], + tet_vertices[v_ids[i]].posf[1], + tet_vertices[v_ids[i]].posf[2]})); + } + out_tets.reserve(std::count(t_is_removed.begin(), t_is_removed.end(), false)); + for (int i = 0; i < tets.size(); i++) { + if (t_is_removed[i]) + continue; + out_tets.push_back(std::array({map_ids[tets[i][0]], map_ids[tets[i][1]], map_ids[tets[i][2]], + map_ids[tets[i][3]]})); } - void gtet_new(const Eigen::MatrixXd& V_in, const Eigen::MatrixXi& F_in, - std::vector>& out_vertices, - std::vector>& out_tets) { - is_using_energy_max = true; - is_use_project = false; - is_using_sampling = true; + if(GArgs::args().is_quiet) + return; - int energy_type = ENERGY_AMIPS; - bool is_sm_single = true; - bool is_preprocess = true; - bool is_check_correctness = false; - bool is_ec_check_quality = true; + outputFinalQuality(tmp_time, tet_vertices, tets, t_is_removed, tet_qualities, v_ids); +} -#ifndef MUTE_COUT - igl::Timer igl_timer; - double tmp_time = 0; - double sum_time = 0; -#endif - - ////pipeline - MeshRefinement MR; - {/// STAGE 1 - //preprocess -#ifndef MUTE_COUT - igl_timer.start(); - cout << "Preprocessing..." << endl; -#endif - Preprocess pp; - if (!pp.init(V_in, F_in, MR.geo_b_mesh, MR.geo_sf_mesh)) { - cout << "Empty!" << endl; - //todo: output a empty tetmesh - PyMesh::MshSaver mSaver(g_working_dir + g_postfix + ".msh", true); - Eigen::VectorXd oV; - Eigen::VectorXi oT; - oV.resize(0); - oT.resize(0); - mSaver.save_mesh(oV, oT, 3, mSaver.TET); - exit(250); - } -#ifndef MUTE_COUT - addRecord(MeshRecord(MeshRecord::OpType::OP_INIT, 0, MR.geo_sf_mesh.vertices.nb(), MR.geo_sf_mesh.facets.nb())); -#endif - - std::vector m_vertices; - std::vector> m_faces; - pp.process(MR.geo_sf_mesh, m_vertices, m_faces); -#ifndef MUTE_COUT - tmp_time = igl_timer.getElapsedTime(); - addRecord(MeshRecord(MeshRecord::OpType::OP_PREPROCESSING, tmp_time, m_vertices.size(), m_faces.size())); - sum_time += tmp_time; - cout << "time = " << tmp_time << "s" << endl; - cout << endl; -#endif - - //delaunay tetrahedralization -#ifndef MUTE_COUT - igl_timer.start(); - cout << "Delaunay tetrahedralizing..." << endl; -#endif - DelaunayTetrahedralization DT; - std::vector raw_e_tags; - std::vector> raw_conn_e4v; - std::vector m_f_tags;//need to use it in ST - DT.init(m_vertices, m_faces, m_f_tags, raw_e_tags, raw_conn_e4v); - std::vector bsp_vertices; - std::vector bsp_edges; - std::vector bsp_faces; - std::vector bsp_nodes; - DT.tetra(m_vertices, MR.geo_sf_mesh, bsp_vertices, bsp_edges, bsp_faces, bsp_nodes); -#ifndef MUTE_COUT - cout << "# bsp_vertices = " << bsp_vertices.size() << endl; - cout << "# bsp_edges = " << bsp_edges.size() << endl; - cout << "# bsp_faces = " << bsp_faces.size() << endl; - cout << "# bsp_nodes = " << bsp_nodes.size() << endl; - cout << "Delaunay tetrahedralization done!" << endl; - tmp_time = igl_timer.getElapsedTime(); - addRecord(MeshRecord(MeshRecord::OpType::OP_DELAUNEY_TETRA, tmp_time, bsp_vertices.size(), bsp_nodes.size())); - sum_time += tmp_time; - cout << "time = " << tmp_time << "s" << endl; - cout << endl; -#endif - - //mesh conforming -#ifndef MUTE_COUT - igl_timer.start(); - cout << "Divfaces matching..." << endl; -#endif - MeshConformer MC(m_vertices, m_faces, bsp_vertices, bsp_edges, bsp_faces, bsp_nodes); - MC.match(); -#ifndef MUTE_COUT - cout << "Divfaces matching done!" << endl; - tmp_time = igl_timer.getElapsedTime(); - addRecord(MeshRecord(MeshRecord::OpType::OP_DIVFACE_MATCH, tmp_time, bsp_vertices.size(), bsp_nodes.size())); - cout << "time = " << tmp_time << "s" << endl; - cout << endl; -#endif - - //bsp subdivision -#ifndef MUTE_COUT - igl_timer.start(); - cout << "BSP subdivision ..." << endl; -#endif - BSPSubdivision BS(MC); - BS.init(); - BS.subdivideBSPNodes(); -#ifndef MUTE_COUT - cout << "Output: " << endl; - cout << "# node = " << MC.bsp_nodes.size() << endl; - cout << "# face = " << MC.bsp_faces.size() << endl; - cout << "# edge = " << MC.bsp_edges.size() << endl; - cout << "# vertex = " << MC.bsp_vertices.size() << endl; - cout << "BSP subdivision done!" << endl; - tmp_time = igl_timer.getElapsedTime(); - addRecord(MeshRecord(MeshRecord::OpType::OP_BSP, tmp_time, bsp_vertices.size(), bsp_nodes.size())); - sum_time += tmp_time; - cout << "time = " << tmp_time << "s" << endl; - cout << endl; -#endif - - //simple tetrahedralization -#ifndef MUTE_COUT - igl_timer.start(); - cout << "Tetrehedralizing ..." << endl; -#endif - SimpleTetrahedralization ST(MC); - ST.tetra(MR.tet_vertices, MR.tets); - ST.labelSurface(m_f_tags, raw_e_tags, raw_conn_e4v, MR.tet_vertices, MR.tets, MR.is_surface_fs); - ST.labelBbox(MR.tet_vertices, MR.tets); - if (!g_is_close)//if input is an open mesh - ST.labelBoundary(MR.tet_vertices, MR.tets, MR.is_surface_fs); -#ifndef MUTE_COUT - cout << "# tet_vertices = " << MR.tet_vertices.size() << endl; - cout << "# tets = " << MR.tets.size() << endl; - cout << "Tetrahedralization done!" << endl; - tmp_time = igl_timer.getElapsedTime(); - addRecord(MeshRecord(MeshRecord::OpType::OP_SIMPLE_TETRA, tmp_time, MR.tet_vertices.size(), MR.tets.size())); - sum_time += tmp_time; - cout << "time = " << tmp_time << "s" << endl; - cout << endl; - - cout << "Total time for the first stage = " << sum_time << endl; -#endif +void gtet_new(const Eigen::MatrixXd& V_in, const Eigen::MatrixXi& F_in, + std::vector>& out_vertices, + std::vector>& out_tets) { + State::state().is_using_energy_max = true; + State::state().is_use_project = false; + State::state().is_using_sampling = true; + + int energy_type = State::state().ENERGY_AMIPS; + bool is_sm_single = true; + bool is_preprocess = true; + bool is_check_correctness = false; + bool is_ec_check_quality = true; + + igl::Timer igl_timer; + double tmp_time = 0; + double sum_time = 0; + + ////pipeline + MeshRefinement MR; + {/// STAGE 1 + //preprocess + igl_timer.start(); + logger().info("Preprocessing..."); + Preprocess pp; + if (!pp.init(V_in, F_in, MR.geo_b_mesh, MR.geo_sf_mesh)) { + logger().debug("Empty!"); + //todo: output a empty tetmesh + PyMesh::MshSaver mSaver(State::state().g_working_dir + State::state().g_postfix + ".msh", true); + Eigen::VectorXd oV; + Eigen::VectorXi oT; + oV.resize(0); + oT.resize(0); + mSaver.save_mesh(oV, oT, 3, mSaver.TET); + exit(250); } + addRecord(MeshRecord(MeshRecord::OpType::OP_INIT, 0, MR.geo_sf_mesh.vertices.nb(), MR.geo_sf_mesh.facets.nb())); + + std::vector m_vertices; + std::vector> m_faces; + pp.process(MR.geo_sf_mesh, m_vertices, m_faces); + tmp_time = igl_timer.getElapsedTime(); + addRecord(MeshRecord(MeshRecord::OpType::OP_PREPROCESSING, tmp_time, m_vertices.size(), m_faces.size())); + sum_time += tmp_time; + logger().debug("time = {}s", tmp_time); + + //delaunay tetrahedralization + igl_timer.start(); + logger().info("Delaunay tetrahedralizing..."); + DelaunayTetrahedralization DT; + std::vector raw_e_tags; + std::vector> raw_conn_e4v; + std::vector m_f_tags;//need to use it in ST + DT.init(m_vertices, m_faces, m_f_tags, raw_e_tags, raw_conn_e4v); + std::vector bsp_vertices; + std::vector bsp_edges; + std::vector bsp_faces; + std::vector bsp_nodes; + DT.tetra(m_vertices, MR.geo_sf_mesh, bsp_vertices, bsp_edges, bsp_faces, bsp_nodes); + logger().debug("# bsp_vertices = {}", bsp_vertices.size()); + logger().debug("# bsp_edges = {}", bsp_edges.size()); + logger().debug("# bsp_faces = {}", bsp_faces.size()); + logger().debug("# bsp_nodes = {}", bsp_nodes.size()); + logger().info("Delaunay tetrahedralization done!"); + tmp_time = igl_timer.getElapsedTime(); + addRecord(MeshRecord(MeshRecord::OpType::OP_DELAUNEY_TETRA, tmp_time, bsp_vertices.size(), bsp_nodes.size())); + sum_time += tmp_time; + logger().info("time = {}s", tmp_time); + + //mesh conforming + igl_timer.start(); + logger().info("Divfaces matching..."); + MeshConformer MC(m_vertices, m_faces, bsp_vertices, bsp_edges, bsp_faces, bsp_nodes); + MC.match(); + logger().info("Divfaces matching done!"); + tmp_time = igl_timer.getElapsedTime(); + addRecord(MeshRecord(MeshRecord::OpType::OP_DIVFACE_MATCH, tmp_time, bsp_vertices.size(), bsp_nodes.size())); + logger().info("time = {}s", tmp_time); + + //bsp subdivision + igl_timer.start(); + logger().info("BSP subdivision ..."); + BSPSubdivision BS(MC); + BS.init(); + BS.subdivideBSPNodes(); + logger().debug("Output: "); + logger().debug("# node = {}", MC.bsp_nodes.size()); + logger().debug("# face = {}", MC.bsp_faces.size()); + logger().debug("# edge = {}", MC.bsp_edges.size()); + logger().debug("# vertex = {}", MC.bsp_vertices.size()); + logger().info("BSP subdivision done!"); + tmp_time = igl_timer.getElapsedTime(); + addRecord(MeshRecord(MeshRecord::OpType::OP_BSP, tmp_time, bsp_vertices.size(), bsp_nodes.size())); + sum_time += tmp_time; + logger().info("time = {}s", tmp_time); + + //simple tetrahedralization + igl_timer.start(); + logger().info("Tetrehedralizing ..."); + SimpleTetrahedralization ST(MC); + ST.tetra(MR.tet_vertices, MR.tets); + ST.labelSurface(m_f_tags, raw_e_tags, raw_conn_e4v, MR.tet_vertices, MR.tets, MR.is_surface_fs); + ST.labelBbox(MR.tet_vertices, MR.tets); + if (!State::state().g_is_close)//if input is an open mesh + ST.labelBoundary(MR.tet_vertices, MR.tets, MR.is_surface_fs); + logger().debug("# tet_vertices = {}", MR.tet_vertices.size()); + logger().debug("# tets = {}", MR.tets.size()); + logger().info("Tetrahedralization done!"); + tmp_time = igl_timer.getElapsedTime(); + addRecord(MeshRecord(MeshRecord::OpType::OP_SIMPLE_TETRA, tmp_time, MR.tet_vertices.size(), MR.tets.size())); + sum_time += tmp_time; + logger().info("time = {}s", tmp_time); - /// STAGE 2 - //init -#ifndef MUTE_COUT - cout << "Refinement initializing..." << endl; -#endif - MR.prepareData(); -#ifndef MUTE_COUT - cout << "Refinement intialization done!" << endl; - cout << endl; -#endif - - //improvement - MR.refine(energy_type); - - outputFinalTetmesh(MR, out_vertices, out_tets); //do winding number and output the tetmesh + logger().info("Total time for the first stage = {}", sum_time); } - void tetrahedralization(const std::vector>& in_vertices, - const std::vector>& in_faces, - std::vector>& out_vertices, - std::vector>& out_tets) { - out_vertices.clear(); - out_tets.clear(); - - args.i_epsilon = parameters.i_epsilon; - args.bg_mesh = parameters.bg_mesh; - args.filter_energy = parameters.filter_energy; - args.i_ideal_edge_length = parameters.i_ideal_edge_length; - args.is_laplacian = parameters.is_laplacian; - args.is_quiet = parameters.is_quiet; - args.max_pass = parameters.max_pass; - args.stage = parameters.stage; - args.targeted_num_v = parameters.targeted_num_v; - - //initalization - GEO::initialize(); - g_postfix = args.postfix; - g_working_dir = args.input.substr(0, args.input.size() - 4); - - if (args.csv_file == "") - g_stat_file = g_working_dir + g_postfix + ".csv"; - else - g_stat_file = args.csv_file; - - if (args.output == "") - g_output_file = g_working_dir + g_postfix + ".msh"; - else - g_output_file = args.output; - - if (args.is_quiet) { - args.is_output_csv = false; - cout.setstate(std::ios_base::failbit);//use std::cout.clear() to get it back - } + /// STAGE 2 + //init + logger().info("Refinement initializing..."); + MR.prepareData(); + logger().info("Refinement intialization done!"); - //do tetrahedralization - Eigen::MatrixXd V; - Eigen::MatrixXi F; - V.resize(in_vertices.size(), 3); - F.resize(in_faces.size(), 3); - for (int i = 0; i < in_vertices.size(); i++) - for (int j = 0; j < 3; j++) - V(i, j) = in_vertices[i][j]; - for (int i = 0; i < in_faces.size(); i++) - for (int j = 0; j < 3; j++) - F(i, j) = in_faces[i][j]; - gtet_new(V, F, out_vertices, out_tets); - - if (args.is_quiet) { - std::cout.clear(); - } + //improvement + MR.refine(energy_type); + + outputFinalTetmesh(MR, out_vertices, out_tets); //do winding number and output the tetmesh +} + +void tetrahedralization(const std::vector>& in_vertices, + const std::vector>& in_faces, + std::vector>& out_vertices, + std::vector>& out_tets, + Args parameters) +{ + out_vertices.clear(); + out_tets.clear(); + + GArgs::args().i_epsilon = parameters.i_epsilon; + GArgs::args().bg_mesh = parameters.bg_mesh; + GArgs::args().filter_energy = parameters.filter_energy; + GArgs::args().i_ideal_edge_length = parameters.i_ideal_edge_length; + GArgs::args().is_laplacian = parameters.is_laplacian; + GArgs::args().is_quiet = parameters.is_quiet; + GArgs::args().max_pass = parameters.max_pass; + GArgs::args().stage = parameters.stage; + GArgs::args().targeted_num_v = parameters.targeted_num_v; + + //initalization + GEO::initialize(); + State::state().g_postfix = GArgs::args().postfix; + State::state().g_working_dir = GArgs::args().input.substr(0, GArgs::args().input.size() - 4); + + if (GArgs::args().csv_file == "") + State::state().g_stat_file = State::state().g_working_dir + State::state().g_postfix + ".csv"; + else + State::state().g_stat_file = GArgs::args().csv_file; + + if (GArgs::args().output == "") + State::state().g_output_file = State::state().g_working_dir + State::state().g_postfix + ".msh"; + else + State::state().g_output_file = GArgs::args().output; + + if (GArgs::args().is_quiet) { + GArgs::args().is_output_csv = false; + std::cout.setstate(std::ios_base::failbit);//use std::cout.clear() to get it back + } + + //do tetrahedralization + Eigen::MatrixXd V; + Eigen::MatrixXi F; + V.resize(in_vertices.size(), 3); + F.resize(in_faces.size(), 3); + for (int i = 0; i < in_vertices.size(); i++) + for (int j = 0; j < 3; j++) + V(i, j) = in_vertices[i][j]; + for (int i = 0; i < in_faces.size(); i++) + for (int j = 0; j < 3; j++) + F(i, j) = in_faces[i][j]; + gtet_new(V, F, out_vertices, out_tets); + + if (GArgs::args().is_quiet) { + std::cout.clear(); } } + +} // namespace tetwild +