MATLAB: Total energy in Time and Frequency Dimain

time domaintotal energy

I truly appreciated if somebody can direct me how to measure total energy in a certain frequency band?

I have a signal and use this code to measure total enery and the result is reasonable.

E1_timedomain=sum(abs(x.^2))

but it measuers total energy of whole signal. I just want to measure it between [0-100] Hz. Appreciared!

Best Answer

If you have R2018a or later, and the Signal Processing Toolbox, use the bandpass function to selectively filter the frequencies-of-interest. Then do whatever calculations you want to do on the output. If you have an earlier version, it is easy to design filters in MATLAB. I will help you with that, if necessary.

Related Solutions

MATLAB: How to Butterworth Filter with Bandpass [10 500] with sampling rate 1000

You have to divide your passband and stopband frequencies by the Nyquist frequency to normalise them:

Fs = 1000;                                          % Sampling Frequency
Fn = Fs/2;                                          % Nyquist Frequency
Wp = [10  495]/Fn;                                  % Norrmalised Passband Frequencies
Ws = [5   499]/Fn;                                  % Normalised Stopband Frequencies
Rp =  1;                                            % Passband Ripple (dB)
Rs = 25;                                            % Stopband Ripple (dB)
[n,Wn] = buttord(Wp, Ws, Rp, Rs);                   % Optimal Filter Order
[b,a] = butter(n, Wn);                              % Calculate Filter Coefficients
[sos,g] = tf2sos(b,a);                              % Convert To Second-Order Sections For Stability
figure(1)
freqz(sos, 4096, Fs)                                % Filter Bode Plot

You cannot have any frequency exactly at zero or exactly at the Nyquist frequency, so I designed your passband to be close enough to them to be stable. In your filter, you can also use a high-pass design with a lower cutoff of 10 Hz, since your upper passband is at the Nyquist frequency. The passband and stopband ripple values are necessary to specify in the design but irrelevant to a Butterworth filter. These are acceptable values for most filter designs (for example Chebyshev) that require them.

You would use the ‘sos’ and ‘g’ variables in your actual filter:

filtered_signal = filtfilt(sos, g, signal);

(This design works well in R2016a and should also work in earlier releases.)

MATLAB: How to transform time series into frequency dimain and obtain amplitude and phase

There's really nothing wrong with the code; you have followed the example from the documentation for FFT and the following lines

    P_amp2=abs(Y/L);
    P_amp(:,jj) = P_amp2(1:L/2+1);
    P_amp(2:end-1,jj) = 2*P_amp(2:end-1,jj);
    f=Fs*(0:(L/2))/L;

return the positive frequency half of the estimated PSD and normalized it to match total power in the double-sided PSD (P_amp2) by the factor of 2 applied to all frequencies except DC and Fmax which are not replicated in the two-sided spectrum.

The frequency vector f is for the one-sided PSD estimate.

If you really want the two-sided PSD, then don't take the positive half elements in P_amp; just use P_amp2 and the corresponding double-sided frequency vector. At that point, for convenience since apparently wouldn't need them, could rename P_amp2 to P_amp and fs to f.

BTW, the FFT and other MATLAB routines are vectorized, you can compute all the results using the original array and dispense with the for loop...

nn_ss=data(1:1800,1:88);           % if size(data)==>1800,88, then this is not needed; use the data array
Ts = mean(diff(time));             % sampling time (by column)
Fs=500;                            %sampling frequency    ** It should be 1/Ts
Fn=Fs/2;                           %Nyquist frequency
L=length(time);
Y=fft(nn_ss);                      % FFT() by column
fts=Y/L;
P_amp=abs(Y/L);                    % 2-sided PSD estimate
phase=angle(Y).*180/pi;
ly=length(Y);
fs=(-ly/2:ly/2-1)/ly*Fs;

The result will be NxM PSD (2-sided) by column same as if had used for...end loop

Best Answer

Related Solutions

MATLAB: How to Butterworth Filter with Bandpass [10 500] with sampling rate 1000

MATLAB: How to transform time series into frequency dimain and obtain amplitude and phase

Related Question