Re: Efficiently compute Fourier coefficients for discrete function
- To: mathgroup at smc.vnet.net
- Subject: [mg105450] Re: Efficiently compute Fourier coefficients for discrete function
- From: Bill Rowe <readnews at sbcglobal.net>
- Date: Thu, 3 Dec 2009 06:14:20 -0500 (EST)
On 12/2/09 at 6:25 AM, kelly.terry.jones at gmail.com (Kelly Jones)
wrote:
>Mathematica-wise, how to efficiently compute the Fourier
>coefficients for a function defined at finite "random" points.
What is it you see as being efficient. Possibly, you are trying
to maximize the amount of computation that can be done in a
given time on a given computer system. Another possibility it
you want to minimize the total time required to obtain a
solution. Note this last would also include time you spend
writing code in addition to execution time for that code.
Possibly you have something else in mind.
I will assume you want to minimize the total time to get a solution.
>Example:
>This defines Sin[3*x] for 100 random values of x:
>l = Table[Random[],{i,1,100}]
>l2 = Table[{l[[i]],Sin[3*l[[i]]]},{i,1,Length[l]}]
Since Sin operates on lists l2 can be computed more efficiently by
l = RandomReal[1,{100}];
l2 = Transpose@{l,Sin[3 l]};
>This calculates the nth Fourier coefficient:
>pfourier[s_,n_] := Sum[s[[i]][[2]]*Exp[-I*n*s[[i]][[1]]],
>{i,1,Length[s]}]
>(I realize I'm off by a factor of pi; I'm only interested in
>relative size for now).
>As expected, the imaginary/sine part of the 3rd coefficient has the
>greatest magnitude:
>Table[Im[pfourier[l2,n]],{n,1,10}]
>I could compute pfourier[l2,n] for all n, but that seems
>inefficient.
>Is there a simpler function (similar to Fourier and
>FourierTransform) that does this?
A few observations. Since you have chosen values from a uniform
distribution, this will fairly closely approximate a uniformly
sampled function as the number of samples gets higher. That is
using Fourier directly should be effective at finding the
maximum frequency component of the data.
For the rest of my remarks, I will use a sample of 128 points
since I know Fourier is a bit more efficient when the number of
points is a power of 2.
Notice the two plots
x = RandomReal[1,{128}];
ListPlot[Transpose@{x, Sin[3 x]}]
and
ListPlot[Sin[3 Sort[x]],DataRange->{0,1}]
are fairly similar.
So the simple straight forward way of finding the index of the
maximum frequency component is:
In[19]:= Ordering[Abs[Fourier[Sin[3 Sort@x]][[2 ;; 64]]], -1]
Out[19]= {1}
Here, I take advantage of knowing the first coefficient is the
dc component and all of the useful frequency information is
contained in the first half of the values returned by Fourier.
The result isn't 3 as you might expect because the data does not
represent sin[3 x]. Instead, it is a half-sine wave since you
restricted x to the range of 0 to 1. If you had sampled over a
full cycle, the result would be 3 as should be expected. That is:
In[20]:= x = RandomReal[2 \[Pi], {128}];
Ordering[Abs[Fourier[Sin[3 Sort@x]][[2 ;; 64]]], -1]
Out[21]= {3}
Unless you have some purpose in mind you haven't made clear, I
see no point to creating your own function to compute the nth coefficient.
Do note, the cheat I am doing here (treating the random sample
intervals as being uniform intervals) will work increasingly
well as the number of sample points increases. But clearly as
the number of sample points decrease, this cheat will become a problem.