最新のリリースでは、このページがまだ翻訳されていません。 このページの最新版は英語でご覧になれます。
rpmtrack
振動信号の RPM プロファイルの追跡および抽出
構文
rpm = rpmtrack(x,fs,order,p)rpm = rpmtrack(xt,order,p)rpm = rpmtrack(___,Name,Value)[rpm,tout] = rpmtrack(___)rpmtrack(___)説明
rpm = rpmtrack(___,Name,Value)
出力引数なしの rpmtrack(___) は、パワー時間-周波数マップおよび推定 RPM プロファイルを対話型 Figure 上にプロットします。
例
振動信号の RPM プロファイル
3 つの高調波成分を含む振動信号を生成します。信号は 1 kHz で 16 秒間サンプリングされます。この信号の瞬時周波数は、エンジンの試運転と惰行試験に似ています。台形則を使用して周波数を積分することにより、瞬間位相を計算します。
fs = 1000; t = 0:1/fs:16; ifq = 20 + t.^6.*exp(-t); phi = 2*pi*cumtrapz(t,ifq);
信号の高調波成分は、1、2 および 3 の次数に対応します。次数 2 の正弦波の振幅は他の 2 倍です。
ol = [1 2 3]; amp = [5 10 5]; vib = amp*cos(ol'.*phi);
次数 2 のリッジ上の点を使用して信号の RPM プロファイルを抽出し、可視化します。
time = 3; order = 2; p = [time order*ifq(t==time)]; rpmtrack(vib,fs,order,p)
空ぶかしエンジンの RPM プロファイル
自動車エンジンの空ぶかしによって発生する振動に似た信号を生成します。この信号は 1 kHz で 30 秒間サンプリングされ、それぞれ振幅 5、4、および 0.5 で次数が 1、2.4、および 3 である 3 つの高調波成分が含まれています。信号を単位分散のホワイト ガウス ノイズに組み込み、MATLAB® timetable に保存します。瞬時周波数を 60 倍にして RPM プロファイルを取得します。RPM プロファイルをプロットします。
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 を使用して振動信号から RPM プロファイルを導きます。5 秒間隔で 4 ポイントを使用して 2.4 の次数に対応するリッジを指定します。出力 timetable の概要を表示します。
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.5
Max 8879.5
再構成した RPM プロファイルと再構成に使用したポイントをプロットします。
hold on plot(seconds(rpmest.tout),rpmest.rpm,'.-') plot(t(ndx),fis(ndx)*60,'ok') hold off legend('Original','Reconstructed','Ridge points','Location','northwest')
抽出した RPM プロファイルを使用して信号の次数-RPM マップを生成します。
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])
ファン スイッチオフ RPM プロファイル
スイッチオフ後ファン ブレードが遅くなっていくときの RPM プロファイルを推定します。
20,000 rpm で回転する工業用天井型扇風機をオフにします。空気抵抗 (ベアリングの摩擦からの寄与は無視できるとします) により約 6 秒でファンのローターが停止します。高速カメラでファン ブレードの 1 つの x 座標を 1 kHz のレートで測定します。
fs = 1000; t = 0:1/fs:6-1/fs; rpm0 = 20000;
このファン ブレードを、ローターを中心として半径 50 cm を回転する質点として理想化します。ブレードには速度に比例する抵抗力がかかり、その結果、位相角は次の式で求められます。
ここで、
は初期周波数、
秒は減衰時間です。
a = 0.5; f0 = rpm0/60; T = 0.75; phi = 2*pi*f0*T*(1-exp(-t/T));
ブレードの x 座標と y 座標を計算し、プロットします。分散 
のホワイト ガウス ノイズを追加します。
x = a*cos(phi) + randn(size(phi))/10; y = a*sin(phi) + randn(size(phi))/10; plot(t,x,t,y)
関数 rpmtrack を使用して、RPM プロファイルを決定します。次のコマンド
rpmtrack(x,fs)
をコマンド ラインで入力し、対話型 Figure を開きます。
スライダーを使用して時間-周波数マップの周波数分解能を約 11 Hz に調節します。信号成分が次数 1 に相当すると仮定し、リッジ抽出の終了時間を 3.0 秒に設定します。時間-周波数マップで十字カーソルを使用し、[追加] ボタンを使用してリッジ上に 3 つのポイントを追加します (または、選択した場所でカーソルをダブルクリックしてポイントを追加します)。[推定] をクリックして RPM プロファイルを追跡および抽出します。
RPM プロファイルが指数関数的に減衰することを確認します。[エクスポート] タブで、[エクスポート] をクリックして [MATLAB スクリプトの生成] を選択します。スクリプトがエディターで開きます。
% 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);
スクリプトを実行します。RPM プロファイルを片対数プロットに表示します。
semilogy(tOut,rpmOut) ylim([500 20000])
入力引数
出力引数
アルゴリズム
rpmtrack は、2 ステップ (おおまか-細かい) の推定手法を使用します。
参照
[1] Urbanek, Jacek, Tomasz Barszcz, and Jerome Antoni. "A Two-Step Procedure for Estimation of Instantaneous Rotational Speed with Large Fluctuations." Mechanical Systems and Signal Processing. Vol. 38, 2013, pp. 96–102.
