1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2010-2011 Hauke Heibel <heibel@gmail.com>
5 //
6 // This Source Code Form is subject to the terms of the Mozilla
7 // Public License v. 2.0. If a copy of the MPL was not distributed
8 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
9
10 #include "main.h"
11
12 #include <unsupported/Eigen/Splines>
13
14 // lets do some explicit instantiations and thus
15 // force the compilation of all spline functions...
16 template class Spline<double, 2, Dynamic>;
17 template class Spline<double, 3, Dynamic>;
18
19 template class Spline<double, 2, 2>;
20 template class Spline<double, 2, 3>;
21 template class Spline<double, 2, 4>;
22 template class Spline<double, 2, 5>;
23
24 template class Spline<float, 2, Dynamic>;
25 template class Spline<float, 3, Dynamic>;
26
27 template class Spline<float, 3, 2>;
28 template class Spline<float, 3, 3>;
29 template class Spline<float, 3, 4>;
30 template class Spline<float, 3, 5>;
31
closed_spline2d()32 Spline<double, 2, Dynamic> closed_spline2d()
33 {
34 RowVectorXd knots(12);
35 knots << 0,
36 0,
37 0,
38 0,
39 0.867193179093898,
40 1.660330955342408,
41 2.605084834823134,
42 3.484154586374428,
43 4.252699478956276,
44 4.252699478956276,
45 4.252699478956276,
46 4.252699478956276;
47
48 MatrixXd ctrls(8,2);
49 ctrls << -0.370967741935484, 0.236842105263158,
50 -0.231401860693277, 0.442245185027632,
51 0.344361228532831, 0.773369994120753,
52 0.828990216203802, 0.106550882647595,
53 0.407270163678382, -1.043452922172848,
54 -0.488467813584053, -0.390098582530090,
55 -0.494657189446427, 0.054804824897884,
56 -0.370967741935484, 0.236842105263158;
57 ctrls.transposeInPlace();
58
59 return Spline<double, 2, Dynamic>(knots, ctrls);
60 }
61
62 /* create a reference spline */
spline3d()63 Spline<double, 3, Dynamic> spline3d()
64 {
65 RowVectorXd knots(11);
66 knots << 0,
67 0,
68 0,
69 0.118997681558377,
70 0.162611735194631,
71 0.498364051982143,
72 0.655098003973841,
73 0.679702676853675,
74 1.000000000000000,
75 1.000000000000000,
76 1.000000000000000;
77
78 MatrixXd ctrls(8,3);
79 ctrls << 0.959743958516081, 0.340385726666133, 0.585267750979777,
80 0.223811939491137, 0.751267059305653, 0.255095115459269,
81 0.505957051665142, 0.699076722656686, 0.890903252535799,
82 0.959291425205444, 0.547215529963803, 0.138624442828679,
83 0.149294005559057, 0.257508254123736, 0.840717255983663,
84 0.254282178971531, 0.814284826068816, 0.243524968724989,
85 0.929263623187228, 0.349983765984809, 0.196595250431208,
86 0.251083857976031, 0.616044676146639, 0.473288848902729;
87 ctrls.transposeInPlace();
88
89 return Spline<double, 3, Dynamic>(knots, ctrls);
90 }
91
92 /* compares evaluations against known results */
eval_spline3d()93 void eval_spline3d()
94 {
95 Spline3d spline = spline3d();
96
97 RowVectorXd u(10);
98 u << 0.351659507062997,
99 0.830828627896291,
100 0.585264091152724,
101 0.549723608291140,
102 0.917193663829810,
103 0.285839018820374,
104 0.757200229110721,
105 0.753729094278495,
106 0.380445846975357,
107 0.567821640725221;
108
109 MatrixXd pts(10,3);
110 pts << 0.707620811535916, 0.510258911240815, 0.417485437023409,
111 0.603422256426978, 0.529498282727551, 0.270351549348981,
112 0.228364197569334, 0.423745615677815, 0.637687289287490,
113 0.275556796335168, 0.350856706427970, 0.684295784598905,
114 0.514519311047655, 0.525077224890754, 0.351628308305896,
115 0.724152914315666, 0.574461155457304, 0.469860285484058,
116 0.529365063753288, 0.613328702656816, 0.237837040141739,
117 0.522469395136878, 0.619099658652895, 0.237139665242069,
118 0.677357023849552, 0.480655768435853, 0.422227610314397,
119 0.247046593173758, 0.380604672404750, 0.670065791405019;
120 pts.transposeInPlace();
121
122 for (int i=0; i<u.size(); ++i)
123 {
124 Vector3d pt = spline(u(i));
125 VERIFY( (pt - pts.col(i)).norm() < 1e-14 );
126 }
127 }
128
129 /* compares evaluations on corner cases */
eval_spline3d_onbrks()130 void eval_spline3d_onbrks()
131 {
132 Spline3d spline = spline3d();
133
134 RowVectorXd u = spline.knots();
135
136 MatrixXd pts(11,3);
137 pts << 0.959743958516081, 0.340385726666133, 0.585267750979777,
138 0.959743958516081, 0.340385726666133, 0.585267750979777,
139 0.959743958516081, 0.340385726666133, 0.585267750979777,
140 0.430282980289940, 0.713074680056118, 0.720373307943349,
141 0.558074875553060, 0.681617921034459, 0.804417124839942,
142 0.407076008291750, 0.349707710518163, 0.617275937419545,
143 0.240037008286602, 0.738739390398014, 0.324554153129411,
144 0.302434111480572, 0.781162443963899, 0.240177089094644,
145 0.251083857976031, 0.616044676146639, 0.473288848902729,
146 0.251083857976031, 0.616044676146639, 0.473288848902729,
147 0.251083857976031, 0.616044676146639, 0.473288848902729;
148 pts.transposeInPlace();
149
150 for (int i=0; i<u.size(); ++i)
151 {
152 Vector3d pt = spline(u(i));
153 VERIFY( (pt - pts.col(i)).norm() < 1e-14 );
154 }
155 }
156
eval_closed_spline2d()157 void eval_closed_spline2d()
158 {
159 Spline2d spline = closed_spline2d();
160
161 RowVectorXd u(12);
162 u << 0,
163 0.332457030395796,
164 0.356467130532952,
165 0.453562180176215,
166 0.648017921874804,
167 0.973770235555003,
168 1.882577647219307,
169 2.289408593930498,
170 3.511951429883045,
171 3.884149321369450,
172 4.236261590369414,
173 4.252699478956276;
174
175 MatrixXd pts(12,2);
176 pts << -0.370967741935484, 0.236842105263158,
177 -0.152576775123250, 0.448975001279334,
178 -0.133417538277668, 0.461615613865667,
179 -0.053199060826740, 0.507630360006299,
180 0.114249591147281, 0.570414135097409,
181 0.377810316891987, 0.560497102875315,
182 0.665052120135908, -0.157557441109611,
183 0.516006487053228, -0.559763292174825,
184 -0.379486035348887, -0.331959640488223,
185 -0.462034726249078, -0.039105670080824,
186 -0.378730600917982, 0.225127015099919,
187 -0.370967741935484, 0.236842105263158;
188 pts.transposeInPlace();
189
190 for (int i=0; i<u.size(); ++i)
191 {
192 Vector2d pt = spline(u(i));
193 VERIFY( (pt - pts.col(i)).norm() < 1e-14 );
194 }
195 }
196
check_global_interpolation2d()197 void check_global_interpolation2d()
198 {
199 typedef Spline2d::PointType PointType;
200 typedef Spline2d::KnotVectorType KnotVectorType;
201 typedef Spline2d::ControlPointVectorType ControlPointVectorType;
202
203 ControlPointVectorType points = ControlPointVectorType::Random(2,100);
204
205 KnotVectorType chord_lengths; // knot parameters
206 Eigen::ChordLengths(points, chord_lengths);
207
208 // interpolation without knot parameters
209 {
210 const Spline2d spline = SplineFitting<Spline2d>::Interpolate(points,3);
211
212 for (Eigen::DenseIndex i=0; i<points.cols(); ++i)
213 {
214 PointType pt = spline( chord_lengths(i) );
215 PointType ref = points.col(i);
216 VERIFY( (pt - ref).matrix().norm() < 1e-14 );
217 }
218 }
219
220 // interpolation with given knot parameters
221 {
222 const Spline2d spline = SplineFitting<Spline2d>::Interpolate(points,3,chord_lengths);
223
224 for (Eigen::DenseIndex i=0; i<points.cols(); ++i)
225 {
226 PointType pt = spline( chord_lengths(i) );
227 PointType ref = points.col(i);
228 VERIFY( (pt - ref).matrix().norm() < 1e-14 );
229 }
230 }
231 }
232
233
test_splines()234 void test_splines()
235 {
236 CALL_SUBTEST( eval_spline3d() );
237 CALL_SUBTEST( eval_spline3d_onbrks() );
238 CALL_SUBTEST( eval_closed_spline2d() );
239 CALL_SUBTEST( check_global_interpolation2d() );
240 }
241