Merge branch 'devel' into filter_developability

This commit is contained in:
lucaceschi 2022-09-20 21:52:01 +02:00 committed by GitHub
commit 2918e9e92c
No known key found for this signature in database
GPG Key ID: 4AEE18F83AFDEB23
24 changed files with 1855 additions and 1 deletions

View File

@ -162,6 +162,7 @@ if(NOT DEFINED MESHLAB_PLUGINS) # it may be already defined in parent directory
meshlabplugins/filter_color_projection
meshlabplugins/filter_colorproc
meshlabplugins/filter_create
meshlabplugins/filter_cubization
meshlabplugins/filter_developability
meshlabplugins/filter_dirt
meshlabplugins/filter_fractal

View File

@ -0,0 +1,34 @@
# Copyright 2019-2020, Collabora, Ltd.
# SPDX-License-Identifier: BSL-1.0
if (TARGET external-libigl AND TARGET Threads::Threads)
set(SOURCES
filter_cubization.cpp
CubicStylizationFiles/src/cube_style_precomputation.cpp
CubicStylizationFiles/src/cube_style_single_iteration.cpp
CubicStylizationFiles/src/fit_rotations_l1.cpp
CubicStylizationFiles/src/normalize_unitbox.cpp
CubicStylizationFiles/src/orthogonal_procrustes.cpp
CubicStylizationFiles/src/shrinkage.cpp)
set(HEADERS
filter_cubization.h
curvdata.h
curvedgeflip.h
CubicStylizationFiles/conversionMeshes.h
CubicStylizationFiles/cubic_stylizing.h
CubicStylizationFiles/src/cube_style_data.h
CubicStylizationFiles/src/cube_style_precomputation.h
CubicStylizationFiles/src/cube_style_single_iteration.h
CubicStylizationFiles/src/fit_rotations_l1.h
CubicStylizationFiles/src/normalize_unitbox.h
CubicStylizationFiles/src/orthogonal_procrustes.h
CubicStylizationFiles/src/shrinkage.h)
add_meshlab_plugin(filter_cubization ${SOURCES} ${HEADERS})
target_link_libraries(filter_cubization PRIVATE external-libigl Threads::Threads)
else()
message(
STATUS "Skipping filter_cubization - don't know about libigl or Threads::Threads on this system.")
endif()

View File

@ -0,0 +1,18 @@
# Cubic Stylization filter
Cubic stylization is a 3D stylization tool. Unlike image stylization (2D to 2D) and non-photorealistic rendering (3D to 2D), cubic stylization is a 3D to 3D stylization algorithm which takes a manifold triangle mesh as the input and outputs a cubified triangle mesh.
This is a MeshLab filter based on "[Cubic Stylization](https://www.dgp.toronto.edu/projects/cubic-stylization/)" by [Hsueh-Ti Derek Liu](https://www.dgp.toronto.edu/~hsuehtil/) and [Alec Jacobson](https://www.cs.toronto.edu/~jacobson/).
The folder `/src` contains the original C++ implementation of Cubic stylization performed by [Hsueh-Ti Derek Liu](https://www.dgp.toronto.edu/~hsuehtil/) and [Alec Jacobson](https://www.cs.toronto.edu/~jacobson/).
Their code is licensed under [MPL2](https://www.mozilla.org/en-US/MPL/2.0/).
## Bibtex
```
@article{Liu:CubicStyle:2019,
title = {Cubic Stylization},
author = {Hsueh-Ti Derek Liu and Alec Jacobson},
year = {2019},
journal = {ACM Transactions on Graphics},
}
```

View File

@ -0,0 +1,48 @@
#ifndef CONVERSIONMESHES_H
#define CONVERSIONMESHES_H
#include <vcg/complex/algorithms/mesh_to_matrix.h>
#include <vcg/complex/allocate.h>
#include <vector>
namespace vcg{
namespace tri{
template<class MeshType >
void Mesh2Matrix(MeshType& convertingMesh, Eigen::MatrixXd &output_vertexes, Eigen::MatrixXi &output_faces){
typedef typename MeshType::ScalarType ScalarType;
typedef typename Eigen::Matrix<ScalarType, Eigen::Dynamic, Eigen::Dynamic> MatrixXm;
MatrixXm temp_verts;
//convert mesh in matrix
vcg::tri::MeshToMatrix< MeshType >::GetTriMeshData( convertingMesh, temp_verts, output_faces);
output_vertexes = temp_verts.template cast<double>();
}
template<class MeshType >
void Matrix2Mesh(MeshType& convertedMesh, Eigen::MatrixXd vertexes, Eigen::MatrixXi faces){
typedef typename MeshType::VertexPointer VertexPointer;
typedef typename MeshType::FaceIterator FaceIterator;
//reconvert V and F matrixes in a mesh
convertedMesh.Clear();
Allocator<MeshType>::AddVertices(convertedMesh,vertexes.rows());
Allocator<MeshType>::AddFaces(convertedMesh,faces.rows());
std::vector<VertexPointer> ivp;
int i;
for (i=0; i < vertexes.rows(); i++){
for (int j = 0; j < 3; j++)
convertedMesh.vert[i].P()[j] = vertexes(i,j);
ivp.push_back(&convertedMesh.vert[i]);
}
FaceIterator fi;
for (i=0,fi=convertedMesh.face.begin();fi!=convertedMesh.face.end();i++,fi++)
for (int j = 0; j < 3; j++)
(*fi).V(j) = ivp[faces(i, j)];
}
}}
#endif // CONVERSIONMESHES_H

View File

@ -0,0 +1,53 @@
#ifndef CUBIC_STYLIZING_H
#define CUBIC_STYLIZING_H
#include <CubicStylizationFiles/src/cube_style_data.h>
#include <CubicStylizationFiles/src/cube_style_precomputation.h>
#include <CubicStylizationFiles/src/cube_style_single_iteration.h>
#include <CubicStylizationFiles/src/normalize_unitbox.h>
#include <vcg/complex/algorithms/mesh_to_matrix.h>
#include <vcg/complex/allocate.h>
#include <common/utilities/eigen_mesh_conversions.h>
#include <CubicStylizationFiles/conversionMeshes.h>
namespace vcg{
namespace tri{
template<class MeshType >
class Cubization{
public:
static void Init(MeshType& m, Eigen::MatrixXd& V, Eigen::MatrixXd& U, Eigen::MatrixXi& F, cube_style_data& data){
// check requirements
vcg::tri::VertexVectorHasPerVertexTexCoord( m.vert );
vcg::tri::VertexVectorHasPerVertexFlags( m.vert );
Mesh2Matrix(m, V, F);
//normalize_unitbox(V);
U = V;
data.bc.resize(1,3);
data.bc << V.row(F(0,0));
data.b.resize(1);
data.b << F(0,0);
// precomputation ARAP and initialize ADMM parameters
cube_style_precomputation(V,F,data);
}
static double Stylize(Eigen::MatrixXd &V, Eigen::MatrixXd &U, Eigen::MatrixXi &F, cube_style_data& data, Eigen::VectorXd& energy_verts, double& energyCubic)
{
// apply cubic stylization
cube_style_single_iteration(V,U,data, energy_verts);
energyCubic = data.objVal;
return data.reldV;
}
};
}}
#endif // CUBIC_STYLIZING_H

View File

@ -0,0 +1,373 @@
Mozilla Public License Version 2.0
==================================
1. Definitions
--------------
1.1. "Contributor"
means each individual or legal entity that creates, contributes to
the creation of, or owns Covered Software.
1.2. "Contributor Version"
means the combination of the Contributions of others (if any) used
by a Contributor and that particular Contributor's Contribution.
1.3. "Contribution"
means Covered Software of a particular Contributor.
1.4. "Covered Software"
means Source Code Form to which the initial Contributor has attached
the notice in Exhibit A, the Executable Form of such Source Code
Form, and Modifications of such Source Code Form, in each case
including portions thereof.
1.5. "Incompatible With Secondary Licenses"
means
(a) that the initial Contributor has attached the notice described
in Exhibit B to the Covered Software; or
(b) that the Covered Software was made available under the terms of
version 1.1 or earlier of the License, but not also under the
terms of a Secondary License.
1.6. "Executable Form"
means any form of the work other than Source Code Form.
1.7. "Larger Work"
means a work that combines Covered Software with other material, in
a separate file or files, that is not Covered Software.
1.8. "License"
means this document.
1.9. "Licensable"
means having the right to grant, to the maximum extent possible,
whether at the time of the initial grant or subsequently, any and
all of the rights conveyed by this License.
1.10. "Modifications"
means any of the following:
(a) any file in Source Code Form that results from an addition to,
deletion from, or modification of the contents of Covered
Software; or
(b) any new file in Source Code Form that contains any Covered
Software.
1.11. "Patent Claims" of a Contributor
means any patent claim(s), including without limitation, method,
process, and apparatus claims, in any patent Licensable by such
Contributor that would be infringed, but for the grant of the
License, by the making, using, selling, offering for sale, having
made, import, or transfer of either its Contributions or its
Contributor Version.
1.12. "Secondary License"
means either the GNU General Public License, Version 2.0, the GNU
Lesser General Public License, Version 2.1, the GNU Affero General
Public License, Version 3.0, or any later versions of those
licenses.
1.13. "Source Code Form"
means the form of the work preferred for making modifications.
1.14. "You" (or "Your")
means an individual or a legal entity exercising rights under this
License. For legal entities, "You" includes any entity that
controls, is controlled by, or is under common control with You. For
purposes of this definition, "control" means (a) the power, direct
or indirect, to cause the direction or management of such entity,
whether by contract or otherwise, or (b) ownership of more than
fifty percent (50%) of the outstanding shares or beneficial
ownership of such entity.
2. License Grants and Conditions
--------------------------------
2.1. Grants
Each Contributor hereby grants You a world-wide, royalty-free,
non-exclusive license:
(a) under intellectual property rights (other than patent or trademark)
Licensable by such Contributor to use, reproduce, make available,
modify, display, perform, distribute, and otherwise exploit its
Contributions, either on an unmodified basis, with Modifications, or
as part of a Larger Work; and
(b) under Patent Claims of such Contributor to make, use, sell, offer
for sale, have made, import, and otherwise transfer either its
Contributions or its Contributor Version.
2.2. Effective Date
The licenses granted in Section 2.1 with respect to any Contribution
become effective for each Contribution on the date the Contributor first
distributes such Contribution.
2.3. Limitations on Grant Scope
The licenses granted in this Section 2 are the only rights granted under
this License. No additional rights or licenses will be implied from the
distribution or licensing of Covered Software under this License.
Notwithstanding Section 2.1(b) above, no patent license is granted by a
Contributor:
(a) for any code that a Contributor has removed from Covered Software;
or
(b) for infringements caused by: (i) Your and any other third party's
modifications of Covered Software, or (ii) the combination of its
Contributions with other software (except as part of its Contributor
Version); or
(c) under Patent Claims infringed by Covered Software in the absence of
its Contributions.
This License does not grant any rights in the trademarks, service marks,
or logos of any Contributor (except as may be necessary to comply with
the notice requirements in Section 3.4).
2.4. Subsequent Licenses
No Contributor makes additional grants as a result of Your choice to
distribute the Covered Software under a subsequent version of this
License (see Section 10.2) or under the terms of a Secondary License (if
permitted under the terms of Section 3.3).
2.5. Representation
Each Contributor represents that the Contributor believes its
Contributions are its original creation(s) or it has sufficient rights
to grant the rights to its Contributions conveyed by this License.
2.6. Fair Use
This License is not intended to limit any rights You have under
applicable copyright doctrines of fair use, fair dealing, or other
equivalents.
2.7. Conditions
Sections 3.1, 3.2, 3.3, and 3.4 are conditions of the licenses granted
in Section 2.1.
3. Responsibilities
-------------------
3.1. Distribution of Source Form
All distribution of Covered Software in Source Code Form, including any
Modifications that You create or to which You contribute, must be under
the terms of this License. You must inform recipients that the Source
Code Form of the Covered Software is governed by the terms of this
License, and how they can obtain a copy of this License. You may not
attempt to alter or restrict the recipients' rights in the Source Code
Form.
3.2. Distribution of Executable Form
If You distribute Covered Software in Executable Form then:
(a) such Covered Software must also be made available in Source Code
Form, as described in Section 3.1, and You must inform recipients of
the Executable Form how they can obtain a copy of such Source Code
Form by reasonable means in a timely manner, at a charge no more
than the cost of distribution to the recipient; and
(b) You may distribute such Executable Form under the terms of this
License, or sublicense it under different terms, provided that the
license for the Executable Form does not attempt to limit or alter
the recipients' rights in the Source Code Form under this License.
3.3. Distribution of a Larger Work
You may create and distribute a Larger Work under terms of Your choice,
provided that You also comply with the requirements of this License for
the Covered Software. If the Larger Work is a combination of Covered
Software with a work governed by one or more Secondary Licenses, and the
Covered Software is not Incompatible With Secondary Licenses, this
License permits You to additionally distribute such Covered Software
under the terms of such Secondary License(s), so that the recipient of
the Larger Work may, at their option, further distribute the Covered
Software under the terms of either this License or such Secondary
License(s).
3.4. Notices
You may not remove or alter the substance of any license notices
(including copyright notices, patent notices, disclaimers of warranty,
or limitations of liability) contained within the Source Code Form of
the Covered Software, except that You may alter any license notices to
the extent required to remedy known factual inaccuracies.
3.5. Application of Additional Terms
You may choose to offer, and to charge a fee for, warranty, support,
indemnity or liability obligations to one or more recipients of Covered
Software. However, You may do so only on Your own behalf, and not on
behalf of any Contributor. You must make it absolutely clear that any
such warranty, support, indemnity, or liability obligation is offered by
You alone, and You hereby agree to indemnify every Contributor for any
liability incurred by such Contributor as a result of warranty, support,
indemnity or liability terms You offer. You may include additional
disclaimers of warranty and limitations of liability specific to any
jurisdiction.
4. Inability to Comply Due to Statute or Regulation
---------------------------------------------------
If it is impossible for You to comply with any of the terms of this
License with respect to some or all of the Covered Software due to
statute, judicial order, or regulation then You must: (a) comply with
the terms of this License to the maximum extent possible; and (b)
describe the limitations and the code they affect. Such description must
be placed in a text file included with all distributions of the Covered
Software under this License. Except to the extent prohibited by statute
or regulation, such description must be sufficiently detailed for a
recipient of ordinary skill to be able to understand it.
5. Termination
--------------
5.1. The rights granted under this License will terminate automatically
if You fail to comply with any of its terms. However, if You become
compliant, then the rights granted under this License from a particular
Contributor are reinstated (a) provisionally, unless and until such
Contributor explicitly and finally terminates Your grants, and (b) on an
ongoing basis, if such Contributor fails to notify You of the
non-compliance by some reasonable means prior to 60 days after You have
come back into compliance. Moreover, Your grants from a particular
Contributor are reinstated on an ongoing basis if such Contributor
notifies You of the non-compliance by some reasonable means, this is the
first time You have received notice of non-compliance with this License
from such Contributor, and You become compliant prior to 30 days after
Your receipt of the notice.
5.2. If You initiate litigation against any entity by asserting a patent
infringement claim (excluding declaratory judgment actions,
counter-claims, and cross-claims) alleging that a Contributor Version
directly or indirectly infringes any patent, then the rights granted to
You by any and all Contributors for the Covered Software under Section
2.1 of this License shall terminate.
5.3. In the event of termination under Sections 5.1 or 5.2 above, all
end user license agreements (excluding distributors and resellers) which
have been validly granted by You or Your distributors under this License
prior to termination shall survive termination.
************************************************************************
* *
* 6. Disclaimer of Warranty *
* ------------------------- *
* *
* Covered Software is provided under this License on an "as is" *
* basis, without warranty of any kind, either expressed, implied, or *
* statutory, including, without limitation, warranties that the *
* Covered Software is free of defects, merchantable, fit for a *
* particular purpose or non-infringing. The entire risk as to the *
* quality and performance of the Covered Software is with You. *
* Should any Covered Software prove defective in any respect, You *
* (not any Contributor) assume the cost of any necessary servicing, *
* repair, or correction. This disclaimer of warranty constitutes an *
* essential part of this License. No use of any Covered Software is *
* authorized under this License except under this disclaimer. *
* *
************************************************************************
************************************************************************
* *
* 7. Limitation of Liability *
* -------------------------- *
* *
* Under no circumstances and under no legal theory, whether tort *
* (including negligence), contract, or otherwise, shall any *
* Contributor, or anyone who distributes Covered Software as *
* permitted above, be liable to You for any direct, indirect, *
* special, incidental, or consequential damages of any character *
* including, without limitation, damages for lost profits, loss of *
* goodwill, work stoppage, computer failure or malfunction, or any *
* and all other commercial damages or losses, even if such party *
* shall have been informed of the possibility of such damages. This *
* limitation of liability shall not apply to liability for death or *
* personal injury resulting from such party's negligence to the *
* extent applicable law prohibits such limitation. Some *
* jurisdictions do not allow the exclusion or limitation of *
* incidental or consequential damages, so this exclusion and *
* limitation may not apply to You. *
* *
************************************************************************
8. Litigation
-------------
Any litigation relating to this License may be brought only in the
courts of a jurisdiction where the defendant maintains its principal
place of business and such litigation shall be governed by laws of that
jurisdiction, without reference to its conflict-of-law provisions.
Nothing in this Section shall prevent a party's ability to bring
cross-claims or counter-claims.
9. Miscellaneous
----------------
This License represents the complete agreement concerning the subject
matter hereof. If any provision of this License is held to be
unenforceable, such provision shall be reformed only to the extent
necessary to make it enforceable. Any law or regulation which provides
that the language of a contract shall be construed against the drafter
shall not be used to construe this License against a Contributor.
10. Versions of the License
---------------------------
10.1. New Versions
Mozilla Foundation is the license steward. Except as provided in Section
10.3, no one other than the license steward has the right to modify or
publish new versions of this License. Each version will be given a
distinguishing version number.
10.2. Effect of New Versions
You may distribute the Covered Software under the terms of the version
of the License under which You originally received the Covered Software,
or under the terms of any subsequent version published by the license
steward.
10.3. Modified Versions
If you create software not governed by this License, and you want to
create a new license for such software, you may create and use a
modified version of this License if you rename the license and remove
any references to the name of the license steward (except to note that
such modified license differs from this License).
10.4. Distributing Source Code Form that is Incompatible With Secondary
Licenses
If You choose to distribute Source Code Form that is Incompatible With
Secondary Licenses under the terms of this version of the License, the
notice described in Exhibit B of this License must be attached.
Exhibit A - Source Code Form License Notice
-------------------------------------------
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/.
If it is not possible or desirable to put the notice in a particular
file, then You may include the notice in a location (such as a LICENSE
file in a relevant directory) where a recipient would be likely to look
for such a notice.
You may add additional accurate notices of copyright ownership.
Exhibit B - "Incompatible With Secondary Licenses" Notice
---------------------------------------------------------
This Source Code Form is "Incompatible With Secondary Licenses", as
defined by the Mozilla Public License, v. 2.0.

View File

@ -0,0 +1,78 @@
#ifndef CUBE_STYLE_DATA_H
#define CUBE_STYLE_DATA_H
#include <Eigen/Core>
#include <Eigen/Sparse>
#include <limits>
#include <igl/min_quad_with_fixed.h>
struct cube_style_data
{
// user should tune these
double lambda = 0.0;
double rhoInit = 1e-3;
double ABSTOL = 1e-6;
double RELTOL = 1e-3;
// usually these don't need to tune
double mu = 10;
double tao = 2;
double maxIter_ADMM = 100;
double objVal = 0;
double reldV = std::numeric_limits<float>::max();
std::vector<double> objHis;
std::vector<Eigen::MatrixXi> hEList;
std::vector<Eigen::MatrixXd> dVList, UHis;
std::vector<Eigen::VectorXd> WVecList;
Eigen::SparseMatrix<double> K, L;
Eigen::MatrixXd N, VA, zAll, uAll;
Eigen::VectorXd rhoAll, objValVec;
// bool useBc = false;
Eigen::MatrixXd bc;
Eigen::VectorXi b;
igl::min_quad_with_fixed_data<double> solver_data;
// for plane constraints
Eigen::VectorXi bx, by, bz;
igl::min_quad_with_fixed_data<double> solver_data_x, solver_data_y, solver_data_z;
double xPlane = 0.0;
double yPlane = 0.0;
double zPlane = 0.0;
// for plane constraints
void reset()
{
// user should tune these
ABSTOL = 1e-5;
rhoInit = 1e-3;
RELTOL = 1e-3;
// usually these don't need to tune
double mu = 10;
double tao = 2;
double maxIter_ADMM = 100;
double objVal = 0;
double reldV = std::numeric_limits<float>::max();
objHis.clear();
hEList.clear();
dVList.clear();
UHis.clear();
WVecList.clear();
K = Eigen::SparseMatrix<double>();
L = Eigen::SparseMatrix<double>();
N = Eigen::MatrixXd();
VA = Eigen::MatrixXd();
zAll =Eigen::MatrixXd();
uAll = Eigen::MatrixXd();
rhoAll = Eigen::VectorXd();
objValVec = Eigen::VectorXd();
igl::min_quad_with_fixed_data<double> solver_data;
}
};
#endif // CUBE_STYLE_DATA_H

View File

@ -0,0 +1,70 @@
#include "cube_style_precomputation.h"
void cube_style_precomputation(
const Eigen::MatrixXd & V,
const Eigen::MatrixXi & F,
cube_style_data & data)
{
using namespace Eigen;
using namespace std;
data.reset();
igl::per_vertex_normals(V,F, data.N);
igl::cotmatrix(V,F,data.L);
SparseMatrix<double> M;
igl::massmatrix(V,F,igl::MASSMATRIX_TYPE_BARYCENTRIC,M);
data.VA = M.diagonal();
vector<vector<int>> adjFList, VI;
igl::vertex_triangle_adjacency(V.rows(),F,adjFList,VI);
igl::arap_rhs(V,F,V.cols(),igl::ARAP_ENERGY_TYPE_SPOKES_AND_RIMS, data.K);
data.hEList.resize(V.rows());
data.WVecList.resize(V.rows());
data.dVList.resize(V.rows());
vector<int> adjF;
for (int ii=0; ii<V.rows(); ii++)
{
adjF = adjFList[ii];
data.hEList[ii].resize(adjF.size()*3, 2);
data.WVecList[ii].resize(adjF.size()*3);
for (int jj=0; jj<adjF.size(); jj++)
{
int v0 = F(adjF[jj],0);
int v1 = F(adjF[jj],1);
int v2 = F(adjF[jj],2);
// compute adjacent half-edge indices of a vertex
data.hEList[ii](3*jj ,0) = v0;
data.hEList[ii](3*jj ,1) = v1;
data.hEList[ii](3*jj+1,0) = v1;
data.hEList[ii](3*jj+1,1) = v2;
data.hEList[ii](3*jj+2,0) = v2;
data.hEList[ii](3*jj+2,1) = v0;
// compute WVec = vec(W)
data.WVecList[ii](3*jj ) = data.L.coeff(v0,v1);
data.WVecList[ii](3*jj+1) = data.L.coeff(v1,v2);
data.WVecList[ii](3*jj+2) = data.L.coeff(v2,v0);
}
// compute [dV] matrix for each vertex
data.dVList[ii].resize(3, adjF.size()*3);
MatrixXd V_hE0, V_hE1;
igl::slice(V,data.hEList[ii].col(0),1,V_hE0);
igl::slice(V,data.hEList[ii].col(1),1,V_hE1);
// data.dVList[ii].block(0,0,3,data.hEList[ii].rows()) = (V_hE1 - V_hE0).transpose();
data.dVList[ii] = (V_hE1 - V_hE0).transpose();
}
igl::min_quad_with_fixed_precompute(data.L,data.b,SparseMatrix<double>(),false,data.solver_data);
data.zAll.resize(3, V.rows()); data.zAll.setRandom();
data.uAll.resize(3, V.rows()); data.uAll.setRandom();
data.rhoAll.resize(V.rows()); data.rhoAll.setConstant(data.rhoInit);
}

View File

@ -0,0 +1,24 @@
#ifndef CUBE_STYLE_PRECOMPUTATION_H
#define CUBE_STYLE_PRECOMPUTATION_H
#include <iostream>
#include <ctime>
#include <vector>
#include <Eigen/Core>
// include libigl functions
#include <igl/cotmatrix.h>
#include <igl/massmatrix.h>
#include <igl/per_vertex_normals.h>
#include <igl/vertex_triangle_adjacency.h>
#include <igl/arap_rhs.h>
#include <igl/columnize.h>
#include <igl/slice.h>
#include <igl/min_quad_with_fixed.h>
#include <CubicStylizationFiles/src/cube_style_data.h>
void cube_style_precomputation(
const Eigen::MatrixXd & V,
const Eigen::MatrixXi & F,
cube_style_data & data);
#endif

View File

@ -0,0 +1,41 @@
#include "cube_style_single_iteration.h"
void cube_style_single_iteration(
const Eigen::MatrixXd & V,
Eigen::MatrixXd & U,
cube_style_data & data,
Eigen::VectorXd & energyVects)
{
using namespace Eigen;
using namespace std;
// local step
MatrixXd RAll(3,V.rows()*3);
{
fit_rotations_l1(V, U, RAll, data, energyVects);
}
//R.setZero();
//R = RAll;
// global step
MatrixXd Upre = U;
{
VectorXd Rcol;
igl::columnize(RAll, V.rows(), 2, Rcol);
VectorXd Bcol = data.K * Rcol;
for(int dim=0; dim<V.cols(); dim++)
{
VectorXd Uc,Bc,bcc;
Bc = Bcol.block(dim*V.rows(),0,V.rows(),1);
bcc = data.bc.col(dim);
min_quad_with_fixed_solve(
data.solver_data,Bc,bcc,VectorXd(),Uc);
U.col(dim) = Uc;
}
}
// print optimization date
data.reldV = (U-Upre).cwiseAbs().maxCoeff() / (U-V).cwiseAbs().maxCoeff();
//cout << "reldV:" << scientific << data.reldV << ", obj:" << data.objVal << endl;
}

View File

@ -0,0 +1,23 @@
#ifndef CUBE_STYLE_SINGLE_ITERATION_H
#define CUBE_STYLE_SINGLE_ITERATION_H
#include <iostream>
#include <ctime>
#include <vector>
#include <Eigen/Core>
// include libigl functions
#include <igl/columnize.h>
#include <igl/slice.h>
#include <igl/min_quad_with_fixed.h>
// include cube flow functions
#include <CubicStylizationFiles/src/fit_rotations_l1.h>
#include <CubicStylizationFiles/src/cube_style_data.h>
void cube_style_single_iteration(
const Eigen::MatrixXd & V,
Eigen::MatrixXd & U,
cube_style_data & data,
Eigen::VectorXd & energyVects);
#endif // CUBE_STYLE_SINGLE_ITERATION_H

View File

@ -0,0 +1,100 @@
#include "fit_rotations_l1.h"
void fit_rotations_l1(
const Eigen::MatrixXd & V,
Eigen::MatrixXd & U,
Eigen::MatrixXd & RAll,
cube_style_data & data,
Eigen::VectorXd & energy_vects)
{
using namespace Eigen;
using namespace std;
data.objValVec.setZero(V.rows());
igl::parallel_for(
V.rows(),
[&data, &RAll, &U](const int ii)
{
// warm start parameters
VectorXd z = data.zAll.col(ii);
VectorXd u = data.uAll.col(ii);
VectorXd n = data.N.row(ii).transpose();
double rho = data.rhoAll(ii);
Matrix3d R;
// get energy parameters
// Note: dVn = [dV n], dUn = [dU z-u]
MatrixXi hE = data.hEList[ii];
MatrixXd dU(3,hE.rows());
{
MatrixXd U_hE0, U_hE1;
igl::slice(U,hE.col(0),1,U_hE0);
igl::slice(U,hE.col(1),1,U_hE1);
dU = (U_hE1 - U_hE0).transpose();
}
// Note:
// S = [dV n] * [W 0; 0 rho] * [dU (z-u)]'
// = dV * W * dU' + n * rho * (z-u)'
// = Spre + n * rho * (z-u)'
MatrixXd dV = data.dVList[ii];
VectorXd WVec = data.WVecList[ii];
Matrix3d Spre = dV * WVec.asDiagonal() * dU.transpose();
// ADMM
for (int k=0; k<data.maxIter_ADMM; k++)
{
// R step
Matrix3d S = Spre + (rho * n * (z-u).transpose());
// S /= S.norm();
orthogonal_procrustes(S, R);
// z step
VectorXd zOld = z;
shrinkage(R*n+u, data.lambda* data.VA(ii)/rho, z);
// u step
u.noalias() += R*n - z;
// compute residual
double r_norm = (z - R*n).norm();
double s_norm = (-rho * (z - zOld)).norm();
// rho step
if (r_norm > data.mu * s_norm)
{
rho = data.tao * rho;
u = u / data.tao;
}
else if (s_norm > data.mu * r_norm)
{
rho = rho / data.tao;
u = u * data.tao;
}
// stopping criteria
double nz = double(z.size());
double eps_pri = sqrt(2.0*nz)*data.ABSTOL + data.RELTOL*max( (R*n).norm(),z.norm() );
double eps_dual = sqrt(1.0*nz)*data.ABSTOL + data.RELTOL* ((rho*u).norm());
if ( (r_norm<eps_pri) && (s_norm<eps_dual) )
{
// save parameters
data.zAll.col(ii) = z;
data.uAll.col(ii) = u;
data.rhoAll(ii) = rho;
RAll.block(0,3*ii,3,3) = R;
// save objective
double objVal =
0.5*((R*dV-dU)*WVec.asDiagonal()*(R*dV-dU).transpose()).trace()
+ data.lambda * data.VA(ii) * (R*n).cwiseAbs().sum();
data.objValVec(ii) = objVal;
break;
}
} // ADMM end
}
,1000);
energy_vects = data.objValVec;
data.objVal = data.objValVec.sum();
}

View File

@ -0,0 +1,20 @@
#ifndef FIT_ROTATIONS_L1_H
#define FIT_ROTATIONS_L1_H
#include <iostream>
#include <Eigen/Core>
#include <Eigen/Sparse>
#include <CubicStylizationFiles/src/cube_style_data.h>
#include <CubicStylizationFiles/src/orthogonal_procrustes.h>
#include <CubicStylizationFiles/src/shrinkage.h>
#include <igl/slice.h>
#include <igl/parallel_for.h>
#include <math.h>
void fit_rotations_l1(
const Eigen::MatrixXd & V,
Eigen::MatrixXd & U,
Eigen::MatrixXd & RAll,
cube_style_data & data,
Eigen::VectorXd & energy_vects);
#endif // FIT_ROTATIONS_L1_H

View File

@ -0,0 +1,8 @@
#include "normalize_unitbox.h"
void normalize_unitbox(
Eigen::MatrixXd & V)
{
V = V.rowwise() - V.colwise().minCoeff();
V /= V.maxCoeff();
}

View File

@ -0,0 +1,15 @@
#ifndef NORMALIZE_UNITBOX_H
#define NORMALIZE_UNITBOX_H
#include <iostream>
#include <Eigen/Core>
// Inputs:
// V a matrix of vertex positions
// Outputs:
// V a matrix of vertex positions (in a unit box)
void normalize_unitbox(
Eigen::MatrixXd & V);
#endif // NORMALIZE_UNITBOX_H

View File

@ -0,0 +1,34 @@
#include "orthogonal_procrustes.h"
void orthogonal_procrustes(
const Eigen::Matrix3d & S,
Eigen::Matrix3d & R)
{
using namespace Eigen;
using namespace std;
// using svd3x3 which has problems for double precision
// Matrix3d SU, SV;
// Matrix<double,3,1> SS;
// igl::svd3x3(S, SU, SS, SV);
// R = SV * SU.transpose();
// if (R.determinant() < 0)
// {
// SU.col(2) = -SU.col(2);
// R = SV * SU.transpose();
// }
// using Eigen svd
JacobiSVD<Matrix3d> svd;
svd.compute(S, Eigen::ComputeFullU | Eigen::ComputeFullV );
Matrix3d SU = svd.matrixU();
Matrix3d SV = svd.matrixV();
R = SV * SU.transpose();
if (R.determinant() < 0)
{
SU.col(2) = -SU.col(2);
R = SV * SU.transpose();
}
assert(R.determinant() > 0);
}

View File

@ -0,0 +1,12 @@
#ifndef ORTHOGONAL_PROCRUSTES_H
#define ORTHOGONAL_PROCRUSTES_H
#include <iostream>
#include <Eigen/Core>
#include <igl/svd3x3.h>
#include <igl/polar_svd.h>
void orthogonal_procrustes(
const Eigen::Matrix3d & S,
Eigen::Matrix3d & R);
#endif // ORTHOGONAL_PROCRUSTES_H

View File

@ -0,0 +1,32 @@
#include "shrinkage.h"
void shrinkage(
const Eigen::VectorXd & x,
const double & k,
Eigen::VectorXd & z)
{
using namespace Eigen;
using namespace std;
VectorXd tmp1 = x.array() - k;
VectorXd posMax = tmp1.array().max(0.0);
VectorXd tmp2 = -x.array() - k;
VectorXd negMax = tmp2.array().max(0.0);
// VectorXd posMax = x.array() - k;
// for (int ii=0; ii<posMax.size(); ii++)
// {
// if (posMax(ii)<0)
// posMax(ii) = 0.0;
// }
// VectorXd negMax = -x.array() - k;
// for (int ii=0; ii<negMax.size(); ii++)
// {
// if (negMax(ii)<0)
// negMax(ii) = 0.0;
// }
z = posMax - negMax;
}

View File

@ -0,0 +1,11 @@
#ifndef SHRINKAGE_H
#define SHRINKAGE_H
#include <Eigen/Core>
#include <iostream>
void shrinkage(
const Eigen::VectorXd & x,
const double & k,
Eigen::VectorXd & z);
#endif // SHRINKAGE_H

View File

@ -0,0 +1,93 @@
/****************************************************************************
* VCGLib o o *
* Visual and Computer Graphics Library o o *
* _ O _ *
* Copyright(C) 2004 \/)\/ *
* Visual Computing Lab /\/| *
* ISTI - Italian National Research Council | *
* \ *
* All rights reserved. *
* *
* This program is free software; you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation; either version 2 of the License, or *
* (at your option) any later version. *
* *
* This program is distributed in the hope that it will be useful, *
* but WITHOUT ANY WARRANTY; without even the implied warranty of *
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
* GNU General Public License (http://www.gnu.org/licenses/gpl.txt) *
* for more details. *
* *
****************************************************************************/
#ifndef __CURVDATA
#define __CURVDATA
namespace vcg
{
class CurvData
{
public:
CurvData()
{
A = 0; // area
H = 0; // absolute mean curvature
K = 0; // gaussian curvature
}
virtual ~CurvData() {}
friend const CurvData operator+(const CurvData &lhs, const CurvData &rhs)
{
CurvData result;
result.A = lhs.A + rhs.A;
result.H = lhs.H + rhs.H;
result.K = lhs.K + rhs.K;
return result;
}
friend CurvData &operator+=(CurvData &lhs, const CurvData &rhs)
{
lhs.A += rhs.A;
lhs.H += rhs.H;
lhs.K += rhs.K;
return lhs;
}
float A;
float H;
float K;
};
class NSMCEval
{
public:
float operator() (const CurvData& c)
{
return (powf(c.H / 4.0, 2.0f) / c.A);
}
};
class MeanCEval
{
public:
float operator() (const CurvData& c)
{
return (c.H / 4.0);
}
};
class AbsCEval
{
public:
float operator() (const CurvData& c)
{
float k = 2 * M_PI - c.K;
if(k > 0) return 2.0 * (c.H / 4.0);
else return 2.0 * math::Sqrt(powf(c.H / 4.0, 2.0) - (c.A * k));
}
};
} //end namespace vcg
#endif // __CURVDATA

View File

@ -0,0 +1,339 @@
/****************************************************************************
* MeshLab o o *
* A versatile mesh processing toolbox o o *
* _ O _ *
* Copyright(C) 2005 \/)\/ *
* Visual Computing Lab /\/| *
* ISTI - Italian National Research Council | *
* \ *
* All rights reserved. *
* *
* This program is free software; you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation; either version 2 of the License, or *
* (at your option) any later version. *
* *
* This program is distributed in the hope that it will be useful, *
* but WITHOUT ANY WARRANTY; without even the implied warranty of *
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
* GNU General Public License (http://www.gnu.org/licenses/gpl.txt) *
* for more details. *
* *
****************************************************************************/
#ifndef __CURVEDGEFLIP
#define __CURVEDGEFLIP
#include <vcg/container/simple_temporary_data.h>
#include <vcg/complex/algorithms/local_optimization/tri_edge_flip.h>
#include <vcg/space/triangle3.h>
#include <vcg/space/point3.h>
#include <vcg/complex/algorithms/update/normal.h>
#include "curvdata.h"
namespace vcg
{
namespace tri
{
/*
* This flip happens only if decreases the curvature of the surface
* Edge flip optimization based on the paper
* "optimizing 3d triangulations using discrete curvature analysis"
* http://www2.in.tu-clausthal.de/~hormann/papers/Dyn.2001.OTU.pdf
*/
template <class TRIMESH_TYPE, class MYTYPE, class CURVEVAL>
class CurvEdgeFlip : public TopoEdgeFlip<TRIMESH_TYPE, MYTYPE>
{
protected:
typedef typename TRIMESH_TYPE::FaceType FaceType;
typedef typename TRIMESH_TYPE::FacePointer FacePointer;
typedef typename TRIMESH_TYPE::FaceIterator FaceIterator;
typedef typename TRIMESH_TYPE::VertexIterator VertexIterator;
typedef typename TRIMESH_TYPE::VertexType VertexType;
typedef typename TRIMESH_TYPE::ScalarType ScalarType;
typedef typename TRIMESH_TYPE::VertexPointer VertexPointer;
typedef typename TRIMESH_TYPE::CoordType CoordType;
typedef vcg::face::Pos<FaceType> PosType;
typedef vcg::face::VFIterator<FaceType> VFIteratorType;
typedef typename LocalOptimization<TRIMESH_TYPE>::HeapElem HeapElem;
typedef typename LocalOptimization<TRIMESH_TYPE>::HeapType HeapType;
typedef TriEdgeFlip<TRIMESH_TYPE, MYTYPE> Parent;
typedef typename vcg::Triangle3<ScalarType> TriangleType;
typedef typename TRIMESH_TYPE::VertContainer VertContainer;
// New curvature precomputed for the vertices of the faces
// adjacent to edge to be flipped
ScalarType _cv0, _cv1, _cv2, _cv3;
static CurvData FaceCurv(VertexPointer v0,
VertexPointer v1,
VertexPointer v2,
CoordType fNormal)
{
CurvData res;
float ang0 = math::Abs(Angle(v1->P() - v0->P(), v2->P() - v0->P()));
float ang1 = math::Abs(Angle(v0->P() - v1->P(), v2->P() - v1->P()));
float ang2 = M_PI - ang0 - ang1;
float s01 = SquaredDistance(v1->P(), v0->P());
float s02 = SquaredDistance(v2->P(), v0->P());
// voronoi cell of vertex
if (ang0 >= M_PI/2) {
TriangleType triangle(v0->P(), v1->P(), v2->P());
res.A += (0.5f * DoubleArea(triangle) -
(s01 * tan(ang1) + s02 * tan(ang2)) / 8.0 );
} else if (ang1 >= M_PI/2)
res.A += (s01 * tan(ang0)) / 8.0;
else if (ang2 >= M_PI/2)
res.A += (s02 * tan(ang0) / 8.0);
else // non obctuse triangle
res.A += ((s02 / tan(ang1)) + (s01 / tan(ang2))) / 8.0;
res.K += ang0;
ang1 = math::Abs(Angle(fNormal, v1->N()));
ang2 = math::Abs(Angle(fNormal, v2->N()));
res.H += ( (math::Sqrt(s01) / 2.0) * ang1 +
(math::Sqrt(s02) / 2.0) * ang2 );
return res;
}
/* Compute vertex curvature walking around a verte with VF adjacency
* f1, f2 --> this faces are to be ignored
* */
static CurvData Curvature(VertexPointer v, FacePointer f1 = NULL, FacePointer f2 = NULL)
{
CurvData curv;
VFIteratorType vfi(v);
while (!vfi.End()) {
if (vfi.F() != f1 && vfi.F() != f2 && !vfi.F()->IsD()) {
int i = vfi.I();
curv += FaceCurv(vfi.F()->V0(i),
vfi.F()->V1(i),
vfi.F()->V2(i),
vfi.F()->N());
}
++vfi;
}
return curv;
}
public:
CurvEdgeFlip() {}
CurvEdgeFlip(PosType pos, int mark,BaseParameterClass *pp)
{
this->_pos = pos;
this->_localMark = mark;
this->_priority = ComputePriority(pp);
}
CurvEdgeFlip(CurvEdgeFlip &par)
{
this->_pos = par.GetPos();
this->_localMark = par.GetMark();
this->_priority = par.Priority();
}
// temporary empty (flip is already done in constructor)
void Execute(TRIMESH_TYPE& m, BaseParameterClass *)
{
int i = this->_pos.E();
FacePointer f1 = this->_pos.F();
int j = f1->FFi(i);
FacePointer f2 = f1->FFp(i);
VertexPointer v0 = f1->V0(i), v1 = f1->V1(i),
v2 = f1->V2(i), v3 = f2->V2(j);
// save precomputed curvature in vertex quality
v0->Q() = _cv0;
v1->Q() = _cv1;
v2->Q() = _cv2;
v3->Q() = _cv3;
CoordType n1 = Normal(v0->P(), v3->P(), v2->P());
CoordType n2 = Normal(v1->P(), v2->P(), v3->P());
// compute and update normals on vertices after flip
v0->N() = v0->N() - f1->N() - f2->N() + n1;
v1->N() = v1->N() - f1->N() - f2->N() + n2;
v2->N() = v2->N() - f1->N() + n1 + n2;
v3->N() = v3->N() - f2->N() + n1 + n2;
// fix VF adjacency - detach
assert(f1->V1(i) == v1);
vcg::face::VFDetach(*f1, (i + 1) % 3); // detach v1 from f1
assert(f2->V1(j) == v0);
vcg::face::VFDetach(*f2, (j + 1) % 3); // detach v0 from f2
// do the flip
vcg::face::FlipEdge(*this->_pos.F(), this->_pos.E());
// fix VF adjacency - append
assert(f2->V1(j) == v2);
vcg::face::VFAppend(f2, (j + 1) % 3); // attach v2 with f2
assert(f1->V1(i) == v3);
vcg::face::VFAppend(f1, (i + 1) % 3); // attach v3 with f1
// update face normals after the flip
f1->N() = n1;
f2->N() = n2;
// avoid texture coordinates swap after flip
if (tri::HasPerWedgeTexCoord(m)) {
f2->WT((j + 1) % 3) = f1->WT((i + 2) % 3);
f1->WT((i + 1) % 3) = f2->WT((j + 2) % 3);
}
}
virtual bool IsFeasible(BaseParameterClass *_pp)
{
PlanarEdgeFlipParameter *pp=(PlanarEdgeFlipParameter *)_pp;
// First the flip must be topologically correct.
if(!vcg::face::CheckFlipEdge(*this->_pos.F(), this->_pos.E()))
return false;
// then the angle between the involved normals must be greater???
if (math::ToDeg(Angle(this->_pos.FFlip()->cN(), this->_pos.F()->cN()) ) <= pp->CoplanarAngleThresholdDeg )
return false;
CoordType v0, v1, v2, v3;
int i = this->_pos.E();
FacePointer f = this->_pos.F();
v0 = f->P0(i);
v1 = f->P1(i);
v2 = f->P2(i);
v3 = f->FFp(i)->P2(f->FFi(i));
// Take the parallelogram formed by the adjacent faces of edge
// If a corner of the parallelogram on extreme of edge to flip is >= 180
// the flip will produce two identical faces - avoid this
if( (Angle(v2 - v0, v1 - v0) + Angle(v3 - v0, v1 - v0) >= M_PI) ||
(Angle(v2 - v1, v0 - v1) + Angle(v3 - v1, v0 - v1) >= M_PI))
return false;
// if any of two faces adj to edge in non writable, the flip is unfeasible
if (!this->_pos.F()->IsW() || !this->_pos.F()->FFp(i)->IsW())
return false;
return true;
}
ScalarType ComputePriority(BaseParameterClass *pp)
{
/*
1
/|\
/ | \
2 | 3
\ | /
\|/
0
*/
if(!this->IsFeasible(pp))
return std::numeric_limits<ScalarType>::infinity();
VertexPointer v0, v1, v2, v3;
int i = this->_pos.E();
FacePointer f1 = this->_pos.F();
v0 = f1->V0(i);
v1 = f1->V1(i);
v2 = f1->V2(i);
FacePointer f2 = this->_pos.F()->FFp(i);
v3 = f2->V2(f1->FFi(i));
// save sum of curvatures of vertices
float cbefore = v0->Q() + v1->Q() + v2->Q() + v3->Q();
// saving current vertex normals
CoordType nv0orig = v0->N();
CoordType nv1orig = v1->N();
CoordType nv2orig = v2->N();
CoordType nv3orig = v3->N();
CurvData cd0, cd1, cd2, cd3;
CoordType n1 = Normal(v0->P(), v3->P(), v2->P());
CoordType n2 = Normal(v1->P(), v2->P(), v3->P());
// new vertex normal after the flip
// set this now to evaluate new curvatures
v0->N() = nv0orig - f1->N() - f2->N() + n1;
v1->N() = nv1orig - f1->N() - f2->N() + n2;
v2->N() = nv2orig - f1->N() + n1 + n2;
v3->N() = nv3orig - f2->N() + n1 + n2;
cd0 = FaceCurv(v0, v3, v2, n1) + Curvature(v0, f1, f2);
cd1 = FaceCurv(v1, v2, v3, n2) + Curvature(v1, f1, f2);
cd2 = FaceCurv(v2, v0, v3, n1) + FaceCurv(v2, v3, v1, n2) + Curvature(v2, f1, f2);
cd3 = FaceCurv(v3, v2, v0, n1) + FaceCurv(v3, v1, v2, n2) + Curvature(v3, f1, f2);
// restore original pervertex normals after computation
v0->N() = nv0orig;
v1->N() = nv1orig;
v2->N() = nv2orig;
v3->N() = nv3orig;
CURVEVAL curveval;
_cv0 = curveval(cd0);
_cv1 = curveval(cd1);
_cv2 = curveval(cd2);
_cv3 = curveval(cd3);
float cafter = _cv0 + _cv1 + _cv2 + _cv3;
// The priority of an edge flip is **how much we lower the overall curvature**.
// If after the flip the sum of the curvature is decreased it is a good move;
// good flips have: cafter < cbefore
// Since the local optimization is designed to make the minimum cost move we put inside
// negative values (the more negative the better).
this->_priority = (cafter - cbefore);
//qDebug("computed curvature change, %f->%f (priority = %f)", cbefore,cafter,this->_priority);
return this->_priority;
}
static void Init(TRIMESH_TYPE &m, HeapType &heap, BaseParameterClass *pp)
{
CURVEVAL curveval;
heap.clear();
// comuputing edge flip priority require non normalized vertex normals AND non normalized face normals.
vcg::tri::UpdateNormal<TRIMESH_TYPE>::PerVertexPerFace(m);
VertexIterator vi;
for (vi = m.vert.begin(); vi != m.vert.end(); ++vi)
if (!(*vi).IsD() && (*vi).IsW())
(*vi).Q() = curveval(Curvature(&(*vi)));
FaceIterator fi;
for (fi = m.face.begin(); fi != m.face.end(); ++fi)
if (!(*fi).IsD())
for (unsigned int i = 0; i < 3; i++)
if ((*fi).V1(i) - (*fi).V0(i) > 0) {
PosType newpos(&*fi, i);
TopoEdgeFlip<TRIMESH_TYPE, MYTYPE>::Insert(heap, newpos, tri::IMark(m),pp);
}
}
}; // end CurvEdgeFlip class
}; // namespace tri
}; // namespace vcg
#endif // __CURVEDGEFLIP

View File

@ -0,0 +1,354 @@
/****************************************************************************
* MeshLab o o *
* A versatile mesh processing toolbox o o *
* _ O _ *
* Copyright(C) 2005 \/)\/ *
* Visual Computing Lab /\/| *
* ISTI - Italian National Research Council | *
* \ *
* All rights reserved. *
* *
* This program is free software; you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation; either version 2 of the License, or *
* (at your option) any later version. *
* *
* This program is distributed in the hope that it will be useful, *
* but WITHOUT ANY WARRANTY; without even the implied warranty of *
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
* GNU General Public License (http://www.gnu.org/licenses/gpl.txt) *
* for more details. *
* *
****************************************************************************/
#include <math.h>
#include <stdlib.h>
#include <time.h>
#include "filter_cubization.h"
#include <vcg/complex/algorithms/clean.h>
#include <vcg/complex/algorithms/smooth.h>
#include <CubicStylizationFiles/cubic_stylizing.h>
#include <CubicStylizationFiles/src/cube_style_data.h>
#include <vcg/complex/algorithms/local_optimization.h>
#include <vcg/complex/algorithms/local_optimization/tri_edge_flip.h>
#include "curvedgeflip.h"
#include "curvdata.h"
using namespace std;
using namespace vcg;
// forward declarations
class MeanCEFlip;
typedef Histogram<Scalarm> Histogramm;
class MeanCEFlip : public vcg::tri::CurvEdgeFlip<CMeshO, MeanCEFlip, MeanCEval >
{
public:
MeanCEFlip(PosType pos, int mark,BaseParameterClass *_pp) :
vcg::tri::CurvEdgeFlip<CMeshO, MeanCEFlip, MeanCEval >(pos, mark,_pp) {}
};
// 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
CubizationPlugin::CubizationPlugin()
{
typeList = {
FP_CUBIZATION
};
for(ActionIDType tt: types())
actionList.push_back(new QAction(filterName(tt), this));
cubic_ApplyEdgeFlip = false;
cubic_ApplyColorize = false;
}
QString CubizationPlugin::pluginName() const
{
return "FilterCubization";
}
QString CubizationPlugin::pythonFilterName(ActionIDType f) const
{
switch (f) {
case FP_CUBIZATION: return tr("apply_coord_cubic_stylization");
default: assert(0); return QString();
}
}
QString CubizationPlugin::filterName(ActionIDType filterId) const
{
switch (filterId) {
case FP_CUBIZATION: return tr("Cubic stylization");
default: assert(0); return QString();
}
}
FilterPlugin::FilterArity CubizationPlugin::filterArity(const QAction*) const
{
return SINGLE_MESH;
}
int CubizationPlugin::getPreConditions(const QAction*) const
{
return MeshModel::MM_VERTCOORD | MeshModel::MM_FACEVERT;
}
int CubizationPlugin::getRequirements(const QAction *action)
{
switch (ID(action)) {
case FP_CUBIZATION:
return MeshModel::MM_VERTCOORD |
MeshModel::MM_FACEFACETOPO |
MeshModel::MM_VERTFACETOPO |
MeshModel::MM_VERTMARK;
}
return 0;
}
QString CubizationPlugin::filterInfo(ActionIDType filterId) const
{
switch(filterId) {
case FP_CUBIZATION:
return tr("Turn a mesh into a cube's style maintaining its original shape."
"<br>For all detailed about cubic stylization see:<br>"
"<br><i> Hsueh-Ti Derek Liu and Alec Jacobson.</i><br>"
"<b>Cubic Stylization</b> (<a href='https://www.dgp.toronto.edu/projects/cubic-stylization/cubicStyle_high.pdf'>pdf</a>)<br>"
"in ACM Transactions on Graphics, 2019<br/><br/> ");
default : assert(0);
}
return {};
}
CubizationPlugin::FilterClass CubizationPlugin::getClass(const QAction *action) const
{
switch(ID(action)) {
case FP_CUBIZATION: return FilterPlugin::Remeshing;
}
return FilterPlugin::Generic;
}
int CubizationPlugin::postCondition(const QAction *a) const
{
switch(ID(a))
{
case FP_CUBIZATION : return MeshModel::MM_ALL |
MeshModel::MM_VERTQUALITY;
default : assert(0);
}
return {};
}
// 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 string shown in the dialog
// - the default value
// - a possibly long string describing the meaning of that parameter (shown as a popup help in the dialog)
RichParameterList CubizationPlugin::initParameterList(const QAction *action, const MeshModel &m)
{
RichParameterList parlst;
if (ID(action) == FP_CUBIZATION) {
parlst.addParam(RichFloat("lcubeness", 0.2f,
tr("Cubeness parameter (λ)"),
tr("Control the cubeness of the mesh. Generally, the higher the cubeness parameter, the more cubic the mesh is. λ ∈ [0, 1] ")));
parlst.addParam(RichBool("applyef", cubic_ApplyEdgeFlip, tr("Apply edge flipping"), tr("Apply edge flip optimization on cubic stylization.")));
parlst.addParam(RichBool("applycol", cubic_ApplyColorize, tr("Colorize by vertex Quality"), tr("Color vertices depending on their cubization energy.")));
}
return parlst;
}
void Freeze(MeshModel &m)
{
tri::UpdatePosition<CMeshO>::Matrix(m.cm, m.cm.Tr,true);
tri::UpdateBounding<CMeshO>::Box(m.cm);
m.cm.shot.ApplyRigidTransformation(m.cm.Tr);
m.cm.Tr.SetIdentity();
}
void ApplyTransform(MeshModel &m, const Matrix44m &tr, bool freeze,
bool invertFlag=false, bool composeFlage=true)
{
if(invertFlag) m.cm.Tr = Inverse(m.cm.Tr);
if(composeFlage) m.cm.Tr = tr * m.cm.Tr;
else m.cm.Tr=tr;
if(freeze) Freeze(m);
}
std::map<std::string, QVariant> CubizationPlugin::applyFilter(
const QAction* filter,
const RichParameterList& par,
MeshDocument& md,
unsigned int&,
vcg::CallBackPos*)
{
if (ID(filter) == FP_CUBIZATION) {
// get bounding box
Box3m bbox = md.mm()->cm.bbox;
Matrix44m transfM, prevScaleTran, scaleTran, trTran, trTranInv, prevTransfM;
Point3m tranVec = Point3m(0, 0, 0);
float maxSide = max(bbox.DimX(), max(bbox.DimY(), bbox.DimZ()));
if (maxSide >= 1.0) {
// compute unit box scale
scaleTran.SetScale(1.0 / maxSide, 1.0 / maxSide, 1.0 / maxSide);
// reverse scale computation
prevScaleTran.SetScale(
1.0 / (1.0 / maxSide), 1.0 / (1.0 / maxSide), 1.0 / (1.0 / maxSide));
// compute translate
trTran.SetTranslate(tranVec);
trTranInv.SetTranslate(-tranVec);
transfM = trTran * scaleTran * trTranInv;
prevTransfM = trTran * prevScaleTran * trTranInv;
}
MeshModel& m = *(md.mm());
if (maxSide >= 1.0) {
// apply transform
ApplyTransform(m, transfM, true);
}
double energyTotal = 0.f;
time_t start = clock();
bool isColorizing = par.getBool("applycol");
m = ComputeCubicStylization(md, par, energyTotal, isColorizing);
m.updateDataMask(MeshModel::MM_VERTQUALITY);
if (isColorizing) {
m.updateDataMask(MeshModel::MM_FACEFACETOPO);
m.updateDataMask(MeshModel::MM_VERTCOLOR);
Histogramm H;
vcg::tri::Stat<CMeshO>::ComputePerVertexQualityHistogram(m.cm, H);
vcg::tri::UpdateColor<CMeshO>::PerVertexQualityRamp(
m.cm, H.Percentile(0.01f), H.Percentile(0.99f));
}
if (maxSide >= 1.0) {
// re-apply transform to get previous mesh scale
ApplyTransform(m, prevTransfM, true);
}
m.updateBoxAndNormals();
log("cubic stylization performed in %.2f sec. with cubic energy equal to %.5f",
(clock() - start) / (float) CLOCKS_PER_SEC,
energyTotal);
}
else {
wrongActionCalled(filter);
}
return std::map<std::string, QVariant>();
}
MeshModel CubizationPlugin::ComputeCubicStylization(
MeshDocument& md,
const RichParameterList& par,
double& totalEnergy,
bool isColorizing){
MeshModel &m=*(md.mm());
float limit = -std::numeric_limits<float>::epsilon();
int delvert = tri::Clean<CMeshO>::RemoveUnreferencedVertex(m.cm);
if (delvert)
log( "Pre-Curvature Cleaning: Removed %d unreferenced vertices",
delvert);
tri::Allocator<CMeshO>::CompactVertexVector(m.cm);
tri::Allocator<CMeshO>::CompactFaceVector(m.cm);
m.updateDataMask(MeshModel::MM_FACEFACETOPO);
vcg::tri::UpdateFlags<CMeshO>::FaceBorderFromFF(m.cm);
bool isApplyEdgeFlip = par.getBool("applyef");
if ( tri::Clean<CMeshO>::CountNonManifoldEdgeFF(m.cm) >0 && isApplyEdgeFlip) {
throw MLException("Mesh has some not 2-manifold faces, edge flips requires manifoldness");
}
Eigen::MatrixXd verts;
Eigen::MatrixXi faces;
Eigen::MatrixXd u_verts;
Eigen::VectorXd energy_verts;
float lambda = par.getFloat("lcubeness");
cube_style_data data;
data.lambda = lambda;
vcg::tri::Cubization<CMeshO>::Init(m.cm, verts, u_verts, faces, data);
// apply cubic stylization
int maxIter = 1000;
double stopReldV = 1e-3; // stopping criteria for relative displacement
double reldV = 0;
for (int iter=0; iter<maxIter; iter++)
{
reldV = vcg::tri::Cubization<CMeshO>::Stylize(verts, u_verts, faces, data, energy_verts, totalEnergy);
//apply Edge Flips
if((iter%30 == 0 || reldV < stopReldV) && isApplyEdgeFlip){
Matrix2Mesh(m.cm, u_verts, faces);
vcg::tri::PlanarEdgeFlipParameter pp;
float limit = -std::numeric_limits<float>::epsilon();
vcg::LocalOptimization<CMeshO> optimiz(m.cm,&pp);
float pthr = 10;
// VF adjacency needed for edge flips based on vertex curvature
vcg::tri::UpdateTopology<CMeshO>::VertexFace(m.cm);
vcg::tri::UpdateTopology<CMeshO>::TestVertexFace(m.cm);
pp.CoplanarAngleThresholdDeg = pthr;
optimiz.Init<MeanCEFlip>();
// stop when flips become harmful
optimiz.SetTargetMetric(limit);
optimiz.DoOptimization();
optimiz.h.clear();
Mesh2Matrix(m.cm, u_verts, faces);
log( "Iteration %d: %d curvature edge flips performed", iter, optimiz.nPerformedOps);
}
if (reldV < stopReldV) break;
}
Matrix2Mesh(m.cm, u_verts, faces);
if(isColorizing){
assert(energy_verts.size() == m.cm.vert.size());
for(int i = 0; i < energy_verts.size(); i++){
m.cm.vert[i].Q() = energy_verts[i];
}
}
return m;
}
MESHLAB_PLUGIN_NAME_EXPORTER(CubizationPlugin)

View File

@ -0,0 +1,73 @@
/****************************************************************************
* MeshLab o o *
* A versatile mesh processing toolbox o o *
* _ O _ *
* Copyright(C) 2005 \/)\/ *
* Visual Computing Lab /\/| *
* ISTI - Italian National Research Council | *
* \ *
* All rights reserved. *
* *
* This program is free software; you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation; either version 2 of the License, or *
* (at your option) any later version. *
* *
* This program is distributed in the hope that it will be useful, *
* but WITHOUT ANY WARRANTY; without even the implied warranty of *
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
* GNU General Public License (http://www.gnu.org/licenses/gpl.txt) *
* for more details. *
* *
****************************************************************************/
#ifndef FILTER_CUBIZATION_H
#define FILTER_CUBIZATION_H
#include <QObject>
#include <common/plugins/interfaces/filter_plugin.h>
class CubizationPlugin : public QObject, public FilterPlugin
{
Q_OBJECT
MESHLAB_PLUGIN_IID_EXPORTER(FILTER_PLUGIN_IID)
Q_INTERFACES(FilterPlugin)
public:
enum {
// mesh improvement by edge flipping
FP_CUBIZATION
};
CubizationPlugin();
QString pythonFilterName(ActionIDType f) const;
QString pluginName() const;
QString filterName(ActionIDType filter) const;
QString filterInfo(ActionIDType filter) 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);
int getPreConditions(const QAction *) const;
int getRequirements(const QAction*);
FilterClass getClass(const QAction *) const;
int postCondition(const QAction* ) const;
FilterArity filterArity(const QAction*) const;
private:
MeshModel ComputeCubicStylization(
MeshDocument& md,
const RichParameterList& par,
double& totalEnergy,
bool isColorizing);
protected:
bool cubic_ApplyEdgeFlip;
bool cubic_ApplyColorize;
};
#endif // FILTER_CUBIZATION_H

@ -1 +1 @@
Subproject commit fa707e78ba38ee93aeb0f3edbfe5957852dff61b
Subproject commit 77db079863e5279625b088af35f98eb148f1d8fd