Biphasic Pulse signal with FFT

4 views (last 30 days)
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')

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')
Good luck with your FES research!
.