I have a polynomial function which I need to approximate by a sum of sine waves with constant amplitude along a given domain.
From what I hear, this might be a good time to make use of Fourier Transforms. However, my function is a simple low-degree polynomial (typically just degree 2) and I need this to be computed very quickly on an everyday computer. I'm not familiar with Fourier Transforms, so I don't know if they're good for this or not.
So here's my approach as of yet:
Let's take the function $f(x) = x – \frac16 x^3 + \frac1{180} x^5$. looking around $-\pi$ to $+\pi$, it appears to have an angular frequency of approximately 1 (This is just a slightly modified expansion of the taylor series for sin(x))
Taking advantage of the fact that sin(x) is approximately equal to x around 0, we can estimate f(x) with $\sin(x) – \frac16 \sin(x)^3 + \frac1{180} \sin(x)^5$.
But this is not the best estimation.
The closer we get to 0, the more sin(x) looks like x.
So if we instead use
$L \sin(x/L) – \frac16 (L \sin(x/L))^3 + \frac1{180} (L \sin(x/L))^5$
and make L a large constant, say 1000, the approximation looks very accurate on the given domain.
This simplifies to:
$L \sin(x/L) – \frac14\frac{L^3}6(3 \sin(x/L) – \sin(3x/L)) + \frac1{16}\frac{L^5}{180} (10 \sin(x/L) – 5 \sin(3x/L) + \sin(5x/L))$,
or $(L – \frac{L^3}8 + \frac{L^5}{288}) \sin(x/L) + (\frac{L^3}{24} – \frac{L^5}{576}) \sin(3x/L)+\frac{L^5}{2880} \sin(5x/L)$
Now, the problem: This is an extremely good approximation on the domain, but the frequency of these sign waves carries very little meaning because they're mostly arbitrary based on my value of L.
There is no unique mapping from the polynomial to the frequencies encoded in it. For instance, I could claim my polynomial is composed of waves of frequencies [0.1, 0.3, and 0.5], or [0.001, 0.003, and 0.005], based on what I set L to.
And with sufficiently large L, all the frequencies are effectively 0.
However, most signal processing algorithms all approximately agree with eachother on what frequencies reside in a signal, unlike mine. And it's clear that f(x) is approximately sin(x) on the domain. So it seems like there should be a way to determine that the frequencies in this signal are all around 1.
What are your recommendations for estimating a polynomial function with a sum of sine waves?
As a final note, I did try to simplify my application for the sake of making a concise question. However, if more information helps, my actual function is $g(x) = P(x) \sin(Bx+c)$ where P(x) is a polynomial function and B and c are constants. My goal is to re-express g as a sum of sine waves. If I can express P as a sum of sine waves (which is is what the body of my post here is about), then I can multiply each sine wave by sin(Bx+c) and make use of the fact that $\sin(a) \sin(b) = \frac12 (\cos(a-b)-\cos(a+b))$ to write g as a sum of sine waves. Also, the resultant sine waves do not need to have constant frequencies.
Best Answer
Representing a signal as a sum of sine waves is exactly what the Fourier transform does and it can be computed extremely fast with the Fast Fourier Transform (FFT) algorithm.
Evaluate your function on a set of $N$ evenly spaced points. For instance, let $x_n = -\pi + n \Delta x$ for $n=0,1,\dots,N-1$ and $\Delta x = 2\pi / (N-1)$ similar to your example. Then take the the FFT of $f(\mathbf{x})$. The result of will be the weights of a set of complex exponentials such that their weigthed sum is equal to $f(\mathbf{x})$. That is, if $\mathbf{w} = \operatorname{FFT}(f(\mathbf{x}))$, then
$$ f(x_n) = \sum_{m=0}^{N-1} w_m e^{j 2\pi n m/N} $$
If you want to stay away from complex exponentials, you can use the Discrete Cosine Transform instead, which as the name suggests will only use cosines in the weighted sum and will look something like
$$ f(x_n) = \sum_{m=0}^{N-1} w_m \cos \left( {\pi\over N}\left(m+{1\over 2}\right) n \right). $$