1 // This file is part of libigl, a simple c++ geometry processing library.
2 //
3 // Copyright (C) 2014 Daniele Panozzo <daniele.panozzo@gmail.com>, Olga Diamanti <olga.diam@gmail.com>
4 //
5 // This Source Code Form is subject to the terms of the Mozilla Public License
6 // v. 2.0. If a copy of the MPL was not distributed with this file, You can
7 // obtain one at http://mozilla.org/MPL/2.0/.
8
9 #include "comb_cross_field.h"
10
11 #include <vector>
12 #include <deque>
13 #include <Eigen/Geometry>
14 #include "per_face_normals.h"
15 #include "is_border_vertex.h"
16 #include "rotation_matrix_from_directions.h"
17
18 #include "triangle_triangle_adjacency.h"
19
20 namespace igl {
21 template <typename DerivedV, typename DerivedF>
22 class Comb
23 {
24 public:
25
26 const Eigen::PlainObjectBase<DerivedV> &V;
27 const Eigen::PlainObjectBase<DerivedF> &F;
28 const Eigen::PlainObjectBase<DerivedV> &PD1;
29 const Eigen::PlainObjectBase<DerivedV> &PD2;
30 DerivedV N;
31
32 private:
33 // internal
34 DerivedF TT;
35 DerivedF TTi;
36
37
38 private:
39
40
Sign(double a)41 static inline double Sign(double a){return (double)((a>0)?+1:-1);}
42
43
44 private:
45
46 // returns the 90 deg rotation of a (around n) most similar to target b
47 /// a and b should be in the same plane orthogonal to N
K_PI_new(const Eigen::Matrix<typename DerivedV::Scalar,3,1> & a,const Eigen::Matrix<typename DerivedV::Scalar,3,1> & b,const Eigen::Matrix<typename DerivedV::Scalar,3,1> & n)48 static inline Eigen::Matrix<typename DerivedV::Scalar, 3, 1> K_PI_new(const Eigen::Matrix<typename DerivedV::Scalar, 3, 1>& a,
49 const Eigen::Matrix<typename DerivedV::Scalar, 3, 1>& b,
50 const Eigen::Matrix<typename DerivedV::Scalar, 3, 1>& n)
51 {
52 Eigen::Matrix<typename DerivedV::Scalar, 3, 1> c = (a.cross(n)).normalized();
53 typename DerivedV::Scalar scorea = a.dot(b);
54 typename DerivedV::Scalar scorec = c.dot(b);
55 if (fabs(scorea)>=fabs(scorec))
56 return a*Sign(scorea);
57 else
58 return c*Sign(scorec);
59 }
60
61
62
63 public:
Comb(const Eigen::PlainObjectBase<DerivedV> & _V,const Eigen::PlainObjectBase<DerivedF> & _F,const Eigen::PlainObjectBase<DerivedV> & _PD1,const Eigen::PlainObjectBase<DerivedV> & _PD2)64 inline Comb(const Eigen::PlainObjectBase<DerivedV> &_V,
65 const Eigen::PlainObjectBase<DerivedF> &_F,
66 const Eigen::PlainObjectBase<DerivedV> &_PD1,
67 const Eigen::PlainObjectBase<DerivedV> &_PD2
68 ):
69 V(_V),
70 F(_F),
71 PD1(_PD1),
72 PD2(_PD2)
73 {
74 igl::per_face_normals(V,F,N);
75 igl::triangle_triangle_adjacency(F,TT,TTi);
76 }
comb(Eigen::PlainObjectBase<DerivedV> & PD1out,Eigen::PlainObjectBase<DerivedV> & PD2out)77 inline void comb(Eigen::PlainObjectBase<DerivedV> &PD1out,
78 Eigen::PlainObjectBase<DerivedV> &PD2out)
79 {
80 // PD1out = PD1;
81 // PD2out = PD2;
82 PD1out.setZero(F.rows(),3);PD1out<<PD1;
83 PD2out.setZero(F.rows(),3);PD2out<<PD2;
84
85 Eigen::VectorXi mark = Eigen::VectorXi::Constant(F.rows(),false);
86
87 std::deque<int> d;
88
89 d.push_back(0);
90 mark(0) = true;
91
92 while (!d.empty())
93 {
94 int f0 = d.at(0);
95 d.pop_front();
96 for (int k=0; k<3; k++)
97 {
98 int f1 = TT(f0,k);
99 if (f1==-1) continue;
100 if (mark(f1)) continue;
101
102 Eigen::Matrix<typename DerivedV::Scalar, 3, 1> dir0 = PD1out.row(f0);
103 Eigen::Matrix<typename DerivedV::Scalar, 3, 1> dir1 = PD1out.row(f1);
104 Eigen::Matrix<typename DerivedV::Scalar, 3, 1> n0 = N.row(f0);
105 Eigen::Matrix<typename DerivedV::Scalar, 3, 1> n1 = N.row(f1);
106
107
108 Eigen::Matrix<typename DerivedV::Scalar, 3, 1> dir0Rot = igl::rotation_matrix_from_directions(n0, n1)*dir0;
109 dir0Rot.normalize();
110 Eigen::Matrix<typename DerivedV::Scalar, 3, 1> targD = K_PI_new(dir1,dir0Rot,n1);
111
112 PD1out.row(f1) = targD;
113 PD2out.row(f1) = n1.cross(targD).normalized();
114
115 mark(f1) = true;
116 d.push_back(f1);
117
118 }
119 }
120
121 // everything should be marked
122 for (int i=0; i<F.rows(); i++)
123 {
124 assert(mark(i));
125 }
126 }
127
128
129
130 };
131 }
132 template <typename DerivedV, typename DerivedF>
comb_cross_field(const Eigen::PlainObjectBase<DerivedV> & V,const Eigen::PlainObjectBase<DerivedF> & F,const Eigen::PlainObjectBase<DerivedV> & PD1,const Eigen::PlainObjectBase<DerivedV> & PD2,Eigen::PlainObjectBase<DerivedV> & PD1out,Eigen::PlainObjectBase<DerivedV> & PD2out)133 IGL_INLINE void igl::comb_cross_field(const Eigen::PlainObjectBase<DerivedV> &V,
134 const Eigen::PlainObjectBase<DerivedF> &F,
135 const Eigen::PlainObjectBase<DerivedV> &PD1,
136 const Eigen::PlainObjectBase<DerivedV> &PD2,
137 Eigen::PlainObjectBase<DerivedV> &PD1out,
138 Eigen::PlainObjectBase<DerivedV> &PD2out)
139 {
140 igl::Comb<DerivedV, DerivedF> cmb(V, F, PD1, PD2);
141 cmb.comb(PD1out, PD2out);
142 }
143
144 #ifdef IGL_STATIC_LIBRARY
145 // Explicit template instantiation
146 template void igl::comb_cross_field<Eigen::Matrix<double, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, 3, 0, -1, 3> >(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> >&);
147 template void igl::comb_cross_field<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&);
148 #endif
149