library(pracma)
library(imager)

# Importing the libraries
for (f in list.files(path="nt_toolbox/toolbox_general/", pattern="*.R")) {
    source(paste("nt_toolbox/toolbox_general/", f, sep=""))
}
for (f in list.files(path="nt_toolbox/toolbox_signal/", pattern="*.R")) {
    source(paste("nt_toolbox/toolbox_signal/", f, sep=""))
}

n = 512
f = load_image("nt_toolbox/data/lena.png", n)

options(repr.plot.width=4, repr.plot.height=4)
imageplot(f, 'Image f')

imageplot(f[c(((n/2) - 32):((n/2) + 32)), c(((n/2) - 32):((n/2) + 32))], 'Zoom')

options(repr.plot.width=5, repr.plot.height=5)
imageplot(-f, '-f', c(1, 2, 1))
imageplot(f[dim(f)[1]:1,], 'Flipped', c(1, 2, 2))

k = 9
h = matrix(1, k, k)
h = h / sum(h) #normalize

fh = convolution(f[,], h)

imageplot(fh, 'Blurred image')

F = fft(f[,]) / n

print(paste("Energy of Image: ", toString(norm(f[,]))))
print(paste("Energy of Fourier: ", toString(norm(Mod(F[,])))))

L = fftshift(log(Mod(F) + 1e-1))

imageplot(L, 'Log(Fourier transform)')

M = n**2 / 64

source("nt_solutions/introduction_4_fourier_wavelets/exo1.R")

# Insert your code here.

# Insert your code here.

T = 0.2

F = fft(f) / n

FT = F * (abs(F) > T)

L = fftshift(log(abs(FT[,]) + 1e-1))
imageplot(L, 'thresholded Log(Fourier transform)')

fM = Re(fft(FT, inverse=TRUE) * n)
# Normalize to 0-1
fM = (fM-min(fM))/(max(fM)-min(fM))

imageplot(fM, paste("Non - Linear, Fourier, SNR = ",  round(snr(f, fM), 1), "dB"))

m = sum(FT != 0)
print(paste("M/N = 1 /",as.integer(n**2/m)))

source("nt_solutions/introduction_4_fourier_wavelets/exo2.R")

# Insert your code here.

Jmin = 0

fw = perform_wavelet_transf(f, Jmin + 1, 1)

options(repr.plot.width=6, repr.plot.height=6)
plot_wavelet(fw)
title(main="Wavelet coefficients")

source("nt_solutions/introduction_4_fourier_wavelets/exo3.R")

# Insert your code here.

T = .15

fwT = fw * (abs(fw) > T)

plot_wavelet(fw)
title(main='Original coefficients')

plot_wavelet(fwT)
title(main='Thresholded coefficients')

fM = perform_wavelet_transf(fwT, Jmin + 1, -1)

# Normalize to 0-1
fM = (fM-min(fM))/(max(fM)-min(fM))

imageplot(fM, paste("Approximation, SNR = ",  round(snr(f, fM), 1), "dB"))

source("nt_solutions/introduction_4_fourier_wavelets/exo4.R")

# Insert your code here.

# Insert your code here.