Labview and Matlab for PT

Just a place to discuss with Matlab and Labview

ad

EMG頻域(frequency domain)之中位頻率(Median frequency)計算


  當我們將肌電訊號轉換成頻譜(Power spectrum)後,我們的確可
以看到頻率和其強度的分佈,不過這代表什麼意義呢?

  我們知道我們的肌肉由不同型態組成,簡單的說可以分為快肌
(fast twitch)跟慢肌(slow twitch),快肌肌肉纖維較粗,動作速
度快,不過較容易疲勞,慢肌動作速度慢,但是可以作用持續時間長,
就傳遞肌肉表面膜電位而言,快肌就像比較粗的電線,傳遞訊號的速度
較快,因此由電極片所接收到單位時間內的電位差訊號相對而言就會比
較多,所以頻率會落在比較高的位置,反觀慢肌,則為相反地現象。



  而我們的肌肉,則由不同比例的快慢肌所組成,例如我們強而有力
的肱二頭肌(biceps),快肌的組成較多,所以也比較容易酸。

  不過因為頻譜看出來的是各個頻率強度的變化,如果要相互比較的
話,則必須要像一個常態分佈圖形一樣,有平均數(mean),等數值,
比較常用來作為計算的為中位頻率(Median frequency)。
 
  中位頻率(Median frequency),其定義為頻譜分析圖中,將其強
度由低頻積分到高頻,積分面積為全部面積一半時的頻率。因此就可以
中位頻率來作為訓練或是疲勞測試後,改變肌肉組成比例(訓練)或是
影響肌肉放電(疲勞)等改變觀察。



  因此利用傅立葉轉換以及if指令來求得中位頻率,當積分到某一點
時,大於全部等積分的一半時的頻率即是中位頻率。

  由於肌電訊號是時間性的,我們比較常做的方式會像之前所用的以
window的方式,每多少點計算一次,重複多少點,然而取樣頻率
(sampling rate)會影響最後每一個點所代表的頻率,所以也必須考慮



function [med_fq]=fft_i(data,n,op,sr)
% fft for n point
% op: the point of overlapping
% med_fq is the median frequency of the power spectrum
% sr is the sampling frequency

% calculate the times of the loop (-1 due to the initial loop is 0)
l=fix(length(data)/op-1);
fq=(1:n/2).*sr/n; % calculate the frequency of each point
for i = 0:(l-1)
int=1+op*i; % the initial value of the loop
f=fft(data(int:(int+511))); % fft for the data
mag=abs(f); % Amplitude of the fft
mag=mag(1:n/2);
for k=1:n/2
if sum(mag(1:k)) >= sum(mag)/2
med=k;
break
end
end
med_fq(i+1,1)=k*sr/n;
end

% calculate the final part of the sample
int=1+op*l;
f=fft(data(int:end)); % fft for the data
mag=abs(f); % Amplitude of the fft
mag=mag(1:n/2);
for k=1:n/2
if sum(mag(1:k)) >= sum(mag)/2
med=k;
break
end
end
med_fq(l+1,1)=k*sr/n;


  首先必須先計算迴圈的次數,以所有的點數除上重複的點數(通常
會以window點數的一半)減一(種樹法則),然後就每多少點作一次傅
立葉轉換並且求得中位頻率,值得注意的是,因為常常所選取的資料數
並不能剛好夠最後一筆,舉例而言,也就是說原本應該要以512個點作計
算,但是迴圈到最後一次時可能只剩下500點,由於fft這個函數當不夠
點的時候會比剩下的點以零帶入,所以最後一個窗口則以另外方式計算


執行結果
>> data=textread('D:/EMG_fft.txt');
>> plot(data);
>> med_fq=fft_i(data,512,256,1000);
>> figure(2);plot(med_fq)



這是原始資料圖,資料可以由此下載

經由計算後為


一開始很高可能是因為在還沒有任何訊號時就已經存在的高頻雜訊,
此時就要利用之前提得濾波器給予濾波。

0 意見: