[_pv 52]

вместо быстрого преобразования, посчитайте обычные интегралы Фурье s=\int{f(x)*sin(x)}; c=\int{f(x)*cos(x)}, на частоте 941 +- пара Гц с мелким шагом в 0.1Гц например чтобы не 3 в десяток-другой точек получился, на результат (s^2+c^2) натяните параболу ax^2+bx+c максимум будет в -b/2a.

ну или то же самое с автокорреляционной функцией можно попробовать.

 

исходные данные сюда выложить можете?

[serjj1333 0]

вместо быстрого преобразования, посчитайте обычные интегралы Фурье s=\int{f(x)*sin(x)}; c=\int{f(x)*cos(x)}, на частоте 941 +- пара Гц с мелким шагом в 0.1Гц например чтобы не 3 в десяток-другой точек получился, на результат (s^2+c^2) натяните параболу ax^2+bx+c максимум будет в -b/2a.

ну или то же самое с автокорреляционной функцией можно попробовать.

Подойдёт если всегда а приори известно в каком узком диапазоне будут лежать интересующие частоты. Если нет, то можно также дополнить выборку нулями и сделать БПФ, получится аналог интерполяции между частотными отсчётами. Но все эти операции размазывают пики. Есть гарантия что при таком подходе не будет ложного выделения пика частоты например в условиях шума? Ведь придётся выделять пик из группы отсчётов с отличием по амплитуде << 3 дБ. Как ни крути а разрешающая способность Фурье (любого, в т.ч. обычного ДПФ) ограничена только размерностью входной выборки и никакие операции и ухищрения её не поднимут. Сверхразрешение же даёт всегда выраженные однозначные пики, но они могут сместиться из-за шума/малого размера выборок/криворукости. Кароче смотреть надо :laughing: Но до третьего знака после запятой это сильно. Обычно речь идёт о Гц, ну максимум десятых долях Герца. Для таких точностей потребуется большой сигнал-шум и выборки порядка не менее 1е4 отсчётов. Это даже для сверхразрешения. Для Фурье будут нужны миллионы отсчётов.

Изменено 8 мая, 2015 пользователем serjj

[Alex11 3]

Кроме всего прочего, при Гауссовском окне чтобы использовать параболическую интерполяцию нужно взять логарифм амплитуды, и только его аппроксимировать параболой. И полезно еще ввести весовые коэффициенты пропорциональные амплитуде или ее степени.

[Hose 0]

В самом начале рассуждений вы назначили что будет бпф. Почему?

Интереса ради, попробуйте взять интервал для анализа кратный не 2^н, а кратный периоду искомой частоты: результат вас порадует )

[Alex11 3]

Это отличный совет, но там две частоты не кратные друг другу, да и, скорее всего, поменять частоту дискретизации в готовом устройстве не удастся.

[[email protected] 0]

Размер выборки кстати ограничен? Или можно копить до посинения?

размер выборки ограниченг самим сигналом, DTMF по стандарту идут от 40 мс сам сигнал и от 40 мс пауза между другим сигналом. В реальности же одна метка - 40-50 мс+перерыв 40-50мс.

 

вместо быстрого преобразования, посчитайте обычные интегралы Фурье s=\int{f(x)*sin(x)}; c=\int{f(x)*cos(x)}, на частоте 941 +- пара Гц с мелким шагом в 0.1Гц например чтобы не 3 в десяток-другой точек получился, на результат (s^2+c^2) натяните параболу ax^2+bx+c максимум будет в -b/2a.

ну или то же самое с автокорреляционной функцией можно попробовать.

исходные данные сюда выложить можете?

 

исходные данные: DTMF.txt , fd=44100 Гц

 

а на счет обычных интегралов Фурье, непонятно что туда подставлять, если у меня есть только отсчеты, что за f(x) и sin (x), что за x?

я параболическую интерполяцию применила, результат близок к истине, но не истина.

 

на счет автокоореляционной функции - можно поподробнее?

мне предлагают корреляционную функцию использовать, посчитать Rs=СУМ(sin(2pifn)*s(n), n=0..N и Rc=СУМ(cos(2pifn)*s(n). не непонятно становится, первое, у меня 2 частоты, значит не sin надо брать а сумму двух sin от ВЧ и НЧ. Я посчитала Rs и Rc, далее R=sqr(Rs^2+Rc^2), и ничего у меня не получается, да и не понятно, даже если я найду этот коэф корреляции, как мне это поможет найти искомую частоту? по идее я должна КАК-ТО построить график R от f и f взять близлежащие к искомой частоте, но! в формулах Rs и Rc присутсвует N, а оно определено и опять я кручусь на том, что выборка ограничена.

 

Подойдёт если всегда а приори известно в каком узком диапазоне будут лежать интересующие частоты. Если нет, то можно также дополнить выборку нулями и сделать БПФ, получится аналог интерполяции между частотными отсчётами. Но все эти операции размазывают пики. Есть гарантия что при таком подходе не будет ложного выделения пика частоты например в условиях шума? Ведь придётся выделять пик из группы отсчётов с отличием по амплитуде << 3 дБ. Как ни крути а разрешающая способность Фурье (любого, в т.ч. обычного ДПФ) ограничена только размерностью входной выборки и никакие операции и ухищрения её не поднимут. Сверхразрешение же даёт всегда выраженные однозначные пики, но они могут сместиться из-за шума/малого размера выборок/криворукости. Кароче смотреть надо :laughing: Но до третьего знака после запятой это сильно. Обычно речь идёт о Гц, ну максимум десятых долях Герца. Для таких точностей потребуется большой сигнал-шум и выборки порядка не менее 1е4 отсчётов. Это даже для сверхразрешения. Для Фурье будут нужны миллионы отсчётов.

 

впринципе известно в каком диапазоне лежат частоты, можно вначале определить грубо частоты (по стандарту есть 4 НЧ и 4 ВЧ) и уже после грубого определения начинать ухищряться точность определить.

на счет шума - это следующая задача, помехоустойчивость метода определения частот проверить.

 

Кроме всего прочего, при Гауссовском окне чтобы использовать параболическую интерполяцию нужно взять логарифм амплитуды, и только его аппроксимировать параболой. И полезно еще ввести весовые коэффициенты пропорциональные амплитуде или ее степени.

 

я делала: рассчет окна, перемножение окна с сигналом, БПФ от получившегося сигнала. и уже параболическая интерполяция. Вы предлагаете взять логарифм от |БПФ| и после этого аппрокимировать? правильно?

 

В самом начале рассуждений вы назначили что будет бпф. Почему?

Интереса ради, попробуйте взять интервал для анализа кратный не 2^н, а кратный периоду искомой частоты: результат вас порадует )

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

 

Кроме всего прочего, при Гауссовском окне чтобы использовать параболическую интерполяцию нужно взять логарифм амплитуды, и только его аппроксимировать параболой. И полезно еще ввести весовые коэффициенты пропорциональные амплитуде или ее степени.

 

прологарифмировала, результат порадовал, спасибо.

на сечт коэф, у нас же амплитуда меняется, как я могу пропорционально ей или ее степени выбрать? я поперебирала, и нашла что 0.1 - ближе всего результат дает.

[_pv 52]

размер выборки ограниченг самим сигналом, DTMF по стандарту идут от 40 мс сам сигнал и от 40 мс пауза между другим сигналом. В реальности же одна метка - 40-50 мс+перерыв 40-50мс.

исходные данные: DTMF.txt , fd=44100 Гц

 

а на счет обычных интегралов Фурье, непонятно что туда подставлять, если у меня есть только отсчеты, что за f(x) и sin (x), что за x?

я параболическую интерполяцию применила, результат близок к истине, но не истина.

 

на счет автокоореляционной функции - можно поподробнее?

 

впринципе известно в каком диапазоне лежат частоты, можно вначале определить грубо частоты (по стандарту есть 4 НЧ и 4 ВЧ) и уже после грубого определения начинать ухищряться точность определить.

на счет шума - это следующая задача, помехоустойчивость метода определения частот проверить.

про автокорреляцию это я погорячился, не поможет она тут.

а для фурье выборка пожалуй коротковата чтобы 0.001Гц вытянуть за 40мс.

 

похоже надо грубо определять частоту через БПФ или как обычно в декодерах герцелем чтобы просто понять какая цифра передаётся и знать частоту +-10Гц, а потом просто подгонять коэффициенты у суммы двух синусов наименьшими квадратами например.

d = Import["d:\\downloads\\dtmf.dat"];
fit = FindFit[d, a0*Sin[2*Pi*f0*t + p0] + a1*Sin[2*Pi*f1*t + p1], {{a0, 1}, {f0, 940}, {p0, 0}, {a1, 1}, {f1, 1200}, {p1, 0}}, t]
{a0 -> 0.305163, f0 -> 941., p0 -> 1.76526*10^-6, a1 -> 0.305164, f1 -> 1209., p1 -> 7.98036*10^-8}

 

[[email protected] 0]

про автокорреляцию это я погорячился, не поможет она тут.

а для фурье выборка пожалуй коротковата чтобы 0.001Гц вытянуть за 40мс.

 

похоже надо грубо определять частоту через БПФ или как обычно в декодерах герцелем чтобы просто понять какая цифра передаётся и знать частоту +-10Гц, а потом просто подгонять коэффициенты у суммы двух синусов наименьшими квадратами например.

d = Import["d:\\downloads\\dtmf.dat"];
fit = FindFit[d, a0*Sin[2*Pi*f0*t + p0] + a1*Sin[2*Pi*f1*t + p1], {{a0, 1}, {f0, 940}, {p0, 0}, {a1, 1}, {f1, 1200}, {p1, 0}}, t]
{a0 -> 0.305163, f0 -> 941., p0 -> 1.76526*10^-6, a1 -> 0.305164, f1 -> 1209., p1 -> 7.98036*10^-8}

 

Подскажите, это код на чем? попробовала подобную функцию в математике, не получилось, dat файл создавала так же в математике,а не просто переименованием тхт.

[petrov 6]

именно измерить.

 

немного про DTMF:

в DTMF передается 2 частоты: НЧ и ВЧ, то мало сказать что, например, 941 Гц и 1209 Гц = "*", по стандарту отклонение от передаваемой частоты может составлять 1.8% (то есть при 924.062<НЧ<957.938 Гц и 1187.238<ВЧ<1230.762 Гц примет как "*" - все верно/

 

Моя задача определить что передавалось, с точностью до 3го знака.

 

А зачем определять точность генератора на соответствие стандарту в режиме коротких посылок? Определяйте в непрерывном режиме, в принципе достаточно частотомером опорный генератор померить.

[_pv 52]

я там руками столбец со временем добавил. Чтобы с исходным txt работало надо

d = Import["d:\\downloads\\dtmf.txt", "Data"][[;; , 1]];
d = Transpose[{Table[i/44100, {i, Length[d]}], d}];
fit = FindFit[d, a0*Sin[2*Pi*f0*t + p0] + a1*Sin[2*Pi*f1*t + p1], {{a0, 1}, {f0, 940}, {p0, 0}, {a1, 1}, {f1, 1200}, {p1, 0}}, t]

[[email protected] 0]

А зачем определять точность генератора на соответствие стандарту в режиме коротких посылок? Определяйте в непрерывном режиме, в принципе достаточно частотомером опорный генератор померить.

 

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

[petrov 6]

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

 

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

 

[[email protected] 0]

я там руками столбец со временем добавил. Чтобы с исходным txt работало надо

d = Import["d:\\downloads\\dtmf.txt", "Data"][[;; , 1]];
d = Transpose[{Table[i/44100, {i, Length[d]}], d}];
fit = FindFit[d, a0*Sin[2*Pi*f0*t + p0] + a1*Sin[2*Pi*f1*t + p1], {{a0, 1}, {f0, 940}, {p0, 0}, {a1, 1}, {f1, 1200}, {p1, 0}}, t]

 

спасибо большое, заработало с этим файлом, завтра буду разбираться с другими "подопытными"...

 

 

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

 

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

[blackfin 16]

именно измерить.

 

в DTMF передается 2 частоты: НЧ и ВЧ, то мало сказать что, например, 941 Гц и 1209 Гц..

 

Моя задача определить что передавалось, с точностью до 3го знака.

Касательно "измерить".. и именно "с точностью до 3го знака".

 

Что как-бы намекает на необходимую точность измерения: 0,001/1209 = 0,8 ppm.

 

Недавно как раз обсуждали похожий случай:

Да дело даже не в алгоритме вычисления.

 

Спектр, который мы видим после АЦП, равен произведению АЧХ всего тракта (включая АЦП) на спектр измеряемого сигнала: S_ацп(f)=K(f)*S_вх(f).

 

Но спектр самого входного сигнала S_вх(f) из-за того, что импульс короткий, оказывается достаточно широким ~0,5 МГц и, как следствие, будет иметь достаточно пологий максимум.

Ну и немного классики:

Этот вопрос был детально изучен в классической статье

Впервые была дана оценка предельной точности измерений параметров (комплексной) синусоиды и предложен алгоритм их измерения, достигающий этих пределов

 

Single-Tone parameter estimation from Discret-Time observations

То есть, цифру то получить можно.. и даже "до 3го знака"..

 

Но вот будет ли эта цифра иметь какое-то отношение к действительности? Сомнительно..

 

[petrov 6]

Ну и немного классики:

 

Причём это для чистой комплексной синусоиды, а тут фильтры для разделения комплексных тонов друг от друга дадут ошибки и смещения больше 3-го знака, не считая нестабильности опорных генераторов.