00001
00025
00026
00027 #ifndef __ETL__BSPLINE_H
00028 #define __ETL__BSPLINE_H
00029
00030
00031
00032 #include <vector>
00033 #include <functional>
00034 #include "_curve_func.h"
00035
00036
00037
00038
00039
00040
00041
00042 _ETL_BEGIN_NAMESPACE
00043
00044 template <class T, class K=float, class C=affine_combo<T,K>, class D=distance_func<T> >
00045 class bspline : public std::unary_function<K,T>
00046 {
00047 public:
00048 typedef T value_type;
00049 typedef K knot_type;
00050 typedef std::vector<knot_type> knot_container;
00051 typedef std::vector<value_type> cpoint_container;
00052 typedef typename knot_container::iterator knot_iterator;
00053 typedef typename cpoint_container::iterator cpoint_iterator;
00054
00055 typedef C affine_func_type;
00056 typedef D distance_func_type;
00057
00058 protected:
00059 affine_func_type affine_func;
00060 distance_func_type distance_func;
00061
00062 private:
00063 int m;
00064 knot_container _knots;
00065 cpoint_container _cpoints;
00066 bool _loop;
00067
00068 public:
00069 bspline():m(2),_loop(false) { }
00070
00071 int get_m()const { return m-1; };
00072 int set_m(int new_m) { m=new_m+1; return m-1; };
00073
00074 bool set_loop(bool x) { _loop=x; reset_knots(); return _loop; }
00075
00076 knot_container & knots() { return _knots; };
00077 cpoint_container & cpoints() { return _cpoints; };
00078
00079 const knot_container & knots()const { return _knots; };
00080 const cpoint_container & cpoints()const { return _cpoints; };
00081
00082 void reset_knots()
00083 {
00084 int i;
00085
00086 if(!_loop)
00087 {
00088 _knots.clear();
00089 if(!_cpoints.size())
00090 return;
00091 while(m>(signed)_cpoints.size())
00092 m--;
00093 for(i=0;i<m;i++)
00094 *_knots.insert(_knots.end())=0;
00095 for(i=1;i<(signed)_cpoints.size()-m+1;i++)
00096 *_knots.insert(_knots.end())=i;
00097 for(i=0;i<m;i++)
00098 *_knots.insert(_knots.end())=_cpoints.size()-m+1;
00099 }
00100 else
00101 {
00102 _knots.clear();
00103 if(!_cpoints.size())
00104 return;
00105 while(m>(signed)_cpoints.size())
00106 m--;
00107 for(i=0;i<=(signed)_cpoints.size()-m+1;i++)
00108 *_knots.insert(_knots.end())=i;
00109 }
00110 }
00111
00112 int calc_curve_segment(knot_type t)const
00113 {
00114 int k;
00115 if(t<0)
00116 t=0;
00117 if(t>=_knots.back())
00118 t=_knots.back()-0.0001;
00119 for(k=0;_knots[k]>t || _knots[k+1]<=t;k++);
00120
00121 return k;
00122 }
00123
00124 knot_container get_segment_knots(int i)const
00125 {
00126 if(i+1<m)
00127 {
00128 knot_container ret(_knots.begin(),_knots.begin()+i+m+1);
00129
00130 return ret;
00131 }
00132 if(i+1>=(signed)_knots.size())
00133 {
00134 knot_container ret(_knots.begin()+i-m+1,_knots.end());
00135 return ret;
00136 }
00137 return knot_container(_knots.begin()+i-m+1, _knots.begin()+i+m);
00138 }
00139
00140 cpoint_container get_segment_cpoints(int i)const
00141 {
00142 if(i+1<m)
00143 {
00144 return cpoint_container();
00145 }
00146 if(i+1>=(signed)_knots.size())
00147 {
00148 return cpoint_container();
00149 }
00150 return cpoint_container(_cpoints.begin()+i-m+1, _cpoints.begin()+i+1);
00151 }
00152
00153 cpoint_container calc_shell(knot_type t, int level)const
00154 {
00155 int
00156 i=calc_curve_segment(t),
00157 j,k;
00158
00159 knot_container u=get_segment_knots(i);
00160
00161 cpoint_container d=get_segment_cpoints(i);
00162
00163 if(!d.size())
00164 return cpoint_container();
00165
00166 for(j=0;d.size()>1 && j<level;d.pop_back(),j++)
00167 {
00168 for(k=0;k<d.size()-1;k++)
00169 {
00170 d[k]=affine_func(d[k],d[k+1],((t-u[j+k+1])/(u[m+k]-u[j+k+1])));
00171 }
00172 }
00173 return d;
00174 }
00175
00176 value_type operator()(knot_type t)const
00177 {
00178 return get_curve_val(calc_curve_segment(t),t);
00179 }
00180
00181 value_type get_curve_val(int i,knot_type t)const
00182 {
00183 int
00184 j,k;
00185
00186 knot_container u=get_segment_knots(i);
00187
00188 cpoint_container d=get_segment_cpoints(i);
00189
00190 if(!d.size())
00191 return value_type();
00192
00193 for(j=0;d.size()>1;d.pop_back(),j++)
00194 {
00195 for(k=0;k<(signed)d.size()-1;k++)
00196 {
00197 d[k]=affine_func(d[k],d[k+1],((t-u[j+k+1])/(u[m+k]-u[j+k+1])));
00198 }
00199 }
00200 return d.front();
00201 }
00202
00203 cpoint_iterator find_closest_cpoint(const value_type &point, typename distance_func_type::result_type max)
00204 {
00205 cpoint_iterator i=_cpoints.begin();
00206 cpoint_iterator ret=i;
00207 typename distance_func_type::result_type dist=distance_func(point,_cpoints[0]);
00208
00209
00210
00211
00212 max=distance_func.cook(max);
00213
00214 for(++i;i<_cpoints.end();i++)
00215 {
00216 typename distance_func_type::result_type thisdist=distance_func(point,*i);
00217
00218 if(thisdist<dist)
00219 {
00220 dist=thisdist;
00221 ret=i;
00222 }
00223 }
00224 if(dist<max)
00225 return ret;
00226 return _cpoints.end();
00227 }
00228 };
00229
00230 _ETL_END_NAMESPACE
00231
00232
00233
00234
00235
00236 #endif