Сгенерируйте сигнал, который напоминает колебания, вызванные путем газования на автомобильном механизме. Сигнал производится на уровне 1 кГц в течение 30 секунд и содержит три гармонических компонента порядков 1, 2.4, и 3, с амплитудами 5, 4, и 0.5, соответственно. Встройте сигнал в белый Гауссов шум модульного отклонения и сохраните его в расписании MATLAB®. Умножьте мгновенную частоту на 60, чтобы получить профиль об/мин. Постройте профиль об/мин.

fs = 1000;
t = (0:1/fs:30)';

fit = @(a,x) (t-x).^6.*exp(-(t-x)).*((t-x)>=0)*a';

fis = fit([0.4 1 0.6 1],[0 6 13 17]);
phi = 2*pi*cumtrapz(t,fis);

ol = [1 2.4 3];
amp = [5 4 0.5]';
vib = cos(phi.*ol)*amp + randn(size(t));

xt = timetable(seconds(t),vib);

plot(t,fis*60)

Используйте rpmtrack функция, чтобы вывести профиль об/мин из сигнала вибрации. Используйте четыре точки в 5-секундных интервалах, чтобы задать гребенчатое соответствие к порядку 2.4. Отобразите сводные данные выходного расписания.

ndx = (5:5:20)*fs;
order = ol(2);

p = [t(ndx) order*fis(ndx)];

rpmest = rpmtrack(xt,order,p);

summary(rpmest)

RowTimes:

    tout: 30001x1 duration
        Values:
            Min           0 sec     
            Median        15 sec    
            Max           30 sec    
            TimeStep      0.001 sec 

Variables:

    rpm: 30001x1 double

        Values:

            Min       2.2204e-16
            Median        4327.2
            Max           8879.8

Постройте восстановленный профиль об/мин и точки, используемые в реконструкции.

hold on
plot(seconds(rpmest.tout),rpmest.rpm,'.-')
plot(t(ndx),fis(ndx)*60,'ok')
hold off
legend('Original','Reconstructed','Ridge points','Location','northwest')

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

rpmordermap(vib,fs,rpmest.rpm)

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

xrc = orderwaveform(vib,fs,rpmest.rpm,ol);

figure
plot(t,xrc)
legend([repmat('Order = ',[3 1]) num2str(ol')])
xlim([5 20])

Оцените профиль об/мин лопасти вентилятора, когда он замедляется после switchoff.

Выключен промышленный вентилятор крыши, вращающийся на уровне 20 000 об/мин. Сопротивление воздуха (с незначительным вкладом от подшипника трения) заставляет ротор вентилятора остановиться приблизительно через 6 секунд. Высокоскоростная камера измеряет x-координату одной из лопастей вентилятора на уровне 1 кГц.

fs = 1000;
t = 0:1/fs:6-1/fs;

rpm0 = 20000;

Идеализируйте лопасть вентилятора как массу точки окружение центра ротора в радиусе 50 см. Лопатка испытывает силу сопротивления, пропорциональную, чтобы ускориться, приводя к следующему выражению для угла фазы:

где ${\mathit{f}}_0$ начальная частота и $\mathit{T}=0.75$ второй время затухания.

a = 0.5;
f0 = rpm0/60;
T = 0.75;

phi = 2*pi*f0*T*(1-exp(-t/T));

Вычислите и постройте x-и y-координаты лопатки. Добавьте белый Гауссов шум отклонения $0.1^2$.

x = a*cos(phi) + randn(size(phi))/10;
y = a*sin(phi) + randn(size(phi))/10;

plot(t,x,t,y)

Используйте rpmtrack функция, чтобы определить профиль об/мин. Ввод

rpmtrack(x,fs)

в командной строке, чтобы открыть интерактивную фигуру.

Используйте ползунок, чтобы настроить, разрешение частоты частоты времени сопоставляют приблизительно с 11 Гц. Примите, что компонент сигнала соответствует порядку 1 и установил время окончания для гребенчатой экстракции к 3,0 секундам. Используйте курсор перекрестия в карте частоты времени и кнопке Add, чтобы добавить три точки, лежащие на гребне. (В качестве альтернативы дважды кликните курсор, чтобы добавить точки в местоположениях, которые вы выбираете.) Нажимают Estimate, чтобы отследить и извлечь профиль об/мин.

Проверьте, что профиль об/мин затухает экспоненциально. На вкладке Export нажмите Export и выберите Generate MATLAB Script. Скрипт появляется в Редакторе.

% MATLAB Code from rpmtrack GUI

% Generated by MATLAB 9.4 and Signal Processing Toolbox 8.0

% Generated on 18-Dec-2017 19:00:35

% Set sample rate
fs = 1000.0000;
% Set order of ridge of interest
order = 1.0000;
% Set ridge points on ridge of interest
ridgePoints = [...
    0.4501 179.8246;...
    0.9944 88.5965;...
    2.4161 11.4035];
% Estimate RPM
[rpmOut,tOut] = rpmtrack(x,fs,order,ridgePoints,...
    'Method','stft',...
    'FrequencyResolution',11.1612,...
    'PowerPenalty',Inf,...
    'FrequencyPenalty',0.0000,...
    'StartTime',0.0000,...
    'EndTime',3.0000);

Запустите скрипт. Отобразите профиль об/мин в полулогарифмическом графике.

semilogy(tOut,rpmOut)
ylim([500 20000])