Embedded HandbookПриложения117. Данные с датчика: окна, признаки, CMSIS-DSP
Модель из раздела 116 получала готовые числа — СКЗ вибрации и частоту пика. На приборе их надо посчитать самому из потока отсчётов акселерометра: нарезать поток на окна, убрать смещение, взять спектр, найти пик. От того, как это сделано, зависит больше, чем от выбора модели: неудачный признак портит любой классификатор. Раздел считает признаки библиотекой CMSIS-DSP — той самой, что идёт на F407, — только на компьютере, где её можно сверить с расчётом по формулам, и показывает, как один правильно выбранный признак поднимает распознавание подшипника с 80 до 100 случаев из 100.
Модели — раздел 116, АЦП — раздел 14, DMA — раздел 9, акселерометр — раздел 99.
Как это называется в документации
Оригинал Перевод sampling rate частота выборки Nyquist frequency частота Найквиста: половина частоты выборки aliasing наложение частот anti-aliasing filter фильтр перед АЦП против наложения frame, window окно: кусок сигнала из N отсчётов fast Fourier transform (FFT) быстрое преобразование Фурье frequency bin бин: полоса частот одного отсчёта спектра spectral leakage утечка: энергия тона растекается по соседним бинам window function (Hann) оконная функция (окно Ханна) parabolic interpolation уточнение пика параболой DC offset постоянная составляющая, смещение нуля crest factor пик-фактор: наибольшее значение / СКЗ band energy энергия в полосе частот
117.1 От датчика до признаков
Цепочка на приборе: акселерометр → АЦП по таймеру, отсчёты кладёт DMA → окно из 256 отсчётов → признаки → решение модели. При частоте выборки 1024 Гц окно длится 250 мс, и решение готово четыре раза в секунду; окна с перекрытием наполовину дают решение каждые 125 мс. Всё, что ниже, — о том, как из этих 256 чисел получить два-три признака, которым можно верить.
117.2 Частота выборки и наложение
Выборка различает частоты только до половины своей: при 1024 Гц — до 512 Гц. Всё, что выше, не пропадает, а притворяется низкой частотой. Опыт: тон 160 Гц, снятый с частотой 256 Гц, лёг в спектре на 96 Гц — 256 − 160. Отличить такой ложный пик от настоящего по спектру нельзя, поэтому частоту выборки берут с запасом над самой высокой интересной частотой, а перед АЦП ставят фильтр, который гасит всё выше половины частоты выборки.
117.3 Окно: время против разрешения
Спектр окна из N отсчётов делится на бины шириной FS / N: при 1024 Гц и 256 отсчётах — 4 Гц. Длиннее окно — уже бины, но дольше ждать решения и хуже видны быстрые изменения. Тон 24 Гц ложится ровно на бин 6, и пик получается 24,0000 Гц, а СКЗ — 0,3536 при амплитуде 0,5: ровно 0,5 / √2, потому что в окно поместилось целое число периодов.
117.4 Утечка и окно Ханна
Тон 25 Гц на бин не ложится, и его энергия растекается по соседним бинам — утечка. Без окна в пяти бинах от пика остаётся −24,89 дБ: этого хватает, чтобы спрятать слабый тон подшипника рядом с сильным. Окно Ханна плавно гасит края отрезка, и там же остаётся −52,81 дБ. Окно накладывают перед преобразованием, одним умножением:
if (hann)
arm_mult_f32(buf, win, buf, N_WIN);
117.5 Пик между бинами
Номер бина даёт частоту с шагом 4 Гц. Точнее — вершина параболы через модули трёх соседних бинов:
static float refine(uint32_t k)
{
if (k == 0 || k >= N_WIN / 2 - 1)
return (float)k;
float a = mag[k - 1], b = mag[k], c = mag[k + 1];
float den = a - 2.0f * b + c;
return den != 0.0f ? (float)k + 0.5f * (a - c) / den : (float)k;
}
Для тона 25 Гц без окна парабола даёт 24,19 Гц, с окном Ханна — 24,80 Гц. Окно нужно и здесь: у его пика форма ближе к параболе.
117.6 Смещение, СКЗ и пик-фактор
Акселерометр с питанием от одного источника выдаёт вибрацию вокруг середины шкалы — например, 1,65 В. Эту постоянную составляющую вычитают до всего остального, иначе она войдёт в СКЗ и займёт нулевой бин. Опыт: тон 160 Гц амплитудой 0,5 поверх смещения 1,65 — СКЗ 0,3536, как без смещения, пик ровно на 160,00 Гц. Пик-фактор — наибольшее значение, делённое на СКЗ: у чистого синуса 1,4142, то есть √2, а удары выкрошенного подшипника поднимают его заметно выше — ещё один признак, почти даром.
117.7 Признак под физику дефекта
Три состояния мотора — по 100 окон с разбросом амплитуды ±30 %, частоты и шумом — и порог из раздела 116 на признаках, посчитанных из сигнала:
| Норма | Дисбаланс | Подшипник | |
|---|---|---|---|
| По частоте пика: верно | 100 | 82 | 80 |
| По энергии в полосе 100–248 Гц: верно | 100 | 82 | 100 |
По частоте самого высокого пика подшипник узнан лишь в 80 окнах из 100: когда тон вращения 25 Гц оказывается сильнее тона дефекта, пик уходит на 25 Гц, и в 17 окнах подшипник принят за норму. Средняя частота пика у подшипника — 131,7 Гц вместо 160. Признак «доля энергии в полосе 100–248 Гц больше 0,1» спрашивает не «где главный пик», а «есть ли энергия там, где звучит подшипник», — и узнаёт его во всех 100 окнах. Дисбаланс по-прежнему путается с нормой в 18 окнах: у них одна частота и близкое СКЗ, и это уже вопрос данных и модели — разделы 118–120.
arm_power_f32(mag + 1, N_WIN / 2 - 1, &e_all);
arm_power_f32(mag + BAND_LO, BAND_HI - BAND_LO + 1, &e_band);
f->band = e_all > 0.0f ? e_band / e_all : 0.0f;
117.8 CMSIS-DSP на F407
Весь расчёт признаков — вызовы CMSIS-DSP:
arm_mean_f32(x, N_WIN, &mean); /* смещение нуля датчика и АЦП */
arm_offset_f32(x, -mean, buf, N_WIN);
arm_rms_f32(buf, N_WIN, &f->rms);
arm_absmax_f32(buf, N_WIN, &peak, &k);
Образ для F407 со стартовым файлом Keil занял 8964 байта кода и 2920 байт
таблиц БПФ во флеш-памяти и 6272 байта ОЗУ — буферы окна, спектра
и отсчётов плюс стек и куча стартового файла. Одна строка решает, сколько
займут таблицы: общая arm_rfft_fast_init_f32 принимает любой размер
и потому тянет в образ таблицы всех размеров — 79 504 байта, в 27 раз
больше. Инициализация под один размер оставляет только нужные:
return arm_rfft_fast_init_256_f32(&fft) == ARM_MATH_SUCCESS;
На компьютере та же библиотека собирается gcc в режиме __GNUC_PYTHON__,
и все опыты раздела посчитаны ею; проверка сверила их с расчётом по
формулам.
117.9 Грабли
В спектре пик, которого в механизме нет. Наложение: частота выборки мала или нет фильтра перед АЦП. 160 Гц при 256 Гц выглядят как 96 Гц.
Пик «гуляет» между соседними частотами. Тон не ложится на бин. Окно Ханна и уточнение параболой.
СКЗ не падает до нуля на стоящем моторе. Не вычтена постоянная составляющая датчика.
Прошивка выросла на 75 КБ от одного БПФ. Общая инициализация тянет
таблицы всех размеров; arm_rfft_fast_init_256_f32.
Подшипник то виден, то нет. Признак «частота главного пика»: при сильном тоне вращения пик уходит на него. Энергия в полосе дефекта.
117.10 Итоги
- 1024 Гц и 256 отсчётов: окно 250 мс, бин 4 Гц, различимо до 512 Гц.
- Наложение не лечится после АЦП: 160 Гц при 256 Гц — это 96 Гц.
- Окно Ханна: утечка −52,81 дБ вместо −24,89; пик 25 Гц уточняется до 24,80 Гц.
- Признак выбирают под физику: энергия в полосе подшипника — 100 из 100 против 80 у частоты пика.
- CMSIS-DSP на F407: 2920 байт таблиц при инициализации под размер против 79 504 у общей.
Полные исходники
feat117.c — признаки вибрации из окна отсчётов на CMSIS-DSP.
/* Раздел 117: признаки вибрации из окна отсчётов — на CMSIS-DSP.
* Окно N_WIN отсчётов при FS: постоянная составляющая долой, СКЗ и пик-
* фактор, окно Ханна, быстрое преобразование Фурье, модуль спектра, пик
* и его уточнение параболой между бинами. Тот же файл собирается armclang
* для F407 и gcc на компьютере (CMSIS-DSP в режиме __GNUC_PYTHON__) —
* там feat117_host.c ставит опыты. Проверка — _verify/ml117_check.py. */
#include "arm_math.h"
#define FS 1024.0f /* Гц: частота выборки */
#define N_WIN 256u /* отсчётов в окне: 250 мс, 4 Гц на бин */
#define BAND_LO 25u /* бин 25 — 100 Гц: полоса дефекта подшипника */
#define BAND_HI 62u /* бин 62 — 248 Гц */
typedef struct {
float rms; /* СКЗ без постоянной составляющей */
float peak_hz; /* частота пика, уточнённая параболой */
float crest; /* пик-фактор: наибольший модуль / СКЗ */
float band; /* доля энергии спектра в полосе 100–248 Гц */
uint32_t bin; /* номер бина пика до уточнения */
} features_t;
static float win[N_WIN], buf[N_WIN], spec[N_WIN], mag[N_WIN / 2];
static arm_rfft_fast_instance_f32 fft;
/* инициализация под один размер: общая arm_rfft_fast_init_f32 тянет в образ таблицы всех размеров */
_Static_assert(N_WIN == 256u, "размер окна совпадает с arm_rfft_fast_init_256_f32");
int features_init(void)
{
arm_hanning_f32(win, N_WIN);
#ifdef FFT_ALL_SIZES /* для сравнения размеров образа */
return arm_rfft_fast_init_f32(&fft, N_WIN) == ARM_MATH_SUCCESS;
#else
return arm_rfft_fast_init_256_f32(&fft) == ARM_MATH_SUCCESS;
#endif
}
/* пик между бинами: вершина параболы через три соседних модуля */
static float refine(uint32_t k)
{
if (k == 0 || k >= N_WIN / 2 - 1)
return (float)k;
float a = mag[k - 1], b = mag[k], c = mag[k + 1];
float den = a - 2.0f * b + c;
return den != 0.0f ? (float)k + 0.5f * (a - c) / den : (float)k;
}
void features(const float *x, int hann, features_t *f)
{
float mean, peak, m;
uint32_t k;
arm_mean_f32(x, N_WIN, &mean); /* смещение нуля датчика и АЦП */
arm_offset_f32(x, -mean, buf, N_WIN);
arm_rms_f32(buf, N_WIN, &f->rms);
arm_absmax_f32(buf, N_WIN, &peak, &k);
f->crest = f->rms > 0.0f ? peak / f->rms : 0.0f;
if (hann)
arm_mult_f32(buf, win, buf, N_WIN);
arm_rfft_fast_f32(&fft, buf, spec, 0); /* buf служит рабочим — он больше не нужен */
arm_cmplx_mag_f32(spec, mag, N_WIN / 2);
mag[0] = 0.0f; /* в первой паре — ноль и Найквист вместе: отбросить */
arm_max_f32(mag, N_WIN / 2, &m, &f->bin);
f->peak_hz = refine(f->bin) * FS / N_WIN;
float e_all, e_band; /* энергия — сумма квадратов модулей */
arm_power_f32(mag + 1, N_WIN / 2 - 1, &e_all);
arm_power_f32(mag + BAND_LO, BAND_HI - BAND_LO + 1, &e_band);
f->band = e_all > 0.0f ? e_band / e_all : 0.0f;
}
/* модуль спектра последнего окна — для опытов и рисунка */
const float *features_spectrum(void)
{
return mag;
}
Дальше: раздел 118 «Обучение маленькой модели в контейнере» — от окон с признаками к весам в прошивке.