From 7befff7becc7364fd6ed40135ef848012ce0eb07 Mon Sep 17 00:00:00 2001 From: cnr-isti-vclab Date: Tue, 28 Oct 2008 00:59:46 +0000 Subject: [PATCH] make point2 derived Eigen's Matrix, and a set of minimal fixes to make meshlab compile with both old and new version. The fixes include: - dot product: vec0 * vec1 => vec0.dot(vec1) (I added .dot() to the old Point classes too) - Transpose: Transpose is an Eigen type, so we cannot keep it if Eigen is used. Therefore I added a .tranpose() to old matrix classes, and modified most of the Transpose() to transpose() both in vcg and meshlab. In fact, transpose() are free with Eigen, it simply returns a transpose expression without copies. On the other be carefull: m = m.transpose() won't work as expected, here me must evaluate to a temporary: m = m.transpose().eval(); However, this operation in very rarely needed: you transpose at the same sime you set m, or you use m.transpose() directly. - the last issue is Normalize which both modifies *this and return a ref to it. This behavior don't make sense anymore when using expression template, e.g., in (a+b).Normalize(), the type of a+b if not a Point (or whatever Vector types), it an expression of the addition of 2 points, so we cannot modify the value of *this, since there is no value. Therefore I've already changed all those .Normalize() of expressions to the Eigen's version .normalized(). - Finally I've changed the Zero to SetZero in the old Point classes too. --- .../tri_edge_collapse_quadric.h | 10 +- vcg/complex/trimesh/autoalign_4pcs.h | 28 +- vcg/complex/trimesh/clean.h | 2 +- vcg/complex/trimesh/clustering.h | 56 +-- vcg/complex/trimesh/create/resampler.h | 2 +- vcg/complex/trimesh/hole.h | 2 +- vcg/complex/trimesh/inertia.h | 60 +-- vcg/complex/trimesh/refine.h | 2 +- vcg/complex/trimesh/update/curvature.h | 146 +++---- vcg/math/base.h | 20 +- vcg/math/eigen.h | 2 + vcg/math/eigen_vcgaddons.h | 37 +- vcg/math/lin_algebra.h | 360 +++++++++--------- vcg/math/linear.h | 6 +- vcg/math/matrix.h | 17 +- vcg/math/matrix33.h | 2 +- vcg/math/quadric.h | 24 +- vcg/simplex/vertex/distance.h | 2 +- vcg/space/color4.h | 84 ++-- vcg/space/deprecated_point2.h | 359 +++++++++++++++++ vcg/space/deprecated_point4.h | 42 +- vcg/space/fitting3.h | 2 +- vcg/space/intersection3.h | 22 +- vcg/space/normal_extrapolation.h | 58 +-- vcg/space/point.h | 122 +++--- vcg/space/point2.h | 331 +++++----------- vcg/space/point3.h | 11 + wrap/gl/addons.h | 13 +- wrap/gl/space.h | 51 ++- wrap/gl/trimesh.h | 186 ++++----- wrap/gui/coordinateframe.cpp | 2 +- wrap/io_trimesh/import_ptx.h | 2 +- 32 files changed, 1150 insertions(+), 913 deletions(-) create mode 100644 vcg/space/deprecated_point2.h diff --git a/vcg/complex/local_optimization/tri_edge_collapse_quadric.h b/vcg/complex/local_optimization/tri_edge_collapse_quadric.h index f333a91b..1b22241e 100644 --- a/vcg/complex/local_optimization/tri_edge_collapse_quadric.h +++ b/vcg/complex/local_optimization/tri_edge_collapse_quadric.h @@ -403,7 +403,7 @@ public: { if(Params().NormalCheck){ nn=NormalizedNormal(*x.F()); - ndiff=nn*on[i++]; + ndiff=nn.dot(on[i++]); if(ndiffcP()); + p.SetOffset( p.Direction().dot((*pf).V(0)->cP())); // Calcolo quadrica delle facce q.ByPlane(p); @@ -559,10 +559,10 @@ static void InitQuadric(TriMeshType &m) // Nota che la lunghezza dell'edge DEVE essere Normalizzata // poiche' la pesatura in funzione dell'area e'gia fatta in p.Direction() // Senza la normalize il bordo e' pesato in funzione della grandezza della mesh (mesh grandi non decimano sul bordo) - pb.SetDirection(p.Direction() ^ ( (*pf).V1(j)->cP() - (*pf).V(j)->cP() ).Normalize()); + pb.SetDirection(p.Direction() ^ ( (*pf).V1(j)->cP() - (*pf).V(j)->cP() ).normalized()); if( (*pf).IsB(j) ) pb.SetDirection(pb.Direction()* (ScalarType)Params().BoundaryWeight); // amplify border planes else pb.SetDirection(pb.Direction()* (ScalarType)(Params().BoundaryWeight/100.0)); // and consider much less quadric for quality - pb.SetOffset(pb.Direction() * (*pf).V(j)->cP()); + pb.SetOffset(pb.Direction().dot((*pf).V(j)->cP())); q.ByPlane(pb); if( (*pf).V (j)->IsW() ) QH::Qd((*pf).V (j)) += q; // Sommo le quadriche diff --git a/vcg/complex/trimesh/autoalign_4pcs.h b/vcg/complex/trimesh/autoalign_4pcs.h index 566ceb21..4ae5c4ca 100644 --- a/vcg/complex/trimesh/autoalign_4pcs.h +++ b/vcg/complex/trimesh/autoalign_4pcs.h @@ -70,7 +70,7 @@ public: typedef typename MeshType::CoordType CoordType; typedef typename MeshType::VertexIterator VertexIterator; typedef typename MeshType::VertexType VertexType; - typedef vcg::Point4< vcg::Point3 > FourPoints ; + typedef vcg::Point4< vcg::Point3 > FourPoints; typedef vcg::GridStaticPtr GridType; /* class for Parameters */ @@ -110,11 +110,11 @@ private: /* returns the closest point between to segments x1-x2 and x3-x4. */ - void IntersectionLineLine( const CoordType & x1,const CoordType & x2,const CoordType & x3,const CoordType & x4, - CoordType&x){ - CoordType a = x2-x1, b = x4-x3, c = x3-x1; - x = x1 + a * ( (c^b)*(a^b)) / (a^b).SquaredNorm(); - } + void IntersectionLineLine(const CoordType & x1,const CoordType & x2,const CoordType & x3,const CoordType & x4, CoordType&x) + { + CoordType a = x2-x1, b = x4-x3, c = x3-x1; + x = x1 + a * ((c^b).dot(a^b)) / (a^b).SquaredNorm(); + } @@ -289,7 +289,7 @@ FourPCS::SelectCoplanarBase(){ int id = rand()/(float)RAND_MAX * (P->vert.size()-1); ScalarType dd = (P->vert[id].P() - B[1]).Norm(); if( ( dd < side + dtol) && (dd > side - dtol)){ - ScalarType angle = fabs( ( P->vert[id].P()-B[1]).Normalize() * (B[1]-B[0]).Normalize()); + ScalarType angle = fabs( ( P->vert[id].P()-B[1]).normalized().dot((B[1]-B[0]).normalized())); if( angle < bestv){ bestv = angle; best = id; @@ -301,7 +301,7 @@ FourPCS::SelectCoplanarBase(){ B[2] = P->vert[best].P(); //printf("B[2] %d\n",best); - CoordType n = ((B[0]-B[1]).Normalize() ^ (B[2]-B[1]).Normalize()).Normalize(); + CoordType n = ((B[0]-B[1]).normalized() ^ (B[2]-B[1]).normalized()).normalized(); CoordType B4 = B[1] + (B[0]-B[1]) + (B[2]-B[1]); VertexType * v =0; ScalarType dist,radius = dtol*4.0; @@ -322,7 +322,7 @@ FourPCS::SelectCoplanarBase(){ return false; best = -1; bestv=std::numeric_limits::max(); for(i = 0; i P() - B[1]).Normalize() * n); + ScalarType angle = fabs((closests[i]->P() - B[1]).normalized().dot(n)); if( angle < bestv){ bestv = angle; best = i; @@ -337,8 +337,8 @@ FourPCS::SelectCoplanarBase(){ std::swap(B[1],B[2]); IntersectionLineLine(B[0],B[1],B[2],B[3],x); - r1 = (x - B[0])*(B[1]-B[0]) / (B[1]-B[0]).SquaredNorm(); - r2 = (x - B[2])*(B[3]-B[2]) / (B[3]-B[2]).SquaredNorm(); + r1 = (x - B[0]).dot(B[1]-B[0]) / (B[1]-B[0]).SquaredNorm(); + r2 = (x - B[2]).dot(B[3]-B[2]) / (B[3]-B[2]).SquaredNorm(); if( ((B[0]+(B[1]-B[0])*r1)-(B[2]+(B[3]-B[2])*r2)).Norm() > prs.delta ) return false; @@ -374,10 +374,10 @@ FourPCS::IsTransfCongruent(FourPoints fp,vcg::Matrix44 & m for(int i = 0 ; i < 4; ++i) fix.push_back(fp[i]); vcg::Point3 n,p; - n = (( B[1]-B[0]).Normalize() ^ ( B[2]- B[0]).Normalize())*( B[1]- B[0]).Norm(); + n = (( B[1]-B[0]).normalized() ^ ( B[2]- B[0]).normalized())*( B[1]- B[0]).Norm(); p = B[0] + n; mov.push_back(p); - n = (( fp[1]-fp[0]).Normalize() ^ (fp[2]- fp[0]).Normalize())*( fp[1]- fp[0]).Norm(); + n = (( fp[1]-fp[0]).normalized() ^ (fp[2]- fp[0]).normalized())*( fp[1]- fp[0]).Norm(); p = fp[0] + n; fix.push_back(p); @@ -565,7 +565,7 @@ int FourPCS::EvaluateSample(CandiType & fp, CoordType & tp, CoordType >(*Q,ugridQ,vq,radius, dist ); if(v!=0) - if( v->N() * np -angle >0) return 1; else return -1; + if( v->N().dot(np) -angle >0) return 1; else return -1; } diff --git a/vcg/complex/trimesh/clean.h b/vcg/complex/trimesh/clean.h index 27c0393f..a3f4f30a 100644 --- a/vcg/complex/trimesh/clean.h +++ b/vcg/complex/trimesh/clean.h @@ -1152,7 +1152,7 @@ static bool TestIntersection(FaceType *f0,FaceType *f1) //no adiacent faces if ( (f0!=f1) && (!ShareEdge(f0,f1)) && (!ShareVertex(f0,f1)) ) - return (vcg::Intersection((*f0),(*f1))); + return (vcg::Intersection_((*f0),(*f1))); return false; } diff --git a/vcg/complex/trimesh/clustering.h b/vcg/complex/trimesh/clustering.h index 7becb097..f7e0226d 100644 --- a/vcg/complex/trimesh/clustering.h +++ b/vcg/complex/trimesh/clustering.h @@ -8,7 +8,7 @@ * \ * * All rights reserved. * * * - * This program is free software; you can redistribute it and/or modify * + * 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. * @@ -101,7 +101,7 @@ class HashedPoint3i : public Point3i { public: - const size_t Hash() const + const size_t Hash() const { return (V(0)*HASH_P0 ^ V(1)*HASH_P1 ^ V(2)*HASH_P2); } @@ -129,9 +129,9 @@ class AverageCell inline void Add(MeshType &m, FaceType &f, int i) { p+=f.cV(i)->cP(); - // we prefer to use the un-normalized face normal so small faces facing away are dropped out + // we prefer to use the un-normalized face normal so small faces facing away are dropped out // and the resulting average is weighed with the size of the faces falling here. - n+=f.cN(); + n+=f.cN(); cnt++; } AverageCell(): p(0,0,0), n(0,0,0),cnt(0){} @@ -139,7 +139,7 @@ class AverageCell CoordType n; int cnt; int id; - CoordType Pos() const + CoordType Pos() const { return p/cnt; } @@ -159,9 +159,9 @@ class AverageColorCell p+=f.cV(i)->cP(); c+=CoordType(f.cV(i)->C()[0],f.cV(i)->C()[1],f.cV(i)->C()[2]); - // we prefer to use the un-normalized face normal so small faces facing away are dropped out + // we prefer to use the un-normalized face normal so small faces facing away are dropped out // and the resulting average is weighed with the size of the faces falling here. - n+=f.cN(); + n+=f.cN(); cnt++; } AverageColorCell(): p(0,0,0), n(0,0,0), c(0,0,0),cnt(0){} @@ -170,12 +170,12 @@ class AverageColorCell CoordType c; int cnt; int id; - Color4b Col() const + Color4b Col() const { return Color4b(c[0]/cnt,c[1]/cnt,c[2]/cnt,255); } - CoordType Pos() const + CoordType Pos() const { return p/cnt; } @@ -187,7 +187,7 @@ class AverageColorCell */ template class Clustering -{ +{ public: typedef typename MeshType::ScalarType ScalarType; typedef typename MeshType::CoordType CoordType; @@ -198,7 +198,7 @@ class Clustering typedef typename MeshType::FaceIterator FaceIterator; // DuplicateFace == bool means that during the clustering doublesided surface (like a thin shell) that would be clustered to a single surface - // will be merged into two identical but opposite faces. + // will be merged into two identical but opposite faces. // So in practice: // DuplicateFace=true a model with looks ok if you enable backface culling // DuplicateFace=false a model with looks ok if you enable doublesided lighting and disable backfaceculling @@ -207,7 +207,7 @@ class Clustering class SimpleTri { - public: + public: CellType *v[3]; const int ii(int i) const {return *((int *)(&(v[i])));} bool operator < ( const SimpleTri &p) const { @@ -218,7 +218,7 @@ class Clustering // Sort the vertex of the face maintaining the original face orientation (it only ensure that v0 is the minimum) void sortOrient() - { + { if(v[1] < v[0] && v[1] < v[2] ) { std::swap(v[0],v[1]); std::swap(v[1],v[2]); return; } // v1 was the minimum if(v[2] < v[0] && v[2] < v[1] ) { std::swap(v[0],v[2]); std::swap(v[1],v[2]); return; } // v2 was the minimum return; // v0 was the minimum; @@ -238,13 +238,13 @@ class Clustering // The init function Take two parameters - // _size is the approximate total number of cells composing the grid surrounding the objects (usually a large number) + // _size is the approximate total number of cells composing the grid surrounding the objects (usually a large number) // eg _size==1.000.000 means a 100x100x100 grid // _cellsize is the absolute lenght of the edge of the grid cell. // eg _cellsize==2.0 means that all the vertexes in a 2.0x2.0x2.0 cell are clustered togheter // Notes: - // _size is used only if the cell edge IS zero. + // _size is used only if the cell edge IS zero. // _cellsize gives you an absolute measure of the maximum error introduced // during the simplification (e.g. half of the cell edge lenght) @@ -267,8 +267,8 @@ class Clustering Grid.voxel[1] = Grid.dim[1]/Grid.siz[1]; Grid.voxel[2] = Grid.dim[2]/Grid.siz[2]; } - - + + BasicGrid Grid; #ifdef _MSC_VER @@ -283,7 +283,7 @@ class Clustering #endif STDEXT::hash_map GridCell; - + void Add(MeshType &m) { FaceIterator fi; @@ -300,17 +300,17 @@ class Clustering Grid.PToIP((*fi).cV(i)->cP(), pi ); st.v[i]=&(GridCell[pi]); st.v[i]->Add(m,*(fi),i); - } + } if( (st.v[0]!=st.v[1]) && (st.v[0]!=st.v[2]) && (st.v[1]!=st.v[2]) ) { // if we allow the duplication of faces we sort the vertex only partially (to maintain the original face orientation) - if(DuplicateFaceParam) st.sortOrient(); + if(DuplicateFaceParam) st.sortOrient(); else st.sort(); TriSet.insert(st); } // printf("Inserted %8i triangles, clustered to %8i tri and %i cells\n",distance(m.face.begin(),fi),TriSet.size(),GridCell.size()); } } - + void Extract(MeshType &m) { m.Clear(); @@ -339,23 +339,23 @@ class Clustering m.face[i].V(0)=&(m.vert[(*ti).v[0]->id]); m.face[i].V(1)=&(m.vert[(*ti).v[1]->id]); m.face[i].V(2)=&(m.vert[(*ti).v[2]->id]); - // if we are merging faces even when opposite we choose + // if we are merging faces even when opposite we choose // the best orientation according to the averaged normal - if(!DuplicateFaceParam) + if(!DuplicateFaceParam) { CoordType N=vcg::Normal(m.face[i]); int badOrient=0; - if( N*(*ti).v[0]->n <0) ++badOrient; - if( N*(*ti).v[1]->n <0) ++badOrient; - if( N*(*ti).v[2]->n <0) ++badOrient; + if( N.dot((*ti).v[0]->n) <0) ++badOrient; + if( N.dot((*ti).v[1]->n) <0) ++badOrient; + if( N.dot((*ti).v[2]->n) <0) ++badOrient; if(badOrient>2) std::swap(m.face[i].V(0),m.face[i].V(1)); } i++; } - + } -}; //end class clustering +}; //end class clustering } // namespace tri } // namespace vcg diff --git a/vcg/complex/trimesh/create/resampler.h b/vcg/complex/trimesh/create/resampler.h index a2cebc1e..4d6f09a2 100644 --- a/vcg/complex/trimesh/create/resampler.h +++ b/vcg/complex/trimesh/create/resampler.h @@ -176,7 +176,7 @@ template dir.Normalize(); //direction of normal inside the mesh - if ((dir*closestNorm)<0) + if ((dir.dot(closestNorm))<0) dist=-dist; //the intersection exist return true; diff --git a/vcg/complex/trimesh/hole.h b/vcg/complex/trimesh/hole.h index 19ba3bae..286b4183 100644 --- a/vcg/complex/trimesh/hole.h +++ b/vcg/complex/trimesh/hole.h @@ -212,7 +212,7 @@ namespace vcg { void ComputeAngle() { angle=Angle(cP(2)-cP(0), cP(1)-cP(0)); - ScalarType flipAngle = n * e0.v->N(); + ScalarType flipAngle = n.dot(e0.v->N()); if(flipAngle<0) angle = (2.0 *(float)M_PI) - angle; } diff --git a/vcg/complex/trimesh/inertia.h b/vcg/complex/trimesh/inertia.h index 97183b3d..0e279156 100644 --- a/vcg/complex/trimesh/inertia.h +++ b/vcg/complex/trimesh/inertia.h @@ -8,7 +8,7 @@ * \ * * All rights reserved. * * * -* This program is free software; you can redistribute it and/or modify * +* 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. * @@ -40,17 +40,17 @@ Revision 1.13 2005/11/17 00:42:03 cignoni #define _VCG_INERTIA_ /* -The algorithm is based on a three step reduction of the volume integrals -to successively simpler integrals. The algorithm is designed to minimize -the numerical errors that can result from poorly conditioned alignment of -polyhedral faces. It is also designed for efficiency. All required volume -integrals of a polyhedron are computed together during a single walk over -the boundary of the polyhedron; exploiting common subexpressions reduces +The algorithm is based on a three step reduction of the volume integrals +to successively simpler integrals. The algorithm is designed to minimize +the numerical errors that can result from poorly conditioned alignment of +polyhedral faces. It is also designed for efficiency. All required volume +integrals of a polyhedron are computed together during a single walk over +the boundary of the polyhedron; exploiting common subexpressions reduces floating point operations. For more information, check out: -Brian Mirtich, +Brian Mirtich, ``Fast and Accurate Computation of Polyhedral Mass Properties,'' journal of graphics tools, volume 1, number 2, 1996 @@ -67,7 +67,7 @@ namespace vcg template class Inertia { - typedef InertiaMeshType MeshType; + typedef InertiaMeshType MeshType; typedef typename MeshType::VertexType VertexType; typedef typename MeshType::VertexPointer VertexPointer; typedef typename MeshType::VertexIterator VertexIterator; @@ -121,7 +121,7 @@ public: db = b1 - b0; a0_2 = a0 * a0; a0_3 = a0_2 * a0; a0_4 = a0_3 * a0; b0_2 = b0 * b0; b0_3 = b0_2 * b0; b0_4 = b0_3 * b0; - a1_2 = a1 * a1; a1_3 = a1_2 * a1; + a1_2 = a1 * a1; a1_3 = a1_2 * a1; b1_2 = b1 * b1; b1_3 = b1_2 * b1; C1 = a1 + a0; @@ -180,7 +180,7 @@ void CompFaceIntegrals(FaceType &f) Faaa = k1 * Paaa; Fbbb = k1 * Pbbb; - Fccc = -k4 * (CUBE(n[A])*Paaa + 3*SQR(n[A])*n[B]*Paab + Fccc = -k4 * (CUBE(n[A])*Paaa + 3*SQR(n[A])*n[B]*Paab + 3*n[A]*SQR(n[B])*Pabb + CUBE(n[B])*Pbbb + 3*w*(SQR(n[A])*Paa + 2*n[A]*n[B]*Pab + SQR(n[B])*Pbb) + w*w*(3*(n[A]*Pa + n[B]*Pb) + w*P1)); @@ -197,13 +197,13 @@ void Compute(MeshType &m) tri::UpdateNormals::PerFaceNormalized(m); double nx, ny, nz; - T0 = T1[X] = T1[Y] = T1[Z] - = T2[X] = T2[Y] = T2[Z] + T0 = T1[X] = T1[Y] = T1[Z] + = T2[X] = T2[Y] = T2[Z] = TP[X] = TP[Y] = TP[Z] = 0; FaceIterator fi; for (fi=m.face.begin(); fi!=m.face.end();++fi) if(!(*fi).IsD()) { FaceType &f=(*fi); - + nx = fabs(f.N()[0]); ny = fabs(f.N()[1]); nz = fabs(f.N()[2]); @@ -233,12 +233,12 @@ void Compute(MeshType &m) } ScalarType Mass() -{ +{ return static_cast(T0); } Point3 CenterOfMass() -{ +{ Point3 r; r[X] = T1[X] / T0; r[Y] = T1[Y] / T0; @@ -261,9 +261,9 @@ void InertiaTensor(Matrix33 &J ){ J[X][X] -= T0 * (r[Y]*r[Y] + r[Z]*r[Z]); J[Y][Y] -= T0 * (r[Z]*r[Z] + r[X]*r[X]); J[Z][Z] -= T0 * (r[X]*r[X] + r[Y]*r[Y]); - J[X][Y] = J[Y][X] += T0 * r[X] * r[Y]; - J[Y][Z] = J[Z][Y] += T0 * r[Y] * r[Z]; - J[Z][X] = J[X][Z] += T0 * r[Z] * r[X]; + J[X][Y] = J[Y][X] += T0 * r[X] * r[Y]; + J[Y][Z] = J[Z][Y] += T0 * r[Y] * r[Z]; + J[Z][X] = J[X][Z] += T0 * r[Z] * r[X]; } void InertiaTensor(Matrix44 &J ) @@ -284,9 +284,9 @@ void InertiaTensor(Matrix44 &J ) J[X][X] -= T0 * (r[Y]*r[Y] + r[Z]*r[Z]); J[Y][Y] -= T0 * (r[Z]*r[Z] + r[X]*r[X]); J[Z][Z] -= T0 * (r[X]*r[X] + r[Y]*r[Y]); - J[X][Y] = J[Y][X] += T0 * r[X] * r[Y]; - J[Y][Z] = J[Z][Y] += T0 * r[Y] * r[Z]; - J[Z][X] = J[X][Z] += T0 * r[Z] * r[X]; + J[X][Y] = J[Y][X] += T0 * r[X] * r[Y]; + J[Y][Z] = J[Z][Y] += T0 * r[Y] * r[Z]; + J[Z][X] = J[X][Z] += T0 * r[Z] * r[X]; } @@ -321,8 +321,8 @@ static void Covariance(const MeshType & m, vcg::Point3 & bary, vcg:: area+=vcg::DoubleArea(*fi); } bary/=area; - - + + C.SetZero(); // C as covariance of triangle (0,0,0)(1,0,0)(0,1,0) vcg::Matrix33 C0; @@ -345,9 +345,15 @@ static void Covariance(const MeshType & m, vcg::Point3 & bary, vcg:: const float da = n.Norm(); n/=da*da; - A.SetColumn(0,P1-P0); - A.SetColumn(1,P2-P0); - A.SetColumn(2,n); + #ifndef VCG_USE_EIGEN + A.SetColumn(0, P1-P0); + A.SetColumn(1, P2-P0); + A.SetColumn(2, n); + #else + A.col(0) = P1-P0; + A.col(1) = P2-P0; + A.col(2) = n; + #endif CoordType delta = P0 - bary; /* DC is calculated as integral of (A*x+delta) * (A*x+delta)^T over the triangle, diff --git a/vcg/complex/trimesh/refine.h b/vcg/complex/trimesh/refine.h index 4cb934dc..156ef8ce 100644 --- a/vcg/complex/trimesh/refine.h +++ b/vcg/complex/trimesh/refine.h @@ -118,7 +118,7 @@ struct MidPoint : public std::unary_functionV(ep.z)->P()+ep.f->V1(ep.z)->P())/2.0; if( MESH_TYPE::HasPerVertexNormal()) - nv.N()= (ep.f->V(ep.z)->N()+ep.f->V1(ep.z)->N()).Normalize(); + nv.N()= (ep.f->V(ep.z)->N()+ep.f->V1(ep.z)->N()).normalized(); if( MESH_TYPE::HasPerVertexColor()) nv.C().lerp(ep.f->V(ep.z)->C(),ep.f->V1(ep.z)->C(),.5f); diff --git a/vcg/complex/trimesh/update/curvature.h b/vcg/complex/trimesh/update/curvature.h index ac0778a2..471c884c 100644 --- a/vcg/complex/trimesh/update/curvature.h +++ b/vcg/complex/trimesh/update/curvature.h @@ -8,7 +8,7 @@ * \ * * All rights reserved. * * * -* This program is free software; you can redistribute it and/or modify * +* 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. * @@ -73,12 +73,12 @@ the vertex namespace vcg { namespace tri { -/// \ingroup trimesh +/// \ingroup trimesh /// \headerfile curvature.h vcg/complex/trimesh/update/curvature.h /// \brief Management, updating and computation of per-vertex and per-face normals. -/** +/** This class is used to compute or update the normals that can be stored in the vertex or face component of a mesh. */ @@ -98,18 +98,18 @@ public: typedef typename MeshType::CoordType CoordType; typedef typename CoordType::ScalarType ScalarType; - + private: typedef struct AdjVertex { VertexType * vert; float doubleArea; bool isBorder; }; - + public: /// \brief Compute principal direction and magniuto of curvature. - + /* Compute principal direction and magniuto of curvature as describe in the paper: @InProceedings{bb33922, @@ -144,10 +144,10 @@ public: VertexType* tempV; float totalDoubleAreaSize = 0.0f; - // compute the area of each triangle around the central vertex as well as their total area - do + // compute the area of each triangle around the central vertex as well as their total area + do { - // this bring the pos to the next triangle counterclock-wise + // this bring the pos to the next triangle counterclock-wise pos.FlipF(); pos.FlipE(); @@ -157,13 +157,13 @@ public: AdjVertex v; v.isBorder = pos.IsBorder(); - v.vert = tempV; + v.vert = tempV; v.doubleArea = vcg::DoubleArea(*pos.F()); totalDoubleAreaSize += v.doubleArea; - vertices.push_back(v); - } - while(tempV != firstV); + vertices.push_back(v); + } + while(tempV != firstV); // compute the weights for the formula computing matrix M for (int i = 0; i < vertices.size(); ++i) { @@ -190,8 +190,8 @@ public: M.SetZero(); for (int i = 0; i < vertices.size(); ++i) { CoordType edge = (central_vertex->cP() - vertices[i].vert->cP()); - float curvature = (2.0f * (central_vertex->cN() * edge) ) / edge.SquaredNorm(); - CoordType T = (Tp*edge).Normalize(); + float curvature = (2.0f * (central_vertex->cN().dot(edge)) ) / edge.SquaredNorm(); + CoordType T = (Tp*edge).normalized(); tempMatrix.ExternalProduct(T,T); M += tempMatrix * weights[i] * curvature ; } @@ -201,7 +201,7 @@ public: CoordType e1(1.0f,0.0f,0.0f); if ((e1 - central_vertex->cN()).SquaredNorm() > (e1 + central_vertex->cN()).SquaredNorm()) W = e1 - central_vertex->cN(); - else + else W = e1 + central_vertex->cN(); W.Normalize(); @@ -282,7 +282,7 @@ public: float Principal_Curvature2 = (3.0f * StMS[1][1]) - StMS[0][0]; CoordType Principal_Direction1 = T1 * c - T2 * s; - CoordType Principal_Direction2 = T1 * s + T2 * c; + CoordType Principal_Direction2 = T1 * s + T2 * c; (*vi).PD1() = Principal_Direction1; (*vi).PD2() = Principal_Direction2; @@ -291,8 +291,8 @@ public: } } } - - + + class AreaData @@ -301,7 +301,7 @@ public: float A; }; - /** Curvature meseaure as described in the paper: + /** Curvature meseaure as described in the paper: Robust principal curvatures on Multiple Scales, Yong-Liang Yang, Yu-Kun Lai, Shi-Min Hu Helmut Pottmann SGP 2004 If pointVSfaceInt==true the covariance is computed by montecarlo sampling on the mesh (faster) @@ -325,11 +325,11 @@ public: PointsGridType pGrid; // Fill the grid used - if(pointVSfaceInt){ + if(pointVSfaceInt){ area = Stat::ComputeMeshArea(m); - vcg::tri::SurfaceSampling >::Montecarlo(m,vs,1000 * area / (2*M_PI*r*r )); + vcg::tri::SurfaceSampling >::Montecarlo(m,vs,1000 * area / (2*M_PI*r*r )); vi = vcg::tri::Allocator::AddVertices(tmpM,m.vert.size()); - for(int y = 0; y < m.vert.size(); ++y,++vi) (*vi).P() = m.vert[y].P(); + for(int y = 0; y < m.vert.size(); ++y,++vi) (*vi).P() = m.vert[y].P(); pGrid.Set(tmpM.vert.begin(),tmpM.vert.end()); } else{ mGrid.Set(m.face.begin(),m.face.end()); } @@ -340,7 +340,7 @@ public: // sample the neighborhood if(pointVSfaceInt) - { + { vcg::tri::GetInSphereVertex< MeshType, PointsGridType,std::vector, @@ -356,25 +356,25 @@ public: vcg::tri::Inertia::Covariance(tmpM,_bary,A); } - Jacobi(A, eigenvalues , eigenvectors, nrot); + Jacobi(A, eigenvalues , eigenvectors, nrot); // get the estimate of curvatures from eigenvalues and eigenvectors // find the 2 most tangent eigenvectors (by finding the one closest to the normal) - int best = 0; ScalarType bestv = fabs( (*vi).cN() * eigenvectors.GetColumn(0).Normalize()); + int best = 0; ScalarType bestv = fabs( (*vi).cN().dot(eigenvectors.GetColumn(0).normalized()) ); for(int i = 1 ; i < 3; ++i){ - ScalarType prod = fabs((*vi).cN() * eigenvectors.GetColumn(i).Normalize()); + ScalarType prod = fabs((*vi).cN().dot(eigenvectors.GetColumn(i).normalized())); if( prod > bestv){bestv = prod; best = i;} } - (*vi).PD1() = eigenvectors.GetColumn( (best+1)%3).Normalize(); - (*vi).PD2() = eigenvectors.GetColumn( (best+2)%3).Normalize(); + (*vi).PD1() = eigenvectors.GetColumn( (best+1)%3).normalized(); + (*vi).PD2() = eigenvectors.GetColumn( (best+2)%3).normalized(); // project them to the plane identified by the normal vcg::Matrix33 rot; - ScalarType angle = acos((*vi).PD1()*(*vi).N()); + ScalarType angle = acos((*vi).PD1().dot((*vi).N())); rot.SetRotateRad( - (M_PI*0.5 - angle),(*vi).PD1()^(*vi).N()); (*vi).PD1() = rot*(*vi).PD1(); - angle = acos((*vi).PD2()*(*vi).N()); + angle = acos((*vi).PD2().dot((*vi).N())); rot.SetRotateRad( - (M_PI*0.5 - angle),(*vi).PD2()^(*vi).N()); (*vi).PD2() = rot*(*vi).PD2(); @@ -388,26 +388,26 @@ public: std::swap((*vi).PD1(),(*vi).PD2()); } } - + } - /// \brief Computes the discrete gaussian curvature. - + /// \brief Computes the discrete gaussian curvature. + /** For further details, please, refer to: \n - Discrete Differential-Geometry Operators for Triangulated 2-Manifolds Mark Meyer, Mathieu Desbrun, Peter Schroder, Alan H. Barr VisMath '02, Berlin -*/ - static void MeanAndGaussian(MeshType & m) +*/ + static void MeanAndGaussian(MeshType & m) { assert(HasFFAdjacency(m)); float area0, area1, area2, angle0, angle1, angle2, e01, e12, e20; FaceIterator fi; - VertexIterator vi; + VertexIterator vi; typename MeshType::CoordType e01v ,e12v ,e20v; - SimpleTempData TDAreaPtr(m.vert); - SimpleTempData TDContr(m.vert); + SimpleTempData TDAreaPtr(m.vert); + SimpleTempData TDContr(m.vert); vcg::tri::UpdateNormals::PerVertexNormalized(m); //Compute AreaMix in H (vale anche per K) @@ -425,49 +425,49 @@ public: angle0 = math::Abs(Angle( (*fi).P(1)-(*fi).P(0),(*fi).P(2)-(*fi).P(0) )); angle1 = math::Abs(Angle( (*fi).P(0)-(*fi).P(1),(*fi).P(2)-(*fi).P(1) )); angle2 = M_PI-(angle0+angle1); - + if((angle0 < M_PI/2) && (angle1 < M_PI/2) && (angle2 < M_PI/2)) // triangolo non ottuso - { + { float e01 = SquaredDistance( (*fi).V(1)->cP() , (*fi).V(0)->cP() ); float e12 = SquaredDistance( (*fi).V(2)->cP() , (*fi).V(1)->cP() ); float e20 = SquaredDistance( (*fi).V(0)->cP() , (*fi).V(2)->cP() ); - + area0 = ( e20*(1.0/tan(angle1)) + e01*(1.0/tan(angle2)) ) / 8.0; area1 = ( e01*(1.0/tan(angle2)) + e12*(1.0/tan(angle0)) ) / 8.0; area2 = ( e12*(1.0/tan(angle0)) + e20*(1.0/tan(angle1)) ) / 8.0; - + (TDAreaPtr)[(*fi).V(0)].A += area0; (TDAreaPtr)[(*fi).V(1)].A += area1; (TDAreaPtr)[(*fi).V(2)].A += area2; } else // obtuse - { + { (TDAreaPtr)[(*fi).V(0)].A += vcg::DoubleArea((*fi)) / 6.0; (TDAreaPtr)[(*fi).V(1)].A += vcg::DoubleArea((*fi)) / 6.0; - (TDAreaPtr)[(*fi).V(2)].A += vcg::DoubleArea((*fi)) / 6.0; + (TDAreaPtr)[(*fi).V(2)].A += vcg::DoubleArea((*fi)) / 6.0; } - } - + } + for(fi=m.face.begin();fi!=m.face.end();++fi) if( !(*fi).IsD() ) - { + { angle0 = math::Abs(Angle( (*fi).P(1)-(*fi).P(0),(*fi).P(2)-(*fi).P(0) )); angle1 = math::Abs(Angle( (*fi).P(0)-(*fi).P(1),(*fi).P(2)-(*fi).P(1) )); angle2 = M_PI-(angle0+angle1); - + e01v = ( (*fi).V(1)->cP() - (*fi).V(0)->cP() ) ; e12v = ( (*fi).V(2)->cP() - (*fi).V(1)->cP() ) ; e20v = ( (*fi).V(0)->cP() - (*fi).V(2)->cP() ) ; - + TDContr[(*fi).V(0)] += ( e20v * (1.0/tan(angle1)) - e01v * (1.0/tan(angle2)) ) / 4.0; TDContr[(*fi).V(1)] += ( e01v * (1.0/tan(angle2)) - e12v * (1.0/tan(angle0)) ) / 4.0; TDContr[(*fi).V(2)] += ( e12v * (1.0/tan(angle0)) - e20v * (1.0/tan(angle1)) ) / 4.0; - + (*fi).V(0)->Kg() -= angle0; (*fi).V(1)->Kg() -= angle1; (*fi).V(2)->Kg() -= angle2; - + for(int i=0;i<3;i++) { if(vcg::face::IsBorder((*fi), i)) @@ -475,7 +475,7 @@ public: CoordType e1,e2; vcg::face::Pos hp(&*fi, i, (*fi).V(i)); vcg::face::Pos hp1=hp; - + hp1.FlipV(); e1=hp1.v->cP() - hp.v->cP(); hp1.FlipV(); @@ -485,7 +485,7 @@ public: } } } - + for(vi=m.vert.begin(); vi!=m.vert.end(); ++vi) if(!(*vi).IsD() /*&& !(*vi).IsB()*/) { if((TDAreaPtr)[*vi].A<=std::numeric_limits::epsilon()) @@ -495,30 +495,30 @@ public: } else { - (*vi).Kh() = (((TDContr)[*vi]* (*vi).cN()>0)?1.0:-1.0)*((TDContr)[*vi] / (TDAreaPtr) [*vi].A).Norm(); + (*vi).Kh() = (((TDContr)[*vi].dot((*vi).cN())>0)?1.0:-1.0)*((TDContr)[*vi] / (TDAreaPtr) [*vi].A).Norm(); (*vi).Kg() /= (TDAreaPtr)[*vi].A; } } } - - + + /// \brief Update the mean and the gaussian curvature of a vertex. - + /** - The function uses the VF adiacency to walk around the vertex. + The function uses the VF adiacency to walk around the vertex. \return It will return the voronoi area around the vertex. If (norm == true) the mean and the gaussian curvature are normalized. Based on the paper "Optimizing 3d triangulations using discrete curvature analysis" */ - + static float VertexCurvature(VertexPointer v, bool norm = true) { // VFAdjacency required! assert(FaceType::HasVFAdjacency()); assert(VertexType::HasVFAdjacency()); - + VFIteratorType vfi(v); float A = 0; - + v->Kh() = 0; v->Kg() = 2 * M_PI; @@ -527,7 +527,7 @@ public: FacePointer f = vfi.F(); int i = vfi.I(); VertexPointer v0 = f->V0(i), v1 = f->V1(i), v2 = f->V2(i); - + float ang0 = math::Abs(Angle(v1->P() - v0->P(), v2->P() - v0->P() )); float ang1 = math::Abs(Angle(v0->P() - v1->P(), v2->P() - v1->P() )); float ang2 = M_PI - ang0 - ang1; @@ -544,22 +544,22 @@ public: A += (s02 * tan(ang0)) / 8.0; else // non obctuse triangle A += ((s02 / tan(ang1)) + (s01 / tan(ang2))) / 8.0; - + // gaussian curvature update v->Kg() -= ang0; // mean curvature update ang1 = math::Abs(Angle(f->N(), v1->N())); ang2 = math::Abs(Angle(f->N(), v2->N())); - v->Kh() += ( (math::Sqrt(s01) / 2.0) * ang1 + + v->Kh() += ( (math::Sqrt(s01) / 2.0) * ang1 + (math::Sqrt(s02) / 2.0) * ang2 ); } - + ++vfi; } - + v->Kh() /= 4.0f; - + if(norm) { if(A <= std::numeric_limits::epsilon()) { v->Kh() = 0; @@ -597,7 +597,7 @@ public: assert(FaceType::HasFFAdjacency()); assert(FaceType::HasFaceNormal()); - + typename MeshType::VertexIterator vi; for(vi = m.vert.begin(); vi != m.vert.end(); ++vi) @@ -616,7 +616,7 @@ public: normalized_edge/=edge_length; Point3 n1 = p.F()->cN();n1.Normalize(); Point3 n2 = p.FFlip()->cN();n2.Normalize(); - ScalarType n1n2 = (n1 ^ n2)* normalized_edge; + ScalarType n1n2 = (n1 ^ n2).dot(normalized_edge); n1n2 = math::Max(math::Min ( 1.0,n1n2),-1.0); ScalarType beta = math::Asin(n1n2); m33[0][0] += beta*edge_length*normalized_edge[0]*normalized_edge[0]; @@ -645,7 +645,11 @@ public: ScalarType normal = std::numeric_limits::min(); int normI = 0; for(int i = 0; i < 3; ++i) - if( fabs((*vi).N().Normalize() * vect.GetRow(i)) > normal ){normal= fabs((*vi).N().Normalize() * vect.GetRow(i)); normI = i;} + if( fabs((*vi).N().Normalize().dot(vect.GetRow(i))) > normal ) + { + normal= fabs((*vi).N().Normalize().dot(vect.GetRow(i))); + normI = i; + } int maxI = (normI+2)%3; int minI = (normI+1)%3; if(fabs(lambda[maxI]) < fabs(lambda[minI])) std::swap(maxI,minI); diff --git a/vcg/math/base.h b/vcg/math/base.h index c51f398c..8327d976 100644 --- a/vcg/math/base.h +++ b/vcg/math/base.h @@ -8,7 +8,7 @@ * \ * * All rights reserved. * * * -* This program is free software; you can redistribute it and/or modify * +* 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. * @@ -109,16 +109,16 @@ namespace vcg { namespace math { - template + template class MagnitudoComparer { public: inline bool operator() ( const SCALAR a, const SCALAR b ) { return fabs(a)>fabs(b); } }; - + inline float Sqrt(const short v) { return sqrtf(v); } inline float Sqrt(const int v) { return sqrtf((float)v); } - + inline float Sqrt(const float v) { return sqrtf(v); } inline float Abs(const float v) { return fabsf(v); } inline float Cos(const float v) { return cosf(v); } @@ -134,7 +134,9 @@ namespace math { inline double Acos(const double v) { return acos(v); } inline double Asin(const double v) { return asin(v); } inline double Atan2(const double v0,const double v1) { return atan2(v0,v1); } - + + template inline static T Sqr(T a) { return a*a; } + template inline const T & Min(const T &a, const T &b){ if (a inline void Sort(T &a, T &b, T &c){ if (a>b) Swap(a,b); if (b>c) {Swap(b,c); if (a>b) Swap(a,b);} - } + } /* Some files do not define M_PI... */ #ifndef M_PI @@ -161,8 +163,8 @@ namespace math { #ifndef SQRT_TWO #define SQRT_TWO 1.4142135623730950488 #endif - -template + +template inline SCALAR Clamp( const SCALAR & val, const SCALAR& minval, const SCALAR& maxval) { if(val < minval) return minval; @@ -191,7 +193,7 @@ template int IsNAN(T t) return t != t; } -#endif +#endif } // End math namespace /// a type that stands for "void". Useful for Parameter type of a point. diff --git a/vcg/math/eigen.h b/vcg/math/eigen.h index 34a81b9b..642616df 100644 --- a/vcg/math/eigen.h +++ b/vcg/math/eigen.h @@ -24,6 +24,8 @@ #ifndef EIGEN_VCGLIB #define EIGEN_VCGLIB +// TODO enable the vectorization +#define EIGEN_DONT_VECTORIZE #define EIGEN_MATRIXBASE_PLUGIN #include "../Eigen/LU" diff --git a/vcg/math/eigen_vcgaddons.h b/vcg/math/eigen_vcgaddons.h index f9211b25..9be8b7e6 100644 --- a/vcg/math/eigen_vcgaddons.h +++ b/vcg/math/eigen_vcgaddons.h @@ -197,23 +197,6 @@ EIGEN_DEPRECATED Derived& operator-=(const Scalar k) // } // }; -/*! -* \deprecated use (*this) * vec.asDiagonal() or (*this) * mat.mark() -* Matrix multiplication by a diagonal matrix -*/ -// EIGEN_DEPRECATED Matrix operator*(const MatrixDiagBase &m) const -// { -// assert(_columns == _rows); -// assert(_columns == m.Dimension()); -// int i,j; -// Matrix result(_rows, _columns); -// -// for (i=0; i >(v,cols(),1); }; +/** \deprecated use *this.col(i) = other */ +template +EIGEN_DEPRECATED void SetColumn(unsigned int j, const MatrixBase& other) +{ + col(j) = other; +}; + /*! * \deprecated use *this.row(i) = expression * Set the elements of the i-th row to v[j] * \param i the row index * \param v ... */ -EIGEN_DEPRECATED void SetRow(const unsigned int i, Scalar* v) +EIGEN_DEPRECATED void SetRow(unsigned int i, Scalar* v) { row(i) = Map >(v,1,rows()); }; +/** \deprecated use *this.row(i) = other */ +template +EIGEN_DEPRECATED void SetRow(unsigned int j, const MatrixBase& other) +{ + row(j) = other; +}; + /*! * \deprecated use *this.diagonal() = expression * Set the diagonal elements vi,i to v[i] @@ -333,5 +330,5 @@ EIGEN_DEPRECATED inline Scalar SquaredNorm() const { return norm2(); } EIGEN_DEPRECATED inline Derived& Normalize() { normalize(); return derived(); } /** \deprecated use .cross(p) */ -inline EvalType operator ^ (const Derived& p ) const { return this->cross(p); } +EIGEN_DEPRECATED inline EvalType operator ^ (const Derived& p ) const { return this->cross(p); } diff --git a/vcg/math/lin_algebra.h b/vcg/math/lin_algebra.h index b9e07ecc..b4032571 100644 --- a/vcg/math/lin_algebra.h +++ b/vcg/math/lin_algebra.h @@ -66,109 +66,109 @@ namespace vcg typename MATRIX_TYPE::ScalarType g=A[i][j]; typename MATRIX_TYPE::ScalarType h=A[k][l]; A[i][j]=g-s*(h+g*tau); - A[k][l]=h+s*(g-h*tau); + A[k][l]=h+s*(g-h*tau); }; /*! * Computes all eigenvalues and eigenvectors of a real symmetric matrix . - * On output, elements of the input matrix above the diagonal are destroyed. - * \param d returns the eigenvalues of a. - * \param v is a matrix whose columns contain, the normalized eigenvectors - * \param nrot returns the number of Jacobi rotations that were required. + * On output, elements of the input matrix above the diagonal are destroyed. + * \param d returns the eigenvalues of a. + * \param v is a matrix whose columns contain, the normalized eigenvectors + * \param nrot returns the number of Jacobi rotations that were required. */ template - static void Jacobi(MATRIX_TYPE &w, POINT_TYPE &d, MATRIX_TYPE &v, int &nrot) - { + static void Jacobi(MATRIX_TYPE &w, POINT_TYPE &d, MATRIX_TYPE &v, int &nrot) + { typedef typename MATRIX_TYPE::ScalarType ScalarType; assert(w.RowsNumber()==w.ColumnsNumber()); int dimension = w.RowsNumber(); - int j,iq,ip,i; + int j,iq,ip,i; //assert(w.IsSymmetric()); - typename MATRIX_TYPE::ScalarType tresh, theta, tau, t, sm, s, h, g, c; - POINT_TYPE b, z; + typename MATRIX_TYPE::ScalarType tresh, theta, tau, t, sm, s, h, g, c; + POINT_TYPE b, z; v.SetIdentity(); - for (ip=0;ip4 && (float)(fabs(d[ip])+g) == (float)fabs(d[ip]) && (float)(fabs(d[iq])+g) == (float)fabs(d[iq])) - w[ip][iq]=ScalarType(0.0); - else if (fabs(w[ip][iq]) > tresh) - { - h=d[iq]-d[ip]; - if ((float)(fabs(h)+g) == (float)fabs(h)) - t=(w[ip][iq])/h; //t =1/(2#) - else - { - theta=ScalarType(0.5)*h/(w[ip][iq]); //Equation (11.1.10). - t=ScalarType(1.0)/(fabs(theta)+sqrt(ScalarType(1.0)+theta*theta)); - if (theta < ScalarType(0.0)) t = -t; - } - c=ScalarType(1.0)/sqrt(ScalarType(1.0)+t*t); - s=t*c; - tau=s/(ScalarType(1.0)+c); - h=t*w[ip][iq]; - z[ip] -= h; - z[iq] += h; - d[ip] -= h; - d[iq] += h; - w[ip][iq]=ScalarType(0.0); - for (j=0;j<=ip-1;j++) { //Case of rotations 1 <= j < p. + for (iq=ip+1;iq4 && (float)(fabs(d[ip])+g) == (float)fabs(d[ip]) && (float)(fabs(d[iq])+g) == (float)fabs(d[iq])) + w[ip][iq]=ScalarType(0.0); + else if (fabs(w[ip][iq]) > tresh) + { + h=d[iq]-d[ip]; + if ((float)(fabs(h)+g) == (float)fabs(h)) + t=(w[ip][iq])/h; //t =1/(2#) + else + { + theta=ScalarType(0.5)*h/(w[ip][iq]); //Equation (11.1.10). + t=ScalarType(1.0)/(fabs(theta)+sqrt(ScalarType(1.0)+theta*theta)); + if (theta < ScalarType(0.0)) t = -t; + } + c=ScalarType(1.0)/sqrt(ScalarType(1.0)+t*t); + s=t*c; + tau=s/(ScalarType(1.0)+c); + h=t*w[ip][iq]; + z[ip] -= h; + z[iq] += h; + d[ip] -= h; + d[iq] += h; + w[ip][iq]=ScalarType(0.0); + for (j=0;j<=ip-1;j++) { //Case of rotations 1 <= j < p. JacobiRotate(w,s,tau,j,ip,j,iq) ; - } - for (j=ip+1;j<=iq-1;j++) { //Case of rotations p < j < q. + } + for (j=ip+1;j<=iq-1;j++) { //Case of rotations p < j < q. JacobiRotate(w,s,tau,ip,j,j,iq); - } - for (j=iq+1;j(w,s,tau,ip,j,iq,j); - } - for (j=0;j(v,s,tau,j,ip,j,iq); - } - ++nrot; - } - } - } - for (ip=0;ip= p) + if (eigenvalues[j] >= p) p = eigenvalues[ k=j ]; } - - if (k != i) + + if (k != i) { - eigenvalues[k] = eigenvalues[i]; // i.e. - eigenvalues[i] = p; // swaps the value of the elements i-th and k-th - - for (j=0; j inline static TYPE pythagora(TYPE a, TYPE b) { TYPE abs_a = fabs(a); TYPE abs_b = fabs(b); - if (abs_a > abs_b) + if (abs_a > abs_b) return abs_a*sqrt((TYPE)1.0+sqr(abs_b/abs_a)); - else + else return (abs_b == (TYPE)0.0 ? (TYPE)0.0 : abs_b*sqrt((TYPE)1.0+sqr(abs_a/abs_b))); }; @@ -243,12 +243,12 @@ namespace vcg } /*! - * + * */ enum SortingStrategy {LeaveUnsorted=0, SortAscending=1, SortDescending=2}; template< typename MATRIX_TYPE > void Sort(MATRIX_TYPE &U, typename MATRIX_TYPE::ScalarType W[], MATRIX_TYPE &V, const SortingStrategy sorting) ; - + /*! * Given a matrix Amxn, this routine computes its singular value decomposition, @@ -258,7 +258,7 @@ namespace vcg * \param W the diagonal matrix of singular values W, stored as a vector W[1...N] * \param V the matrix V (not the transpose VT) * \param max_iters max iteration number (default = 30). - * \return + * \return */ template static bool SingularValueDecomposition(MATRIX_TYPE &A, typename MATRIX_TYPE::ScalarType *W, MATRIX_TYPE &V, const SortingStrategy sorting=LeaveUnsorted, const int max_iters=30) @@ -271,20 +271,20 @@ namespace vcg bool convergence = true; rv1 = new ScalarType[n]; - g = scale = anorm = 0; + g = scale = anorm = 0; // Householder reduction to bidiagonal form. - for (i=0; i( sqrt(s), f ); h = f*g - s; A[i][i]=f-g; - for (j=l; j(sqrt(s),f); h = f*g - s; A[i][l] = f-g; - for (k=l; k=0; i--) - { + for (i=(n-1); i>=0; i--) + { //Accumulation of right-hand transformations. - if (i < (n-1)) + if (i < (n-1)) { - if (g) + if (g) { for (j=l; j=0; i--) + for (i=math::Min(m,n)-1; i>=0; i--) { l = i+1; g = W[i]; - for (j=l; j=0; k--) - { - for (its=1; its<=max_iters; its++) + for (k=(n-1); k>=0; k--) + { + for (its=1; its<=max_iters; its++) { flag=1; - for (l=k; l>=0; l--) - { + for (l=k; l>=0; l--) + { // Test for splitting. - nm=l-1; + nm=l-1; // Note that rv1[1] is always zero. - if ((double)(fabs(rv1[l])+anorm) == anorm) + if ((double)(fabs(rv1[l])+anorm) == anorm) { flag=0; break; } - if ((double)(fabs(W[nm])+anorm) == anorm) + if ((double)(fabs(W[nm])+anorm) == anorm) break; } - if (flag) + if (flag) { c=0.0; //Cancellation of rv1[l], if l > 1. s=1.0; - for (i=l ;i<=k; i++) + for (i=l ;i<=k; i++) { f = s*rv1[i]; rv1[i] = c*rv1[i]; - if ((double)(fabs(f)+anorm) == anorm) + if ((double)(fabs(f)+anorm) == anorm) break; g = W[i]; h = pythagora(f,g); @@ -424,7 +424,7 @@ namespace vcg h = (ScalarType)1.0/h; c = g*h; s = -f*h; - for (j=0; j(f,1.0); f=((x-z)*(x+z) + h*((y/(f+sign(g,f)))-h))/x; - c=s=1.0; + c=s=1.0; //Next QR transformation: - for (j=l; j<= nm;j++) + for (j=l; j<= nm;j++) { i = j+1; g = rv1[i]; @@ -472,7 +472,7 @@ namespace vcg g = g*c - x*s; h = y*s; y *= c; - for (jj=0; jj(f,h); W[j] = z; // Rotation can be arbitrary if z = 0. - if (z) + if (z) { z = (ScalarType)1.0/z; c = f*z; @@ -490,7 +490,7 @@ namespace vcg } f = c*g + s*y; x = c*y - s*g; - for (jj=0; jj - void Sort(MATRIX_TYPE &U, typename MATRIX_TYPE::ScalarType W[], MATRIX_TYPE &V, const SortingStrategy sorting) + void Sort(MATRIX_TYPE &U, typename MATRIX_TYPE::ScalarType W[], MATRIX_TYPE &V, const SortingStrategy sorting) { typedef typename MATRIX_TYPE::ScalarType ScalarType; assert(U.ColumnsNumber()==V.ColumnsNumber()); - int mu = U.RowsNumber(); - int mv = V.RowsNumber(); + int mu = U.RowsNumber(); + int mv = V.RowsNumber(); int n = U.ColumnsNumber(); - - //ScalarType* u = &U[0][0]; + + //ScalarType* u = &U[0][0]; //ScalarType* v = &V[0][0]; for (int i=0; i p) - { - k = j; - p = W[j]; - } + { + if (W[j] > p) + { + k = j; + p = W[j]; + } } break; } @@ -572,24 +572,24 @@ namespace vcg //ScalarType* ujk = u + k; // ujk = &U[0][k] //ScalarType* vji = v + i; // vji = &V[0][i] //ScalarType* vjk = v + k; // vjk = &V[0][k] - //if (j) - //{ + //if (j) + //{ // for(;;) for( ; j!=0; --j, uji+=n, ujk+=n) // { { // p = *uji; p = *uji; // i.e. // *uji = *ujk; *uji = *ujk; // swap( U[s][i], U[s][k] ) // *ujk = p; *ujk = p; // // if (!(--j)) } - // break; - // uji += n; - // ujk += n; - // } - //} + // break; + // uji += n; + // ujk += n; + // } + //} for(int s=0; j!=0; ++s, --j) std::swap(U[s][i], U[s][k]); j = mv; - //if (j!=0) + //if (j!=0) //{ // for(;;) for ( ; j!=0; --j, vji+=n, ujk+=n) // { { @@ -598,7 +598,7 @@ namespace vcg // *vjk = p; *vjk = p; // // if (!(--j)) } // break; - // vji += n; + // vji += n; // vjk += n; // } //} @@ -610,7 +610,7 @@ namespace vcg /*! - * Solves AxX = B for a vector X, where A is specified by the matrices Umxn, + * Solves AxX = B for a vector X, where A is specified by the matrices Umxn, * Wnx1 and Vnxn as returned by SingularValueDecomposition. * No input quantities are destroyed, so the routine may be called sequentially with different bxs. * \param x is the output solution vector (xnx1) @@ -630,20 +630,20 @@ namespace vcg ScalarType s; ScalarType *tmp = new ScalarType[columns_number]; for (j=0; jm�m matrix. - */ -class MatrixDiagBase{public: -virtual const int & Dimension()const =0; -virtual const float operator[](const int & i) const = 0; -}; - -template -class MatrixDiag: public Point, public MatrixDiagBase{ -public: - const int & Dimension() const {return N;} - MatrixDiag(const Point&p):Point(p){} -}; - /*! * This class represent a generic m�n matrix. The class is templated over the scalar type field. * @param Scalar (Templete Parameter) Specifies the ScalarType field. @@ -72,7 +57,7 @@ template class Matrix : public Eigen::Matrix<_Scalar,Eigen::Dynamic,Eigen::Dynamic,Eigen::RowMajor> // FIXME col or row major ? { typedef Eigen::Matrix<_Scalar,Eigen::Dynamic,Eigen::Dynamic,Eigen::RowMajor> _Base; - + public: _EIGEN_GENERIC_PUBLIC_INTERFACE(Matrix,_Base); diff --git a/vcg/math/matrix33.h b/vcg/math/matrix33.h index 013fdd4e..cbbf5b05 100644 --- a/vcg/math/matrix33.h +++ b/vcg/math/matrix33.h @@ -243,7 +243,7 @@ Matrix33 RotationMatrix(vcg::Point3 v0,vcg::Point3 v1,bool normalized=t v0.Normalize(); v1.Normalize(); } - S dot=(v0*v1); + S dot=v0.dot(v1); ///control if there is no rotation if (dot>((S)1-epsilon)) { diff --git a/vcg/math/quadric.h b/vcg/math/quadric.h index 117db161..e4454a97 100644 --- a/vcg/math/quadric.h +++ b/vcg/math/quadric.h @@ -8,7 +8,7 @@ * \ * * All rights reserved. * * * -* This program is free software; you can redistribute it and/or modify * +* 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. * @@ -63,7 +63,7 @@ public: ScalarType a[6]; // Matrice 3x3 simmetrica: a11 a12 a13 a22 a23 a33 ScalarType b[3]; // Vettore r3 ScalarType c; // Fattore scalare (se -1 quadrica nulla) - + inline Quadric() { c = -1; } // Necessari se si utilizza stl microsoft @@ -155,14 +155,14 @@ template /* Versione veloce */ return ResultScalarType ( - p[0]*p[0]*a[0] + 2*p[0]*p[1]*a[1] + 2*p[0]*p[2]*a[2] + p[0]*b[0] + p[0]*p[0]*a[0] + 2*p[0]*p[1]*a[1] + 2*p[0]*p[2]*a[2] + p[0]*b[0] + p[1]*p[1]*a[3] + 2*p[1]*p[2]*a[4] + p[1]*b[1] + p[2]*p[2]*a[5] + p[2]*b[2] + c); } // spostare..risolve un sistema 3x3 -template +template bool Gauss33( FLTYPE x[], FLTYPE C[3][3+1] ) { const FLTYPE keps = (FLTYPE)1e-3; @@ -227,7 +227,7 @@ bool Gauss33( FLTYPE x[], FLTYPE C[3][3+1] ) // determina il punto di errore minimo template bool Minimum(Point3 &x) -{ +{ ReturnScalarType C[3][4]; C[0][0]=a[0]; C[0][1]=a[1]; C[0][2]=a[2]; C[1][0]=a[1]; C[1][1]=a[3]; C[1][2]=a[4]; @@ -240,10 +240,10 @@ bool Minimum(Point3 &x) } // determina il punto di errore minimo usando le fun di inversione di matrice che usano svd -// Molto + lento +// Molto + lento template bool MinimumSVD(Point3 &x) -{ +{ Matrix33 C; C[0][0]=a[0]; C[0][1]=a[1]; C[0][2]=a[2]; C[1][0]=a[1]; C[1][1]=a[3]; C[1][2]=a[4]; @@ -259,7 +259,7 @@ bool MinimumSVD(Point3 &x) bool MinimumNew(Point3 &x) const -{ +{ ScalarType c0=-b[0]/2; ScalarType c1=-b[1]/2; ScalarType c2=-b[2]/2; @@ -279,12 +279,12 @@ bool MinimumNew(Point3 &x) const } // determina il punto di errore minimo vincolato nel segmento (a,b) bool Minimum(Point3 &x,Point3 &pa,Point3 &pb){ -ScalarType t1,t2, t4, t5, t8, t9, +ScalarType t1,t2, t4, t5, t8, t9, t11,t12,t14,t15,t17,t18,t25,t26,t30,t34,t35, t41,t42,t44,t45,t50,t52,t54, t56,t21,t23,t37,t64,lambda; - t1 = a[4]*pb.z(); + t1 = a[4]*pb.z(); t2 = t1*pa.y(); t4 = a[1]*pb.y(); t5 = t4*pa.x(); @@ -320,14 +320,14 @@ ScalarType t1,t2, t4, t5, t8, t9, if(lambda<0) lambda=0; else if(lambda>1) lambda = 1; - x = pa*(1.0-lambda)+pb*lambda; + x = pa*(1.0-lambda)+pb*lambda; return true; } void operator *= ( const ScalarType & w ) // Amplifica una quadirca { assert( IsValid() ); - + a[0] *= w; a[1] *= w; a[2] *= w; diff --git a/vcg/simplex/vertex/distance.h b/vcg/simplex/vertex/distance.h index 9a972b69..1867be72 100644 --- a/vcg/simplex/vertex/distance.h +++ b/vcg/simplex/vertex/distance.h @@ -92,7 +92,7 @@ class PointNormalDistanceFunctor { inline bool operator () (const VERTEXTYPE & v, const VERTEXTYPE & vp, SCALARTYPE & minDist, Point3 & q) { float h = vcg::Distance(v.cP(),vp.P()) ; - float dev = InterPoint ()* ( pow((ScalarType) (1-v.cN()*vp.cN()),(ScalarType) Beta()) / (Gamma()*h +0.1)); + float dev = InterPoint() * ( pow((ScalarType) (1-v.cN().dot(vp.cN())), (ScalarType)Beta()) / (Gamma()*h +0.1)); if(h+dev < minDist){ minDist = h+dev; q = v.P(); diff --git a/vcg/space/color4.h b/vcg/space/color4.h index 3ec05cf0..49b2acd6 100644 --- a/vcg/space/color4.h +++ b/vcg/space/color4.h @@ -8,7 +8,7 @@ * \ * * All rights reserved. * * * -* This program is free software; you can redistribute it and/or modify * +* 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. * @@ -84,28 +84,28 @@ namespace vcg { /** The templated class for representing 4 entity color. The class is templated over the ScalarType. class that is used to represent color with float or with unsigned chars. All the usual - operator overloading (* + - ...) is present. + operator overloading (* + - ...) is present. */ -template +template class Color4 : public Point4 { typedef Point4 Base; public: - /// Constant for storing standard colors. + /// Constant for storing standard colors. /// Each color is stored in a simple in so that the bit pattern match with the one of Color4b. enum ColorConstant { Black =0xff000000, Gray =0xff808080, White =0xffffffff, - + Red =0xff0000ff, Green =0xff00ff00, Blue =0xffff0000, - + Cyan =0xffffff00, Yellow =0xff00ffff, Magenta=0xffff00ff, - + LightGray =0xffc0c0c0, LightRed =0xff8080ff, LightGreen =0xff80ff80, @@ -116,7 +116,7 @@ public: DarkGreen =0xff004000, DarkBlue =0xff400000 }; - + inline Color4 ( const T nx, const T ny, const T nz , const T nw ) :Point4(nx,ny,nz,nw) {}; // inline Color4 ( Color4 &c) :Point4(c) {}; inline Color4 ( const Point4 &c) :Point4(c) {}; @@ -129,26 +129,26 @@ public: // TODO make sure the types are the same } #endif - - template + + template inline void Import(const Color4 & b ) { - Point4::V()[0] = T(b[0]); - Point4::V()[1] = T(b[1]); - Point4::V()[2] = T(b[2]); - Point4::V()[3] = T(b[3]); + (*this)[0] = T(b[0]); + (*this)[1] = T(b[1]); + (*this)[2] = T(b[2]); + (*this)[3] = T(b[3]); } - template + template inline void Import(const Point4 & b ) - { - Point4::V()[0] = T(b[0]); - Point4::V()[1] = T(b[1]); - Point4::V()[2] = T(b[2]); - Point4::V()[3] = T(b[3]); + { + (*this)[0] = T(b[0]); + (*this)[1] = T(b[1]); + (*this)[2] = T(b[2]); + (*this)[3] = T(b[3]); } - template + template static inline Color4 Construct( const Color4 & b ) { return Color4(T(b[0]),T(b[1]),T(b[2]),T(b[3])); @@ -156,13 +156,13 @@ public: //inline void Import(const Color4 &b); //inline void Import(const Color4 &b); - + inline Color4 operator + ( const Color4 & p) const - { + { return Color4( (*this)[0]+p.V()[0], (*this)[1]+p.V()[1], (*this)[2]+p.V()[2], (*this)[3]+p.V()[3] ); } - - + + inline void lerp(const Color4 &c0, const Color4 &c1, const float x); inline void lerp(const Color4 &c0, const Color4 &c1, const Color4 &c2, const Point3f &ip); /// given a float and a range set the corresponding color in the well known red->green->blue color ramp. To reverse the direction of the ramp just swap minf and maxf. @@ -203,10 +203,10 @@ public: case 4: r=t; g=p; b=v; break; case 5: r=v; g=p; b=q; break; } - Point4::V()[0]=(unsigned char)(255*r); - Point4::V()[1]=(unsigned char)(255*g); - Point4::V()[2]=(unsigned char)(255*b); - Point4::V()[3]=255; + (*this)[0]=(unsigned char)(255*r); + (*this)[1]=(unsigned char)(255*g); + (*this)[2]=(unsigned char)(255*b); + (*this)[3]=255; // V()[0]=r*256;V()[1]=g*256;V()[2]=b*256; } @@ -221,9 +221,9 @@ inline void SetGrayShade(float f) } -/** Given an integer returns a well ordering of colors +/** Given an integer returns a well ordering of colors // so that every color differs as much as possible form the previous one -// params: +// params: // n is the maximum expected value (max of the range) // v is the requested position */ @@ -239,7 +239,7 @@ inline static Color4 Scatter(int n, int a,float Sat=.3f,float Val=.9f) a -= (m+1)>>1; m >>= 1; } - else m = (m+1)>>1; + else m = (m+1)>>1; if (r>n-b) r = n-b; //TRACE("Scatter range 0..%i, in %i out %i\n",n,a,b); @@ -265,7 +265,7 @@ template inline void Color4::lerp(const Color4 &c0, const Color4 &c1, const Color4 &c2, const Point3f &ip) { assert(fabs(ip[0]+ip[1]+ip[2]-1)<0.00001); - + (*this)[0]=(T)(c0[0]*ip[0] + c1[0]*ip[1]+ c2[0]*ip[2]); (*this)[1]=(T)(c0[1]*ip[0] + c1[1]*ip[1]+ c2[1]*ip[2]); (*this)[2]=(T)(c0[2]*ip[0] + c1[2]*ip[1]+ c2[2]*ip[2]); @@ -276,7 +276,7 @@ inline void Color4::lerp(const Color4 &c0, const Color4 &c1, const Colo template inline void Color4::ColorRamp(const float &minf,const float &maxf ,float v ) { - if(minf>maxf) { ColorRamp(maxf,minf,maxf+(minf-v)); return; } + if(minf>maxf) { ColorRamp(maxf,minf,maxf+(minf-v)); return; } if(v < minf ) { *this=Color4(Color4::Red); return; } //the case v > maxf is handled automatically at the end of the function @@ -298,7 +298,7 @@ inline void Color4::ColorRamp(const float &minf,const float &maxf ,float v ) #if !defined(__GNUC__) || (__GNUC__ > 3) template <> #endif -template <> +template <> inline void Color4::Import(const Color4 &b) { (*this)[0]=b[0]/255.0f; @@ -334,7 +334,7 @@ inline void Color4::Import(const Point4 &b) #if !defined(__GNUC__) || (__GNUC__ > 3) template <> #endif -template <> +template <> inline Color4 Color4::Construct( const Color4 & b ) { return Color4( @@ -347,7 +347,7 @@ inline Color4 Color4::Construct( const Color4 3) template <> #endif -template <> +template <> inline Color4 Color4::Construct( const Color4 & b ) { return Color4( @@ -357,7 +357,7 @@ inline Color4 Color4::Construct( const Color4 & b ) (float)(b[3])/255.0f); } -//template +//template //inline void Color4::Import(const Color4 &b) //{ // V()[0] = T(b[0]); @@ -369,13 +369,13 @@ inline Color4 Color4::Construct( const Color4 & b ) template<> inline Color4::Color4(Color4::ColorConstant cc) { - *((int *)this )= cc; + *((int *)this )= cc; } template<> inline Color4::Color4(Color4::ColorConstant cc) { - Import(Color4((Color4::ColorConstant)cc)); + Import(Color4((Color4::ColorConstant)cc)); } inline Color4 Clamp(Color4 &c) @@ -389,8 +389,8 @@ inline Color4 Clamp(Color4 &c) template<> inline Color4 Color4::operator + ( const Color4 & p) const -{ - return Color4( +{ + return Color4( math::Clamp(int((*this)[0])+int(p[0]),0,255), math::Clamp(int((*this)[1])+int(p[1]),0,255), math::Clamp(int((*this)[2])+int(p[2]),0,255), diff --git a/vcg/space/deprecated_point2.h b/vcg/space/deprecated_point2.h new file mode 100644 index 00000000..eab330b5 --- /dev/null +++ b/vcg/space/deprecated_point2.h @@ -0,0 +1,359 @@ +/**************************************************************************** +* 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. * +* * +****************************************************************************/ +/**************************************************************************** + History + +$Log: not supported by cvs2svn $ +Revision 1.9 2006/10/07 16:51:43 m_di_benedetto +Implemented Scale() method (was only declared). + +Revision 1.8 2006/01/19 13:53:19 m_di_benedetto +Fixed product by scalar and SquaredNorm() + +Revision 1.7 2005/10/15 19:11:49 m_di_benedetto +Corrected return type in Angle() and protected member access in unary operator - + +Revision 1.6 2005/03/18 16:34:42 fiorin +minor changes to comply gcc compiler + +Revision 1.5 2004/05/10 13:22:25 cignoni +small syntax error Math -> math in Angle + +Revision 1.4 2004/04/05 11:57:32 cignoni +Add V() access function + +Revision 1.3 2004/03/10 17:42:40 tarini +Added comments (Dox) ! +Added Import(). Costruct(), ScalarType... Corrected cross prod (sign). Added Angle. Now using Math:: stuff for trigon. etc. + +Revision 1.2 2004/03/03 15:07:40 cignoni +renamed protected member v -> _v + +Revision 1.1 2004/02/13 00:44:53 cignoni +First commit... + + +****************************************************************************/ + +#ifndef __VCGLIB_POINT2 +#define __VCGLIB_POINT2 + +#include +#include + +namespace vcg { + +/** \addtogroup space */ +/*@{*/ + /** + The templated class for representing a point in 2D space. + The class is templated over the ScalarType class that is used to represent coordinates. + All the usual operator overloading (* + - ...) is present. + */ +template class Point2 +{ +protected: + /// The only data member. Hidden to user. + P2ScalarType _v[2]; +public: + /// the scalar type + typedef P2ScalarType ScalarType; + enum {Dimension = 2}; + +//@{ + + /** @name Access to Coords. + access to coords is done by overloading of [] or explicit naming of coords (X,Y,) + ("p[0]" or "p.X()" are equivalent) **/ + inline const ScalarType &X() const {return _v[0];} + inline const ScalarType &Y() const {return _v[1];} + inline ScalarType &X() {return _v[0];} + inline ScalarType &Y() {return _v[1];} + inline const ScalarType * V() const + { + return _v; + } + inline ScalarType & V( const int i ) + { + assert(i>=0 && i<2); + return _v[i]; + } + inline const ScalarType & V( const int i ) const + { + assert(i>=0 && i<2); + return _v[i]; + } + inline const ScalarType & operator [] ( const int i ) const + { + assert(i>=0 && i<2); + return _v[i]; + } + inline ScalarType & operator [] ( const int i ) + { + assert(i>=0 && i<2); + return _v[i]; + } +//@} + /// empty constructor (does nothing) + inline Point2 () { } + /// x,y constructor + inline Point2 ( const ScalarType nx, const ScalarType ny ) + { + _v[0] = nx; _v[1] = ny; + } + /// copy constructor + inline Point2 ( Point2 const & p) + { + _v[0]= p._v[0]; _v[1]= p._v[1]; + } + /// copy + inline Point2 & operator =( Point2 const & p) + { + _v[0]= p._v[0]; _v[1]= p._v[1]; + return *this; + } + /// sets the point to (0,0) + inline void SetZero() + { _v[0] = 0;_v[1] = 0;} + /// dot product + inline ScalarType operator * ( Point2 const & p ) const + { + return ( _v[0]*p._v[0] + _v[1]*p._v[1] ); + } + /// cross product + inline ScalarType operator ^ ( Point2 const & p ) const + { + return _v[0]*p._v[1] - _v[1]*p._v[0]; + } +//@{ + /** @name Linearity for 2d points (operators +, -, *, /, *= ...) **/ + inline Point2 operator + ( Point2 const & p) const + { + return Point2( _v[0]+p._v[0], _v[1]+p._v[1] ); + } + inline Point2 operator - ( Point2 const & p) const + { + return Point2( _v[0]-p._v[0], _v[1]-p._v[1] ); + } + inline Point2 operator * ( const ScalarType s ) const + { + return Point2( _v[0] * s, _v[1] * s ); + } + inline Point2 operator / ( const ScalarType s ) const + { + return Point2( _v[0] / s, _v[1] / s ); + } + inline Point2 & operator += ( Point2 const & p) + { + _v[0] += p._v[0]; _v[1] += p._v[1]; + return *this; + } + inline Point2 & operator -= ( Point2 const & p) + { + _v[0] -= p._v[0]; _v[1] -= p._v[1]; + return *this; + } + inline Point2 & operator *= ( const ScalarType s ) + { + _v[0] *= s; _v[1] *= s; + return *this; + } + inline Point2 & operator /= ( const ScalarType s ) + { + _v[0] /= s; _v[1] /= s; + return *this; + } + //@} + /// returns the norm (Euclidian) + inline ScalarType Norm( void ) const + { + return math::Sqrt( _v[0]*_v[0] + _v[1]*_v[1] ); + } + /// returns the squared norm (Euclidian) + inline ScalarType SquaredNorm( void ) const + { + return ( _v[0]*_v[0] + _v[1]*_v[1] ); + } + inline Point2 & Scale( const ScalarType sx, const ScalarType sy ) + { + _v[0] *= sx; + _v[1] *= sy; + return * this; + } + /// normalizes, and returns itself as result + inline Point2 & Normalize( void ) + { + ScalarType n = math::Sqrt(_v[0]*_v[0] + _v[1]*_v[1]); + if(n>0.0) { _v[0] /= n; _v[1] /= n; } + return *this; + } + /// points equality + inline bool operator == ( Point2 const & p ) const + { + return (_v[0]==p._v[0] && _v[1]==p._v[1]); + } + /// disparity between points + inline bool operator != ( Point2 const & p ) const + { + return ( (_v[0]!=p._v[0]) || (_v[1]!=p._v[1]) ); + } + /// lexical ordering + inline bool operator < ( Point2 const & p ) const + { + return (_v[1]!=p._v[1])?(_v[1] ( Point2 const & p ) const + { + return (_v[1]!=p._v[1])?(_v[1]>p._v[1]): + (_v[0]>p._v[0]); + } + /// lexical ordering + inline bool operator <= ( Point2 const & p ) const + { + return (_v[1]!=p._v[1])?(_v[1]< p._v[1]): + (_v[0]<=p._v[0]); + } + /// lexical ordering + inline bool operator >= ( Point2 const & p ) const + { + return (_v[1]!=p._v[1])?(_v[1]> p._v[1]): + (_v[0]>=p._v[0]); + } + /// returns the distance to another point p + inline ScalarType Distance( Point2 const & p ) const + { + return Norm(*this-p); + } + /// returns the suqared distance to another point p + inline ScalarType SquaredDistance( Point2 const & p ) const + { + return Norm2(*this-p); + } + /// returns the angle with X axis (radiants, in [-PI, +PI] ) + inline ScalarType Angle() const { + return math::Atan2(_v[1],_v[0]); + } + /// transform the point in cartesian coords into polar coords + inline Point2 & Cartesian2Polar() + { + ScalarType t = Angle(); + _v[0] = Norm(); + _v[1] = t; + return *this; + } + /// transform the point in polar coords into cartesian coords + inline Point2 & Polar2Cartesian() + { + ScalarType l = _v[0]; + _v[0] = (ScalarType)(l*math::Cos(_v[1])); + _v[1] = (ScalarType)(l*math::Sin(_v[1])); + return *this; + } + /// rotates the point of an angle (radiants, counterclockwise) + inline Point2 & Rotate( const ScalarType rad ) + { + ScalarType t = _v[0]; + ScalarType s = math::Sin(rad); + ScalarType c = math::Cos(rad); + + _v[0] = _v[0]*c - _v[1]*s; + _v[1] = t *s + _v[1]*c; + + return *this; + } + + /// Questa funzione estende il vettore ad un qualsiasi numero di dimensioni + /// paddando gli elementi estesi con zeri + inline ScalarType Ext( const int i ) const + { + if(i>=0 && i<2) return _v[i]; + else return 0; + } + /// imports from 2D points of different types + template + inline void Import( const Point2 & b ) + { + _v[0] = b.X(); _v[1] = b.Y(); + } + /// constructs a 2D points from an existing one of different type + template + static Point2 Construct( const Point2 & b ) + { + return Point2(b.X(),b.Y()); + } + + +}; // end class definition + + +template +inline T Angle( Point2 const & p0, Point2 const & p1 ) +{ + return p1.Angle() - p0.Angle(); +} + +template +inline Point2 operator - ( Point2 const & p ){ + return Point2( -p[0], -p[1] ); +} + +template +inline Point2 operator * ( const T s, Point2 const & p ){ + return Point2( p[0] * s, p[1] * s ); +} + +template +inline T Norm( Point2 const & p ){ + return p.Norm(); +} + +template +inline T SquaredNorm( Point2 const & p ){ + return p.SquaredNorm(); +} + +template +inline Point2 & Normalize( Point2 & p ){ + return p.Normalize(); +} + +template +inline T Distance( Point2 const & p1,Point2 const & p2 ){ + return Norm(p1-p2); +} + +template +inline T SquaredDistance( Point2 const & p1,Point2 const & p2 ){ + return SquaredNorm(p1-p2); +} + +typedef Point2 Point2s; +typedef Point2 Point2i; +typedef Point2 Point2f; +typedef Point2 Point2d; + +/*@}*/ +} // end namespace +#endif diff --git a/vcg/space/deprecated_point4.h b/vcg/space/deprecated_point4.h index 73375ade..e95182ec 100644 --- a/vcg/space/deprecated_point4.h +++ b/vcg/space/deprecated_point4.h @@ -8,7 +8,7 @@ * \ * * All rights reserved. * * * -* This program is free software; you can redistribute it and/or modify * +* 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. * @@ -68,8 +68,8 @@ namespace vcg { /*@{*/ /** The templated class for representing a point in 4D space. - The class is templated over the ScalarType class that is used to represent coordinates. - All the usual operator (* + - ...) are defined. + The class is templated over the ScalarType class that is used to represent coordinates. + All the usual operator (* + - ...) are defined. */ template class Point4 @@ -84,7 +84,7 @@ public: //@{ - /** @name Standard Constructors and Initializers + /** @name Standard Constructors and Initializers No casting operators have been introduced to avoid automatic unattended (and costly) conversion between different point types **/ @@ -94,15 +94,15 @@ public: _v[0] = nx; _v[1] = ny; _v[2] = nz; _v[3] = nw; } inline Point4 ( const T p[4] ) - { + { _v[0] = p[0]; _v[1]= p[1]; _v[2] = p[2]; _v[3]= p[3]; } inline Point4 ( const Point4 & p ) - { + { _v[0]= p._v[0]; _v[1]= p._v[1]; _v[2]= p._v[2]; _v[3]= p._v[3]; } - inline void Zero() - { + inline void SetZero() + { _v[0] = _v[1] = _v[2] = _v[3]= 0; } template @@ -114,7 +114,7 @@ public: _v[3] = T(b[3]); } /// constuctor that imports from different Point4 types - template + template static inline Point4 Construct( const Point4 & b ) { return Point4(T(b[0]),T(b[1]),T(b[2]),T(b[3])); @@ -124,7 +124,7 @@ public: //@{ - /** @name Data Access. + /** @name Data Access. access to data is done by overloading of [] or explicit naming of coords (x,y,z,w) **/ inline const T & operator [] ( const int i ) const @@ -155,7 +155,7 @@ public: assert(i>=0 && i<4); return _v[i]; } - /// Padding function: give a default 0 value to all the elements that are not in the [0..2] range. + /// Padding function: give a default 0 value to all the elements that are not in the [0..2] range. /// Useful for managing in a consistent way object that could have point2 / point3 / point4 inline T Ext( const int i ) const { @@ -163,12 +163,12 @@ public: else return 0; } //@} - + //@{ /** @name Linear operators and the likes **/ inline Point4 operator + ( const Point4 & p) const - { + { return Point4( _v[0]+p._v[0], _v[1]+p._v[1], _v[2]+p._v[2], _v[3]+p._v[3] ); } inline Point4 operator - ( const Point4 & p) const @@ -208,7 +208,7 @@ public: return Point4( -_v[0], -_v[1], -_v[2], -_v[3] ); } inline Point4 VectProd ( const Point4 &x, const Point4 &z ) const - { + { Point4 res; const Point4 &y = *this; @@ -223,7 +223,7 @@ public: return res; } //@} - + //@{ /** @name Norms and normalizations **/ @@ -237,7 +237,7 @@ public: { return _v[0]*_v[0] + _v[1]*_v[1] + _v[2]*_v[2] + _v[3]*_v[3]; } - /// Euclidian normalization + /// Euclidian normalization inline Point4 & Normalize() { T n = sqrt(_v[0]*_v[0] + _v[1]*_v[1] + _v[2]*_v[2] + _v[3]*_v[3] ); @@ -251,14 +251,14 @@ public: }; //@} - + //@{ /** @name Comparison operators (lexicographical order) **/ inline bool operator == ( const Point4& p ) const { return _v[0]==p._v[0] && _v[1]==p._v[1] && _v[2]==p._v[2] && _v[3]==p._v[3]; - } + } inline bool operator != ( const Point4 & p ) const { return _v[0]!=p._v[0] || _v[1]!=p._v[1] || _v[2]!=p._v[2] || _v[3]!=p._v[3]; @@ -292,7 +292,7 @@ public: (_v[0]>=p._v[0]); } //@} - + //@{ /** @name Dot products **/ @@ -326,7 +326,7 @@ public: if (exp2>exp3) { math::Swap(k2,k3); math::Swap(exp2,exp3); } return ( (k0 + k1) + k2 ) +k3; - } + } //@} @@ -371,7 +371,7 @@ template double StableDot ( Point4 const & p0, Point4 const & p1 ) { return p0.StableDot(p1); -} +} typedef Point4 Point4s; typedef Point4 Point4i; diff --git a/vcg/space/fitting3.h b/vcg/space/fitting3.h index e9c34922..939a93d8 100644 --- a/vcg/space/fitting3.h +++ b/vcg/space/fitting3.h @@ -88,7 +88,7 @@ Point3 PlaneFittingPoints( std::vector< Point3 > & samples,Plane3 &p){ d[1]=res[1][mini]; d[2]=res[2][mini]; - p.SetOffset(c*d/d.Norm()); + p.SetOffset(c.dot(d)/d.Norm()); p.SetDirection(d/d.Norm()); return eval; diff --git a/vcg/space/intersection3.h b/vcg/space/intersection3.h index 6e1ef6f4..9890759f 100644 --- a/vcg/space/intersection3.h +++ b/vcg/space/intersection3.h @@ -271,12 +271,12 @@ namespace vcg { Point3t p21 = p2-p1; Point3t p20 = p2-p0; - ScalarType delta0_p01 = p10*p1; - ScalarType delta1_p01 = -p10*p0; - ScalarType delta0_p02 = p20*p2; - ScalarType delta2_p02 = -p20*p0; - ScalarType delta1_p12 = p21*p2; - ScalarType delta2_p12 = -p21*p1; + ScalarType delta0_p01 = p10.dot(p1); + ScalarType delta1_p01 = -p10.dot(p0); + ScalarType delta0_p02 = p20.dot(p2); + ScalarType delta2_p02 = -p20.dot(p0); + ScalarType delta1_p12 = p21.dot(p2); + ScalarType delta2_p12 = -p21.dot(p1); // the closest point can be one of the vertices of the triangle if (delta1_p01<=ScalarType(0.0) && delta2_p02<=ScalarType(0.0)) { witness = p0; } @@ -284,10 +284,10 @@ namespace vcg { else if (delta0_p02<=ScalarType(0.0) && delta1_p12<=ScalarType(0.0)) { witness = p2; } else { - ScalarType temp = p10*p2; + ScalarType temp = p10.dot(p2); ScalarType delta0_p012 = delta0_p01*delta1_p12 + delta2_p12*temp; ScalarType delta1_p012 = delta1_p01*delta0_p02 - delta2_p02*temp; - ScalarType delta2_p012 = delta2_p02*delta0_p01 - delta1_p01*(p20*p1); + ScalarType delta2_p012 = delta2_p02*delta0_p01 - delta1_p01*(p20.dot(p1)); // otherwise, can be a point lying on same edge of the triangle if (delta0_p012<=ScalarType(0.0)) @@ -366,10 +366,10 @@ namespace vcg { typedef typename SEGMENTTYPE::ScalarType T; const T epsilon = T(1e-8); - T k = pl.Direction() * (sg.P1()-sg.P0()); + T k = pl.Direction().dot((sg.P1()-sg.P0())); if( (k > -epsilon) && (k < epsilon)) return false; - T r = (pl.Offset() - pl.Direction()*sg.P0())/k; // Compute ray distance + T r = (pl.Offset() - pl.Direction().dot(sg.P0()))/k; // Compute ray distance if( (r<0) || (r > 1.0)) return false; po = sg.P0()*(1-r)+sg.P1() * r; @@ -404,7 +404,7 @@ namespace vcg { /// intersection between two triangles template - inline bool Intersection(const TRIANGLETYPE & t0,const TRIANGLETYPE & t1){ + inline bool Intersection_(const TRIANGLETYPE & t0,const TRIANGLETYPE & t1){ return NoDivTriTriIsect(t0.P0(0),t0.P0(1),t0.P0(2), t1.P0(0),t1.P0(1),t1.P0(2)); } diff --git a/vcg/space/normal_extrapolation.h b/vcg/space/normal_extrapolation.h index 172b5e71..ccc1a2da 100644 --- a/vcg/space/normal_extrapolation.h +++ b/vcg/space/normal_extrapolation.h @@ -8,7 +8,7 @@ * \ * * All rights reserved. * * * -* This program is free software; you can redistribute it and/or modify * +* 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. * @@ -25,7 +25,7 @@ #define VCG_SPACE_NORMAL_EXTRAPOLATION_H #include -#include +#include #include #include #include @@ -57,7 +57,7 @@ namespace vcg typedef typename vcg::Matrix33 MatrixType; enum NormalOrientation {IsCorrect=0, MustBeFlipped=1}; - + private: /************************************************* * Inner class definitions @@ -68,10 +68,10 @@ namespace vcg // Object functor: compute the distance between a vertex and a point struct VertPointDistanceFunctor { - inline bool operator()(const VertexType &v, const CoordType &p, ScalarType &d, CoordType &q) const - { - ScalarType distance = vcg::Distance(p, v.P()); - if (distance>d) + inline bool operator()(const VertexType &v, const CoordType &p, ScalarType &d, CoordType &q) const + { + ScalarType distance = vcg::Distance(p, v.P()); + if (distance>d) return false; d = distance; @@ -95,10 +95,10 @@ namespace vcg // Object functor: compute the distance between a point and the plane struct PlanePointDistanceFunctor { - inline bool operator()(const Plane &plane, const CoordType &p, ScalarType &d, CoordType &q) const - { - ScalarType distance = vcg::Distance(p, plane.center); - if (distance>d) + inline bool operator()(const Plane &plane, const CoordType &p, ScalarType &d, CoordType &q) const + { + ScalarType distance = vcg::Distance(p, plane.center); + if (distance>d) return false; d = distance; @@ -134,7 +134,7 @@ namespace vcg VertexPointer vertex; std::vector< MSTNode* > sons; }; - + typedef std::vector< Plane > PlaneContainer; typedef typename PlaneContainer::iterator PlaneIterator; @@ -155,7 +155,7 @@ namespace vcg int percentage; char message[128]; sprintf(message, "Locating tangent planes..."); - std::vector< Plane > tangent_planes(vertex_count); + std::vector< Plane > tangent_planes(vertex_count); vcg::Octree< VertexType, ScalarType > octree_for_planes; octree_for_planes.Set( begin, end ); @@ -182,8 +182,8 @@ namespace vcg for (unsigned int n=0; ncenter; - for (int i=0; i<3; i++) - for (int j=0; j<3; j++) + for (int i=0; i<3; i++) + for (int j=0; j<3; j++) covariance_matrix[i][j]+=diff[i]*diff[j]; } @@ -195,10 +195,10 @@ namespace vcg for (int d=0; d<3; d++) plane->normal[d] = eigenvectors[d][2]; plane->normal.Normalize(); - iter->N() = plane->normal; + iter->N() = plane->normal; plane->index = int( std::distance(begin, iter) ); } - + // Step 2: build the Riemannian graph, i.e. the graph where each point is connected to the k-nearest neigbours. dataset_bb.SetNull(); PlaneIterator ePlane = tangent_planes.end(); @@ -215,7 +215,7 @@ namespace vcg for (PlaneIterator iPlane=tangent_planes.begin(); iPlane!=ePlane; iPlane++) { if (callback!=NULL && (++progress%step)==0 && (percentage=int((progress*100)/vertex_count))<100) (callback)(percentage, message); - + unsigned int kk = k; PlanePointDistanceFunctor ppdf; DummyObjectMarker dom; @@ -224,7 +224,7 @@ namespace vcg for (unsigned int n=0; nindexindex) - riemannian_graph[iPlane->index].push_back( RiemannianEdge( nearest_planes[n], 1.0f - fabs(iPlane->normal * nearest_planes[n]->normal)) ); + riemannian_graph[iPlane->index].push_back(RiemannianEdge(nearest_planes[n], 1.0f - fabs(iPlane->normal .dot(nearest_planes[n]->normal)))); } // Step 3: compute the minimum spanning tree (MST) over the Riemannian graph (we use the Kruskal algorithm) @@ -237,7 +237,7 @@ namespace vcg std::sort( E.begin(), E.end() ); vcg::DisjointSet planeset; - + for (typename std::vector< Plane >::iterator iPlane=tangent_planes.begin(); iPlane!=ePlane; iPlane++) planeset.MakeSet( &*iPlane ); @@ -261,7 +261,7 @@ namespace vcg for ( ; iMSTEdge!=eMSTEdge; iMSTEdge++) { if (callback!=NULL && (++progress%step)==0 && (percentage=int((progress*100)/mst_size))<100) (callback)(percentage, message); - + int u_index = int(iMSTEdge->u->index); int v_index = int(iMSTEdge->v->index); incident_edges[ u_index ].push_back( v_index ), @@ -271,15 +271,15 @@ namespace vcg // Traverse the incident_edges vector and build the MST VertexIterator iCurrentVertex, iSonVertex; std::vector< MSTNode > MST(vertex_count); - + typename std::vector< Plane >::iterator iFirstPlane = tangent_planes.begin(); typename std::vector< Plane >::iterator iCurrentPlane, iSonPlane; - + MSTNode *mst_root; int r_index = (root_index!=-1)? root_index : rand()*vertex_count/RAND_MAX; mst_root = &MST[ r_index ]; mst_root->parent = mst_root; //the parent of the root is the root itself - + if (orientation==MustBeFlipped) { iCurrentVertex = begin; @@ -304,7 +304,7 @@ namespace vcg std::advance((iCurrentVertex=begin), current_node_index); //retrieve the pointer to the correspective vertex current_node->vertex = &*iCurrentVertex; //and associate it to the MST node - std::vector< int >::iterator iSon = incident_edges[ current_node_index ].begin(); + std::vector< int >::iterator iSon = incident_edges[ current_node_index ].begin(); std::vector< int >::iterator eSon = incident_edges[ current_node_index ].end(); for ( ; iSon!=eSon; iSon++) { @@ -316,7 +316,7 @@ namespace vcg //std::advance((iSonVertex=begin), *iSon);//retrieve the pointer to the Vertex associated to son border.push( *iSon ); } - maxSize = vcg::math::Max(maxSize, queueSize); + maxSize = vcg::math::Max(maxSize, queueSize); } } } @@ -337,11 +337,11 @@ namespace vcg for (int s=0; ssons.size()); s++) { if (callback!=NULL && ((++progress%step)==0) && (percentage=int((maxSize-queueSize)*100/maxSize))<100) (callback)(percentage, message); - - if (current_node->vertex->N()*current_node->sons[s]->vertex->N()vertex->N().dot(current_node->sons[s]->vertex->N())sons[s]->vertex->N() *= ScalarType(-1.0f); border.push( current_node->sons[s] ); - maxSize = vcg::math::Max(maxSize, queueSize); + maxSize = vcg::math::Max(maxSize, queueSize); } } } diff --git a/vcg/space/point.h b/vcg/space/point.h index fb2f87a2..9fb74eff 100644 --- a/vcg/space/point.h +++ b/vcg/space/point.h @@ -8,7 +8,7 @@ * \ * * All rights reserved. * * * -* This program is free software; you can redistribute it and/or modify * +* 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. * @@ -62,21 +62,21 @@ Revision 1.1 2004/03/16 03:07:38 tarini namespace vcg { namespace ndim{ - -//template + +//template //class Point; /** \addtogroup space */ /*@{*/ /** The templated class for representing a point in R^N space. - The class is templated over the ScalarType class that is used to represent coordinates. + The class is templated over the ScalarType class that is used to represent coordinates. PointBase provides the interface and the common operators for points - of any dimensionality. + of any dimensionality. */ -template +template class Point { public: @@ -94,14 +94,14 @@ public: //@{ - /** @name Standard Constructors and Initializers + /** @name Standard Constructors and Initializers No casting operators have been introduced to avoid automatic unattended (and costly) conversion between different PointType types **/ inline Point () { }; // inline Point ( const S nv[N] ); - - /// Padding function: give a default 0 value to all the elements that are not in the [0..2] range. + + /// Padding function: give a default 0 value to all the elements that are not in the [0..2] range. /// Useful for managing in a consistent way object that could have point2 / point3 / point4 inline S Ext( const int i ) const { @@ -127,7 +127,7 @@ public: return p; } - /// importer for homogeneous points + /// importer for homogeneous points template inline void ImportHomo( const Point & b ) { @@ -137,7 +137,7 @@ public: _v[N-1] = 1.0; } - /// constructor for homogeneus point. + /// constructor for homogeneus point. template static inline PointType Construct( const Point & b ) { @@ -149,7 +149,7 @@ public: //@{ - /** @name Data Access. + /** @name Data Access. access to data is done by overloading of [] or explicit naming of coords (x,y,z)**/ inline S & operator [] ( const int i ) @@ -162,10 +162,10 @@ public: assert(i>=0 && i2); return _v[2]; } - /// W is in any case the last coordinate. + /// W is in any case the last coordinate. /// (in a 2D point, W() == Y(). In a 3D point, W()==Z() /// in a 4D point, W() is a separate component) inline const S &W() const { return _v[N-1]; } @@ -190,11 +190,11 @@ public: //@} //@{ - /** @name Linearity for points + /** @name Linearity for points **/ /// sets a PointType to Zero - inline void Zero() + inline void SetZero() { for(unsigned int ii = 0; ii < Dimension;++ii) V(ii) = S(); @@ -259,7 +259,7 @@ public: V(ii) *= s; return *this; } - + inline PointType operator - () const { PointType res; @@ -336,7 +336,7 @@ public: //@{ - /** @name Comparison Operators. + /** @name Comparison Operators. Lexicographic order. **/ @@ -350,10 +350,10 @@ public: //@{ - /** @name - Glocal to Local and viceversa + /** @name + Glocal to Local and viceversa (provided for uniformity with other spatial classes. trivial for points) - **/ + **/ inline PointType LocalToGlobal(ParamType p) const{ return *this; } @@ -372,9 +372,9 @@ public: -template +template class Point2 : public Point<2,S> { -public: +public: typedef S ScalarType; typedef Point2 PointType; using Point<2,S>::_v; @@ -435,13 +435,13 @@ public: inline Point2 ( const S nv[2] ){ _v[0]=nv[0]; _v[1]=nv[1]; }; - inline Point2 operator + ( Point2 const & p) const { + inline Point2 operator + ( Point2 const & p) const { return Point2( _v[0]+p._v[0], _v[1]+p._v[1]); } - inline Point2 operator - ( Point2 const & p) const { + inline Point2 operator - ( Point2 const & p) const { return Point2( _v[0]-p._v[0], _v[1]-p._v[1]); } - inline Point2 operator * ( const S s ) const { + inline Point2 operator * ( const S s ) const { return Point2( _v[0]*s, _v[1]*s ); } inline Point2 operator / ( const S s ) const { @@ -499,7 +499,7 @@ public: inline Point2 & Normalize() { PointType n = Norm(); if(n!=0.0) { n=1.0/n; _v[0]*=n; _v[1]*=n;} return *this;}; - template Point2 & Normalize(const PT &p){ + template Point2 & Normalize(const PT &p){ PointType n = Norm(); if(n!=0.0) { n=1.0/n; V(0)*=n; V(1)*=n; } return *this;}; @@ -507,10 +507,10 @@ public: if (_v[2]!=0.0) { _v[0] /= W(); W()=1.0; } return *this;}; inline S NormInfinity() const { - return math::Max( math::Abs(_v[0]), math::Abs(_v[1]) ); } + return math::Max( math::Abs(_v[0]), math::Abs(_v[1]) ); } inline S NormOne() const { - return math::Abs(_v[0])+ math::Abs(_v[1]);} + return math::Abs(_v[0])+ math::Abs(_v[1]);} inline S operator % ( Point2 const & p ) const { return _v[0] * p._v[1] - _v[1] * p._v[0]; } @@ -519,10 +519,10 @@ public: return _v[0]+_v[1];} inline S Max() const { - return math::Max( _v[0], _v[1] ); } + return math::Max( _v[0], _v[1] ); } inline S Min() const { - return math::Min( _v[0], _v[1] ); } + return math::Min( _v[0], _v[1] ); } inline int MaxI() const { return (_v[0] < _v[1]) ? 1:0; }; @@ -539,9 +539,9 @@ public: }; -template +template class Point3 : public Point<3,S> { -public: +public: typedef S ScalarType; typedef Point3 PointType; using Point<3,S>::_v; @@ -576,13 +576,13 @@ public: inline Point3 ( const S nv[3] ){ _v[0]=nv[0]; _v[1]=nv[1]; _v[2]=nv[2]; }; - inline Point3 operator + ( Point3 const & p) const{ + inline Point3 operator + ( Point3 const & p) const{ return Point3( _v[0]+p._v[0], _v[1]+p._v[1], _v[2]+p._v[2]); } - inline Point3 operator - ( Point3 const & p) const { + inline Point3 operator - ( Point3 const & p) const { return Point3( _v[0]-p._v[0], _v[1]-p._v[1], _v[2]-p._v[2]); } - inline Point3 operator * ( const S s ) const { + inline Point3 operator * ( const S s ) const { return Point3( _v[0]*s, _v[1]*s , _v[2]*s ); } inline Point3 operator / ( const S s ) const { @@ -642,13 +642,13 @@ public: (_v[1]!=p._v[1])?(_v[1]> p._v[1]) : (_v[0]>=p._v[0]); } inline PointType & Normalize() { - S n = Norm(); + S n = Norm(); if(n!=0.0) { n=S(1.0)/n; _v[0]*=n; _v[1]*=n; _v[2]*=n; } return *this;}; - template PointType & Normalize(const PT &p){ + template PointType & Normalize(const PT &p){ S n = Norm(); if(n!=0.0) { n=1.0/n; V(0)*=n; V(1)*=n; V(2)*=n; } return *this;}; @@ -658,10 +658,10 @@ public: inline S NormInfinity() const { return math::Max( math::Max( math::Abs(_v[0]), math::Abs(_v[1]) ), - math::Abs(_v[3]) ); } + math::Abs(_v[3]) ); } inline S NormOne() const { - return math::Abs(_v[0])+ math::Abs(_v[1])+math::Max(math::Abs(_v[2]));} + return math::Abs(_v[0])+ math::Abs(_v[1])+math::Max(math::Abs(_v[2]));} inline S operator % ( PointType const & p ) const { S t = (*this)*p; /* Area, general formula */ @@ -671,10 +671,10 @@ public: return _v[0]+_v[1]+_v[2];} inline S Max() const { - return math::Max( math::Max( _v[0], _v[1] ), _v[2] ); } + return math::Max( math::Max( _v[0], _v[1] ), _v[2] ); } inline S Min() const { - return math::Min( math::Min( _v[0], _v[1] ), _v[2] ); } + return math::Min( math::Min( _v[0], _v[1] ), _v[2] ); } inline int MaxI() const { int i= (_v[0] < _v[1]) ? 1:0; if (_v[i] < _v[2]) i=2; return i;}; @@ -691,16 +691,16 @@ public: frexp( double(k0), &exp0 ); frexp( double(k1), &exp1 ); frexp( double(k2), &exp2 ); - if( exp0 +template class Point4 : public Point<4,S> { public: typedef S ScalarType; @@ -727,13 +727,13 @@ public: inline Point4 ( const S nv[4] ){ _v[0]=nv[0]; _v[1]=nv[1]; _v[2]=nv[2]; _v[3]=nv[3]; }; - inline Point4 operator + ( Point4 const & p) const { + inline Point4 operator + ( Point4 const & p) const { return Point4( _v[0]+p._v[0], _v[1]+p._v[1], _v[2]+p._v[2], _v[3]+p._v[3] ); } - inline Point4 operator - ( Point4 const & p) const { + inline Point4 operator - ( Point4 const & p) const { return Point4( _v[0]-p._v[0], _v[1]-p._v[1], _v[2]-p._v[2], _v[3]-p._v[3] ); } - inline Point4 operator * ( const S s ) const { + inline Point4 operator * ( const S s ) const { return Point4( _v[0]*s, _v[1]*s , _v[2]*s , _v[3]*s ); } inline PointType operator ^ ( PointType const & p ) const { @@ -801,7 +801,7 @@ public: PointType n = Norm(); if(n!=0.0) { n=1.0/n; _v[0]*=n; _v[1]*=n; _v[2]*=n; _v[3]*=n; } return *this;}; - template PointType & Normalize(const PT &p){ + template PointType & Normalize(const PT &p){ PointType n = Norm(); if(n!=0.0) { n=1.0/n; V(0)*=n; V(1)*=n; V(2)*=n; V(3)*=n; } return *this;}; @@ -811,10 +811,10 @@ public: inline S NormInfinity() const { return math::Max( math::Max( math::Abs(_v[0]), math::Abs(_v[1]) ), - math::Max( math::Abs(_v[2]), math::Abs(_v[3]) ) ); } + math::Max( math::Abs(_v[2]), math::Abs(_v[3]) ) ); } inline S NormOne() const { - return math::Abs(_v[0])+ math::Abs(_v[1])+math::Max(math::Abs(_v[2]),math::Abs(_v[3]));} + return math::Abs(_v[0])+ math::Abs(_v[1])+math::Max(math::Abs(_v[2]),math::Abs(_v[3]));} inline S operator % ( PointType const & p ) const { S t = (*this)*p; /* Area, general formula */ @@ -824,10 +824,10 @@ public: return _v[0]+_v[1]+_v[2]+_v[3];} inline S Max() const { - return math::Max( math::Max( _v[0], _v[1] ), math::Max( _v[2], _v[3] )); } + return math::Max( math::Max( _v[0], _v[1] ), math::Max( _v[2], _v[3] )); } inline S Min() const { - return math::Min( math::Min( _v[0], _v[1] ), math::Min( _v[2], _v[3] )); } + return math::Min( math::Min( _v[0], _v[1] ), math::Min( _v[2], _v[3] )); } inline int MaxI() const { int i= (_v[0] < _v[1]) ? 1:0; if (_v[i] < _v[2]) i=2; if (_v[i] < _v[3]) i=3; @@ -872,9 +872,9 @@ template inline S AngleN( Point3 const & p1, Point3 const & p2 ) { S w = p1*p2; - if(w>1) + if(w>1) w = 1; - else if(w<-1) + else if(w<-1) w=-1; return (S) acos(w); } @@ -914,14 +914,14 @@ inline S SquaredDistance( Point const & p1,Point const & p2 ) //template //struct Point2:public Point<2,S>{ -// inline Point2(){}; +// inline Point2(){}; // inline Point2(Point<2,S> const & p):Point<2,S>(p){} ; // inline Point2( const S a, const S b):Point<2,S>(a,b){}; //}; // //template //struct Point3:public Point3 { -// inline Point3(){}; +// inline Point3(){}; // inline Point3(Point3 const & p):Point3 (p){} // inline Point3( const S a, const S b, const S c):Point3 (a,b,c){}; //}; @@ -929,7 +929,7 @@ inline S SquaredDistance( Point const & p1,Point const & p2 ) // //template //struct Point4:public Point4{ -// inline Point4(){}; +// inline Point4(){}; // inline Point4(Point4 const & p):Point4(p){} // inline Point4( const S a, const S b, const S c, const S d):Point4(a,b,c,d){}; //}; @@ -942,7 +942,7 @@ typedef Point2 Vector2s; typedef Point2 Vector2i; typedef Point2 Vector2f; typedef Point2 Vector2d; - + typedef Point3 Point3s; typedef Point3 Point3i; typedef Point3 Point3f; @@ -951,8 +951,8 @@ typedef Point3 Vector3s; typedef Point3 Vector3i; typedef Point3 Vector3f; typedef Point3 Vector3d; - - + + typedef Point4 Point4s; typedef Point4 Point4i; typedef Point4 Point4f; diff --git a/vcg/space/point2.h b/vcg/space/point2.h index ff79fbbc..1fdfe20b 100644 --- a/vcg/space/point2.h +++ b/vcg/space/point2.h @@ -8,7 +8,7 @@ * \ * * All rights reserved. * * * -* This program is free software; you can redistribute it and/or modify * +* 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. * @@ -20,282 +20,166 @@ * for more details. * * * ****************************************************************************/ -/**************************************************************************** - History -$Log: not supported by cvs2svn $ -Revision 1.9 2006/10/07 16:51:43 m_di_benedetto -Implemented Scale() method (was only declared). - -Revision 1.8 2006/01/19 13:53:19 m_di_benedetto -Fixed product by scalar and SquaredNorm() - -Revision 1.7 2005/10/15 19:11:49 m_di_benedetto -Corrected return type in Angle() and protected member access in unary operator - - -Revision 1.6 2005/03/18 16:34:42 fiorin -minor changes to comply gcc compiler - -Revision 1.5 2004/05/10 13:22:25 cignoni -small syntax error Math -> math in Angle - -Revision 1.4 2004/04/05 11:57:32 cignoni -Add V() access function - -Revision 1.3 2004/03/10 17:42:40 tarini -Added comments (Dox) ! -Added Import(). Costruct(), ScalarType... Corrected cross prod (sign). Added Angle. Now using Math:: stuff for trigon. etc. - -Revision 1.2 2004/03/03 15:07:40 cignoni -renamed protected member v -> _v - -Revision 1.1 2004/02/13 00:44:53 cignoni -First commit... - - -****************************************************************************/ +#ifndef VCG_USE_EIGEN +#include "deprecated_point2.h" +#else #ifndef __VCGLIB_POINT2 #define __VCGLIB_POINT2 -#include +#include "../math/eigen.h" #include +namespace vcg{ +template class Point2; +} + +namespace Eigen{ +template +struct ei_traits > : ei_traits > {}; +} + namespace vcg { /** \addtogroup space */ /*@{*/ /** The templated class for representing a point in 2D space. - The class is templated over the ScalarType class that is used to represent coordinates. - All the usual operator overloading (* + - ...) is present. + The class is templated over the Scalar class that is used to represent coordinates. + All the usual operator overloading (* + - ...) is present. */ -template class Point2 +template class Point2 : public Eigen::Matrix<_Scalar,2,1> { -protected: - /// The only data member. Hidden to user. - P2ScalarType _v[2]; + typedef Eigen::Matrix<_Scalar,2,1> _Base; + using _Base::coeff; + using _Base::coeffRef; + using _Base::setZero; + using _Base::data; + using _Base::V; + public: - /// the scalar type - typedef P2ScalarType ScalarType; + + _EIGEN_GENERIC_PUBLIC_INTERFACE(Point2,_Base); + typedef Scalar ScalarType; + + VCG_EIGEN_INHERIT_ASSIGNMENT_OPERATORS(Point2) + enum {Dimension = 2}; //@{ - - /** @name Access to Coords. + /** @name Access to Coords. access to coords is done by overloading of [] or explicit naming of coords (X,Y,) ("p[0]" or "p.X()" are equivalent) **/ - inline const ScalarType &X() const {return _v[0];} - inline const ScalarType &Y() const {return _v[1];} - inline ScalarType &X() {return _v[0];} - inline ScalarType &Y() {return _v[1];} - inline const ScalarType * V() const - { - return _v; - } - inline ScalarType & V( const int i ) + inline const Scalar &X() const {return data()[0];} + inline const Scalar &Y() const {return data()[1];} + inline Scalar &X() {return data()[0];} + inline Scalar &Y() {return data()[1];} + + inline Scalar & V( const int i ) { assert(i>=0 && i<2); - return _v[i]; + return data()[i]; } - inline const ScalarType & V( const int i ) const + inline const Scalar & V( const int i ) const { assert(i>=0 && i<2); - return _v[i]; - } - inline const ScalarType & operator [] ( const int i ) const - { - assert(i>=0 && i<2); - return _v[i]; - } - inline ScalarType & operator [] ( const int i ) - { - assert(i>=0 && i<2); - return _v[i]; + return data()[i]; } //@} - /// empty constructor (does nothing) + + /// empty constructor (does nothing) inline Point2 () { } - /// x,y constructor - inline Point2 ( const ScalarType nx, const ScalarType ny ) + /// x,y constructor + inline Point2 ( const Scalar nx, const Scalar ny ) : Base(nx,ny) {} + /// copy constructor + inline Point2(Point2 const & p) : Base(p) {} + template + inline Point2(const Eigen::MatrixBase& other) : Base(other) {} + + /// cross product + inline Scalar operator ^ ( Point2 const & p ) const { - _v[0] = nx; _v[1] = ny; + return data()[0]*p.data()[1] - data()[1]*p.data()[0]; } - /// copy constructor - inline Point2 ( Point2 const & p) - { - _v[0]= p._v[0]; _v[1]= p._v[1]; - } - /// copy - inline Point2 & operator =( Point2 const & p) + + inline Point2 & Scale( const Scalar sx, const Scalar sy ) { - _v[0]= p._v[0]; _v[1]= p._v[1]; - return *this; - } - /// sets the point to (0,0) - inline void Zero() - { _v[0] = 0;_v[1] = 0;} - /// dot product - inline ScalarType operator * ( Point2 const & p ) const - { - return ( _v[0]*p._v[0] + _v[1]*p._v[1] ); - } - /// cross product - inline ScalarType operator ^ ( Point2 const & p ) const - { - return _v[0]*p._v[1] - _v[1]*p._v[0]; - } -//@{ - /** @name Linearity for 2d points (operators +, -, *, /, *= ...) **/ - inline Point2 operator + ( Point2 const & p) const - { - return Point2( _v[0]+p._v[0], _v[1]+p._v[1] ); - } - inline Point2 operator - ( Point2 const & p) const - { - return Point2( _v[0]-p._v[0], _v[1]-p._v[1] ); - } - inline Point2 operator * ( const ScalarType s ) const - { - return Point2( _v[0] * s, _v[1] * s ); - } - inline Point2 operator / ( const ScalarType s ) const - { - return Point2( _v[0] / s, _v[1] / s ); - } - inline Point2 & operator += ( Point2 const & p) - { - _v[0] += p._v[0]; _v[1] += p._v[1]; - return *this; - } - inline Point2 & operator -= ( Point2 const & p) - { - _v[0] -= p._v[0]; _v[1] -= p._v[1]; - return *this; - } - inline Point2 & operator *= ( const ScalarType s ) - { - _v[0] *= s; _v[1] *= s; - return *this; - } - inline Point2 & operator /= ( const ScalarType s ) - { - _v[0] /= s; _v[1] /= s; - return *this; - } - //@} - /// returns the norm (Euclidian) - inline ScalarType Norm( void ) const - { - return math::Sqrt( _v[0]*_v[0] + _v[1]*_v[1] ); - } - /// returns the squared norm (Euclidian) - inline ScalarType SquaredNorm( void ) const - { - return ( _v[0]*_v[0] + _v[1]*_v[1] ); - } - inline Point2 & Scale( const ScalarType sx, const ScalarType sy ) - { - _v[0] *= sx; - _v[1] *= sy; + data()[0] *= sx; + data()[1] *= sy; return * this; } - /// normalizes, and returns itself as result - inline Point2 & Normalize( void ) - { - ScalarType n = math::Sqrt(_v[0]*_v[0] + _v[1]*_v[1]); - if(n>0.0) { _v[0] /= n; _v[1] /= n; } - return *this; - } - /// points equality - inline bool operator == ( Point2 const & p ) const - { - return (_v[0]==p._v[0] && _v[1]==p._v[1]); - } - /// disparity between points - inline bool operator != ( Point2 const & p ) const - { - return ( (_v[0]!=p._v[0]) || (_v[1]!=p._v[1]) ); - } - /// lexical ordering + + /// lexical ordering inline bool operator < ( Point2 const & p ) const { - return (_v[1]!=p._v[1])?(_v[1] ( Point2 const & p ) const { - return (_v[1]!=p._v[1])?(_v[1]>p._v[1]): - (_v[0]>p._v[0]); + return (data()[1]!=p.data()[1])?(data()[1]>p.data()[1]): + (data()[0]>p.data()[0]); } - /// lexical ordering + /// lexical ordering inline bool operator <= ( Point2 const & p ) const { - return (_v[1]!=p._v[1])?(_v[1]< p._v[1]): - (_v[0]<=p._v[0]); + return (data()[1]!=p.data()[1])?(data()[1]< p.data()[1]): + (data()[0]<=p.data()[0]); } - /// lexical ordering + /// lexical ordering inline bool operator >= ( Point2 const & p ) const { - return (_v[1]!=p._v[1])?(_v[1]> p._v[1]): - (_v[0]>=p._v[0]); + return (data()[1]!=p.data()[1])?(data()[1]> p.data()[1]): + (data()[0]>=p.data()[0]); } - /// returns the distance to another point p - inline ScalarType Distance( Point2 const & p ) const - { - return Norm(*this-p); - } - /// returns the suqared distance to another point p - inline ScalarType SquaredDistance( Point2 const & p ) const - { - return Norm2(*this-p); - } - /// returns the angle with X axis (radiants, in [-PI, +PI] ) - inline ScalarType Angle() const { - return math::Atan2(_v[1],_v[0]); + + /// returns the angle with X axis (radiants, in [-PI, +PI] ) + inline Scalar Angle() const { + return math::Atan2(data()[1],data()[0]); } /// transform the point in cartesian coords into polar coords inline Point2 & Cartesian2Polar() { - ScalarType t = Angle(); - _v[0] = Norm(); - _v[1] = t; + Scalar t = Angle(); + data()[0] = this->norm(); + data()[1] = t; return *this; } /// transform the point in polar coords into cartesian coords inline Point2 & Polar2Cartesian() { - ScalarType l = _v[0]; - _v[0] = (ScalarType)(l*math::Cos(_v[1])); - _v[1] = (ScalarType)(l*math::Sin(_v[1])); + Scalar l = data()[0]; + data()[0] = (Scalar)(l*math::Cos(data()[1])); + data()[1] = (Scalar)(l*math::Sin(data()[1])); return *this; } /// rotates the point of an angle (radiants, counterclockwise) - inline Point2 & Rotate( const ScalarType rad ) + inline Point2 & Rotate( const Scalar rad ) { - ScalarType t = _v[0]; - ScalarType s = math::Sin(rad); - ScalarType c = math::Cos(rad); + Scalar t = data()[0]; + Scalar s = math::Sin(rad); + Scalar c = math::Cos(rad); - _v[0] = _v[0]*c - _v[1]*s; - _v[1] = t *s + _v[1]*c; + data()[0] = data()[0]*c - data()[1]*s; + data()[1] = t *s + data()[1]*c; return *this; } /// Questa funzione estende il vettore ad un qualsiasi numero di dimensioni /// paddando gli elementi estesi con zeri - inline ScalarType Ext( const int i ) const + inline Scalar Ext( const int i ) const { - if(i>=0 && i<2) return _v[i]; + if(i>=0 && i<2) return data()[i]; else return 0; } /// imports from 2D points of different types template inline void Import( const Point2 & b ) { - _v[0] = b.X(); _v[1] = b.Y(); + data()[0] = b.X(); data()[1] = b.Y(); } /// constructs a 2D points from an existing one of different type template @@ -303,8 +187,6 @@ public: { return Point2(b.X(),b.Y()); } - - }; // end class definition @@ -314,41 +196,6 @@ inline T Angle( Point2 const & p0, Point2 const & p1 ) return p1.Angle() - p0.Angle(); } -template -inline Point2 operator - ( Point2 const & p ){ - return Point2( -p[0], -p[1] ); -} - -template -inline Point2 operator * ( const T s, Point2 const & p ){ - return Point2( p[0] * s, p[1] * s ); -} - -template -inline T Norm( Point2 const & p ){ - return p.Norm(); -} - -template -inline T SquaredNorm( Point2 const & p ){ - return p.SquaredNorm(); -} - -template -inline Point2 & Normalize( Point2 & p ){ - return p.Normalize(); -} - -template -inline T Distance( Point2 const & p1,Point2 const & p2 ){ - return Norm(p1-p2); -} - -template -inline T SquaredDistance( Point2 const & p1,Point2 const & p2 ){ - return SquaredNorm(p1-p2); -} - typedef Point2 Point2s; typedef Point2 Point2i; typedef Point2 Point2f; @@ -357,3 +204,5 @@ typedef Point2 Point2d; /*@}*/ } // end namespace #endif + +#endif diff --git a/vcg/space/point3.h b/vcg/space/point3.h index 1e338daf..a1299fcd 100644 --- a/vcg/space/point3.h +++ b/vcg/space/point3.h @@ -38,6 +38,17 @@ template class Point3; namespace Eigen{ template struct ei_traits > : ei_traits > {}; + +template +struct NumTraits > : NumTraits +{ + enum { + ReadCost = 3, + AddCost = 3, + MulCost = 3 + }; +}; + } namespace vcg { diff --git a/wrap/gl/addons.h b/wrap/gl/addons.h index d1b95425..002ff1ba 100644 --- a/wrap/gl/addons.h +++ b/wrap/gl/addons.h @@ -51,13 +51,22 @@ namespace vcg enum DrawMode {DMUser,DMWire,DMSolid} ; private: - ///used to find right trasformation in case of rotation - static void XAxis( vcg::Point3f zero, vcg::Point3f uno, Matrix44f & tr){ + ///used to find right transformation in case of rotation + static void XAxis(vcg::Point3f zero, vcg::Point3f uno, Matrix44f & tr) + { + #ifndef VCG_USE_EIGEN tr.SetZero(); *((vcg::Point3f*)&tr[0][0]) = uno-zero; GetUV(*((vcg::Point3f*)tr[0]),*((vcg::Point3f*)tr[1]),*((vcg::Point3f*)tr[2])); tr[3][3] = 1.0; *((vcg::Point3f*)&tr[3][0]) = zero; + #else + tr.col(0).start<3>().setZero(); + tr.row(0).start<3>() = (uno-zero).normalized(); // n + tr.row(1).start<3>() = tr.row(0).start<3>().unitOrthogonal(); // u + tr.row(2).start<3>() = tr.row(0).start<3>().cross(tr.row(1).start<3>()).normalized(); // v + tr.row(3) << zero.transpose(), 1.; + #endif } //set drawingmode parameters diff --git a/wrap/gl/space.h b/wrap/gl/space.h index dea392f2..d412a82d 100644 --- a/wrap/gl/space.h +++ b/wrap/gl/space.h @@ -8,7 +8,7 @@ * \ * * All rights reserved. * * * -* This program is free software; you can redistribute it and/or modify * +* 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. * @@ -60,7 +60,7 @@ First working version! #ifndef VCG_GL_SPACE_H #define VCG_GL_SPACE_H -// Please note that this file assume that you have already included your +// Please note that this file assume that you have already included your // gl-extension wrapping utility, and that therefore all the extension symbol are already defined. #include @@ -108,16 +108,16 @@ namespace vcg { inline void glColor(Color4b const & c) { glColor4ubv(c.V());} inline void glClearColor(Color4b const &c) { ::glClearColor(float(c[0])/255.0f,float(c[1])/255.0f,float(c[2])/255.0f,1.0f);} - inline void glLight(GLenum light, GLenum pname, Color4b const & c) { - static float cf[4]; - cf[0]=float(cf[0]/255.0); cf[1]=float(c[1]/255.0); cf[2]=float(c[2]/255.0); cf[3]=float(c[3]/255.0); + inline void glLight(GLenum light, GLenum pname, Color4b const & c) { + static float cf[4]; + cf[0]=float(cf[0]/255.0); cf[1]=float(c[1]/255.0); cf[2]=float(c[2]/255.0); cf[3]=float(c[3]/255.0); glLightfv(light,pname,cf); } template - inline void glBoxWire(Box3 const & b) -{ + inline void glBoxWire(Box3 const & b) +{ glPushAttrib(GL_ENABLE_BIT); glDisable(GL_LIGHTING); glBegin(GL_LINE_STRIP); @@ -137,13 +137,13 @@ namespace vcg { glBegin(GL_LINES); glVertex3f((float)b.min[0],(float)b.min[1],(float)b.min[2]); glVertex3f((float)b.min[0],(float)b.min[1],(float)b.max[2]); - + glVertex3f((float)b.max[0],(float)b.min[1],(float)b.min[2]); glVertex3f((float)b.max[0],(float)b.min[1],(float)b.max[2]); - + glVertex3f((float)b.max[0],(float)b.max[1],(float)b.min[2]); glVertex3f((float)b.max[0],(float)b.max[1],(float)b.max[2]); - + glVertex3f((float)b.min[0],(float)b.max[1],(float)b.min[2]); glVertex3f((float)b.min[0],(float)b.max[1],(float)b.max[2]); glEnd(); @@ -151,7 +151,7 @@ namespace vcg { }; template /// Funzione di utilita' per la visualizzazione in OpenGL (flat shaded) -inline void glBoxFlat(Box3 const & b) +inline void glBoxFlat(Box3 const & b) { glPushAttrib(GL_SHADE_MODEL); glShadeModel(GL_FLAT); @@ -190,10 +190,10 @@ inline void glBoxFlat(Box3 const & b) template - /// Setta i sei clip planes di opengl a far vedere solo l'interno del box -inline void glBoxClip(const Box3 & b) + /// Setta i sei clip planes di opengl a far vedere solo l'interno del box +inline void glBoxClip(const Box3 & b) { - double eq[4]; + double eq[4]; eq[0]= 1; eq[1]= 0; eq[2]= 0; eq[3]=(double)-b.min[0]; glClipPlane(GL_CLIP_PLANE0,eq); eq[0]=-1; eq[1]= 0; eq[2]= 0; eq[3]=(double) b.max[0]; @@ -211,8 +211,8 @@ inline void glBoxClip(const Box3 & b) glClipPlane(GL_CLIP_PLANE5,eq); } template - inline void glBoxWire(const Box2 & b) -{ + inline void glBoxWire(const Box2 & b) +{ glPushAttrib(GL_ENABLE_BIT); glDisable(GL_LIGHTING); glBegin(GL_LINE_LOOP); @@ -222,7 +222,7 @@ inline void glBoxClip(const Box3 & b) glVertex2f((float)b.max[0],(float)b.max[1]); glVertex2f((float)b.min[0],(float)b.max[1]); glEnd(); - + glPopAttrib(); }; template @@ -280,8 +280,21 @@ template #ifdef VCG_USE_EIGEN -#define _WRAP_EIGEN_XPR(FUNC) template \ - inline void FUNC(const Eigen::MatrixBase& p) { FUNC(p.eval()); } +template +struct EvalToKnownPointType; + +template struct EvalToKnownPointType +{ typedef Point2 Type; }; + +template struct EvalToKnownPointType +{ typedef Point3 Type; }; + +template struct EvalToKnownPointType +{ typedef Point4 Type; }; + +#define _WRAP_EIGEN_XPR(FUNC) template \ + inline void FUNC(const Eigen::MatrixBase& p) { \ + FUNC(typename EvalToKnownPointType::Type(p)); } _WRAP_EIGEN_XPR(glVertex) _WRAP_EIGEN_XPR(glNormal) diff --git a/wrap/gl/trimesh.h b/wrap/gl/trimesh.h index 1b547b99..836b1a66 100644 --- a/wrap/gl/trimesh.h +++ b/wrap/gl/trimesh.h @@ -8,7 +8,7 @@ * \ * * All rights reserved. * * * -* This program is free software; you can redistribute it and/or modify * +* 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. * @@ -123,8 +123,8 @@ first draft: it includes glew ! namespace vcg { -// classe base di glwrap usata solo per poter usare i vari drawmode, normalmode senza dover -// specificare tutto il tipo (a volte lunghissimo) +// classe base di glwrap usata solo per poter usare i vari drawmode, normalmode senza dover +// specificare tutto il tipo (a volte lunghissimo) // della particolare classe glwrap usata. class GLW { @@ -135,17 +135,17 @@ public: enum TextureMode{TMNone, TMPerVert, TMPerWedge, TMPerWedgeMulti}; enum Hint { HNUseTriStrip = 0x0001, // ha bisogno che ci sia la fftopology gia calcolata! -// HNUseEdgeStrip = 0x0002, // - HNUseDisplayList = 0x0004, +// HNUseEdgeStrip = 0x0002, // + HNUseDisplayList = 0x0004, HNCacheDisplayList = 0x0008, // Each mode has its dl; - HNLazyDisplayList = 0x0010, // Display list are generated only when requested - HNIsTwoManifold = 0x0020, // There is no need to make DetachComplex before . - HNUsePerWedgeNormal = 0x0040, // + HNLazyDisplayList = 0x0010, // Display list are generated only when requested + HNIsTwoManifold = 0x0020, // There is no need to make DetachComplex before . + HNUsePerWedgeNormal = 0x0040, // HNHasFFTopology = 0x0080, // E' l'utente che si preoccupa di tenere aggiornata la topologia ff HNHasVFTopology = 0x0100, // E' l'utente che si preoccupa di tenere aggiornata la topologia vf HNHasVertNormal = 0x0200, // E' l'utente che si preoccupa di tenere aggiornata le normali per faccia HNHasFaceNormal = 0x0400, // E' l'utente che si preoccupa di tenere aggiornata le normali per vertice - HNUseVArray = 0x0800, + HNUseVArray = 0x0800, HNUseLazyEdgeStrip = 0x1000, // Edge Strip are generated only when requested HNUseVBO = 0x2000 // Use Vertex Buffer Object }; @@ -160,7 +160,7 @@ public: CHAll = 0xff }; enum HintParami { - HNPDisplayListSize =0 + HNPDisplayListSize =0 }; enum HintParamf { HNPCreaseAngle =0, // crease angle in radians @@ -181,16 +181,16 @@ public: // GL Array Elemet class GLAElem { - public : + public : int glmode; int len; int start; }; - + }; -template > +template > class GlTrimesh : public GLW { public: @@ -207,7 +207,7 @@ public: // The parameters of hints int HNParami[8]; float HNParamf[8]; - + MESH_TYPE *m; GlTrimesh() { @@ -217,7 +217,7 @@ public: cdm=DMNone; ccm=CMNone; cnm=NMNone; - + SetHintParamf(HNPCreaseAngle,float(M_PI/5)); SetHintParamf(HNPZTwist,0.00005f); SetHintParamf(HNPPointSize,1.0f); @@ -228,12 +228,12 @@ public: //Delete the VBOs if(curr_hints&HNUseVBO) { - for(int i=0;i<3;++i) + for(int i=0;i<3;++i) if(glIsBuffer(array_buffers[i])) glDeleteBuffersARB(1, (GLuint *)(array_buffers+i)); } } - + void SetHintParami(const HintParami hip, const int value) { HNParami[hip]=value; @@ -250,16 +250,16 @@ public: { return HNParamf[hip]; } - void SetHint(Hint hn) + void SetHint(Hint hn) { curr_hints |= hn; } - void ClearHint(Hint hn) + void ClearHint(Hint hn) { curr_hints&=(~hn); } - unsigned int dl; + unsigned int dl; std::vector indices; DrawMode cdm; // Current DrawMode @@ -269,7 +269,7 @@ public: void Update(/*Change c=CHAll*/) { if(m==0) return; - + if(curr_hints&HNUseVArray || curr_hints&HNUseVBO) { typename MESH_TYPE::FaceIterator fi; @@ -285,15 +285,15 @@ void Update(/*Change c=CHAll*/) { if(!glIsBuffer(array_buffers[1])) glGenBuffers(2,(GLuint*)array_buffers); - glBindBuffer(GL_ARRAY_BUFFER,array_buffers[0]); - glBufferData(GL_ARRAY_BUFFER_ARB, m->vn * sizeof(typename MESH_TYPE::VertexType), - (char *)&(m->vert[0].P()), GL_STATIC_DRAW_ARB); + glBindBuffer(GL_ARRAY_BUFFER,array_buffers[0]); + glBufferData(GL_ARRAY_BUFFER_ARB, m->vn * sizeof(typename MESH_TYPE::VertexType), + (char *)&(m->vert[0].P()), GL_STATIC_DRAW_ARB); - glBindBuffer(GL_ARRAY_BUFFER,array_buffers[1]); - glBufferData(GL_ARRAY_BUFFER_ARB, m->vn * sizeof(typename MESH_TYPE::VertexType), + glBindBuffer(GL_ARRAY_BUFFER,array_buffers[1]); + glBufferData(GL_ARRAY_BUFFER_ARB, m->vn * sizeof(typename MESH_TYPE::VertexType), (char *)&(m->vert[0].N()), GL_STATIC_DRAW_ARB); } - + glVertexPointer(3,GL_FLOAT,sizeof(typename MESH_TYPE::VertexType),0); glNormalPointer(GL_FLOAT,sizeof(typename MESH_TYPE::VertexType),0); } @@ -316,12 +316,12 @@ void Update(/*Change c=CHAll*/) // if(!(curr_hints&HNHasVFTopology)) m->VFTopology(); // CreaseWN(*m,MESH_TYPE::scalar_type(GetHintParamf(HNPCreaseAngle))); //} - //if(C!=0) { // force the recomputation of display list + //if(C!=0) { // force the recomputation of display list // cdm=DMNone; // ccm=CMNone; // cnm=NMNone; //} - //if((curr_hints&HNUseVArray) && (curr_hints&HNUseTriStrip)) + //if((curr_hints&HNUseVArray) && (curr_hints&HNUseTriStrip)) // { // ConvertTriStrip(*m,TStrip,TStripF,TStripVED,TStripVEI); // } @@ -420,14 +420,14 @@ void DrawFill() typename FACE_POINTER_CONTAINER::iterator fp; typename MESH_TYPE::FaceIterator fi; - + typename std::vector::iterator fip; short curtexname=-1; if(cm == CMPerMesh) glColor(m->C()); - + if(tm == TMPerWedge || tm == TMPerWedgeMulti ) glDisable(GL_TEXTURE_2D); @@ -441,10 +441,10 @@ void DrawFill() if (nm==NMPerVert) { - glBindBuffer(GL_ARRAY_BUFFER,array_buffers[1]); + glBindBuffer(GL_ARRAY_BUFFER,array_buffers[1]); glNormalPointer(GL_FLOAT,sizeof(typename MESH_TYPE::VertexType),0); } - glBindBuffer(GL_ARRAY_BUFFER,array_buffers[0]); + glBindBuffer(GL_ARRAY_BUFFER,array_buffers[0]); glVertexPointer(3,GL_FLOAT,sizeof(typename MESH_TYPE::VertexType),0); glDrawElements(GL_TRIANGLES ,m->fn*3,GL_UNSIGNED_INT, &(*indices.begin()) ); @@ -458,8 +458,8 @@ void DrawFill() } } - - if(curr_hints&HNUseVArray) + + if(curr_hints&HNUseVArray) { if( (cm==CMNone) || (cm==CMPerMesh) ) { @@ -469,7 +469,7 @@ void DrawFill() if (nm==NMPerVert) glNormalPointer(GL_FLOAT,sizeof(typename MESH_TYPE::VertexType),&(m->vert.begin()->N()[0])); - glVertexPointer(3,GL_FLOAT,sizeof(typename MESH_TYPE::VertexType),&(m->vert.begin()->P()[0])); + glVertexPointer(3,GL_FLOAT,sizeof(typename MESH_TYPE::VertexType),&(m->vert.begin()->P()[0])); glDrawElements(GL_TRIANGLES ,m->fn*3,GL_UNSIGNED_INT, &(*indices.begin()) ); glDisableClientState (GL_VERTEX_ARRAY); @@ -480,8 +480,8 @@ void DrawFill() } } else - - if(curr_hints&HNUseTriStrip) + + if(curr_hints&HNUseTriStrip) { //if( (nm==NMPerVert) && ((cm==CMNone) || (cm==CMPerMesh))) // if(curr_hints&HNUseVArray){ @@ -492,7 +492,7 @@ void DrawFill() // std::vector::iterator vi; // for(vi=TStripVED.begin();vi!=TStripVED.end();++vi) // glDrawElements(vi->glmode ,vi->len,GL_UNSIGNED_SHORT,&TStripVEI[vi->start] ); - // + // // glDisableClientState (GL_NORMAL_ARRAY ); // glDisableClientState (GL_VERTEX_ARRAY); // return; @@ -519,7 +519,7 @@ void DrawFill() } else { - if(partial) + if(partial) fp = face_pointers.begin(); else fi = m->face.begin(); @@ -537,9 +537,9 @@ void DrawFill() glDisable(GL_TEXTURE_2D); } } - + glBegin(GL_TRIANGLES); - + while( (partial)?(fp!=face_pointers.end()):(fi!=m->face.end())) { typename MESH_TYPE::FaceType & f = (partial)?(*(*fp)): *fi; @@ -551,7 +551,7 @@ void DrawFill() { curtexname=(*fi).WT(0).n(); glEnd(); - + if (curtexname >= 0) { glEnable(GL_TEXTURE_2D); @@ -561,10 +561,10 @@ void DrawFill() { glDisable(GL_TEXTURE_2D); } - + glBegin(GL_TRIANGLES); } - + if(nm == NMPerFace) glNormal(f.cN()); if(nm == NMPerVert) glNormal(f.V(0)->cN()); if(nm == NMPerWedge)glNormal(f.WN(0)); @@ -595,14 +595,14 @@ void DrawFill() else ++fi; } - + glEnd(); - + } } -/// Basic Point drawing fucntion +/// Basic Point drawing fucntion // works also for mesh with deleted vertices template void DrawPointsBase() @@ -610,12 +610,12 @@ void DrawPointsBase() typename MESH_TYPE::VertexIterator vi; glBegin(GL_POINTS); if(cm==CMPerMesh) glColor(m->C()); - + for(vi=m->vert.begin();vi!=m->vert.end();++vi)if(!(*vi).IsD()) { - if(nm==NMPerVert) glNormal((*vi).cN()); - if(cm==CMPerVert) glColor((*vi).C()); - glVertex((*vi).P()); + if(nm==NMPerVert) glNormal((*vi).cN()); + if(cm==CMPerVert) glColor((*vi).C()); + glVertex((*vi).P()); } glEnd(); } @@ -628,27 +628,27 @@ double CameraDistance(){ Point3 c=m->bbox.Center(); res=mm*c; return Norm(res); -} +} template void DrawPoints() { glPointSize(GetHintParamf(HNPPointSize)); - + float camDist=CameraDistance(); float quadratic[] = { 0.0f, 0.0f, 1.0f/(camDist*camDist) }; glPointParameterfv( GL_POINT_DISTANCE_ATTENUATION, quadratic ); glPointParameterf( GL_POINT_SIZE_MAX, 16.0f ); - glPointParameterf( GL_POINT_SIZE_MIN, 1.0f ); - - if(m->vn!=(int)m->vert.size()) + glPointParameterf( GL_POINT_SIZE_MIN, 1.0f ); + + if(m->vn!=(int)m->vert.size()) { DrawPointsBase(); return; } - + // Perfect case, no deleted stuff, // draw the vertices using vertex arrays - if (nm==NMPerVert) + if (nm==NMPerVert) { glEnableClientState (GL_NORMAL_ARRAY); glNormalPointer(GL_FLOAT,sizeof(typename MESH_TYPE::VertexType),&(m->vert.begin()->N()[0])); @@ -658,17 +658,17 @@ void DrawPoints() glEnableClientState (GL_COLOR_ARRAY); glColorPointer(4,GL_UNSIGNED_BYTE,sizeof(typename MESH_TYPE::VertexType),&(m->vert.begin()->C()[0])); } - + glEnableClientState (GL_VERTEX_ARRAY); - glVertexPointer(3,GL_FLOAT,sizeof(typename MESH_TYPE::VertexType),&(m->vert.begin()->P()[0])); - + glVertexPointer(3,GL_FLOAT,sizeof(typename MESH_TYPE::VertexType),&(m->vert.begin()->P()[0])); + //glDrawElements(GL_POINTS ,m->vn,GL_UNSIGNED_INT, &(*indices.begin()) ); glDrawArrays(GL_POINTS,0,m->vn); - + glDisableClientState (GL_VERTEX_ARRAY); if (nm==NMPerVert) glDisableClientState (GL_NORMAL_ARRAY); if (cm==CMPerVert) glDisableClientState (GL_COLOR_ARRAY); - + return; } @@ -709,7 +709,7 @@ void DrawFlatWire() DrawWire(); glPopAttrib(); //glDepthRange(0,1.0f); -} +} template void DrawRadar() @@ -720,7 +720,7 @@ void DrawRadar() glDepthMask(0); glDepthRange(ZTWIST,1.0f); - if (cm == CMNone) + if (cm == CMNone) glColor4f(0.2f, 1.0f, 0.4f, 0.2f); // DrawFill(); Draw(); @@ -760,12 +760,12 @@ void DrawTexture_NPV_TPW2() glMultiTexCoordARB(GL_TEXTURE1_ARB, (*fi).WT(0).t(0)); glNormal((*fi).V(0)->N()); glVertex((*fi).V(0)->P()); - + glMultiTexCoordARB(GL_TEXTURE0_ARB, (*fi).WT(1).t(0)); glMultiTexCoordARB(GL_TEXTURE1_ARB, (*fi).WT(1).t(0)); glNormal((*fi).V(1)->N()); glVertex((*fi).V(1)->P()); - + glMultiTexCoordARB(GL_TEXTURE0_ARB, (*fi).WT(2).t(0)); glMultiTexCoordARB(GL_TEXTURE1_ARB, (*fi).WT(2).t(0)); glNormal((*fi).V(2)->N()); @@ -799,10 +799,10 @@ void DrawWire() DrawFill(); glPopAttrib(); // } - //else + //else // { // if(!HasEdges()) ComputeEdges(); - + //if(cm==CMPerMesh) glColor(m->C()); //std::vector< MESH_TYPE::VertexType *>::iterator vi; //glBegin(GL_LINE_STRIP); @@ -818,7 +818,7 @@ void DrawWire() // } //} //glEnd(); - // } + // } } void DrawBBox(ColorMode cm) @@ -838,7 +838,7 @@ la mesh non abbia complex (o se li aveva fossero stati detached) Abbia le normali per faccia normalizzate!! -Prende una mesh e duplica tutti gli edge le cui normali nelle facce incidenti formano un angolo maggiore +Prende una mesh e duplica tutti gli edge le cui normali nelle facce incidenti formano un angolo maggiore di (espresso in rad). foreach face foreach unvisited vert vi @@ -859,17 +859,17 @@ void Crease(MESH_TYPE &m, typename MESH_TYPE::scalar_type angleRad) std::vector > SPL; std::vector newvert; newvert.reserve(m.fn*3); - // indica se un il vertice z della faccia e' stato processato - enum {VISITED_0= MESH_TYPE::FaceType::USER0, + // indica se un il vertice z della faccia e' stato processato + enum {VISITED_0= MESH_TYPE::FaceType::USER0, VISITED_1= MESH_TYPE::FaceType::USER0<<1, VISITED_2= MESH_TYPE::FaceType::USER0<<2} ; - int vis[3]={VISITED_0,VISITED_1,VISITED_2}; - + int vis[3]={VISITED_0,VISITED_1,VISITED_2}; + int _t2=clock(); typename MESH_TYPE::FaceIterator fi; for(fi=m.face.begin();fi!=m.face.end();++fi) if(!(*fi).IsD()) (*fi).Supervisor_Flags()&= (~(VISITED_0 | VISITED_1 | VISITED_2)); - + for(fi=m.face.begin();fi!=m.face.end();++fi) if(!(*fi).IsD()) for(int j=0;j<3;++j) @@ -883,7 +883,7 @@ void Crease(MESH_TYPE &m, typename MESH_TYPE::scalar_type angleRad) GLW::VertToSplit spl; spl.newp=false; spl.edge=-1; - + //Primo giro per trovare un bordo da cui partire do { he.FlipF(); @@ -897,15 +897,15 @@ void Crease(MESH_TYPE &m, typename MESH_TYPE::scalar_type angleRad) he.FlipE(); nextf=he.f->F(he.z); typename MESH_TYPE::scalar_type ps=nextf->N()*he.f->N(); - if(psV(he.z)) vz=he.z; if(he.v == he.f->V((he.z+1)%3)) vz=(he.z+1)%3; assert((he.f->Supervisor_Flags() & vis[vz] )==0); } while(he!=she); } - he.FlipE(); - + he.FlipE(); + she=he; newvert.push_back(*(*fi).V(j)); typename MESH_TYPE::vertex_pointer curvert=&newvert.back(); @@ -920,14 +920,14 @@ void Crease(MESH_TYPE &m, typename MESH_TYPE::scalar_type angleRad) if(he.v == he.f->V(he.z)) spl.z=he.z; if(he.v == he.f->V((he.z+1)%3)) spl.z=(he.z+1)%3; assert(spl.z>=0); - //VCTRACE(" -- spinning face vert %i Adding spl face %i vert %i\n", + //VCTRACE(" -- spinning face vert %i Adding spl face %i vert %i\n", // he.v-m.vert.begin(), spl.f-m.face.begin(), spl.z ); assert((spl.f->Supervisor_Flags() & vis[spl.z] )==0); spl.f->Supervisor_Flags() |= vis[spl.z]; SPL.push_back(spl); spl.newp=false; spl.edge=-1; - if(he.IsBorder()) break; + if(he.IsBorder()) break; nextf=he.f->F(he.z); if(nextf==she.f) break; typename MESH_TYPE::scalar_type ps=nextf->N()*he.f->N(); @@ -941,7 +941,7 @@ void Crease(MESH_TYPE &m, typename MESH_TYPE::scalar_type angleRad) he.FlipF(); if(spl.newp) spl.edge=he.z; he.FlipE(); - + }while(he!=she); } assert(SPL.size()==m.fn*3); @@ -952,9 +952,9 @@ void Crease(MESH_TYPE &m, typename MESH_TYPE::scalar_type angleRad) (*vsi).f->V((*vsi).z)=(*vsi).v; if((*vsi).newp){ assert((*vsi).edge>=0 && (*vsi).edge<3); - if(!(*vsi).f->IsBorder( (*vsi).edge) ) + if(!(*vsi).f->IsBorder( (*vsi).edge) ) (*vsi).f->Detach((*vsi).edge); - + } } @@ -978,20 +978,20 @@ void CreaseWN(MESH_TYPE &m, typename MESH_TYPE::scalar_type angle) } typename MESH_TYPE::scalar_type cosangle=Cos(angle); - + typename MESH_TYPE::FaceIterator fi; // Clear the per wedge normals - for(fi=m.face.begin();fi!=m.face.end();++fi) if(!(*fi).IsD()) - { + for(fi=m.face.begin();fi!=m.face.end();++fi) if(!(*fi).IsD()) + { (*fi).WN(0)=MESH_TYPE::vectorial_type(0,0,0); (*fi).WN(1)=MESH_TYPE::vectorial_type(0,0,0); (*fi).WN(2)=MESH_TYPE::vectorial_type(0,0,0); } typename MESH_TYPE::FaceType::vectorial_type nn; - - for(fi=m.face.begin();fi!=m.face.end();++fi) if(!(*fi).IsD()) + + for(fi=m.face.begin();fi!=m.face.end();++fi) if(!(*fi).IsD()) { nn=(*fi).cN(); for(int i=0;i<3;++i) @@ -1003,9 +1003,9 @@ void CreaseWN(MESH_TYPE &m, typename MESH_TYPE::scalar_type angle) } } } - + }*/ -} // end namespace +} // end namespace #endif diff --git a/wrap/gui/coordinateframe.cpp b/wrap/gui/coordinateframe.cpp index fde0caa4..ab0477a9 100644 --- a/wrap/gui/coordinateframe.cpp +++ b/wrap/gui/coordinateframe.cpp @@ -352,7 +352,7 @@ void MovableCoordinateFrame::RotateToAlign(const Point3f source, const Point3f d Point3f axis = dest ^ source; float sinangle = axis.Norm(); - float cosangle = dest * source; + float cosangle = dest.dot(source); float angle = math::Atan2(sinangle,cosangle); if( math::Abs(angle) < EPSILON ) diff --git a/wrap/io_trimesh/import_ptx.h b/wrap/io_trimesh/import_ptx.h index e955ba36..1baeee2e 100644 --- a/wrap/io_trimesh/import_ptx.h +++ b/wrap/io_trimesh/import_ptx.h @@ -417,7 +417,7 @@ namespace io { { raggio = -((*fi).V(0)->P() + (*fi).V(1)->P() + (*fi).V(2)->P()) / 3.0; raggio.Normalize(); - if((raggio * (*fi).N()) < limit) + if((raggio.dot((*fi).N())) < limit) Allocator::DeleteFace(m,*fi); }