Как найти частоту звукового сигнала

Как найти частоту звуковых колебаний

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

Вам понадобится

• — калькулятор;
• — физический справочник;
• — тахометр;
• — звуковой датчик;
• — осциллограф.

Инструкция

Найдите частоту звуковых колебаний, если известна длина их волн и скорость звука в среде, где они распространяются. Вычисления следует производить по формуле F=V/L. Здесь V — скорость звука в среде, а L — длина волны (известная величина). Значения скоростей звука для разных сред можно узнать из физических справочников. Так, для воздуха при нормальных условиях (температуре в районе 20°C и давлении, близком к атмосферному) это значение составляет 341 м/с. Поэтому, например, звуковые колебания в воздухе с длиной волны в 0,25 м будут иметь частоту 341/0,25=1364 Гц.

Найти частоту звуковых колебаний, зная их период, можно по простой формуле: F=1/T. Обратите внимание на то, что для получения корректных значений частоты, представленных в герцах, период T должен быть выражен системе СИ, то есть иметь размерность в секундах.

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

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

Обратите внимание

При вычислениях, связанных с получением значений частоты звука, всегда переводите все известные значения в систему СИ.

Полезный совет

Если дополнительно требуется найти циклическую частоту звуковых колебаний, рассчитайте ее по формуле w=2*PI*F, где F — частота, выраженная в герцах, полученная одним из описанных способов.

Войти на сайт

или

Забыли пароль?
Еще не зарегистрированы?

Звуковая волна – период, длина, частота и скорость распространения

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

Звуковая волна – это механические колебания, которые в результате колебаний молекул вещества распространяются в какой-либо
среде (в газе, жидкости или твёрдом теле) и, достигнув органов слуха человека, воспринимаются им как звук. Источник, создающий
возмущение (колебания воздуха), называется источником звука.
Как уже было сказано, для распространения звука необходима какая-либо упругая среда. Поэтому в вакууме ори, не ори – тебя никто не
услышит, по причине того, что звуковые волны распространяться не смогут, так как там нечему колебаться.., да и слушать там, по большому
счёту, тоже некому.

Так же, как и в случае с электромагнитными волнами, соотношение, связывающее длину звуковой волны с частотой колебаний,
в общем случае выглядит следующим образом:
λ (м) = V (м/сек) / F (Гц), где V (м/сек) — это скорость распространения
звука в среде.

Период колебаний также не претерпел никаких изменений и по-прежнему равен:

T(сек) = 1 / F (Гц) = λ (м) / V (м/сек).

Частота колебаний звукового сигнала F (Гц) – это параметр стабильный, практически не зависящий от среды распространения.

А вот скорость звука V (м/сек), а соответственно и длина звуковой волны – это величины, которые зависят
не только от плотности вещества, но и от его упругости, а в случае с жидкостями и газами ещё – и от температуры, и атмосферного
давления.

Зависимость скорости звуковой волны от свойств упругой среды легко прослеживается по следующей формуле:
V (м/сек) = √Eупр (паскаль) / ρ (кг/м³)
,
где Eупр представляет собой модуль объёмной упругости среды, а ρ – плотность среды.
Модуль упругости, так же как и плотность – это справочные величины, прописанные для конкретных материалов.

В качестве примера, ниже приведена таблица величины скорости распространения звука в различных средах:

Среда	Скорость звука, м/сек
Воздух при 0°	331
Воздух при 30°	350
Вода	1450
Медь	3800
Дерево	4800
Железо	4900
Сталь	5600

Для газов параметры модуля объёмной упругости и плотности имеют ярко выраженную зависимость от температуры и атмосферного давления.
Если углубиться, то скорость звука в газах можно вычислить по следующей формуле:

V (м/сек) = √γ*Ратм / ρ ,
где

γ = cp/сv – это отношение удельной теплоёмкости при постоянном давлении
к удельной теплоёмкости при постоянном объёме, а Pатм – атмосферное давление,
которое связано с температурой газообразной среды.

Поэтому, чтобы никого сильно не грузить, приведу и приближённую зависимость скорости звука (при нормальном
атмосферном давлении) от температуры среды:
V (м/сек) = (331 + 0,6 * T°), где 331 м/сек – это скорость звука при 0°С,
а T° – температура в градусах Цельсия.

Теперь можно совместить формулы и получить простое соотношение, связывающее длину звуковой волны с частотой колебаний с учётом
температуры среды:

λ (м) = (331 + 0,6 * T°) / F (Гц).

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

Онлайн калькулятор расчёта длины звуковой волны по частоте

Частота звуковых колебаний f
Температура Т(°С) (по умолчанию 20°)
Длина волны

Онлайн калькулятор расчёта частоты по длине звуковой волны

Длина волны λ при заданной Т
Температура Т(°С) (по умолчанию 20°)
Частота колебаний

Полный диапазон звуковых частот условно находится в пределах:
16…20 000 Гц.
Ниже ( 0,001…16Гц ) – инфразвук.
Выше ( 20…100кГц ) – низкочастотный ультразвук,
ещё выше (100кГц…1МГц) – высокочастотный ультразвук.

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

Частота (Гц)
Октава	Нота
До	До — диез	Ре	Ми — бемоль	Ми	Фа	Фа — диез	Си	Си- диез	Ля	Соль-бемоль	Соль
C	C#	D	Eb	E	F	F#	G	G#	A	Bb	B
0	16.35	17.32	18.35	19.45	20.60	21.83	23.12	24.50	25.96	27.50	29.14	30.87
1	32.70	34.65	36.71	38.89	41.20	43.65	46.25	49.00	51.91	55.00	58.27	61.74
2	65.41	69.30	73.42	77.78	82.41	87.31	92.50	98.00	103.8	110.0	116.5	123.5
3	130.8	138.6	146.8	155.6	164.8	174.6	185.0	196.0	207.7	220.0	233.1	246.9
4	261.6	277.2	293.7	311.1	329.6	349.2	370.0	392.0	415.3	440.0	466.2	493.9
5	523.3	554.4	587.3	622.3	659.3	698.5	740.0	784.0	830.6	880.0	932.3	987.8
6	1047	1109	1175	1245	1319	1397	1480	1568	1661	1760	1865	1976
7	2093	2217	2349	2489	2637	2794	2960	3136	3322	3520	3729	3951
8	4186	4435	4699	4978	5274	5588	5920	6272	6645	7040	7459	7902

Преобразование Фурье в действии: точное определение частоты сигнала и выделение нот

последняя редакция статьи доступна на сайте makeloft.xyz

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

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

Наглядная иллюстрация нотного рисунка

Определение частоты (режим гитарного тюнера)

Обратная задача состоит в том, чтобы разобрать звучащую музыкальную композицию на ноты. То есть разложить суммарный акустический сигнал, улавливаемый ухом, на исходные синусоиды. По сути, этот процесс и представляет собой прямое преобразование Фурье. А нажатие на клавиши и извлечение звука есть процесс обратного преобразования Фурье.

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

Ортогональность, в некотором роде, обозначает несмешиваемость функций. Например, если мы возьмём несколько кусочков цветного пластилина и склеим их, то потом всё же сможем разобрать, какие цвета были изначально, но если хорошенько перемешаем несколько баночек гуашевых красок, то точно восстановить исходные цвета без дополнительной информации уже будет невозможно.

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

Теперь опишем весь процесс анализа более подробно:

1. Всё начинается с того, что звуковые волны колеблют мембрану микрофона, который преобразует их в аналоговые колебания электрического тока.

2. Затем происходит дискретизация аналогового электрического сигнала в цифровую форму. На этом моменте стоит остановиться подробно.

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

Как правило, значения-отсчёты берутся через небольшие равные временные промежутки, то есть с определённой частотой, например, 16000 или 22000 Гц. Однако в общем случае дискретные отсчёты могут идти и неравномерно, но это усложняет математический аппарат анализа, поэтому на практике обычно не применяется.

Существует важная теорема Котельникова-Найквиста-Шеннона, которая гласит, что аналоговый периодический сигнал, имеющий конечный (ограниченный по ширине) спектр, может быть однозначно восстановлен без искажений и потерь по своим отсчётам, взятым с частотой, большей или равной удвоенной верхней частоте спектра (называемой частотой дискретизации или Найквиста).

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

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

3. На следующем этапе происходит само дискретное прямое преобразование Фурье.

Мы выделяем короткий кадр (интервал) композиции, состоящий из дискретных отсчётов, который условно считаем периодическим и применяем к нему преобразование Фурье. В результате преобразования получаем массив комплексных чисел, содержащий информацию об амплитудном и фазовом спектрах анализируемого кадра. Причём спектры также являются дискретными с шагом равным (частота дискретизации)/(количество отсчётов). То есть чем больше мы берём отсчётов, тем более точное разрешение получаем по частоте. Однако при постоянной частоте дискретизации увеличивая число отсчётов, мы увеличиваем анализируемый временной интервал, а поскольку в реальных музыкальных произведениях ноты имеют различную длительность звучания и могут быстро сменять друг друга, происходит их наложение, поэтому амплитуда длительных нот «затмевает» собой амплитуду коротких. С другой стороны для гитарных тюнеров такой способ увеличения разрешения по частоте подходит хорошо, поскольку нота, как правило, звучит долго и одна.

Существует также довольно простой трюк для увеличения разрешения по частоте — нужно исходный дискретный сигнал заполнить нулями между отсчётами. Однако в результате такого заполнения сильно искажается фазовый спектр, но зато увеличивается разрешение амплитудного. Также возможно применение фильтров Фарроу и искусственное увеличение частоты дискретизации, однако и оно вносит искажения в спектры.

Длительность кадра обычно составляет приблизительно от 30 мс до 1 с. Чем он короче, тем лучшее разрешение мы получаем по времени, но худшее по частоте, чем сэмпл длиннее, тем лучшее по частоте, но худшее по времени. Это очень напоминает принцип неопределённости Гейзенберга из квантовой механики..и не с проста, как гласит Википедия, соотношение неопределенностей в квантовой механике в математическом смысле есть прямое следствие свойств преобразования Фурье…

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

Синусоидальный сигнал, ограниченный прямоугольным окном, и его «дифракция»

Дифракция световых волн

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

Применяется оконная функция ко входному кадру очень просто:

for (var i = 0; i < frameSize; i++)
{
    frame[i] *= Window.Gausse(i, frameSize);
}

using System;
using System.Numerics;

namespace Rainbow
{
    public class Window
    {
        private const double Q = 0.5;

        public static double Rectangle(double n, double frameSize)
        {
            return 1;
        }

        public static double Gausse(double n, double frameSize)
        {
            var a = (frameSize - 1)/2;
            var t = (n - a)/(Q*a);
            t = t*t;
            return Math.Exp(-t/2);
        }

        public static double Hamming(double n, double frameSize)
        {
            return 0.54 - 0.46*Math.Cos((2*Math.PI*n)/(frameSize - 1));
        }

        public static double Hann(double n, double frameSize)
        {
            return 0.5*(1 - Math.Cos((2*Math.PI*n)/(frameSize - 1)));
        }

        public static double BlackmannHarris(double n, double frameSize)
        {
            return 0.35875 - (0.48829*Math.Cos((2*Math.PI*n)/(frameSize - 1))) +
                   (0.14128*Math.Cos((4*Math.PI*n)/(frameSize - 1))) - (0.01168*Math.Cos((4*Math.PI*n)/(frameSize - 1)));
        }
    }
}

Что касается компьютеров, в своё время был разработан алгоритм быстрого преобразования Фурье, который минимизирует число математических операций, необходимых для его вычисления. Единственное требование алгоритма состоит в том, чтобы число отсчётов было кратно степени двойки (256, 512, 1024 и так далее).

Ниже его классическая рекурсивная реализация на языке C#.

using System;
using System.Numerics;

namespace Rainbow
{
    public static class Butterfly
    {
        public const double SinglePi = Math.PI;
        public const double DoublePi = 2*Math.PI;

        public static Complex[] DecimationInTime(Complex[] frame, bool direct)
        {
            if (frame.Length == 1) return frame;
            var frameHalfSize = frame.Length >> 1; // frame.Length/2
            var frameFullSize = frame.Length;

            var frameOdd = new Complex[frameHalfSize];
            var frameEven = new Complex[frameHalfSize];
            for (var i = 0; i < frameHalfSize; i++)
            {
                var j = i << 1; // i = 2*j;
                frameOdd[i] = frame[j + 1];
                frameEven[i] = frame[j];
            }

            var spectrumOdd = DecimationInTime(frameOdd, direct);
            var spectrumEven = DecimationInTime(frameEven, direct);

            var arg = direct ? -DoublePi/frameFullSize : DoublePi/frameFullSize;
            var omegaPowBase = new Complex(Math.Cos(arg), Math.Sin(arg));
            var omega = Complex.One;
            var spectrum = new Complex[frameFullSize];

            for (var j = 0; j < frameHalfSize; j++)
            {
                spectrum[j] = spectrumEven[j] + omega*spectrumOdd[j];
                spectrum[j + frameHalfSize] = spectrumEven[j] - omega*spectrumOdd[j];
                omega *= omegaPowBase;
            }

            return spectrum;
        }

        public static Complex[] DecimationInFrequency(Complex[] frame, bool direct)
        {
            if (frame.Length == 1) return frame;
            var halfSampleSize = frame.Length >> 1; // frame.Length/2
            var fullSampleSize = frame.Length;

            var arg = direct ? -DoublePi/fullSampleSize : DoublePi/fullSampleSize;
            var omegaPowBase = new Complex(Math.Cos(arg), Math.Sin(arg));
            var omega = Complex.One;
            var spectrum = new Complex[fullSampleSize];

            for (var j = 0; j < halfSampleSize; j++)
            {
                spectrum[j] = frame[j] + frame[j + halfSampleSize];
                spectrum[j + halfSampleSize] = omega*(frame[j] - frame[j + halfSampleSize]);
                omega *= omegaPowBase;
            }

            var yTop = new Complex[halfSampleSize];
            var yBottom = new Complex[halfSampleSize];
            for (var i = 0; i < halfSampleSize; i++)
            {
                yTop[i] = spectrum[i];
                yBottom[i] = spectrum[i + halfSampleSize];
            }

            yTop = DecimationInFrequency(yTop, direct);
            yBottom = DecimationInFrequency(yBottom, direct);
            for (var i = 0; i < halfSampleSize; i++)
            {
                var j = i << 1; // i = 2*j;
                spectrum[j] = yTop[i];
                spectrum[j + 1] = yBottom[i];
            }

            return spectrum;
        }
    }
}

Существует две разновидности алгоритма БПФ — с прореживанием по времени и по частоте, но оба дают идентичный результат. Функции принимают массив комплексных чисел, заполненный реальными значениями амплитуд сигнала во временной области, а после своего выполнения возвращают массив комплексных чисел, содержащий информацию об амплитудном и фазовом спектрах. Стоит помнить, что реальная и мнимая части комплексного числа — это далеко не то же самое, что его амплитуда и фаза!

magnitude = Math.Sqrt(x.Real*x.Real + x.Imaginary*x.Imaginary)
phase = Math.Atan2(x.Imaginary, x.Real)

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

Также существует разновидность алгоритма БПФ без рекурсии по Кули-Тьюки, которая часто применяется на практике, но она чуть более сложна для восприятия.

Сразу после вычисления преобразования Фурье удобно нормализовать амплитудный спектр:

var spectrum = Butterfly.DecimationInTime(frame, true);
for (var i = 0; i < frameSize; i++)  
{
	spectrum[i] /= frameSize;
}

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

Вычислив амплитудный и частотный спектры, легко производить обработку сигнала, например, применять частотную фильтрацию или производить сжатие. По сути, таким образом можно сделать эквалайзер: выполнив прямое преобразование Фурье, легко увеличить или уменьшить амплитуду определённой области частот, после чего выполнить обратное преобразование Фурье (хотя работа настоящих эквалайзеров обычно основана на другом принципе — фазовом сдвиге сигнала). Да и сжать сигнал очень просто — нужно всего лишь сделать словарь, где ключом является частота, а значением соответствующее комплексное число. В словарь нужно занести лишь те частоты, амплитуда сигнала на которых превышает какой-то минимальный порог. Информация о «тихих» частотах, не слышимых ухом, будет потеряна, но получится ощутимое сжатие при сохранении приемлемого качества звучания. Отчасти этот принцип лежит в основе многих кодеков.

4. Точное определение частоты

Дискретное преобразование Фурье даёт нам дискретный спектр, где каждое значение амплитуды отстоит от соседних на равные промежутки по частоте. И если частота в сигнале кратна шагу равному (частота дискретизации)/(количество отсчётов), то мы получим выраженный остроконечный пик, но если частота сигнала лежит где-то между границами шага ближе к середине у нас выйдет пик со «срезанной» вершиной и нам будет затруднительно сказать, что же там за частота. Очень может быть что в сигнале присутствуют две частоты лежащие рядом друг с другом. В этом и заключается ограничение разрешения по частоте. Так же как на фотоснимке с низким разрешением мелкие предметы склеиваются и становятся неразличимы, так же и тонкие детали спектра могут теряться.

Но частоты музыкальных нот лежат далеко не на сетке шагов преобразования Фурье, а для повседневных задач настройки музыкальных инструментов и распознавания нот необходимо знать именно точную частоту. Более того, на низких октавах при разрешении от 1024 отсчётов и ниже сетка частот Фурье становится настолько редкой, что попросту на одном шаге начинают умещаться несколько нот и определить какая же на самом деле из них играет становится фактически невозможно.

Чтобы как-то обойти это ограничение иногда применяют аппроксимирующие функции, например, параболические.
www.ingelec.uns.edu.ar/pds2803/Materiales/Articulos/AnalisisFrecuencial/04205098.pdf
mgasior.web.cern.ch/mgasior/pap/biw2004_poster.pdf
Но всё это искусственные меры, которые улучшая одни показатели могут давать искажения в других.

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

Подробно прочитать о нём можно по ссылкам:
www.guitarpitchshifter.com/algorithm.html
www.dspdimension.com/admin/pitch-shifting-using-the-ft (+ примеры кода)
eudl.eu/pdf/10.1007/978-3-642-29157-9_43
ctuner.googlecode.com (примеры использования алгоритма на C++ и Java)

На C# реализация метода выглядит довольно просто:

using System;
using System.Collections.Generic;
using System.Linq;
using System.Numerics;

namespace Rainbow
{
    //  Δ∂ωπ
    public static class Filters
    {
        public const double SinglePi = Math.PI;
        public const double DoublePi = 2*Math.PI;

        public static Dictionary<double, double> GetJoinedSpectrum(
            IList<Complex> spectrum0, IList<Complex> spectrum1,
            double shiftsPerFrame, double sampleRate)
        {
            var frameSize = spectrum0.Count;
            var frameTime = frameSize/sampleRate;
            var shiftTime = frameTime/shiftsPerFrame;
            var binToFrequancy = sampleRate/frameSize;
            var dictionary = new Dictionary<double, double>();

            for (var bin = 0; bin < frameSize; bin++)
            {
                var omegaExpected = DoublePi*(bin*binToFrequancy); // ω=2πf
                var omegaActual = (spectrum1[bin].Phase - spectrum0[bin].Phase)/shiftTime; // ω=∂φ/∂t
                var omegaDelta = Align(omegaActual - omegaExpected, DoublePi); // Δω=(∂ω + π)%2π - π
                var binDelta = omegaDelta/(DoublePi*binToFrequancy);
                var frequancyActual = (bin + binDelta)*binToFrequancy;
                var magnitude = spectrum1[bin].Magnitude + spectrum0[bin].Magnitude;
                dictionary.Add(frequancyActual, magnitude*(0.5 + Math.Abs(binDelta)));
            }

            return dictionary;
        }

        public static double Align(double angle, double period)
        {
            var qpd = (int) (angle/period);
            if (qpd >= 0) qpd += qpd & 1;
            else qpd -= qpd & 1;
            angle -= period*qpd;
            return angle;
        }
	}
}

Применение также несложное:

                        var spectrum0 = Butterfly.DecimationInTime(frame0, true);
                        var spectrum1 = Butterfly.DecimationInTime(frame1, true);
                        for (var i = 0; i < frameSize; i++)
                        {
                            spectrum0[i] /= frameSize;
                            spectrum1[i] /= frameSize;
                        }

                        var spectrum = Filters.GetJoinedSpectrum(spectrum0, spectrum1, ShiftsPerFrame, Device.SampleRate);

Обычно исходные кадры сдвинуты на 1/16 или 1/32 своей длины, то есть ShiftsPerFrame равно 16 или 32.

В результате мы получим словарь частота-амплитуда, где значения частот будут довольно близки к реальным. Однако «срезанные пики» всё ещё будут наблюдаться, хоть и менее выражено. Чтобы устранить этот недостаток, можно просто «дорисовать» их.

using System;
using System.Collections.Generic;
using System.Linq;
using System.Numerics;

namespace Rainbow
{
    public static class Filters
    {
        public static Dictionary<double, double> Antialiasing(Dictionary<double, double> spectrum)
        {
            var result = new Dictionary<double, double>();
            var data = spectrum.ToList();
            for (var j = 0; j < spectrum.Count - 4; j++)
            {
                var i = j;
                var x0 = data[i].Key;
                var x1 = data[i + 1].Key;
                var y0 = data[i].Value;
                var y1 = data[i + 1].Value;

                var a = (y1 - y0)/(x1 - x0);
                var b = y0 - a*x0;

                i += 2;
                var u0 = data[i].Key;
                var u1 = data[i + 1].Key;
                var v0 = data[i].Value;
                var v1 = data[i + 1].Value;

                var c = (v1 - v0)/(u1 - u0);
                var d = v0 - c*u0;

                var x = (d - b)/(a - c);
                var y = (a*d - b*c)/(a - c);

                if (y > y0 && y > y1 && y > v0 && y > v1 &&
                    x > x0 && x > x1 && x < u0 && x < u1)
                {
                    result.Add(x1, y1);
                    result.Add(x, y);
                }
                else
                {
                    result.Add(x1, y1);
                }
            }

            return result;
        }
    }
}

Перспективы

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

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

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

Итоги

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

P.S. Радужный фреймворк (Rainbow Framework) содержащий все вышеизложенные примеры кода можно скачать с Кодплекса.

P.P.S. Возможно, эта статья когда-нибудь серьёзно пригодится вам в профессиональной деятельности и поможет сохранить уйму времени и сил, поэтому если вы желаете отблагодарить автора за труд, то можете сделать пожертвование, купить приложение (бесплатная версия с рекламой), благодаря которому статья возникла, либо просто выразить благодарность добрым словом.

Литература

1. Основы спектрального анализа звуков pandia.org/text/77/481/644.php

2. Алгоритм Кули-Тьюки www.codeproject.com/Articles/32172/FFT-Guitar-Tuner

3. Передискретизация (ресэмплинг)
www.dsplib.ru/forum/viewtopic.php?f=5&t=11
www.dsplib.ru/content/farrow/farrow.html

4. Коррекция частоты по смещению фазы
www.guitarpitchshifter.com/algorithm.html
www.dspdimension.com/admin/pitch-shifting-using-the-ft (+ примеры кода)
eudl.eu/pdf/10.1007/978-3-642-29157-9_43
ctuner.googlecode.com (примеры использования алгоритма на C++ и Java)

    Очень часто походя употребляют такие вроде бы понятные термины, как спектр, фаза, частота и прочие. Но зачастую мы до конца не понимаем, что же это на самом деле такое.  Что значат эти термины на самом деле, как можно «пощупать» их истинное значение? Можно пойти в библиотеку и почитать там книги по теории радиотехники и цифровой обработке сигналов, но времени постоянно не хватает даже на более важные дела. Поэтому автор попытался дать читателю выжимки из радиотехнических учебников, объясненные «на пальцах» и самый минимум формул (если кто-то заинтересовался более «математическим» изложением материала).

Волновая форма сигнала (звука). Период. Частота

     Что такое звук? Это переменное звуковое (воздушное) давление на барабанную перепонку. Ухо воспринимает как звук только изменение давления. Когда звучит отдельная нота давление периодически то нарастает, то убывает и этот процесс циклически повторяется.

Период (T, сек) — длительность этого цикла. 

Частота (f, Гц, Герц) — количество периодов, помещающихся в одной секунде. 1 Герц — это 1 период за секунду.

f = 1 / T (формула частоты)

    Причем закон (форма) изменения звукового давления не изменяется от периода к периоду.

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

    А что же такое шум?

    Шум — это сигнал (волновая форма не имеет периода), который в любой момент времени имеет случайное значение звукового давления. Шум не имеет периода.

   Звук, как известно распространяется с задержкой, которая зависит от расстояния от источника до человеческого уха. Как это происходит?

Длина волны

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

     Расстояние, между областями одинакового сжатия воздуха называется длиной звуковой волны.

L = M / f (формула длины волны),

где

L — длина волны в метрах;

M — скорость звука (331,46 м/с) в метрах в секунду;

f — частота звука в Герцах.

   Длина волны для:

    20 Гц L20 = (331,46 м/с) / (20 Гц) = 16,5 м.

    100 Гц L100 = (331,46 м/с) / (100 Гц) = 3,3 м.

    1000 Гц L1000 = (331,46 м/с) / (1000 Гц) = 0,33 м = 33 см.

    10000 Гц L10000 = (331,46 м/с) / (10000 Гц) = 0,033 м = 3,3 см.

    20000 Гц L10000 = (331,46 м/с) / (20000 Гц) = 0,017 м = 1,7 см.

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

     Расстояние  вносит  задержку распространения звука не зависящую от частоты, так как скорость звука на разных частотах одинакова.

   Dt = l / M (формула задержки распространения звука),

где

Dt — задержка в секундах;

l — расстояние в метрах;

M — скорость звука (331,46 м/с) в метрах в секунду.

      1 метр вносит задержку распространения звука

      Dt= (1 м) / (331,46 м/с) = 0,003 секунды или 3 миллисекунды (мс).

Автор: Юрий Корзунов (2010)

ПРОДОЛЖЕНИЕ…>>