Detection of sinusoids in the discrete sequence

10 views (last 30 days)
Good day! Please excuse me for my English. I have the following problem: There is a period of a discrete signal. This signal is the sum of sinusoidal signals (sinusoids). These sinusoids not always begin and end at the beginning and end of the period. As an example, I suggest the following (see figure) Period length = 300. There are two sinusoidal signals: one signal all over the period (300), the second - only half of the period (150). Question: how to detect both signals in this example? That is, how to determine the frequency of these two signals. If both signals were all over the period, the usual Fourier (fft) would perform this task. But in this example, Fourier did not determine the frequency correctly. Thank you in advance! http://i064.radikal.ru/1210/90/ed1741e392a6.jpg
  2 Comments
Greg Heath
Greg Heath on 5 Oct 2012
The link doesn't work. I get the error message
The page you were looking for doesn't exist.
You may have mistyped the address or the page may have moved.
Will you please submit a short piece of code that will generate the signal and produce the plot?
TIA
Andrey
Andrey on 8 Oct 2012
Edited: Andrey on 8 Oct 2012
Sorry. I corrected

Sign in to comment.

Answers (4)

Greg Heath
Greg Heath on 6 Oct 2012
If both sampling frequencies are equal to Fs = 1/dt, then the folowing conditions should be met
1. The signal amplitudes are comparable so that neither spectrum peak is comparable to the other signal spectrum's peak sidelobe.
2. Each signal length is longer than a fundamental period:
t1max > 1/f1 and t2max = 2*t1max > 1/f2
3. If x1 is zeropadded to length t2max, the peak spectra can be resolved when
df = 1/(t2max+dt) < abs(f1-f2)
Thank you for formally accepting my answer.

Andrey
Andrey on 8 Oct 2012
Edited: Andrey on 8 Oct 2012
Thanks for the answers! This example was very simple. In this example I would draw your attention to the following: If the sinusoid is not defined on the whole window and when it intersects with other unknown sinusoid, the Fourier is not exactly determine the frequency. Of course, we can take a smaller window, but .... How about this example: Note that I am defining period of sinusoid (not the frequency). That is, how many points should be sinusoid to make the full period. Sorry, but I feel so comfortable. In my example, sinusoids with three different lengths, with different periods (frequencies). It are put together and create a new signal. Objective: to determine the period (frequency) of the source sinusoids.
function [ ] = Test()
x = 1:100;
y1 = CreateSin(x, 20, 78, 30);
y2 = CreateSin(x, 13, 60, 45);
y3 = CreateSin(x, 3, 100, 60);
y = y1 + y2 + y3;
subplot(4, 1, 1);
plot(x, y1);
axis([0, 100, -2, 2]);
subplot(4, 1, 2);
plot(x, y2);
axis([0, 100, -2, 2]);
subplot(4, 1, 3);
plot(x, y3);
axis([0, 100, -2, 2]);
subplot(4, 1, 4);
plot(x, y);
axis([0, 100, -2, 2]);
end
% function create sinusoid
% x - domain
% begin - point from which to begin a sinusoid
% finish - the point at which ends sinusoid
% period - sinusoid period
function [ y ] = CreateSin(x, begin, finish, period )
y = zeros(1, length(x));
y(begin:finish) = y(begin:finish) + sin(x(begin:finish) * 2 * pi / period);
end
This example is not so simple. Because there is sinusoids is almost completely overlap. And sinusoids very short. This makes the task very difficult. I would be grateful for any ideas

Andrey
Andrey on 24 Oct 2012
Edited: Andrey on 24 Oct 2012
Here is a fantastic example:
Imagine a planet. The water level at the beach depends on the moon. We all know that the moon pulls the water. And on the last day of our planet 100 Earth hours. Here's the plot: http://s019.radikal.ru/i611/1210/2a/7ba1b19ff15d.jpg
Now imagine that we have a second moon, which rotates twice as fast around the planet (the period is 50 Earth hours) Moons can be at one point, or can on the different sides of the moon. Therefore, the gravitational force of the two moons are summarized as two vectors, so we get the following picture. http://s005.radikal.ru/i210/1210/48/37f102edf9bd.jpg
And that's not all. Second Moon permanently disappear somewhere. No one knows where. Just disappears. So when the second moon is absent the effect of gravity on the water only to the first moon. Our second moon appears stable at 27 o'clock, and disappears by 64 o'clock http://s14.radikal.ru/i187/1210/d2/86c0e6daae87.jpg
Now imagine that we want to know about when the moon disappear. We also need to know how fast to spin the moon around the planet. And all this we can know, having the only the graph (the water level):
Who has ideas how to do this algorithmically?
function [ ] = Run()
close all;
clear all;
clc;
moon1 = CreateMoon(100, 1, 100, 100);
moon2 = CreateMoon(100, 27, 64, 50);
moon3 = CreateMoon(100, 50, 70, 75);
subplot(4, 1, 1);
plot(moon1, '.');
axis([0, 100, -2, 2]);
subplot(4, 1, 2);
plot(moon2, '.');
axis([0, 100, -2, 2]);
subplot(4, 1, 3);
plot(moon3, '.');
axis([0, 100, -2, 2]);
subplot(4, 1, 4);
plot(moon1+moon2+moon3, '.');
axis([0, 100, -2, 2]);
end
% function create moon
% dayLength - the number of hours in a day
% appearance - time of appearance
% disappearance - time of disappearance
% period - time for a complete orbit around the planet
function [ y ] = CreateMoon(dayLength, appearance, disappearance, period )
x = 1:dayLength
y = zeros(1, length(x));
y(appearance:disappearance) = y(appearance:disappearance) + sin(x(appearance:disappearance) * 2 * pi / period);
end

Andrey
Andrey on 24 Oct 2012
Edited: Andrey on 25 Oct 2012
I am interested in a method that, in my view, better helps to solve my problem: Empirical Mode Decomposition (EMD) Look what results achieved by the example with moons http://s56.radikal.ru/i154/1210/cd/322b7fd4904b.jpg
Unfortunately, we were able to identify only one moon (which is the full width of the window). But we were able to algorithmically find at least one of the signal composes. I think this is a good start to finding a solution.
moon1 = CreateMoon(1000, 1, 1000, 100);
moon2 = CreateMoon(1000, 340, 1000, 33);
moon3 = CreateMoon(1000, 1, 700, 70);
Look at the results and compare them to the sine wave. Of course, there is no 100% accurate. But we're just getting started. I marked in red pencil, which would have placed the boundary of the sine wave and their periods.
Can yCan you tell me something about EMD (or other similar methods) that would help me move on?ou tell us anything about the EMD, that would help me move on?

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!