From 07f2983e87cb8a6ce927d3f543f7eb43d077bc79 Mon Sep 17 00:00:00 2001 From: Gael Guennebaud guennebaud Date: Tue, 21 Oct 2008 17:32:06 +0000 Subject: [PATCH] mlsfilter: add color from curvature filters --- src/meshlabplugins/mlsfilters/apss.h | 5 +- src/meshlabplugins/mlsfilters/apss.tpp | 19 ++ src/meshlabplugins/mlsfilters/mlsplugin.cpp | 234 +++++++------------- src/meshlabplugins/mlsfilters/mlsplugin.h | 21 +- 4 files changed, 120 insertions(+), 159 deletions(-) diff --git a/src/meshlabplugins/mlsfilters/apss.h b/src/meshlabplugins/mlsfilters/apss.h index 2573d4a3f..59ebd7fb5 100644 --- a/src/meshlabplugins/mlsfilters/apss.h +++ b/src/meshlabplugins/mlsfilters/apss.h @@ -63,9 +63,12 @@ class APSS : public MlsSurface<_MeshType> virtual Scalar potential(const VectorType& x, int* errorMask = 0) const; virtual VectorType gradient(const VectorType& x, int* errorMask = 0) const; - virtual MatrixType hessian(const VectorType& x, int* errorMask) const; + virtual MatrixType hessian(const VectorType& x, int* errorMask) const; virtual VectorType project(const VectorType& x, VectorType* pNormal = 0, int* errorMask = 0) const; + /** \returns the approximation of the mean curvature obtained from the radius of the fitted sphere */ + virtual Scalar approxMeanCurvature(const VectorType& x, int* errorMask = 0) const; + void setSphericalParameter(Scalar v); protected: diff --git a/src/meshlabplugins/mlsfilters/apss.tpp b/src/meshlabplugins/mlsfilters/apss.tpp index cf6539592..02b0cad6b 100644 --- a/src/meshlabplugins/mlsfilters/apss.tpp +++ b/src/meshlabplugins/mlsfilters/apss.tpp @@ -65,6 +65,25 @@ typename APSS<_MeshType>::Scalar APSS<_MeshType>::potential(const VectorType& x, } } +template +typename APSS<_MeshType>::Scalar APSS<_MeshType>::approxMeanCurvature(const VectorType& x, int* errorMask) const +{ + if ((!mCachedQueryPointIsOK) || mCachedQueryPoint!=x) + { + if (!fit(x)) + { + if (errorMask) + *errorMask = MLS_TOO_FAR; + return Base::InvalidValue(); + } + } + + if (mStatus==ASS_SPHERE) + return (uQuad>0.?1.0:-1.0)/mRadius; + else + return 0.; +} + template typename APSS<_MeshType>::VectorType APSS<_MeshType>::gradient(const VectorType& x, int* errorMask) const { diff --git a/src/meshlabplugins/mlsfilters/mlsplugin.cpp b/src/meshlabplugins/mlsfilters/mlsplugin.cpp index 3da610116..763f02566 100644 --- a/src/meshlabplugins/mlsfilters/mlsplugin.cpp +++ b/src/meshlabplugins/mlsfilters/mlsplugin.cpp @@ -41,8 +41,7 @@ #include #include #include -// #include -// #include +#include #include "mlsmarchingcube.h" #include "mlsplugin.h" @@ -54,6 +53,7 @@ #include "smallcomponentselection.h" using namespace GaelMls; +using namespace vcg; // Constructor usually performs only two simple tasks of filling the two lists // - typeList: with all the possible id of the filtering actions @@ -65,6 +65,7 @@ MlsPlugin::MlsPlugin() << FP_RIMLS_PROJECTION << FP_APSS_PROJECTION // << FP_RIMLS_AFRONT << FP_APSS_AFRONT << FP_RIMLS_MCUBE << FP_APSS_MCUBE + << FP_RIMLS_COLORIZE << FP_APSS_COLORIZE << FP_RADIUS_FROM_DENSITY << FP_SELECT_SMALL_COMPONENTS; @@ -84,6 +85,8 @@ const QString MlsPlugin::filterName(FilterIDType filterId) case FP_RIMLS_AFRONT : return QString("MLS meshing/##### Advancing Front"); case FP_APSS_MCUBE : return QString("Marching Cubes (APSS)"); case FP_RIMLS_MCUBE : return QString("Marching Cubes (#####)"); + case FP_APSS_COLORIZE : return QString("Colorize curvature (APSS)"); + case FP_RIMLS_COLORIZE : return QString("Colorize curvature (#####)"); case FP_RADIUS_FROM_DENSITY : return QString("Estimate radius from density"); case FP_SELECT_SMALL_COMPONENTS : return QString("Small component selection"); default : assert(0); @@ -119,6 +122,11 @@ const QString MlsPlugin::filterInfo(FilterIDType filterId) "step onto the MLS, and an extra zero removal procedure.
"; } + if (filterId & _COLORIZE_) + { + str += "Colorize the vertices of a mesh or point set using the curfvature of the underlying surface.
"; + } + if (filterId & _APSS_) { str += @@ -179,7 +187,7 @@ void MlsPlugin::initParameterSet(QAction* action, MeshDocument& md, FilterParame return; } - if (id & _PROJECTION_) + if ((id & _PROJECTION_) || (id & _COLORIZE_)) { parlst.addMesh( "ControlMesh", target, "Point set", "The point set (or mesh) which defines the MLS surface."); @@ -218,7 +226,8 @@ void MlsPlugin::initParameterSet(QAction* action, MeshDocument& md, FilterParame "1 to a pure spherical fit, values between 0 and 1 gives intermediate results," "while others real values might give interresting results, but take care with extreme" "settings !"); - parlst.addBool( "AccurateNormal", + if (!(id & _COLORIZE_)) + parlst.addBool( "AccurateNormal", true, "Accurate normals", "If checked, use the accurate MLS gradient instead of the local approximation" @@ -255,6 +264,14 @@ void MlsPlugin::initParameterSet(QAction* action, MeshDocument& md, FilterParame { } + if ((id & _COLORIZE_) && (id & _APSS_)) + { + parlst.addBool( "ApproxCurvature", + false, + "Approx curvature", + "If checked, use the radius of the fitted sphere as an approximation of the mean curvature."); + } + if (id & _MCUBE_) { parlst.addInt( "Resolution", @@ -291,144 +308,6 @@ struct EdgeAnglePredicate } }; -#if 0 -namespace vcg { -namespace tri { - -template class AdvancingMLS : public AdvancingFront { -public: - typedef typename MeshType::VertexType VertexType; - typedef typename MeshType::VertexIterator VertexIterator; - typedef typename MeshType::FaceType FaceType; - typedef typename MeshType::FaceIterator FaceIterator; - - typedef typename MeshType::ScalarType ScalarType; - typedef typename MeshType::VertexType::CoordType VectorType; - - typedef MlsSurface MlsType; - - AdvancingMLS(MeshType &_mesh, MlsType& mls) - : AdvancingFront(_mesh), mMls(mls) - {} - bool _SeedFace() {return this->SeedFace();} - bool Seed(int &v0, int &v1, int &v2) - { - // pick the first poin and build a face around - VectorType p0 = mMls.points()[0]; - VectorType n0 = mMls.normals()[0]; - ScalarType r0 = mMls.radii()[0]; - - p0 = mMls.project(p0, &n0); - VectorType vu = vcg::Cross(n0,VectorType(0,1,0)); - if (Norm(vu)<0.1) - vu = Cross(n0,VectorType(0,0,1)); - vu.Normalize(); - VectorType vv = Cross(n0, vu); - VertexType v[3]; - v[0].P() = VectorType(p0[0], p0[1], p0[2]); - v[0].N() = n0; - v[1].P() = mMls.project(v[0].P() + vu * r0, &v[1].N()); - v[2].P() = mMls.project(v[0].P() + vu * (r0*0.5) + vv * (r0*sqrt(3)*0.5), &v[2].N() ); - - std::cout << r0 << " " - << vcg::Distance(v[0].P(), v[1].P()) << " " - << vcg::Distance(v[1].P(), v[2].P()) << " " - << vcg::Distance(v[2].P(), v[0].P()) << "\n"; - - v[0].ClearFlags(); - v[1].ClearFlags(); - v[2].ClearFlags(); - - v0 = this->mesh.vert.size(); - AddVertex(v[0]); - v1 = this->mesh.vert.size(); - AddVertex(v[1]); - v2 = this->mesh.vert.size(); - AddVertex(v[2]); - - std::cout << "seed created\n"; - - return true; - } - - int Place(FrontEdge &e, std::list::iterator &touch) - { - VectorType p[3]; - p[0] = this->mesh.vert[e.v0].P(); - p[1] = this->mesh.vert[e.v1].P(); - p[2] = this->mesh.vert[e.v2].P(); - VectorType point = p[0] + p[1] - p[2]; - ScalarType scale = 0.33 * (vcg::Distance(p[0], p[1]) + vcg::Distance(p[1], p[2]) + vcg::Distance(p[0], p[2])); - std::cout - << vcg::Distance(p[0], p[1]) << " " - << vcg::Distance(p[1], p[2]) << " " - << vcg::Distance(p[2], p[0]) << "\n"; - - // project p2 on the tangent plane - VectorType m = (p[0]+p[1])*0.5; - VectorType n0 = this->mesh.vert[e.v0].N() + this->mesh.vert[e.v1].N(); - n0.Normalize(); -// VectorType h = vcg::Cross(p[0]-p[1],n0); -// h.Normalize(); -// VectorType pp2 = p[2] - n0 * (vcg::Dot(p[2]-m,n0)); -// if (vcg::Dot(pp2-m, h)>0) -// h = -h; -// VectorType point = m + h * (scale*sqrt(3)*0.5); - - point = mMls.project(point, &n0); - - int vn = this->mesh.vert.size(); - // find closest - int minId = -1; - ScalarType minD2 = 0.5*scale; - ScalarType bm = 1.4*scale; - bm = bm * bm; - minD2 = minD2*minD2; - for(int i = 0; i < this->mesh.vert.size(); i++) - { - ScalarType d2 = (this->mesh.vert[i].P() - point).SquaredNorm(); - if (d2mesh.vert[i].IsB() - && (this->mesh.vert[i].P() - m).SquaredNorm() < bm ) - { - minId = i; - minD2 = d2; - } - } - if (minId>=0) - { -// if((this->mesh.vert[i].P() - point).Norm() < 0.5*scale && i!= e.v0 && i!=e.v1 && i!=e.v2) { - std::cout << "reuse existing vertex\n"; - vn = minId; - //find the border - //assert(this->mesh.vert[i].IsB()); - for(std::list::iterator k = this->front.begin(); k != this->front.end(); k++) - if((*k).v0 == vn) touch = k; - for(std::list::iterator k = this->deads.begin(); k != this->deads.end(); k++) - if((*k).v0 == vn) touch = k; -// break; -// } - } - else //if(vn == this->mesh.vert.size()) - { - std::cout << "create new vertex\n"; - VertexType v; - v.P() = point;//mMls.project(point, &v.N()); - v.N() = n0; - v.ClearFlags(); - AddVertex(v); - } - return vn; - } - - MlsType& mMls; -// ScalarType r; -}; - -} -} - -#endif - /** compute the normal of a face as the average of its vertices */ template void UpdateFaceNormalFromVertex(MeshType& m) @@ -453,7 +332,6 @@ bool MlsPlugin::applyFilter(QAction* filter, MeshDocument& md, FilterParameterSe { int id = ID(filter); - // get source and dest meshes if (id == FP_RADIUS_FROM_DENSITY) { CMeshO::VertContainer& points = md.mm()->cm.vert; @@ -535,7 +413,8 @@ bool MlsPlugin::applyFilter(QAction* filter, MeshDocument& md, FilterParameterSe if (apss) { apss->setSphericalParameter(par.getFloat("SphericalParameter")); - apss->setGradientHint(par.getBool("AccurateNormal") ? GaelMls::MLS_DERIVATIVE_ACCURATE : GaelMls::MLS_DERIVATIVE_APPROX); + if (!(id & _COLORIZE_)) + apss->setGradientHint(par.getBool("AccurateNormal") ? GaelMls::MLS_DERIVATIVE_ACCURATE : GaelMls::MLS_DERIVATIVE_APPROX); } MeshModel * mesh = 0; @@ -576,6 +455,61 @@ bool MlsPlugin::applyFilter(QAction* filter, MeshDocument& md, FilterParameterSe Log(0, "Successfully projected %i vertices", mesh->cm.vn); } + else if (id & _COLORIZE_) + { + mesh = par.getMesh("ProxyMesh"); + + bool selectionOnly = par.getBool("SelectionOnly"); + bool approx = apss && par.getBool("ApproxCurvature"); + + int size = mesh->cm.vert.size(); + //std::vector curvatures(size); + float minc=1e9, maxc=-1e9, minabsc=1e9; + vcg::Point3f grad; + vcg::Matrix33f hess; + // pass 1: computes curvatures and statistics + for (unsigned int i = 0; i< size; i++) + { + cb(1+98*i/pPoints->cm.vert.size(), "MLS colorization..."); + + if ( (!selectionOnly) || (pPoints->cm.vert[i].IsS()) ) + { +// grad = mls->gradient(mesh->cm.vert[i].P()); +// hess = mls->hessian(mesh->cm.vert[i].P()); +// curvatures[i] = mls->meanCurvature(grad,hess); + //float c = mesh->cm.vert[i].Q() = mls->meanCurvature(grad,hess); + float c = 0; + if (approx) + c = mesh->cm.vert[i].Q() = apss->approxMeanCurvature(mesh->cm.vert[i].P()); + else + { + grad = mls->gradient(mesh->cm.vert[i].P()); + hess = mls->hessian(mesh->cm.vert[i].P()); + c = mesh->cm.vert[i].Q() = mls->meanCurvature(grad,hess); + } + minc = std::min(c,minc); + maxc = std::max(c,maxc); + minabsc = std::min(fabsf(c),minabsc); + } + } + // pass 2: convert the curvature to color + cb(99, "Curvature to color..."); + float d = maxc-minc; + minc += 0.1*d; + maxc -= 0.1*d; + + vcg::Histogramf H; + vcg::tri::Stat::ComputePerVertexQualityHistogram(mesh->cm,H); + vcg::tri::UpdateColor::VertexQualityRamp(mesh->cm,H.Percentile(0.01f),H.Percentile(0.99f)); +// for (unsigned int i = 0; i< size; i++) +// { +// if ( (!selectionOnly) || (mesh->cm.vert[i].IsS()) ) +// { +// // float c = curvatures[i]; +// mesh->cm.vert[i].C().ColorRamp(minc,maxc, mesh->cm.vert[i].Q()); +// } +// } + } // else if (id & _AFRONT_) // { // // create a new mesh @@ -589,7 +523,6 @@ bool MlsPlugin::applyFilter(QAction* filter, MeshDocument& md, FilterParameterSe // } else if (id & _MCUBE_) { - using namespace vcg; // create a new mesh mesh = md.addNewMesh("mc_mesh"); @@ -628,10 +561,13 @@ bool MlsPlugin::applyFilter(QAction* filter, MeshDocument& md, FilterParameterSe delete pPoints; } - cb(99, "Update face normals..."); - vcg::tri::UpdateNormals::PerFace(mesh->cm); - cb(100, "Update box..."); - vcg::tri::UpdateBounding::Box(mesh->cm); + if (mesh) + { + cb(99, "Update face normals..."); + vcg::tri::UpdateNormals::PerFace(mesh->cm); + cb(100, "Update box..."); + vcg::tri::UpdateBounding::Box(mesh->cm); + } } // end MLS based stuff diff --git a/src/meshlabplugins/mlsfilters/mlsplugin.h b/src/meshlabplugins/mlsfilters/mlsplugin.h index 19755a95b..01793d4f8 100644 --- a/src/meshlabplugins/mlsfilters/mlsplugin.h +++ b/src/meshlabplugins/mlsfilters/mlsplugin.h @@ -37,20 +37,24 @@ class MlsPlugin : public QObject, public MeshFilterInterface public: enum { - _RIMLS_ = 0x1, - _APSS_ = 0x2, - _PROJECTION_ = 0x1000, - _AFRONT_ = 0x2000, - _MCUBE_ = 0x4000, + _RIMLS_ = 0x1, + _APSS_ = 0x2, + _PROJECTION_ = 0x1000, + _AFRONT_ = 0x2000, + _MCUBE_ = 0x4000, + _COLORIZE_ = 0x8000, FP_RIMLS_PROJECTION = _RIMLS_ | _PROJECTION_, FP_APSS_PROJECTION = _APSS_ | _PROJECTION_, FP_RIMLS_AFRONT = _RIMLS_ | _AFRONT_, - FP_APSS_AFRONT = _APSS_ | _AFRONT_, + FP_APSS_AFRONT = _APSS_ | _AFRONT_, - FP_RIMLS_MCUBE = _RIMLS_ | _MCUBE_, - FP_APSS_MCUBE = _APSS_ | _MCUBE_, + FP_RIMLS_MCUBE = _RIMLS_ | _MCUBE_, + FP_APSS_MCUBE = _APSS_ | _MCUBE_, + + FP_RIMLS_COLORIZE = _RIMLS_ | _COLORIZE_, + FP_APSS_COLORIZE = _APSS_ | _COLORIZE_, FP_RADIUS_FROM_DENSITY = 0x10000, FP_SELECT_SMALL_COMPONENTS = 0x20000 @@ -58,7 +62,6 @@ public: MlsPlugin(); - //virtual const QString filterName(FilterIDType filter) {return filterNames(filter);} virtual const QString filterName(FilterIDType filter); virtual const QString filterInfo(FilterIDType filter); const FilterClass getClass(QAction *a);