How to remove unwanted region in a graph

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

Answer 1

[Mathieu NOE]

hello

see filloutliers , it was made for you

Detect and replace outliers in data - MATLAB filloutliers - MathWorks France

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)')

[Mathieu NOE]

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');