Fourier Series

Share icon

Hello,

I am trying to implement the Fourier Series (attachment1: "fourier_series") on Igor as a programming exercise. My biggest issue is integrating. While calculating the integral, I'm not sure how the variable of integration should be declared. Because "integrate1d" uses another function (usr function) to calculate the definite integral, should I declare it in the usrf?. The code I have written is

<br />
<br />
#pragma rtGlobals=1     // Use modern global access method.<br />
<br />
function integrand_ak (func)<br />
<br />
wave func<br />
<br />
nvar i, t<br />
<br />
return func*cos(i*pi/t)<br />
<br />
end<br />
<br />
function integrand_bk (func)<br />
<br />
wave func<br />
<br />
nvar i, t<br />
<br />
return func*sin(i*pi/t)<br />
<br />
end<br />
<br />
function fseries(func, L, harmonic)<br />
<br />
Wave func<br />
Variable L<br />
Variable harmonic<br />
variable a0<br />
make/n=(harmonic) ak<br />
make/n=(harmonic) bk<br />
make/n=(harmonic) f<br />
variable/G i, t=L<br />
<br />
a0=(1/L)*integrate1d(func, -L, L)<br />
for(i=1; i<(harmonic+1); i+=1)<br />
    ak[i]=(1/L)*integrate1d(integrand_ak, -L, L)<br />
    bk[i]=(1/L)*integrate1d(integrand_bk, -L, L)<br />
    f[i]=a0/2+ak[i]*cos(i*pi/L)+bk[i]*sin(i*pi/L)<br />
endfor<br />
<br />
end<br />
<br />
When I run it I get the following message: "While executing Integrate_1D, the following error ococured: Bad user function format"<br />
<br />
All help is appreciated.<br />
<br />
Thanks in advance,<br />
<br />
Gabriel </span>
 fourier_series.jpg (53.58 KB)
gzvitiello wrote:

I am trying to implement the Fourier Series (attachment1: "fourier_series") on Igor as a programming exercise.


For various reasons it may not be a good idea to try and compute the integrals directly but that is not the central issue of your question.

Your user functions for Integrate1D should return the value of the integrands for an input parameter (say inx). Since the functions take only a single input parameter you need to rely on global variables to access 'm', 'L' and the data in a wave 'func'. I would use a couple of functions. The first is used to interpolate the data in the wave:
Function interpF(inWave,inx)
    Wave inWave
    Variable inx
   
    NVAR startValue,endValue    // faster than computing this.
    if(inx<startValue || inx>endValue)
        return 0
    else
        return inWave(inx)
    endif
End

Next come a function for expressing one of the integrands:
Function integrand1(inx)
Variable inx

    NVAR mmLL   // globals that define m*pi/L
        Wave func         // fill this wave reference as needed.
    return interpF(func,inx)*sin(mmLL*inx)
End


Somewhere in your calling routine you need to define the globals, e.g.,
Variable/G startValue=...
Variable/G endValue=...
Variable/G mmLL

and then in your loop set
mmLL=m*pi/L


I hope this helps you understand the use of Integrate1D. Before you spend too much time trying to evaluate these integrals (especially for large mmLL), you may want to take a break and reconsider if integration is the the appropriate approach.

A.G.
WaveMetrics, Inc.

November 3, 2014 at 12:55 pm - Permalink

I suspect that AG was politely trying to steer the user to a more efficient FFT analysis. However, I want to mention one case where the Integrate1D( ) method may be appropriate. In many microwave or optical waveguide problems, the even-symmetry basis modes are cosines, but with zero-value at the guide boundaries. These modes do not have complete 2*pi repetition periods, and are not the same as the cosine functions of Fourier analysis. To find the decomposition of an arbitrary even-symmetry input function into these modes, Integrate1D( ) seems like the simplest approach.

A possible FFT solution might be to artificially double the window size, and use only the FFT sine results, with an input function restricted to one (half) side of the waveguide. The added complications of such a procedure might not be worth the effort if speed is not an issue.

November 5, 2014 at 09:43 am - Permalink