%%%%%%%%%% crmz_guess.m %%%%%%%%%% 
function  [fext,iext] = crmz_guess( edges, grid, nfcns )
% REMGRID   generate initial guess of "EXTREMAL FREQS"
%            for remez exchange algorithm [jMc 11/92]
%   usage:
%        [Fext,Iext] = rmz_guess( Edges, Grid, Nfcns )
%
%     Edges:  band-edges moved onto the grid (Revised edges)
%     Grid:   frequencies on the fine grid
%     Nfcns:  number of basis functions in the approx
%     Fext:   initial guess of "extremal frequencies"
%     Iext:   indices for the "extremal frequencies"
%
TOL = 5*eps;
next = nfcns+1;
Nbands = length(edges(:,2));
tt = edges;
merged = tt;
if( Nbands > 1 )
  jkl = find( abs(tt(1:(Nbands-1),2)-tt(2:Nbands,1)) > TOL );
  merged = [ [tt(1,1); tt(jkl+1,1)] , [tt(jkl,2); tt(Nbands,2)] ];
end
Nbands = length(merged(:,1));
bw = merged(:,2)-merged(:,1);
if( any(bw<0) )
   edges
   error('REM_GUESS: negative bandwidth!')
end
percent_bw = bw / sum(bw);
fext = zeros(next,1);
n = 0;  i = 1;
while( n<next & i<=Nbands )
   nfreqs_i = min( next-n, ceil( percent_bw(i)*next ) );
   if( nfreqs_i == 0 )
      n=n+1;   fext(n) = merged(i,1);
   else
      fext(n+[1:nfreqs_i]) = ...
           linspace(merged(i,1),merged(i,2),nfreqs_i);
      n = n+nfreqs_i;
   end
   i = i+1;
end
iext = 0*fext;
for i=1:next
   [tt,iext(i)] = min( abs(fext(i)-grid) );
end
if( any(diff(iext)==0) )
   error('RMZ_GUESS: two extremal freqs at same grid pt, NOT!')
end
fext = grid(iext);
%%%%%%%%%% End crmz_guess.m %%%%%%%%%% 
