function [ generator_new, change_l2, seed ] = cvt_iteration ( m, n, ...
  generator, sample_num_cvt, sample_function_cvt, seed )

%*****************************************************************************80
%
%%  CVT_ITERATION takes one step of the CVT iteration.
%
%  Discussion:
%
%    The routine is given a set of points, called "generators", which
%    define a tessellation of the region into Voronoi cells.  Each point
%    defines a cell.  Each cell, in turn, has a centroid, but it is
%    unlikely that the centroid and the generator coincide.
%
%    Each time this CVT iteration is carried out, an attempt is made
%    to modify the generators in such a way that they are closer and
%    closer to being the centroids of the Voronoi cells they generate.
%
%    A large number of sample points are generated, and the nearest generator
%    is determined.  A count is kept of how many points were nearest to each
%    generator.  Once the sampling is completed, the location of all the
%    generators is adjusted.  This step should decrease the discrepancy
%    between the generators and the centroids.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    10 February 2015
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, integer M, the spatial dimension.
%
%    Input, integer N, the number of Voronoi cells.
%
%    Input, real GENERATOR(M,N), the Voronoi cell generators.
%
%    Input, integer SAMPLE_NUM_CVT, the number of sample points.
%
%    Input, integer SAMPLE_FUNCTION_CVT, region sampling function:
%    -1, sampling function is RAND (MATLAB intrinsic);
%    0, sampling function is UNIFORM;
%    1, sampling function is HALTON;
%    2, sampling function is GRID;
%
%    Input, integer SEED, the random number seed.
%
%    Output, real GENERATOR_NEW(M,N), the new Voronoi cell generators.  
%
%    Output, integer SEED, the new random number seed.
%
%    Output, real CHANGE_L2, the L2 norm of the difference between
%    the input and output data.
%
  generator_new(1:m,1:n) = 0.0;
  tally(1:n) = 0;
  reset = 1;

  for j = 1 : sample_num_cvt
%
%  Generate a sampling point X.
%
    [ x(1:m), seed ] = region_sampler ( m, 1, sample_num_cvt, ...
       sample_function_cvt, reset, seed );
%
%  Ensure that X is a column vector.
%
    x = x(:);

    reset = 0;
%
%  Find the nearest cell generator.
%
    nearest = find_closest ( m, n, 1, x, generator );
%
%  Add X to the averaging data for GENERATOR(*,NEAREST).
%
    for i = 1 : m
      generator_new(i,nearest) = generator_new(i,nearest) + x(i);
    end
    tally(nearest) = tally(nearest) + 1;

  end
%
%  Compute the new generators.
%
  for j = 1 : n
    if ( tally(j) ~= 0 )
      generator_new(1:m,j) = generator_new(1:m,j) / tally(j);
    end
  end
%
%  Determine the change.
%
  change_l2 = 0.0;
  for j = 1 : n
    for i = 1 : m
      change_l2 = change_l2 + ( generator_new(i,j) - generator(i,j) )^2;
    end
  end
  change_l2 = sqrt ( change_l2 );

  return
end
