cutoff = 0.1;
nr = size(d,1);
mfstep = round(nr/100);
%mfstep = 1;
clear mfs
OPT = foptions; OPT(14) = 1000;

tic;
for i = 1 : NI,
  [ds,fs]=smooth(d(1:mfstep:nr,i),f(1:mfstep:nr,:),cutoff); figure(i); clf; plot(ds,fs); drawnow;
  mf = mffit(ds,fs,MFTYPE,OPT,[1 0 0]);
  lim = mf(:,3) == min(mf(:,3)); mf(lim,2:3) = ones(sum(lim),1)*[rlimits(i,1) rlimits(i,1)];
  if mf(lim,4) < rlimits(i,1), mf(lim,4) = rlimits(i,1); end;
  lim = mf(:,4) == max(mf(:,4)); mf(lim,4:5) = ones(sum(lim),1)*[rlimits(i,2) rlimits(i,2)];
  if mf(lim,3) > rlimits(i,2), mf(lim,3) = rlimits(i,2); end;
  mfplot(ds',mf); 
  mfs = [mfs; mf];
end;

mftime = toc
clear ds fs mf

rls = [ones(NC,1) reshape(1:NI*NC,NC,NI)];
