Minimum working example (for version 0.1.1):
using EllipseSampling
import Distributions
using Plots
gr()
N=8
a, b = 2.0, 1.0
α = 0.25*π
Cx, Cy= 2.0, 2.0
points=generate_N_equally_spaced_points(N, a, b, α, Cx, Cy, start_point_shift=0.0)
plt=scatter(points[1,:], points[2,:], label=nothing, mw=0, ms=8, aspect_ratio=:equal)
Hw11 = (cos(α)^2 / a^2 + sin(α)^2 / b^2)
Hw22 = (sin(α)^2 / a^2 + cos(α)^2 / b^2)
Hw12 = cos(α)*sin(α)*(1/a^2 - 1/b^2)
Hw_norm = [Hw11 Hw12; Hw12 Hw22]
confidence_level=0.95
Hw = Hw_norm ./ (0.5 ./ (Distributions.quantile(Distributions.Chisq(2), confidence_level)*0.5))
Γ = inv(Hw)
points2 = generate_N_equally_spaced_points(N, Γ, [Cx, Cy], 1, 2, confidence_level=confidence_level, start_point_shift=0.0)
scatter!(points2[1, :], points2[2, :], label=nothing, markeralpha=0.5, mw=0, ms=4)
display(plt)
![precision_issue](https://private-user-images.githubusercontent.com/79883375/245961236-c2b344b1-cbfc-42c4-83f3-d19e31ba9c2c.png?jwt=eyJhbGciOiJIUzI1NiIsInR5cCI6IkpXVCJ9.eyJpc3MiOiJnaXRodWIuY29tIiwiYXVkIjoicmF3LmdpdGh1YnVzZXJjb250ZW50LmNvbSIsImtleSI6ImtleTUiLCJleHAiOjE3MjEzMTQzNzEsIm5iZiI6MTcyMTMxNDA3MSwicGF0aCI6Ii83OTg4MzM3NS8yNDU5NjEyMzYtYzJiMzQ0YjEtY2JmYy00MmM0LTgzZjMtZDE5ZTMxYmE5YzJjLnBuZz9YLUFtei1BbGdvcml0aG09QVdTNC1ITUFDLVNIQTI1NiZYLUFtei1DcmVkZW50aWFsPUFLSUFWQ09EWUxTQTUzUFFLNFpBJTJGMjAyNDA3MTglMkZ1cy1lYXN0LTElMkZzMyUyRmF3czRfcmVxdWVzdCZYLUFtei1EYXRlPTIwMjQwNzE4VDE0NDc1MVomWC1BbXotRXhwaXJlcz0zMDAmWC1BbXotU2lnbmF0dXJlPTJiNDU5NWExM2IyYjE1NjRlNWQ1ZGEzMzk4OThiYTgxMjM0YWFmMjFiMWI4OTU1NWE3MTNiNzI3NTg0ODhiZTYmWC1BbXotU2lnbmVkSGVhZGVycz1ob3N0JmFjdG9yX2lkPTAma2V5X2lkPTAmcmVwb19pZD0wIn0.vBh8pQNQKzzl2FPcUyHe5t0G7eylzvhwJRBkUPfl_20)
If we modify the rotation slightly away from 0.25pi the problem goes away:
a, b = 2.0, 1.0
α = 0.251*π
Cx, Cy= 2.0, 2.0
points=generate_N_equally_spaced_points(N, a, b, α, Cx, Cy, start_point_shift=0.0)
plt=scatter(points[1,:], points[2,:], label=nothing, mw=0, ms=8, aspect_ratio=:equal)
Hw11 = (cos(α)^2 / a^2 + sin(α)^2 / b^2)
Hw22 = (sin(α)^2 / a^2 + cos(α)^2 / b^2)
Hw12 = cos(α)*sin(α)*(1/a^2 - 1/b^2)
Hw_norm = [Hw11 Hw12; Hw12 Hw22]
confidence_level=0.95
Hw = Hw_norm ./ (0.5 ./ (Distributions.quantile(Distributions.Chisq(2), confidence_level)*0.5))
Γ = inv(Hw)
points2 = generate_N_equally_spaced_points(N, Γ, [Cx, Cy], 1, 2, confidence_level=confidence_level, start_point_shift=0.0)
scatter!(points2[1, :], points2[2, :], label=nothing, markeralpha=0.5, mw=0, ms=4)
display(plt)
![precision_issue_0 251pi](https://private-user-images.githubusercontent.com/79883375/245961544-5f9b76c7-8092-4adf-9711-9a78b31d1770.png?jwt=eyJhbGciOiJIUzI1NiIsInR5cCI6IkpXVCJ9.eyJpc3MiOiJnaXRodWIuY29tIiwiYXVkIjoicmF3LmdpdGh1YnVzZXJjb250ZW50LmNvbSIsImtleSI6ImtleTUiLCJleHAiOjE3MjEzMTQzNzEsIm5iZiI6MTcyMTMxNDA3MSwicGF0aCI6Ii83OTg4MzM3NS8yNDU5NjE1NDQtNWY5Yjc2YzctODA5Mi00YWRmLTk3MTEtOWE3OGIzMWQxNzcwLnBuZz9YLUFtei1BbGdvcml0aG09QVdTNC1ITUFDLVNIQTI1NiZYLUFtei1DcmVkZW50aWFsPUFLSUFWQ09EWUxTQTUzUFFLNFpBJTJGMjAyNDA3MTglMkZ1cy1lYXN0LTElMkZzMyUyRmF3czRfcmVxdWVzdCZYLUFtei1EYXRlPTIwMjQwNzE4VDE0NDc1MVomWC1BbXotRXhwaXJlcz0zMDAmWC1BbXotU2lnbmF0dXJlPTc5YzNiMTYyYmQwYWQ1MzcyYWU3ODQ1NThjN2NlYWJiMWY3MDU0MDhkOWI3NWVhY2U2YTBjNzU0MDViNzFiNjQmWC1BbXotU2lnbmVkSGVhZGVycz1ob3N0JmFjdG9yX2lkPTAma2V5X2lkPTAmcmVwb19pZD0wIn0.x0VX5VVMyh5CQyptHAGk9mVVc896945wmfcS3tm_cFc)
We can confirm this issue by inspecting the ellipse values returned by calculate_ellipse_parameters
:
# for α = 0.25*π
a_calc, b_calc, _, _, α_calc = EllipseSampling.calculate_ellipse_parameters(Γ, 1, 2, confidence_level)
println("a=", a_calc); println("b=", b_calc); println("α/pi=", α_calc/pi)
a=1.7994708216848752
b=1.0307764064044151
α=0.24999999999999992