March 18, 2006

Fourier fitting with gnuplot

gnuplot can be used to fit a fourier series to a series of data points. In my case it was an observation log for the star BE Lyncis from the Bossom Observatory, curtsey of Geir Klingenberg.



Firstly one needs to define a time function t(x) such that as x progresses from its minimum value to its maximum value t(x) ranges from 0 to 2PI. In my case minimum value of x was 2453416.20172 and maximum value was 2453416.42447. So my t(x) was constructed as follows:



t(x)=(x-2453416.20172)/0.22275*2*pi


We now need to define a fourier series, f(x). It takes the generic form of:



f(t)=c+a0*cos(t)+a1*sin(t)+a2*cos(2*t)+a3*sin(2*t)+
...
+{a#*cos((#+1)*t)+a(#+1)*sin((#+1)*t)}


Where #=0,1,2,3... infinity, and is the degree of the fourier series. Now we can use gnuplot's fit function to determine the values of all co-efficients and the constant c:



fit f(t(x)) 'be_lyn-2005.txt' using 1:2 via c,a0,a1,a2,...,a#,a(#+1)


A quick explanation: we are trying to fit f(t) where t=t(x) using data from the file be_lyn-2005.txt with data from the first and second column. The order of 1:2 is significant as it is in the form of x:y. Thus 2:1 means the first column is y value, second column is x value. (This is from experience, if any one knows better, tell me!)



Co-efficient values and the constant will be computed and shown, along with various statistic data. An plot of the fourier series and the observed data is shown below.




The plain text file be_lyn-2005.txt has 3 columns: JD, instrumental magnitude, and error estimate, the later of which was not used.


Cheers,

Steve