Effacer les filtres
Effacer les filtres

How to remove unwanted region in a graph

6 vues (au cours des 30 derniers jours)
okoth ochola
okoth ochola le 1 Mai 2024
Commenté : Mathieu NOE le 3 Mai 2024
Hi, I have this slight problem. consider the following attached code and the resultant graph also inserted as a picture. My problem is, the shoots on top of each wave is unwanted. I know they arise from points when sin(pi*t/T)=0. From my simulation, these are not important, since at these points are physically meaningless. Now how can I remove these values so the that the wave curves smoothly and peaks true maximum around the region circled(red). So that the curve smoothly peaks as illustrated with yellow marker. I thought I coluld use average around the vicinity wehere the peak appear and use it to replace the vauus that are shooting. How can I go about this?
clear
clc
r= 1e-7; %input('Enter the radius of each electron\n');
R=0.02; %input('Enter the the radius of the electron path\n');
T= 5; %input('Input the desire period\n');
t=0:0.01:30;
n=R/r;
h=6.626e-34;
me=9.109e-31;
acc1=900;
c=3e8;
f=3e15;
e=1.602e-19;
beta=1;
alpha=(h*f)/(me*c^2);
Keis=(me*c^2/(h*f))*(alpha-(alpha/(1+alpha*(1-cos(6*pi/6)))));
Z_i=beta.*sqrt(Keis);
A=Z_i*sqrt((h*f*e^2)/(R*acc1*me));
for i=1:1:numel(t)
m1=n*(abs(sin(pi*t(i)/T)));
m3=n*sqrt(abs(sin(pi*t(i)/T)));
m2=n*abs(cos(pi*t(i)/T));
if m1<2*cos(pi/6)
I(i,1)=-1000;
else
m=1:1:fix((m1)/(2*cos(pi/6)));
for j=1:1:numel(m)
A01=abs(sin(pi*t(i)/T));
A02=sqrt(A01);
I1(j,1)=A*(fix(m2/(2*A02))+fix((((m2*6)/(n*T)*180*(m(j))^2)/(A01^2*sqrt((n^2*A01^2)-(3*(m(j))^2))))+fix((m2*2*sin(acos(m(j)*sqrt(3)/(m1))))/(2*A02))));
if j==m(end);
mmax=m(end);
f10=m3*2*sin(acos((mmax*sqrt(3)/m1)));
f11=(fix(f10/2))/2;
f12=f11-fix(f11);
if f12==0
Na=1:1:((fix(f10/2)-2)/2);
for a=1:1:((fix(f10/2)-2)/2)
kon=2*180*R/(r*T);
kon1=4*180*r/(R*T);% define new constant
h=2*abs(csc(pi*t(i))/T)*Na(a)*r/R;
f13(a,1)=2*fix((2*(m2/n)*((kon*sin(acos(2*Na(a)*r/(m1*r))))+(kon1*Na(a)^2/(A01^3*sqrt(1-h^2))))));
%go to the next line
if a==numel(Na)
q=abs(csc(pi*t(i))/T)*r/R;
f14=2*fix(2*(m2/n)*((kon*sin(acos(r/(m1*r))))+(0.5*kon1*Na(a)^2/(A01^3*sqrt(1-q^2)))));
f15=A*fix((f14+sum(f13(1:a))));
I1(j,1)=f15;
end
end
else
Nb=1:1:((fix(f10/2)-1)/2);
for b=1:1:((fix(f10/2)-1)/2)
kon=2*180*R/(r*T);
kon1=4*180*r/(R*T);
h=2*abs(csc(pi*t(i))/T)*Nb(b)*r/R;
f16(b,1)=2*fix(2*(m2/n)*((kon*sin(acos(2*Nb(b)*r/(m1*r))))+(kon1*Nb(b)^2/(A01^3*sqrt(1-h^2)))));
if Nb==numel(Nb)
f17=fix(kon*(m2/n)*cos(pi*t(i)/T));
f18=A*fix(f17+sum(f16(1:b)));
I1(j,1)=f18;
end
end
end
f2=(sum(I1(1:j)));
I(i,1)=f2;
end
end
end
if i==numel(t)
for k=1:1:numel(t)
if I(k)==-1000
I(k)=max(I);%(I(2))*1;
if k==numel(t)
figure()
plot(t,I)
xlabel('time (s)')
ylabel('I(A)')
end
end
end
end
end

Réponse acceptée

Mathieu NOE
Mathieu NOE le 3 Mai 2024
Modifié(e) : Mathieu NOE le 3 Mai 2024
hello
see filloutliers , it was made for you
I changed a bit your code , I prefer to store the data (array tt and II) then plot it , instead of calling plot multiple times inside a for loop
makes also post processing easier
clear
clc
r= 1e-7; %input('Enter the radius of each electron\n');
R=0.02; %input('Enter the the radius of the electron path\n');
T= 5; %input('Input the desire period\n');
t=0:0.01:30;
n=R/r;
h=6.626e-34;
me=9.109e-31;
acc1=900;
c=3e8;
f=3e15;
e=1.602e-19;
beta=1;
alpha=(h*f)/(me*c^2);
Keis=(me*c^2/(h*f))*(alpha-(alpha/(1+alpha*(1-cos(6*pi/6)))));
Z_i=beta.*sqrt(Keis);
A=Z_i*sqrt((h*f*e^2)/(R*acc1*me));
tt = [];
II = [];
for i=1:1:numel(t)
m1=n*(abs(sin(pi*t(i)/T)));
m3=n*sqrt(abs(sin(pi*t(i)/T)));
m2=n*abs(cos(pi*t(i)/T));
if m1<2*cos(pi/6)
I(i,1)=-1000;
else
m=1:1:fix((m1)/(2*cos(pi/6)));
for j=1:1:numel(m)
A01=abs(sin(pi*t(i)/T));
A02=sqrt(A01);
I1(j,1)=A*(fix(m2/(2*A02))+fix((((m2*6)/(n*T)*180*(m(j))^2)/(A01^2*sqrt((n^2*A01^2)-(3*(m(j))^2))))+fix((m2*2*sin(acos(m(j)*sqrt(3)/(m1))))/(2*A02))));
if j==m(end);
mmax=m(end);
f10=m3*2*sin(acos((mmax*sqrt(3)/m1)));
f11=(fix(f10/2))/2;
f12=f11-fix(f11);
if f12==0
Na=1:1:((fix(f10/2)-2)/2);
for a=1:1:((fix(f10/2)-2)/2)
kon=2*180*R/(r*T);
kon1=4*180*r/(R*T);% define new constant
h=2*abs(csc(pi*t(i))/T)*Na(a)*r/R;
f13(a,1)=2*fix((2*(m2/n)*((kon*sin(acos(2*Na(a)*r/(m1*r))))+(kon1*Na(a)^2/(A01^3*sqrt(1-h^2))))));
%go to the next line
if a==numel(Na)
q=abs(csc(pi*t(i))/T)*r/R;
f14=2*fix(2*(m2/n)*((kon*sin(acos(r/(m1*r))))+(0.5*kon1*Na(a)^2/(A01^3*sqrt(1-q^2)))));
f15=A*fix((f14+sum(f13(1:a))));
I1(j,1)=f15;
end
end
else
Nb=1:1:((fix(f10/2)-1)/2);
for b=1:1:((fix(f10/2)-1)/2)
kon=2*180*R/(r*T);
kon1=4*180*r/(R*T);
h=2*abs(csc(pi*t(i))/T)*Nb(b)*r/R;
f16(b,1)=2*fix(2*(m2/n)*((kon*sin(acos(2*Nb(b)*r/(m1*r))))+(kon1*Nb(b)^2/(A01^3*sqrt(1-h^2)))));
if Nb==numel(Nb)
f17=fix(kon*(m2/n)*cos(pi*t(i)/T));
f18=A*fix(f17+sum(f16(1:b)));
I1(j,1)=f18;
end
end
end
f2=(sum(I1(1:j)));
I(i,1)=f2;
end
end
end
if i==numel(t)
for k=1:1:numel(t)
if I(k)==-1000
I(k)=max(I);%(I(2))*1;
if k==numel(t)
tt = [tt; t];
II = [II;I];
% figure()
% plot(t,I)
% xlabel('time (s)')
% ylabel('I(A)')
end
end
end
end
end
IIs = filloutliers(II,'linear','movmedian',25);
plot(tt,II,'b',tt,IIs,'r')
xlabel('time (s)')
ylabel('I(A)')
  1 commentaire
Mathieu NOE
Mathieu NOE le 3 Mai 2024
if you have an older matlab release , you can achieve the same result with interp1
load data.mat
% filloutliers METHOD
IIout = filloutliers(II,'linear','movmedian',25);
% OLD MATLAB METHOD (with interp1)
all_idx = 1:length(tt);
IIs = smoothdata(II,'movmedian',25);
outlier_idx = (abs(II./IIs) > 1.2); % Find outlier idx
IIold = II;
IIold(outlier_idx) = []; % remove outlier from dataset
IIold = interp1(all_idx(~outlier_idx), IIold, all_idx)'; % replace outliers by interpolated data
plot(tt,II,'b',tt,IIout,'*r',tt,IIold,'g')
xlabel('time (s)')
ylabel('I(A)')
legend('raw','filloutliers','interp1');

Connectez-vous pour commenter.

Plus de réponses (0)

Catégories

En savoir plus sur Interpolation dans Help Center et File Exchange

Produits


Version

R2018a

Community Treasure Hunt

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

Start Hunting!

Translated by