0 引言
    傳統化石能源的日益匱乏，使得以太陽能、風能、海浪能等新型清潔可再生能源獲得前所未有的發展機遇。中國擁有狹長的海岸線以及不計其數大小島嶼，擁有豐富的波浪能儲備。據粗略預估，我國波浪能的儲量為7 000萬千瓦，故而開發的前景十分明朗。波浪能發電(wave power generation)是以海洋波浪的能量為動力生產出電能。雖然目前海浪發電的裝置有很多種，但是成本投入大，并且總轉換效率比較低。一種有效提高轉換效率的思路即是使發電裝置的機械振動固有頻率和海浪的顯著頻率實時保持一致，以實現機械裝置對海浪主動共振，便可提取最大的能量[1]。這類發電裝置迫切需要實時獲取海浪頻率等特征參數。本文則采用了一種基于MEMS運動傳感器的方法，通過捕捉海浪的運動過程來建立模型，搭建軟硬件平臺采集加速度和陀螺儀等傳感器捕捉到的運動信號并使用時頻域分析方法，完成海浪特征參數的提取。
1 系統整體設計
    整套系統的設計包括監測終端的硬件設計部分和軟件算法實現兩個部分。硬件監測終端如圖1所示，主要是由微控制器、加速度計、陀螺儀、GSM無線模塊、SD數據存儲模塊以及電源電路組成。

圖1  硬件框架圖

    軟件部分框架如圖2所示，包括MEMS模塊運動數據的采集、基于互補濾波的狀態融合、對數據進行時頻域分析提取特征參數，并通過GPRS將數據無線傳送至遠程服務器。

圖2  軟件流程框圖

2 硬件設計
　　整套硬件系統采用ST公司推出的基于ARM Cortex-M4F內核STM32F407VG作為主控制芯片。它擁有符點運算能力和增強的DSP處理指令，在168 MHz的主頻下可以達到210 DMIPS的處理能力,其高達1 MB的閃存（Flash）、196 KB的內嵌RAM以及豐富的IO口和外設可以滿足系統設計的需要。
    MEMS運動傳感器則采用InvenSense公司MPU6050和Honeywell公司的HMC5883L。前者可以提供加速度計和陀螺儀數據，后者提供三軸磁場強度。MPU6050為全球首例整合性6軸運動處理組件，整合了16位高精度加速度計和陀螺儀，并提供擴展I2C接口外接電子羅盤以便構成9軸運動系統。單片機使用I2C接口與傳感器通信。其電路圖如圖3所示。

圖3  MEMS傳感器

    通信采用SIM900A工業級的GSM模塊，以便于將海浪特征數據通過GPRS形式使用TCP協議發送至遠程服務器端。單片機可以通過串口與GSM模塊相連接，使用AT指令集來操作GSM模塊。其原理圖如圖4所示。

圖4  GSM模塊外部電路

    為了進行算法研究，需要將實時采集的數據保存在SD卡，然后在PC平臺用MATLAB等工具處理數據。單片機可以通過SPI接口和SDIO接口來實現SD卡的讀寫，此處選擇SDIO接口可以獲得更快的讀寫速度。單片機系統移植FatFS文件系統來完成SD卡的文件讀寫。FatFs是一個可移植的通用的FAT文件系統，適用于小型嵌入式系統。SD卡與單片機的接口原理圖如圖5所示。
3 算法分析
3.1 海浪理論
    根據水波理論，不規則的長峰波海浪可以看作是無數不同波幅、波長和初始相位的規則正弦波的疊加^[2]。假設t時刻波面某點距離平衡位置的高度為h(t)，即有：

圖5  SD卡接口電路

    垂直方向的加速度a(t)可以通過MEMS傳感器獲取，每個波的ai(t)和?棕i也可通過對a(t)進行頻域分析來得到。通過這種方法得到主波峰幅值，比直接用加速度二次積分得到計算量更小，更準確。
3.2 MEMS傳感器數據處理
    MEMS傳感器的浮標在水中不斷運動，加速度計的軸不能和水平面垂直方向時刻保持一致。要得到垂直方向的真實加速度值，需要對海浪浮標做姿態解算，求解出運動過程中浮標與海平面的縱搖、橫搖角度，并以此計算出所需的真實加速度值。
    四元數是簡單的超復數，由實數加上3個虛數單位i、j、k組成。由于i、j、k本身的幾何意義可以理解為對應空間坐標系x、y、z 3個軸的一種旋轉，所以四元數常用在描述物體的空間旋轉過程，且易于和傳統表示姿態的歐拉角相互轉換。只要能求解出當前的四元數，就可以計算出當前的姿態角。根據四元數的定義，當向量V1繞著軸n旋轉?茲角度到達V₂，可以構造一個四元數q(q₀，q₁，q₂，q₃)：
   

    根據四元數的微分方程，對任意時刻有：
   

    陀螺儀的作用是測量角速度信號，通過對角速度積分便能得到角度值。由于陀螺儀本身極易受噪聲干擾，并且會產生漂移誤差，通過積分產生較大的累積誤差。而加速度計不存在累積誤差，可以使用互補濾波對陀螺儀進行傾角補償[3，4]。互補濾波流程如圖6所示。

圖6  互補濾波流程

    先對加速度計和陀螺儀除去零漂，并對原始數據使用低通濾波器除去傳感器自身的機械和電子噪聲。假定測得t時刻的加速度計值為a_x、a_y、a_z，陀螺儀值為g_x、g_y、g_z。使用數據融合，完成加速度對陀螺儀進行補償。將加速度計近似看作只有重力加速度的作用，其四元數為a_g=(0，0，g)^T。
   

    如圖7，定義為繞z軸、y軸、x軸的旋轉角度，即艏搖、縱搖和橫搖。垂直方向的加速度也只與橫搖及縱搖有關，只需要求出角即可。

圖7  旋轉角的定義

    由系統此刻的四元數q=(q₀，q₁，q₂，q₃)，則可以將四元數轉換為歐拉角求出浮標的橫搖和縱搖：

3.3 時頻域分析
3.3.1 快速傅里葉變換
    海浪特征參數的變化是一個慢過程，在短時間內可以認為不變。由于海平面在短時間內不會發生明顯變化，故可以短時間內將海浪垂直方向的加速度看作是一個平穩信號，可以直接對其進行FFT分析，求解出主頻率。使用STM32F4系列單片機自帶的DSP庫中的FFT算法可以快速實現。
3.3.2 短時傅里葉變換
    從長時間的角度來看，海浪是一個時變系統，海浪的頻率和幅值在時刻發生變化。純粹地采用頻域分析方法完全將時域和頻域割裂開來，不能得到頻率與時間之間的關系，也不能捕捉到偶爾某時刻海浪突變的信號。要想描述信號的頻譜含量與時間之間的關系，就需要采用時頻域的方法來分析信號。
    根據垂直方向的加速度az(t)，進行頻域分析。需要用到時頻域分析的方法，這里使用短時傅里葉變換[5]。給定一個時間寬度很短的窗函數g(t)，沿著時間軸滑動，不斷去截取信號，并對截取信號作傅里葉變換[3]。則信號x(t)的STFT變換的定義為：
   

    根據式(15)可以確定如圖8所示的短時傅里葉變換的流程，容易在單片機系統內實現。

圖8  短時傅里葉變換流程

    短時傅里葉里存在不確定原理，時間分辨率和頻率分辨率不能同時達到最優。當窗函數確定后，其時間和頻率分辨率便確定。由于海浪變化的速度有限，在短時間內可以看作頻率不變，所以可以適當犧牲時域的分辨率以提高對頻率的分辨率。
    實際應用中將其離散化處理。以x(k)表示為信號x(t)的離散形式，則其離散變化的形式為：
   

    從式(16)來看，其變換結果為兩維的復數矩陣，對其求幅值可得到在給定時刻某一頻率的幅值大小。由于其本質還是傅里葉變換，在使用窗函數截取后，依舊可以使用FFT對信號進行快速處理。
4 實際驗證
4.1 數據采集并修正垂直方向的加速度
    完成硬件系統并制作成浮標，采用密封裝置保證了系統的防水性。在武漢東湖進行水試。為了方便陀螺儀的積分運算，將采樣頻率提高到100 Hz，并將原始運動數據存儲到SD卡中。
    根據采集的數據進行姿態融合，實時求出姿態角，并推算出垂直海平面的加速度，如圖9所示是修正后的一段加速度數據曲線，波浪相對較平穩。

圖9  時域曲線

4.2 FFT分析
    對圖9中采集到的數據作FFT分析，FFT分析結果如圖10所示?？梢缘玫斤@著主頻率大約為0.62 Hz，比較接近目測觀察值。

圖10  FFT分析

    根據主波的加速度變化頻率和其幅值，按照式(4)可以計算出主波的幅值約為9.88 cm。
4.3 短時傅里葉變換
    在FFT基礎上，對加速度值進行短時傅里葉變換。為突顯主瓣，直接使用矩形窗。取較長時間使用短時傅里葉變換來分析數據，得可如圖11所示的時頻圖。

圖11  時頻曲線

    從圖11中可以明顯地看到這100 s內海浪的頻譜主要能量集中在低頻部分，主頻率在0.6～1 Hz，也隨時間略微變化。在離散短時傅里葉變換生成矩陣中尋找每一時刻對應極值，可得到如圖12所示的頻率變化圖。

圖12  頻率變化曲線

5 結語
    本文借助STM32微控制器及MEMS傳感器，實現了對海浪頻率和幅值的監測，并通過GSM模塊將數據傳輸到云服務器端，完成了對海浪特征參數提取的目的。
    該設計具有低成本特點以及良好的實時性和準確性。其創新點在于將四元數姿態解析算法用于分析浮標的運動過程，并通過短時傅里葉的算法分析加速度信號，得到海浪頻率隨時間變化的趨勢。同時STM32單片機控制系統方便擴展硬件電路和接口，實現海浪發電系統的其他的控制算法的要求。實驗表明，本方案適用于海浪特性的測量，能達到較理想的效果。
參考文獻
[1] 蔡元奇.共振波力發電裝置：中國，CN201110233308.7[P].2011-08-16.
[2] 文圣常.海浪理論與計算原理[M].北京：科學出版社，1984.
[3] 張榮輝，賈宏光，陳濤，等．基于四元數法的捷聯式慣性導航系統的姿態解算[J]．光學精密工程，2008，16(10)：1963-1970.
[4] 周金金，王小英.基于MEMS的海洋浮標云臺穩定控制算法[J].制造業自動化，2015(4)：76-79.
[5] 趙鳳展，楊仁剛.基于短時傅里葉變換的電壓暫降擾動檢測[J].中國電機工程學報，2007(4)：28-34.