1. 程式人生 > >[Matlab] 濾波器filter函數造輪子及使用代碼生成進行速度優化

[Matlab] 濾波器filter函數造輪子及使用代碼生成進行速度優化

之前 環境 支持 測試數據 。。 內置 說明 cell isp

之前做腦機接口上位機的時候需要對數據進行實時濾波,也就是需要對數據進行分段濾波,保存濾波器前一次的歷史狀態。翻了翻MATLAB官方文檔的filter函數發現有這個功能,不過他們的函數說明是用相位及延遲進行設置,看了半天沒理解什麽意思,為了保險起見我自己造了個輪子實現簡單的IIR濾波器。filter函數的官方文檔:https://ww2.mathworks.cn/help/matlab/ref/filter.html?s_tid=doc_ta 在官方文檔中IIR濾波器的詳細原理都已經講的很明白了不再煞述。自己造的輪子也很簡單,關於濾波器的參數都被保存在myIIR的結構體中,包括濾波器系數a,b,歷史數據x,y,以及濾波後的數據FilteredData。在實時檢測系統中一般對速度由較強的要求,所以我們可以將MATLAB代碼轉成mex文件進行提速。不是所有MATLAB代碼都是可以轉成mex或C的,總結來說如果要想順利轉成C的話還是得用C的思想寫MATLAB代碼,MATLAB那些奇技淫巧能少用就少用,要用某些函數的手先到官方文檔上查一下有沒有相應的C/C++ Coder的支持。這樣這個簡單的濾波函數代碼如下:

function myIIR = myFilter(data,myIIR)
    for iSample = 1:size(data,2)
        myIIR.x = circshift(myIIR.x,1,2);
        myIIR.x(:,1) = data(:,iSample);
%         myIIR.x = [data(:,iSample),myIIR.x];    % PUSH
%         myIIR.x(:,end) = [];    % POP
        y = sum(myIIR.b.*myIIR.x,2)-sum(myIIR.a.*myIIR.y,2);
        myIIR.y = circshift(myIIR.y,1,2);
        myIIR.y(:,1) = y;
%         myIIR.y = [y,myIIR.y];  % PUSH
%         myIIR.y(:,end) = [];    % POP
        myIIR.FilteredData(:,iSample) = y;
    end
end

代碼註釋掉的部分就是代碼生成不支持的語句,用循環位移函數circshift函數代替了。為了對比,內置函數也寫成這種結構體的形式,代碼如下:

function IIR = Filter(data,IIR)
    [IIR.FilteredData,IIR.zf] = filter(IIR.b,IIR.a,data,IIR.zf,2);
end

進行這樣設置以後,可以進行mex文件生成,生成後在本機上進行測試。測試數據是一個cell元胞數組,包含10000個8x256矩陣,模擬10000次8通道濾波。本機測試環境:

技術分享圖片

程序運行的MATLAB版本為MATLAB-2018b,我分別在MATLAB-2016a, -2018a,-2018b三個不同版本的MATLAB上生成了3種mex文件對比不同版本MATLAB的差異,結果如下:

技術分享圖片

技術分享圖片

技術分享圖片

可以看出來自己用MATLAB造的filter函數的輪子遠遠沒有內置函數的快,但轉為mex文件以後可以提高至少10倍的速度。而內置函數一直是最快的,甚至生成mex文件後速度還有一定的下降。不同版本之間生成的代碼也是有差異的,2016a中生成的輪子的mex明顯比2018a和2018b中生成的要快,但內置函數轉為mex後反而慢了。。這個現象也是挺奇葩

綜上,我感覺除非真的對速度有較為苛刻的要求,否則還是用MATLAB語法裏的奇技淫巧寫代碼效率比較高,若要轉為mex的話需要參考不少C 的語法,那MATLAB代碼寫著就沒啥意思了

[Matlab] 濾波器filter函數造輪子及使用代碼生成進行速度優化