Simbody  3.5
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
SplineFitter.h
Go to the documentation of this file.
1 #ifndef SimTK_SIMMATH_SPLINE_FITTER_H_
2 #define SimTK_SIMMATH_SPLINE_FITTER_H_
3 
4 /* -------------------------------------------------------------------------- *
5  * Simbody(tm): SimTKmath *
6  * -------------------------------------------------------------------------- *
7  * This is part of the SimTK biosimulation toolkit originating from *
8  * Simbios, the NIH National Center for Physics-Based Simulation of *
9  * Biological Structures at Stanford, funded under the NIH Roadmap for *
10  * Medical Research, grant U54 GM072970. See https://simtk.org/home/simbody. *
11  * *
12  * Portions copyright (c) 2008-12 Stanford University and the Authors. *
13  * Authors: Peter Eastman *
14  * Contributors: *
15  * *
16  * Licensed under the Apache License, Version 2.0 (the "License"); you may *
17  * not use this file except in compliance with the License. You may obtain a *
18  * copy of the License at http://www.apache.org/licenses/LICENSE-2.0. *
19  * *
20  * Unless required by applicable law or agreed to in writing, software *
21  * distributed under the License is distributed on an "AS IS" BASIS, *
22  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. *
23  * See the License for the specific language governing permissions and *
24  * limitations under the License. *
25  * -------------------------------------------------------------------------- */
26 
27 #include "SimTKcommon.h"
31 
32 namespace SimTK {
33 
68 template <class T>
69 class SplineFitter {
70 public:
71  SplineFitter(const SplineFitter& copy) : impl(copy.impl) {
72  impl->referenceCount++;
73  }
75  impl = copy.impl;
76  impl->referenceCount++;
77  return *this;
78  }
80  impl->referenceCount--;
81  if (impl->referenceCount == 0)
82  delete impl;
83  }
84 
93  (int degree, const Vector& x, const Vector_<T>& y) {
94  Vector_<T> coeff;
95  Vector wk;
96  int ier;
97  GCVSPLUtil::gcvspl(x, y, Vector(x.size(), 1.0), T(1),
98  degree, 2, 0, coeff, wk, ier);
100  (degree, Spline_<T>(degree, x, coeff), wk[3], wk[4], wk[2]));
101  }
102 
112  (int degree, const Vector& x, const Vector_<T>& y, Real error) {
113  Vector_<T> coeff;
114  Vector wk;
115  int ier;
116  GCVSPLUtil::gcvspl(x, y, Vector(x.size(), 1.0), T(1),
117  degree, 3, error, coeff, wk, ier);
119  (degree, Spline_<T>(degree, x, coeff), wk[3], wk[4], wk[2]));
120  }
121 
130  static SplineFitter fitFromDOF
131  (int degree, const Vector& x, const Vector_<T>& y, Real dof) {
132  Vector_<T> coeff;
133  Vector wk;
134  int ier;
135  GCVSPLUtil::gcvspl(x, y, Vector(x.size(), 1.0), T(1),
136  degree, 4, dof, coeff, wk, ier);
138  (degree, Spline_<T>(degree, x, coeff), wk[3], wk[4], wk[2]));
139  }
140 
150  (int degree, const Vector& x, const Vector_<T>& y, Real p) {
151  Vector_<T> coeff;
152  Vector wk;
153  int ier;
154  GCVSPLUtil::gcvspl(x, y, Vector(x.size(), 1.0), T(1),
155  degree, 1, p, coeff, wk, ier);
157  (degree, Spline_<T>(degree, x, coeff), wk[3], wk[4], wk[2]));
158  }
163  return impl->spline;
164  }
169  return impl->p;
170  }
175  return impl->error;
176  }
181  return impl->dof;
182  }
183 private:
184  class SplineFitterImpl;
185  SplineFitter(SplineFitterImpl *impl) : impl(impl) {
186  }
187  SplineFitterImpl* impl;
188 };
189 
190 template <class T>
192 public:
193  SplineFitterImpl(int degree, const Spline_<T>& spline, Real p, Real error, Real dof) : referenceCount(1), degree(degree), spline(spline), p(p), error(error), dof(dof) {
194  }
196  assert(referenceCount == 0);
197  }
199  int degree;
201  Real p, error, dof;
202 };
203 
204 } // namespace SimTK
205 
206 #endif // SimTK_SIMMATH_SPLINE_FITTER_H_
207 
208 
Real error
Definition: SplineFitter.h:201
Real dof
Definition: SplineFitter.h:201
static SplineFitter fitFromErrorVariance(int degree, const Vector &x, const Vector_< T > &y, Real error)
Perform a fit, choosing a value of the smoothing parameter based on the known error variance in the d...
Definition: SplineFitter.h:112
static SplineFitter fitFromDOF(int degree, const Vector &x, const Vector_< T > &y, Real dof)
Perform a fit, choosing a value of the smoothing parameter based on the expect number of degrees of f...
Definition: SplineFitter.h:131
Real p
Definition: SplineFitter.h:201
This is the top-level SimTK namespace into which all SimTK names are placed to avoid collision with o...
Definition: Assembler.h:37
int size() const
Definition: VectorBase.h:396
SplineFitterImpl(int degree, const Spline_< T > &spline, Real p, Real error, Real dof)
Definition: SplineFitter.h:193
~SplineFitter()
Definition: SplineFitter.h:79
static SplineFitter fitFromGCV(int degree, const Vector &x, const Vector_< T > &y)
Perform a fit, choosing a value of the smoothing parameter that minimizes the Generalized Cross Valid...
Definition: SplineFitter.h:93
SimTK_Real Real
This is the default compiled-in floating point type for SimTK, either float or double.
Definition: SimTKcommon/include/SimTKcommon/internal/common.h:577
SplineFitter operator=(const SplineFitter &copy)
Definition: SplineFitter.h:74
Real getDegreesOfFreedom()
Get the estimate of the number of degrees of freedom of the residual that was determined by the fitti...
Definition: SplineFitter.h:180
SplineFitter(const SplineFitter &copy)
Definition: SplineFitter.h:71
Includes internal headers providing declarations for the basic SimTK Core classes, including Simmatrix.
Given a set of data points, this class creates a Spline_ which interpolates or approximates them...
Definition: SplineFitter.h:69
~SplineFitterImpl()
Definition: SplineFitter.h:195
Definition: SplineFitter.h:191
const Spline_< T > & getSpline()
Get the Spline_ that was generated by the fitting.
Definition: SplineFitter.h:162
Spline_< T > spline
Definition: SplineFitter.h:200
int degree
Definition: SplineFitter.h:199
This is the header file that every Simmath compilation unit should include first. ...
int referenceCount
Definition: SplineFitter.h:198
This class implements a non-uniform Bezier curve.
Definition: Spline.h:52
static SplineFitter fitForSmoothingParameter(int degree, const Vector &x, const Vector_< T > &y, Real p)
Perform a fit, using a specified fixed value for the smoothing parameter.
Definition: SplineFitter.h:150
Vector_< Real > Vector
Variable-size column vector of Real elements; abbreviation for Vector_.
Definition: BigMatrix.h:1473
void copy(Row< 1, E1, S1 > &r1, const Row< 1, E2, S2 > &r2)
Definition: Row.h:105
static void gcvspl(const Vector &x, const Vector &y, const Vector &wx, Real wy, int m, int md, Real val, Vector &c, Vector &wk, int &ier)
Real getSmoothingParameter()
Get the smoothing parameter that was used for the fitting.
Definition: SplineFitter.h:168
Real getMeanSquaredError()
Get the estimate of the true mean squared error in the data that was determined by the fitting...
Definition: SplineFitter.h:174