mirror of
https://github.com/lucaspalomodevelop/meshlab.git
synced 2026-03-13 08:09:39 +00:00
candidate final
- removed heatgeodesic filter - merged heatgeodesic filter into geodesic filter - moved algorithmic part to vcglib
This commit is contained in:
davidboening
parent
f14dca3693
commit
3faad7bd73
@ -152,8 +152,7 @@ if(NOT DEFINED MESHLAB_PLUGINS) # it may be already defined in parent directory
|
||||
meshlabplugins/filter_developability
|
||||
meshlabplugins/filter_dirt
|
||||
meshlabplugins/filter_fractal
|
||||
meshlabplugins/filter_func
|
||||
meshlabplugins/filter_heatgeodesic # added
|
||||
meshlabplugins/filter_func
|
||||
meshlabplugins/filter_img_patch_param
|
||||
meshlabplugins/filter_icp
|
||||
meshlabplugins/filter_io_nxs
|
||||
|
||||
@ -44,7 +44,8 @@ FilterGeodesic::FilterGeodesic()
|
||||
typeList = {
|
||||
FP_QUALITY_BORDER_GEODESIC,
|
||||
FP_QUALITY_POINT_GEODESIC,
|
||||
FP_QUALITY_SELECTED_GEODESIC
|
||||
FP_QUALITY_SELECTED_GEODESIC,
|
||||
FP_QUALITY_SELECTED_GEODESIC_HEAT
|
||||
};
|
||||
|
||||
for(ActionIDType tt : types())
|
||||
@ -67,6 +68,8 @@ QString FilterGeodesic::filterName(ActionIDType filter) const
|
||||
return QString("Colorize by geodesic distance from a given point");
|
||||
case FP_QUALITY_SELECTED_GEODESIC:
|
||||
return QString("Colorize by geodesic distance from the selected points");
|
||||
case FP_QUALITY_SELECTED_GEODESIC_HEAT:
|
||||
return QString("Colorize by approximated geodesic distance from the selected points");
|
||||
default: assert(0); return QString();
|
||||
}
|
||||
}
|
||||
@ -79,6 +82,8 @@ QString FilterGeodesic::pythonFilterName(ActionIDType f) const
|
||||
return QString("compute_scalar_by_geodesic_distance_from_given_point_per_vertex");
|
||||
case FP_QUALITY_SELECTED_GEODESIC:
|
||||
return QString("compute_scalar_by_geodesic_distance_from_selection_per_vertex");
|
||||
case FP_QUALITY_SELECTED_GEODESIC_HEAT:
|
||||
return QString("compute_scalar_by_heat_geodesic_distance_from_selection_per_vertex");
|
||||
default: assert(0); return QString();
|
||||
}
|
||||
}
|
||||
@ -87,9 +92,10 @@ QString FilterGeodesic::filterInfo(ActionIDType filterId) const
|
||||
{
|
||||
switch(filterId)
|
||||
{
|
||||
case FP_QUALITY_BORDER_GEODESIC : return tr("Store in the quality field the geodesic distance from borders and color the mesh accordingly.");
|
||||
case FP_QUALITY_POINT_GEODESIC : return tr("Store in the quality field the geodesic distance from a given point on the mesh surface and color the mesh accordingly.");
|
||||
case FP_QUALITY_SELECTED_GEODESIC : return tr("Store in the quality field the geodesic distance from the selected points on the mesh surface and color the mesh accordingly.");
|
||||
case FP_QUALITY_BORDER_GEODESIC : return tr("Store in the quality field the geodesic distance from borders and color the mesh accordingly.");
|
||||
case FP_QUALITY_POINT_GEODESIC : return tr("Store in the quality field the geodesic distance from a given point on the mesh surface and color the mesh accordingly.");
|
||||
case FP_QUALITY_SELECTED_GEODESIC : return tr("Store in the quality field the geodesic distance from the selected points on the mesh surface and color the mesh accordingly.");
|
||||
case FP_QUALITY_SELECTED_GEODESIC_HEAT : return tr("Store in the quality field the approximated geodesic distance, computed via heat method, from the selected points on the mesh surface and color the mesh accordingly. This method is very sensitive to triangulation and currently does not support disconnected components.");
|
||||
default : assert(0);
|
||||
}
|
||||
return QString("error!");
|
||||
@ -99,9 +105,10 @@ FilterGeodesic::FilterClass FilterGeodesic::getClass(const QAction *a) const
|
||||
{
|
||||
switch(ID(a))
|
||||
{
|
||||
case FP_QUALITY_BORDER_GEODESIC :
|
||||
case FP_QUALITY_SELECTED_GEODESIC :
|
||||
case FP_QUALITY_POINT_GEODESIC : return FilterGeodesic::FilterClass(FilterPlugin::VertexColoring + FilterPlugin::Quality);
|
||||
case FP_QUALITY_BORDER_GEODESIC :
|
||||
case FP_QUALITY_SELECTED_GEODESIC :
|
||||
case FP_QUALITY_POINT_GEODESIC : return FilterGeodesic::FilterClass(FilterPlugin::VertexColoring + FilterPlugin::Quality);
|
||||
case FP_QUALITY_SELECTED_GEODESIC_HEAT : return FilterGeodesic::FilterClass(FilterPlugin::VertexColoring + FilterPlugin::Quality);
|
||||
default : assert(0);
|
||||
}
|
||||
return FilterPlugin::Generic;
|
||||
@ -111,9 +118,10 @@ int FilterGeodesic::getRequirements(const QAction *action)
|
||||
{
|
||||
switch(ID(action))
|
||||
{
|
||||
case FP_QUALITY_BORDER_GEODESIC :
|
||||
case FP_QUALITY_SELECTED_GEODESIC:
|
||||
case FP_QUALITY_POINT_GEODESIC : return MeshModel::MM_VERTFACETOPO;
|
||||
case FP_QUALITY_BORDER_GEODESIC :
|
||||
case FP_QUALITY_SELECTED_GEODESIC :
|
||||
case FP_QUALITY_POINT_GEODESIC : return MeshModel::MM_VERTFACETOPO;
|
||||
case FP_QUALITY_SELECTED_GEODESIC_HEAT : return MeshModel::MM_VERTFACETOPO + MeshModel::MM_FACEFACETOPO;
|
||||
default: assert(0);
|
||||
}
|
||||
return 0;
|
||||
@ -236,7 +244,56 @@ std::map<std::string, QVariant> FilterGeodesic::applyFilter(const QAction *filte
|
||||
else
|
||||
log("Warning: no vertices are selected! aborting geodesic computation.");
|
||||
}
|
||||
break;
|
||||
break;
|
||||
case FP_QUALITY_SELECTED_GEODESIC_HEAT:
|
||||
{
|
||||
m.updateDataMask(MeshModel::MM_FACEFACETOPO);
|
||||
m.updateDataMask(MeshModel::MM_VERTFACETOPO);
|
||||
m.updateDataMask(MeshModel::MM_FACEQUALITY);
|
||||
m.updateDataMask(MeshModel::MM_FACENORMAL);
|
||||
m.updateDataMask(MeshModel::MM_VERTQUALITY);
|
||||
m.updateDataMask(MeshModel::MM_VERTCOLOR);
|
||||
|
||||
std::vector<CMeshO::VertexPointer> seedVec;
|
||||
ForEachVertex(m.cm, [&seedVec] (CMeshO::VertexType & v) {
|
||||
if (v.IsS()) seedVec.push_back(&v);
|
||||
});
|
||||
|
||||
float param_m = par.getFloat("m");
|
||||
|
||||
if (seedVec.size() > 0){
|
||||
// cache factorizations to reduce runtime on successive computations
|
||||
std::pair<float, tri::GeodesicHeatCache> cache;
|
||||
if (!vcg::tri::HasPerMeshAttribute(m.cm, "GeodesicHeatCache")){
|
||||
cache = std::make_pair(param_m, tri::GeodesicHeat<CMeshO>::BuildCache(m.cm, param_m));
|
||||
// save cache for next compute
|
||||
auto GeodesicHeatCacheHandle = tri::Allocator<CMeshO>::GetPerMeshAttribute<std::pair<float, tri::GeodesicHeatCache>>(m.cm, std::string("GeodesicHeatCache"));
|
||||
GeodesicHeatCacheHandle() = cache;
|
||||
}
|
||||
else {
|
||||
// recover cache
|
||||
auto GeodesicHeatCacheHandle = vcg::tri::Allocator<CMeshO>::GetPerMeshAttribute<std::pair<float, tri::GeodesicHeatCache>>(m.cm, std::string("GeodesicHeatCache"));
|
||||
// if m has changed rebuild everything
|
||||
if (std::get<0>(GeodesicHeatCacheHandle()) != param_m){
|
||||
GeodesicHeatCacheHandle() = std::make_pair(param_m, tri::GeodesicHeat<CMeshO>::BuildCache(m.cm, param_m));
|
||||
}
|
||||
cache = GeodesicHeatCacheHandle();
|
||||
}
|
||||
|
||||
if (tri::GeodesicHeat<CMeshO>::ComputeFromCache(m.cm, seedVec, std::get<1>(cache))){
|
||||
tri::UpdateColor<CMeshO>::PerVertexQualityRamp(m.cm);
|
||||
}
|
||||
else{
|
||||
log("Warning: heat method has failed. The mesh is most likely badly conditioned (e.g. angles ~ 0deg) or has disconnected components");
|
||||
// delete cache as its most likely useless after failure
|
||||
auto GeodesicHeatCacheHandle = vcg::tri::Allocator<CMeshO>::GetPerMeshAttribute<std::pair<float, tri::GeodesicHeatCache>>(m.cm, std::string("GeodesicHeatCache"));
|
||||
tri::Allocator<CMeshO>::DeletePerMeshAttribute<std::pair<float, tri::GeodesicHeatCache>>(m.cm, GeodesicHeatCacheHandle);
|
||||
}
|
||||
}
|
||||
else
|
||||
log("Warning: no vertices are selected! aborting geodesic computation.");
|
||||
}
|
||||
break;
|
||||
default:
|
||||
wrongActionCalled(filter);
|
||||
break;
|
||||
@ -255,7 +312,10 @@ RichParameterList FilterGeodesic::initParameterList(const QAction *action, const
|
||||
break;
|
||||
case FP_QUALITY_SELECTED_GEODESIC :
|
||||
parlst.addParam(RichPercentage("maxDistance",m.cm.bbox.Diag(),0,m.cm.bbox.Diag()*2,"Max Distance","If not zero it indicates a cut off value to be used during geodesic distance computation."));
|
||||
break;
|
||||
break;
|
||||
case FP_QUALITY_SELECTED_GEODESIC_HEAT :
|
||||
parlst.addParam(RichFloat("m", 1, "m", "Multiplier used in backward Euler timestep."));
|
||||
break;
|
||||
default: break; // do not add any parameter for the other filters
|
||||
}
|
||||
return parlst;
|
||||
@ -265,10 +325,11 @@ int FilterGeodesic::postCondition(const QAction * filter) const
|
||||
{
|
||||
switch (ID(filter))
|
||||
{
|
||||
case FP_QUALITY_BORDER_GEODESIC :
|
||||
case FP_QUALITY_SELECTED_GEODESIC :
|
||||
case FP_QUALITY_POINT_GEODESIC : return MeshModel::MM_VERTCOLOR + MeshModel::MM_VERTQUALITY;
|
||||
default : return MeshModel::MM_ALL;
|
||||
case FP_QUALITY_BORDER_GEODESIC :
|
||||
case FP_QUALITY_SELECTED_GEODESIC :
|
||||
case FP_QUALITY_POINT_GEODESIC : return MeshModel::MM_VERTCOLOR + MeshModel::MM_VERTQUALITY;
|
||||
case FP_QUALITY_SELECTED_GEODESIC_HEAT : return MeshModel::MM_VERTCOLOR + MeshModel::MM_VERTQUALITY;
|
||||
default : return MeshModel::MM_ALL;
|
||||
}
|
||||
}
|
||||
MESHLAB_PLUGIN_NAME_EXPORTER(FilterGeodesic)
|
||||
|
||||
@ -26,6 +26,7 @@
|
||||
#include <QObject>
|
||||
#include <common/plugins/interfaces/filter_plugin.h>
|
||||
#include <vcg/complex/algorithms/geodesic.h>
|
||||
#include <vcg/complex/algorithms/geodesic_heat.h>
|
||||
|
||||
|
||||
class FilterGeodesic : public QObject, public FilterPlugin
|
||||
@ -42,7 +43,8 @@ class FilterGeodesic : public QObject, public FilterPlugin
|
||||
enum {
|
||||
FP_QUALITY_BORDER_GEODESIC,
|
||||
FP_QUALITY_POINT_GEODESIC,
|
||||
FP_QUALITY_SELECTED_GEODESIC
|
||||
FP_QUALITY_SELECTED_GEODESIC,
|
||||
FP_QUALITY_SELECTED_GEODESIC_HEAT
|
||||
|
||||
} ;
|
||||
|
||||
|
||||
@ -1,5 +0,0 @@
|
||||
set(SOURCES filter_heatgeodesic.cpp)
|
||||
|
||||
set(HEADERS filter_heatgeodesic.h trimesh_heatmethod.h)
|
||||
|
||||
add_meshlab_plugin(filter_heatgeodesic ${SOURCES} ${HEADERS})
|
||||
@ -1,249 +0,0 @@
|
||||
#include "trimesh_heatmethod.h"
|
||||
|
||||
#include "filter_heatgeodesic.h"
|
||||
|
||||
#include <common/plugins/interfaces/meshlab_plugin_logger.h>
|
||||
|
||||
/**
|
||||
* @brief Constructor usually performs only two simple tasks of filling the two lists
|
||||
* - typeList: with all the possible id of the filtering actions
|
||||
* - actionList with the corresponding actions. If you want to add icons to
|
||||
* your filtering actions you can do here by construction the QActions accordingly
|
||||
*/
|
||||
FilterHeatGeodesicPlugin::FilterHeatGeodesicPlugin()
|
||||
{
|
||||
typeList = {FP_COMPUTE_HEATGEODESIC_FROM_SELECTION};
|
||||
|
||||
for(ActionIDType tt : types())
|
||||
actionList.push_back(new QAction(filterName(tt), this));
|
||||
}
|
||||
|
||||
FilterHeatGeodesicPlugin::~FilterHeatGeodesicPlugin()
|
||||
{
|
||||
}
|
||||
|
||||
QString FilterHeatGeodesicPlugin::pluginName() const
|
||||
{
|
||||
return "FilterHeatGeodesic";
|
||||
}
|
||||
|
||||
/**
|
||||
* @brief ST() must return the very short string describing each filtering action
|
||||
* (this string is used also to define the menu entry)
|
||||
* @param filterId: the id of the filter
|
||||
* @return the name of the filter
|
||||
*/
|
||||
QString FilterHeatGeodesicPlugin::filterName(ActionIDType filterId) const
|
||||
{
|
||||
switch(filterId) {
|
||||
case FP_COMPUTE_HEATGEODESIC_FROM_SELECTION :
|
||||
return QString("Compute Geodesic From Selection (Heat)");
|
||||
default :
|
||||
assert(0);
|
||||
return QString();
|
||||
}
|
||||
}
|
||||
|
||||
/**
|
||||
* @brief FilterHeatGeodesicPlugin::pythonFilterName if you want that your filter should have a different
|
||||
* name on pymeshlab, use this function to return its python name.
|
||||
* @param f
|
||||
* @return
|
||||
*/
|
||||
QString FilterHeatGeodesicPlugin::pythonFilterName(ActionIDType f) const
|
||||
{
|
||||
switch(f) {
|
||||
case FP_COMPUTE_HEATGEODESIC_FROM_SELECTION :
|
||||
return "compute_approximate_geodesic_from_selection";
|
||||
default :
|
||||
assert(0);
|
||||
return QString();
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
/**
|
||||
* @brief // Info() must return the longer string describing each filtering action
|
||||
* (this string is used in the About plugin dialog)
|
||||
* @param filterId: the id of the filter
|
||||
* @return an info string of the filter
|
||||
*/
|
||||
QString FilterHeatGeodesicPlugin::filterInfo(ActionIDType filterId) const
|
||||
{
|
||||
switch(filterId) {
|
||||
case FP_COMPUTE_HEATGEODESIC_FROM_SELECTION :
|
||||
return QString(
|
||||
"Computes an approximated geodesic distance from the selected vertices to all others. "
|
||||
"This algorithm implements the heat method. "
|
||||
"If the mesh has changed, update m to a different value to reflect these changes as"
|
||||
" this rebuilds the existing cache."
|
||||
);
|
||||
default :
|
||||
assert(0);
|
||||
return QString("Unknown Filter");
|
||||
}
|
||||
}
|
||||
|
||||
/**
|
||||
* @brief The FilterClass describes in which generic class of filters it fits.
|
||||
* This choice affect the submenu in which each filter will be placed
|
||||
* More than a single class can be chosen.
|
||||
* @param a: the action of the filter
|
||||
* @return the class od the filter
|
||||
*/
|
||||
FilterHeatGeodesicPlugin::FilterClass FilterHeatGeodesicPlugin::getClass(const QAction *a) const
|
||||
{
|
||||
switch(ID(a)) {
|
||||
case FP_COMPUTE_HEATGEODESIC_FROM_SELECTION :
|
||||
return FilterPlugin::Measure; // Note: not sure
|
||||
// return FilterGeodesic::FilterClass(FilterPlugin::VertexColoring + FilterPlugin::Quality);
|
||||
default :
|
||||
assert(0);
|
||||
return FilterPlugin::Generic;
|
||||
}
|
||||
}
|
||||
|
||||
/**
|
||||
* @brief FilterSamplePlugin::filterArity
|
||||
* @return
|
||||
*/
|
||||
FilterPlugin::FilterArity FilterHeatGeodesicPlugin::filterArity(const QAction*) const
|
||||
{
|
||||
return SINGLE_MESH;
|
||||
}
|
||||
|
||||
/**
|
||||
* @brief FilterSamplePlugin::getPreConditions
|
||||
* @return
|
||||
*/
|
||||
int FilterHeatGeodesicPlugin::getPreConditions(const QAction*) const
|
||||
{
|
||||
// NOTE: note sure about these
|
||||
return MeshModel::MM_FACEFACETOPO | MeshModel::MM_VERTFACETOPO;
|
||||
}
|
||||
|
||||
/**
|
||||
* @brief FilterSamplePlugin::postCondition
|
||||
* @return
|
||||
*/
|
||||
int FilterHeatGeodesicPlugin::postCondition(const QAction*) const
|
||||
{
|
||||
// NOTE: note sure about these
|
||||
return MeshModel::MM_VERTQUALITY | MeshModel::MM_VERTCOLOR;
|
||||
}
|
||||
|
||||
/**
|
||||
* @brief This function define the needed parameters for each filter. Return true if the filter has some parameters
|
||||
* it is called every time, so you can set the default value of parameters according to the mesh
|
||||
* For each parameter you need to define,
|
||||
* - the name of the parameter,
|
||||
* - the default value
|
||||
* - the string shown in the dialog
|
||||
* - a possibly long string describing the meaning of that parameter (shown as a popup help in the dialog)
|
||||
* @param action
|
||||
* @param m
|
||||
* @param parlst
|
||||
*/
|
||||
RichParameterList FilterHeatGeodesicPlugin::initParameterList(const QAction *action,const MeshModel &m)
|
||||
{
|
||||
RichParameterList parlst;
|
||||
switch(ID(action)) {
|
||||
case FP_COMPUTE_HEATGEODESIC_FROM_SELECTION :
|
||||
parlst.addParam(RichFloat("m", 1, "m", "Multiplier used in backward Euler timestep."));
|
||||
break;
|
||||
default :
|
||||
assert(0);
|
||||
}
|
||||
return parlst;
|
||||
}
|
||||
|
||||
/**
|
||||
* @brief The Real Core Function doing the actual mesh processing.
|
||||
* @param action
|
||||
* @param md: an object containing all the meshes and rasters of MeshLab
|
||||
* @param par: the set of parameters of each filter
|
||||
* @param cb: callback object to tell MeshLab the percentage of execution of the filter
|
||||
* @return true if the filter has been applied correctly, false otherwise
|
||||
*/
|
||||
std::map<std::string, QVariant> FilterHeatGeodesicPlugin::applyFilter(const QAction * action, const RichParameterList & parameters, MeshDocument &md, unsigned int& /*postConditionMask*/, vcg::CallBackPos *cb)
|
||||
{
|
||||
MeshModel &mm=*(md.mm());
|
||||
switch(ID(action)) {
|
||||
case FP_COMPUTE_HEATGEODESIC_FROM_SELECTION :
|
||||
{
|
||||
mm.updateDataMask(MeshModel::MM_FACEFACETOPO);
|
||||
mm.updateDataMask(MeshModel::MM_VERTFACETOPO);
|
||||
mm.updateDataMask(MeshModel::MM_FACEQUALITY);
|
||||
mm.updateDataMask(MeshModel::MM_FACENORMAL);
|
||||
mm.updateDataMask(MeshModel::MM_VERTQUALITY);
|
||||
mm.updateDataMask(MeshModel::MM_VERTCOLOR);
|
||||
CMeshO &mesh = mm.cm;
|
||||
// TODO: compact all vectors
|
||||
computeHeatGeodesicFromSelection(mesh, cb, parameters.getFloat("m"));
|
||||
}
|
||||
break;
|
||||
default :
|
||||
wrongActionCalled(action);
|
||||
}
|
||||
return std::map<std::string, QVariant>();
|
||||
}
|
||||
|
||||
inline void FilterHeatGeodesicPlugin::computeHeatGeodesicFromSelection(CMeshO& mesh, vcg::CallBackPos* cb, float m){
|
||||
vcg::tri::Allocator<CMeshO>::CompactEveryVector(mesh);
|
||||
|
||||
cb(1, "Checking Selection...");
|
||||
// build and check source vertices
|
||||
Eigen::VectorXd sourcePoints(mesh.VN());
|
||||
int selection_count = 0;
|
||||
for(int i=0; i < mesh.VN(); i++){
|
||||
sourcePoints(i) = mesh.vert[i].IsS() ? (++selection_count, 1) : 0;
|
||||
}
|
||||
if (selection_count < 1){
|
||||
log("Warning: no vertices are selected! aborting computation.");
|
||||
return;
|
||||
}
|
||||
|
||||
// TODO detect if the mesh has changed and if so delete state
|
||||
// NOW changing the value of m rebuilds the cache
|
||||
//
|
||||
// if (meshHasChanged) {
|
||||
// auto handle = vcg::tri::Allocator<CMeshO>::GetPerMeshAttribute<std::shared_ptr<HeatMethodSolver>>(mesh, std::string("HeatMethodSolver"));
|
||||
// vcg::tri::Allocator<CMeshO>::DeletePerMeshAttribute<std::shared_ptr<HeatMethodSolver>>(mesh, handle);
|
||||
// }
|
||||
|
||||
std::shared_ptr<HeatMethodSolver> HMSolver;
|
||||
// if no state can be recovered, initialize it
|
||||
if (!vcg::tri::HasPerMeshAttribute(mesh, "HeatMethodSolver")){
|
||||
HMSolver = std::make_shared<HeatMethodSolver>(HeatMethodSolver(mesh, m, cb));
|
||||
if (HMSolver->factorization_failed())
|
||||
log("Warning: factorization has failed. The mesh is non-manifold or badly conditioned (e.g. angles ~ 0deg) or disconnected components");
|
||||
auto HMSolverHandle = vcg::tri::Allocator<CMeshO>::GetPerMeshAttribute<std::shared_ptr<HeatMethodSolver>>(mesh, std::string("HeatMethodSolver"));
|
||||
HMSolverHandle() = HMSolver;
|
||||
}
|
||||
else {
|
||||
cb(5, "Attempting to Recover State...");
|
||||
// recover solver
|
||||
auto HMSolverHandle = vcg::tri::Allocator<CMeshO>::GetPerMeshAttribute<std::shared_ptr<HeatMethodSolver>>(mesh, std::string("HeatMethodSolver"));
|
||||
// if m has changed rebuild everything
|
||||
if (HMSolverHandle()->get_m() != m){
|
||||
HMSolverHandle() = std::make_shared<HeatMethodSolver>(HeatMethodSolver(mesh, m, cb));
|
||||
if (HMSolverHandle()->factorization_failed())
|
||||
log("Warning: factorization has failed. The changed mesh is non-manifold or badly conditioned (e.g. angles ~ 0deg) or disconnected components");
|
||||
}
|
||||
HMSolver = HMSolverHandle();
|
||||
}
|
||||
|
||||
cb(85, "Computing Geodesic Distance...");
|
||||
Eigen::VectorXd geodesicDistance = HMSolver->solve(mesh, sourcePoints);
|
||||
if (HMSolver->solver_failed())
|
||||
log("Warning: solver has failed.");
|
||||
|
||||
cb(99, "Updating Mesh...");
|
||||
// set geodesic distance as quality and color
|
||||
for(int i=0; i < mesh.VN(); i++){
|
||||
mesh.vert[i].Q() = geodesicDistance(i);
|
||||
}
|
||||
vcg::tri::UpdateColor<CMeshO>::PerVertexQualityRamp(mesh);
|
||||
}
|
||||
|
||||
MESHLAB_PLUGIN_NAME_EXPORTER(FilterSamplePlugin)
|
||||
@ -1,39 +0,0 @@
|
||||
#ifndef FILTERSAMPLE_PLUGIN_H
|
||||
#define FILTERSAMPLE_PLUGIN_H
|
||||
|
||||
#include <common/plugins/interfaces/filter_plugin.h>
|
||||
|
||||
class FilterHeatGeodesicPlugin : public QObject, public FilterPlugin
|
||||
{
|
||||
Q_OBJECT
|
||||
MESHLAB_PLUGIN_IID_EXPORTER(FILTER_PLUGIN_IID)
|
||||
Q_INTERFACES(FilterPlugin)
|
||||
|
||||
public:
|
||||
enum { FP_COMPUTE_HEATGEODESIC_FROM_SELECTION } ;
|
||||
|
||||
FilterHeatGeodesicPlugin();
|
||||
virtual ~FilterHeatGeodesicPlugin();
|
||||
|
||||
QString pluginName() const;
|
||||
|
||||
QString filterName(ActionIDType filter) const;
|
||||
QString pythonFilterName(ActionIDType f) const;
|
||||
QString filterInfo(ActionIDType filter) const;
|
||||
FilterClass getClass(const QAction* a) const;
|
||||
FilterArity filterArity(const QAction*) const;
|
||||
int getPreConditions(const QAction *) const;
|
||||
int postCondition(const QAction* ) const;
|
||||
RichParameterList initParameterList(const QAction*, const MeshModel &/*m*/);
|
||||
std::map<std::string, QVariant> applyFilter(
|
||||
const QAction* action,
|
||||
const RichParameterList & parameters,
|
||||
MeshDocument &md,
|
||||
unsigned int& postConditionMask,
|
||||
vcg::CallBackPos * cb);
|
||||
|
||||
private:
|
||||
void computeHeatGeodesicFromSelection(CMeshO& mesh, vcg::CallBackPos* cb, float m);
|
||||
};
|
||||
|
||||
#endif
|
||||
@ -1,329 +0,0 @@
|
||||
#ifndef TRIMESH_HEAT_METHOD
|
||||
#define TRIMESH_HEAT_METHOD
|
||||
|
||||
#include <Eigen/Sparse>
|
||||
#include <Eigen/SparseCholesky>
|
||||
|
||||
#include <vcg/complex/algorithms/update/quality.h>
|
||||
|
||||
#include <common/plugins/interfaces/filter_plugin.h>
|
||||
|
||||
#include <memory>
|
||||
|
||||
|
||||
// LOCAL FUNCTIONS
|
||||
|
||||
namespace {
|
||||
inline Eigen::Vector3d toEigen(const vcg::Point3f& p)
|
||||
{
|
||||
return Eigen::Vector3d(p.X(), p.Y(), p.Z());
|
||||
};
|
||||
|
||||
|
||||
inline double cotan(const Eigen::Vector3d& v0, const Eigen::Vector3d& v1)
|
||||
{
|
||||
// cos(theta) / sin(theta)
|
||||
return v0.dot(v1) / v0.cross(v1).norm();
|
||||
};
|
||||
}
|
||||
|
||||
// simple struct to track progress bar status
|
||||
struct ProgressBarCallbackT {
|
||||
vcg::CallBackPos* callback;
|
||||
int current_progress;
|
||||
int allocated_progress;
|
||||
};
|
||||
|
||||
|
||||
inline void buildMassMatrix(CMeshO &mesh, Eigen::SparseMatrix<double> &mass){
|
||||
mass.resize(mesh.VN(), mesh.VN());
|
||||
mass.reserve(Eigen::VectorXi::Constant(mesh.VN(),1));
|
||||
|
||||
// compute area of all faces
|
||||
for (CMeshO::FaceIterator fi = mesh.face.begin(); fi != mesh.face.end(); ++fi)
|
||||
{
|
||||
vcg::Point3f p0 = fi->V(0)->P();
|
||||
vcg::Point3f p1 = fi->V(1)->P();
|
||||
vcg::Point3f p2 = fi->V(2)->P();
|
||||
double e0 = toEigen(p1 - p0).norm();
|
||||
double e1 = toEigen(p2 - p0).norm();
|
||||
double e2 = toEigen(p2 - p1).norm();
|
||||
double s = (e0 + e1 + e2) / 2;
|
||||
double area = std::sqrt(s * (s - e0) * (s - e1) * (s - e2));
|
||||
// we store face areas in quality field to avoid a hash table
|
||||
// this will also be useful for the gradient computation
|
||||
fi->Q() = area;
|
||||
}
|
||||
// compute area of the dual cell for each vertex
|
||||
for (int i = 0; i < mesh.VN(); ++i){
|
||||
CMeshO::VertexType *vp = &mesh.vert[i];
|
||||
|
||||
std::vector<CMeshO::FaceType*> faces;
|
||||
std::vector<int> indices;
|
||||
vcg::face::VFStarVF<CMeshO::FaceType>(vp, faces, indices);
|
||||
|
||||
double area = 0;
|
||||
for (int j = 0; j < faces.size(); ++j)
|
||||
{
|
||||
area += faces[j]->Q();
|
||||
}
|
||||
area /= 3;
|
||||
mass.insert(i, i) = area;
|
||||
}
|
||||
|
||||
mass.makeCompressed();
|
||||
}
|
||||
|
||||
|
||||
inline void buildCotanLowerTriMatrix(CMeshO &mesh, Eigen::SparseMatrix<double> &cotanOperator){
|
||||
cotanOperator.resize(mesh.VN(), mesh.VN());
|
||||
|
||||
typedef Eigen::Triplet<double> T;
|
||||
std::vector<T> tripletList;
|
||||
tripletList.reserve(3*mesh.VN() + 3*mesh.FN());
|
||||
|
||||
// compute cotan weights
|
||||
for (int i = 0; i < mesh.FN(); ++i){
|
||||
CMeshO::FaceType* fi = &mesh.face[i];
|
||||
|
||||
CMeshO::VertexType* v0 = fi->V(0);
|
||||
CMeshO::VertexType* v1 = fi->V(1);
|
||||
CMeshO::VertexType* v2 = fi->V(2);
|
||||
|
||||
vcg::Point3f p0 = v0->P();
|
||||
vcg::Point3f p1 = v1->P();
|
||||
vcg::Point3f p2 = v2->P();
|
||||
|
||||
Eigen::Vector3d e0 = toEigen(p2 - p1);
|
||||
Eigen::Vector3d e1 = toEigen(p0 - p2);
|
||||
Eigen::Vector3d e2 = toEigen(p1 - p0);
|
||||
|
||||
// first edge is inverted to get correct orientation
|
||||
double alpha0 = cotan(-e1, e2) / 2;
|
||||
double alpha1 = cotan(-e2, e0) / 2;
|
||||
double alpha2 = cotan(-e0, e1) / 2;
|
||||
|
||||
int i0 = vcg::tri::Index(mesh,v0);
|
||||
int i1 = vcg::tri::Index(mesh,v1);
|
||||
int i2 = vcg::tri::Index(mesh,v2);
|
||||
|
||||
// save only lower triangular part
|
||||
if (i0 > i1)
|
||||
tripletList.push_back(T(i0,i1,alpha2));
|
||||
else
|
||||
tripletList.push_back(T(i1,i0,alpha2));
|
||||
if (i0 > i2)
|
||||
tripletList.push_back(T(i0,i2,alpha1));
|
||||
else
|
||||
tripletList.push_back(T(i2,i0,alpha1));
|
||||
if (i1 > i2)
|
||||
tripletList.push_back(T(i1,i2,alpha0));
|
||||
else
|
||||
tripletList.push_back(T(i2,i1,alpha0));
|
||||
|
||||
tripletList.push_back(T(i0,i0,-(alpha1 + alpha2)));
|
||||
tripletList.push_back(T(i1,i1,-(alpha0 + alpha2)));
|
||||
tripletList.push_back(T(i2,i2,-(alpha0 + alpha1)));
|
||||
}
|
||||
|
||||
cotanOperator.setFromTriplets(tripletList.begin(), tripletList.end());
|
||||
cotanOperator.makeCompressed();
|
||||
}
|
||||
|
||||
|
||||
inline double computeAverageEdgeLength(CMeshO &mesh){
|
||||
// compute total length of all edges
|
||||
double total_length = 0;
|
||||
for (CMeshO::FaceIterator fi = mesh.face.begin(); fi != mesh.face.end(); ++fi)
|
||||
{
|
||||
vcg::Point3f p0 = fi->V(0)->P();
|
||||
vcg::Point3f p1 = fi->V(1)->P();
|
||||
vcg::Point3f p2 = fi->V(2)->P();
|
||||
double e2 = toEigen(p1 - p0).norm();
|
||||
double e1 = toEigen(p2 - p0).norm();
|
||||
double e0 = toEigen(p2 - p1).norm();
|
||||
total_length += (e0 + e1 + e2) / 2;
|
||||
}
|
||||
return total_length / (3./2. * mesh.FN());
|
||||
}
|
||||
|
||||
|
||||
inline Eigen::MatrixX3d computeVertexGradient(CMeshO &mesh, const Eigen::VectorXd &heat){
|
||||
Eigen::MatrixX3d heatGradientField(mesh.FN(), 3);
|
||||
// compute gradient of heat function at each vertex
|
||||
for (int i = 0; i < mesh.FN(); ++i){
|
||||
CMeshO::FaceType *fp = &mesh.face[i];
|
||||
|
||||
CMeshO::VertexType* v0 = fp->V(0);
|
||||
CMeshO::VertexType* v1 = fp->V(1);
|
||||
CMeshO::VertexType* v2 = fp->V(2);
|
||||
|
||||
vcg::Point3f p0 = v0->P();
|
||||
vcg::Point3f p1 = v1->P();
|
||||
vcg::Point3f p2 = v2->P();
|
||||
|
||||
int i0 = vcg::tri::Index(mesh,v0);
|
||||
int i1 = vcg::tri::Index(mesh,v1);
|
||||
int i2 = vcg::tri::Index(mesh,v2);
|
||||
|
||||
// normal unit vector
|
||||
Eigen::Vector3d n = toEigen(fp->N());
|
||||
n /= n.norm();
|
||||
// face area
|
||||
double faceArea = fp->Q();
|
||||
// edge unit vectors (counter-clockwise)
|
||||
Eigen::Vector3d e0 = toEigen(p2 - p1);
|
||||
e0 /= e0.norm();
|
||||
Eigen::Vector3d e1 = toEigen(p0 - p2);
|
||||
e1 /= e1.norm();
|
||||
Eigen::Vector3d e2 = toEigen(p1 - p0);
|
||||
e2 /= e2.norm();
|
||||
// gradient unit vectors
|
||||
Eigen::Vector3d g0 = n.cross(e0); //v0 grad
|
||||
Eigen::Vector3d g1 = n.cross(e1); //v1 grad
|
||||
Eigen::Vector3d g2 = n.cross(e2); //v2 grad
|
||||
|
||||
// add vertex gradient contributions
|
||||
Eigen::Vector3d tri_grad = (g0 * heat(i0) + g1 * heat(i1) + g2 * heat(i2)) / (2 * faceArea);
|
||||
|
||||
heatGradientField.row(i) = tri_grad;
|
||||
}
|
||||
return heatGradientField;
|
||||
}
|
||||
|
||||
|
||||
inline Eigen::MatrixX3d normalizeVectorField(const Eigen::MatrixX3d &field){
|
||||
Eigen::MatrixX3d normalizedField(field.rows(), 3);
|
||||
normalizedField.setZero();
|
||||
// normalize vector field at each vertex
|
||||
for (int i = 0; i < field.rows(); ++i){
|
||||
Eigen::Vector3d v = field.row(i);
|
||||
normalizedField.row(i) = v / v.norm();
|
||||
}
|
||||
return normalizedField;
|
||||
}
|
||||
|
||||
|
||||
inline Eigen::VectorXd computeVertexDivergence(CMeshO &mesh, const Eigen::MatrixX3d &field){
|
||||
Eigen::VectorXd divergence(mesh.VN());
|
||||
divergence.setZero();
|
||||
|
||||
// compute divergence of vector field at each vertex
|
||||
for (int i = 0; i < mesh.VN(); ++i){
|
||||
CMeshO::VertexType *vp = &mesh.vert[i];
|
||||
|
||||
std::vector<CMeshO::FaceType*> faces;
|
||||
std::vector<int> indices;
|
||||
vcg::face::VFStarVF<CMeshO::FaceType>(vp, faces, indices);
|
||||
for (int j = 0; j < faces.size(); ++j)
|
||||
{
|
||||
CMeshO::FaceType *fp = faces[j];
|
||||
int index = indices[j];
|
||||
|
||||
vcg::Point3f p0 = fp->V(0)->P();
|
||||
vcg::Point3f p1 = fp->V(1)->P();
|
||||
vcg::Point3f p2 = fp->V(2)->P();
|
||||
|
||||
// (ORDERING) edge vectors
|
||||
Eigen::Vector3d el, er, eo; //left, right, opposite to vp
|
||||
if (index == 0){
|
||||
el = toEigen(p2 - p0); //-e1
|
||||
er = toEigen(p1 - p0); //e2
|
||||
eo = toEigen(p2 - p1); //e0
|
||||
} else if (index == 1){
|
||||
el = toEigen(p0 - p1); //-e2
|
||||
er = toEigen(p2 - p1); //e0
|
||||
eo = toEigen(p0 - p2); //e1
|
||||
} else if (index == 2){
|
||||
el = toEigen(p1 - p2); //-e0
|
||||
er = toEigen(p0 - p2); //e1
|
||||
eo = toEigen(p1 - p0); //e2
|
||||
}
|
||||
// compute left and right cotangents
|
||||
double cotl = cotan(-el, -eo);
|
||||
double cotr = cotan(-er, eo);
|
||||
// normalize edge vectors after cotangent computation
|
||||
el /= el.norm();
|
||||
er /= er.norm();
|
||||
// add divergence contribution of given face
|
||||
Eigen::Vector3d x = field.row(vcg::tri::Index(mesh,fp));
|
||||
divergence(i) += (cotl * er.dot(x) + cotr * el.dot(x)) / 2;
|
||||
}
|
||||
}
|
||||
return divergence;
|
||||
}
|
||||
|
||||
|
||||
class HeatMethodSolver {
|
||||
public:
|
||||
HeatMethodSolver(CMeshO& mesh, float m, vcg::CallBackPos* cb = NULL)
|
||||
: m(m)
|
||||
{
|
||||
if(cb != NULL) cb(10, "Updating Topology...");
|
||||
// update topology and face normals
|
||||
vcg::tri::UpdateTopology<CMeshO>::VertexFace(mesh);
|
||||
vcg::tri::UpdateTopology<CMeshO>::FaceFace(mesh);
|
||||
vcg::tri::UpdateNormal<CMeshO>::PerFaceNormalized(mesh);
|
||||
|
||||
// state variables
|
||||
Eigen::SparseMatrix<double> massMatrix;
|
||||
Eigen::SparseMatrix<double> cotanMatrix;
|
||||
double averageEdgeLength;
|
||||
|
||||
buildMassMatrix(mesh, massMatrix);
|
||||
// this step is very slow (95% of total time) hence we pass the callback
|
||||
// to update progress and avoid freezes
|
||||
if(cb != NULL) cb(20, "Building Cotan Matrix...");
|
||||
buildCotanLowerTriMatrix(mesh, cotanMatrix);
|
||||
averageEdgeLength = computeAverageEdgeLength(mesh);
|
||||
|
||||
// initialize pointers
|
||||
solver1 = std::shared_ptr<Eigen::SimplicialLDLT<Eigen::SparseMatrix<double>>>(new Eigen::SimplicialLDLT<Eigen::SparseMatrix<double>>);
|
||||
solver2 = std::shared_ptr<Eigen::SimplicialLDLT<Eigen::SparseMatrix<double>>>(new Eigen::SimplicialLDLT<Eigen::SparseMatrix<double>>);
|
||||
|
||||
// build factorization
|
||||
double timestep = m * averageEdgeLength * averageEdgeLength;
|
||||
if(cb != NULL) cb(50, "Matrix Factorization...");
|
||||
solver1->compute(massMatrix - timestep * cotanMatrix);
|
||||
if(cb != NULL) cb(75, "Matrix Factorization...");
|
||||
solver2->compute(cotanMatrix);
|
||||
if((solver1->info() != Eigen::Success) || (solver2->info() != Eigen::Success))
|
||||
factorizationFailed = true;
|
||||
else
|
||||
factorizationFailed = false;
|
||||
}
|
||||
|
||||
HeatMethodSolver() = delete;
|
||||
|
||||
HeatMethodSolver& operator=(const HeatMethodSolver&) = delete;
|
||||
|
||||
Eigen::VectorXd solve(CMeshO &mesh, Eigen::VectorXd &sourcePoints){
|
||||
Eigen::VectorXd heatflow = solver1->solve(sourcePoints); // (VN)
|
||||
Eigen::MatrixX3d heatGradient = computeVertexGradient(mesh, heatflow); // (FN, 3)
|
||||
Eigen::MatrixX3d unitVectorField = normalizeVectorField(-heatGradient); // (FN, 3)
|
||||
Eigen::VectorXd divergence = computeVertexDivergence(mesh, unitVectorField); // (VN)
|
||||
Eigen::VectorXd geodesicDistance = solver2->solve(divergence); // (VN)
|
||||
if((solver1->info() != Eigen::Success) || (solver2->info() != Eigen::Success))
|
||||
solverFailed = true;
|
||||
else
|
||||
solverFailed = false;
|
||||
// shift to impose dist(source) = 0
|
||||
geodesicDistance.array() -= geodesicDistance.minCoeff();
|
||||
return geodesicDistance;
|
||||
}
|
||||
|
||||
double get_m(){return m;}
|
||||
|
||||
bool factorization_failed(){return factorizationFailed;}
|
||||
|
||||
bool solver_failed(){return solverFailed;}
|
||||
private:
|
||||
// Eigen::SimplicialLDLT has no copy/move constructor
|
||||
std::shared_ptr<Eigen::SimplicialLDLT<Eigen::SparseMatrix<double>>> solver1;
|
||||
std::shared_ptr<Eigen::SimplicialLDLT<Eigen::SparseMatrix<double>>> solver2;
|
||||
double m;
|
||||
bool factorizationFailed;
|
||||
bool solverFailed;
|
||||
};
|
||||
|
||||
#endif
|
||||
Loading…
x
Reference in New Issue
Block a user