|
|
We designed an API to automatically fit a polynomial from a finite set of samples. After specifying a basis to the finite-dimensional linear space (over the field of real numbers) of n-dimensional polynomials with degree less than, or equal to, a non-negative integer d, we are able to both: fit (point,value)-pairs and fit the zero-set of a polynomial from a given set of points. In both applications, weights may be assigned (both to samples and to the elements in the basis) to bias the fitting process.
 |
 |
| Degree 25 approximation of the sinc function | Degree 3 fitting of 8 samples from S1 |
The polyfit library was designed to be of simple use and easy to extend. For simplicity, the user just needs to provide enough data to perform the desired fitting (e.g. basis, samples etc.) by instancing a few classes and calling other few methods. The design is based on the main concepts involved in the fitting process: the specification of an arbitrary basis to a finite-dimensional function space and the evaluation of the fitted polynomial. This is the rationale which motivated our architecture (below follows a simplified class diagram containing the public interface avaliable to the user):
To use the polyfit library, the user just needs to include the main header polyfit.h in his code, instantiate an object of type PolynomialBasis<Real>, one from the class Fitter<Real> and provide the data required to perform the fitting. In these templates, Real is assumed to be float or double, providing flexibility to work with both of these data types (if memory is not a huge problem, we recommend to use double whenever possible).
For a more detailed description of the parameters in the above classes and their methods and how to provide them, we refer the reader to polyfit's API documentation, to the demo program in the source code and to the examples in the following section.
Observations:
0. The polyfit library requires both netlib's version of clapack [1] and the ATLAS' implementation of the cblas routines [3,2].
1. When compiling a project which uses the polyfit library (with GNU's g++), don't forget to add the libraries at the end of the command-line in the same order as done with the demo program (change the library names according to those in your system):
g++ -Wall -O2 -o polyfit polyfit.cpp -I../include -lclapack_netlib -lcblaswr -lcblas -latlas -lF77 -lI77 -lm2. For numerical stability reasons, we do not recommend the use of MonomialBasis in calculations involving large dimensions and large polynomial degrees, i.e. if you really need to work with high-order polynomials in large dimensions, you should extend the library with a better behaved polynomial basis.
3. The fitting of zero-sets by the implemented method is an unstable problem and, therefore, its use in noisy point clouds (with moderate-order polynomials) is prone to errors. Although some works try to overcome these problems using more information about the desired fit (e.g. tangents and normals), we only implemented the basic (unstable) method. For more information on these methods, see [6] and references cited therein.
ToDo List:
1. Implement exception handling to improve robustness and debugging information.
2. Optimize the polynomial evaluation code both for efficiency and numerical stability reasons [7].
3. Extend the library with more (and numerically better conditioned) bases.
4. Inspect the least squares solver in clapack's {s|d}gels_ routines and reimplement a solution adapted to our applications.
5. Profile a real application and optimize the library's hot spots.
Below, we provide simple examples of the two main fitting procedures available through the polyfit library.
 |
 |
|
Fitting a parabola: parabola.cpp
|
Fitting an ellipse: ellipse.cpp
|
Next we show some results of applying the polynomial fitting process to approximate a halftone image.
Observation: Notice that this was done for illustration purposes, the discontinuities in (and the complexity of) the original image make the approximation process very hard to accomplish with the use of a single bivariate polynomial.
 |
 |
|
Original
|
Degree 16
|
 |
 |
|
Degree 32
|
Degree 48
|
Source Code: [gzipped tar archive - 25K] (use it by your own risk... :D)
Documentation (generated by Doxygen): [HTML] [HTML tarball - 327K] [PDF file - 726K]
References:
[1] LAPACK -- Linear Algebra PACKage
[2] BLAS (Basic Linear Algebra Subprograms)
[3] Automatically Tuned Linear Algebra Software (ATLAS)
[4] Trefethen, Lloyd N.; Bau, David, III, Numerical Linear Algebra, SIAM, Philadelphia, PA, 1997
[5] Demmel, James W., Applied Numerical Linear Algebra, SIAM, Philadelphia, PA, 1997
[6] Tasdizen, Tolga; Tarel, Jean-Philippe; Cooper, David B., Improving the Stability of Algebraic Curves for Applications, IEEE Transactions on Image Processing, Volume 9, Number 3, March 2000
[7] Peña, J. M., On the Multivariate Horner Scheme, SIAM Journal on Numerical Analysis, Volume 37, Issue 4, Pages 1186-1197, 2000
[8] Nijenhuis, Albert; Wilf, Herbert S., Combinatorial Algorithms - For Computers and Calculators, Second Edition, Academic Press, 1978
|
|