how to solve repeated y-axis problem

1 view (last 30 days)
common fernando
common fernando on 5 Dec 2020
Commented: common fernando on 6 Dec 2020
function [t_est,f_est]=threeDFT(v,fs,tmax,N0)
% v : volt as function of time
% fs : sampling frequency (Hz)
% tmax : time of final estimation
% N0 : number of samples in the window
% to test: [t,f]=threeDFT(@(t)(220*sin(2*pi*50.1*t+pi/2)),50*512,1,512)
% to test: [t,f]=threeDFT(@(t)(220*sin(2*pi*50.1*t+pi/2)+randn(size(t))*.1),50*512,1,512)
fs=200; %sampling freq.
dt =1/fs;
N0=fs/50; %number of samples/cycle
m=3; %no. of cycles
t = dt*(0:200); %data window
fi=50; %Frequency test
v = @(t)0.5*exp(-1.*t/0.3)+ 2*sin(2*pi*fi*t + pi/6);
figure()
plot(t,v(t))
tmax=1;
n=N0-1:-1:0;
f0=50;
f=50.88;
Hc=2/N0*cos(2*pi*n/N0+pi/N0);
Hs=-2/N0*sin(2*pi*n/N0+pi/N0);
t_est=[];
f_est=[];
j_max=tmax*fs;
for j=1:j_max+1
x=v((j-1:j+N0-2)*dt);
c(j)=x*Hc';
s(j)=x*Hs';
if(j>N0)
Ac(j-N0)=sqrt(sum(c(end-N0+1:end).^2)/N0);
As(j-N0)=sqrt(sum(s(end-N0+1:end).^2)/N0);
cc(j-N0)=c(end-N0+1:end)*Hc';
ss(j-N0)=c(end-N0+1:end)*Hs';
if(j>2*N0)
Acc(j-2*N0)=sqrt(sum(cc(end-N0+1:end).^2)/N0);
Ass(j-2*N0)=sqrt(sum(ss(end-N0+1:end).^2)/N0);
ccc(j-2*N0)=cc(end-N0+1:end)*Hc';
ccs(j-2*N0)=cc(end-N0+1:end)*Hs';
ssc(j-2*N0)=ss(end-N0+1:end)*Hc';
sss(j-2*N0)=ss(end-N0+1:end)*Hs';
ff=f0*N0/pi*atan(tan(pi/N0)*((ccc(j-2*N0).^2+ccs(j-2*N0).^2)./(ssc(j-2*N0).^2+sss(j-2*N0).^2)).^.25);
t_est=[t_est;(j-1)*dt];
f_est=[f_est;ff];
end
end
end
t_est;
f_est
plot(t_est, f_est,'red')
hold on
plot (t,fi)
hold off
  3 Comments
common fernando
common fernando on 5 Dec 2020
the error is in the above fig y-axis repeated value

Sign in to comment.

Answers (1)

common fernando
common fernando on 5 Dec 2020
the error in this fig y-axis
  2 Comments
Cris LaPierre
Cris LaPierre on 5 Dec 2020
I get something slightly different but equivalent.
fs=200; %sampling freq.
dt =1/fs;
N0=fs/50; %number of samples/cycle
m=3; %no. of cycles
t = dt*(0:200); %data window
fi=50; %Frequency test
v = @(t)0.5*exp(-1.*t/0.3)+ 2*sin(2*pi*fi*t + pi/6);
% figure()
% plot(t,v(t))
tmax=1;
n=N0-1:-1:0;
f0=50;
f=50.88;
Hc=2/N0*cos(2*pi*n/N0+pi/N0);
Hs=-2/N0*sin(2*pi*n/N0+pi/N0);
t_est=[];
f_est=[];
j_max=tmax*fs;
for j=1:j_max+1
x=v((j-1:j+N0-2)*dt);
c(j)=x*Hc';
s(j)=x*Hs';
if(j>N0)
Ac(j-N0)=sqrt(sum(c(end-N0+1:end).^2)/N0);
As(j-N0)=sqrt(sum(s(end-N0+1:end).^2)/N0);
cc(j-N0)=c(end-N0+1:end)*Hc';
ss(j-N0)=c(end-N0+1:end)*Hs';
if(j>2*N0)
Acc(j-2*N0)=sqrt(sum(cc(end-N0+1:end).^2)/N0);
Ass(j-2*N0)=sqrt(sum(ss(end-N0+1:end).^2)/N0);
ccc(j-2*N0)=cc(end-N0+1:end)*Hc';
ccs(j-2*N0)=cc(end-N0+1:end)*Hs';
ssc(j-2*N0)=ss(end-N0+1:end)*Hc';
sss(j-2*N0)=ss(end-N0+1:end)*Hs';
ff=f0*N0/pi*atan(tan(pi/N0)*((ccc(j-2*N0).^2+ccs(j-2*N0).^2)./(ssc(j-2*N0).^2+sss(j-2*N0).^2)).^.25);
t_est=[t_est;(j-1)*dt];
f_est=[f_est;ff];
end
end
end
t_est;
f_est
f_est = 193×1
50.0000 50.0000 50.0000 50.0000 50.0000 50.0000 50.0000 50.0000 50.0000 50.0000
plot(t_est, f_est,'red')
% hold on
% plot (t,fi)
% hold off
This is because your values of f_est are all essencially 50
format long
min(f_est)
ans =
49.999999635905169
max(f_est)
ans =
50.000000376435871

Sign in to comment.

Categories

Find more on Argument Definitions 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!