
簡介面向設備健康監測與故障預測場景的MATLAB振動故障診斷源碼包適合機械工程、自動化及故障診斷方向的學生與工程師學習也可作為課程設計或畢業設計的參考。資源圍繞振動信號時域、頻域及復頻域分析框架覆蓋數據預處理、FFT、小波變換、希爾伯特變換等算法實現以及峭度、峰值因子等統計特征提取和SVM、神經網絡等模式識別建模流程。共26個文件以25個.m腳本為主輔以1個說明文檔壓縮包僅49KB輕量且便于直接運行驗證。目前已有337人學習下載。通過運行源碼可系統掌握從振動信號處理到故障類型識別的完整鏈路理解頻譜異常、幅值增大等特征背后的機械故障機理還可在示例基礎上擴展自己的診斷方法是理論與實踐結合的優質參考。1. 振動故障診斷為什么繞不開MATLAB把加速度傳感器貼在電機軸承座上采一段振動數據時域波形看起來只是雜亂無章的“毛刺”。但同樣的信號經過傅里葉變換展開到頻域轉頻、倍頻、邊頻帶依次浮現軸承外圈故障的特征頻率就藏在某個峰值下面——這就是振動故障診斷的基本盤。MATLAB在這個領域里幾乎是默認工作臺原因不復雜它同時具備信號處理工具箱、統計工具箱和深度學習工具箱從原始波形到故障判據可以在同一個腳本里走完。對于旋轉機械、齒輪箱、泵和風機這類設備的預測性維護振動故障診斷的價值在于把“壞了再修”變成“有趨勢再停機”。這篇內容面向做設備健康管理、狀態監測或工業數據分析的工程師目標是讓你拿到振動數據后能自己跑通從頻譜分析到故障判據的整條鏈路。2. 振動故障診斷技術的原理從時域到頻域的信號解構2.1 時域特征為什么容易誤判振動故障診斷技術最早期的做法是看時域指標比如峰值、均方根值RMS、峭度Kurtosis和波形因子。RMS 反映振動能量總體水平峭度對沖擊型故障敏感這套指標在設備劣化趨勢監測里至今仍在使用。但它的局限性很明確時域指標是全局統計量無法區分“轉子不平衡”和“軸承外圈點蝕”這兩種完全不同性質的故障它們可能表現出相近的 RMS 增量。換句話說時域指標能告訴你“設備狀態變了”但說不清“哪里變了”。2.2 FFT頻譜分析把振動信號拆成頻率分量頻域分析是振動故障診斷技術的核心支撐。對一段離散振動信號 x[n] 做離散傅里葉變換得到的是各頻率分量的幅值分布轉頻 f_r 及其諧波、齒輪嚙合頻率 f_m f_r × ZZ 為齒數、軸承故障特征頻率都會以譜峰形式出現。MATLAB 里一次fft(x)就能完成變換但有三個參數必須自己確認采樣率 Fs、采樣點數 N、窗函數類型。頻率分辨率 Δf Fs / N采樣時長 T N / Fs 決定你能分辨多近的兩根譜線。例如 Fs 12800 Hz采樣 1 秒Δf 就是 1 Hz這通常夠用來分離轉頻和軸承故障頻率如果采樣只采 0.1 秒Δf 變成 10 Hz兩根相距 5 Hz 的譜峰就會糊成一根故障特征直接丟失。Fs 12800; % 采樣率 12.8 kHz N Fs * 1; % 采樣 1 秒 t (0:N-1) / Fs; % 模擬一個轉頻 25 Hz、帶 2 倍頻成分的信號 x 1.0 * sin(2*pi*25*t) 0.4 * sin(2*pi*50*t) 0.3 * randn(1, N); X fft(x); f (0:N-1) * Fs / N; % 頻率軸單邊譜取前 N/2 個點 X_mag abs(X(1:N/2)) / (N/2); f_half f(1:N/2); plot(f_half, X_mag); xlim([0 200]);這段代碼里X_mag的歸一化是除以 N/2原因是 MATLAB 的fft結果是對稱的單邊譜幅值要乘 2 再除以 N 才能還原真實幅值。頻率軸f (0:N-1) * Fs / N的單位是 Hz取前 N/2 個點對應奈奎斯特頻率以內的正頻率部分。如果不做歸一化譜峰幅值會隨采樣點數變化后續做幅值對比就失去意義。2.3 包絡譜與階次分析定位軸承和齒輪故障頻率滾動軸承的故障特征頻率不是直接出現在原始頻譜上的。外圈點蝕產生的沖擊會激發軸承座固有頻率形成高頻調制載波是結構共振頻率包絡是故障沖擊頻率。直接看原始頻譜看到的是一簇高頻共振峰故障頻率藏在邊頻帶里。處理辦法是先用帶通濾波器把共振頻帶切出來再做 Hilbert 變換得到包絡信號最后對包絡信號做 FFT得到包絡譜。故障特征頻率在這個譜圖上以明顯峰值出現。各部件特征頻率與轉頻 f_r 的關系如下表Z 為滾動體個數d 為滾動體直徑D 為節圓直徑α 為接觸角故障位置特征頻率計算公式外圈故障 (BPFO)f_o Z × f_r / 2 × (1 - d/D × cosα)內圈故障 (BPFI)f_i Z × f_r / 2 × (1 d/D × cosα)滾動體故障 (BSF)f_b D × f_r / (2d) × (1 - (d/D × cosα)^2)保持架故障 (FTF)f_c f_r / 2 × (1 - d/D × cosα)齒輪嚙合頻率f_m f_r × Z包絡譜分析對軸承早期故障尤其敏感因為沖擊能量分散在寬帶共振里包絡解調把分散的能量重新集中到故障特征頻率上。內圈故障調制現象比外圈更明顯包絡譜上會出現 f_i 以及 f_i ± f_r 的邊頻帶這個邊頻特征本身就是判據。齒輪故障則看嚙合頻率 f_m 兩側的邊頻帶間隔齒面磨損時 f_m 幅值上升斷齒時 f_m ± f_r 邊頻帶顯著增強。2.4 時頻分析非平穩信號的補充手段FFT 假設信號是平穩的但變轉速工況下的振動信號頻率隨時間變化直接用 FFT 會把瞬態成分平均掉。短時傅里葉變換STFT把信號切成短窗口再做 FFT得到時間-頻率二維譜圖能看出頻率成分隨時間的變化軌跡。MATLAB 里spectrogram(x, window, noverlap, nfft, Fs)一行就能出圖。但這個工具箱函數的輸出是復數矩陣調用時要留意參數順序window是窗函數向量noverlap是重疊點數nfft是 FFT 點數Fs用于頻率軸標注。變轉速場景下階次跟蹤比 STFT 更實用恒定轉速設備用 STFT 或小波變換已經足夠。3. 用MATLAB實現振動故障診斷的最小流程3.1 從CSV導入振動數據到MATLAB工業現場采集的振動數據大多以 CSV 格式保存常見布局是每行一個采樣點列包含時間戳和加速度值也可能一列存一個測點序列。readtable讀入后先確認變量類型時間列是datetime類型振動數據必須是double。如果 CSV 里時間戳是字符串用datetime(str, InputFormat, yyyy-MM-dd HH:mm:ss)轉換。導入后第一件事不是做頻譜而是檢查數據質量有沒有 NaN、有沒有毛刺、傳感器有沒有飽和削頂。data readtable(vibration_data.csv); t_raw data.timestamp; x_raw data.acc_x; % 假設加速度列名為 acc_x % 檢查缺失值 fprintf(NaN count: %d\n, sum(isnan(x_raw))); % 去除趨勢項消除傳感器零漂 x_detrend detrend(x_raw, linear);detrend默認消除線性趨勢對零漂和緩慢熱漂移有效。注意detrend不改變采樣率只改幅值基線。如果數據里有明顯的離群尖峰幅值超過均值加 10 倍標準差建議先做中值濾波或直接剔除該段數據再做后續分析否則 FFT 后會出現全頻段抬高。3.2 用濾波器把關注頻帶切出來包絡譜分析的第一步是帶通濾波頻帶選擇直接決定分析質量。如果已知軸承座的共振頻帶可通過敲擊測試或歷史數據分析獲得直接在該頻帶內濾波。未知時常見做法是用傅里葉變換先看共振峰集中在哪個頻段再選擇覆蓋它的頻帶。MATLAB 里推薦用designfilt而不是老式的butterfiltfilt前者的參數校驗更友好濾波器階數和波紋指標可以在一條語句里配完。bpFilt designfilt(bandpassiir, FilterOrder, 4, ... HalfPowerFrequency1, 2000, ... HalfPowerFrequency2, 6000, ... SampleRate, Fs); x_band filtfilt(bpFilt, x_detrend);這里HalfPowerFrequency1和HalfPowerFrequency2是通帶邊界頻率單位 Hz對應 -3 dB 點。FilterOrder為 4 時濾波器的陡降特性比較溫和適合現場數據如果頻帶很窄或者干擾很強可以提高到 8但要留意filtfilt是零相位濾波不產生相位偏移這是包絡分析必須的條件——Hilbert 變換對相位敏感任何非線性相位誤差都會污染包絡波形。3.3 計算包絡譜并標記特征頻率濾波完成后做 Hilbert 變換取包絡再對包絡信號做 FFT得到包絡譜。現在的關鍵操作是在包絡譜上畫出特征頻率的標記線并讀取目標頻點的幅值。env abs(hilbert(x_band)); % Hilbert 變換取包絡 % 去掉包絡中的直流分量再做 FFT env_ac env - mean(env); N_env length(env_ac); EnvSpec abs(fft(env_ac)) / (N_env/2); f_env (0:N_env/2-1) * Fs / N_env; EnvSpec_half EnvSpec(1:N_env/2); % 定義特征頻率轉頻已知 25 Hz軸承外圈故障系數 3.05 為例 fr 25; coef_BPFO 3.05; % 外圈故障系數 f_BPFO fr * coef_BPFO; hold on; xline(f_BPFO, --r, BPFO); [~, idx] min(abs(f_env - f_BPFO)); fprintf(BPFO amplitude: %.3f\n, EnvSpec_half(idx));這段代碼的最后一步用最小絕對值誤差找到離理論特征頻率最近的譜線讀取該點幅值。xline在圖上畫豎線標記--r是紅色虛線樣式。注意hilbert返回的是解析信號abs(hilbert(x))得到的才是包絡。包絡信號的幅值單位與原始信號相同但物理意義已經變成沖擊幅度的時變軌跡譜峰幅值代表該頻率沖擊成分的強度。3.4 輸出特征矩陣用于趨勢記錄單次診斷結果不足以支撐決策振動故障診斷的工程價值在于趨勢同一個特征頻率的幅值隨時間的增長率比它的絕對大小更重要。每次分析后把特征頻率幅值寫入一張表追加保存后續可以做劣化曲線。result_row table(datetime(now), fr, f_BPFO, EnvSpec_half(idx), ... VariableNames, {Time, RotSpeedHz, BPFO_Hz, BPFO_Amp}); if isfile(diagnosis_history.csv) writetable(result_row, diagnosis_history.csv, WriteMode, Append); else writetable(result_row, diagnosis_history.csv); endWriteMode參數在writetable中的取值是Append或Overwrite。第一次運行時文件不存在直接writetable之后再運行時加WriteMode, Append追加行。這樣積累兩周數據后用plot畫 BPFO 幅值隨時間的曲線劣化趨勢比任何單次閾值判斷都可靠。4. 從特征到故障識別閾值判據與機器學習4.1 基于特征頻率比對的閾值判據振動故障診斷最常見的落地方式不是訓練模型而是“頻率比對 幅值判據”把包絡譜或頻譜里的顯著峰值頻率與理論特征頻率對比誤差在 ±0.5% 以內視為匹配匹配后看幅值是否超過報警閾值。閾值設定有兩種來源一是 ISO 10816 這類標準給出的振動速度均方根值分段二是設備自身歷史數據的 90 分位數。前者偏保守適合通用設備后者更貼合實際但要保證歷史數據覆蓋正常運行工況。幅值判據有一個容易出錯的地方不同測點、不同工況下同一故障的絕對幅值差異巨大。軸承座垂直方向和水平方向的幅值可能差一倍以上負載 50% 和 100% 時特征頻率幅值也沒有可比性。因此閾值必須按測點和工況分組設定不能全局共用一個值。實際做法是在特征矩陣里增加工況列轉速、負載按組別統計閾值實時數據只跟同組歷史數據比較。4.2 用MATLAB訓練一個簡單的故障分類器當設備種類多、故障模式多頻率比對腳本會變成一大串 if-else維護成本高。這時可以用機器學習方法替代人肉判據。常見做法是提取每組數據的特征向量——包括時域的 RMS、峭度、峰值因子頻域的 1 倍頻幅值、2 倍頻幅值、特征頻率幅值、邊頻帶能量——構成特征矩陣用帶標簽的歷史數據訓練分類器。MATLAB 的fitcecoc多類 SVM和fitcknnK 近鄰是兩類典型的入門選項。% 假設 feat_mat 是 n x m 的特征矩陣label_vec 是 n x 1 的分類標簽 % 數據先劃分訓練集和測試集 cv cvpartition(label_vec, HoldOut, 0.3); X_train feat_mat(training(cv), :); Y_train label_vec(training(cv), :); X_test feat_mat(test(cv), :); Y_test label_vec(test(cv), :); % 特征標準化對 SVM 和 KNN 至關重要 [X_train_s, mu, sigma] zscore(X_train); X_test_s (X_test - mu) ./ sigma; mdl fitcecoc(X_train_s, Y_train); pred predict(mdl, X_test_s); acc sum(pred Y_test) / length(Y_test); fprintf(Classification accuracy: %.2f%%\n, acc * 100);這段代碼里的特征標準化是決定分類效果的關鍵點。SVM 和 KNN 都依賴樣本間距離如果某個特征如 RMS量綱遠大于其他特征模型會忽略小量綱特征。zscore把每個特征變成零均值單位方差保證各特征在距離計算中權重一致。cvpartition的HoldOut參數表示隨機抽取 30% 樣本做測試集分類器在剩余 70% 上訓練。這類模型的輸出結果是“故障類別”適合故障模式明確的場景如果故障樣本很少類別不均衡會嚴重降低準確率需要配合Prior參數設置先驗概率或使用fitcensemble集成方法。4.3 深度學習的邊界樣本量與可解釋性深度學習近年被引入振動故障診斷領域MATLAB 的深度學習工具箱也提供了從 LSTM 到 Transformer 的完整工具鏈。基于 BiLSTM 的時序分類模型確實能直接從原始振動波形學習特征省去人工特征工程在學術數據集上準確率可以達到 95% 以上。但工程現場的現實是故障樣本極其稀少一個正常運行的工廠很難采集到“外圈點蝕”“保持架斷裂”“齒輪斷齒”各幾百段標注數據。深度學習模型在樣本不足時容易過擬合預測結果又缺乏可解釋性——它告訴你“這是內圈故障”但不告訴你依據哪條特征頻率維修人員沒法核驗。因此在工程部署中深度學習更適合作為“第二層判斷”先用頻譜比對和包絡譜判據篩出可疑樣本再用深度模型在疑似樣本里做細分類或者反過來用深度模型做異常檢測發現偏離正常分布的數據后轉人工用頻譜分析確認故障類型。特征頻率計算是物理依據機器學習是統計規律兩者結合是目前工業界落地最穩的路徑。5. 調參、驗證與MATLAB診斷里的易踩坑5.1 采樣率與頻譜分辨率先定這兩個參數現場采集參數如果定錯了后續所有分析都是白費。采樣率至少要覆蓋到你需要觀察的最高頻率的 2.56 倍這是工程上考慮抗混疊濾波器過渡帶后的常用系數。比如要分析到 5 kHz采樣率至少 12.8 kHz。頻率分辨率取決于采樣時長Δf Fs / N 1 / T要分辨 1 Hz 以內的邊頻帶采樣時長至少 1 秒。高速設備轉頻 100 Hz 以上建議采樣時長覆蓋至少 20 轉低速設備轉頻 1 Hz則需要 30-60 秒數據才能把故障特征頻率分辨清楚。MATLAB 里可以用N 2^nextpow2(Fs * T)來取 2 的冪次長度配合fft計算效率更高。5.2 驗證診斷結果的三條路徑診斷結論不能只靠一次分析就下。第一個驗證方法是仿真信號驗證用正弦疊加沖擊脈沖構造已知特征的信號跑同一套分析代碼確認算法本身沒有 bug。第二個方法是多測點交叉驗證同一個軸承故障在軸向和徑向測點上的特征頻率相同但幅值分布不同如果只在一個測點看到特征峰另一個測點沒有可能是傳感器安裝位置距離故障源太遠。第三個方法是趨勢驗證連續監測數據里特征頻率幅值呈單調上升趨勢比單次高幅值更可信啟動、停機、變載瞬間的瞬態沖擊也可能產生偽譜峰但這些峰不會穩定出現在每次測量的同一頻率位置。5.3 MATLAB實現中的5個常見坑第一個坑是fft結果的幅值歸一化很多人直接把abs(fft(x))當幅值用導致同一信號不同采樣點數時幅值不一致。正確做法是單邊譜乘 2 除 N。第二個坑是頻率軸計算錯誤f (0:N-1) * Fs / N得到的是從 0 到 Fs-Fs/N 的雙邊頻率軸畫單邊譜時只取前 N/2 個點去掉后面一半如果直接畫完整頻率軸會看到對稱的鏡像頻譜誤以為存在兩倍數量的頻率成分。第三個坑是濾波器的相位污染用filter而不是filtfilt做帶通濾波會讓包絡產生幾十個采樣點的延遲導致后續頻譜分析的相位信息失真。第四個坑是 CSV 導入時數據類型的隱式轉換某些 CSV 的加速度列是整數類型readtable讀進來后double()轉換不及時后續fft計算時出現意外截斷。第五個坑是抗混疊濾波器的前提被忽略如果數據采集時沒有硬件抗混疊濾波高于采樣率一半的信號會折疊進低頻段頻譜上出現虛假峰值對這種數據先用lowpass(x, 0.48*Fs, Fs)做一次軟件抗混疊是補救手段但效果取決于原始頻譜中高頻成分的強度。驗證仿真代碼是對抗這些坑最直接的訓練。用已知頻率和幅值的正弦信號注入分析流程對比輸出值與理論值的誤差任何歸一化、頻率軸或濾波器的錯誤都會在誤差里暴露。振動故障診斷做久了會發現算法本身不難難的是每一步都確認“這個數的單位是什么、物理意義是什么”把 MATLAB 當計算器用還是當分析工具用區別就在這些細節里。本文還有配套的精品資源點擊獲取