using Random
using Distributions
using Plots

# ===============================================================

function max_conv_operator(samples, f_samples, input, L) # max convolution operator, with random uniform sampling
    return maximum(f_samples .- L .* abs.(input .- samples))
end

function sharktooth(f, domain, number_of_sharkteeth, L, plot_arg)
    samples = rand(Uniform(domain[1], domain[2]), number_of_sharkteeth) #samples are uniformly taken from domain
    f_samples = f.(samples)        #y_i are f applied to samples componentwise
    x = LinRange(domain[1], domain[2], 10000) #approximant is plotted over a mesh of resolution 10000
    approximant = broadcast(x -> max_conv_operator(samples, f_samples, x, L), x) 
    error = maximum(abs.(f.(x) - approximant))
    if plot_arg == 1
        plot(x, approximant)
    else
        return error
    end
end


function sharktooth_error_plot(L_values, sharktooth_upper_limit, f, domain)
    # Initialize an empty plot
    p = plot(legend=false, xlabel="Number of Sharkteeth", ylabel="Error")

    for L in L_values
        errors = Float64[]  # Initialize an empty array to store errors
        for i in 1:sharktooth_upper_limit
            push!(errors, sharktooth(f, domain[1], domain[2], i, L, 0))
        end
        plot!(1:sharktooth_upper_limit, errors, label="L = $L")  # Add the plot to the existing figure
    end

    # Show the legend
    plot!(legend=true)

    # Display the final plot
    display(p)
end

sharktooth_error_plot([1, 2, 3, 4], 150, sin, [0, 30])

sharktooth(tanh, [-5, 5], 1500, 1, 1)