From 1b99810b042b1fea2a4b64a2616e8fda4c940f94 Mon Sep 17 00:00:00 2001 From: Greg Burri Date: Thu, 8 Jul 2021 11:41:10 +0200 Subject: [PATCH] Initial commit. --- eth_staking.jl | 61 +++++ fourier.jl | 661 +++++++++++++++++++++++++++++++++++++++++++++++++ 2 files changed, 722 insertions(+) create mode 100644 eth_staking.jl create mode 100644 fourier.jl diff --git a/eth_staking.jl b/eth_staking.jl new file mode 100644 index 0000000..7a42f56 --- /dev/null +++ b/eth_staking.jl @@ -0,0 +1,61 @@ +### A Pluto.jl notebook ### +# v0.14.8 + +using Markdown +using InteractiveUtils + +# ╔═╡ 701a0e50-df22-11eb-2eea-a56013a4ffdd +html""" +""" + +# ╔═╡ 6605e130-e014-4de3-8837-baf4eb77418c +md""" +**Ethereum staking project** + +The goal is to set up a passive income for the beginning of year 2025. +""" + +# ╔═╡ 9a96995c-89be-4fa9-97f6-81b7bd5898b8 +price_eth_today = 2380 + +# ╔═╡ ff389da3-6761-4d78-873e-94f1b5319a05 +price_btc_today = 34_800 + +# ╔═╡ 52becaeb-d325-4b17-a4ba-de19d6eb2cf8 +yield = 0.05 # 5% per year. + +# ╔═╡ 5bd92f00-0b9d-4217-92dc-ed3f7234666f +usd_per_month_to_live = 9_000 + +# ╔═╡ 3a5fdd52-6445-4b10-ab77-de5b2f833009 +#= +This is a prediction of eth growth by the year 2025. +=# +grow_factor_of_eth = 6 + +# ╔═╡ 8d251250-a9be-4689-9318-106eb240f4e6 +let + usd_per_year_per_eth = (price_eth_today * grow_factor_of_eth) * yield + eth_needed = usd_per_month_to_live * 12 / usd_per_year_per_eth + bitcoin_needed = eth_needed * price_eth_today / price_btc_today + md""" + Number of ETH needed today: $eth_needed + + Number of bitcoin needed today: $bitcoin_needed + """ +end + +# ╔═╡ Cell order: +# ╟─701a0e50-df22-11eb-2eea-a56013a4ffdd +# ╟─6605e130-e014-4de3-8837-baf4eb77418c +# ╟─9a96995c-89be-4fa9-97f6-81b7bd5898b8 +# ╟─ff389da3-6761-4d78-873e-94f1b5319a05 +# ╟─52becaeb-d325-4b17-a4ba-de19d6eb2cf8 +# ╟─5bd92f00-0b9d-4217-92dc-ed3f7234666f +# ╟─3a5fdd52-6445-4b10-ab77-de5b2f833009 +# ╟─8d251250-a9be-4689-9318-106eb240f4e6 diff --git a/fourier.jl b/fourier.jl new file mode 100644 index 0000000..6530094 --- /dev/null +++ b/fourier.jl @@ -0,0 +1,661 @@ +### A Pluto.jl notebook ### +# v0.14.8 + +using Markdown +using InteractiveUtils + +# This Pluto notebook uses @bind for interactivity. When running this notebook outside of Pluto, the following 'mock version' of @bind gives bound variables a default value (instead of an error). +macro bind(def, element) + quote + local el = $(esc(element)) + global $(esc(def)) = Core.applicable(Base.get, el) ? Base.get(el) : missing + el + end +end + +# ╔═╡ 799b8d42-e5ff-11ea-37b7-ede54edaf352 +begin + using PlutoUI, Plots, LinearAlgebra, Statistics, CSV, DataFrames, StatsPlots + plotly() + default(legend = false, minorgrid = true) +end + +# ╔═╡ 60408320-3f7b-11eb-1773-e706096dc457 +html""" +""" + +# ╔═╡ 7f284800-e606-11ea-0e36-a5e6934dbcf1 +md"## Fourier series" + +# ╔═╡ 404bafa0-e7aa-11ea-0158-55902780a54e +md""" +The following code will parse a path from a SVG file (see ) +Then will compute the fourier series of the path with a number of coefficients to compute as a parameter. + +Great explanation about the Fourier series: +""" + +# ╔═╡ 34294d10-e785-11ea-373a-8f17483f4468 +# An SVG path command. +struct Command + command # Can be :move_to_rel :move_to_abs :vertical_rel :vertical_abs :horizontal_rel :horizontal_abs :line_rel :quadratic_bézier_rel :cubic_bézier_rel :cubic_bézier_abs :close_path + points :: Vector{Vector{Float64}} +end + +# ╔═╡ c90fa10e-e7ab-11ea-1dd1-7b5d60129245 +md"**`parse_parth`** take a path as a string and return a list of commands, see the struct `Command`" + +# ╔═╡ 2ae11130-e6df-11ea-3a92-55d4d69ef9a8 +function parse_path(path :: String) :: Vector{Command} + commands_str = split(path, " ") + commands = [] + for word ∈ commands_str + if word == "m" + push!(commands, Command(:move_to_rel, [])) + elseif word == "M" + push!(commands, Command(:move_to_abs, [])) + elseif word == "v" + push!(commands, Command(:vertical_rel, [])) + elseif word == "V" + push!(commands, Command(:vertical_abs, [])) + elseif word == "h" + push!(commands, Command(:horizontal_rel, [])) + elseif word == "H" + push!(commands, Command(:horizontal_abs, [])) + elseif word == "l" + push!(commands, Command(:line_rel, [])) + elseif word == "q" + push!(commands, Command(:quadratic_bézier_rel, [])) + elseif word == "c" + push!(commands, Command(:cubic_bézier_rel, [])) + elseif word == "C" + push!(commands, Command(:cubic_bézier_abs, [])) + elseif word == "z" || word == "Z" + push!(commands, Command(:close_path, [])) + elseif word[1] >= 'A' && word[1] <= 'z' + throw(ErrorException("Command not supported: $(word[1])")) + else + current = last(commands) + if current.command == :vertical_rel || current.command == :vertical_abs + push!(current.points, [0.0, parse(Float64, word)]) + elseif current.command == :horizontal_rel || current.command == :horizontal_abs + push!(current.points, [parse(Float64, word), 0.0]) + else + push!(current.points, parse.(Float64, split(word, ","))) + end + end + end + commands +end + +# ╔═╡ 077ff940-e7ac-11ea-3461-5f0c5a2b38b8 +md"Each bézier curve is converted in a number of points given by **`points_per_bézier`**" + +# ╔═╡ 28287c70-e785-11ea-37e0-8b6c4f8d31ef +points_per_bézier = 20 + +# ╔═╡ 1cde9fc0-e785-11ea-0644-7d58bfd03a5e +function commands_to_points(commands, reverse_y = true) :: Vector{Vector{Float64}} + points = [] + + B2(t, P₀, P₁, P₂) = (1 - t) * ((1 - t) * P₀ + t * P₁) + t * ((1 - t) * P₁ + t * P₂) + + B3(t, P₀, P₁, P₂, P₃) = (1 - t)^3 * P₀ + 3 * (1 - t)^2 * t * P₁ + 3 * (1 - t) * t^2 * P₂ + t^3 * P₃ + + for c ∈ commands + if c.command == :move_to_rel + push!(points, c.points[1] + (length(points) > 0 ? last(points) : [0.0, 0.0])) + elseif c.command == :move_to_abs + push!(points, c.points[1]) + elseif c.command == :vertical_rel || c.command == :horizontal_rel || c.command == :line_rel + for p ∈ c.points + push!(points, last(points) + p) + end + elseif c.command == :vertical_abs + push!(points, [last(points)[1], c.points[1][2]]) + elseif c.command == :horizontal_abs + push!(points, [c.points[1][1], last(points)[2]]) + elseif c.command == :quadratic_bézier_rel + for i ∈ 1:2:length(c.points)-1 + P₀ = last(points) + P₁ = P₀ + c.points[i] + P₂ = P₀ + c.points[i + 1] + + for t ∈ 0:1/points_per_bézier:1 + push!(points, B2(t, P₀, P₁, P₂)) + end + end + elseif c.command == :cubic_bézier_rel + for i ∈ 1:3:length(c.points)-2 + P₀ = last(points) + P₁ = P₀ + c.points[i] + P₂ = P₀ + c.points[i + 1] + P₃ = P₀ + c.points[i + 2] + for t ∈ 0:1/points_per_bézier:1 + push!(points, B3(t, P₀, P₁, P₂, P₃)) + end + end + + elseif c.command == :cubic_bézier_abs + for i ∈ 1:3:length(c.points)-2 + P₀ = last(points) + P₁ = c.points[i] + P₂ = c.points[i + 1] + P₃ = c.points[i + 2] + for t ∈ 0:1/points_per_bézier:1 + push!(points, B3(t, P₀, P₁, P₂, P₃)) + end + end + elseif c.command == :close_path + push!(points, copy(points[1])) + else + throw(ErrorException("command not supported: $c")) + end + end + + # Reverse y coordinate of each points because the y axis is to the bottom in . + if reverse_y + for p ∈ points + p[2] = -p[2] + end + end + + points +end + +# ╔═╡ 140a7e00-e785-11ea-2599-6fa3e0e8191e +function position_on_points(points, k :: Float64) :: Vector{Float64} + if !(0 <= k <= 1) + throw(ArgumentError("k must be in [0, 1]")) + end + + total_len = 0 + lengths = [] + for i ∈ 2:length(points) + l = norm(points[i] - points[i-1]) + total_len += l + push!(lengths, l) + end + + current_len = 0 + for i ∈ 1:length(lengths) + next_len = current_len + lengths[i] + + if current_len == next_len + continue + end + + if k >= current_len / total_len && k <= next_len / total_len + len_in_segment = k * total_len - current_len + return normalize(points[i+1] - points[i]) * len_in_segment + points[i] + end + current_len = next_len + end +end + +# ╔═╡ 50395b50-e82d-11ea-2da6-3baecb6100de +md"Some vecorial drawings (path):" + +# ╔═╡ a59ac5ce-e782-11ea-09ac-295f9e2292b4 +begin + path_box = "m 40.82143,73.327377 h 71.05952 V 136.07142 H 40.82143 Z" + commands_box = parse_path(path_box) + points_box = commands_to_points(commands_box) + + plot(points_box .|> Tuple, aspect_ratio = true) +end + +# ╔═╡ 62f97bb0-e762-11ea-0561-c75618e97061 +begin + path_letter_g = "m 65.054135,67.883872 v -1.007821 l 3.638702,-0.0059 v 3.187527 q -0.837898,0.667974 -1.72853,1.007821 -0.890633,0.333987 -1.82814,0.333987 -1.265636,0 -2.302754,-0.539067 -1.031259,-0.544886 -1.558607,-1.570286 -0.527348,-1.025399 -0.527348,-2.291034 0,-1.253917 0.521489,-2.33791 0.527348,-1.089853 1.511731,-1.617201 0.984383,-0.527349 2.267597,-0.527349 0.931649,0 1.681655,0.30469 0.755866,0.298831 1.183604,0.837898 0.427738,0.539067 0.650396,1.406262 l -1.0254,0.281252 q -0.193361,-0.656255 -0.480472,-1.031258 -0.287112,-0.375004 -0.82032,-0.597662 -0.533207,-0.228517 -1.183603,-0.228517 -0.779304,0 -1.347668,0.240236 -0.568364,0.234377 -0.919929,0.621099 -0.345706,0.386722 -0.539067,0.849617 -0.328128,0.796881 -0.328128,1.72853 0,1.148447 0.392581,1.921891 0.398441,0.773444 1.154307,1.148447 0.755866,0.375003 1.605482,0.375003 0.738287,0 1.441418,-0.281252 0.703131,-0.287112 1.066415,-0.609381 v -1.599622 z" + commands_letter_g = parse_path(path_letter_g) + points_letter_g = commands_to_points(commands_letter_g) + + plot(points_letter_g .|> Tuple, aspect_ratio = true) +end + +# ╔═╡ ede054d0-e86e-11ea-1002-6f4e3d3c517c +md"## 🐱" + +# ╔═╡ 5ddbafc0-e7a5-11ea-113a-232748354da7 +begin + path_cat = "M 1607.2,4071.3 C 1196.4,3961.8 1036.6,3831.7 865.4,3480.2 680.5,3099 612.1,2742.9 657.7,2416.6 726.2,1957.8 892.8,1542.5 1242,955.9 l 255.6,-429.1 -16,-189.5 c -9.1,-105 -43.4,-310.4 -75.3,-454.2 -34.2,-143.8 -68.5,-362.9 -82.2,-483.9 -13.7,-130.1 -41.1,-264.7 -68.5,-328.7 -79.9,-180.3 -239.6,-378.9 -449.6,-561.4 l -203.1,-175.8 -146.1,-360.6 c -226,-561.5 -253.3,-641.3 -301.3,-901.5 -22.8,-132.4 -47.9,-255.6 -54.8,-273.9 -13.7,-34.2 146.1,-230.5 216.8,-269.3 57.1,-29.7 146.1,11.4 198.6,89 27.4,43.3 38.8,111.8 41.1,253.3 0,175.8 11.4,216.8 107.3,454.2 59.3,143.8 132.4,287.6 164.3,319.5 125.5,130.1 289.9,203.1 698.4,308.1 219.1,54.8 429.1,114.1 467.9,127.8 57.1,20.5 70.8,18.3 95.9,-13.7 41.1,-54.8 36.5,-116.4 -16,-392.5 -25.1,-127.8 -45.6,-278.4 -45.6,-333.2 0,-205.4 255.6,-595.7 680.1,-1040.7 l 175.7,-184.9 271.6,6.8 c 262.5,6.8 271.6,9.1 317.3,66.2 38.8,54.8 41.1,68.5 13.7,123.3 -45.6,93.6 -136.9,150.6 -267,166.6 -111.8,16 -125.5,25.1 -310.4,226 -203.1,219.1 -235.1,278.4 -210,401.7 18.3,95.9 235.1,413.1 417.7,614 102.7,111.8 175.7,219.1 235.1,351.5 79.9,171.2 95.9,191.7 152.9,194 728.1,20.5 938,38.8 1063.6,84.5 l 130.1,47.9 200.8,-57.1 c 538.7,-148.4 495.3,-127.8 502.1,-225.9 13.7,-157.5 -86.7,-787.4 -157.5,-995.1 -59.4,-175.7 -66.2,-228.2 -66.2,-470.2 0,-244.2 4.6,-280.7 47.9,-335.5 45.7,-57.1 54.8,-59.3 280.7,-59.3 262.5,0 287.6,11.4 287.6,139.2 0,66.2 -16,98.2 -68.5,141.5 -36.5,32 -89,66.2 -114.1,73 -93.6,29.7 -38.8,189.4 248.8,739.5 367.5,707.5 415.4,755.5 725.8,755.5 109.6,0 219.1,13.7 283,38.8 l 107.3,36.5 146.1,-105 c 159.8,-111.8 292.1,-248.8 292.1,-299 0,-43.3 152.9,-178 730.4,-643.6 269.3,-219.1 547.8,-445 616.2,-502.1 68.5,-59.4 166.6,-125.5 219.1,-148.3 50.2,-22.8 139.2,-84.5 196.3,-134.7 l 105,-91.3 H 9527 c 237.4,0 244.2,2.3 310.4,66.2 86.7,89 84.5,141.5 -13.7,228.2 -73,63.9 -95.9,70.8 -226,70.8 H 9454 l -182.6,189.4 c -399.4,410.8 -582,634.5 -732.6,899.3 -237.4,413.1 -228.2,365.2 -237.4,1116.1 l -6.8,661.9 77.6,152.9 c 84.4,166.6 225.9,312.7 397.1,413.1 86.7,52.5 148.3,66.2 335.5,82.2 125.5,11.4 237.4,25.1 248.8,31.9 9.1,6.9 66.2,95.9 125.5,198.6 77.6,130.1 105,200.8 95.8,237.4 -6.8,29.7 -79.9,134.6 -159.8,235.1 -82.2,100.4 -155.2,210 -164.3,241.9 -20.6,84.5 -210,248.8 -362.9,312.7 l -127.8,54.8 6.8,228.2 c 6.8,223.7 4.6,228.2 -43.4,223.7 -31.9,-4.6 -52.5,9.1 -57.1,34.2 -6.8,29.7 -77.6,-27.4 -271.6,-216.8 -230.5,-226 -283,-267 -477,-356 -166.6,-79.9 -264.7,-146.1 -415.4,-283 -109.4,-98.4 -271.5,-217.1 -360.5,-265 -148.3,-82.2 -173.4,-89 -344.6,-86.7 -100.4,0 -303.6,20.5 -449.6,45.6 -146.1,25.1 -369.7,57.1 -495.3,70.8 -125.5,13.7 -511.2,68.5 -855.9,125.5 -951.7,150.6 -1024.8,159.8 -1517.8,157.5 -326.4,0 -575.2,-16 -876.4,-52.5 l -420,-50.2 -123.3,63.9 c -152.9,77.6 -415.4,315 -493,447.3 -130.1,221.4 -312.7,887.8 -312.7,1136.6 2.3,248.8 118.7,620.8 235.1,746.3 107.3,116.4 374.3,198.6 529.5,164.3 82.2,-18.3 178,-79.9 248.8,-157.5 66.2,-73 164.3,-63.9 207.7,20.6 43.4,84.5 43.4,235.1 2.3,317.2 -43.4,84.5 -216.8,239.6 -330.9,296.7 -109.5,57.1 -358.3,64 -538.6,13.8 z" + + commands_cat = parse_path(path_cat) + points_cat = commands_to_points(commands_cat, false) .* 1/100 # Here we down scale the cat to have more readable coordinates. + + plot(points_cat .|> Tuple, aspect_ratio = true) +end + +# ╔═╡ 71520a10-e785-11ea-3aea-19575bbadbfc +function to_complexes(points, n :: Integer) :: Vector{Complex{Float64}} + f(t) = let + pos = position_on_points(points, t) + complex(pos[1], pos[2]) + end + [f(t) for t ∈ LinRange(0, 1-1/n, n)] +end + +# ╔═╡ 3d763860-e771-11ea-0edc-094a80f44f74 +function plot_scatter(points, nb_points = 100) + scatter(to_complexes(points, nb_points), aspect_ratio = true) +end + +# ╔═╡ 762ceed0-e82d-11ea-34ec-e92e32a51815 +md"Same drawings sampled and transformed in complex numbers:" + +# ╔═╡ 88f3eb30-e79d-11ea-0df4-cd4ab90ae8cb +plot_scatter(points_box, 100) + +# ╔═╡ 724535f2-e780-11ea-1a67-69aceba56c81 +plot_scatter(points_letter_g, 100) + +# ╔═╡ 45da7c50-e7a8-11ea-2e48-35182730df7a +plot_scatter(points_cat, 200) + +# ╔═╡ 0b84ca5e-e785-11ea-3c6d-f75c72b7b825 +# Also known as 'Fourier descriptor'. +function nth_coefficient(complex_points, nth) :: Complex{Float64} + n = length(complex_points) + (complex_points .* (exp(-nth * 2π * im * t) for t ∈ LinRange(0, 1-1/n, n)) |> sum) / n +end + +# ╔═╡ 869d9bd0-eb92-11ea-2e0e-05306004cd45 +function fft(complex_points, nb_coef_positive) + [nth_coefficient(complex_points, k) for k ∈ -nb_coef_positive:nb_coef_positive] +end + +# ╔═╡ c1e999ee-e7aa-11ea-2a8a-83ae8f10dffc +begin + drawings = + Dict( + "box" => points_box, + "g_letter" => points_letter_g, + "cat" => points_cat + ) + + @bind drawing Select(["box" => "Box", "g_letter" => "G letter", "cat" => "Cat"]) +end + +# ╔═╡ 27aaf490-e770-11ea-3579-61a1d791051d +@bind nb_coef Slider(3:2:101; default = 25) + +# ╔═╡ 866c0140-e770-11ea-0e97-eb2d15bf29d6 +md"Number of coefficients: **$nb_coef**" + +# ╔═╡ 2593e000-eb84-11ea-3019-2775c11bb903 +begin + N = 500 + md"Number of complex numbers: **$N**" +end + +# ╔═╡ a9b906e0-eb83-11ea-21eb-f763fbf6ec1d +complexes = to_complexes(drawings[drawing], N) + +# ╔═╡ 47fe4880-e76d-11ea-307e-69f5f9f4cb3e +let + nb_coef_positive = nb_coef ÷ 2 + + coefficients = fft(complexes, nb_coef_positive) + + coefficient(freq) = let + i = freq + nb_coef_positive + 1 + if i >= 1 && i <= nb_coef + coefficients[i] + else + 0im + end + end + + t = LinRange(0, 1, N) + f(t) = (coefficient(freq) * exp(freq * 2π * im * t) for freq in -nb_coef_positive:nb_coef_positive) |> sum + + plot(f.(t), aspect_ratio = true) +end + +# ╔═╡ 111a7f50-eb7d-11ea-3a2e-8328258daddf + md"## Matching scores" + +# ╔═╡ b7175880-ec2f-11ea-3275-9b2c99a431e9 +begin + limits_real = real(complexes) |> extrema + limits_imag = imag(complexes) |> extrema + + real_size = abs(limits_real[2] - limits_real[1]) + imag_size = abs(limits_imag[2] - limits_imag[1]) + + md"""**real limits**: $limits_real + + **imaginary limits**: $limits_imag""" +end + +# ╔═╡ b8db0c70-3f7c-11eb-2441-9b0bf3b848fe +@bind reset Button("Reset to default values") + +# ╔═╡ 42062500-eb88-11ea-3443-c32b192a16a5 +begin + reset + md"θ: $(@bind θ Slider(range(0, 2π; length = 100); default = 0, show_value = true))" +end + +# ╔═╡ b0e10d00-eb88-11ea-379e-e1f03cf9d8c0 +begin + reset + md"Δx: $(@bind Δx Slider(range(-real_size, real_size; length = 100); default = 0, show_value = true))" +end + +# ╔═╡ c4ca66e0-eb88-11ea-1b3b-279a2ce5cabb +begin + reset + md"Δy: $(@bind Δy Slider(range(-imag_size, imag_size; length = 100); default = 0, show_value = true))" +end + +# ╔═╡ ca636560-eb89-11ea-0906-d1892d1d514c +begin + reset + md"Scale: $(@bind α Slider(range(5^-1, 5; length = 25); default = 1, show_value = true))" +end + +# ╔═╡ 664a3a20-eb90-11ea-099a-2fd33e5daf89 +begin + reset + md"noise: $(@bind noise Slider(range(0, .1; length = 100); default = 0, show_value = true))" +end + +# ╔═╡ cfb46240-92ea-11eb-06ee-5d71c9bfc4f5 +let + c = complexes[1:100] + nth_coefficient(c, 1), nth_coefficient(reverse(c), 1) +end + +# ╔═╡ 96980942-eb87-11ea-25e9-f710357f69c3 +function rotate(points, θ) + points * exp(θ * im) +end + +# ╔═╡ 1941bf40-eb7d-11ea-0624-ffe3a941cfd2 +function translate(points, Δxy) + points .+ Δxy +end + +# ╔═╡ 9eaf1c90-eb87-11ea-16af-4d85ca6c35eb +function scale(points, α) + points .* α +end + +# ╔═╡ a7cb1cc0-eb87-11ea-2f53-9339f8433bcc +function add_noise(points) + noise_factor = mean([real_size, imag_size]) * noise + points + noise_factor * randn(Complex{Float64}, length(points)) +end + +# ╔═╡ d59fbf70-eb87-11ea-0183-8de6a0137ee0 +begin + using Random + Random.seed!(1234) + + noised = add_noise(complexes) + translated = translate(noised, complex(Δx, Δy)) + rotated = rotate(translated, θ) + scaled = scale(rotated, α) + + scatter(scaled, aspect_ratio = true) +end + +# ╔═╡ 76ccc310-9228-11eb-063f-432008968caf +function relative_fft(complex_points, n) + coefficients = fft(complex_points, n ÷ 2) + + # Remove the frequency 0 to avoid to take in account the position. + coefficients[n ÷ 2 + 1] = 0 + + coefficients +end + +# ╔═╡ 0ff985fe-92e1-11eb-266d-e777af1ab6fd +a = norm.(relative_fft(complexes, 7)) + +# ╔═╡ 25cf94b2-92e1-11eb-2c86-6f980f7a5d29 +b = norm.(relative_fft(reverse(complexes), 7)) + +# ╔═╡ db88b910-eb90-11ea-3b93-054a0e9d1a10 +function score(coefficients₁, coefficients₂) + nb_coef = length(coefficients₁) + sum(abs.(norm.(coefficients₁) .- norm.(coefficients₂))) / (nb_coef - 1) +end + +# ╔═╡ 0abdd950-eb95-11ea-16a6-d19923afaabb +let + coefficients_ref = relative_fft(complexes, nb_coef) + coefficients_measured = relative_fft(scaled, nb_coef) + + s = score(coefficients_ref, coefficients_measured) + if abs(s) < 0.00001 + s = 0.0 + end + md"Score = **$s**" +end + +# ╔═╡ 02fac600-92e4-11eb-16e1-9f2a7322b8cb +score(a, b) + +# ╔═╡ 07c48fe0-ec2c-11ea-33be-1f9d4a6a2721 +md""" +## Finding a part of the image + +Input: +* A pattern to find (points) +* A curve: list of consecutive points + +Output: +* A score for each given point +""" + +# ╔═╡ 44923ef0-920c-11eb-3306-e51c06ba41df +begin + profile_filename = "ProfileRail.csv" + profile_nominal_filename = "NominalHead.csv" + + md""" + ### Example with real data + + * Rail profile: $profile_filename + * Nominal profile: $profile_nominal_filename""" +end + +# ╔═╡ 5280ed20-921d-11eb-148b-d582dabbfb45 +begin + profile = CSV.File(profile_filename; header=["X", "Y"]) |> DataFrame + transform!(profile, :, :Y => (-) => :Y, :X => (-) => :X) + nominal = CSV.File(profile_nominal_filename; header=["X", "Y"]) |> DataFrame +end + +# ╔═╡ 3381b020-9223-11eb-1198-434e98bcdf31 +@df nominal scatter(:X, :Y, aspect_ratio = true) + +# ╔═╡ 91f5dbb0-922b-11eb-05d2-e5799199c053 +begin + # Parameters. + search_box_enlarge_factor = 2 + nb_of_coefficients = 61 + md""" + #### Parameters + + *Search box enlarge factor* = $search_box_enlarge_factor + """ +end + +# ╔═╡ 3972e7e0-92c0-11eb-3df0-67b830fbc023 +rect((x1, y1), (x2, y2)) = Shape([(x1, y1), (x1, y2), (x2, y2), (x2, y1)]) + +# ╔═╡ 55f9daa0-92f6-11eb-1c48-6587d126d413 +md"Current point: $(@bind current_point NumberField(0:200, default = 0))" + +# ╔═╡ c45c2dc0-9226-11eb-2ac1-252cb75d7c12 +begin + local x_min, x_max = extrema(nominal[!, :X]) + local y_min, y_max = extrema(nominal[!, :Y]) + local x_mean = mean(nominal[!, :X]) + local y_mean = mean(nominal[!, :Y]) + + local nearest_point = 400 + #= + nearest_d = typemax(Float64) + for (i, (x, y)) ∈ enumerate(eachrow(nominal)) + d = norm([x - x_mean, y - y_mean]) + if d < nearest_d + nearest_point = i + nearest_d = d + end + end + =# + + @debug nominal[nearest_point, [:X, :Y]] + + search_box_width = (x_max - x_min) * search_box_enlarge_factor + search_box_height = (y_max - y_min) * search_box_enlarge_factor + + # Shift from the center. + δx = x_min + (x_max - x_min) / 2 - nominal[nearest_point, :X] + δy = y_min + (y_max - y_min) / 2 - nominal[nearest_point, :Y] + + @debug "------" + + @debug "width = $search_box_width, height = $search_box_height" + + #fft_profile = select(profile, :, :Y => ((y) -> -y) => :Y, :X => ((x) -> -x) => :X) + fft_nominal = relative_fft([Complex(x, y) for (x, y) ∈ eachrow(nominal)], nb_of_coefficients) + + local best_point = 1 + local best_score = typemax(Float64) + local windows = Vector{Shape}() + + scores = Vector{Tuple{Int64, Float64}}() + + for (i, (x, y)) ∈ enumerate(eachrow(profile)) + x_min = x + δx - search_box_width / 2 + x_max = x + δx + search_box_width / 2 + y_min = y + δy - search_box_height / 2 + y_max = y + δy + search_box_height / 2 + + surrounding_points = [(x, y) for (x, y) ∈ eachrow(profile) if x_min ≤ x ≤ x_max && y_min ≤ y ≤ y_max] + + push!(windows, rect((x_min, y_min), (x_max, y_max))) + + if length(surrounding_points) ≥ nb_of_coefficients + + fft_profile = + relative_fft( + map(surrounding_points) do (x, y) + Complex(x, y) + end, + nb_of_coefficients + ) + + s = score(fft_nominal, fft_profile) + + push!(scores, (i, s)) + + if s < best_score + best_score = s + best_point = i + end + + #@debug "score: $s" + + #= + @debug "fft_profile: $fft_profile" + @debug "surrounding_points: $surrounding_points" + @debug "current point: ($x, $y)" + @debug "x limits: ($x_min, $x_max)" + @debug "y limits: ($y_min, $y_max)" + =# + + #break + end + end + + function highlight_point_with_window(i, color = :red) + point_value = profile[i, [:X, :Y]] + plot!((point_value.X, point_value.Y), aspect_ratio = true, seriestype = :scatter, markershape = :hexagon, markersize = 6, seriescolor = color) + plot!(windows[i], fillcolor = :transparent, linecolor = color, linealpha = 0.5) + end + + @df profile plot(:X, :Y, aspect_ratio = true, seriestype = :scatter, size = (1200, 800)) + + highlight_point_with_window(best_point, :orange) + highlight_point_with_window(current_point) + #highlight_point_with_window(537) + +end + +# ╔═╡ 4e116eb0-92c5-11eb-1449-ad3ef2e3cbde +let + scatter(scores, size = (1200, 800)) +end + +# ╔═╡ Cell order: +# ╟─60408320-3f7b-11eb-1773-e706096dc457 +# ╠═799b8d42-e5ff-11ea-37b7-ede54edaf352 +# ╟─7f284800-e606-11ea-0e36-a5e6934dbcf1 +# ╟─404bafa0-e7aa-11ea-0158-55902780a54e +# ╠═34294d10-e785-11ea-373a-8f17483f4468 +# ╟─c90fa10e-e7ab-11ea-1dd1-7b5d60129245 +# ╟─2ae11130-e6df-11ea-3a92-55d4d69ef9a8 +# ╟─077ff940-e7ac-11ea-3461-5f0c5a2b38b8 +# ╟─28287c70-e785-11ea-37e0-8b6c4f8d31ef +# ╟─1cde9fc0-e785-11ea-0644-7d58bfd03a5e +# ╟─140a7e00-e785-11ea-2599-6fa3e0e8191e +# ╟─50395b50-e82d-11ea-2da6-3baecb6100de +# ╠═a59ac5ce-e782-11ea-09ac-295f9e2292b4 +# ╠═62f97bb0-e762-11ea-0561-c75618e97061 +# ╟─ede054d0-e86e-11ea-1002-6f4e3d3c517c +# ╟─5ddbafc0-e7a5-11ea-113a-232748354da7 +# ╠═71520a10-e785-11ea-3aea-19575bbadbfc +# ╠═3d763860-e771-11ea-0edc-094a80f44f74 +# ╟─762ceed0-e82d-11ea-34ec-e92e32a51815 +# ╟─88f3eb30-e79d-11ea-0df4-cd4ab90ae8cb +# ╟─724535f2-e780-11ea-1a67-69aceba56c81 +# ╟─45da7c50-e7a8-11ea-2e48-35182730df7a +# ╠═0b84ca5e-e785-11ea-3c6d-f75c72b7b825 +# ╠═869d9bd0-eb92-11ea-2e0e-05306004cd45 +# ╠═c1e999ee-e7aa-11ea-2a8a-83ae8f10dffc +# ╟─27aaf490-e770-11ea-3579-61a1d791051d +# ╟─866c0140-e770-11ea-0e97-eb2d15bf29d6 +# ╟─2593e000-eb84-11ea-3019-2775c11bb903 +# ╠═a9b906e0-eb83-11ea-21eb-f763fbf6ec1d +# ╠═47fe4880-e76d-11ea-307e-69f5f9f4cb3e +# ╟─111a7f50-eb7d-11ea-3a2e-8328258daddf +# ╟─b7175880-ec2f-11ea-3275-9b2c99a431e9 +# ╟─42062500-eb88-11ea-3443-c32b192a16a5 +# ╟─b0e10d00-eb88-11ea-379e-e1f03cf9d8c0 +# ╟─c4ca66e0-eb88-11ea-1b3b-279a2ce5cabb +# ╟─ca636560-eb89-11ea-0906-d1892d1d514c +# ╟─664a3a20-eb90-11ea-099a-2fd33e5daf89 +# ╠═b8db0c70-3f7c-11eb-2441-9b0bf3b848fe +# ╟─d59fbf70-eb87-11ea-0183-8de6a0137ee0 +# ╠═0abdd950-eb95-11ea-16a6-d19923afaabb +# ╠═0ff985fe-92e1-11eb-266d-e777af1ab6fd +# ╠═25cf94b2-92e1-11eb-2c86-6f980f7a5d29 +# ╠═cfb46240-92ea-11eb-06ee-5d71c9bfc4f5 +# ╠═02fac600-92e4-11eb-16e1-9f2a7322b8cb +# ╟─96980942-eb87-11ea-25e9-f710357f69c3 +# ╟─1941bf40-eb7d-11ea-0624-ffe3a941cfd2 +# ╟─9eaf1c90-eb87-11ea-16af-4d85ca6c35eb +# ╟─a7cb1cc0-eb87-11ea-2f53-9339f8433bcc +# ╠═76ccc310-9228-11eb-063f-432008968caf +# ╠═db88b910-eb90-11ea-3b93-054a0e9d1a10 +# ╟─07c48fe0-ec2c-11ea-33be-1f9d4a6a2721 +# ╠═44923ef0-920c-11eb-3306-e51c06ba41df +# ╠═5280ed20-921d-11eb-148b-d582dabbfb45 +# ╠═3381b020-9223-11eb-1198-434e98bcdf31 +# ╠═91f5dbb0-922b-11eb-05d2-e5799199c053 +# ╠═3972e7e0-92c0-11eb-3df0-67b830fbc023 +# ╟─55f9daa0-92f6-11eb-1c48-6587d126d413 +# ╠═c45c2dc0-9226-11eb-2ac1-252cb75d7c12 +# ╠═4e116eb0-92c5-11eb-1449-ad3ef2e3cbde -- 2.45.2