function Fu = fray( u, t, rcfn, a, eps)
% Apply the fraying operator defined by a 
% rising cutoff function rcfn() to a signal
% function u(), at the point a, with reach eps.
%
% Inputs:
%     u     signal function
%     t     vector of time values
%     rcfn  rising cutoff function
%     a     point at which to fray
%     eps   reach of the fraying operator
%
% Output:
%     Fu    vector of Fu(t) for all times t.
%
% Example:
%     r=@(x) rcf(x,1); u=@(x)cos(x/2); t=-2:0.1:2; a=0; eps=1;
%     plot(t,fray(u,t,r,a,eps)); % plot the frayed signal
%     hold on; plot(t,u(t));    % overlay the original signal
%
% M.V.Wickerhauser, 2017-02-21
% Example added 2023-02-22
%
  Fu=t;   % initialize output array
  for(i=1:length(t))
    s=t(i);
    rplus = rcfn( (s-a)/eps );
    rminus = rcfn( (a-s)/eps );
    uplus = u(s);
    uminus = u(2*a-s);
    if(s>a)
      Fu(i)= uplus*rplus + uminus*rminus;
    else
      if(s<a)
        Fu(i)= uplus*rminus - uminus*rplus;
      else
        Fu(i)=u(a);
      endif
    endif
  end
end

