Использование MS Excel и Python для извлечения вибрационных мод из видео

Аннотация

Сочетание высокоскоростных камер и скриптов в MS Excel и Python предоставляет мощный инструмент для исследования колебательных систем и извлечения их вибрационных характеристик. В этой статье анализируется записанное видео колеблющегося конца стержня, чтобы извлечь траекторию колебаний во времени, а затем на её основе вычислить информацию о вибрационных модах с помощью дискретного преобразования Фурье (DFT). Это достигается за счёт использования доступных встроенных библиотек и инструментов, позволяющих исследователям просто импортировать своё высокоскоростное видео и/или данные трекинга и получать ценную информацию о механическом поведении материалов и конструкций. Здесь мы приводим реальные видеоданные колеблющегося объекта, снятые высокоскоростной камерой Phantom, формируем данные трекинга и показываем, как Excel и Python можно использовать для извлечения вибрационных мод, присутствующих в видео.

Введение

Используя программное обеспечение, реализующее алгоритмы оптического потока, сопоставления признаков или сопоставления по шаблону, пользователь может выбрать объект, который нужно автоматически отслеживать во времени, то есть от кадра к кадру. Эти данные трекинга, состоящие из пространственных координат объекта во времени, служат основой для вычисления кинематических параметров, таких как скорость и ускорение. Поскольку видео обычно имеет фиксированную частоту дискретизации, то есть частоту кадров, и заданное число кадров, данные трекинга можно использовать для извлечения вибрационной информации из движения объекта. Выполняя направленное DFT по данным смещения, можно выделить доминирующие частоты, раскрывающие вибрационные моды отслеживаемого объекта. Интеграция трекинга объекта с постобработкой позволяет всесторонне анализировать как кинематическое, так и динамическое поведение, что делает этот подход особенно ценным для изучения механических вибраций, структурной динамики и других зависящих от времени явлений. Особенно важно, что это весьма доступный инструмент, поскольку для получения такой ключевой информации нужны только камера и компьютер с MS Excel или Python.

ЭТАПЫ ВЫПОЛНЕНИЯ DFT В MS EXCEL

Шаг 1. Сбор экспериментальных данных (запись высокоскоростного видео)

Запишите видео высокоскоростного события, обеспечив съёмку с частотой кадров, значительно превышающей частоту исследуемого объекта, например 10 кГц для объекта с частотой 1 кГц. В идеале такие измерения должны выполняться через равные интервалы времени, то есть при фиксированной частоте кадров. Частота кадров обычно записывается в метаданных файла .CINE.

Шаг 2. Использование ПО трекинга для извлечения координат X, Y во времени
Используйте программу трекинга объектов, чтобы отслеживать характерную точку на вибрирующем объекте во всех кадрах серии. В примере справа PCC 3.9 использовался для отслеживания аморфного признака, показанного квадратной вставкой, во времени на протяжении примерно 1200 кадров. Несмотря на сложные условия слабой освещённости, программа смогла отслеживать этот признак на протяжении всего видео с помощью сопоставления по шаблону. График справа показывает колебательное изменение координаты Y во времени для данной серии и построен с помощью графической функции PCC.

Шаг 3. Экспорт данных трекинга

Разместите данные отслеживаемого положения в формате, удобном для анализа, обычно в виде столбца значений в таблице, например Excel, или массива numpy при использовании Python. При экспорте из PCC 3.9 вы получите номер кадра и время от момента триггера вместе с отслеженными координатами X и Y. Эти столбцы данных, а именно время от триггера и смещение, будут служить основными наборами данных для анализа DFT. В некоторых случаях может быть предпочтительно использовать субпиксельную интерполяцию, однако в данном случае никакого апскейлинга или интерполяции не применялось.

Шаг 4. Удаление смещения (DC-компоненты)
Перед выполнением анализа DFT по данным трекинга удалите постоянную составляющую, то есть среднее значение, из данных X и Y, чтобы сосредоточиться на колебательном поведении. При работе в MS Excel это можно сделать, создав новый столбец, в котором из исходных данных вычитается AVERAGE(data). При использовании Python это удобнее всего делать в Numpy, где среднее значение вычисляется как numpy.mean(array). На рисунке слева показаны данные X и Y до и после вычитания среднего.

Шаг 5. Выполнение DFT-анализа данных смещения

В MS Excel убедитесь, что у вас включена надстройка Analysis ToolPak. Это делается через File > Options > Add-ins > Manage > Excel Add-ins > Go > Check Analysis ToolPak. После этого на вкладке Data при нажатии Data Analysis, крайняя правая функция, становится доступным инструмент Fourier Analysis.

Этапы выполнения DFT в MS Excel

Организуйте данные по столбцам:
• Столбец C: данные времени, то есть временной шаг для каждого измерения.
• Столбцы E и F: ваши данные во временной области, например смещение, амплитуда сигнала и т. д. В данном случае это координаты X и Y.
• Примечание: убедитесь, что число точек данных N является степенью двойки, например 16, 32, 64, 128, 256 и т. д., поскольку это оптимизирует вычисление DFT. Если число точек не является степенью двойки, можно дополнить данные нулями, чтобы общее количество стало степенью двойки. Если выбрать количество точек, не являющееся степенью двойки, MS Excel запросит корректное значение.

Шаг 3. Выполнение преобразования Фурье (FT)

1. Перейдите на вкладку Data на ленте.
2. В группе Analysis нажмите Data Analysis.
3. В окне Data Analysis прокрутите список вниз, выберите Fourier Analysis и нажмите OK.
4. В окне Fourier Analysis:
   o Input Range: выберите диапазон ваших данных во временной области, например A1:A4096.
   o Output Range: укажите пустой диапазон, куда будут выведены результаты FFT, например B1.
   o Select Complex Output: установите флажок Complex.
5. Нажмите OK, чтобы запустить DFT.

Шаг 4. Интерпретация результата
Результаты FT будут представлять собой комплексные числа, то есть действительную и мнимую части, в выходном диапазоне.
• Действительная часть: представляет косинусную составляющую каждой частоты.
• Мнимая часть: представляет синусную составляющую каждой частоты.
• Модуль: вычислите модуль в новом столбце, используя формулу SQRT(G1^2 + H1^2), где G — действительная часть, а H — мнимая часть.
• Частотная область: частота, соответствующая каждому результату FT, определяется частотой дискретизации, то есть частотой кадров, и количеством обработанных точек. Если частота дискретизации равна Fs, а число отсчётов равно N, то шаг по частоте для каждого бина FT равен Fs / N. Таким образом, если камера была настроена на 1000 fps и вы обрабатываете 4096 кадров, шаг по частоте составит 1000 / 4096 Гц, а первая строка будет соответствовать 0 Гц, то есть DC-компоненте.

Шаг 6. Построение частотного спектра

1. Постройте график модулей в зависимости от частоты: выберите значения частоты, вычисленные выше, и соответствующие им модули, вычисленные ранее. Справа показан пример такого результата.
2. Подпишите оси:
   o Ось X: частота, Гц
   o Ось Y: модуль

ЭТАПЫ ВЫПОЛНЕНИЯ DFT В PYTHON

Тот же подход, что и выше, можно использовать и в среде Python.

Краткое описание процедуры:

1. Импорт библиотек: импортируйте необходимые библиотеки Python, такие как numpy, scipy, matplotlib и pandas.
2. Загрузка данных трекинга: импортируйте данные трекинга с помощью pandas, например из CSV или txt.
3. Извлечение нужных данных: выберите из данных столбцы времени и положения.
4. Коррекция по среднему: вычтите среднее значение из данных положения, чтобы убрать постоянную составляющую.
5. Применение DFT: используйте быстрое преобразование Фурье, FFT, чтобы перевести данные из временной области в частотную.
6. Формирование частотных данных: определите шаг частотной области, разделив частоту кадров на общее количество анализируемых отсчётов.
7. Выделение пиков: найдите пики в частотном спектре, соответствующие доминирующим вибрационным модам.
8. Визуализация: постройте графики сигнала во временной области и частотного спектра, отметив пики.

Пример кода на Python

import numpy as np
import matplotlib.pyplot as plt
from scipy.signal import find_peaks
# Function to perform DFT, plot the results, and label the largest three peaks with different colors
def plot_dft_analysis_with_top_peaks(data, sampling_rate):
# Perform the Discrete Fourier Transform (DFT)
n = len(data)
dft_result = np.fft.fft(data)
# Compute the corresponding frequency values
frequencies = np.fft.fftfreq(n, d=1/sampling_rate)
# Get the magnitude of the DFT (to represent the frequency content)
dft_magnitude = np.abs(dft_result)
# Only consider the positive frequencies (first half of the DFT result)
positive_frequencies = frequencies[:n//2]
positive_dft_magnitude = dft_magnitude[:n//2]
# Find the peaks in the frequency spectrum (using scipy's find_peaks)
peaks, _ = find_peaks(positive_dft_magnitude)
# Get the indices of the largest three peaks
largest_peaks_indices = np.argsort(positive_dft_magnitude[peaks])[-3:]
largest_peaks = peaks[largest_peaks_indices]
# Plot the original temporal data
P a g e | 5 of 5 Using MS Excel & Python to Extract Vibrational Modes from Video
plt.subplot(1, 2, 1)
plt.plot(np.arange(n) / sampling_rate, data, label="Original Data")
plt.xlabel('Time (s)', weight='bold')
plt.ylabel('Amplitude', weight='bold')
plt.title('Original Temporal Data')
plt.grid(False)
# Plot the frequency spectrum (Magnitude of DFT)
plt.subplot(1, 2, 2)
plt.plot(positive_frequencies, positive_dft_magnitude, label="Frequency Spectrum")
# Mark and label the largest three peaks in the frequency plot with different colors
colors = ['red', 'green', 'blue']
for i, peak in enumerate(largest_peaks):
plt.scatter(positive_frequencies[peak], positive_dft_magnitude[peak], color=colors[i], label=f"Peak {i+1} ({positive_frequencies[peak]:.2f} Hz)")
plt.text(positive_frequencies[peak]+85, positive_dft_magnitude[peak] - 5,
f'{positive_frequencies[peak]:.2f} Hz', fontsize=8, color=colors[i], ha='center')
plt.xlabel('Frequency (Hz)', weight='bold')
plt.ylabel('Magnitude', weight='bold')
plt.title('Frequency Spectrum')
plt.grid(False)
#plt.legend()
plt.tight_layout()
plt.show()
# Example usage
if __name__ == "__main__":
# Example temporal data (e.g., three sinusoidal signals with noise)
sampling_rate = 1000 # in Hz
t = np.linspace(0, 1, sampling_rate, endpoint=False) # 1 second of data
# Frequencies of the three sine waves
freq_signal1 = 50 # Frequency of the first sine wave (Hz)
freq_signal2 = 150 # Frequency of the second sine wave (Hz)
freq_signal3 = 325 # Frequency of the third sine wave (Hz)
# Generate three sine waves with different frequencies and some noise
signal1 = .75*np.sin(2 * np.pi * freq_signal1 * t)
signal2 = .5*np.sin(2 * np.pi * freq_signal2 * t)
signal3 = .6*np.sin(2 * np.pi * freq_signal3 * t)
# Combine all signals and add some noise
signal = signal1 + signal2 + signal3 + 0.5 * np.random.randn(sampling_rate)
# Call the function to plot the DFT analysis with the top 3 peaks

Ожидаемый результат