We need functions from the notebook L13 Sparse + Low-Rank Splitting.
In [1]:
using LinearAlgebra
# Shrinkage
function Shr(x::Array{T},τ::T) where T
sign.(x).*max.(abs.(x).-τ,zero(T))
end
# Singular value thresholding
function D(A::Array{T},τ::T) where T
# This can be replaced by a faster approach
V=svd(A)
S=Shr(V.S,τ)
k=count(!iszero,S)
return (V.U[:,1:k]*Diagonal(S[1:k]))*V.Vt[1:k,:]
end
Out[1]:
In [2]:
function PCPAD(A::Array{T}) where T
# Initialize
δ=1.0e-7
tol=δ*norm(A)
m,n=size(A)
S=zero(A)
Y=zero(A)
L=zero(A)
T₁=zero(A)
μ=(m*n)/(4*(norm(A[:],1)))
μ₁=one(T)/μ
λ=one(T)/sqrt(max(map(T,m),n))
λμ₁=λ*μ₁
ν=1e20
maxiter=1000
iterations=0
# Iterate
while (ν>tol) && iterations<maxiter
# println(iterations," ",ν)
iterations+=1
L.=D(A-S+μ₁.*Y,μ₁)
S.=Shr(A-L+μ₁.*Y,λμ₁)
T₁.=(A-L-S)
Y.+=(μ.*T₁)
ν=norm(T₁)
end
L,S, iterations
end
Out[2]:
In [3]:
# For compilation
A0=rand(3,3)
Out[3]:
In [4]:
L,S,iter=PCPAD(A0)
Out[4]:
In [5]:
rank(L), L+S-A0
Out[5]:
In [6]:
L
Out[6]:
In [7]:
# Packages for video manipulation
# import Pkg; Pkg.add("VideoIO")
# import Pkg; Pkg.add("Makie")
In [8]:
using Images
using Makie, VideoIO
In [9]:
varinfo(VideoIO)
Out[9]:
In [11]:
video=VideoIO.open("files/visor.avi")
Out[11]:
In [12]:
playvideo(video)
In [10]:
VideoIO.get_duration("files/visor.avi")
Out[10]:
In [13]:
# Read the frames in Gray
f = VideoIO.openvideo("files/visor.avi",target_format=VideoIO.AV_PIX_FMT_GRAY8)
frames=Array{Array{Gray{Normed{UInt8,8}},2},1}(undef,0)
while !eof(f)
push!(frames,read(f))
end
close(f)
In [14]:
# Number of frames
length(frames)
Out[14]:
In [25]:
# See particular frame
frames[1401]
Out[25]:
In [18]:
# Make shorter video clip
clip=frames[1401:1450]
Out[18]:
In [19]:
size(clip[1])
Out[19]:
In [20]:
# Turn frames into tall matrix
mi,ni=size(clip[1])
m=mi*ni
n=length(clip)
A=Array{Float64}(undef,m,n)
for i=1:n
A[:,i]=vec(float(clip[i]))
end
In [21]:
size(A)
Out[21]:
In [22]:
norm(A)
Out[22]:
In [24]:
# For orientation
@time Q=qr(A);
@time D(A,0.5);
In [22]:
# Compute the splitting - 4 minutes
@time L,S,iters=PCPAD(A)
Out[22]:
In [24]:
# Reconstruct the low-rank video component
rank(L), norm(A-L-S), norm(A)
Out[24]:
In [25]:
svdvals(L)[1:20]
Out[25]:
In [26]:
# How to restore the video?
# First frame of the low-rank part
v1=L[:,1]
v2=reshape(v1,mi,ni)
map(Gray{Normed{UInt8,8}},v2)
Out[26]:
In [27]:
# First frame of the sparse part
s1=S[:,1]
s2=reshape(s1,mi,ni)
map(Gray{Normed{UInt8,8}},clamp01.(s2.+0.5))
Out[27]:
In [28]:
LowRank=similar(clip)
Sparse=similar(clip)
for i=1:n
LowRank[i]=reshape(L[:,i],mi,ni)
Sparse[i]=clamp01.(reshape(S[:,i],mi,ni).+0.5)
end
In [29]:
# import Pkg; Pkg.add("JLD")
In [26]:
# Let us save the results
using JLD
In [30]:
@save "files/visor_results.jld" LowRank Sparse
In [27]:
@load "files/visor_results.jld"
Out[27]:
In [29]:
LowRank[1]
Out[29]:
In [30]:
Sparse[1]
Out[30]:
In [28]:
# Play the LowRank part a bit slower (framerate=10)
props = [:priv_data => ("crf"=>"22","preset"=>"medium")]
encodevideo("files/LowRank.mp4",LowRank,framerate=10,AVCodecContextProperties=props)
In [32]:
videoLowRank=VideoIO.open("files/LowRank.mp4")
playvideo(videoLowRank)
In [33]:
# Play the Sparse part
encodevideo("files/Sparse.mp4",Sparse,framerate=10,AVCodecContextProperties=props)
In [36]:
videoSparse=VideoIO.open("files/Sparse.mp4")
playvideo(videoSparse)
In [ ]: