В предыдущей главе мы рассматривали ситуацию, когда частота
колебания равнялась m/T, где m - целое.
Теперь рассмотрим ситуацию, когда это не так. Положим, что
m = n + q, где n - целое и 0 < q < 1.
Воспользуемся формулой (35). Поскольку первые несколько
условий для нецелого m не выполняются, то
остается последняя, самая сложная формула, помеченная словами
"Для остальных k":
Подставим в эту формулу n + q
вместо m и выполним упрощения,
воспользовавшись формулой (34) и введя обозначение ρ = 2πj.
Итого:
, ρ = 2πj (37)
Теперь построим график функции, чтобы понять, как
она себя ведет. Ниже показана трехмерная поверхность.
По горизонтальной оси отложено k,
по вертикальной |Xk| и по
оси, уходящей вглубь плоскости, отложено q
от 0.01 до 0.99.
Рис. 1.
На рисунке видно два ярко выраженных ребра. Первое из них
всегда приходится на k = n и k = n + 1.
Второе ребро получается в результате зеркального эффекта.
Высота пика наименьшая в окрестности q = 0.5.
А наибольшая в окрестности q = 1 и q = 0
- то есть при целочисленном m.
К сожалению, пик не является единственным ненулевым коэффициентом
Фурье. Рядом с ним есть множество меньших, но не нулевых величин.
Если при целочисленном m можно наблюдать единственную
полоску, то при нецелом m = n + q эта полоска
размазывается:
Рис. 2.
На рисунке приведена практическая ситуация. Это - ДПФ для звука,
содержащегося в обычном WAV-файле. Высота синего штриха цвет соответствует
|Xk|. Исходный сигнал содержал
ноту "ля" второй октавы с частотой 440 гц и фазой в 90 градусов. ДПФ было
выполнено для N = 1000. Однако частота дискретизации
звука в WAV-файле составляла 44100 Гц, так что период дискретизации
был равен T = 1000/44100 секунд и из формулы m/T = 440
получим m = 440*(1000/44100) = 9.97,
то есть, не целое. В результате ярко выраженный пик окружают дополнительные
ненулевые значения.
На следующем рисунке:
Рис. 3.
показана "хорошая" ситуация, когда частота исходного звука составляла 441 Гц,
и m = 441*1000/44100 = 10, то есть целое. Вы видите только один
ненулевой отсчет.
Этот эффект будем называть эффектом размазывания. Вы видите, что он
определяет погрешность, с которой можно найти частоту
исходного колебания. Погрешность равна 1/T. При достаточно
большом отклонении m от целого эффект может быть очень заметен.
Например ниже вы видите ДПФ для сигнала, соответствующего ноте "ля-бемоль":
Рис. 4.
Точнее можно попытаться определить параметры m,
A и φ численными методами.
Для поиска φ следует учесть, что изменение
A не повлияет на комплексную фазу (аргумент)
коэффициентов Xk.
В самом деле, мы можем представить коэффициенты в виде:
Xk = (A/2)Z(m,&phi)k,
где Z(m,&phi)k - комплексное число, не зависящее
от действительного числа A, но зависящее
от m и φ:
Z(m,&phi)k =
Также не зависит от A отношение коэффициентов
Xk/Xl = Z(m,&phi)k/Z(m,&phi)l.
Это значит, что у нас есть две целевые функции, с помощью которых
мы можем найти частоту m/T и фазу φ.
Возьмем Xk, максимальное по модулю. Если соседние
отсчеты Xk-1 и Xk+1
равны нулю, то у нас нет эффекта размазывания и параметры
восстанавливаются так, как описано в предыдущей главе. На самом деле
нам придется сравнивать не с нулем, а с некоторым малым числом, поскольку
некоторая погрешность при вычислении ДПФ неизбежна.
Теперь, когда мы убедились в наличии эффекта размазывания, попробуем найти
m и φ после чего восстановим A
по формуле: A = |2Xk / Z(m,&phi)k|.
Сначала найдем два максимальных
отсчета Xk и Xk+1.
Теперь мы знаем, что искомое m лежит на
интервале (k, k+1).
Для нахождения m и φ нужно численно решить
систему уравнений:
|Z(m,&phi)k/Z(m,&phi)k+1| = |Xk/Xk+1|
Arg(Z(m,&phi)k) = Arg(Xk)
Здесь только две неизвестных - m и φ, неизвестная A исключена.
Систему можно решить итерационно. Для этого сначала фиксируем φ и подбираем m,
которое удовлетворяет первому уравнению системы. Потом - наоборот фиксируем найденное m и
подбираем φ, которое удовлетворяет второму уравнению системы. Потом снова возавращаемся
к первому уравнения - до тех пор, пока оба уравнения не окажутся сбалансированы с достаточной точностью.
Также существует метод вычисления частоты m, основанный на сопоставлении фаз. Стефан Бернси использует его
в своем алгоритме транспонирования.
Данный метод быстрее (не требует нескольких итераций для решения системы уравнений), однако требует двух
преобразований Фурье для одного фрагмента с небольшим сдвигом на время ΔT. Этот
алгоритм предполагает, что спектр не меняется за это короткое время - хотя бы в отношении частот.
Допустим, у нас есть фрагмент длиной T. Мы делаем преобразование Фурье не для всего этого
фрагмента, а для двух перекрывающихся фрагментов [0, T - ΔT] и [ΔT, T].
Затем рассматриваем каждую получившуюся гармонику и ее фазу в первом и втором преобразовании. По величине
сдвига фазы вносится поправка q и вычисляется точная частота (n + q)/T.
Наконец, эффект размазывания можно уменьшить, если применять оконные функции.
Хотя в результате применения
оконной функции, спектр может сильно исказиться, но все-таки форма пиков изменяется, сужается их "подошва".
Вот пример для T = 1 сек, N = 44100, m/T = 440.2 Гц, φ = 0.
Слева - спектр с "размазыванием", полученный обычным ДПФ. Справа - с применением оконной функции Хамминга.
Знаете ли Вы, что такое "усталость света"? Усталость света, анг. tired light - это явление потери энергии квантом электромагнитного излучения при прохождении космических расстояний, то же самое, что эффект красного смещения спектра далеких галактик, обнаруженный Эдвином Хабблом в 1926 г. На самом деле кванты света, проходя миллиарды световых лет, отдают свою энергию эфиру, "пустому пространству", так как он является реальной физической средой - носителем электромагнитных колебаний с ненулевой вязкостью или трением, и, следовательно, колебания в этой среде должны затухать с расходом энергии на трение. Трение это чрезвычайно мало, а потому эффект "старения света" или "красное смещение Хаббла" обнаруживается лишь на межгалактических расстояниях. Таким образом, свет далеких звезд не суммируется со светом ближних. Далекие звезды становятся красными, а совсем далекие уходят в радиодиапазон и перестают быть видимыми вообще. Это реально наблюдаемое явление астрономии глубокого космоса. Подробнее читайте в FAQ по эфирной физике.
, 
