function xy = Isomap(D,ndims,options); 

% isomap - computes Isomap embedding using the algorithm of 
%             Tenenbaum, de Silva, and Langford (2000). 
%
%   xy = Isomap(X,ndims,options); 
%       or 
%   xy = Isomap(D,ndims,options); 
%
%   'X' is a D x N matrix (D = dimensionality, N = #points)
%   'D' is either a local distance matrix (with Inf when no connection),
%       or a global one (i.e. it contains already computed geodesic
%       distances between pair of points).
%   OPTIONAL:
%   'ndims' is the number of output dimensions (default =2).
%   'options.distance_mode' is either 'local' (set if D is a local distance matrix)
%       or 'global' (set if D is a geodesic distance).
%   'options.nn_epsilon' : set it if you want to compute the nearest
%       neighbor by thresholding.
%   'options.nn_nbr' : set it if you want to compute the nearest
%       neighbor by using a fixed number of neighbors.
%       
%   Modified by Gabriel Peyr