function s = discrete_cdf_to_xy ( m1, m2, cdf_mat, xb, yb, n, seed, u )

%*****************************************************************************80
%
%% DISCRETE_CDF_TO_XY finds sample points S corresponding to discrete CDF values.
%
%  Discussion:
%
%    This program is given a discrete CDF function and a set of N random
%    values U.  Each value of U corresponds to a particular (I,J) subregion
%    whose CDF value just exceeds the value of U.  Inside that subregion,
%    we pick a point at random - this is equivalent to assuming the PDF
%    is constant over the subregion.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    23 August 2016
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, integer M1, M2, the number of rows and columns
%    of PDF data.
%
%    Input, real CDF_MAT(M1,M2), the CDF values associated with each subcell.
%    A particular ordering has been given to the subcells so that the
%    CDF is a monotonoe function when the subcells are listed in that order.
%
%    Input, real XB(M1+1), YB(M2+1), the break points.
%
%    Input, integer N, the number of sample points.
%
%    Input, real U(N), N random values.
%
%    Input/output, integer SEED, a seed for the random number generator.
%
%    Output, real S(2,N), the sample points.
%
  s(1:2,1:n) = 0.0;

  low = 0.0;
  for j = 1 : m2
    for i = 1 : m1
      high = cdf_mat(i,j);

      for k = 1 : n
        if ( low <= u(k) && u(k) <= high )
          [ r, seed ] = r8vec_uniform_01 ( 2, seed );
          s(1,k) = ( 1.0 - r(1) ) * xb(i) ... 
                 +         r(1)   * xb(i+1);
          s(2,k) = ( 1.0 - r(2) ) * yb(j) ...
                 +         r(2)   * yb(j+1);
        end
      end

      low = high;
    end
  end

  return
end
