library(trackViewer)
library(InteractionSet)
# load the interaction data.
# to import your own data, please refer trackViewer help documents
gi_nij <- readRDS(system.file("extdata", "nij.chr6.51120000.53200000.gi.rds",
  package = "geomeTriD"
))
# the data is a GInteractions object with metadata score
head(gi_nij, n = 2)

# define a range to plot
range_chr6 <- GRanges("chr6", IRanges(51120000, 53200000))
# plot it if interactive
if (interactive()) { ## not run to reduce the package size.
  mdsPlot(gi_nij, range = range_chr6, k = 3, render = "threejs")
}

# create gene annotations
library(TxDb.Hsapiens.UCSC.hg19.knownGene)
library(org.Hs.eg.db)
## get all genes
feature.gr <- genes(TxDb.Hsapiens.UCSC.hg19.knownGene)
## subset the data by target viewer region
feature.gr <- subsetByOverlaps(feature.gr, range(regions(gi_nij)))
## assign symbols for each gene
symbols <- mget(feature.gr$gene_id, org.Hs.egSYMBOL, ifnotfound = NA)
feature.gr$label[lengths(symbols) == 1] <- 
  unlist(symbols[lengths(symbols) == 1])
## assign colors for each gene
feature.gr$col <- sample(1:7, length(feature.gr), replace = TRUE)
## Here we replaced the some gene feature as cis-regulatory elements
## to validate how it will be treated by the plot
feature.gr$type <- sample(c("cRE", "gene"),
  length(feature.gr),
  replace = TRUE,
  prob = c(0.1, 0.9)
)
## set the size the cRE in the plot
feature.gr$pch <- rep(NA, length(feature.gr))
feature.gr$pch[feature.gr$type == "cRE"] <- 0.5
## features header
head(feature.gr, n = 2)

# plot it if interactive
if (interactive()) { ## not run to reduce the package size.
  mdsPlot(gi_nij,
    range = range_chr6, feature.gr = feature.gr,
    k = 3, render = "threejs"
  )
}

# one layer of signals. It is a `track` object
ctcf <- readRDS(system.file("extdata", "ctcf.sample.rds",
  package = "trackViewer"
))
# plot it if interactive
if (interactive()) { ## not run to reduce the package size.
  mdsPlot(gi_nij,
    range = range_chr6, feature.gr = feature.gr,
    genomicSigs = ctcf,
    ## we reverse the ATAC-seq signale to show where is open
    reverseGenomicSigs = TRUE,
    k = 3, render = "threejs"
  )
}

# create a random signal for demonstration
set.seed(1)
randomSig <- ctcf$dat
randomSig <- randomSig[sort(sample(seq_along(randomSig), 50))]
randomSig$score <- randomSig$score * 2 * runif(50)
head(randomSig, n = 2)

objs <- mdsPlot(gi_nij,
  range = range_chr6, feature.gr = feature.gr,
  genomicSigs = list(ctcf = ctcf, test = randomSig),
  reverseGenomicSigs = c(TRUE, FALSE),
  k = 3, render = "none"
)

threeJsViewer(objs)

## not run to reduce the package size.
library(rgl)
library(manipulateWidget)
rgl::clear3d() # Remove the earlier display
rglViewer(objs, background = "white")
rgl::rglwidget() %>%
  rgl::toggleWidget(tags = "tick_minor") %>%
  toggleWidget(tags = "tick_major") %>%
  toggleWidget(tags = "tick_labels") %>%
  toggleWidget(tags = "ctcf") %>%
  toggleWidget(tags = "test") %>%
  toggleWidget(tags = "backbone") %>%
  toggleWidget(tags = "gene_body") %>%
  toggleWidget(tags = "tss_labels") %>%
  toggleWidget(tags = "gene_labels") %>%
  toggleWidget(tags = "cRE") %>%
  rgl::asRow(last = 10)

## get the backbone 3d positions,
## it will be used as a position reference for the interaction data.
xyz <- extractBackbonePositions(objs)
# make the example easy to see
gi.subset <- gi_nij[distance(first(gi_nij), second(gi_nij))>50000]
hic <- create3dGenomicSignals(
  gi.subset, 
  xyz,
  name='segment', # name prefix for the geometry
  tag='hic', # name for the layer in the scene
  color = c('green', 'darkred'),
  topN=3, # only plot the top 3 events ordered by the scores.
  lwd.maxGenomicSigs = 3,
  alpha=0.5
  )
hic2 <- create3dGenomicSignals(
  gi.subset,
  xyz,
  name='polygon', # name prefix for the geometry
  tag='hic', # name for the layer in the scene
  color = c('blue', 'darkred'),
  topN=seq(4, 10), # only plot the top 4-10 events ordered by the scores.
  type = 'polygon'
  )
width(regions(gi.subset)) <- 101#for high resolution input (such as reads pairs)
hic3 <- create3dGenomicSignals(
  gi.subset, 
  xyz,
  name='segment', # name prefix for the geometry
  tag='hic', # name for the layer in the scene
  color = rainbow(40),
  topN=length(gi.subset),
  lwd.maxGenomicSigs = 5,
  alpha = 0.1
)
# plot it if interactive
if (interactive()) { ## not run to reduce the package size.
  threeJsViewer(c(objs, hic, hic2, hic3))
}

for(i in seq_along(objs)){
  if(grepl('tick|arrow|cRE|tss', names(objs)[i])){
    objs[[i]]$layer <- 'bottom'
  }
}
if(interactive()){
  threeJsViewer(objs)
}

objs2 <- lapply(objs, function(.ele){
  .ele$side <- 'right'
  .ele
})
if(interactive()){
  threeJsViewer(objs, objs2)
}

## simulate a backbone
library(GenomicRanges)
library(geomeTriD)
backbone <- GRanges(1, IRanges(seq(1, 10000, by=1000), width = 1000))
## set x, y, z position, we will use the torus knots equation
t <- seq(0, 2*pi, length.out=length(backbone))
p <- 3
q <- 2
r <- cos(p*t) + 2
backbone$x <- 2*r*cos(q*t)
backbone$y <- 2*r*sin(q*t)
backbone$z <- -2*sin(p*t)

## simulate the signals to be add
N1 <- 30
sig1 <- GRanges(1, IRanges(sample(seq.int(10000), N1), width = 10),
                score=runif(N1, min=1, max=10))
N2 <- 8
sig2 <- GRanges(1, IRanges(sample(seq.int(10000), N2), width = 1),
                score=runif(N2, min=0.05, max=1),
                shape=sample(c('dodecahedron', 'icosahedron',
                               'octahedron', 'tetrahedron'),
                             N2, replace = TRUE),
                color=sample(seq.int(7)[-1], N2, replace = TRUE))

## do smooth for the coordinates
backbone <- smooth3dPoints(backbone, resolution=10)

## create the line geometry
geo_backbone <- threeJsGeometry(
  x=c(backbone$x0, backbone$x1[length(backbone)]),
  y=c(backbone$y0, backbone$y1[length(backbone)]),
  z=c(backbone$z0, backbone$y1[length(backbone)]),
  colors = "black",
  type = "line",
  tag = "tagNameHere",
  properties= list(size = 2)
)

## create the sphere geometry for sig1
geo_sig1_sphere <- create3dGenomicSignals(
  GenoSig = sig1,
  targetObj = backbone,
  positionTransformFun = function(x) {
      x$y0 <- x$y0 + 0.5 # offset from y
      x$y1 <- x$y1 + 0.5
      x
    }, reverseGenomicSigs = FALSE,
    type = 'sphere',
    tag = 'head', name = 'head',
  color = c('blue', 'red'),
  radius = .25)
## create the cylinder geometry for sig1
geo_sig1_body <- create3dGenomicSignals(
  GenoSig = sig1,
  targetObj = backbone,
  positionTransformFun = function(x) {
      x$y0 <- x$y0 + 0.25 # offset from y
      x$y1 <- x$y1 + 0.25
      x
    }, reverseGenomicSigs = FALSE,
    type = 'cylinder',
    tag = 'body', name = 'body',
  height = .5, radiusTop = .05, radiusBottom = .15)
## create two small spheres
geo_sig1_ear1 <- create3dGenomicSignals(
  GenoSig = sig1,
  targetObj = backbone,
  positionTransformFun = function(x) {
      x$x0 <- x$x0 - 0.25 # offset from x
      x$x1 <- x$x1 - 0.25
      x$y0 <- x$y0 + 0.75 # offset from y
      x$y1 <- x$y1 + 0.75
      x
    }, reverseGenomicSigs = FALSE,
    type = 'sphere',
    tag = 'ear', name = 'ear1',
  color = c('blue', 'red'),
  radius = .1)
geo_sig1_ear2 <- create3dGenomicSignals(
  GenoSig = sig1,
  targetObj = backbone,
  positionTransformFun = function(x) {
      x$x0 <- x$x0 + 0.25 # offset from x
      x$x1 <- x$x1 + 0.25
      x$y0 <- x$y0 + 0.75 # offset from y
      x$y1 <- x$y1 + 0.75
      x
    }, reverseGenomicSigs = FALSE,
    type = 'sphere',
    tag = 'ear', name = 'ear2',
  color = c('blue', 'red'),
  radius = .1)
geo_sig1_nose <- create3dGenomicSignals(
  GenoSig = sig1,
  targetObj = backbone,
  positionTransformFun = function(x) {
      x$y0 <- x$y0 + 0.5 # offset from y
      x$y1 <- x$y1 + 0.5
      x$z0 <- x$z0 - 0.25 # offset from z
      x$z1 <- x$z1 - 0.25
      x
    }, reverseGenomicSigs = FALSE,
    type = 'capsule',
    tag = 'nose', name = 'nose',
  color = c('red', 'blue'),
  radius = .05, height=.05)

## create the icosahedron geometry for sig2
geo_sig2 <- lapply(seq_along(sig2), function(id){
  create3dGenomicSignals(
    GenoSig = sig2[id],
    targetObj = backbone,
    reverseGenomicSigs = FALSE,
    type=sig2[id]$shape,
    tag = 'sig2', name = as.character(id),
    color = sig2[id]$color,
    radius = sig2[id]$score
  )})

## create the a customized geometry by a json file
geo_json <- create3dGenomicSignals(
  GenoSig = GRanges(1, IRanges(5000, width=1), score=1),
  targetObj = backbone,
  type = 'json',
  tag = 'json', name='monkey',
  color = 'green',
  rotation = c(pi, 0, pi), ## we need to rotate the model to change the facing side
  path = system.file('extdata', 'suzanne_buffergeometry.json',
                     package = 'geomeTriD')
)
threeJsViewer(c(backbone=geo_backbone,
                geo_sig1_sphere, geo_sig1_body,
                geo_sig1_ear1, geo_sig1_ear2, geo_sig1_nose,
                geo_sig2, geo_json))