function [output]=quzao1(input,M1,lg1,M2,lg2)
%%%%%%%%%%%%%%%%%%%%%%%
N=length(input);
for b=1:1:M1;
   if lg1==1;
       tg1(b)=0.0001*(b-1);
      g1(b)=((0.0001*(M1-1)/2)^2-(tg1(b)-(0.0001*(M1-1)/2))^2)^0.5;
   else 
       g1(b)=1;
   end;
 end;
 for b=1:1:M2;
   if lg2==1;
       tg(b)=0.0001*(b-1);
      g2(b)=((0.0001*(M2-1)/2)^2-(tg(b)-(0.0001*(M2-1)/2))^2)^0.5;
   else 
       g2(b)=1;
   end;
 end;
%補充首尾數據
for n=1:1:N+400
     if n<201
        f1(n)=0;
        else if n>N+200
        f1(n)=0;
        else
       f1(n)=input(n-200);
            end;
     end;
end;
%開閉運算**********************************8k,b
%對此波形進行開運算去噪
N1=N+400;
for i=1:1:N1-M1+1;
    y0=f1(i:i+M1-1);
    for j=1:1:M1;
        kma(j)=y0(j)-g1(j);
    end;
    kma1(i)=min(kma);
    end;
for n=M1:1:N1-M1+1
    y1=kma1(n-M1+1:n);
   for m=1:1:M1
        kma2(m)=y1(m)+g1(m);
   end
    kma3(n)=max(kma2);
end;
%對jieguo進行閉運算去噪
for i=M2:1:N1-M1+1;
    y2=kma3(i-M2+1:i);
    for j=1:1:M2;
        bma(j)=y2(j)+g2(j);
    end;
    bma1(i)=max(bma);
end;
for n=1:1:N1-M1-M2+2
    y3=bma1(n:n+M2-1);
   for m=1:1:M2
        bma2(m)=y3(m)-g2(m);
   end
    bma3(n)=min(bma2);
end;
yb=bma3(1:N1-M1-M2+2);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%閉開運算d q
%對此波形進行閉運算去噪
for i=M1:1:N1;
    y4=f1(i-M1+1:i);
    for j=1:1:M1;
        dma(j)=y4(j)+g1(j);
    end;
    dma1(i)=max(dma);
end;
for n=1:1:N1-M1+1
    y5=dma1(n:n+M1-1);
   for m=1:1:M1
        dma2(m)=y5(m)-g1(m);
   end
    dma3(n)=min(dma2);
 end;
%對閉運算結果進行開運算去噪
for i=1:1:N1-M1-M2+2;
    y6=dma3(i:i+M2-1);
    for j=1:1:M2;
        qma(j)=y6(j)-g2(j);
    end;
    qma1(i)=min(qma);
    
end;
for n=M2:1:N1-M1-M2+2
    y7=qma1(n-M2+1:n);
   for m=1:1:M2
        qma2(m)=y7(m)+g2(m);
   end
    qma3(n)=max(qma2);
end;
yk=qma3(1:N1-M1-M2+2);
output1=0.5*(yk+yb);
output=output1(201:N+200);