diff -r a64766c44642 -r 843e7611b068 src/AlgorithmNoPrediction.jl --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/src/AlgorithmNoPrediction.jl Thu Apr 18 10:47:16 2024 +0300 @@ -0,0 +1,124 @@ +#################################################################### +# Predictive online PDPS for optical flow with known velocity field +#################################################################### + +__precompile__() + +module AlgorithmNoPrediction + +identifier = "pdps_known_noprediction" + +using Printf + +using AlgTools.Util +import AlgTools.Iterate +using ImageTools.Gradient + +using ..OpticalFlow: ImageSize, + Image, + pdflow! + +######################### +# Iterate initialisation +######################### + +function init_rest(x::Image) + imdim=size(x) + + y = zeros(2, imdim...) + Δx = copy(x) + Δy = copy(y) + x̄ = copy(x) + + return x, y, Δx, Δy, x̄ +end + +function init_iterates(xinit::Image) + return init_rest(copy(xinit)) +end + +function init_iterates(dim::ImageSize) + return init_rest(zeros(dim...)) +end + +############ +# Algorithm +############ + +function solve( :: Type{DisplacementT}; + dim :: ImageSize, + iterate = AlgTools.simple_iterate, + params::NamedTuple) where DisplacementT + + ################################ + # Extract and set up parameters + ################################ + + α, ρ = params.α, params.ρ + R_K² = ∇₂_norm₂₂_est² + γ = 1.0 + Λ = params.Λ + τ₀, σ₀ = params.τ₀, params.σ₀ + + τ = τ₀/γ + @assert(1+γ*τ ≥ Λ) + σ = σ₀*1/(τ*R_K²) + + println("Step length parameters: τ=$(τ), σ=$(σ)") + + ###################### + # Initialise iterates + ###################### + + x, y, Δx, Δy, x̄ = init_iterates(dim) + init_data = (params.init == :data) + + #################### + # Run the algorithm + #################### + + v = iterate(params) do verbose :: Function, + b :: Image, + v_known :: DisplacementT, + 🚫unused_b_next :: Image + + ################## + # Prediction step + ################## + if init_data + x .= b + init_data = false + end + + + + ############ + # PDPS step + ############ + + ∇₂ᵀ!(Δx, y) # primal step: + @. x̄ = x # | save old x for over-relax + @. x = (x-τ*(Δx-b))/(1+τ) # | prox + @. x̄ = 2x - x̄ # over-relax + ∇₂!(Δy, x̄) # dual step: y + @. y = (y + σ*Δy)/(1 + σ*ρ/α) # | + proj_norm₂₁ball!(y, α) # | prox + + ################################ + # Give function value if needed + ################################ + v = verbose() do + ∇₂!(Δy, x) + value = norm₂²(b-x)/2 + params.α*γnorm₂₁(Δy, params.ρ) + value, x, [NaN, NaN], nothing + end + + v + end + + return x, y, v +end + +end # Module + +