include("even_sampling.jl")

"""
Returns a matrix of convex coefficients to construct points contained within
a `dim`-dimensional simplex. Coefficients can either be random (generated
according to a uniform distribution) or be constructed such that the resulting
points are uniformly distributed within the simplex (generated by shape-preserving
splitting a generic simplex). The `sample_randomly` argument controls this
behaviour.

The default is to not sample randomly (`sample_randomly = false`). In this case,
you won't get the exact number of convex coefficient sets you want, but the
smallest number of set of coefficients such that you have *at least* `n_pts`
convex coefficient sets.

If `sample_randomly = true`, you get as many convex coefficient combinations as
`n_pts`.
"""
function subsample_coeffs(dim::Int, n_randpts::Int, sample_randomly::Bool)
    # Create a set of convex coefficients to be used for all simplices
    if sample_randomly
        convex_coeffs = rand(Distributions.Uniform(0, 1), n_randpts, dim + 1)
        convex_coeffs .= convex_coeffs ./ sum(convex_coeffs, dims=2)
    else
        minimum_split_factor = ceil(Int, n_randpts^(1 / dim))[1]
        convex_coeffs = copy(transpose(even_sampling_rules(dim, minimum_split_factor)))
    end
     copy(transpose(convex_coeffs))
end


"""
    childpoint(parentsimplex::AbstractArray{Float64, 2})

Generates a random simplex which is entirely contained within `parentsimplex`,
which is a (dim+1)-by-dim array.
"""
function childpoint(parentsimplex::AbstractArray{Float64, 2})
    dim = size(parentsimplex, 2)
    # Random linear combination coefficients
    R = rand(Float64, 1, dim + 1)

    # Normalise the coefficients so that they sum to one. We can then create
	# the new point as a convex linear combination of the vertices of the parent
	# simplex.
    R = (1 ./ sum(R, dims=2)) .* R

    R * parentsimplex
end

"""
    childpoint(parentsimplex::AbstractArray{Float64, 2})

Generates a random simplex which is entirely contained within `parentsimplex`,
which is a (dim+1)-by-dim array.
"""
function childpoint(parentsimplex::AbstractArray{Float64, 2}, dim::Int)
    # Random linear combination coefficients
    R = rand(Float64, 1, dim + 1)

    # Normalise the coefficients so that they sum to one. We can then create
	# the new point as a convex linear combination of the vertices of the parent
	# simplex.
    R = (1 ./ sum(R, dims=2)) .* R

    R * parentsimplex
end


"""
Create `n` (default: n = 100) points that lies inside `parent_simplex`,
which is a (dim+1)-by-dim array.
"""
function childpoints(parent_simplex::AbstractArray{Float64, 2}, n::Int)
    dim = size(parent_simplex, 2)

    # Random linear combination coefficients
    R = rand(Float64, n, dim + 1)

    # Normalise the coefficients so that they sum to one. We can then create
    # the new point
    # as a convex linear combination of the vertices of the parent simplex.
    normalised_coeffs = (1 ./ sum(R, dims=2)) .* R

    normalised_coeffs * parent_simplex
end


"""
Create a point that lies inside `parent_simplex`, which is a (dim+1)-by-dim
array.
"""
function outsidepoint(parentsimplex::AbstractArray{Float64, 2})
    dim = size(parentsimplex, 2)
    # Random linear combination coefficients
    R = rand(Float64, 1, dim + 1)

    # Normalise the coefficients so that they sum to one. We can then create the new point
    # as a convex linear combination of the vertices of the parent simplex.
    normalised_coeffs = (1 ./ sum(R, dims=2)) .* R
    normalised_coeffs[1] += 1

    normalised_coeffs * parentsimplex
end
