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