Why does Fourier give a shifted frequency?

Question

I have a signal that I want to identify the frequencies in it, I used the Fourier function but I can't get the frequency correctly. Here is a simplified example:

dt = 1/100;
ls = Table[0.1 Cos[30 x] + 2 Sin[x]^2, {x, 0, 200 dt, dt}];
ListPlot[ls, Mesh -> All, MeshStyle -> Red]

and the Fourier transform

ListPlot[Abs[Fourier[ls]]^2, PlotRange -> {{0, 10}, {0, 1}}, 
 DataRange -> {0, 1/dt}, 
 FrameLabel -> {"Frequency", "Intensity"}, 
 Mesh -> All, MeshStyle -> Red, GridLines -> {{30/(2 π)}, None}]

Why do I get the peak not at the original frequency 30/(2Pi), but with a frequency shift? Did I make a terrible mistake? What's the correct way to recover the original frequency using Mathematica's signal processing features?

I tried to padding zeros but still have a frequency shift.

ListPlot[Abs[Fourier[PadRight[ls, 2000]]]^2, 
 PlotRange -> {{0, 10}, {0, .1}}, DataRange -> {0, 1/dt}, 
 FrameLabel -> {"Frequency", "Intensity"}, 
 Mesh -> All, MeshStyle -> Red, GridLines -> {{30/(2 π)}, None}]

Edit

I'm still not convinced that this problem is due to that there are too few periods contains in the signal, nor that number of periods contained in the signal is not an integer. Consider this signal

ls = Table[0.1 Cos[30 x], {x, 0, 200 dt, dt}];

it contains the same number of periods and the number of periods is not an integer, it gives a peak that are not centered,

ListPlot[Abs[Fourier[ls]]^2, PlotRange -> {{0, 10}, {0, 1}}, 
 DataRange -> {0, 1/dt}, 
 FrameLabel -> {"Frequency", "Intensity"}, 
 Mesh -> All, MeshStyle -> Red, GridLines -> {{30/(2 π)}, None}]

but padding zero helps

ListPlot[Abs[Fourier[PadRight[ls, 2000]]]^2, 
 PlotRange -> {{0, 10}, {0, .1}}, DataRange -> {0, 1/dt}, 
 FrameLabel -> {"Frequency", "Intensity"}, 
 Mesh -> All, MeshStyle -> Red, GridLines -> {{30/(2 π)}, None}]

Answer 1

I believe the frequency mismatch arises because the endpoints of your 200 point series are offset. The first point has amplitude 0.1, the last 1.5584. As others mention, the Fourier transform assumes periodicity. So the signal you are transforming has a sine component, a cosine component, and a step function offset of the first and last points. The Fourier transform of a step function is

FourierTransform[UnitStep[t], t, x]

which evaluates to roughly $1/x+\delta[x]$, where $\delta$ is the delta function, and $x$ is the frequency variable. This spectrum is peaked at $x=0$ and drops off with increasing $|x|$. Thus, the time-domain step function contributes a sloping frequency spectrum to your frequency-domain delta-functions from the sine and cosine. The sloped spectrum increases the amplitudes of frequencies less than $30/(2\pi)$ more than the amplitudes of frequencies greater than $30/(2\pi)$.

If your signal were sampled with 297 points, the first and last points would be 0.1 and 0.100489. With the endpoints roughly matched, the spectrum peak is not shifted from the expected $30/(2\pi)$.

EDIT

1) You comment that the partial period only gives the finite peak width. With respect, you cannot ignore the frequency response of the step function inherent in your data. It is in the spectrum, period. The "Edit" example you give has a smaller step-function offset, so its influence is reduced.

Plot[Abs[FourierTransform[UnitStep[t], t, w]]^2, {w, -1, 1}]

2) As you say, in a "real situation" signals cannot always contain an integer number of periods. In practice, signals are tapered on both ends with, for example, a cosine-squared bell function to eliminate step-function artefacts. Thus, an integral number of periods is not required. When I taper your 200 point function with the following very rough function, the DFT amplitudes at the sampled frequencies on either side of $f=30/(2\pi)$ are much more equal.

Block[{dt=1/100, ls, a=0.2, nn=201},
   ls = Table[0.1 Cos[30 x], {x, 0, 200 dt, dt}];
   ls = ls * Table[(1 - a)/2 - 1/2 Cos[2 \[Pi] n/(nn - 1)] + 
                   a/2 Cos[4 \[Pi] n/(nn - 1)], {n, 0, 200}];
   ListPlot[Abs[Fourier[PadRight[ls,2000]]]^2, 
      PlotRange->{{0,10},{0,0.015}}, DataRange->{0, 1/dt}, 
      Joined->True, FrameLabel->{"Frequency","Intensity"}, 
      Mesh->All, MeshStyle->Red, GridLines->{{30/(2 \[Pi])},None}] ]

3) The DFT does what it does, and its answers are not "wrong"; however, its answers require some interpretation in light of limited data, non-integral periods, and mismatched endpoints. (Thank you @LeoFang)

Answer 2

Fourier Transform is based on assumptions of periodicity related to the duration of the data. If you choose a misleading duration, you will get misleading results.

The duration (200 dt) of your array ls is not a multiple of the periods of your waveforms; this introduces artefacts arising from the Fourier transform of the 'top hat' function of width (200 dt). No amount of padding will remove this completely.

Try with the duration of ls much closer to a multiple of the duration of the fundamental; results will be closer to expected.

A dramatic improvement is seen with 314 points, i.e. about 100 times \Pi

Answer 3

A windowing function should help:

ls2 = (ls - Mean[ls]) Array[TukeyWindow, Length@ls, {{-0.5, 0.5}}];

ListLinePlot[ls2]

ListLinePlot[Abs[Fourier[ls2]]^2, PlotRange -> {{0, 10}, {0, 1}}, DataRange -> {0, 1/dt}, 
 FrameLabel -> {"Frequency", "Intensity"}, Mesh -> All, MeshStyle -> Red, 
 GridLines -> {{30/(2 π)}, None}]

Answer 4

This is just a long supplement comment to @James Cunnane's answer which is correct. Try

dt = 1/100;
T=Pi;
ls = Table[0.1 Cos[30 x] + 2 Sin[x]^2, {x, 0, T, dt}];
ListPlot[ls, Mesh -> All, MeshStyle -> Red]
ListPlot[Abs[Fourier[ls]]^2, PlotRange -> {{0, 10}, {0, 1}}, 
DataRange -> {0, 1/dt}, FrameLabel -> {"Frequency", "Intensity"}, 
Mesh -> All, MeshStyle -> Red, GridLines -> {{30/(2 \[Pi])}, None}, 
Frame -> True, Joined -> True]

Note the only thing changed is the duration T whose inverse (1/Pi) defines the grid size in the Fourier spectrum. As a result, the desired frequency 30/Pi is now an integer multiple of the grid size and therefore can be captured correctly.

EDIT

Because the grid size 1/T=1/Pi is not fine enough, it is less easier to observe the lower frequency 2/2Pi (not 1/Pi because it's a sine square). Try to increase T to a larger value would resolve this issue if you like:

dt = 1/100;
T = 5 Pi;
ls = Table[0.1 Cos[30 x] + 2 Sin[x]^2, {x, 0, T, dt}];
ListPlot[Abs[Fourier[ls]]^2, PlotRange -> {{0, 1}, All}, 
DataRange -> {0, 1/dt}, FrameLabel -> {"Frequency", "Intensity"}, 
Mesh -> All, MeshStyle -> Red, GridLines -> {{2/(2 \[Pi])}, None}, 
Frame -> True, Joined -> True]

Last comment: The reason you see a very high peak at zero frequency is due to the constant term when you write Sin[x]^2 as (1 - Cos[2 x])/2.