From 83d559efe63aafcaff2b6d57d45a11549293c9aa Mon Sep 17 00:00:00 2001 From: Paolo Cignoni cignoni Date: Tue, 9 Dec 2008 21:29:03 +0000 Subject: [PATCH] renamed mlsfilters to the standard naming scheme filter_mls --- src/meshlabplugins/mlsfilters/apss.cpp | 11 - src/meshlabplugins/mlsfilters/apss.h | 120 ---- src/meshlabplugins/mlsfilters/apss.tpp | 499 --------------- src/meshlabplugins/mlsfilters/balltree.cpp | 159 ----- src/meshlabplugins/mlsfilters/balltree.h | 119 ---- src/meshlabplugins/mlsfilters/implicits.h | 187 ------ src/meshlabplugins/mlsfilters/kdtree.cpp | 211 ------ src/meshlabplugins/mlsfilters/kdtree.h | 103 --- src/meshlabplugins/mlsfilters/mlsfilters.pro | 18 - .../mlsfilters/mlsmarchingcube.h | 280 -------- src/meshlabplugins/mlsfilters/mlsplugin.cpp | 605 ------------------ src/meshlabplugins/mlsfilters/mlsplugin.h | 76 --- src/meshlabplugins/mlsfilters/mlssurface.h | 203 ------ src/meshlabplugins/mlsfilters/mlssurface.tpp | 262 -------- src/meshlabplugins/mlsfilters/mlsutils.h | 140 ---- src/meshlabplugins/mlsfilters/priorityqueue.h | 122 ---- src/meshlabplugins/mlsfilters/rimls.cpp | 11 - src/meshlabplugins/mlsfilters/rimls.h | 104 --- src/meshlabplugins/mlsfilters/rimls.tpp | 307 --------- .../mlsfilters/smallcomponentselection.h | 150 ----- 20 files changed, 3687 deletions(-) delete mode 100644 src/meshlabplugins/mlsfilters/apss.cpp delete mode 100644 src/meshlabplugins/mlsfilters/apss.h delete mode 100644 src/meshlabplugins/mlsfilters/apss.tpp delete mode 100644 src/meshlabplugins/mlsfilters/balltree.cpp delete mode 100644 src/meshlabplugins/mlsfilters/balltree.h delete mode 100644 src/meshlabplugins/mlsfilters/implicits.h delete mode 100644 src/meshlabplugins/mlsfilters/kdtree.cpp delete mode 100644 src/meshlabplugins/mlsfilters/kdtree.h delete mode 100644 src/meshlabplugins/mlsfilters/mlsfilters.pro delete mode 100644 src/meshlabplugins/mlsfilters/mlsmarchingcube.h delete mode 100644 src/meshlabplugins/mlsfilters/mlsplugin.cpp delete mode 100644 src/meshlabplugins/mlsfilters/mlsplugin.h delete mode 100644 src/meshlabplugins/mlsfilters/mlssurface.h delete mode 100644 src/meshlabplugins/mlsfilters/mlssurface.tpp delete mode 100644 src/meshlabplugins/mlsfilters/mlsutils.h delete mode 100644 src/meshlabplugins/mlsfilters/priorityqueue.h delete mode 100644 src/meshlabplugins/mlsfilters/rimls.cpp delete mode 100644 src/meshlabplugins/mlsfilters/rimls.h delete mode 100644 src/meshlabplugins/mlsfilters/rimls.tpp delete mode 100644 src/meshlabplugins/mlsfilters/smallcomponentselection.h diff --git a/src/meshlabplugins/mlsfilters/apss.cpp b/src/meshlabplugins/mlsfilters/apss.cpp deleted file mode 100644 index e8c4fe816..000000000 --- a/src/meshlabplugins/mlsfilters/apss.cpp +++ /dev/null @@ -1,11 +0,0 @@ - -#include -#include "apss.h" -#include "apss.tpp" -#include - -namespace GaelMls { - -template class APSS; - -} diff --git a/src/meshlabplugins/mlsfilters/apss.h b/src/meshlabplugins/mlsfilters/apss.h deleted file mode 100644 index d4884fe09..000000000 --- a/src/meshlabplugins/mlsfilters/apss.h +++ /dev/null @@ -1,120 +0,0 @@ -/**************************************************************************** -* 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 APSS_H -#define APSS_H - -#include "mlssurface.h" - -namespace GaelMls { - -template -class APSS : public MlsSurface<_MeshType> -{ - typedef MlsSurface<_MeshType> Base; - - typedef typename Base::Scalar Scalar; - typedef typename Base::VectorType VectorType; - typedef typename Base::MatrixType MatrixType; - typedef _MeshType MeshType; - using Base::mCachedQueryPointIsOK; - using Base::mCachedQueryPoint; - using Base::mNeighborhood; - using Base::mCachedWeights; - using Base::mCachedWeightDerivatives; - using Base::mCachedWeightGradients; - using Base::mCachedWeightSecondDerivatives; - using Base::mBallTree; - using Base::mPoints; - using Base::mFilterScale; - using Base::mMaxNofProjectionIterations; - using Base::mAveragePointSpacing; - using Base::mProjectionAccuracy; - using Base::mGradientHint; - - enum Status {ASS_SPHERE, ASS_PLANE, ASS_UNDETERMINED}; - - public: - - APSS(const MeshType& m) - : Base(m) - { - mSphericalParameter = 1; - } - - 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 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: - bool fit(const VectorType& x) const; - bool mlsGradient(const VectorType& x, VectorType& grad) const; - bool mlsHessian(const VectorType& x, MatrixType& hessian) const; - - protected: - Scalar mSphericalParameter; - - // use double precision anyway - typedef double LScalar; - typedef vcg::Point3 LVector; - - // cached algebraic sphere coefficients - mutable LScalar uConstant; - mutable LVector uLinear; - mutable LScalar uQuad; - - mutable LVector mCenter; - mutable LScalar mRadius; - mutable Status mStatus; - - mutable LVector mCachedSumP; - mutable LVector mCachedSumN; - mutable LScalar mCachedSumDotPP; - mutable LScalar mCachedSumDotPN; - mutable LScalar mCachedSumW; - - mutable LVector mCachedGradSumP[3]; - mutable LVector mCachedGradSumN[3]; - mutable LScalar mCachedGradSumDotPN[3]; - mutable LScalar mCachedGradSumDotPP[3]; - mutable LScalar mCachedGradSumW[3]; - - mutable LScalar mCachedGradNume[3]; - mutable LScalar mCachedGradDeno[3]; - - mutable LScalar mCachedGradUConstant[3]; - mutable LVector mCachedGradULinear[3]; - mutable LScalar mCachedGradUQuad[3]; -}; - -} - -//#include "apss.tpp" - -#endif diff --git a/src/meshlabplugins/mlsfilters/apss.tpp b/src/meshlabplugins/mlsfilters/apss.tpp deleted file mode 100644 index b86267e50..000000000 --- a/src/meshlabplugins/mlsfilters/apss.tpp +++ /dev/null @@ -1,499 +0,0 @@ -/**************************************************************************** -* 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 "apss.h" -#include "kdtree.h" -#include "mlsutils.h" -#include - -namespace GaelMls { - -template -void APSS<_MeshType>::setSphericalParameter(Scalar v) -{ - mSphericalParameter = v; - mCachedQueryPointIsOK = false; -} - -template -typename APSS<_MeshType>::Scalar APSS<_MeshType>::potential(const VectorType& x, int* errorMask) const -{ - if ((!mCachedQueryPointIsOK) || mCachedQueryPoint!=x) - { - if (!fit(x)) - { - if (errorMask) - *errorMask = MLS_TOO_FAR; - return Base::InvalidValue(); - } - } - - LVector lx(x.X(), x.Y(), x.Z()); - - if (mStatus==ASS_SPHERE) - { - Scalar aux = vcg::Norm(lx - mCenter) - mRadius; - if (uQuad<0.) - aux = -aux; - return aux; - } - else if (mStatus==ASS_PLANE) - return vcg::Dot(lx,uLinear) + uConstant; - else - { - return uConstant + vcg::Dot(lx,uLinear) + uQuad * vcg::SquaredNorm(lx); - } -} - -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 -{ - if ((!mCachedQueryPointIsOK) || mCachedQueryPoint!=x) - { - if (!fit(x)) - { - if (errorMask) - *errorMask = MLS_TOO_FAR; - return VectorType(0,0,0); - } - } - - if (mGradientHint==MLS_DERIVATIVE_ACCURATE) - { - VectorType grad; - mlsGradient(x,grad); - return grad; - } - else - { - LVector lx(x.X(), x.Y(), x.Z()); - if (mStatus==ASS_PLANE) - return VectorType(uLinear.X(), uLinear.Y(), uLinear.Z()); - else - { - LVector g = uLinear + lx * (Scalar(2) * uQuad); - return VectorType(g.X(), g.Y(), g.Z()); - } - } -} - -template -typename APSS<_MeshType>::MatrixType APSS<_MeshType>::hessian(const VectorType& x, int* errorMask) const -{ - if ((!mCachedQueryPointIsOK) || mCachedQueryPoint!=x) - { - if (!fit(x)) - { - if (errorMask) - *errorMask = MLS_TOO_FAR; - return MatrixType(); - } - } - - MatrixType hessian; - if (Base::mHessianHint==MLS_DERIVATIVE_ACCURATE) - { - mlsHessian(x, hessian); - } - else - { - // this is very approximate !! - Scalar c = Scalar(2) * uQuad; - for (int i=0; i<3; ++i) - for (int j=0; j<3; ++j) - { - if (i==j) - hessian[i][j] = c; - else - hessian[i][j] = 0; - } - } - return hessian; -} - -template -typename APSS<_MeshType>::VectorType APSS<_MeshType>::project(const VectorType& x, VectorType* pNormal, int* errorMask) const -{ - int iterationCount = 0; - LVector lx(x.X(), x.Y(), x.Z()); - LVector position = lx; - LVector normal; - LVector previousPosition; - LScalar delta2; - LScalar epsilon2 = mAveragePointSpacing * mProjectionAccuracy; - epsilon2 = epsilon2 * epsilon2; - do { - if (!fit(VectorType(position.X(), position.Y(), position.Z()))) - { - if (errorMask) - *errorMask = MLS_TOO_FAR; - //std::cout << " proj failed\n"; - return x; - } - - previousPosition = position; - // local projection - if (mStatus==ASS_SPHERE) - { - normal = lx - mCenter; - normal.Normalize(); - position = mCenter + normal * mRadius; - - normal = uLinear + position * (LScalar(2) * uQuad); - normal.Normalize(); - } - else if (mStatus==ASS_PLANE) - { - normal = uLinear; - position = lx - uLinear * (vcg::Dot(lx,uLinear) + uConstant); - } - else - { - // Newton iterations - LVector grad; - LVector dir = uLinear+lx*(2.*uQuad); - LScalar ilg = 1./vcg::Norm(dir); - dir *= ilg; - LScalar ad = uConstant + vcg::Dot(uLinear,lx) + uQuad * vcg::SquaredNorm(lx); - LScalar delta = -ad*std::min(ilg,1.); - LVector p = lx + dir*delta; - for (int i=0 ; i<2 ; ++i) - { - grad = uLinear+p*(2.*uQuad); - ilg = 1./vcg::Norm(grad); - delta = -(uConstant + vcg::Dot(uLinear,p) + uQuad * vcg::SquaredNorm(p))*std::min(ilg,1.); - p += dir*delta; - } - position = p; - - normal = uLinear + position * (Scalar(2) * uQuad); - normal.Normalize(); - } - - delta2 = vcg::SquaredNorm(previousPosition - position); - } while ( delta2>epsilon2 && ++iterationCount(position), grad); - grad.Normalize(); - *pNormal = grad; - } - else - { - for (int k=0; k<3; ++k) - (*pNormal)[k] = normal[k]; - } - } - - if (iterationCount>=mMaxNofProjectionIterations && errorMask) - *errorMask = MLS_TOO_MANY_ITERS; - - return VectorType(position.X(), position.Y(), position.Z()); -} - -template -bool APSS<_MeshType>::fit(const VectorType& x) const -{ - Base::computeNeighborhood(x, true); - unsigned int nofSamples = mNeighborhood.size(); - - if (nofSamples==0) - { - mCachedQueryPointIsOK = false; - return false; - } - else if (nofSamples==1) - { - int id = mNeighborhood.index(0); - LVector p = vcg::Point3Cast(mPoints[id].cP()); - LVector n = vcg::Point3Cast(mPoints[id].cN()); - - uLinear = n; - uConstant = -vcg::Dot(p, uLinear); - uQuad = 0; - mStatus = ASS_PLANE; - return true; - } - - LVector sumP; sumP.SetZero(); - LVector sumN; sumN.SetZero(); - LScalar sumDotPN = 0.; - LScalar sumDotPP = 0.; - LScalar sumW = 0.; - for (unsigned int i=0; i(mPoints[id].cP()); - LVector n = vcg::Point3Cast(mPoints[id].cN()); - LScalar w = mCachedWeights.at(i); - - sumP += p * w; - sumN += n * w; - sumDotPN += w * vcg::Dot(n,p); - sumDotPP += w * vcg::SquaredNorm(p); - sumW += w; - } - - LScalar invSumW = Scalar(1)/sumW; - LScalar aux4 = mSphericalParameter * LScalar(0.5) * - (sumDotPN - invSumW*vcg::Dot(sumP,sumN)) - /(sumDotPP - invSumW*vcg::SquaredNorm(sumP)); - uLinear = (sumN-sumP*(Scalar(2)*aux4))*invSumW; - uConstant = -invSumW*(Dot(uLinear,sumP) + sumDotPP*aux4); - uQuad = aux4; - - // finalize - if (fabs(uQuad)>1e-7) - { - mStatus = ASS_SPHERE; - LScalar b = 1./uQuad; - mCenter = uLinear*(-0.5*b); - mRadius = sqrt( vcg::SquaredNorm(mCenter) - b*uConstant ); - } - else if (uQuad==0.) - { - mStatus = ASS_PLANE; - LScalar s = LScalar(1)/vcg::Norm(uLinear); - uLinear *= s; - uConstant *= s; - } - else - { - mStatus = ASS_UNDETERMINED; - // normalize the gradient - LScalar f = 1./sqrt(vcg::SquaredNorm(uLinear) - Scalar(4)*uConstant*uQuad); - uConstant *= f; - uLinear *= f; - uQuad *= f; - } - - // cache some values to be used by the mls gradient - mCachedSumP = sumP; - mCachedSumN = sumN; - mCachedSumW = sumW; - mCachedSumDotPP = sumDotPP; - mCachedSumDotPN = sumDotPN; - - mCachedQueryPoint = x; - mCachedQueryPointIsOK = true; - return true; - } - - template - bool APSS<_MeshType>::mlsGradient(const VectorType& x, VectorType& grad) const - { - unsigned int nofSamples = mNeighborhood.size(); - - const LVector& sumP = mCachedSumP; - const LVector& sumN = mCachedSumN; - const LScalar& sumDotPN = mCachedSumDotPN; - const LScalar& sumDotPP = mCachedSumDotPP; - const LScalar& sumW = mCachedSumW; - const LScalar invSumW = 1.f/sumW; - - const LScalar nume = sumDotPN - invSumW * vcg::Dot(sumP, sumN); - const LScalar deno = sumDotPP - invSumW * vcg::Dot(sumP, sumP); - - for (uint k=0 ; k<3 ; ++k) - { - LVector dSumP; dSumP.SetZero(); - LVector dSumN; dSumN.SetZero(); - LScalar dSumDotPN = 0.; - LScalar dSumDotPP = 0.; - LScalar dSumW = 0.; - for (unsigned int i=0; i(mPoints[id].cP()); - LVector n = vcg::Point3Cast(mPoints[id].cN()); - LScalar dw = mCachedWeightGradients.at(i)[k]; - - dSumW += dw; - dSumP += p*dw; - dSumN += n*dw; - dSumDotPN += dw * vcg::Dot(n,p); - dSumDotPP += dw * vcg::SquaredNorm(p); - } - - mCachedGradSumP[k] = dSumP; - mCachedGradSumN[k] = dSumN; - mCachedGradSumDotPN[k] = dSumDotPN; - mCachedGradSumDotPP[k] = dSumDotPP; - mCachedGradSumW[k] = dSumW; - - LScalar dVecU0; - LVector dVecU13; - LScalar dVecU4; - - LScalar dNume = dSumDotPN - invSumW*invSumW*( sumW*(vcg::Dot(dSumP,sumN) + vcg::Dot(sumP,dSumN)) - - dSumW*vcg::Dot(sumP,sumN)); - LScalar dDeno = dSumDotPP - invSumW*invSumW*( 2.*sumW*vcg::Dot(dSumP,sumP) - - dSumW*vcg::Dot(sumP,sumP)); - - dVecU4 = mSphericalParameter * 0.5 * (deno * dNume - dDeno * nume)/(deno*deno); - dVecU13 = ((dSumN - (dSumP*uQuad + sumP*dVecU4)*2.0) - uLinear * dSumW) * invSumW; - dVecU0 = -invSumW*( vcg::Dot(dVecU13,sumP) + dVecU4*sumDotPP + vcg::Dot(uLinear,dSumP) + uQuad*dSumDotPP + dSumW*uConstant); - - grad[k] = dVecU0 + vcg::Dot(dVecU13,vcg::Point3Cast(x)) + dVecU4*vcg::SquaredNorm(x) + uLinear[k] + 2.*x[k]*uQuad; - - mCachedGradDeno[k] = dDeno; - mCachedGradNume[k] = dNume; - mCachedGradUConstant[k] = dVecU0; - mCachedGradULinear[k] = dVecU13; - mCachedGradUQuad[k] = dVecU4; - } - - return true; -} - -template -bool APSS<_MeshType>::mlsHessian(const VectorType& x, MatrixType& hessian) const -{ - this->requestSecondDerivatives(); - - // TODO call mlsGradient first - VectorType grad; - mlsGradient(x,grad); - - uint nofSamples = mNeighborhood.size(); - - const LVector& sumP = mCachedSumP; - const LVector& sumN = mCachedSumN; - const LScalar& sumDotPN = mCachedSumDotPN; - const LScalar& sumDotPP = mCachedSumDotPP; - const LScalar& sumW = mCachedSumW; - const LScalar invSumW = 1.f/sumW; - - const LScalar nume = sumDotPN - invSumW * vcg::Dot(sumP, sumN); - const LScalar deno = sumDotPP - invSumW * vcg::Dot(sumP, sumP); - - for (uint k=0 ; k<3 ; ++k) - { - const LVector& dSumP = mCachedGradSumP[k]; - const LVector& dSumN = mCachedGradSumN[k]; - const LScalar& dSumDotPN = mCachedGradSumDotPN[k]; - const LScalar& dSumDotPP = mCachedGradSumDotPP[k]; - const LScalar& dSumW = mCachedGradSumW[k]; - - LScalar dVecU0 = mCachedGradUConstant[k]; - LVector dVecU13 = mCachedGradULinear[k]; - LScalar dVecU4 = mCachedGradUQuad[k]; - - LScalar dNume = mCachedGradNume[k]; - LScalar dDeno = mCachedGradDeno[k]; - - // second order derivatives - for (uint j=0 ; j<3 ; ++j) - { - LVector d2SumP; d2SumP.SetZero(); - LVector d2SumN; d2SumN.SetZero(); - LScalar d2SumDotPN = 0.; - LScalar d2SumDotPP = 0.; - LScalar d2SumW = 0.; - for (unsigned int i=0; i(mPoints[id].cP()); - LVector n = vcg::Point3Cast(mPoints[id].cN()); - LScalar dw = mCachedWeightGradients.at(i)[j]; - LScalar d2w = ((x[k]-p[k]))*((x[j]-p[j])) * mCachedWeightSecondDerivatives.at(i); - - if (j==k) - d2w += mCachedWeightDerivatives.at(i); - - d2SumW += d2w; - d2SumP += p*d2w; - d2SumN += n*d2w; - d2SumDotPN += d2w * vcg::Dot(n,p); - d2SumDotPP += d2w * vcg::Dot(p,p); - } - - LScalar d2u0; - LVector d2u13; - LScalar d2u4; - - LScalar d2Nume = d2SumDotPN - invSumW*invSumW*invSumW*invSumW*( - - 2.*sumW*mCachedGradSumW[j]*( sumW*(vcg::Dot(dSumP,sumN)+vcg::Dot(sumP,dSumN)) - dSumW* vcg::Dot(sumP,sumN)) - + sumW*sumW*( mCachedGradSumW[j]*(vcg::Dot(dSumP,sumN)+vcg::Dot(sumP,dSumN)) - + sumW*(vcg::Dot(d2SumP,sumN) + vcg::Dot(sumP,d2SumN) - + vcg::Dot(mCachedGradSumP[j],dSumN) + vcg::Dot(dSumP,mCachedGradSumN[j])) - - d2SumW*vcg::Dot(sumP,sumN) - - dSumW*(vcg::Dot(mCachedGradSumP[j],sumN)+vcg::Dot(sumP,mCachedGradSumN[j])) )); - - LScalar d2Deno = d2SumDotPP - invSumW*invSumW*invSumW*invSumW*( - - 2.*sumW*mCachedGradSumW[j] * ( 2.*sumW*vcg::Dot(dSumP,sumP) - dSumW*vcg::Dot(sumP,sumP)) - + sumW*sumW*( 2.*mCachedGradSumW[j]*(vcg::Dot(dSumP,sumP)) - + 2.*sumW*(vcg::Dot(mCachedGradSumP[j],dSumP)+vcg::Dot(d2SumP,sumP)) - - d2SumW*vcg::Dot(sumP,sumP) - dSumW*(2.*vcg::Dot(mCachedGradSumP[j],sumP))) ); - - LScalar deno2 = deno*deno; - d2u4 = mSphericalParameter * 0.5 * (deno2*(d2Nume*deno + mCachedGradDeno[j] * dNume - - d2Deno*nume - dDeno * mCachedGradNume[j]) - - 2.*deno*mCachedGradDeno[j]*(deno * dNume - dDeno * nume))/(deno2*deno2); - - d2u13 = ( -dVecU13 * mCachedGradSumW[j] - + (d2SumN - (dSumP*mCachedGradUQuad[j] + d2SumP*uQuad + sumP*d2u4 + mCachedGradSumP[j]*dVecU4)*2.0 ) - - uLinear*d2SumW - mCachedGradULinear[j]*dSumW ) * invSumW; - - d2u0 = ( -dVecU0 * mCachedGradSumW[j] - - ( vcg::Dot(dVecU13,mCachedGradSumP[j]) + vcg::Dot(d2u13,sumP) - + d2u4*sumDotPP + dVecU4*mCachedGradSumDotPP[j] - + vcg::Dot(uLinear,d2SumP) + vcg::Dot(mCachedGradULinear[j],dSumP) - + dSumDotPP*mCachedGradUQuad[j] + d2SumDotPP*uQuad - + d2SumW*uConstant + dSumW*mCachedGradUConstant[j]) ) * invSumW; - - hessian[j][k] = - dVecU13[j] + 2.*dVecU4*x[j] - + d2u0 + vcg::Dot(d2u13,vcg::Point3Cast(x)) + d2u4*vcg::Dot(x,x) - + mCachedGradULinear[j][k] + (j==k ? 2.*uQuad : 0.) + 2.*x[k]*mCachedGradUQuad[j]; - - } - } - - return true; -} - -} diff --git a/src/meshlabplugins/mlsfilters/balltree.cpp b/src/meshlabplugins/mlsfilters/balltree.cpp deleted file mode 100644 index 7ebf79677..000000000 --- a/src/meshlabplugins/mlsfilters/balltree.cpp +++ /dev/null @@ -1,159 +0,0 @@ -/**************************************************************************** -* 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 "balltree.h" - -namespace GaelMls { - -template -BallTree<_Scalar>::BallTree(const ConstDataWrapper& points, const ConstDataWrapper& radii) - : mPoints(points), mRadii(radii), mRadiusScale(1.), mTreeIsUptodate(false) -{ - mRootNode = 0; - mMaxTreeDepth = 12; - mTargetCellSize = 24; -} - -template -void BallTree<_Scalar>::computeNeighbors(const VectorType& x, Neighborhood* pNei) const -{ - if (!mTreeIsUptodate) - const_cast(this)->rebuild(); - - pNei->clear(); - mQueryPosition = x; - queryNode(*mRootNode, pNei); -} - -template -void BallTree<_Scalar>::queryNode(Node& node, Neighborhood* pNei) const -{ - if (node.leaf) - { - for (unsigned int i=0 ; iinsert(id, d2); - } - } - else - { - if (mQueryPosition[node.dim] - node.splitValue < 0) - queryNode(*node.children[0], pNei); - else - queryNode(*node.children[1], pNei); - } -} - -template -void BallTree<_Scalar>::rebuild(void) -{ - delete mRootNode; - - mRootNode = new Node(); - IndexArray indices(mPoints.size()); - AxisAlignedBoxType aabb; - aabb.Set(mPoints[0]); - for (unsigned int i=0 ; i -void BallTree<_Scalar>::split(const IndexArray& indices, const AxisAlignedBoxType& aabbLeft, const AxisAlignedBoxType& aabbRight, IndexArray& iLeft, IndexArray& iRight) -{ - for (std::vector::const_iterator it=indices.begin(), end=indices.end() ; it!=end ; ++it) - { - unsigned int i = *it; - if (vcg::Distance(mPoints[i], aabbLeft) < mRadii[i]*mRadiusScale) - iLeft.push_back(i); - - if (vcg::Distance(mPoints[i], aabbRight) < mRadii[i]*mRadiusScale) - iRight.push_back(i); - } -} - -template -void BallTree<_Scalar>::buildNode(Node& node, std::vector& indices, AxisAlignedBoxType aabb, int level) -{ - Scalar avgradius = 0.; - for (std::vector::const_iterator it=indices.begin(), end=indices.end() ; it!=end ; ++it) - avgradius += mRadii[*it]; - avgradius = mRadiusScale * avgradius / Scalar(indices.size()); - VectorType diag = aabb.max - aabb.min; - if (int(indices.size()) std::max(std::max(diag.X(), diag.Y()), diag.Z()) - || int(level)>=mMaxTreeDepth) - { - node.leaf = true; - node.size = indices.size(); - node.indices = new unsigned int[node.size]; - for (unsigned int i=0 ; i iLeft, iRight; - split(indices, aabbLeft, aabbRight, iLeft,iRight); - - // we don't need the index list anymore - indices.clear(); - - { - // left child - //mNodes.resize(mNodes.size()+1); - Node* pChild = new Node(); - node.children[0] = pChild; - buildNode(*pChild, iLeft, aabbLeft, level+1); - } - - { - // right child - //mNodes.resize(mNodes.size()+1); - Node* pChild = new Node(); - node.children[1] = pChild; - buildNode(*pChild, iRight, aabbRight, level+1); - } -} - -template class BallTree; -template class BallTree; - -} diff --git a/src/meshlabplugins/mlsfilters/balltree.h b/src/meshlabplugins/mlsfilters/balltree.h deleted file mode 100644 index d3c20424b..000000000 --- a/src/meshlabplugins/mlsfilters/balltree.h +++ /dev/null @@ -1,119 +0,0 @@ -/**************************************************************************** -* 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 BALLTREE_H -#define BALLTREE_H - -#include -#include -#include "mlsutils.h" - -namespace GaelMls { - -template -class Neighborhood -{ - public: - typedef _Scalar Scalar; - - int index(int i) const { return mIndices.at(i); } - Scalar squaredDistance(int i) const { return mSqDists.at(i); } - - void clear() { mIndices.clear(); mSqDists.clear(); } - void resize(int size) { mIndices.resize(size); mSqDists.resize(size); } - void reserve(int size) { mIndices.reserve(size); mSqDists.reserve(size); } - int size() { return mIndices.size(); } - - void insert(int id, Scalar d2) { mIndices.push_back(id); mSqDists.push_back(d2); } - - protected: - std::vector mIndices; - std::vector mSqDists; -}; - -template -class BallTree -{ - public: - typedef _Scalar Scalar; - typedef vcg::Point3 VectorType; - - BallTree(const ConstDataWrapper& points, const ConstDataWrapper& radii); - - void computeNeighbors(const VectorType& x, Neighborhood* pNei) const; - - void setRadiusScale(Scalar v) { mRadiusScale = v; mTreeIsUptodate = false; } - - protected: - - struct Node - { - ~Node() - { - if (!leaf) - { - delete children[0]; - delete children[1]; - } - else - { - delete[] indices; - } - } - Scalar splitValue; - unsigned char dim:2; - unsigned char leaf:1; - union { - Node* children[2]; - struct { - unsigned int* indices; - unsigned int size; - }; - }; - }; - - typedef std::vector IndexArray; - typedef vcg::Box3 AxisAlignedBoxType; - - void rebuild(); - void split(const IndexArray& indices, const AxisAlignedBoxType& aabbLeft, const AxisAlignedBoxType& aabbRight, - IndexArray& iLeft, IndexArray& iRight); - void buildNode(Node& node, std::vector& indices, AxisAlignedBoxType aabb, int level); - void queryNode(Node& node, Neighborhood* pNei) const; - - protected: - ConstDataWrapper mPoints; - ConstDataWrapper mRadii; - Scalar mRadiusScale; - - int mMaxTreeDepth; - int mTargetCellSize; - mutable bool mTreeIsUptodate; - mutable VectorType mQueryPosition; - - Node* mRootNode; -}; - -} - -#endif diff --git a/src/meshlabplugins/mlsfilters/implicits.h b/src/meshlabplugins/mlsfilters/implicits.h deleted file mode 100644 index bd09ac36c..000000000 --- a/src/meshlabplugins/mlsfilters/implicits.h +++ /dev/null @@ -1,187 +0,0 @@ -/**************************************************************************** -* 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 __IMPLICITS -#define __IMPLICITS - -#include -#include -#include - -namespace vcg -{ -namespace implicits -{ - -/** \returns the Gauss curvature directly from the gradient and hessian */ -template -Scalar GaussCurvature(const Point3& gradient, const Matrix33& hessian) -{ - Scalar l2 = gradient.SquaredNorm(); - Matrix33 adjugate; - adjugate[0][0] = hessian[1][1]*hessian[2][2] - hessian[1][2]*hessian[2][1]; - adjugate[1][0] = hessian[0][2]*hessian[2][1] - hessian[0][1]*hessian[2][2]; - adjugate[2][0] = hessian[0][1]*hessian[1][2] - hessian[0][2]*hessian[1][1]; - - adjugate[0][1] = hessian[1][2]*hessian[2][0] - hessian[1][0]*hessian[2][2]; - adjugate[1][1] = hessian[0][0]*hessian[2][2] - hessian[0][2]*hessian[2][0]; - adjugate[2][1] = hessian[1][0]*hessian[0][2] - hessian[0][0]*hessian[1][2]; - - adjugate[0][2] = hessian[1][0]*hessian[2][1] - hessian[1][1]*hessian[2][0]; - adjugate[1][2] = hessian[0][1]*hessian[2][0] - hessian[0][0]*hessian[2][1]; - adjugate[2][2] = hessian[0][0]*hessian[1][1] - hessian[0][1]*hessian[1][0]; - return (vcg::Dot(gradient, adjugate * gradient)) / (l2*l2); -} - -/** \returns the mean curvature directly from the gradient and hessian */ -template -Scalar MeanCurvature(const Point3& gradient, const Matrix33& hessian) -{ - Scalar l = gradient.Norm(); - return (l*l*hessian.Trace() - vcg::Dot(gradient, hessian * gradient)) / (2.*l*l*l); -} - -/** This class computes the Weingarten map of a scalar field and provides - * methods to extract curvatures from it. - * - * The Weingarten map is equal to gradient of the normal vector: - * \f$ W = \nabla \mathbf(n) - * = \nabla \frac{\mathbf{g}}{\Vert \mathbf{g} \Vert} - * = \frac{(I - n n^T) H}{\Vert \mathbf{g} \Vert} \f$ - * This matrix can also be seen as the projection of the hessian - * matrix onto the tangent plane of normal n. - */ -template class WeingartenMap -{ -public: - typedef Point3 VectorType; - typedef Matrix33 MatrixType; - - /** Default constructor computing the Weingarten map from the - * first and second derivatives, namely the gradient vector - * and hessian matrix of the scalar field. - */ - WeingartenMap(const VectorType& grad, const MatrixType& hess) - { - Scalar invL = 1.0/grad.Norm(); - m_normal = grad * invL; - - Matrix33 I; I.SetIdentity(); - m_nnT.ExternalProduct(m_normal,m_normal); - - m_w = (I-m_nnT) * hess * invL; - - m_kgIsDirty = true; - m_kmIsDirty = true; - m_kpAreDirty = true; - m_kdirAreDirty = true; - } - - /** \returns the Weingarten map matrix */ - const MatrixType& W() const { return m_w; } - - /** \returns the Gauss curvature = k1 * k2 */ - Scalar GaussCurvature() const - { - if (m_kgIsDirty) - { - // we add nn^T to W such that the third eigenvalue becomes 1 - // then det(W) = k1 * k2 * 1 = Gauss curvature ! - m_kg = (m_w + m_nnT).Determinant(); - m_kgIsDirty = false; - } - return m_kg; - } - - /** \returns the mean curvature = (k1 + k2)/2 */ - Scalar MeanCurvature() const - { - if (m_kmIsDirty) - { - // the third eigenvalue of W is 0, then tr(W) = k1 + k2 + 0 = 2 k mean ! - m_km = m_w.Trace(); - m_kmIsDirty = false; - } - return m_km; - } - - /** \returns the first principal curvature */ - Scalar K1() const { updateKp(); return m_k1; } - - /** \returns the second principal curvature */ - Scalar K2() const { updateKp(); return m_k2; } - - /** \returns the direction of the first principal curvature */ - const VectorType& K1Dir() const { extractEigenvectors(); return m_k1dir; } - - /** \returns the direction of the second principal curvature */ - const VectorType& K2Dir() const { extractEigenvectors(); return m_k2dir; } - -protected: - - // direct computation of principal curvatures if needed - inline void updateKp() const - { - if (m_kpAreDirty) - { - Scalar delta = sqrt(MeanCurvature()*m_km - 4.0*GaussCurvature()); - m_k1 = 0.5*(m_km + delta); - m_k2 = 0.5*(m_km - delta); - if (fabs(m_k1) -#include - -template -KdTree::KdTree(const ConstDataWrapper& points, unsigned int nofPointsPerCell, unsigned int maxDepth) - : mPoints(points.size()) -{ - // compute the AABB of the input - mPoints[0] = points[0]; - mAABB.Set(mPoints[0]); - for (unsigned int i=1 ; i -KdTree::~KdTree() -{ -} - -template -void KdTree::setMaxNofNeighbors(unsigned int k) -{ - mNeighborQueue.setMaxSize(k); -} - -/** Performs the kNN query. - * - * This algorithm uses the simple distance to the split plane to prune nodes. - * A more elaborated approach consists to track the closest corner of the cell - * relatively to the current query point. This strategy allows to save about 5% - * of the leaves. However, in practice the slight overhead due to this tracking - * reduces the overall performance. - * - * This algorithm also use a simple stack while a priority queue using the squared - * distances to the cells as a priority values allows to save about 10% of the leaves. - * But, again, priority queue insertions and deletions are quite involved, and therefore - * a simple stack is by far much faster. - */ -template -void KdTree::doQueryK(const VectorType& queryPoint) -{ - mNeighborQueue.init(); - mNeighborQueue.insert(0xffffffff, std::numeric_limits::max()); - - mNodeStack[0].nodeId = 0; - mNodeStack[0].sq = 0.f; - unsigned int count = 1; - - while (count) - { - QueryNode& qnode = mNodeStack[count-1]; - Node& node = mNodes[qnode.nodeId]; - - if (qnode.sq < mNeighborQueue.getTopWeight()) - { - if (node.leaf) - { - --count; // pop - unsigned int end = node.start+node.size; - for (unsigned int i=node.start ; i -unsigned int KdTree::split(int start, int end, unsigned int dim, float splitValue) -{ - int l(start), r(end-1); - for ( ; l= start && mPoints[r][dim] >= splitValue) - r--; - if (l > r) - break; - std::swap(mPoints[l],mPoints[r]); - } - return (mPoints[l][dim] < splitValue ? l+1 : l); -} - -/** recursively builds the kdtree - * - * The heuristic is the following: - * - if the number of points in the node is lower than targetCellsize then make a leaf - * - else compute the AABB of the points of the node and split it at the middle of - * the largest AABB dimension. - * - * This strategy might look not optimal because it does not explicitly prune empty space, - * unlike more advanced SAH-like techniques used for RT. On the other hand it leads to a shorter tree, - * faster to traverse and our experience shown that in the special case of kNN queries, - * this strategy is indeed more efficient (and much faster to build). Moreover, for volume data - * (e.g., fluid simulation) pruning the empty space is useless. - * - * Actually, storing at each node the exact AABB (we therefore have a binary BVH) allows - * to prune only about 10% of the leaves, but the overhead of this pruning (ball/ABBB intersection) - * is more expensive than the gain it provides and the memory consumption is x4 higher ! - */ -template -void KdTree::createTree(unsigned int nodeId, unsigned int start, unsigned int end, unsigned int level, unsigned int targetCellSize, unsigned int targetMaxDepth) -{ - Node& node = mNodes[nodeId]; - AxisAlignedBoxType aabb; - aabb.Set(mPoints[start]); - for (unsigned int i=start+1 ; i=targetMaxDepth) - { - child.leaf = 1; - child.start = start; - child.size = midId-start; - } - else - { - child.leaf = 0; - createTree(childId, start, midId, level+1, targetCellSize, targetMaxDepth); - } - } - - { - // right child - unsigned int childId = mNodes[nodeId].firstChildId+1; - Node& child = mNodes[childId]; - if (end-midId <= targetCellSize || level>=targetMaxDepth) - { - child.leaf = 1; - child.start = midId; - child.size = end-midId; - } - else - { - child.leaf = 0; - createTree(childId, midId, end, level+1, targetCellSize, targetMaxDepth); - } - } -} - -template class KdTree; -template class KdTree; diff --git a/src/meshlabplugins/mlsfilters/kdtree.h b/src/meshlabplugins/mlsfilters/kdtree.h deleted file mode 100644 index 4cca5a845..000000000 --- a/src/meshlabplugins/mlsfilters/kdtree.h +++ /dev/null @@ -1,103 +0,0 @@ -/**************************************************************************** -* 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 KDTREE_H -#define KDTREE_H - -#include -#include -#include "mlsutils.h" -#include "priorityqueue.h" -#include - -template -class KdTree -{ -public: - - typedef _Scalar Scalar; - typedef vcg::Point3 VectorType; - typedef vcg::Box3 AxisAlignedBoxType; - - struct Node - { - union { - struct { - Scalar splitValue; - unsigned int firstChildId:24; - unsigned int dim:2; - unsigned int leaf:1; - }; - struct { - unsigned int start; - unsigned short size; - }; - }; - }; - typedef std::vector NodeList; - - inline const NodeList& _getNodes(void) { return mNodes; } - inline const std::vector& _getPoints(void) { return mPoints; } - - void setMaxNofNeighbors(unsigned int k); - inline int getNofFoundNeighbors(void) { return mNeighborQueue.getNofElements(); } - inline const VectorType& getNeighbor(int i) { return mPoints[ mNeighborQueue.getIndex(i) ]; } - inline unsigned int getNeighborId(int i) { return mNeighborQueue.getIndex(i); } - inline float getNeighborSquaredDistance(int i) { return mNeighborQueue.getWeight(i); } - -public: - - KdTree(const ConstDataWrapper& points, unsigned int nofPointsPerCell = 16, unsigned int maxDepth = 64); - - ~KdTree(); - - void doQueryK(const VectorType& p); - -protected: - - // element of the stack - struct QueryNode - { - QueryNode() {} - QueryNode(unsigned int id) : nodeId(id) {} - unsigned int nodeId; // id of the next node - Scalar sq; // squared distance to the next node - }; - - // used to build the tree: split the subset [start..end[ according to dim and splitValue, - // and returns the index of the first element of the second subset - unsigned int split(int start, int end, unsigned int dim, float splitValue); - - void createTree(unsigned int nodeId, unsigned int start, unsigned int end, unsigned int level, unsigned int targetCellsize, unsigned int targetMaxDepth); - -protected: - - AxisAlignedBoxType mAABB; - NodeList mNodes; - std::vector mPoints; - HeapMaxPriorityQueue mNeighborQueue; - QueryNode mNodeStack[64]; -}; - -#endif - diff --git a/src/meshlabplugins/mlsfilters/mlsfilters.pro b/src/meshlabplugins/mlsfilters/mlsfilters.pro deleted file mode 100644 index dba494ca8..000000000 --- a/src/meshlabplugins/mlsfilters/mlsfilters.pro +++ /dev/null @@ -1,18 +0,0 @@ -include(../../shared.pri) - -# CONFIG += debug - -HEADERS = rimls.h mlsplugin.h \ - ../../meshlab/meshmodel.h - -SOURCES = balltree.cpp kdtree.cpp mlsplugin.cpp \ - apss.cpp rimls.cpp \ - ../../meshlab/filterparameter.cpp \ - ../../meshlab/meshmodel.cpp \ - $$GLEWCODE - -TARGET = mlsfilters - -CONFIG += opengl -CONFIG += warn_off - diff --git a/src/meshlabplugins/mlsfilters/mlsmarchingcube.h b/src/meshlabplugins/mlsfilters/mlsmarchingcube.h deleted file mode 100644 index 2f97b2037..000000000 --- a/src/meshlabplugins/mlsfilters/mlsmarchingcube.h +++ /dev/null @@ -1,280 +0,0 @@ -/**************************************************************************** -* 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 MLSMARCHINGCUBES_H -#define MLSMARCHINGCUBES_H - -#include -#include -#include -#include "mlssurface.h" -#include "mlsutils.h" - -namespace vcg { -namespace tri { - -template -class MlsWalker -{ -private: - typedef int VertexIndex; - typedef typename MeshType::ScalarType ScalarType; - typedef vcg::Point3f VectorType; - typedef typename MeshType::VertexPointer VertexPointer; - //typedef GaelMls::MlsSurface SurfaceType; - typedef long long unsigned int Key; - typedef std::map MapType; - - struct GridElement - { - VectorType position; - ScalarType value; - }; - - template - inline bool IsFinite(T value) - { - return (value>=-std::numeric_limits::max()) && (value<=std::numeric_limits::max()); - } - -public: - - int resolution; - - MlsWalker() - { - resolution = 150; - mMaxBlockSize = 64; - mIsoValue = 0; - } - - template - void BuildMesh(MeshType &mesh, SurfaceType &surface, EXTRACTOR_TYPE &extractor, CallBackPos *cb = 0) - { - mpSurface = &surface; - const ScalarType invalidValue = SurfaceType::InvalidValue(); - - // precomputed offsets to access the cell corners - const int offsets[8] = { - 0, - 1, - 1+mMaxBlockSize*mMaxBlockSize, - mMaxBlockSize*mMaxBlockSize, - mMaxBlockSize, - 1+mMaxBlockSize, - 1+mMaxBlockSize+mMaxBlockSize*mMaxBlockSize, - mMaxBlockSize+mMaxBlockSize*mMaxBlockSize - }; - - mAABB = mpSurface->boundingBox(); - VectorType diag = mAABB.max - mAABB.min; - mAABB.min -= diag * 0.1; - mAABB.max += diag * 0.1; - diag = mAABB.max - mAABB.min; - - if ( (diag[0]<=0.) - || (diag[1]<=0.) - || (diag[2]<=0.) - || (resolution==0)) - { - return; - } - - mCache = new GridElement[(mMaxBlockSize)*(mMaxBlockSize)*(mMaxBlockSize)]; - ScalarType step = vcg::MaxCoeff(diag)/ScalarType(resolution); - - unsigned int nofCells[3]; - unsigned int nofBlocks[3]; - - for (uint k=0 ; k<3 ; ++k) - { - nofCells[k] = int(diag[k]/step)+2; - nofBlocks[k] = nofCells[k]/mMaxBlockSize + ( (nofCells[k]%mMaxBlockSize)==0 ? 0 : 1); - } - - _mesh = &mesh; - _mesh->Clear(); - - int countSubSlice = 0; - int totalSubSlices = nofBlocks[2] * nofBlocks[1] * nofCells[0]; - - extractor.Initialize(); - // for each macro block - vcg::Point3i bi; // block id - for (bi[2]=0 ; bi[2](mMaxBlockSize, nofCells[k]-(mMaxBlockSize-1)*bi[k]); - } - VectorType origin = mAABB.min + VectorType(bi[0],bi[1],bi[2]) * (step * (mMaxBlockSize-1)); - - // fill the grid - vcg::Point3i ci; // local cell id - - // for each corners... - for (ci[0]=0 ; ci[0]potential(el.position); - if (!mpSurface->isInDomain(el.position)) - el.value = invalidValue; - } - } - - // polygonize the grid (marching cube) - // for each cell... - for (ci[0]=0 ; ci[0]id2) - std::swap(id1,id2); - Key k = Key(id1) + (Key(id2)<<32); - MapType::iterator it = mVertexMap.find(k); - if (it!=mVertexMap.end()) - { - // a vertex already exist - v = &_mesh->vert[it->second]; - } - else if (create) - { - // let's create a new vertex - VertexIndex vi = (VertexIndex) _mesh->vert.size(); - Allocator::AddVertices( *_mesh, 1 ); - mVertexMap[k] = vi; - v = &_mesh->vert[vi]; - id1 = GetLocalCellIdFromGlobal(p1); - id2 = GetLocalCellIdFromGlobal(p2); - // interpol along the edge - ScalarType epsilon = 1e-5; - const GridElement& c1 = mCache[id1]; - const GridElement& c2 = mCache[id2]; - if (fabs(mIsoValue-c1.value) < epsilon) - v->P() = c1.position; - else if (fabs(mIsoValue-c2.value) < epsilon) - v->P() = c2.position; - else if (fabs(c1.value-c2.value) < epsilon) - v->P() = (c1.position+c1.position)*0.5; - else - { - ScalarType a = (mIsoValue - c1.value) / (c2.value - c1.value); - v->P() = c1.position + (c2.position - c1.position) * a; - } - } - else - { - v = 0; - } - } - - bool Exist(const vcg::Point3i &p0, const vcg::Point3i &p1, VertexPointer &v) - { - GetIntercept(p0, p1, v, false); - return v!=0; - } - - void GetXIntercept(const vcg::Point3i &p1, const vcg::Point3i &p2, VertexPointer &v) - { - GetIntercept(p1, p2, v, true); - } - void GetYIntercept(const vcg::Point3i &p1, const vcg::Point3i &p2, VertexPointer &v) - { - GetIntercept(p1, p2, v, true); - } - void GetZIntercept(const vcg::Point3i &p1, const vcg::Point3i &p2, VertexPointer &v) - { - GetIntercept(p1, p2, v, true); - } - -protected: - Box3f mAABB; - MapType mVertexMap; - MeshType *_mesh; - SurfaceType* mpSurface; - GridElement* mCache; - vcg::Point3i mBlockOrigin; - vcg::Point3i mGridSize; - int mMaxBlockSize; - ScalarType mIsoValue; - -}; -} // end namespace -} // end namespace - -#endif - diff --git a/src/meshlabplugins/mlsfilters/mlsplugin.cpp b/src/meshlabplugins/mlsfilters/mlsplugin.cpp deleted file mode 100644 index cd2f51173..000000000 --- a/src/meshlabplugins/mlsfilters/mlsplugin.cpp +++ /dev/null @@ -1,605 +0,0 @@ -/**************************************************************************** -* 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 - -#define NDEBUG - -#include -#include -#include -#include - -#include -#include - -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include "mlsmarchingcube.h" - -#include "mlsplugin.h" -#include "rimls.h" -#include "apss.h" -#include "mlsutils.h" -#include "implicits.h" -#include "../meshfilter/refine_loop.h" - -#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 -// - actionList with the corresponding actions. If you want to add icons to your filtering actions you can do here by construction the QActions accordingly - -enum { - CT_MEAN = 0, - CT_GAUSS = 1, - CT_K1 = 2, - CT_K2 = 3, - CT_APSS = 4 -}; - -MlsPlugin::MlsPlugin() -{ - typeList - << 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; - -// initFilterList(this); - foreach(FilterIDType tt , types()) - actionList << new QAction(filterName(tt), this); -} - -// ST() must return the very short string describing each filtering action -// (this string is used also to define the menu entry) -const QString MlsPlugin::filterName(FilterIDType filterId) -{ - switch(filterId) { - case FP_APSS_PROJECTION : return QString("MLS projection (APSS)"); - case FP_RIMLS_PROJECTION : return QString("MLS projection (#####)"); - case FP_APSS_AFRONT : return QString("MLS meshing/APSS Advancing Front"); - 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); - } -} - -const MeshFilterInterface::FilterClass MlsPlugin::getClass(QAction *a) -{ - int id = ID(a); - if (id == FP_SELECT_SMALL_COMPONENTS) - { - return MeshFilterInterface::Selection; - } - - return MeshFilterInterface::PointSet; -} - -// Info() must return the longer string describing each filtering action -// (this string is used in the About plugin dialog) -const QString MlsPlugin::filterInfo(FilterIDType filterId) -{ - QString str = ""; - if (filterId & _PROJECTION_) - { - str += "Project a mesh (or a point set) onto the MLS surface defined by itself or another point set.
"; - } - - if (filterId & _MCUBE_) - { - str += - "Extract the iso-surface (as a mesh) of a MLS surface defined by the current point set (or mesh)" - "using the marching cubes algorithm. The coarse extraction is followed by an accurate projection" - "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 += - "
This is the algebraic point set surfaces (APSS) variant which is based on" - "the local fitting of algebraic spheres. It requires points equipped with oriented normals.
" - "See [Guennebaud and Gross, Algebraic Point Set Surfaces, Siggraph 2007]" - "and [Guennebaud et al., Dynamic Sampling and Rendering of APSS, Eurographics 2008]" - "for all the details about APSS."; - } - - if (filterId & _RIMLS_) - { - str += - "
This is the ##### MLS variant."; - } - - if (filterId == FP_RADIUS_FROM_DENSITY) - str = "Estimate the local point spacing (aka radius) around each vertex using a basic estimate of the local density."; - else if (filterId == FP_SELECT_SMALL_COMPONENTS) - str = "Select the small disconnected components of a mesh."; - - return str; -} - - -// 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 parmeter 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) -void MlsPlugin::initParameterSet(QAction* action, MeshDocument& md, FilterParameterSet& parlst) -//void ExtraSamplePlugin::initParList(QAction *action, MeshModel &m, FilterParameterSet &parlst) -{ - int id = ID(action); - MeshModel *target = md.mm(); - - if (id == FP_SELECT_SMALL_COMPONENTS) - { - parlst.addFloat("NbFaceRatio", - 0.1, - "Small component ratio", - "This ratio (between 0 and 1) defines the meaning of small as the threshold ratio between the number of faces" - "of the largest component and the other ones. A larger value will select more components."); - parlst.addBool( "NonClosedOnly", - false, - "Select only non closed components", - ""); - return; - } - else if (id == FP_RADIUS_FROM_DENSITY) - { - parlst.addInt("NbNeighbors", - 16, - "Number of neighbors", - "Number of neighbors used to estimate the local density. Larger values lead to smoother variations."); - return; - } - - if ((id & _PROJECTION_) || (id & _COLORIZE_)) - { - parlst.addMesh( "ControlMesh", target, "Point set", - "The point set (or mesh) which defines the MLS surface."); - parlst.addMesh( "ProxyMesh", target, "Proxy Mesh", - "The mesh that will be projected/resampled onto the MLS surface."); - parlst.addBool( "SelectionOnly", - target->cm.sfn>0, - "Selection only", - "If checked, only selected vertices will be projected."); - } - - if ( (id & _APSS_) || (id & _RIMLS_) ) - { - parlst.addFloat("FilterScale", - 2.0, - "MLS - Filter scale", - "Scale of the spatial low pass filter.\n" - "It is relative to the radius (local point spacing) of the vertices."); - parlst.addFloat("ProjectionAccuracy", - 1e-4, - "Projection - Accuracy (adv)", - "Threshold value used to stop the projections.\n" - "This value is scaled by the mean point spacing to get the actual threshold."); - parlst.addInt( "MaxProjectionIters", - 15, - "Projection - Max iterations (adv)", - "Max number of iterations for the projection."); - } - - if (id & _APSS_) - { - parlst.addFloat("SphericalParameter", - 1, - "MLS - Spherical parameter", - "Control the curvature of the fitted spheres: 0 is equivalent to a pure plane fit," - "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 !"); - if (!(id & _COLORIZE_)) - parlst.addBool( "AccurateNormal", - true, - "Accurate normals", - "If checked, use the accurate MLS gradient instead of the local approximation" - "to compute the normals."); - } - - if (id & _RIMLS_) - { - parlst.addFloat("SigmaN", - 0.75, - "MLS - Sharpness", - "Width of the filter used by the normal refitting weight." - "This weight function is a Gaussian on the distance between two unit vectors:" - "the current gradient and the input normal. Therefore, typical value range between 0.5 (sharp) to 2 (smooth)."); - parlst.addInt( "MaxRefittingIters", - 3, - "MLS - Max fitting iterations", - "Max number of fitting iterations. (0 or 1 is equivalent to the standard IMLS)"); - } - - if (id & _PROJECTION_) - { - parlst.addInt( "MaxSubdivisions", - 0, - "Refinement - Max subdivisions", - "Max number of subdivisions."); - parlst.addFloat("ThAngleInDegree", - 2, - "Refinement - Crease angle (degree)", - "Threshold angle between two faces controlling the refinement."); - } - - if (id & _AFRONT_) - { - } - - if ((id & _COLORIZE_)) - { - QStringList lst; - lst << "Mean" << "Gauss" << "K1" << "K2"; - if (id & _APSS_) - lst << "ApproxMean"; - - parlst.addEnum("CurvatureType", CT_MEAN, - lst, - "Curvature type", - QString("The type of the curvature to plot.") - + ((id & _APSS_) ? "
ApproxMean uses the radius of the fitted sphere as an approximation of the mean curvature." : "")); -// if ((id & _APSS_)) -// parlst.addBool( "ApproxCurvature", -// false, -// "Approx mean curvature", -// "If checked, use the radius of the fitted sphere as an approximation of the mean curvature."); - } - - if (id & _MCUBE_) - { - parlst.addInt( "Resolution", - 200, - "Grid Resolution", - "The resolution of the grid on which we run the marching cubes." - "This marching cube is memory friendly, so you can safely set large values up to 1000 or even more."); - } -} - -const int MlsPlugin::getRequirements(QAction *action) -{ - return 0; -} - -/** Predicate functor for adaptive refinement according to crease angle. - * - */ -template -struct EdgeAnglePredicate -{ - Scalar thCosAngle; - bool operator()(vcg::face::Pos ep) const - { - // FIXME why does the following fails: -// vcg::face::Pos op = ep; -// op.FlipF(); -// if (op.f) -// return vcg::Dot(ep.f->cN(), op.f->cN()) < thCosAngle; -// else -// return true; - - return vcg::Dot(ep.F()->cN(), ep.FFlip()->cN()) < thCosAngle; - } -}; - -/** compute the normal of a face as the average of its vertices */ -template -void UpdateFaceNormalFromVertex(MeshType& m) -{ - typedef typename MeshType::VertexType VertexType; - typedef typename VertexType::NormalType NormalType; - typedef typename VertexType::ScalarType ScalarType; - typedef typename MeshType::FaceIterator FaceIterator; - - for (FaceIterator f=m.face.begin(); f!=m.face.end(); ++f) - { - NormalType n; - n.SetZero(); - for(int j=0; j<3; ++j) - n += f->V(j)->cN(); - n.Normalize(); - f->N() = n; - } -} - -bool MlsPlugin::applyFilter(QAction* filter, MeshDocument& md, FilterParameterSet& par, vcg::CallBackPos* cb) -{ - int id = ID(filter); - - if (id == FP_RADIUS_FROM_DENSITY) - { - CMeshO::VertContainer& points = md.mm()->cm.vert; - if (!points.RadiusEnabled) - { - points.EnableRadius(); - } - int nbNeighbors = par.getInt("NbNeighbors"); - - assert(points.size()>=2); - KdTree knn(ConstDataWrapper(&points[0].cP(), points.size(), - size_t(points[1].cP().V()) - size_t(points[0].cP().V()))); - - knn.setMaxNofNeighbors(nbNeighbors); - for (size_t i = 0; i< points.size(); i++) - { - knn.doQueryK(points[i].cP()); - points[i].R() = 2. * sqrt(knn.getNeighborSquaredDistance(0)/float(knn.getNofFoundNeighbors())); - } - } - else if (id == FP_SELECT_SMALL_COMPONENTS) - { - MeshModel* mesh = md.mm(); - mesh->updateDataMask(MeshModel::MM_FACEFACETOPO | MeshModel::MM_FACEFLAGBORDER); - bool nonClosedOnly = par.getBool("NonClosedOnly"); - float ratio = par.getFloat("NbFaceRatio"); - vcg::tri::SmallComponent::Select(mesh->cm, ratio, nonClosedOnly); - } - else - { - // we are doing some MLS based stuff - - MeshModel* pPoints = 0; - if (id & _PROJECTION_) - { - if (par.getMesh("ControlMesh") == par.getMesh("ProxyMesh")) - { - // clone the control mesh - MeshModel* ref = par.getMesh("ControlMesh"); - pPoints = new MeshModel(); - pPoints->updateDataMask(ref->currentDataMask); - vcg::tri::Append::Mesh(pPoints->cm, ref->cm, false); - vcg::tri::UpdateBounding::Box(pPoints->cm); - pPoints->cm.Tr = ref->cm.Tr; - } - else - pPoints = par.getMesh("ControlMesh"); - } - else - pPoints = md.mm(); - - // create the MLS surface - cb(1, "Create the MLS data structures..."); - MlsSurface* mls = 0; - - RIMLS* rimls = 0; - APSS* apss = 0; - - if (id & _RIMLS_) - mls = rimls = new RIMLS(pPoints->cm); - else if (id & _APSS_) - mls = apss = new APSS(pPoints->cm); - else - { - assert(0); - } - - mls->setFilterScale(par.getFloat("FilterScale")); - mls->setMaxProjectionIters(par.getInt("MaxProjectionIters")); - mls->setProjectionAccuracy(par.getFloat("ProjectionAccuracy")); - - if (rimls) - { - rimls->setMaxRefittingIters(par.getInt("MaxRefittingIters")); - //mls.setMinRefittingIters(par.getFloat("MinRefittingIters")); - rimls->setSigmaN(par.getFloat("SigmaN")); - } - - if (apss) - { - apss->setSphericalParameter(par.getFloat("SphericalParameter")); - if (!(id & _COLORIZE_)) - apss->setGradientHint(par.getBool("AccurateNormal") ? GaelMls::MLS_DERIVATIVE_ACCURATE : GaelMls::MLS_DERIVATIVE_APPROX); - } - - MeshModel * mesh = 0; - - if (id & _PROJECTION_) - { - mesh = par.getMesh("ProxyMesh"); - bool selectionOnly = par.getBool("SelectionOnly"); - - if (selectionOnly) - vcg::tri::UpdateSelection::VertexFromFaceStrict(mesh->cm); - EdgeAnglePredicate edgePred; - edgePred.thCosAngle = cos(M_PI * par.getFloat("ThAngleInDegree")/180.); - - int nbRefinements = par.getInt("MaxSubdivisions"); - for (int k=0; kupdateDataMask(MeshModel::MM_FACEFACETOPO | MeshModel::MM_FACEFLAGBORDER); - - vcg::tri::UpdateNormals::PerFace(mesh->cm); - vcg::tri::UpdateNormals::NormalizeFace(mesh->cm); - //vcg::RefineE >(m.cm, vcg::MidPoint(), edgePred, false, cb); - vcg::RefineOddEvenE, vcg::EvenPointLoop > - (mesh->cm, vcg::OddPointLoop(), vcg::EvenPointLoop(), edgePred, selectionOnly, cb); - } - // project all vertices onto the MLS surface - for (unsigned int i = 0; i< mesh->cm.vert.size(); i++) - { - cb(1+98*i/mesh->cm.vert.size(), "MLS projection..."); - - if ( (!selectionOnly) || (mesh->cm.vert[i].IsS()) ) - mesh->cm.vert[i].P() = mls->project(mesh->cm.vert[i].P(), &mesh->cm.vert[i].N()); - } - } - - 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 ct = par.getEnum("CurvatureType"); - - 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()) ) - { - Point3f p = mls->project(mesh->cm.vert[i].P()); - float c = 0; - - if (ct==CT_APSS) - c = apss->approxMeanCurvature(p); - else - { - grad = mls->gradient(p); - hess = mls->hessian(p); - implicits::WeingartenMap W(grad,hess); - switch(ct) - { - case CT_MEAN: c = W.MeanCurvature(); break; - case CT_GAUSS: c = W.GaussCurvature(); break; - case CT_K1: c = W.K1(); break; - case CT_K2: c = W.K2(); break; - default: assert(0 && "invalid curvature type"); - } - } - mesh->cm.vert[i].Q() = c; - 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.05*d; - maxc -= 0.05*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 - // mesh = md.addNewMesh("afront mesh"); - // vcg::tri::AdvancingMLS afront(mesh->cm, *mls); - // //afront.BuildMesh(cb); - // afront._SeedFace(); - // for (int i=0; i<20120; ++i) - // afront.AddFace(); - // Log(0, "Advancing front MLS meshing done."); - // } - else if (id & _MCUBE_) - { - // create a new mesh - mesh = md.addNewMesh("mc_mesh"); - - typedef vcg::tri::MlsWalker > MlsWalker; - typedef vcg::tri::MarchingCubes MlsMarchingCubes; - MlsWalker walker; - walker.resolution = par.getInt("Resolution"); - - // iso extraction - MlsMarchingCubes mc(mesh->cm, walker); - walker.BuildMesh(mesh->cm, *mls, mc, cb); - - // accurate projection - for (unsigned int i = 0; i< mesh->cm.vert.size(); i++) - { - cb(1+98*i/mesh->cm.vert.size(), "MLS projection..."); - mesh->cm.vert[i].P() = mls->project(mesh->cm.vert[i].P(), &mesh->cm.vert[i].N()); - } - - // extra zero detection and removal - { - mesh->updateDataMask(MeshModel::MM_FACEFACETOPO | MeshModel::MM_FACEFLAGBORDER); - // selection... - vcg::tri::SmallComponent::Select(mesh->cm, 0.1); - // deletion... - vcg::tri::SmallComponent::DeleteFaceVert(mesh->cm); - mesh->clearDataMask(MeshModel::MM_FACEFACETOPO | MeshModel::MM_FACEFLAGBORDER); - } - - Log(0, "Marching cubes MLS meshing done."); - } - - delete mls; - if ( (id & _PROJECTION_) && par.getMesh("ControlMesh")!=pPoints) - { - delete pPoints; - } - - 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 - - return true; -} - -Q_EXPORT_PLUGIN(MlsPlugin) diff --git a/src/meshlabplugins/mlsfilters/mlsplugin.h b/src/meshlabplugins/mlsfilters/mlsplugin.h deleted file mode 100644 index 01793d4f8..000000000 --- a/src/meshlabplugins/mlsfilters/mlsplugin.h +++ /dev/null @@ -1,76 +0,0 @@ -/**************************************************************************** -* 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 MLSPLUGIN_H -#define MLSPLUGIN_H - -#include - -#include -#include - -class MlsPlugin : public QObject, public MeshFilterInterface -{ - Q_OBJECT - Q_INTERFACES(MeshFilterInterface) - -public: - - enum { - _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_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 - }; - - MlsPlugin(); - - virtual const QString filterName(FilterIDType filter); - virtual const QString filterInfo(FilterIDType filter); - const FilterClass getClass(QAction *a); - virtual bool autoDialog(QAction *) {return true;} - virtual void initParameterSet(QAction *,MeshDocument &md, FilterParameterSet &parent); - virtual const int getRequirements(QAction *action); - virtual bool applyFilter(QAction *filter, MeshDocument &m, FilterParameterSet &parent, vcg::CallBackPos *cb) ; - virtual bool applyFilter(QAction *filter, MeshModel &m, FilterParameterSet &parent, vcg::CallBackPos *cb) - {assert(0); return false;} -}; - -#endif diff --git a/src/meshlabplugins/mlsfilters/mlssurface.h b/src/meshlabplugins/mlsfilters/mlssurface.h deleted file mode 100644 index f74f774ab..000000000 --- a/src/meshlabplugins/mlsfilters/mlssurface.h +++ /dev/null @@ -1,203 +0,0 @@ -/**************************************************************************** -* 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 MLSSURFACE_H -#define MLSSURFACE_H - -#include "balltree.h" -#include -#include -#include - -namespace GaelMls { - -enum { - MLS_OK, - MLS_TOO_FAR, - MLS_TOO_MANY_ITERS, - MLS_NOT_SUPPORTED, - - MLS_DERIVATIVE_ACCURATE, - MLS_DERIVATIVE_APPROX, - MLS_DERIVATIVE_FINITEDIFF -}; - -template -class MlsSurface -{ - public: - typedef typename MeshType::ScalarType Scalar; - typedef vcg::Point3 VectorType; - typedef vcg::Matrix33 MatrixType; - typedef typename MeshType::VertContainer PointsType; - - MlsSurface(const MeshType& mesh) - : mMesh(mesh), mPoints(mesh.vert) - { - mCachedQueryPointIsOK = false; - - mAABB = mesh.bbox; - - // compute radii using a basic meshless density estimator - if (!mPoints.RadiusEnabled) - { - const_cast(mPoints).EnableRadius(); - computeVertexRaddi(); - } - - mFilterScale = 4.0; - mMaxNofProjectionIterations = 20; - mProjectionAccuracy = 1e-4; - mBallTree = 0; - mGradientHint = MLS_DERIVATIVE_ACCURATE; - mHessianHint = MLS_DERIVATIVE_ACCURATE; - - mDomainMinNofNeighbors = 4; - mDomainRadiusScale = 2.; - mDomainNormalScale = 1.; - } - - /** \returns the value of the reconstructed scalar field at point \a x */ - virtual Scalar potential(const VectorType& x, int* errorMask = 0) const = 0; - - /** \returns the gradient of the reconstructed scalar field at point \a x - * - * The method used to compute the gradient can be controlled with setGradientHint(). - */ - virtual VectorType gradient(const VectorType& x, int* errorMask = 0) const = 0; - - /** \returns the hessian matrix of the reconstructed scalar field at point \a x - * - * The method used to compute the hessian matrix can be controlled with setHessianHint(). - */ - virtual MatrixType hessian(const VectorType& x, int* errorMask = 0) const - { if (errorMask) *errorMask = MLS_NOT_SUPPORTED; return MatrixType(); } - - /** \returns the projection of point x onto the MLS surface, and optionnaly returns the normal in \a pNormal */ - virtual VectorType project(const VectorType& x, VectorType* pNormal = 0, int* errorMask = 0) const = 0; - - /** \returns whether \a x is inside the restricted surface definition domain */ - virtual bool isInDomain(const VectorType& x) const; - - /** \returns the mean curvature from the gradient vector and Hessian matrix. - */ - Scalar meanCurvature(const VectorType& gradient, const MatrixType& hessian) const; - - /** set the scale of the spatial filter */ - void setFilterScale(Scalar v); - /** set the maximum number of iterations during the projection */ - void setMaxProjectionIters(int n); - /** set the threshold factor to detect convergence of the iterations */ - void setProjectionAccuracy(Scalar v); - - /** set a hint on how to compute the gradient - * - * Possible values are MLS_DERIVATIVE_ACCURATE, MLS_DERIVATIVE_APPROX, MLS_DERIVATIVE_FINITEDIFF - */ - void setGradientHint(int h); - - /** set a hint on how to compute the hessian matrix - * - * Possible values are MLS_DERIVATIVE_ACCURATE, MLS_DERIVATIVE_APPROX, MLS_DERIVATIVE_FINITEDIFF - */ - void setHessianHint(int h); - - inline const MeshType& mesh() const { return mMesh; } - /** a shortcut for mesh().vert */ - inline const PointsType& points() const { return mPoints; } - - inline ConstDataWrapper positions() const - { - return ConstDataWrapper(&mPoints[0].cP(), mPoints.size(), - size_t(mPoints[1].cP().V()) - size_t(mPoints[0].cP().V())); - } - inline ConstDataWrapper normals() const - { - return ConstDataWrapper(&mPoints[0].cN(), mPoints.size(), - size_t(mPoints[1].cN().V()) - size_t(mPoints[0].cN().V())); - } - inline ConstDataWrapper radii() const - { - return ConstDataWrapper(&mPoints[0].cR(), mPoints.size(), - size_t(&mPoints[1].cR()) - size_t(&mPoints[0].cR())); - } - const vcg::Box3& boundingBox() const { return mAABB; } - - static const Scalar InvalidValue() { return 12345679810.11121314151617; } - - protected: - void computeNeighborhood(const VectorType& x, bool computeDerivatives) const; - void computeVertexRaddi(); - void requestSecondDerivatives() const; - - struct PointToPointSqDist - { - inline bool operator()(const VectorType &a, const VectorType &b, Scalar& refD2, VectorType &q) const - { -// std::cout << a.X() << a.Y() << a.Z() << " - " << b.X() << b.Y() << b.Z() << -// " => " << vcg::Distance(a, b) << " < " << refD2 << "\n"; - Scalar d2 = vcg::SquaredDistance(a, b); - if (d2>refD2) - return false; - - refD2 = d2; - q = a; - return true; - } - }; - - class DummyObjectMarker {}; - - protected: - const MeshType& mMesh; - const PointsType& mPoints; - vcg::Box3 mAABB; - int mGradientHint; - int mHessianHint; - - BallTree* mBallTree; - - int mMaxNofProjectionIterations; - Scalar mFilterScale; - Scalar mAveragePointSpacing; - Scalar mProjectionAccuracy; - - int mDomainMinNofNeighbors; - float mDomainRadiusScale; - float mDomainNormalScale; - - // cached values: - mutable bool mCachedQueryPointIsOK; - mutable VectorType mCachedQueryPoint; - mutable Neighborhood mNeighborhood; - mutable std::vector mCachedWeights; - mutable std::vector mCachedWeightDerivatives; - mutable std::vector mCachedWeightGradients; - mutable std::vector mCachedWeightSecondDerivatives; -}; - -} // namespace - -#include "mlssurface.tpp" - -#endif // MLSSURFACE_H diff --git a/src/meshlabplugins/mlsfilters/mlssurface.tpp b/src/meshlabplugins/mlsfilters/mlssurface.tpp deleted file mode 100644 index 62f13b117..000000000 --- a/src/meshlabplugins/mlsfilters/mlssurface.tpp +++ /dev/null @@ -1,262 +0,0 @@ -/**************************************************************************** -* 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 "mlssurface.h" -#include "kdtree.h" -#include "mlsutils.h" -#include -#include -#include - -namespace GaelMls { - -// template -// MlsSurface<_Scalar>::MlsSurface(const MeshModel& m) -// : mMesh(m) -// { -// mCachedQueryPointIsOK = false; -// -// mPoints.resize(m.cm.vert.size()); -// mNormals.resize(m.cm.vert.size()); -// mRadii.resize(m.cm.vert.size()); -// mAABB.Set(mMesh.cm.vert[0].cP()); -// for (uint i = 0; i< m.cm.vert.size(); i++) -// { -// mPoints[i] = /*vcg::vector_cast*/(mMesh.cm.vert[i].cP()); -// mNormals[i] = /*vcg::vector_cast*/(mMesh.cm.vert[i].cN()); -// mAABB.Add(mMesh.cm.vert[i].cP()); -// } -// -// // compute radii using a basic meshless density estimator -// computeVertexRaddi(); -// -// mFilterScale = 4.0; -// mMaxNofProjectionIterations = 20; -// mProjectionAccuracy = 1e-4; -// mBallTree = 0; -// } - -template -void MlsSurface<_MeshType>::setFilterScale(Scalar v) -{ - mFilterScale = v; - mCachedQueryPointIsOK = false; - if (mBallTree) - mBallTree->setRadiusScale(mFilterScale); -} - -template -void MlsSurface<_MeshType>::setMaxProjectionIters(int n) -{ - mMaxNofProjectionIterations = n; - mCachedQueryPointIsOK = false; -} - -template -void MlsSurface<_MeshType>::setProjectionAccuracy(Scalar v) -{ - mProjectionAccuracy = v; - mCachedQueryPointIsOK = false; -} - -template -void MlsSurface<_MeshType>::setGradientHint(int h) -{ - mGradientHint = h; - mCachedQueryPointIsOK = false; -} - -template -void MlsSurface<_MeshType>::setHessianHint(int h) -{ - mHessianHint = h; - mCachedQueryPointIsOK = false; -} - -template -void MlsSurface<_MeshType>::computeVertexRaddi() -{ - #if 0 - int nbNeighbors = 16; - vcg::Octree knn; - knn.Set(mPoints.begin(), mPoints.end()); - std::vector nearest_objects; - std::vector nearest_points; - std::vector sqDistances; - mAveragePointSpacing = 0; - for (uint i = 0; i< mPoints.size(); i++) - { - DummyObjectMarker dom; - PointToPointSqDist dfunc; - Scalar max_dist2 = 1e9;//std::numeric_limits::max(); - knn.GetKClosest(dfunc, dom, nbNeighbors, mPoints[i], - max_dist2, nearest_objects, sqDistances, nearest_points); -// for (int j=0; i=2); - KdTree knn(positions()); - - knn.setMaxNofNeighbors(nbNeighbors); - mAveragePointSpacing = 0; - for (size_t i = 0; i< mPoints.size(); i++) - { - knn.doQueryK(mPoints[i].cP()); - const_cast(mPoints)[i].R() = 2. * sqrt(knn.getNeighborSquaredDistance(0)/Scalar(knn.getNofFoundNeighbors())); - mAveragePointSpacing += mPoints[i].cR(); - } - mAveragePointSpacing /= Scalar(mPoints.size()); - - #endif -} - -template -void MlsSurface<_MeshType>::computeNeighborhood(const VectorType& x, bool computeDerivatives) const -{ - if (!mBallTree) - { - const_cast*&>(mBallTree) = new BallTree(positions(), radii()); - const_cast*>(mBallTree)->setRadiusScale(mFilterScale); - } - mBallTree->computeNeighbors(x, &mNeighborhood); - size_t nofSamples = mNeighborhood.size(); - - // compute spatial weights and partial derivatives - mCachedWeights.resize(nofSamples); - if (computeDerivatives) - { - mCachedWeightDerivatives.resize(nofSamples); - mCachedWeightGradients.resize(nofSamples); - } - else - mCachedWeightGradients.clear(); - - for (size_t i=0; i -void MlsSurface<_MeshType>::requestSecondDerivatives() const -{ - //if (!mSecondDerivativeUptodate) - { - size_t nofSamples = mNeighborhood.size(); - if (nofSamples>mCachedWeightSecondDerivatives.size()) - mCachedWeightSecondDerivatives.resize(nofSamples+10); - - { - for (size_t i=0 ; i -typename MlsSurface<_MeshType>::Scalar -MlsSurface<_MeshType>::meanCurvature(const VectorType& gradient, const MatrixType& hessian) const -{ - Scalar gl = gradient.Norm(); - return (gl*gl*hessian.Trace() - vcg::Dot(gradient, VectorType(hessian * gradient))) / (2.*gl*gl*gl); -} - -template -bool MlsSurface<_MeshType>::isInDomain(const VectorType& x) const -{ - if ((!mCachedQueryPointIsOK) || mCachedQueryPoint!=x) - { - computeNeighborhood(x, false); - } - int nb = mNeighborhood.size(); - if (nb rs2; - ++i; - } - } - else - { - Scalar s = 1./(mDomainNormalScale*mDomainNormalScale) - 1.f; - while (out && i rs2; - ++i; - } - } - return !out; -} - -// template class MlsSurface; -// template class MlsSurface; - -} diff --git a/src/meshlabplugins/mlsfilters/mlsutils.h b/src/meshlabplugins/mlsfilters/mlsutils.h deleted file mode 100644 index 3e4e88053..000000000 --- a/src/meshlabplugins/mlsfilters/mlsutils.h +++ /dev/null @@ -1,140 +0,0 @@ -/**************************************************************************** -* 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 VCGADDONS_H -#define VCGADDONS_H - -template -class ConstDataWrapper -{ -public: - typedef _DataType DataType; - inline ConstDataWrapper() - : mpData(0), mStride(0), mSize(0) - {} - inline ConstDataWrapper(const DataType* pData, int size, int stride = sizeof(DataType)) - : mpData(reinterpret_cast(pData)), mStride(stride), mSize(size) - {} - inline const DataType& operator[] (int i) const - { - return *reinterpret_cast(mpData + i*mStride); - } - inline int size() const { return mSize; } -protected: - const unsigned char* mpData; - int mStride; - int mSize; -}; - -namespace vcg { - -template -inline Point3 CwiseMul(Point3 const & p1, Point3 const & p2) -{ - return Point3(p1.X()*p2.X(), p1.Y()*p2.Y(), p1.Z()*p2.Z()); -} - -template -inline Point3 Min(Point3 const & p1, Point3 const & p2) -{ - return Point3(std::min(p1.X(), p2.X()), std::min(p1.Y(), p2.Y()), std::min(p1.Z(), p2.Z())); -} - -template -inline Point3 Max(Point3 const & p1, Point3 const & p2) -{ - return Point3(std::max(p1.X(), p2.X()), std::max(p1.Y(), p2.Y()), std::max(p1.Z(), p2.Z())); -} - -template -inline Scalar MaxCoeff(Point3 const & p) -{ - return std::max(std::max(p.X(), p.Y()), p.Z()); -} - -template -inline Scalar MinCoeff(Point3 const & p) -{ - return std::min(std::min(p.X(), p.Y()), p.Z()); -} - -template -inline Scalar Dot(Point3 const & p1, Point3 const & p2) -{ - return p1.X() * p2.X() + p1.Y() * p2.Y() + p1.Z() * p2.Z(); -} - -template -inline Point3 Cross(Point3 const & p1, Point3 const & p2) -{ - return p1 ^ p2; -} - -template -inline Point3 CwiseAdd(Point3 const & p1, Scalar s) -{ - return Point3(p1.X() + s, p1.Y() + s, p1.Z() + s); -} - -template -inline int MaxCoeffId(Point3 const & p) -{ - if (p.X()>p.Y()) - return p.X()>p.Z() ? 0 : 2; - else - return p.Y()>p.Z() ? 1 : 2; -} - -template -inline int MinCoeffId(Point3 const & p) -{ - if (p.X() -inline Point3 Point3Cast(const Point3& p) -{ - return Point3(p.X(), p.Y(), p.Z()); -} - -template -Scalar Distance(const Point3 &p, const Box3 &bbox) -{ - Scalar dist2 = 0.; - Scalar aux; - for (int k=0 ; k<3 ; ++k) - { - if ( (aux = (p[k]-bbox.min[k]))<0. ) - dist2 += aux*aux; - else if ( (aux = (bbox.max[k]-p[k]))<0. ) - dist2 += aux*aux; - } - return sqrt(dist2); -} - -} - -#endif diff --git a/src/meshlabplugins/mlsfilters/priorityqueue.h b/src/meshlabplugins/mlsfilters/priorityqueue.h deleted file mode 100644 index 13aae0247..000000000 --- a/src/meshlabplugins/mlsfilters/priorityqueue.h +++ /dev/null @@ -1,122 +0,0 @@ -/**************************************************************************** -* 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 _PriorityQueue_h_ -#define _PriorityQueue_h_ - -/** Implements a bounded-size max priority queue using a heap -*/ -template -class HeapMaxPriorityQueue -{ - struct Element - { - Weight weight; - Index index; - }; - -public: - - HeapMaxPriorityQueue(void) - { - mElements = 0; - mMaxSize = 0; - } - - inline void setMaxSize(int maxSize) - { - if (mMaxSize!=maxSize) - { - mMaxSize = maxSize; - delete[] mElements; - mElements = new Element[mMaxSize]; - mpOffsetedElements = (mElements-1); - } - init(); - } - - inline void init() { mCount = 0; } - - inline bool isFull() const { return mCount == mMaxSize; } - - /** returns number of elements inserted in queue - */ - inline int getNofElements() const { return mCount; } - - inline Weight getWeight(int i) const { return mElements[i].weight; } - inline Index getIndex(int i) const { return mElements[i].index; } - - inline Weight getTopWeight() const { return mElements[0].weight; } - - inline void insert(Index index, Weight weight) - { - if (mCount==mMaxSize) - { - if (weightweight < mpOffsetedElements[k+1].weight)) - z = &(mpOffsetedElements[++k]); - - if(weight >= z->weight) - break; - mpOffsetedElements[j] = *z; - j = k; - k = 2 * j; - } - mpOffsetedElements[j].weight = weight; - mpOffsetedElements[j].index = index; - } - } - else - { - int i, j; - i = ++mCount; - while (i >= 2) - { - j = i >> 1; - Element& y = mpOffsetedElements[j]; - if(weight <= y.weight) - break; - mpOffsetedElements[i] = y; - i = j; - } - mpOffsetedElements[i].index = index; - mpOffsetedElements[i].weight = weight; - } - } - -protected: - - int mCount; - int mMaxSize; - Element* mElements; - Element* mpOffsetedElements; -}; - -#endif diff --git a/src/meshlabplugins/mlsfilters/rimls.cpp b/src/meshlabplugins/mlsfilters/rimls.cpp deleted file mode 100644 index 645f20040..000000000 --- a/src/meshlabplugins/mlsfilters/rimls.cpp +++ /dev/null @@ -1,11 +0,0 @@ - -#include -#include "rimls.h" -#include "rimls.tpp" -#include - -namespace GaelMls { - -template class RIMLS; - -} diff --git a/src/meshlabplugins/mlsfilters/rimls.h b/src/meshlabplugins/mlsfilters/rimls.h deleted file mode 100644 index dbc40f0f8..000000000 --- a/src/meshlabplugins/mlsfilters/rimls.h +++ /dev/null @@ -1,104 +0,0 @@ -/**************************************************************************** -* 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 RIMLS_H -#define RIMLS_H - -#include "mlssurface.h" - -namespace GaelMls { - -template -class RIMLS : public MlsSurface<_MeshType> -{ - typedef _MeshType MeshType; - typedef MlsSurface<_MeshType> Base; - - typedef typename Base::Scalar Scalar; - typedef typename Base::VectorType VectorType; - typedef typename Base::MatrixType MatrixType; - using Base::mCachedQueryPointIsOK; - using Base::mCachedQueryPoint; - using Base::mNeighborhood; - using Base::mCachedWeights; - using Base::mCachedWeightDerivatives; - using Base::mCachedWeightGradients; - using Base::mCachedWeightSecondDerivatives; - using Base::mBallTree; - using Base::mPoints; - using Base::mFilterScale; - using Base::mMaxNofProjectionIterations; - using Base::mAveragePointSpacing; - using Base::mProjectionAccuracy; - - public: - - RIMLS(const MeshType& points) - : Base(points) - { - mSigmaR = 0; - mSigmaN = 0.8; - mRefittingThreshold = 1e-3; - mMinRefittingIters = 1; - mMaxRefittingIters = 3; - } - - 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 = 0) const; - virtual VectorType project(const VectorType& x, VectorType* pNormal = 0, int* errorMask = 0) const; - - void setSigmaR(Scalar v); - void setSigmaN(Scalar v); - void setRefittingThreshold(Scalar v); - void setMinRefittingIters(int n); - void setMaxRefittingIters(int n); - - protected: - bool computePotentialAndGradient(const VectorType& x) const; - bool mlsHessian(const VectorType& x, MatrixType& hessian) const; - - protected: - - int mMinRefittingIters; - int mMaxRefittingIters; - Scalar mRefittingThreshold; - Scalar mSigmaN; - Scalar mSigmaR; - - // cached values: - mutable VectorType mCachedGradient; - mutable Scalar mCachedPotential; - - mutable Scalar mCachedSumW; - mutable std::vector mCachedRefittingWeights;; - mutable VectorType mCachedSumN; - mutable VectorType mCachedSumGradWeight; - mutable VectorType mCachedSumGradPotential; -}; - -} - -//#include "rimls.tpp" - -#endif diff --git a/src/meshlabplugins/mlsfilters/rimls.tpp b/src/meshlabplugins/mlsfilters/rimls.tpp deleted file mode 100644 index 9fbffea0d..000000000 --- a/src/meshlabplugins/mlsfilters/rimls.tpp +++ /dev/null @@ -1,307 +0,0 @@ -/**************************************************************************** -* 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 "rimls.h" -#include "kdtree.h" -#include "mlsutils.h" -#include - -using namespace vcg; - -namespace GaelMls { - -// template -// RIMLS<_Scalar>::RIMLS(const MeshModel& m) -// : Base(m) -// { -// mSigmaR = 0; -// mSigmaN = 0.8; -// mRefittingThreshold = 1e-3; -// mMinRefittingIters = 1; -// mMaxRefittingIters = 3; -// } - -template -void RIMLS<_MeshType>::setSigmaR(Scalar v) -{ - mSigmaR = v; - mCachedQueryPointIsOK = false; -} - -template -void RIMLS<_MeshType>::setSigmaN(Scalar v) -{ - mSigmaN = v; - mCachedQueryPointIsOK = false; -} - -template -void RIMLS<_MeshType>::setRefittingThreshold(Scalar v) -{ - mRefittingThreshold = v; - mCachedQueryPointIsOK = false; -} - -template -void RIMLS<_MeshType>::setMinRefittingIters(int n) -{ - mMinRefittingIters = n; - mCachedQueryPointIsOK = false; -} - -template -void RIMLS<_MeshType>::setMaxRefittingIters(int n) -{ - mMaxRefittingIters = n; - mCachedQueryPointIsOK = false; -} - -template -typename RIMLS<_MeshType>::Scalar RIMLS<_MeshType>::potential(const VectorType& x, int* errorMask) const -{ - if ((!mCachedQueryPointIsOK) || mCachedQueryPoint!=x) - { - if (!computePotentialAndGradient(x)) - { - if (errorMask) - *errorMask = MLS_TOO_FAR; - return Base::InvalidValue(); - } - } - - return mCachedPotential; -} - -template -typename RIMLS<_MeshType>::VectorType RIMLS<_MeshType>::gradient(const VectorType& x, int* errorMask) const -{ - if ((!mCachedQueryPointIsOK) || mCachedQueryPoint!=x) - { - if (!computePotentialAndGradient(x)) - { - if (errorMask) - *errorMask = MLS_TOO_FAR; - return VectorType(0,0,0); - } - } - - return mCachedGradient; -} - -template -typename RIMLS<_MeshType>::MatrixType RIMLS<_MeshType>::hessian(const VectorType& x, int* errorMask) const -{ - if ((!mCachedQueryPointIsOK) || mCachedQueryPoint!=x) - { - if (!computePotentialAndGradient(x)) - { - if (errorMask) - *errorMask = MLS_TOO_FAR; - return MatrixType(); - } - } - - MatrixType hessian; - mlsHessian(x, hessian); - return hessian; -} - -template -typename RIMLS<_MeshType>::VectorType RIMLS<_MeshType>::project(const VectorType& x, VectorType* pNormal, int* errorMask) const -{ - int iterationCount = 0; - VectorType position = x; - VectorType normal; - Scalar delta; - Scalar epsilon = mAveragePointSpacing * mProjectionAccuracy; - do { - if (!computePotentialAndGradient(position)) - { - if (errorMask) - *errorMask = MLS_TOO_FAR; - //std::cerr << " proj failed\n"; - return x; - } - - normal = mCachedGradient; - normal.Normalize(); - delta = mCachedPotential; - position = position - normal*delta; - } while ( fabs(delta)>epsilon && ++iterationCount=mMaxNofProjectionIterations && errorMask) - *errorMask = MLS_TOO_MANY_ITERS; - - if (pNormal) - *pNormal = normal; - - return position; -} - -template -bool RIMLS<_MeshType>::computePotentialAndGradient(const VectorType& x) const -{ - Base::computeNeighborhood(x, true); - unsigned int nofSamples = mNeighborhood.size(); - - if (nofSamples<1) - { - mCachedGradient.SetZero(); - mCachedQueryPoint = x; - mCachedPotential = 1e9; - mCachedQueryPointIsOK = false; - return false; - } - - if (mCachedRefittingWeights.size()0) - invSigmaR2 = Scalar(1) / (mSigmaR*mSigmaR); - VectorType sumGradWeight; - VectorType sumGradWeightPotential; - Scalar sumW; - - int iterationCount = 0; - do - { - previousGrad = grad; - sumGradWeight.SetZero(); - sumGradWeightPotential.SetZero(); - sumN.SetZero(); - potential = 0.; - sumW = 0.; - - for (unsigned int i=0; i 0) - { - refittingWeight = exp(-vcg::SquaredNorm(normal - previousGrad) * invSigma2); -// if (mSigmaR>0) -// { -// Scalar residual = (ddotn - potentialPrev) / mRadii.at(id); -// refittingWeight *= exp(-residual*residual * invSigmaR2); -// } - } - mCachedRefittingWeights.at(i) = refittingWeight; - Scalar w = mCachedWeights.at(i) * refittingWeight; - VectorType gw = mCachedWeightGradients.at(i) * refittingWeight; - - sumGradWeight += gw; - sumGradWeightPotential += gw * f; - sumN += normal * w; - potential += w * f; - sumW += w; - } - - if(sumW==0.) - { - return false; - } - - potential /= sumW; - grad = (-sumGradWeight*potential + sumGradWeightPotential + sumN) * (1./sumW); - - iterationCount++; - - } while ( (iterationCount < mMinRefittingIters) - || ( vcg::SquaredNorm(grad - previousGrad) > mRefittingThreshold && iterationCount < mMaxRefittingIters) ); - - mCachedGradient = grad; - mCachedPotential = potential; - mCachedQueryPoint = x; - mCachedQueryPointIsOK = true; - - mCachedSumGradWeight = sumGradWeight; - mCachedSumN = sumN; - mCachedSumW = sumW; - mCachedSumGradPotential = sumGradWeightPotential; - - return true; -} - -template -bool RIMLS<_MeshType>::mlsHessian(const VectorType& x, MatrixType& hessian) const -{ - this->requestSecondDerivatives(); - // at this point we assume computePotentialAndGradient has been called first - - uint nofSamples = mNeighborhood.size(); - - const VectorType& sumGradWeight = mCachedSumGradWeight; - const VectorType& sumGradWeightPotential = mCachedSumGradPotential ; - const VectorType& sumN = mCachedSumN; - const Scalar& sumW = mCachedSumW; - const Scalar invW = 1.f/sumW; - - for (uint k=0 ; k<3 ; ++k) - { - VectorType sumDGradWeight; sumDGradWeight.SetZero(); - VectorType sumDWeightNormal; sumDWeightNormal.SetZero(); - VectorType sumGradWeightNk; sumGradWeightNk.SetZero(); - VectorType sumDGradWeightPotential; sumDGradWeightPotential.SetZero(); - - for (unsigned int i=0; i -class SmallComponent -{ - -public: - typedef _MeshType MeshType; - typedef typename MeshType::VertexType VertexType; - typedef typename MeshType::VertexPointer VertexPointer; - typedef typename MeshType::VertexIterator VertexIterator; - typedef typename MeshType::FaceType FaceType; - typedef typename MeshType::FacePointer FacePointer; - typedef typename MeshType::FaceIterator FaceIterator; - - static int Select(MeshType &m, float nbFaceRatio = 0.1, bool nonClosedOnly = false) - { - assert(tri::HasFFAdjacency(m) && "The small component selection procedure requires face to face adjacency."); - - // the different components as a list of face pointer - std::vector< std::vector > components; - - for(int faceSeed = 0; faceSeed activefaces; - activefaces.push_back(&m.face[faceSeed]); - while (!activefaces.empty()) - { - FacePointer f = activefaces.back(); - activefaces.pop_back(); - if (f->IsS()) - continue; - - f->SetS(); - components.back().push_back(f); - for (int k=0; k<3; ++k) - { - if (f->IsB(k)) - continue; - FacePointer of = f->FFp(k); - if (!of->IsS()) - activefaces.push_back(of); - } - } - ++faceSeed; - } - UpdateSelection::ClearFace(m); - - // now the segmentation is done, let's compute the absolute face count threshold - int total_selected = 0; - int maxComponent = 0; - for (int i=0; i " << components[i].size() << "\n"; - total_selected += components[i].size(); - maxComponent = std::max(maxComponent,components[i].size()); - } - int remaining = m.face.size() - total_selected; - int th = std::max(maxComponent,remaining) * nbFaceRatio; - - int selCount = 0; - for (int i=0; iSetS(); - } - } - return selCount; - } - - static void DeleteFaceVert(MeshType &m) - { - typename MeshType::FaceIterator fi; - typename MeshType::VertexIterator vi; - UpdateSelection::ClearVertex(m); - UpdateSelection::VertexFromFaceStrict(m); - - for(fi=m.face.begin();fi!=m.face.end();++fi) - if (!(*fi).IsD() && (*fi).IsS() ) - tri::Allocator::DeleteFace(m,*fi); - for(vi=m.vert.begin();vi!=m.vert.end();++vi) - if (!(*vi).IsD() && (*vi).IsS() ) - tri::Allocator::DeleteVertex(m,*vi); - } - -}; // end class - -} // End namespace -} // End namespace - - -#endif