library(depstats)
library(progress)
library(fst)
library(dplyr)
library(imager)

#Additional Tests
#Dependence Structures from Recent Papers

##1
#Correlated Gaussian
t1 <- function(n, rho = 0.5) {
  x <- rnorm(n, 0, 1)
  y <- rho * x + (1 - rho) * rnorm(n, 0, 1)
  X <- cbind(x, y)
  return(X)
}
##2
#Quantile of t distribution with df=1
t2 <- function(n, rho = 0.5) {
  xs <- rnorm(n, 0, 1)
  ys <- rho * xs + (1 - rho) * rnorm(n, 0, 0.2)
  x <- qt(pnorm(xs), df = 4)
  y <- qt(pnorm(ys), df = 4)
  X <- cbind(x, y)
  return(X)
}
##3
#Correlated Laplace
t3 <- function(n, rho) {
  X <- LaplacesDemon::rmvl(n, c(0, 0),
                           matrix(c(1,   rho,
                                    rho, 1  ), nrow = 2, byrow = TRUE))
  return(X)
}
##4
#tanh
t4 <- function(n) {
  x <- rnorm(n)
  y <- tanh(x)
  X <- cbind(x, y)
  return(X)
}
##5
#Ishigami Style
t5 <- function(n, T = 1) {
  U <- runif(n)
  V <- runif(n)
  W <- runif(n)
  x <- U
  y <- sin(T * U) + 4 * sin(V) ^ 2 + 0.5 * W ^ 4 * sin(T * U)
  X <- cbind(x, y)
  return(X)
}
##6
#Tree Ring
t6 <- function(n, sigma = 0.1, rings = 10) {
  L <- sample(1:rings, n, replace = TRUE)
  theta <- runif(n, 0, 2 * pi)
  epsilonx <- rnorm(n, 0, sigma)
  epsilony <- rnorm(n, 0, sigma)
  x <- L * cos(theta) + epsilonx / 4
  y <- L * sin(theta) + epsilony / 4
  X <- cbind(x, y)
  return(X)
}
##7
#Change in variance
t7 <- function(n, p) {
  x <- runif(n)
  e <- rnorm(n)
  y <- abs(x) ^ p * e
  X <- cbind(x, y)
  return(X)
}


sample_sizes <- c(30, 50, 100, 200, 300, 400)


#only generate samples for the remaining ones
sample_sizes_remain <- as.character(sample_sizes) %>% setdiff(list.dirs('./additional_tests/', full.names = FALSE))

I1 <- imager::load.image("I1.jpg") #infinity
I2 <- load.image("I2.jpg") #pi

for (sample_size in sample_sizes_remain) {
  cat('working on sample size', sample_size, '...\n')

  n <- as.integer(sample_size)

  set.seed(n)

  #generate additional test samples
  X <- c()
  X1 <- depgen(1000, n, 't3(n, runif(1, 0.4, 0.5))',
               randrotate = FALSE)
  X2 <- depgen(1000, n, 't3(n, runif(1, 0.5, 0.6))',
               randrotate = FALSE)
  X <- rbind(X, X1, X2)
  X1 <- depgen(1000, n, 't5(n, runif(1, 0.5, 1.5))',
               randrotate = FALSE)
  X2 <- depgen(1000, n, 't5(n, runif(1, pi, 3 * pi))',
               randrotate = FALSE)
  X <- rbind(X, X1, X2)
  X1 <- depgen(1000, n, 't6(n, runif(1, 0, 0.5), sample(2:10, 1))',
               randrotate = FALSE)
  X2 <- depgen(1000, n, 't6(n, runif(1, 0, 1), sample(2:10, 1))',
               randrotate = FALSE)
  X <- rbind(X, X1, X2)
  X1 <- depgen(1000, n, 't7(n, runif(1, 0.4, 0.5))',
               randrotate = FALSE)
  X2 <- depgen(1000, n, 't7(n, runif(1, 0.5, 0.6))',
               randrotate = FALSE)
  X <- rbind(X, X1, X2)

##################################################
##################################################

  X1 <- depgen(1000, n, 'normnoise(randimage(n, I1), runif(1, 0, 0.1))',
                randrotate = TRUE, print = TRUE)
  X2 <- depgen(1000, n, 'normnoise(randimage(n, I1), runif(1, 0, 0.2))',
                randrotate = TRUE, print = TRUE)
  X <- rbind(X, X1, X2)

##################################################
##################################################

  X1 <- depgen(1000, n, 'normnoise(randimage(n, I2), runif(1, 0, 0.5))',
                randrotate = TRUE, print = TRUE)
  X2 <- depgen(1000, n, 'normnoise(randimage(n, I2), runif(1, 0, 1))',
                randrotate = TRUE, print = TRUE)
  X <- rbind(X, X1, X2)

  
  cat('working on dependence measures ...\n')
  Xind <- sampleapply(X, n, 1) %>% as.data.frame()

  cat('working on images ...\n')
  Ximage <- sampleapply(X, n, 0, grid = 25) %>% as.data.frame()

  cat('writing files ...\n')
  dir.create(paste0('./additional_tests/', sample_size)) #create folder

  write.fst(Xind, paste0('./additional_tests/', sample_size, '/', 'add.fst'))
  write.fst(Ximage, paste0('./additional_tests/', sample_size, '/', 'add_image.fst'))
}