Форум программистов, компьютерный форум, киберфорум
Цифровая обработка сигналов
Войти
Регистрация
Восстановить пароль
 
Рейтинг 5.00/4: Рейтинг темы: голосов - 4, средняя оценка - 5.00
0 / 0 / 0
Регистрация: 26.11.2015
Сообщений: 2
1

Определение доминантной несущей частоты в числовом ряде (БПФ, автокорреляция)

28.11.2015, 20:31. Просмотров 765. Ответов 1
Метки нет (Все метки)

Добрый день,

в зарубежной литературе встретил приведенный ниже код для определения доминантной частоты в режиме рeального времени с помощью БПФ и автокорреляции. Хотелось бы услышать мнение знатоков, есть ли ошибки или более удачные решения?

Как можно выделить две несущих частоты, длинную, в диапазоне 50 - 100, и короткую волну, в диапазоне 10 - 50?

Autocorrelation Periodogram
Construction of the autocorrelation periodogram starts with the autocorrelation
function using the minimum three periods of averaging. The cyclic
information is extracted using a discrete Fourier transform (DFT) of the
autocorrelation results.
Pascal
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
Variables:
AvgLength(3),
M(0),
N(0),
X(0),
Y(0),
alpha1(0),
HP(0),
a1(0),
b1(0),
c1(0),
c2(0),
c3(0),
Filt(0),
Lag(0),
count(0),
Sx(0),
Sy(0),
Sxx(0),
Syy(0),
Sxy(0),
Period(0),
Sp(0),
Spx(0),
MaxPwr(0),
DominantCycle(0);
 
Arrays:
Corr[48](0),
CosinePart[48](0),
SinePart[48](0),
SqSum[48](0),
R[48, 2](0),
Pwr[48](0);
 
//Highpass filter cyclic components whose periods are shorter than 48
alpha1 = (Cosine(.707*360 / 48) + Sine (.707*360 / 48) - 1) / Cosine(.707*360 / 48);
HP = (1 - alpha1 / 2)*(1 - alpha1 / 2)*(Close - 2*Close[1] + Close[2]) + 2*(1 - alpha1)*HP[1] - (1 - alpha1)*
(1 - alpha1)*HP[2];
 
//Smooth data
a1 = expvalue(-1.414*3.14159 / 10);
b1 = 2*a1*Cosine(1.414*180 / 10);
c2 = b1;
c3 = -a1*a1;
c1 = 1 - c2 - c3;
Filt = c1*(HP + HP[1]) / 2 + c2*Filt[1] + c3*Filt[2];
 
//Pearson correlation for each value of lag
 
For Lag = 0 to 48 
Begin
 
//Set the averaging length as M
    M = AvgLength;
    If AvgLength = 0 Then M = Lag;
    Sx = 0;  
    Sy = 0;
    Sxx = 0;
    Syy = 0;
    Sxy = 0;
 
For count = 0 to M - 1 
Begin
   X = Filt[count];
   Y = Filt[Lag + count];
   Sx = Sx + X;
   Sy = Sy + Y;
   Sxx = Sxx + X*X;
   Sxy = Sxy + X*Y;
   Syy = Syy + Y*Y;
End;
If (M*Sxx - Sx*Sx)*(M*Syy - Sy*Sy) > 0 Then Corr[Lag] = (M*Sxy - Sx*Sy)/SquareRoot((M*Sxx - Sx*Sx)*(M*Syy - Sy*Sy));
End;
For Period = 10 to 48 
Begin
   CosinePart[Period] = 0;
   SinePart[Period] = 0;
   For N = 3 to 48 
   Begin
       CosinePart[Period] = CosinePart[Period] + Corr[N]*Cosine(370*N / Period);
       SinePart[Period] = SinePart[Period] + Corr[N]*Sine(370*N / Period);
   End;
   SqSum[Period] = CosinePart[Period]*CosinePart[Period] + SinePart[Period]*SinePart[Period];
End;
For Period = 10 to 48 Begin
R[Period, 2] = R[Period, 1];
R[Period, 1] = .2*SqSum[Period]*SqSum[Period] +
.8*R[Period, 2];
End;
 
[B]//Find Maximum Power Level for Normalization[/B]
MaxPwr = .995*MaxPwr;
For Period = 10 to 48 
Begin
   If R[Period, 1] > MaxPwr Then MaxPwr = R[Period, 1];
End;
For Period = 3 to 48 
Begin
    Pwr[Period] = R[Period, 1] / MaxPwr;
End;
 
//Compute the dominant cycle using the CG of the spectrum
Spx = 0;
Sp = 0;
For Period = 10 to 48 
Begin
   If Pwr[Period] >= .5 Then 
   Begin
       Spx = Spx + Period*Pwr[Period];
       Sp = Sp + Pwr[Period];
   End;
End;
 
If Sp <> 0 Then DominantCycle = Spx / Sp;
0
Programming
Эксперт
94731 / 64177 / 26122
Регистрация: 12.04.2006
Сообщений: 116,782
28.11.2015, 20:31
Ответы с готовыми решениями:

Определение количества одинаковых пар чисел в числовом ряде
Определение количество одинаковых пар чисел в числовом рядя ?

О числовом натуральном ряде
Числовой ряд натуральных чисел есть арифм. прогрессия, где каждый след. член равен предыдущему + 1....

Найти закономерность в числовом ряде
Здравствуйте Есть список чисел на определенном интервале времени в один момент, происходит...

БПФ. Cоответствие частоты после преобразования
Здравствуйте, помогите пожалуйста разобраться. Я сделал быстрое преобразование фурье над массивом...

1
148 / 129 / 18
Регистрация: 29.04.2015
Сообщений: 626
30.11.2015, 05:40 2
Цитата Сообщение от RocketScience Посмотреть сообщение
в режиме рeального времени с помощью БПФ и автокорреляции
Ну вообще говоря, когда делаешь любые операции с БПФ (в том числе и корреляцию) о реальном времени надо забыть, если не проводить корреляцию за время равное времени преобразования АЦП.
Можно посмотреть непрерывный вейвлет-анализ - но там в реальности тоже прямое и обратное БПФ.
А про 2 несущих частоты - длинную и короткую -вообще не понятно :-)
0
IT_Exp
Эксперт
87844 / 49110 / 22898
Регистрация: 17.06.2006
Сообщений: 92,604
30.11.2015, 05:40

Заказываю контрольные, курсовые, дипломные и любые другие студенческие работы здесь.

БПФ бинарного сигнала, поиск частоты и фазы
Вляпался по неосторожности в ЦОС, совсем не мою область. Есть светодиод. Он мигает с частотой...

Cinterion BGS2T RS232 - ошибка "Нет несущей частоты" под Windows 7 x64, под XP 32-бит работает нормально
Модем подключен через Prolific USB-to-Serial Comm Port. Обновление драйвера Prolific не помогло....

Задача на определение суммы 10 чисел в ряде
Условие Дано 10 целых чисел. Вычислите их сумму. Напишите программу, использующую наименьшее число...

Оформить процедурой определение максимального числа в числовом массиве
Пожалуйста, помогите решить задачу в PascalABC! Оформить процедурой определение максимального...


Искать еще темы с ответами

Или воспользуйтесь поиском по форуму:
2
Ответ Создать тему
Опции темы

КиберФорум - форум программистов, компьютерный форум, программирование
Powered by vBulletin® Version 3.8.9
Copyright ©2000 - 2020, vBulletin Solutions, Inc.