Biphasic Pulse signal with FFT

4 views (last 30 days)
Negin
Negin on 22 Feb 2024
Answered: Star Strider on 22 Feb 2024
I want to create a biphasic pulse signal as used in Functional electrical stimulation, afterwards I would like to do a FFT to get a better understanding of which frequencies would be good for the stimulation. However my code has some errors. I changed my Input variables and now it says:
"Error using plot
Vectors must be the same length.
Error in BiphasicFFT (line 13)
plot(t,wholeSignal);"
What does it mean and how can I change that? I would like to change my parameters, so that I can analyze different stimulation parameters.
Is my FFT correct?
This is my code:
fs=100000; %sampling freq in Hz
pulsesPerS=60; %frequency in Hz
signalDuration=3; % in seconds
amplitude=22*10^-3; % in mA
pulseWidth=250*10^-6; % in µs
t=0:1/fs:signalDuration-1/fs; %vector over the duration of time and sampled every 1/fs
onePeriod=amplitude*[ones(round(pulseWidth*fs),1); -ones(round(pulseWidth*fs),1); zeros(round(1/pulsesPerS*fs)-2*round(pulseWidth*fs),1)];
wholeSignal=repmat(onePeriod,[pulsesPerS*signalDuration 1]);
figure (1)
plot(t,wholeSignal);
Error using plot
Vectors must be the same length.
xlabel('Time (seconds)');
ylabel('Amplitude');
title('Biphasic Pulse Train')
fourier=abs(fft(wholeSignal));
figure (2)
plot(fourier);
xlabel('Samples');
ylabel('Magnitude');
title('Frequency-Domain plot - FFT');
N=length(wholeSignal); %length of signal
Freq=((0:N-1)/N)*fs; %frequency axis for FFT
Amp=abs(wholeSignal)/(N/2); %amplitude axis for FFT
figure (3)
plot(Freq,fourier);
%stem(Freq,fourier);
xlabel('Frequency in Hz');
ylabel('Amplitude');
title('FFT')

Accepted Answer

Star Strider
Star Strider on 22 Feb 2024
The problem was that ‘wholeSignal’ had 60 more elements than ‘t’ so I artificially shortened it:
wholeSignal = wholeSignal(1:numel(t));
so that they both fit.
After that, your code works. I would suggest doing a one-sided fft rather than two-sided, since it is a bit easier (more intuitive) to interpret. See the fft documentation for help with that.
Windowing the fft may also help, and you can do that with:
fourier=abs(fft(wholeSignal.*hann(numel(wholeSignal))));
That produces a slightly cleaner result. I added it (commented-out) so you can experiment with it.
fs=100000; %sampling freq in Hz
pulsesPerS=60; %frequency in Hz
signalDuration=3; % in seconds
amplitude=22*10^-3; % in mA
pulseWidth=250*10^-6; % in µs
t=0:1/fs:signalDuration-1/fs; %vector over the duration of time and sampled every 1/fs
onePeriod=amplitude*[ones(round(pulseWidth*fs),1); -ones(round(pulseWidth*fs),1); zeros(round(1/pulsesPerS*fs)-2*round(pulseWidth*fs),1)];
wholeSignal=repmat(onePeriod,[pulsesPerS*signalDuration 1]);
wholeSignal = wholeSignal(1:numel(t));
figure (1)
plot(t,wholeSignal);
xlabel('Time (seconds)');
ylabel('Amplitude');
title('Biphasic Pulse Train')
fourier=abs(fft(wholeSignal));
% fourier=abs(fft(wholeSignal.*hann(numel(wholeSignal))));
figure (2)
plot(fourier);
xlabel('Samples');
ylabel('Magnitude');
title('Frequency-Domain plot - FFT');
N=length(wholeSignal); %length of signal
Freq=((0:N-1)/N)*fs; %frequency axis for FFT
Amp=abs(wholeSignal)/(N/2); %amplitude axis for FFT
figure (3)
plot(Freq,fourier);
%stem(Freq,fourier);
xlabel('Frequency in Hz');
ylabel('Amplitude');
title('FFT')
Those were the only changes I made to your code.
Good luck with your FES research!
.

More Answers (0)

Categories

Find more on Fourier Analysis and Filtering in Help Center and File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!