options(warn=-1) # turns off warnings, to turn on: "options(warn=0)"


library(imager)
library(png)

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=""))
}

source("nt_toolbox/toolbox_wavelet_meshes/meshgrid.R")
options(repr.plot.width=3.5, repr.plot.height=3.5)

n <- 256

f0 = load_image("nt_toolbox/data/hibiscus.png",n)

imageplot(f0)

sigma <- 0.04

f <- as.cimg(f0) + sigma*as.cimg(rnorm(n**2))

imageplot(clamp(f))

w <- 3
w1 <- 2*w + 1

grid <- meshgrid_4d(1:n, 1:n, (-w):w, (-w):w)
X <- grid$X ; Y <- grid$Y ; dX <- grid$Z ; dY <- grid$S
X <- X + dX
Y <- Y + dY

X[X < 1] <- 2-X[X < 1] 
Y[Y < 1] <- 2-Y[Y < 1]
X[X > n] <- 2*n-X[X > n]
Y[Y > n] <- 2*n-Y[Y > n]

I <- (X-1) + (Y-1)*n
for (i in 1:(n%/%w)){
    for (j in 1:(n%/%w)){
        I[i,j,,] <- t(I[i,j,,])
    }
}

Pi <- function(f){ array(as.vector(f)[I+1], c(n,n,w1*w1)) }

P <- Pi(f)

options(repr.plot.width=7, repr.plot.height=7)

for (i in 1:16){
    x <- sample(1:n, 1)
    y <- sample(1:n, 1)
    imageplot(array(P[x, y,],c(w1,w1)), '', c(4, 4, i))
}

Pmean <- function(f){ array(apply(Pi(f), c(2,1), mean), c(n,n)) }

options(repr.plot.width=3.5, repr.plot.height=3.5)

imageplot(Pmean(f))

p <- 100
psi <- function(f){ f**(1/p) }
ipsi <- function(f){ f**p }

imageplot(Pmean(abs(f)) - ipsi(Pmean(psi(abs(f)))))

r <- function(beta){ min(ceiling(beta*w1*w1), w1*w1-1) }

subsample <- function(x,s){ x[ , , s+1] }
phi <- function(f, beta){ subsample(aperm(apply(Pi(f), c(2,1), sort), c(2,3,1)), r(beta)) }

options(repr.plot.width=7, repr.plot.height=5)

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

## Insert your code here.

closing <- function(f){ phi(f, 0) }

options(repr.plot.width=3.5, repr.plot.height=3.5)

imageplot(closing(f))

opening <- function(f){ phi(f, 1) }

options(repr.plot.width=3.5, repr.plot.height=3.5)

imageplot(opening(f))

options(repr.plot.width=3.5, repr.plot.height=3.5)

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

## Insert your code here.

options(repr.plot.width=3.5, repr.plot.height=3.5)

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

## Insert your code here.

options(repr.plot.width=14, repr.plot.height=7)

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

## Insert your code here.

medfilt <- function(f){ phi(f, 1/2) }

options(repr.plot.width=3.5, repr.plot.height=3.5)

imageplot(medfilt(f))

options(repr.plot.width=7, repr.plot.height=5)

source("nt_solutions/denoisingadv_7_rankfilters/exo5.R")

## Insert your code here.

options(repr.plot.width=3.5, repr.plot.height=3.5)

imageplot(f1)
<script>
  $(document).ready(function(){
      $('div.prompt').hide();
  });
</script>