Arduino:Библиотеки/FFT

Содержание

Знакомство с Arduino

Продукты

Основы

Справочник языка Arduino

Примеры

Библиотеки

Хакинг

Изменения

Сравнение языков Arduino и Processing

Перевод: Максим Кузьмин

Проверка/Оформление/Редактирование: Мякишев Е.А.

Библиотека FFT

Примечания:

• Если вам нужна более быстрая библиотека для анализа частоты, попробуйте FHT
• Выходные данные, получаемые библиотекой FFT, представлены не в ASCII, а в двоичном формате. То есть данные, которые будут приходить по последовательному порту, удобочитаемыми не будут. Чтобы исправить это, поменяйте serial.write() на serial.print(). Кроме того, если вы хотите видеть каждый выходной отсчет, возможно, понадобится цикл for().

О библиотеке

Библиотека FFT – программное воплощение стандартного алгоритма FFT (от «fast Fourier transform», что значит «быстрое преобразование Фурье»; далее – БПФ), который может оперировать только действительными числами. Функционал библиотеки способен выдавать от 16 до 256 выходных отсчетов разрешением от 8 до 16 бит и минимальной скоростью обновления, равной примерно 7 миллисекунд. Кроме того, выдача результата может быть в четырех разных форматах: линейном (16 и 8 бит), логарифмическом (8 бит) и октавном (8 бит). Подробнее о каждом из этих четырех режимов читайте в разделе «Функции». Поскольку исходные данные для алгоритма имеют 16-битный формат (т.е. являются числом с фиксированной точкой), на низких частотах порог шума будет -72 дБ, а на высоких он будет -76 дБ. Если вы будете использовать встроенный АЦП, то его порог шума рассчитывается по тому же принципу, что и порог шума для БПФ, что дает ОСШ (отношение сигнал/шум) где-то между 9 и 10 битами (-55 дБ).

Скорость

Функция	Запуск	Реорганизация	Функция-окно	Линейный формат (16 бит)	Линейный формат (8 бит)*	Логарифмический формат
N	мс	мкс	мкс	мкс	мкс	мкс
256	6,32	412	608	588	470	608
128	2,59	193	304	286	234	290
64	1,02	97	152	145	114	144
32	0,37	41	76	80	59	74
16	0,12	21	37	46	30	39

Примечание

Значения для линейного формата (8 бит) – приблизительные, т.к. могут немного варьироваться из-за директивы SCALE. Более подробно читайте ниже.

Память

Функция	Запуск	Реорганизация	Функция-окно	Линейный формат (16 бит)	Линейный формат (8 бит)	Логарифмический формат
N	S/F (байты)	F (байты)	F (байты)	S/F (байты)	S/F (байты)	S/F (байты)
256	6,32	412	608	588	470	608
128	2,59	193	304	286	234	290
64	1,02	97	152	145	114	144
32	0,37	41	76	80	59	74
16	0,12	21	37	46	30	39

• S – SRAM-память, F – flash-память.

Загрузка

• Версия для Arduino 1.05 (новая). Исправлена PROGMEM для новой версии AVR-GCC (благодаря форумчанину kirill9617), добавлен пример с передачей данных по последовательному порту.
• Версия для Arduino 1.0 или Arduino-0022 (старая). Должна работать с обеими.

Установка

Чтобы использовать библиотеку, скачайте ZIP-файл, распакуйте его, а затем поместите в папку для библиотек IDE Arduino.

Чтобы узнать, где находится эта папка на PC, откройте IDE Arduino, кликните на Файл > Настройки (File > Preferences) – адрес папки будет указан в самом верху, в поле «Размещение папки скетчей» (Sketchbook location). Затем перезапустите IDE Arduino. Если папки «libraries» нет, ее нужно создать.

На Mac эта папка находится по адресу «/Пользователи/<Имя пользователя>/Документы/Arduino», или «/Пользователи/<Имя пользователя>/Документы/Maple» (если используете Maple). Если папки «libraries» нет, ее нужно создать.

Функции

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

Всего в библиотеке семь функций: две предварительные fft_reorder() и fft_window(), одна главная fft_run() и четыре амплитудные, отвечающие за разные форматы выходных данных: для 16-битного линейного, 8-битного линейного, 8-битного логарифмического и 8-битного октавного.

fft_run()

Функция fft_run() – главная функция библиотеки FFT. Ей не нужны никакие аргументы, и она ничего не возвращает. Если вы запускаете эту функцию, то предполагается, что нужные данные уже реорганизованы и находятся в SRAM-памяти. Эти данные хранятся в архиве под названием fft_input[], который содержит два 16-битных значения – действительное и комплексное. Если вы заполняете этот массив самостоятельно, то действительные значения сохраняйте в четные выходные отсчеты, а комплексные – нечетные.

Таким образом, fft_input[0] – это первое действительное значение, fft_input[1] – это первое комплексное значение, fft_input[2] – это второе действительное значение, fft_input[3] – это второе комплексное значение и т.д.

Следовательно, получается удвоение, поэтому количество отсчетов в массиве будет в два раза больше, чем, собственно, реальных ПБФ-отсчетов. Если вы используете только действительные числа (т.е. значения, полученные от АЦП), то сохраняйте их в четные отсчеты, а в нечетные отсчеты вписывайте «0».

Результат сохраняется в fft_input[], где четные отсчеты отвечают за действительные амплитуды, а нечетные – за комплексные. Сами отсчеты являются диапазонами частоты (по возрастающей).

Таким образом, fft_input[0] и fft_input[1] – это амплитуда для первого отсчета (0hz -> Fs/N), fft_input[2] и fft_input[3] – это амплитуда для второго отсчета (Fs/N -> 2Fs/N) и т.д.

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

fft_reorder()

Функция fft_reorder() реорганизовывает входные данные, подготавливая их к обработке алгоритмом БПФ. Возможно, это функция не понадобится, т.к. эта реорганизация может быть выполнена вручную. Функция fft_reorder() запускается перед fft_run(). Ей не нужны никакие аргументы, и она ничего не возвращает. Кроме того, перед вызовом этой функции нужно заполнить массив fft_input[], т.к. она берет исходные данные именно оттуда.

fft_window()

Функция fft_window() умножает входные данные при помощи функции-окна, чтобы увеличить частотное разрешение данных для ПБФ. Ей не нужны никакие аргументы, и она ничего не возвращает. Кроме того, перед вызовом этой функции нужно заполнить массив fft_input[], т.к. она берет исходные данные именно оттуда. Функция fft_reorder() запускается перед fft_run() и fft_reorder().

fft_mag_lin8()

Функция fft_mag_lin8() дает амплитуду для каждого отсчета БПФ. Она суммирует действительные и комплексные значения, возведенные в квадрат, а потом извлекает квадратный корень, округляя результат до 8-битного значения (она использует таблицу поиска, подгоняя эти значения к полному 8-битному диапазону). Ей не нужны никакие аргументы, и она ничего не возвращает. Исходными данными для нее служат данные из массива fft_input[], а результат сохраняется в массив fft_lin_out8[], после чего эти данные можно использовать для других целей. Амплитуда рассчитывается только для первой половины выходных отсчетов (N/2), поскольку вторая половина БПФ идентична первой для всего диапазона исходных значений. Следовательно, количество 8-битных значений в массиве fft_lin_out8[] будет составлять N/2, где каждый индекс будет соответствовать номеру выходного отсчета (минус «1»).

То есть fft_lin_out8[0] – это амплитуда для первого отсчета (0hz -> Fs/N), fft_lin_out8[1] – это амплитуда для второго отсчета (Fs/N -> 2Fs/N) и т.д.

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

fft_mag_lin()

Функция fft_mag_lin() дает амплитуду для каждого отсчета БПФ. Она суммирует действительные и комплексные значения, возведенные в квадрат, а потом извлекает квадратный корень. Для извлечения квадратного корня используется таблица поиска, поэтому точность имеет некоторые ограничения. Функция охватывает полный 16-битный диапазон, но в любой отдельной точке этого диапазона разрешение все же будет 8-битным.

Это 8-битное значение, где 4 бита отведено под порядок (т.е. под экспоненту). Исходными данными для нее служат данные из массива fft_input[], а результат сохраняется в массив fft_lin_out[]. Итоговые данные расположены в последовательном порядке, а их общее количество составляет N/2, поскольку вторая половина БПФ идентична первой для всего диапазона исходных значений.

fft_mag_log()

Функция fft_mag_log() дает амплитуду для каждого отсчета БПФ. Она суммирует действительные и комплексные значения, возведенные в квадрат, затем извлекает квадратный корень. Далее функция рассчитывает логарифм полученного числа по отношению к двум. Результат будет представлен в логарифмическом формате – фактически, в децибелах. Функции fft_mag_log() не нужны никакие аргументы, и она ничего не возвращает.

Для извлечения квадратного корня используется таблица поиска, а результат подгоняется под полный 8-битный диапазон. То есть уравнение функции будет следующим: 16*(log2((комплексн2 + действ2)1/2)). Исходными данными для функции служат данные из массива fft_input[], а результат сохраняется в массив fft_log_out[]. Итоговые значения имеют тот же порядок, что и выходные отсчеты БПФ, а общее количество отсчетов составляет N/2, поскольку вторая половина БПФ идентична первой для всего диапазона исходных значений.

fft_mag_octave()

Функция fft_mag_octave() дает среднее квадратичное значение всех выходных отсчетов, представляя это все в октавном (удвоение частоты) формате. Этот формат, как правило, более полезен, потому что близок к тому, как люди воспринимают звук. Ей не нужны никакие аргументы, и она ничего не возвращает. Исходными данными для функции служат данные из массива fft_input[], а результат сохраняется в массив fft_oct_out[]. Итоговые данные представляют собой 8-битное значение, рассчитанное по формуле 16*log2(квадр(ампл)). Количество выходных отсчетов рассчитывается следующим образом:

FFT_N = 256 : отсчеты = [0, 1, 2:4, 5:8, 9:16, 17:32, 3:64, 65:128]
FFT_N = 128 : отсчеты = [0, 1, 2:4, 5:8, 9:16, 17:32, 3:64]

Здесь, к примеру, 5:8 – это сумма всех отсчетов, расположенных с 5-го по 8-ой. Данные в каждом отсчете возводятся в квадрат (и действительные, и комплексные значения), а затем суммируются со всеми амплитудами (тоже возведенными в квадрат), находящимися в этом диапазоне. Затем берется количество выходных отсчетов и делится на полученное число (это, впрочем, можно отключить – смотрите следующий раздел) после чего извлекается квадратный корень и рассчитывается логарифм.

Подробнее о принципе работы библиотеки

Хотите «заглянуть под капот»? Окей, попробую объяснить, по каким принципам это все работает. За ускорение вычислений в библиотеке FFT отвечают две вещи.

Во-первых, в каждом БПФ вам нужно умножать множество исходных значений на синусные и косинусные константы. На ATmega это требует больших вычислительных ресурсов, поскольку умножение 16 бит на 16 бит требует 18 тактовых циклов. С другой стороны, сложение 16 бит и 16 бит требует только 2 тактовых цикла. Отсюда вывод, что лучше складывать, чем умножать. Дело в том, что синусные и косинусные константы, используемые в БПФ – это просто «0» и «1», поэтому вам не обязательно использовать операцию умножения, т.к. достаточно и сложения. К примеру, если БПФ имеет 256 выходных отсчетов, то в нем нужно выполнить 1024 сложных умножения, из которых 382 – это «0» или «1». Почти половина!

Библиотека FFT для Arduino ищет условия, где используются эти «0» и «1», и попросту делает там сложение. Дело в том, что эти константы появляются с регулярными интервалами, поэтому искать их довольно просто. Впрочем, с увеличением количества выходных отсчетов преимущества этого метода снижаются. Если за количество выходных отсчетов в БПФ представить N, то экономия времени будет составлять (1.5*N - 2), а общим количеством умножений будет (N/2)*log2(N). Таким образом, формулой, по которой рассчитываются экономия, будет 3/log2(N), и чем выше N, тем меньше времени можно сэкономить.

Второй метод, направленный на ускорение работы БПФ – это таблицы поиска, с помощью которых рассчитываются квадратные корни амплитуд. Сложность этого метода в том, что исходных данных гораздо больше, чем содержимое таблицы поиска. Таким образом, чтобы не тратить зря программную память, исходные данные нужно сжать. К примеру, если извлечь квадратный корень из 16-битного значения, мы получим 64 тысячи входных значений, которые нужно как-то подогнать к 256 выходным значениям. В программной памяти Arduino невозможно провести такой объем вычислений, поэтому в библиотеке ко всему массиву входных данных применяется линейная интерполяция, в результате которой для разных участков будут рассчитаны разные линии. К примеру, для линейного 8-битного формата будет достаточно 3-4 линейных секций, причем без потери разрешения.

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

Версия с 32-битными входными данными немного сложнее, потому что теперь выходные данные нужно подогнать к 16 битам, а линейная интерполяция здесь не поможет, потому что не умеет сжимать значения, разрешение которых больше 16 бит. Поэтому здесь применяется гибридный подход, при котором входное значение конвертируется в 16-битное значение с 8-битным порядком. Это можно сделать очень быстро с основанием «2», что позволяет использовать таблицу поиска квадратного корня и для значений с разрешением выше 16 бит. Если при сжатии входного значения выполнен сдвиг на 2 бита вправо, то выходное значение может быть реконструировано со сдвигом на 1 бит влево. Как правило, после создания порядок (т.е. экспоненту) округляют до четного числа, в результате чего квадратный корень может вернуть целочисленное значение.

Разумеется, 32-битная версия не так точна, как библиотека для вычисления квадратного корня, но требует всего 40 циклов, тогда как эта библиотека требует 500. Таблица поиска возвращает значение, где точными являются только первые 8 бит, но для этого БПФ такой точности вполне достаточно. Общая битовая глубина БПФ не сильно уступает 12 битам, поскольку в нем используется фиксированная точка (во избежание переполнения перед сложением каждое значение должно быть поделено на 2; если БПФ имеет 256 выходных отсчетов, это в конце концов дает деление на 256). Точность зависит от разрешения выходных данных. Если это 8-битное значение, то точность будет максимальной. Если это 9-битное значение, то один младший бит, возможно, будет неправильным. Если это 10-битное значение, то два младших бита, возможно, будут неправильными и т.д. То есть самый худший сценарий – это 16-битное значение с точностью +/- 0,5%.

Директивы

Эти значения позволяют вам модифицировать код FFT под собственные нужды. Большинство из них включают и выключают определенные функции библиотеки. По умолчанию большинство функций выключены, поэтому чтобы использовать, их нужно сначала включить. Эти директивы вписываются при помощи конструкции #define и делать это нужно перед подключением библиотек, что делается конструкцией #include.

• Директива FFT_N задает размер БПФ. Возможные варианты: 16, 32, 64, 128 и 256. По умолчанию выставлен вариант 256.
• Директива SCALE отвечает за масштабирование в функции fft_mag_lin8(). Поскольку 8 бит – это довольно слабое разрешение, пользователь, возможно, захочет масштабировать данные, чтобы максимизировать их до полного диапазона.
• Директива SCALE умножает выходные данные на константу перед извлечением квадратного корня, в результате чего на выходе получается максимизированное разрешение. Это потребует чуть больше ресурсов, чем обычно, но излишний расход будет очень невелик. В директиве SCALE можно указать любое число в диапазоне от 1 до 255. По умолчанию выставлено 1. Наименьшее количество ресурсов расходуется при использовании 1, 2, 4, 128 и 256.
• Директива WINDOW включает и выключает функцию-окно. Если вы не используете fft_window(), то нужно вписать WINDOW 0 (выкл). По умолчанию стоит WINDOW 1 (вкл).
• Директива REORDER включает и выключает функцию реорганизации. Если вы не используете fft_reorder(), то нужно вписать REORDER 0 (выкл). По умолчанию стоит REORDER 1 (вкл).
• Директива LOG_OUT включает и выключает функцию для логарифмического формата. Если вы используете fft_mag_log(), то нужно вписать LOG_OUT 1 (вкл). По умолчанию стоит LOG_OUT 0 (выкл).
• Директива LIN_OUT включает и выключает функцию для линейного (16 бит) формата. Если вы используете fft_mag_lin(), то нужно вписать LIN_OUT 1 (вкл). По умолчанию стоит LIN_OUT 0 (выкл).
• Директива LIN_OUT8 включает и выключает функцию для линейного (8 бит) формата. Если вы используете fft_mag_lin8(), то нужно вписать LIN_OUT8 1 (вкл). По умолчанию стоит LIN_OUT8 0 (выкл).
• Директива OCTAVE включает и выключает функцию для октавного формата. Если вы используете fft_mag_octave(), то нужно вписать OCTAVE 1 (вкл). По умолчанию стоит OCTAVE 0 (выкл).
• Директива OCT_NORM включает и выключает у октавного формата функцию нормализации. Она является частью функции fft_mag_octave(), которая делит каждый сгруппированный отсчет на среднее количество отсчетов. Поскольку многие источники звука издают так называемый «розовый шум» (это падение амплитуды при увеличении частоты), это приводит и к уменьшению масштаба. Будучи выключенной (OCT_NORM 0), эта директива искусственно повышает амплитуду на высоких частотах. По умолчанию нормализация включена (OCT_NORM 1).

Примеры

• fft adc - Применение БПФ на данных от АЦП и отправка результата посредством write()
• fft adc serial - Применение БПФ на данных от АЦП и отправка результата посредством println()
• fft codec - Применение БПФ на данных от модуля Codec Shield