Phase extraction from FFT
13 views (last 30 days)
Show older comments
Hello everyone,
I have written a code that takes as input voltage and current waveforms and is intended to have as output their frequency components and the power consumed. The thing is, by having a really high sampling frequency, tha data files must be very big. So I used zero-padding to fill in for the missing data. My comment is that I am getting wrong phase results and some of the magnitudes differ slightly from the real ones. What is really troubling though is the wrong phase results.
Below is my code:.
- This is a simplified version of the code which is intended to help me find my mistakes. The actual code gets input from data sheets.
Could you please inform me if something is wrong with my code or if I am missing something?
clear all;
clc;
F=13.56e6;
t=-0.5e-5:1/(2.5e9):0.5e-5;
V=-500+230*sin(2*pi*F*t-pi/4)+1*sin(2*pi*2*F*t-pi/6)+0.5*sin(2*pi*3*F*t-pi/3);
I=0.1*sin(2*pi*F*t+pi/8)+0.02*sin(2*pi*2*F*t+pi/10)+0.04*sin(2*pi*3*F*t);
signal1=V;
signal2=I;
N=length(t);
dt=-t(1,1)+t(1,2);
NFFT=5e5;
FS=1/dt;
Fn=FS/2; %nyquist frequency
Y1=fft(signal1,NFFT);
Y2=fft(signal2,NFFT);
F=linspace(0,1,(fix(NFFT/2)+1))*Fn; %freq vector
Iv=1:length(F);
HAmpl1=2*abs(Y1(Iv))/N;
HAmpl1(1,1)=HAmpl1(1,1)/2;
HAmpl2=2*abs(Y2(Iv))/N;
phase1= rad2deg(angle(Y1(Iv)));
phase2= rad2deg(angle(Y2(Iv)));
%% Figuring where the harmonics are in the matrix
C=zeros(1,length(F));
YY=zeros(1,length(Y1));
YYY=zeros(1,length(Y1));
Num=3; %number of harmonics calculated
M=max(HAmpl1(1,max(find(F<1.2e7 & F>1e7)):end));
noharm=ones(1,Num+1);
noharm(2)=find(HAmpl1==M);
temp=noharm(2);
for i=2:Num
noharm(i+1)=i*temp-(i-1);
end
%Values of interest
HAmplinterestV(1,:)=HAmpl1(noharm);
PhaseinterestV(1,:)=phase1(noharm);
HAmplinterestC(1,:)=HAmpl2(noharm);
PhaseinterestC(1,:)=phase2(noharm);
figure(1)
Ampl=subplot(2,1,1)
stem(F(noharm),HAmpl1(noharm))
Phasee=subplot(2,1,2)
stem(F(noharm),phase1(noharm))
5 Comments
Star Strider
on 21 Jan 2022
The waveform is (at the very least) not ‘pure’. Waveforms with multiple frequencies have multiple phase relationships. Most power calculations (that I am aware of) involve single frequencies. Dealing with waveform such as these are analytically difficult. Even though they involve relatively elementary trigonometry, interpreting the results can be a challenge.
% clear all;
% clc;
F=13.56e6;
t=-0.5e-5:1/(2.5e9):0.5e-5;
V=-500+230*sin(2*pi*F*t-pi/4)+1*sin(2*pi*2*F*t-pi/6)+0.5*sin(2*pi*3*F*t-pi/3);
I=0.1*sin(2*pi*F*t+pi/8)+0.02*sin(2*pi*2*F*t+pi/10)+0.04*sin(2*pi*3*F*t);
signal1=V;
signal2=I;
N=length(t);
dt=-t(1,1)+t(1,2);
NFFT=5e5;
FS=1/dt;
Fn=FS/2; %nyquist frequency
Y1=fft(signal1,NFFT);
Y2=fft(signal2,NFFT);
F=linspace(0,1,(fix(NFFT/2)+1))*Fn; %freq vector
Iv=1:length(F);
HAmpl1=2*abs(Y1(Iv))/N;
HAmpl1(1,1)=HAmpl1(1,1)/2;
HAmpl2=2*abs(Y2(Iv))/N;
phase1= rad2deg(angle(Y1(Iv)));
phase2= rad2deg(angle(Y2(Iv)));
%% Figuring where the harmonics are in the matrix
C=zeros(1,length(F));
YY=zeros(1,length(Y1));
YYY=zeros(1,length(Y1));
Num=3; %number of harmonics calculated
M=max(HAmpl1(1,max(find(F<1.2e7 & F>1e7)):end));
noharm=ones(1,Num+1);
noharm(2)=find(HAmpl1==M);
temp=noharm(2);
for i=2:Num
noharm(i+1)=i*temp-(i-1);
end
%Values of interest
HAmplinterestV(1,:)=HAmpl1(noharm);
PhaseinterestV(1,:)=phase1(noharm);
HAmplinterestC(1,:)=HAmpl2(noharm);
PhaseinterestC(1,:)=phase2(noharm);
figure(1)
Ampl=subplot(2,1,1)
stem(F(noharm),HAmpl1(noharm))
Phasee=subplot(2,1,2)
stem(F(noharm),phase1(noharm))
figure
yyaxis left
plot(t,V)
ylabel('V(t)')
yyaxis right
plot(t,I)
ylabel('I(t)')
grid
xlim([-1 1]*2.5E-7)
figure
plot(t, V.*I)
grid
title('Power')
xlim([-1 1]*2.5E-7)
.
Answers (0)
See Also
Categories
Find more on Spectral Measurements 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!