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

[Kseniia]

Всем доброго времени суток!

Задача: есть DTMF метка, необходимо определить частоты которые передавались, с точностью до 3го знака после запятой. (DTMF - передается сумма двух гармоник, НЧ и ВЧ)

Что проделано:
выделяю из сигнала DTMF, fd=44100 Гц, для того, чтобы N была кратна двум, и равнялась 2048, немного обрезаю сигнал до 0.046 секунд (было 0.05 сек), таким образом T=1/fd=2.268*10^5.
Беру БПФ, строю график зависимости уровня(?)(20 log(sqwr(ReБПФ^2+ImБПФ^2)/fd-N), дБ от частоты, Гц (прикрепляю скрин, freqn=n*fd/N, n=0,...,N-1).

Пики приходятся на 947,46Гц и 1205,9Гц, хотя должны быть 941Гц и 1209Гц. такая разница по видимому из-за огромного ШАГА по частоте, в 21,55 Гц, то есть передаваемые частоты не кратны нашему "шагу".

Что делать дальше - не понимаю.
Пробовала инетрполяцию Лагранжа - чушь полная, не получается.
Есть предложение использовать коэффициенты корреляции, относительно синусоидального и косинусоидального сигналов, но результатов пока не дал, да и не до конца понятно, что с этими коэф делать.

вопрос: Как математически можно определить частоты?

[serjj]

Вам тут уже советовали алгоритм MUSIC. Он дает очень хороший результат по точности. Для него не надо генерировать базисных функций, сетка частот заранее не задана. Процедура БПФ тоже не требуется.
В матлабе функция music возвращает оценку частот. Алгоритм основан на разложении оценки автокорреляционной матрицы на собственные вектора, по-моему есть быстрые алгоритмы. Посмотрите в сторону него.
Так же можно попробовать модернизированный ковариационный метод, который есть так же в матлабе - функция armcov. Она возвращает оценку АР-параметров модели сигнала, чтобы посчитать частоты надо вычислить комплексные корни характеристического полинома построенного по АР-параметрам, а потом посчитать их аргумент. Процедура БПФ так же не требуется, сетка частот не задана. Если у вас две комплексные экспоненты, то этот метод может быть эффективнее, чем БПФ и в плане точности и в плане быстродействия.

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

Насчет точности измерения. В книге Кея "Kay S.M. Fundamentals of Statistical Signal Processing: Estimation Theory (v.1). Prentice Hall, Upper Sadle River, NJ, 1998. 595 p." на 542 странице есть выражение границы Крамера-Рао для оценки частоты комплексной экспоненты на фоне шума. Попробуйте подставить в формулу Ваши условия, может быть такая точность и теоретически то не достижима)))

У меня знакомый копал в свое время в сторону АР-моделей и метода Юли-Уолкера, дают хороший результат для малого количества спектральных компонент, это как раз примеры параметрических методов. Из непараметрических еще ESPRIT можно глянуть. Фурье имеет фундаментальное ограничение в точности, связанное с размером выборки и тут как ни крути, это ограничение никуда не денется. Не знаю, зачем себя ограничивать в данном случае, оговаривая именно примение Фурье. Речь ведь идет не о жестком реальном времени с ограниченными ресурсами.
А книжкой не поделитесь?)
ЗЫ: хорошая идея по поводу Крамера-Рао, как раз с этого и стоило бы начинать дискуссию

[blackfin]

.. хорошая идея по поводу Крамера-Рао, как раз с этого и стоило бы начинать дискуссию

Так шума нет! Причем тут граница Крамера-Рао? У Kseniia идеальный сферический DTMF в MATLAB'е!

И тем не менее, даже в отсутствии шума оценки частот неверны!

Но они и не будут верны, пока не придет понимание того, что реально существующий спектр суммы двух гармоник умноженных на любое окно ограниченной длительности,
имеет максимумы отнюдь не на тех частотах из которых эта сумма гармоник была сформирована!

Достаточно наглядно это можно увидеть уже просто сложив спектры двух радиоимпульсов имеющих близкие частоты несущих.
При таком сложении максимум в спектре каждого из радиоимпульсов сместится в сторону максимума соседнего радиоимпульса.

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

Разложив эти спектры в ряд в окрестности точки f₁ и отбросив все члены ряда Тейлора выше второй степени, находим:

S₁(f) = A*(f - f₁)² + C + ...

S₂(f) = B*(f - f₁) + ...

В этих разложениях предполагается, что спектр S₁(f) имеет максимум в точке f₁, а спектр S₂(f) представлен в точке f₁ линейным членом, так как максимум его спектра находится в точке f₂.

Сложив оба спектра, находим суммарный спектр вблизи точки экстремума функции S₁(f):

S_сум(f) = A*(f - f₁)² + B*(f - f₁) + С.

Теперь, чтобы найти максимум функции S_сум(f), мы находим её производную и приравниваем её нулю:

S'_сум(f) = 2*A*(f - f₁) + B = 0.

Откуда:

f_max = f₁ - B/(2*A).

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

Поэтому, если частота несущей из которой был получен сигнал DTMF равнялась f₁, тем не менее любой спектральный метод оценивания даст в результате смещенную оценку f₁±Δf даже и при отсутствии шума вообще.

Точную оценку обеих частот f₁ и f₂ в этом случае можно получить только лишь зная амплитуды и фазы обеих гармоник из которых был получен сигнал DTMF..