function [sx,fm] = smooth(x,f,cutoff)
% SMOOTH smooth a fuzzy set defined point-wise on the universe of
% discourse. Can be used e.g. for obtaining membership functions from
% multidimensional fuzzy clusters projected on the individual dimensions.
%
% [SX,FM] = SMOOTH(X,F,CUTOFF)
%
% Input:   X      ... universe of discourse (column vector)
%          F      ... fuzzy partition matrix N x C, N is number of data
%                     points, C is number of clusters
%          CUTOFF ... cut-off frequency of the 1-st order butterworth
%                     used for smoothing
% Outputs: SX     ... sorted universe of discourse
%          FM     ... smoothed fuzzy partition matrix

% (c) Robert Babuska, 1993-94

if nargin < 3, cutoff = 0.05; end;

[sx,xi] = sort(x);                              % sort X and F accordingly
f = f(xi,:);
[m,n] = size(f);

[b,a] = butter(1,cutoff);                       % design a Butterworth filter

for i = 1 : n,                                  % filter all columns of F
  fi1 = filter(b,a,f(:,i),f(1,i));              % top-bottom
  fi2 = filter(b,a,f(m:-1:1,i),f(m,i));         % bottom-top  fi(:,i) = normal((fi1 + fi2(m:-1:1)) / 2);    % average
  fi(:,i) = normal((fi1 + fi2(m:-1:1)) / 2);    % average
  fi(1,i) = f(1,i);  fi(m,i) = f(m,i);          % remove initial condition
end;

lf = zeros(1,n); rf = lf;
left = zeros(size(f)); right = left;
for i = 1 : m,                                  % peak detector, i.e.
  lf = max(lf,fi(i,:));                         % make F convex by
  rf = max(rf,fi(m-i+1,:));                     % taking its envelope
  right(m-i+1,:) = rf;
  left(i,:) = lf;
end;
fm = min(left,right);

if nargout == 0, plot(sx,f,'o',sx,fm); end;     % plot result

