MATLAB and Simulink resources for Arduino, LEGO, and Raspberry Pi

Learn moreOpportunities for recent engineering grads.

Apply Today
Asked by li on 28 Jun 2013

I am trying to use Matlab to find the median frequency with respect to time of an EMG signal.

B = xlsread('C:/Users//Desktop/fatigue_20.xlsx'); Fe=800; A = B(0.1e4:8.2e5, :); N=length(A); t=0:1/Fe:(N/Fe)-1/Fe; psdest = psd(spectrum.periodogram,A_si,'Fs',800,'NFFT',N); normcumsumpsd = cumsum(psdest.Data)./sum(psdest.Data); Ind = find(normcumsumpsd <=0.5,1,'last'); Median frequency=psdest.Frequencies(Ind)

but now how can i build the relation between the median frequency and the time?

*No products are associated with this question.*

Answer by Wayne King on 28 Jun 2013

Accepted answer

If you want to find the median frequency with respect to time, you'll have to use the short-time Fourier transform (or some other time-frequency analysis technique)

Using the periodogram, you lose all time information.

You can use spectrogram() to obtain a matrix of PSD estimates with some degree of time localization that depends on how long the window (segment) is. You can then compute the median of your segments.

Answer by Wayne King on 28 Jun 2013

Edited by Wayne King on 28 Jun 2013

Here you go:

Fs = 1000; t = 0:1/Fs:5; x = 2*cos(2*pi*100*t).*(t<1)+3*cos(2*pi*200*t).*(t>1 & t<2)+1.5*cos(2*pi*150*t).*(t>2 & t<3)+2*sin(2*pi*50*t-pi/4).*(t>4)+randn(size(t)); window = 100; [S,F,T,P] = spectrogram(x,window,window/2,window,Fs); for nn = 1:size(P,2) normcumsumpsd = cumsum(P(:,nn))./sum(P(:,nn)); Ind = find(normcumsumpsd <=0.5,1,'last'); medianfreqs(nn) = F(Ind); end plot(T,medianfreqs); xlabel('Time (seconds)'); ylabel('Median Frequency (Hz)');

You see that the median frequency is time dependent

li on 28 Jun 2013

yes,thank you,if i remplace the x with a EMG signal, what should i change?beacause i have met a question : Improper assignment with rectangular empty matrix.

Error in Untitled4 (line 9) medianfreqs(nn) = F(Ind);

Answer by Wayne King on 28 Jun 2013

If I just substitute a vector of random noise the same size as your input

Fs = 800; x = randn(841700,1); window = 100; [S,F,T,P] = spectrogram(x,window,window/2,window,Fs); for nn = 1:size(P,2) normcumsumpsd = cumsum(P(:,nn))./sum(P(:,nn)); Ind = find(normcumsumpsd <=0.5,1,'last'); medianfreqs(nn) = F(Ind); end plot(T,medianfreqs); xlabel('Time (seconds)'); ylabel('Median Frequency (Hz)');

The above works fine for me. Of course, you need to make sure the window length makes sense for your data.

li on 28 Jun 2013

you are right,but i haven't succeeded yet,it still gives the following error message

>> Untitled666 Improper assignment with rectangular empty matrix.

Error in Untitled666 (line 8) medianfreqs(nn) = F(Ind);

our matrix looks like : 0.2250 0.3220 0.2250 0.0150 -0.0490 0.0730 0.0830 0.0340 0.0200 -0.0340 0.0630 0.0490 0.0340 -0.0630 -0.0590 0.0490 0.2050 0.2690 0.3120 0.0830 -0.0590 -0.0440 it's a bit of a mistery.

## 0 Comments