I have created two complex exponential functions in MATLAB, both have same frequency but different phase. The MATLAB code used to create complex exponential functions is given below:
% First complex exponential generation
Fs = 5000; % sampling frequency
f1 = 5; % frequency
vector = linspace(0,1,Fs);
exp1 = exp(1i*(2*pi*f1.*vector));
% Second complex exponential generation
Fs = 5000;
f1 = 5;
vector = linspace(0,1,Fs);
exp2 = exp(1i*(2*pi*f1.*vector)).*exp(1i*pi/7); % Giving phase shift of pi/7
I want to calculate the phase difference between these two functions (which is known in this case i.e pi/7) so that in real situation, this phase shift can be used as a feedback to correct the phase difference and both functions can become exactly same.
Since you seem interested in finding the phase for a single frequency, I wrote the code to calculate the fft of each signal, then find the peak frequency, then calculate the phase difference at that frequency. The way you calculate the sine waves generates complex numbers, so I also added the code to use the sin
function to avoid generating complex numbers, but the code works with either real or complex numbers. I have it writing the actual phase and calculated phase on the screen so you can change the actual phase and easily verify that the phase is being calculated properly.
% First sine wave
Fs = 50; % sampling frequency
max_time = 20;
f1 = 5; % frequency
time = (0:(round(max_time * Fs) - 1)) / Fs;
%exp1 = exp(1i*(2*pi*f1*time));
exp1 = sin(2*pi*f1*time);
% Second sine wave
actual_phase_diff = pi/7
%exp2 = exp(1i*(2*pi*f1*time)).*exp(1i*actual_phase_diff); % original sine code
exp2 = sin(2*pi*f1*time + actual_phase_diff); % my sine code
fft_exp1 = fft(exp1);
fft_exp2 = fft(exp2);
if mod(length(fft_exp1), 2) > 0
% length is odd
max_freq_idx = (length(fft_exp1) + 1) / 2;
else
% length is even
max_freq_idx = length(fft_exp1) / 2 + 1;
end
% Find the point in the fft with the highest value
[~, peak_idx] = max(fft_exp1(1:max_freq_idx));
calc_phase_diff = angle(fft_exp2(peak_idx)) - angle(fft_exp1(peak_idx))