Форум разработчиков электроники ELECTRONIX.ru > Определение постоянной составляющей синусоиды...

Помощь - Поиск - Пользователи - Календарь

Полная версия этой страницы: Определение постоянной составляющей синусоиды...

Форум разработчиков электроники ELECTRONIX.ru > Цифровая обработка сигналов - ЦОС (DSP) > Алгоритмы ЦОС (DSP)

Страницы: 1, 2

Stanislav

Nov 27 2008, 00:23

Цитата(mvb @ Nov 27 2008, 03:07)

Для реализации метода, предложенного Станиславом взял Curve Fitting toolbox из матлаба. При бесконечном ОСШ он точно выдает искомую величину. В присутствии шума решение начинает гулять. Только надо выставлять область ограничений для решения не бесконечную, иначе минимумов у целевой функции получается несколько.

На все это я потратил 1 час 7 минут.

Вероятно, Вы неплохо изучили Матлаб, потому, что у меня ушло бы времени больше, хотя работать с ним приходится.

Если не трудно, выложите программу, интересно посмотреть. Постараюсь сделать и свою, если время будет.
Конечно, оценка в шумах будет гулять. Но она будет, и она будет оптимальной, в каком угодно смысле.
Кроме того, в условии о шумах ничего не говорилось. Есть сигнал, нужно определить постоянную составляющую. Формально поставленная задача, имеющая сугубо формальное решение. О ньюансах автор как-то не сообщал.
В условиях больших помех сильно может помочь априорное знание хотя бы одного из параметров - амплитуды или фазы гармонического сигнала, например.
Впрочем, это банальность. Знание вида функции уже даёт многое.

ЗЫ. Кстати, фитнуть можно и не гармонической функцией, а её разложением в ряд. Вероятно, будет менее напряжённо вычислительно.

GetSmart

Nov 27 2008, 04:12

Цитата(mvb)

2. Как здесь можно использовать дпф или акф?

ДПФ более-менее точно измеряет частоту/период, когда в блоке находится много периодов. АКФ достаточно два периода. Хотя имея гарантированно два периода частоты можно "ручками" отсканировать блок (свёртками sin/cos) на дробные частоты. 729 кажется знал ещё какой-то хороший метод определения дробных частот. Только он пропал куда-то.

mvb

Nov 27 2008, 21:10

Вот мой код, без претензий на робастность и скорость:

Код

clc;

N = 100;        % кол-во отсчетов
n = 1:N;
d = .1;         % постоянная составляющая
p = rand()*pi;  % случайная фаза
f = 1/(4*N);    % частота (в N отсчетах 1/4 периода)
SNR = 40;       % ОСШ в дБ

% собсно анализируемые отсчеты
s = sin(2*pi*f*n+p) + d + randn(1, N)*(d+sqrt(2))*10^(-SNR/20);

st_ = [ .5 2*pi*f 0 0 ];
fo_ = fitoptions('method','NonlinearLeastSquares','Lower',[0 0 -1.7 -1],'Upper',[5 0.10000000000000001 1.7 0.10000000000000001]);
set(fo_,'Startpoint',st_);
ft_ = fittype('a*sin(b*x+c)+d',...
     'dependent',{'y'},'independent',{'x'},...
     'coefficients',{'a', 'b', 'c', 'd'});
cf_ = fit(n',s',ft_,fo_);

a = coeffvalues(cf_);

a(4)

Да, у меня тоже была мысль о разложении синуса в ряд, или как-то перейти к полиному Чебышева, но тогда постоянная составляющая будет спрятана в свободном члене полинома, и как ее оттуда вычленить я не знаю.

2GetSmart: ну если в распоряжении есть больше (м.б. больше равно) одного периода, то все становится понятно

Михаил_K

Nov 28 2008, 05:00

Цитата(RadioJunior @ Nov 26 2008, 12:11)

Это кризис возраста или врожденный комплекс Наполеона?...

Нет. Вы ошибаетесь. Нет никакого кризиса возраста и нет никакого Наполеона.
Все оказалось намного банальнее.

Цитата(Stanislav @ Nov 27 2008, 00:29)

Можете считать меня ассенизатором.

Если человек потратил жизнь на борьбу с гав..ом. Сами понимаете....

Stanislav

Nov 28 2008, 06:55

Цитата(mvb @ Nov 28 2008, 00:10)

Вот мой код, без претензий на робастность и скорость:

Код

..................................

Понятно. Собственно, у меня вышло бы так же.
Попробую его "усовершенствовать", когда время будет.

Цитата(mvb @ Nov 28 2008, 00:10)

Тоже подумаю.

Цитата(mvb @ Nov 28 2008, 00:10)

2GetSmart: ну если в распоряжении есть больше (м.б. больше равно) одного периода, то все становится понятно

А мне вот не совсем понятно, что собирается делать GetSmart, если в его распоряжении будет кусок, заведомо больший только одного периода.

Цитата(Михаил_K @ Nov 28 2008, 08:00)

Нет. Вы ошибаетесь. Нет никакого кризиса возраста и нет никакого Наполеона.
Все оказалось намного банальнее.

Если человек потратил жизнь на борьбу с гав..ом. Сами понимаете....

Да не потратил ещё. А трачу.

mvb

Nov 28 2008, 10:03

Если бы у меня был бы фрейм с достаточным количеством отсчетов, я бы сначала попробовал такой способ:

1. Определил бы период синуса (в первом приближении с пом. автоковариационной функции, дальше покопал бы продвинутые методы вычисления спектра)

2. Далее скомпенсировал бы ошибку в ДПФ на 0 частоте от того, что частота синуса не кратна 2*Pi/N
Насколько я помню отклик гармонической функции при ДПФ пропорционален sinc, тогда получается что-то вроде этого:
модуль постоянной составляющей = ДПФ(0) - sinc(0 - fsin)*ДПФ(fsin)
(fsin -- измереная частота синуса)

Наверняка вы найдете здесь ошибки, укажите мне на них пожалуйста!

Stanislav

Nov 28 2008, 11:45

Цитата(mvb @ Nov 28 2008, 13:03)

Щас времени маловато. Вот некоторые соображения:
1. Автоковариационный метод даст большую ошибку в определении частоты, даже если в окне будет несколько периодов. Только по одному периоду, ПМСМ, работать не будет. Собственный простой способ предложу попозже.
2. Не учтена амплитуда гармонического сигнала и его фаза, поэтому ДПФ будет неоднозначной по отношению к этим параметрам. Правда, после оценивания частоты её найти их не так уж сложно. Или я не понял что-то?

Впрочем, способ интересный, есть в нём что-то эдакое...

ЗЫ. А-а, понял, наконец.

Речь идёт, конечно, о модулях спектра сигнала. Что ж, посмотрю на досуге, действительно стало интересно.

Только с фазой, по-моему, всё равно засада...

Лучшим же методом считаю, как и в первом случае, использование приближения входного сигнала известной функцией (фиттинг). Правда, здесь ещё подумать надо: возможно, есть более простой и точный способ.

GetSmart

Nov 28 2008, 13:08

Цитата(mvb)

2. Далее скомпенсировал бы ошибку в ДПФ на 0 частоте от того, что частота синуса не кратна 2*Pi/N
Насколько я помню отклик гармонической функции при ДПФ пропорционален sinc, тогда получается что-то вроде этого:
модуль постоянной составляющей = ДПФ(0) - sinc(0 - fsin)*ДПФ(fsin)
(fsin -- измереная частота синуса)

Я например многое из этого не понял.
ДПФ по какому блоку делается? По одному или по скольки буферам?

sinc(0 - fsin)*ДПФ(fsin) что из себя представляет? Это свёртка?

Кстати, ДПФ(0) - это простое среднее арифметическое блока. Как бы амплитуда постоянки (0 Гц), которая в отличие от остальных спектральных компонентов даже фазы не имеет. У меня сразу мысли возникают - почему такая дискриминация, ведь все частоты должны быть равноправные? Но это к этой теме отношения не имеет

Да, и забыл совсем

Цитата(Stanislav)

...Надеюсь, что попрошайничество на сей раз Вам не поможет снять взыскание...

Бред сивой кобылы

Stanislav

Nov 28 2008, 20:59

Цитата(GetSmart @ Nov 28 2008, 16:08)

Да, и забыл совсем
Бред сивой кобылы

А вот почитайте, позорник. Память, видимо, короткая, да и достоинство на нуле...

mvb

Nov 30 2008, 09:21

Попробовал на досуге определить период синуса в матлабе через автоковариацию. При вышеописаных условиях получается ошибка порядка 10-15%, чего явно не достаточно. Чем больше периодов умещается в окне, тем меньше относителльная погрешность.
Вобщем да, имея один период определить частоту получается плохо.

Станислав: да, под ДПФ(х) я подрузомевал модуль, извините, за то что так криво написал. На самом деле под ДПФ(fsin) я имел ввиду просто амплитуду синусоиды, из ДПФ ее вычленить наверное не получится...

GetSmart: sinc(0 - fsin)*ДПФ(fsin) -- нет, это просто sinc, умноженный на амплитуду синусоиды.
Модуль ДПФ(0) строго равен среднему арифметическому, аргумент же строго равен нулю. Все просто: вычисляется скалярное произведение сигнала с комплексной экспонентной от нуля --
ДПФ(0) = <s, exp(j*0)>/N
exp(j*0) = 1
Т.о. получается что действительная часть ДПФ(0) равна просто сумме сигнала деленного на N, а мнимая нулю.

Для просмотра полной версии этой страницы, пожалуйста, пройдите по ссылке.