Integer Programming Models

Documented in finalize3D getTexture gMOIPTheme ini3D loadView plotCones2D plotCones3D plotCriterion2D plotHull2D plotHull3D plotLines2D plotNDSet2D plotPlane3D plotPoints3D plotPolygon3D plotPolytope plotPolytope2D plotPolytope3D plotRectangle3D pngSize saveView .sizeM texToPng

## Functions for plotting

#' Plot the convex hull of a set of points in 2D.
#'
#' @param pts A matrix with a point in each row.
#' @param drawPoints Draw the points.
#' @param drawLines Draw lines of the facets.
#' @param drawPolygons Fill the hull.
#' @param addText Add text to the points. Currently `coord` (coordinates), `rownames` (rownames)
#'   and `both` supported or a vector with text.
#' @param addRays Add the ray defined by `direction`.
#' @param direction Ray direction. If i'th entry is positive, consider the i'th column of `pts`
#'   plus a value greater than on equal zero (minimize objective $i$). If negative, consider the
#'   i'th column of `pts` minus a value greater than on equal zero (maximize objective $i$).
#' @param drawPlot Draw the `ggplot`. Set to FALSE if you want to combine hulls in a single plot.
#' @param drawBBoxHull If `addRays` then draw the hull areas hitting the bounding box also.
#' @param m Minimum values of the bounding box.
#' @param M Maximum values of the bounding box.
#' @param ... Further arguments passed on the the `ggplot` plotting functions. This must be done as
#'   lists. Currently the following arguments are supported:
#'
#'   * `argsGeom_point`: A list of arguments for [`ggplot2::geom_point`].
#'   * `argsGeom_path`: A list of arguments for [`ggplot2::geom_path`].
#'   * `argsGeom_polygon`: A list of arguments for [`ggplot2::geom_polygon`].
#'   * `argsGeom_label`: A list of arguments for [`ggplot2::geom_label`].
#'
#' @return The `ggplot` object if `drawPlot = TRUE`; otherwise, a list of `ggplot` components.
#' @export
#' @importFrom rlang .data
#' @import ggplot2
#' @examples
#' library(ggplot2)
#' pts<-matrix(c(1,1), ncol = 2, byrow = TRUE)
#' plotHull2D(pts)
#' pts1<-matrix(c(2,2, 3,3), ncol = 2, byrow = TRUE)
#' plotHull2D(pts1, drawPoints = TRUE)
#' plotHull2D(pts1, drawPoints = TRUE, addRays = TRUE, addText = "coord")
#' plotHull2D(pts1, drawPoints = TRUE, addRays = TRUE, addText = "coord", drawBBoxHull = TRUE)
#' plotHull2D(pts1, drawPoints = TRUE, addRays = TRUE, direction = -1, addText = "coord")
#' pts2<-matrix(c(1,1, 2,2, 0,1), ncol = 2, byrow = TRUE)
#' plotHull2D(pts2, drawPoints = TRUE, addText = "coord")
#' plotHull2D(pts2, drawPoints = TRUE, addRays = TRUE, addText = "coord")
#' plotHull2D(pts2, drawPoints = TRUE, addRays = TRUE, direction = -1, addText = "coord")
#' ## Combine hulls
#' ggplot() +
#'   plotHull2D(pts2, drawPoints = TRUE, addText = "coord", drawPlot = FALSE) +
#'   plotHull2D(pts1, drawPoints = TRUE, drawPlot = FALSE) +
#'   gMOIPTheme() +
#'   xlab(expression(x[1])) +
#'   ylab(expression(x[2]))
#'
#' # Plotting an LP
#' A <- matrix(c(-3,2,2,4,9,10), ncol = 2, byrow = TRUE)
#' b <- c(3,27,90)
#' obj <- c(7.75, 10)
#' pts3 <- cornerPoints(A, b)
#' plotHull2D(pts3, drawPoints = TRUE, addText = "coord", argsGeom_polygon = list(fill = "red"))
plotHull2D <- function(pts,
                       drawPoints = FALSE,
                       drawLines = TRUE,
                       drawPolygons = TRUE,
                       addText = FALSE,
                       addRays = FALSE,
                       direction = 1,
                       drawPlot = TRUE,
                       drawBBoxHull = FALSE,
                       m = apply(pts,2,min)-5,
                       M = apply(pts,2,max)+5,
                       ...)
{
   args <- list(...)
   argsGeom_point <- mergeLists(list(size = 1.5, col = "black"), args$argsGeom_point)
   argsGeom_path <- mergeLists(list(size = 0.75, col = "grey40"), args$argsGeom_path)
   argsGeom_polygon <- mergeLists(list(fill = "black", alpha = 0.2), args$argsGeom_polygon)
   argsGeom_label <- mergeLists(list(nudge_x = 0.2), args$argsGeom_label)

   pts <- .checkPts(pts, p = 2)
   hull <- convexHull(pts, addRays = addRays, direction = direction, m = m, M = M)
   set <- hull$pts
   hull <- hull$hull
   d <- dimFace(set[,1:2, drop = FALSE])
   aES <- aes_string(x = colnames(set)[1], y = colnames(set)[2])
   aES1 <- aes_string(x = colnames(set)[1], y = colnames(set)[2], label = 'text')
   aES2 <- aes_string(x = colnames(set)[1], y = colnames(set)[2], label = 'addText')
   lst <- NULL

   if (d > 1) { # a polygon
      if (drawPolygons) {
         ptsT <- set[hull,]
         lst <- c(lst,
            do.call(geom_polygon, args = c(list(aES, data = ptsT), argsGeom_polygon)))
      }
   }
   if (d > 0) { # a line or polygon
      if (drawLines) {
         ptsT <- set[hull,]
         for (i in 2:nrow(ptsT)){
            if (!(ptsT$pt[i-1] == 0 & ptsT$pt[i] == 0) | drawBBoxHull) {
               ptsTT <- ptsT[(i-1):i,]
               lst <- c(lst,
                  do.call(geom_path, args = c(list(aES, data = ptsTT), argsGeom_path)))
            }
         }
         if (!(ptsT$pt[1] == 0 & ptsT$pt[nrow(ptsT)] == 0) | drawBBoxHull) {
            ptsTT <- ptsT[c(1,nrow(ptsT)),]
            lst <- c(lst,
               do.call(geom_path, args = c(list(aES, data = ptsTT), argsGeom_path)))
         }
      }
   }
   if (drawPoints | d == 0) {
      ptsT <- dplyr::filter(set, .data$pt == 1)
      lst <- c(lst,
         do.call(geom_point, args = c(list(aES, data = ptsT), argsGeom_point)))
   }

   ptsT <- dplyr::filter(set, .data$pt == 1)
   if (length(addText) > 1) {
      if (length(addText) == nrow(ptsT)) {
         lst <- c(lst,
            do.call(geom_label, args = c(list(aes_string(label = 'addText'), data = ptsT), argsGeom_label)))
      }
   } else if (length(addText) == 1) {
      if (addText == TRUE | addText == "coord") {
         ptsT$text <- paste0("(",ptsT[,1],",",ptsT[,2],")")
         lst <- c(lst,
            do.call(geom_label, args = c(list(aES1, data = ptsT), argsGeom_label)))
      } else if (addText == "rownames") {
         ptsT$text <- rownames(ptsT)
         lst <- c(lst,
            do.call(geom_label, args = c(list(aES1, data = ptsT), argsGeom_label)))
      } else if (addText == "both") {
         text <- paste0("(",ptsT[,1],",",ptsT[,2],")")
         ptsT$text <- paste(text,rownames(ptsT))
         lst <- c(lst,
            do.call(geom_label, args = c(list(aES1, data = ptsT), argsGeom_label)))
      } else if (addText != FALSE & length(addText) == nrow(ptsT)) {
         lst <- c(lst,
            do.call(geom_label, args = c(list(aES2, data = ptsT), argsGeom_label)))
      }
   }
   # if (latex) lst <- c(lst, xlab("$x_1$"), ylab("$x_2$"))
   # if (!latex) lst <- c(lst, xlab(expression(x[1])), ylab(expression(x[2])))
   if (drawPlot) lst <- ggplot() + lst
   return(lst)
}



#' Plot the polytope (bounded convex set) of a linear mathematical program (Ax <= b)
#'
#' This is a wrapper function calling [plotPolytope2D()] (2D graphics) and
#' [plotPolytope3D()] (3D graphics).
#'
#' @param A The constraint matrix.
#' @param b Right hand side.
#' @param obj A vector with objective coefficients.
#' @param type A character vector of same length as number of variables. If
#'   entry k is 'i' variable \eqn{k} must be integer and if 'c' continuous.
#' @param nonneg A boolean vector of same length as number of variables. If
#'   entry k is TRUE then variable k must be non-negative.
#' @param crit Either max or min (only used if add the iso-profit line)
#' @param faces A character vector of same length as number of variables. If
#'   entry k is 'i' variable \eqn{k} must be integer and if 'c' continuous.
#'   Useful if e.g. want to show the linear relaxation of an IP.
#' @param plotFaces If \code{True} then plot the faces.
#' @param plotFeasible If \code{True} then plot the feasible points/segments
#'   (relevant for IPLP/MILP).
#' @param plotOptimum Show the optimum corner solution point (if alternative solutions
#'   only one is shown) and add the iso-profit line.
#' @param latex If \code{True} make latex math labels for TikZ.
#' @param labels If \code{NULL} don't add any labels. If 'n' no labels but show the points. If equal
#'   `coord` add coordinates to the points. Otherwise number all points from one.
#' @param ... If 2D, further arguments passed on the the `ggplot` plotting functions. This must be
#'   done as lists. Currently the following arguments are supported:
#'
#'   * `argsFaces`: A list of arguments for [`plotHull2D`].
#'   * `argsFeasible`: A list of arguments for `ggplot2` functions:
#'      - `geom_point`: A list of arguments for [`ggplot2::geom_point`].
#'      - `geom_line`: A list of arguments for [`ggplot2::geom_line`].
#'   * `argsLabels`: A list of arguments for `ggplot2` functions:
#'      - `geom_text`: A list of arguments for [`ggplot2::geom_text`].
#'   * `argsOptimum`:
#'      - `geom_point`: A list of arguments for [`ggplot2::geom_point`].
#'      - `geom_abline`: A list of arguments for [`ggplot2::geom_abline`].
#'      - `geom_label`: A list of arguments for [`ggplot2::geom_label`].
#'   * `argsTheme`: A list of arguments for [`ggplot2::theme`].
#'
#' If 3D further arguments passed on the the RGL plotting functions. This must be done as
#'   lists. Currently the following arguments are supported:
#'
#'   * `argsAxes3d`: A list of arguments for [`rgl::axes3d`].
#'   * `argsPlot3d`: A list of arguments for [`rgl::plot3d`] to open the RGL window.
#'   * `argsTitle3d`: A list of arguments for [`rgl::title3d`].
#'   * `argsFaces`: A list of arguments for [`plotHull3D`].
#'   * `argsFeasible`: A list of arguments for RGL functions:
#'      - `points3d`: A list of arguments for [`rgl::points3d`].
#'      - `segments3d`: A list of arguments for [`rgl::segments3d`].
#'      - `triangles3d`: A list of arguments for [`rgl::triangles3d`].
#'   * `argsLabels`: A list of arguments for RGL functions:
#'      - `points3d`: A list of arguments for [`rgl::points3d`].
#'      - `text3d`: A list of arguments for [`rgl::text3d`].
#'   * `argsOptimum`: A list of arguments for RGL functions:
#'      - `points3d`: A list of arguments for [`rgl::points3d`].
#'
#' @note The feasible region defined by the constraints must be bounded (i.e. no extreme rays)
#'   otherwise you may see strange results.
#'
#' @return If 2D a `ggplot` object. If 3D a RGL window with the 3D plot.
#' @author Lars Relund \email{lars@@relund.dk}
#' @export
#' @import rgl ggplot2
#' @example inst/examples/ex_polytope.R
plotPolytope <- function(A,
            b,
            obj = NULL,
            type = rep("c", ncol(A)),
            nonneg = rep(TRUE, ncol(A)),
            crit = "max",
            faces = type,
            plotFaces = TRUE,
            plotFeasible = TRUE,
            plotOptimum = FALSE,
            latex = FALSE,
            labels = NULL, ...)
{
   if (length(type)!=ncol(A)) stop("Arg. 'type' must be same length as columns in A!")
   if (length(faces)!=ncol(A)) stop("Arg. 'faces' must be same length as columns in A!")
   if (!is.null(obj))
      if (length(obj)!=ncol(A))
         stop("Arg. 'obj' must have the same columns as in A and be a single criterion!")
   if (ncol(A)==2) return(plotPolytope2D(A, b, obj, type, nonneg, crit, faces, plotFaces,
                                         plotFeasible, plotOptimum, latex, labels, ...))
   if (ncol(A)==3) return(plotPolytope3D(A, b, obj, type, nonneg, crit, faces, plotFaces,
                                         plotFeasible, plotOptimum, latex, labels, ...))
   stop("Only 2D or 3D variables supported!")
}



#' Plot the polytope (bounded convex set) of a linear mathematical program
#'
#' @param A The constraint matrix.
#' @param b Right hand side.
#' @param obj A vector with objective coefficients.
#' @param type A character vector of same length as number of variables. If
#'   entry k is 'i' variable \eqn{k} must be integer and if 'c' continuous.
#' @param nonneg A boolean vector of same length as number of variables. If
#'   entry k is TRUE then variable k must be non-negative.
#' @param crit Either max or min (only used if add the iso-profit line)
#' @param faces A character vector of same length as number of variables. If
#'   entry k is 'i' variable \eqn{k} must be integer and if 'c' continuous.
#'   Useful if e.g. want to show the linear relaxation of an IP.
#' @param plotFaces If \code{True} then plot the faces.
#' @param plotFeasible If \code{True} then plot the feasible points/segments
#'   (relevant for ILP/MILP).
#' @param plotOptimum Show the optimum corner solution point (if alternative solutions
#'   only one is shown) and add the iso-profit line.
#' @param latex If \code{True} make latex math labels for TikZ.
#' @param labels If \code{NULL} don't add any labels. If 'n' no labels but show the points. If equal
#'   `coord` add coordinates to the points. Otherwise number all points from one.
#' @param ... Further arguments passed on the the `ggplot` plotting functions. This must be done as
#'   lists. Currently the following arguments are supported:
#'
#'   * `argsFaces`: A list of arguments for [`plotHull2D`].
#'   * `argsFeasible`: A list of arguments for `ggplotl2` functions:
#'      - `geom_point`: A list of arguments for [`ggplot2::geom_point`].
#'      - `geom_line`: A list of arguments for [`ggplot2::geom_line`].
#'   * `argsLabels`: A list of arguments for `ggplotl2` functions:
#'      - `geom_text`: A list of arguments for [`ggplot2::geom_text`].
#'   * `argsOptimum`:
#'      - `geom_point`: A list of arguments for [`ggplot2::geom_point`].
#'      - `geom_abline`: A list of arguments for [`ggplot2::geom_abline`].
#'      - `geom_label`: A list of arguments for [`ggplot2::geom_label`].
#'   * `argsTheme`: A list of arguments for [`ggplot2::theme`].
#'
#' @return A `ggplot` object.
#' @note In general use [plotPolytope()] instead of this function. The feasible region defined by the constraints must be bounded otherwise you may see
#'   strange results.
#' @seealso [plotPolytope()] for examples.
#' @author Lars Relund \email{lars@@relund.dk}
#' @import ggplot2
plotPolytope2D <-
   function(A,
            b,
            obj = NULL,
            type = rep("c", ncol(A)),
            nonneg = rep(TRUE, ncol(A)),
            crit = "max",
            faces = rep("c", ncol(A)),
            plotFaces = TRUE,
            plotFeasible = TRUE,
            plotOptimum = FALSE,
            latex = FALSE,
            labels = NULL,
            ...)
{
   args <- list(...)
   argsFaces <- mergeLists(list(argsGeom_polygon = list(fill = "gray90"), argsGeom_path = list(size = 0.5)), args$argsFaces)
   argsFeasible <- list()
   argsFeasible$geom_point <- mergeLists(list(), args$argsFeasible$geom_point)
   argsFeasible$geom_line <- mergeLists(list(), args$argsFeasible$geom_line)
   argsLabels <- list()
   argsLabels$geom_text <- mergeLists(list(size=3, color = "gray50", hjust = 1), args$argsLabels$geom_text)
   argsOptimum <- list()
   argsOptimum$geom_abline <- mergeLists(list(lty="dashed"), args$argsOptimum$geom_abline)
   argsOptimum$geom_point <- mergeLists(list(color = "red"), args$argsOptimum$geom_point)
   argsOptimum$geom_label <- mergeLists(list(nudge_x = 1, nudge_y = 0.5), args$argsOptimum$geom_label)
   argsTheme <- mergeLists(list(), args$argsTheme)

   if (!is.null(obj) & (!is.vector(obj) | !length(obj) == ncol(A)))
     stop("Arg. obj must be a vector of same length as the number of columns in A.")

   # Create solution plot
   p <- ggplot() #+ coord_fixed(ratio = 1)
   if (latex) p <- p + xlab("$x_1$") + ylab("$x_2$")
   if (!latex) p <- p + xlab(expression(x[1])) + ylab(expression(x[2]))

   if (plotFaces) {
      cPoints <- cornerPoints(A, b, faces, nonneg)
      cPoints <- unique(cPoints)
      p <- p +
         do.call(plotHull2D, args = c(list(cPoints, drawPlot = FALSE), argsFaces))
   }

   # find feasible points
   if (all(type == "c")) {
      points <- cornerPoints(A, b, type, nonneg)
   } else if (all(type == "i")) {
      points <- integerPoints(A, b, nonneg)
   } else {
      pl <- slices(A, b, type, nonneg)
      pl <- lapply(pl, unique)
      for (i in 1:length(pl)) pl[[i]] <- cbind(pl[[i]],i)
      pl <- lapply(pl, function(x) {
         colnames(x) <- c("x1", "x2", "g")
         rownames(x) <- NULL
         x<-data.frame(x)
      })
      points <- do.call(rbind, pl)
      #points <- points[,1:2]
   }
   points <- as.data.frame(points)

   if (plotFeasible) {
      if (all(type == "i")) {
         p <- p +
            do.call(geom_point, args = c(list(aes_string(x = 'x1', y = 'x2'), data=points), argsFeasible$geom_point))
      }
      if (length(which(type == "c"))==1) {
         p <- p +
            do.call(geom_line, args = c(list(aes_string(x = 'x1', y = 'x2', group='g'), data=points), argsFeasible$geom_line)) +
            do.call(geom_point, args = c(list(aes_string(x = 'x1', y = 'x2'), data=points), argsFeasible$geom_point))
      }
   }

   if (!is.null(labels)) {
      tmp <- points[,1:ncol(A)]
      tmp <- as.data.frame(tmp)
      if (labels == "coord")
         tmp$lbl <- df2String(tmp)
      else if (labels == "n")
         tmp$lbl <- ""
      else
         tmp$lbl <- 1:nrow(tmp)
      if (length(tmp$lbl)>0) {
         #p <- p + do.call(geom_point, args = c(list(aes_string(x = 'x1', y = 'x2'), data=tmp), argsLabels$geom_point))
         nudgeS=-(max(tmp$x1)-min(tmp$x1))/100
         #if (anyDuplicated(cbind(tmp$x1,tmp$x2), MARGIN = 1) > 0)
            p <- p +
               do.call(ggrepel::geom_text_repel,
                       args = c(list(aes_string(x = 'x1', y = 'x2', label = 'lbl'), data=tmp),
                                argsLabels$geom_text))
      }
   }

   if (plotOptimum) {
      if (!is.null(obj)) {    # add iso-profit line
         tmp <- points
         tmp$lbl <- df2String(tmp)
         tmp$z <- as.matrix(points[,1:2]) %*% obj
         if (crit=="max") i <- which.max(tmp$z)
         if (crit=="min") i <- which.min(tmp$z)
         if (latex) str <- paste0("$z^* = ", round(tmp$z[i], 2) , "$")
         if (!latex) str <- paste0("z* = ", round(tmp$z[i],2) )
         if (obj[2]!=0) {
            p <- p +
               do.call(geom_abline, args = c(list(intercept = tmp$z[i]/obj[2], slope = -obj[1]/obj[2]), argsOptimum$geom_abline))
         } else {
            p <- p +
               do.call(geom_abline, args = c(list(xintercept = tmp$x1[i]), argsOptimum))
         }
         p <- p +
            do.call(geom_point, args = c(list(aes_string(x = 'x1', y = 'x2'), data=tmp[i,]), argsOptimum$geom_point))
         p <- p +
            do.call(geom_label, args = c(list(aes_string(x = 'x1', y = 'x2', label = 'str'), data = tmp[i,]), argsOptimum$geom_label))
      }
   }

   p <- p +
      do.call(gMOIPTheme, args = c(list(), argsTheme))
   return(p)
}


#' Plot the polytope (bounded convex set) of a linear mathematical program
#'
#' @param A The constraint matrix.
#' @param b Right hand side.
#' @param obj A vector with objective coefficients.
#' @param type A character vector of same length as number of variables. If
#'   entry k is 'i' variable \eqn{k} must be integer and if 'c' continuous.
#' @param nonneg A boolean vector of same length as number of variables. If
#'   entry k is TRUE then variable k must be non-negative.
#' @param crit Either max or min (only used if add the iso-profit line)
#' @param faces A character vector of same length as number of variables. If
#'   entry k is 'i' variable \eqn{k} must be integer and if 'c' continuous.
#'   Useful if e.g. want to show the linear relaxation of an IP.
#' @param plotFaces If \code{True} then plot the faces.
#' @param plotFeasible If \code{True} then plot the feasible points/segments
#'   (relevant for ILP/MILP).
#' @param plotOptimum Show the optimum corner solution point (if alternative solutions
#'   only one is shown) and add the iso-profit line.
#' @param latex If \code{True} make latex math labels for TikZ.
#' @param labels If \code{NULL} don't add any labels. If 'n' no labels but show the points. If equal
#'   `coord` add coordinates to the points. Otherwise number all points from one.
#' @param ... Further arguments passed on the the RGL plotting functions. This must be done as
#'   lists. Currently the following arguments are supported:
#'
#'   * `argsAxes3d`: A list of arguments for [`rgl::axes3d`].
#'   * `argsPlot3d`: A list of arguments for [`rgl::plot3d`] to open the RGL window.
#'   * `argsTitle3d`: A list of arguments for [`rgl::title3d`].
#'   * `argsFaces`: A list of arguments for [`plotHull3D`].
#'   * `argsFeasible`: A list of arguments for RGL functions:
#'      - `points3d`: A list of arguments for [`rgl::points3d`].
#'      - `segments3d`: A list of arguments for [`rgl::segments3d`].
#'      - `triangles3d`: A list of arguments for [`rgl::triangles3d`].
#'   * `argsLabels`: A list of arguments for RGL functions:
#'      - `points3d`: A list of arguments for [`rgl::points3d`].
#'      - `text3d`: A list of arguments for [`rgl::text3d`].
#'   * `argsOptimum`: A list of arguments for RGL functions:
#'      - `points3d`: A list of arguments for [`rgl::points3d`].
#'
#' @note In general use [plotPolytope()] instead of this function. The feasible region defined by the constraints must be bounded otherwise you may see
#'   strange results.
#' @seealso [plotPolytope()] for examples.
#' @return A RGL window with 3D plot.
#' @author Lars Relund \email{lars@@relund.dk}
#' @import rgl
plotPolytope3D <-
   function(A,
            b,
            obj = NULL,
            type = rep("c", ncol(A)),
            nonneg = rep(TRUE, ncol(A)),
            crit = "max",
            faces = rep("c", ncol(A)),
            plotFaces = TRUE,
            plotFeasible = TRUE,
            plotOptimum = FALSE,
            latex = FALSE,
            labels = NULL, ...)
   {
      # set plot parameters
      args <- list(...)
      argsPlot3d <- mergeLists(list(xlab = '', box = FALSE, axes = FALSE), args$argsPlot3d)
      argsFaces <- mergeLists(list(drawPoints = FALSE, drawLines = TRUE), args$argsFaces)
      argsFeasible <- list()
      argsFeasible$points3d <- mergeLists(list(size = 3), args$argsFeasible$points3d)
      argsFeasible$segments3d <- mergeLists(list(lwd=1), args$argsFeasible$segments3d)
      argsFeasible$triangles3d <- mergeLists(list(col="grey100", alpha=0.4, smooth = FALSE), args$argsFeasible$triangles3d)
      argsOptimum <- list()
      argsOptimum$points3d <- mergeLists(list(col="red", size = 7), args$argsOptimum$points3d)
      argsLabels <- list()
      argsLabels$points3d <- mergeLists(list(col="grey50", size = 7), args$argsLabels$points3d)
      argsLabels$text3d <- mergeLists(list(), args$argsLabels$text3d)
      argsAxes3d <- mergeLists(list(edges=c('x', 'y', 'z')), args$argsAxes3d)
      argsTitle3d <- mergeLists(list(xlab = 'x1', ylab = 'x2', zlab = 'x3'), args$argsTitle3d)


      do.call(plot3d, args = c(list(x = replicate(2, 1:3), type = 'n'), argsPlot3d))
      aspect3d("iso")

      if (plotFaces) { # plot faces
         vertices <- cornerPoints(A, b, faces, nonneg)
         do.call(plotHull3D, args = c(list(pts = vertices), argsFaces))
      }

      if (plotFeasible) {
         if (all(type == "c")) {
            # don't plot anything
         } else if (all(type == "i")) {
            iPoints <- integerPoints(A, b, nonneg)
            do.call(points3d, args = c(list(iPoints[,1:3]), argsFeasible$points3d))
         } else {
            pl <- slices(A, b, type, nonneg)
            for (i in 1:length(pl)) {
               mat <- pl[[i]]
               # print(mat)
               if (is.null(mat)) next
               if (nrow(mat) == 1) {
                  do.call(points3d, args = c(list(mat), argsFeasible$points3d))
               }
               if (nrow(mat) == 2) {
                  do.call(segments3d, args = c(list(mat), argsFeasible$segments3d))
                  #cyl <- cylinder3d(mat, radius = 0.015)
                  #shade3d(cyl, col = "black")
               }
               if (nrow(mat) == 3) do.call(triangles3d, args = c(list(mat), argsFeasible$triangles3d))
               if (nrow(mat) > 3) {
                  hull <- geometry::convhulln(mat, options = "QJ")
                  tri <- t(hull)
                  do.call(triangles3d, args = c(list(mat[tri,1], mat[tri,2], mat[tri,3]), argsFeasible$triangles3d))
               }
            }
         }
      }

      if (plotOptimum) {
         if (is.null(obj)) stop("You need to specify the objective coefficients when using argument plotOption = TRUE.")
         vertices <- cornerPoints(A, b, type, nonneg)
         val <- vertices[,1:3] %*% as.matrix(obj)
         idx <- which.max(val)
         val <- vertices[idx,1:3]
         do.call(points3d, args = c(list(val[1], val[2], val[3]), argsOptimum$points3d))
      }

      if (!is.null(labels)) {
         if (all(type == "c")) {
            points <- cornerPoints(A, b, type, nonneg)
         } else if (all(type == "i")) {
            points <- integerPoints(A, b, nonneg)
         } else {
            points <- slices(A, b, type, nonneg, collapse = TRUE)
         }
         points <- as.data.frame(points)
         rownames(points) <- NULL
         if ((length(which(type == "c"))<length(type) & length(which(type == "c"))>0) |
             (length(which(type == "c"))==length(type))) {
            #pch3d(iPoints[,1:3], col="black", cex = 0.1, pch = 16)
            do.call(points3d, args = c(list(points[,1:3]), argsLabels$points3d))
         }
         if (labels=="coord")
            points$lbl <- df2String(points)
         else if (labels == "n")
            points$lbl <- ""
         else points$lbl <- 1:nrow(points)
         do.call(text3d, args = c(list(points[,1:3], texts = points$lbl), argsLabels$text3d))
      }

      do.call(axes3d, args = argsAxes3d)
      do.call(title3d, args = argsTitle3d)
      rgl::highlevel()
   }


#' Merge two lists to one
#'
#' @param a First list.
#' @param b Second list.
mergeLists <- function (a,b) {
   c(a[setdiff(names(a), names(b))], b)
}


#' Plot the lines of a linear mathematical program (Ax = b)
#'
#' @param A The constraint matrix.
#' @param b Right hand side.
#' @param nonneg A boolean vector of same length as number of variables. If
#'   entry k is TRUE then variable k must be non-negative and the line is plotted too.
#' @param latex If \code{True} make latex math labels for TikZ.
#' @param ... Further arguments passed on the the `ggplot` plotting functions. This must be done as
#'   lists. Currently the following arguments are supported:
#'
#'   * `argsTheme`: A list of arguments for [`ggplot2::theme`].
#'
#' @return A `ggplot` object.
#' @note In general you will properly use [plotPolytope()] instead of this function.
#' @seealso [plotPolytope()].
#' @author Lars Relund \email{lars@@relund.dk}
#' @import ggplot2
plotLines2D <-
   function(A,
            b,
            nonneg = rep(TRUE, ncol(A)),
            latex = FALSE,
            ...)
{
      args <- list(...)
      argsTheme <- mergeLists(list(), args$argsTheme)

      # Create solution plot
      p <- ggplot() #+ coord_fixed(ratio = 1)
      if (latex) p <- p + xlab("$x_1$") + ylab("$x_2$")
      if (!latex) p <- p + xlab(expression(x[1])) + ylab(expression(x[2]))

      for (i in 1:nrow(A)) {
         if (A[i,1] == 0) {
            p <- p + geom_hline(yintercept = b[i])
            next
         }
         if (A[i,2] == 0) {
            p <- p + geom_vline(xintercept = b[i])
            next
         }
         p <- p + geom_abline(slope = -A[i,1]/A[i,2], intercept = b[i]/A[i,2])
      }
      if (nonneg[1]) p <- p + geom_vline(xintercept = 0)
      if (nonneg[2]) p <- p + geom_hline(yintercept = 0)

      l <- cornerPoints(A, b, nonneg = nonneg)
      p <- p + xlim(min(l[,1]), max(l[,1])) +
         ylim(min(l[,2]), max(l[,2])) +
         do.call(gMOIPTheme, args = c(list(), argsTheme))
      return(p)
}




#' Create a plot of the criterion space of a bi-objective problem
#'
#' @param A The constraint matrix.
#' @param b Right hand side.
#' @param obj A p x n matrix(one row for each criterion).
#' @param type A character vector of same length as number of variables. If
#'   entry k is 'i' variable \eqn{k} must be integer and if 'c' continuous.
#' @param nonneg A boolean vector of same length as number of variables. If
#'   entry k is TRUE then variable k must be non-negative.
#' @param crit Either max or min (only used if add the iso-profit line).
#' @param addTriangles Add search triangles defined by the non-dominated extreme
#'   points.
#' @param addHull Add the convex hull and the rays.
#' @param plotFeasible If \code{True} then plot the criterion points/slices.
#' @param latex If true make latex math labels for TikZ.
#' @param labels If \code{NULL} don't add any labels. If 'n' no labels but show the points. If equal
#'   `coord` add coordinates to the points. Otherwise number all points from one.
#'
#' @note Currently only points are checked for dominance. That is, for MILP
#'   models some nondominated points may in fact be dominated by a segment.
#' @return The `ggplot` object.
#' @author Lars Relund \email{lars@@relund.dk}
#' @export
#' @import ggplot2
#' @example inst/examples/ex_criterion.R
plotCriterion2D <- function(A,
                            b,
                            obj,
                            type = rep("c", ncol(A)),
                            nonneg = rep(TRUE, ncol(A)),
                            crit = "max",
                            addTriangles = FALSE,
                            addHull = TRUE,
                            plotFeasible = TRUE,
                            latex = FALSE,
                            labels = NULL)
{
   # Set Custom theme
   myTheme <- theme_bw()
   myTheme <-
      myTheme + theme(
         panel.grid.major = element_blank(),
         panel.grid.minor = element_blank(),
         panel.border = element_blank(),
         #axis.line = element_blank(),
         axis.line = element_line(
            colour = "black",
            size = 0.5,
            arrow = arrow(length = unit(0.3, "cm"))
         ),
         #axis.ticks = element_blank()
         #axis.text.x = element_text(margin = margin(r = 30))
         # axis.ticks.length = unit(0.5,"mm"),
         #aspect.ratio=4/3,
         legend.position = "none"
      )

   # First find all relevant points
   if (all(type == "i")) {
      points <- integerPoints(A, b, nonneg)
   } else if (all(type == "c")) {
      points <- cornerPoints(A, b, type, nonneg)
   } else {
      points <- slices(A, b, type, nonneg, collapse = TRUE)
   }
   points <- criterionPoints(points, obj, crit, labels)
   if (all(type == "c")) { # if cont then no non-ext
      points$nd[points$nd & !points$se] <- FALSE
      points$se[points$nd & !points$se] <- FALSE
      points$sne[points$nd & !points$se] <- FALSE
   }
   #dat <<- points  # hack to get data as data frame

   if (crit=="max") points <- points[order(-points$z2,-points$z1),]
   if (crit=="min") points <- points[order(points$z2,points$z1),]

   # Initialize plot
   p <- ggplot(points, aes_q(x = quote(z1), y = quote(z2), col = "grey10") )
   if (latex) p <- p + xlab("$z_1$") + ylab("$z_2$")
   if (!latex) p <- p + xlab(expression(z[1])) + ylab(expression(z[2]))

   # Add hull plus rays
   if (addHull) {
      tmp<-points[points$se & !duplicated(cbind(points$z1,points$z2), MARGIN = 1),]
      delta <- max( (max(points$z1)-min(points$z1))/10, (max(points$z2)-min(points$z2))/10 )
      if (crit=="max") {
         tmp<-rbind(tmp[1:2,],tmp,tmp[1,]) # add rows
         tmp$z1[1] <- min(points$z1) - delta
         tmp$z2[1] <- min(points$z2) - delta
         tmp$z1[2] <- min(points$z1) - delta
         tmp$z2[2] <- max(points$z2)
         tmp$z1[length(tmp$z1)] <- max(points$z1)
         tmp$z2[length(tmp$z1)] <- min(points$z2)- delta
      }
      if (crit=="min") {
         tmp<-rbind(tmp[1,],tmp,tmp[1:2,]) # add rows
         tmp$z1[1] <- max(points$z1) + delta
         tmp$z2[1] <- min(points$z2)
         tmp$z1[length(tmp$z1)-1] <- min(points$z1)
         tmp$z2[length(tmp$z1)-1] <- max(points$z2) + delta
         tmp$z1[length(tmp$z1)] <- max(points$z1) + delta
         tmp$z2[length(tmp$z1)] <- max(points$z2) + delta
      }
      p <- p + geom_polygon(fill="gray95", col = NA, data=tmp)
   }

   if (plotFeasible) {
      # plot feasible areas
      if (all(type == "c")) {
         idx <- grDevices::chull(points[,c("z1","z2")])
         p <- p + geom_polygon(data = points[idx,], aes_string(x = 'z1', y = 'z2'),
                               fill=NA, size = 0.5, linetype = 1, col="grey80", alpha=0.6)
         p <- p + geom_point(aes_string(colour = 'nd', shape = 'se'), data = points) +
            scale_colour_grey(start = 0.6, end = 0)
      } else if (all(type == "i")) {
         #iPoints <- integerPoints(A, b, nonneg)
         #iPoints <- criterionPoints(iPoints, obj, crit, labels)
         p <- p + geom_point(aes_string(colour = 'nd', shape = 'se'), data = points) +
            #    #coord_fixed(ratio = 1) +
            scale_colour_grey(start = 0.6, end = 0)
      } else {
         pl <- slices(A, b, type, nonneg)
         pl <- lapply(pl, function(x) {
            x <- x %*% t(obj)
            if (is.vector(x)) x <- matrix(x, nrow=1)
            x[,1:2, drop = FALSE]
         })
         idx <- lapply(pl, grDevices::chull)
         #idx <- lapply(idx, function(x) c(x, x[1]))
         pl <- mapply(function(x,y){
            x[y,,drop = FALSE]
         },
         x = pl, y = idx, SIMPLIFY = FALSE)
         for (i in 1:length(pl)) pl[[i]] <- cbind(pl[[i]],i)
         # pl <- lapply(pl, function(x) {
         #    colnames(x) <- c("z1", "z2", "g")
         #    rownames(x) <- NULL
         #    x<-data.frame(x)
         # })
         pl <- do.call(rbind, pl)
         pl <- as.data.frame(pl)
         colnames(pl) <- c("z1", "z2", "g")
         rownames(pl) <- NULL
         p <- p + geom_polygon(
            data = pl,
            alpha = 0.6, fill = "grey",
            aes_string(
               x = "z1",
               y = "z2",
               group = "g"
            )
         ) + scale_fill_identity()
         points$co <- points$se | points$sne
         p <- p + geom_point(aes_string(colour = 'co', shape = 'se'), data = points) +
            scale_colour_grey(start = 0.8, end = 0)
      }
   }
   # Add triangles
   if (addTriangles) {
      tmp<-points[points$se | points$sne,]
      if (length(tmp$z1)>1) { # triangles
         for (r in 1:(dim(tmp)[1] - 1)) {
            p <- p +
               geom_segment(
                  x = tmp$z1[r],
                  y = tmp$z2[r],
                  xend = tmp$z1[r + 1],
                  yend = tmp$z2[r + 1],
                  colour = "gray"
               ) +
               geom_segment(
                  x = tmp$z1[r],
                  y = tmp$z2[r],
                  xend = tmp$z1[r],
                  yend = tmp$z2[r + 1],
                  colour = "gray"
               ) +
               geom_segment(
                  x = tmp$z1[r],
                  y = tmp$z2[r + 1],
                  xend = tmp$z1[r + 1],
                  yend = tmp$z2[r + 1],
                  colour = "gray"
               )
         }
      }
   }

   nudgeC=-(max(points$z1)-min(points$z1))/100
   if (!is.null(labels) & anyDuplicated(round(cbind(points$z1,points$z2),10), MARGIN = 1) > 0)
      p <- p + ggrepel::geom_text_repel(aes_string(label = 'lbl'), size=3,
                                        colour = "gray50", data=points)
   if (!is.null(labels) & anyDuplicated(round(cbind(points$z1,points$z2),10), MARGIN = 1) == 0)
      p <- p + geom_text(aes_string(label = 'lbl'), nudge_x = nudgeC, nudge_y = nudgeC,
                         hjust=1, size=3, colour = "gray50", data=points)
   p <- p + myTheme
   return(p)
}


#' Help function to save the view angle for the RGL 3D plot
#'
#' @param fname The file name of the view.
#' @param overwrite Overwrite existing file.
#' @param print Print the view so can be copied to R code (no file is saved).
#'
#' @note Only save if the file name don't exists.
#' @return NULL (invisible).
#' @author Lars Relund \email{lars@@relund.dk}
#' @export
#' @examples
#' \donttest{
#' view <- matrix( c(-0.412063330411911, -0.228006735444069, 0.882166087627411, 0,
#' 0.910147845745087, -0.0574885793030262, 0.410274744033813, 0, -0.042830865830183,
#' 0.97196090221405, 0.231208890676498, 0, 0, 0, 0, 1), nc = 4)
#'
#' loadView(v = view)
#' A <- matrix( c(3, 2, 5, 2, 1, 1, 1, 1, 3, 5, 2, 4), nc = 3, byrow = TRUE)
#' b <- c(55, 26, 30, 57)
#' obj <- c(20, 10, 15)
#' plotPolytope(A, b, plotOptimum = TRUE, obj = obj, labels = "coord")
#'
#' # Try to modify the angle in the RGL window
#' saveView(print = TRUE)  # get the view angle to insert into R code
#' }
saveView <- function(fname = "view.RData", overwrite = FALSE, print = FALSE) {
   if (print) {
      view <- rgl::par3d()$userMatrix
      cat(paste0("view <- matrix( c(", paste0(view, collapse = ", "), "), nc = 4)"))
   } else if (!file.exists(fname) | overwrite) {
      view <- rgl::par3d()$userMatrix
      save(view, file = fname)
      message(paste0("RGL view saved to RData file ", fname, "."))
   }
   return(invisible(NULL))
}


#' Help function to load the view angle for the RGL 3D plot from a file or matrix
#'
#' @param fname The file name of the view.
#' @param v The view matrix.
#' @param clear Call [rgl::clear3d()].
#' @param close Call [rgl::close3d()].
#' @param zoom Zoom level.
#' @param ... Additional parameters passed to [rgl::view3d()].
#'
#' @return NULL
#' @author Lars Relund \email{lars@@relund.dk}
#' @export
#' @examples
#' \donttest{
#' view <- matrix( c(-0.412063330411911, -0.228006735444069, 0.882166087627411, 0,
#' 0.910147845745087, -0.0574885793030262, 0.410274744033813, 0, -0.042830865830183,
#' 0.97196090221405, 0.231208890676498, 0, 0, 0, 0, 1), nc = 4)
#'
#' loadView(v = view)
#' A <- matrix( c(3, 2, 5, 2, 1, 1, 1, 1, 3, 5, 2, 4), nc = 3, byrow = TRUE)
#' b <- c(55, 26, 30, 57)
#' obj <- c(20, 10, 15)
#' plotPolytope(A, b, plotOptimum = TRUE, obj = obj, labels = "coord")
#'
#' # Try to modify the angle in the RGL window
#' saveView(print = TRUE)  # get the view angle to insert into R code
#' }
loadView <- function(fname = "view.RData", v = NULL, clear = TRUE, close = FALSE, zoom = 1, ...) {
   if (clear) rgl::clear3d()
   if (close) rgl::close3d()
   if (!is.null(v)) {
      rgl::view3d(userMatrix = v, zoom = zoom, ...)
   } else {
      if (file.exists(fname)) {
         load(fname)
         rgl::view3d(userMatrix = v)
      } else {
         warning(paste0("Can't load view in file ", fname, "!"))
      }
   }
}



#' Create a plot of a discrete non-dominated set.
#'
#' @param points Data frame with non-dominated points.
#' @param crit Either max or min (only used if add the iso-profit line). A vector is currently not
#'   supported.
#' @param addTriangles Add search triangles defined by the non-dominated extreme points.
#' @param addHull Add the convex hull and the rays.
#' @param latex If true make latex math labels for TikZ.
#' @param labels If \code{NULL} don't add any labels. If 'n' no labels but show the points. If equal
#'   `coord` add coordinates to the points. Otherwise number all points from one.
#'
#' @note Currently only points are checked for dominance. That is, for MILP models some
#'   nondominated points may in fact be dominated by a segment.
#' @return The `ggplot` object.
#' @author Lars Relund \email{lars@@relund.dk}
#' @export
#' @import ggplot2
#' @examples
#' dat <- data.frame(z1=c(12,14,16,18,18,18,14,15,15), z2=c(18,16,12,4,2,6,14,14,16))
#' points <- addNDSet(dat, crit = "min", keepDom = TRUE)
#' plotNDSet2D(points, crit = "min", addTriangles = TRUE)
#' plotNDSet2D(points, crit = "min", addTriangles = FALSE)
#' plotNDSet2D(points, crit = "min", addTriangles = TRUE, addHull = FALSE)
#' points <- addNDSet(dat, crit = "max", keepDom = TRUE)
#' plotNDSet2D(points, crit = "max", addTriangles = TRUE)
#' plotNDSet2D(points, crit = "max", addHull = FALSE)
plotNDSet2D <- function(points,
                        crit,
                        addTriangles = FALSE,
                        addHull = TRUE,
                        latex = FALSE,
                        labels = NULL)
{
   # Set Custom theme
   myTheme <- theme_bw()
   myTheme <-
      myTheme + theme(
         panel.grid.major = element_blank(),
         panel.grid.minor = element_blank(),
         panel.border = element_blank(),
         #axis.line = element_blank(),
         axis.line = element_line(
            colour = "black",
            size = 0.5,
            arrow = arrow(length = unit(0.3, "cm"))
         ),
         #axis.ticks = element_blank()
         #axis.text.x = element_text(margin = margin(r = 30))
         # axis.ticks.length = unit(0.5,"mm"),
         #aspect.ratio=4/3,
         legend.position = "none"
      )
   if (crit=="max") points <- points[order(-points$z2,-points$z1),]
   if (crit=="min") points <- points[order(points$z2,points$z1),]
   # Initialize plot
   p <- ggplot(points, aes_q(x = quote(z1), y = quote(z2)) )
   if (latex) p <- p + xlab("$z_1$") + ylab("$z_2$")
   if (!latex) p <- p + xlab(expression(z[1])) + ylab(expression(z[2]))

   # Add hull plus rays
   if (addHull) {
      di <- .mToDirection(crit, 2)
      p <- p +
         plotHull2D(points[points$nd, c("z1", "z2")],
                    addRays = TRUE, drawLines = TRUE, drawBBoxHull = FALSE,
                    direction = di, drawPlot = FALSE)
   }

   # Add triangles
   if (addTriangles) {
      tmp<-points[points$se | points$sne,]
      if (length(tmp$z1)>1) { # triangles
         for (r in 1:(dim(tmp)[1] - 1)) {
            p <- p +
               geom_segment(
                  x = tmp$z1[r],
                  y = tmp$z2[r],
                  xend = tmp$z1[r + 1],
                  yend = tmp$z2[r + 1],
                  colour = "gray"
               ) +
               geom_segment(
                  x = tmp$z1[r],
                  y = tmp$z2[r],
                  xend = tmp$z1[r],
                  yend = tmp$z2[r + 1],
                  colour = "gray"
               ) +
               geom_segment(
                  x = tmp$z1[r],
                  y = tmp$z2[r + 1],
                  xend = tmp$z1[r + 1],
                  yend = tmp$z2[r + 1],
                  colour = "gray"
               )
         }
      }
   }

   p <- p + geom_point(aes_string(colour = 'nd', shape = 'se'), data = points) +
      #    #coord_fixed(ratio = 1) +
      scale_colour_grey(start = 0.6, end = 0)

   nudgeC=-(max(points$z1)-min(points$z1))/100
   if (!is.null(labels) &
       anyDuplicated(round(cbind(points$z1, points$z2), 10), MARGIN = 1) > 0)
      p <- p + ggrepel::geom_text_repel(
         aes_string(label = 'lbl'),
         size = 3,
         colour = "gray50",
         data = points
      )
   if (!is.null(labels) &
       anyDuplicated(round(cbind(points$z1, points$z2), 10), MARGIN = 1) == 0)
      p <-
      p + geom_text(
         aes_string(label = 'lbl'),
         nudge_x = nudgeC,
         nudge_y = nudgeC,
         hjust = 1,
         size = 3,
         colour = "gray50",
         data = points
      )
   p <- p + myTheme
   return(p)
}


#' Plot a rectangle defined by two corner points.
#'
#' The rectangle is defined by {x|a <= x <= b} where a is the minimum values and
#' b is the maximum values.
#'
#' @param a A vector of length 3.
#' @param b A vector of length 3.
#' @param ... Further arguments passed on the the RGL plotting functions. This must be done as
#'   lists (see examples). Currently the following arguments are supported:
#'
#'   * `argsPlot3d`: A list of arguments for [`rgl::plot3d`].
#'   * `argsSegments3d`: A list of arguments for [`rgl::segments3d`][rgl::points3d()].
#'   * `argsPolygon3d`: A list of arguments for [`rgl::polygon3d`].
#'   * `argsShade3d`: A list of arguments for [`rgl::shade3d`][rgl::mesh3d()].
#'
#' @return Object ids (invisible).
#' @export
#'
#' @examples
#' \donttest{
#' ini3D()
#' plotRectangle3D(c(0,0,0), c(1,1,1))
#' plotRectangle3D(c(1,1,1), c(4,4,3), drawPoints = TRUE, drawLines = FALSE,
#'            argsPlot3d = list(size=2, type="s", alpha=0.3))
#' ids <- plotRectangle3D(c(2,2,2), c(3,3,2.5), argsPolygon3d = list(alpha = 1) )
#' finalize3D()
#' # pop3d(id = ids) remove last object
#' }
plotRectangle3D <- function(a, b, ...) {
   args <- list(...)
   argsSegments3d <- mergeLists(list(), args$argsSegments3d)
   argsPlot3d <- mergeLists(list(), args$argsPlot3d)
   argsShade3d <- mergeLists(list(), args$argsShade3d)
   argsPolygon3d <- mergeLists(list(), args$argsPolygon3d)

   if (any(a == b)) stop("Define different min and max values!")
   if (is.matrix(a)) a <- as.vector(a)
   if (is.matrix(b)) b <- as.vector(b)
   x<-matrix(c(a,b), ncol=3, byrow = TRUE)
   x<-expand.grid(x=c(x[1,1],x[2,1]), y=c(x[1,2],x[2,2]), z=c(x[1,3],x[2,3]))
   argsP <- list()
   if (!is.null(args$drawPoints)) argsP$drawPoints <- args$drawPoints
   if (!is.null(args$drawLines)) argsP$drawLines <- args$drawLines
   lst <- do.call(plotHull3D, args = c(list(x), list(argsSegments3d = argsSegments3d,
                                              argsPlot3d = argsPlot3d, argsShade3d = argsShade3d,
                                              argsPolygon3d = argsPolygon3d), argsP))
   return(invisible(lst$ids))
}


#' Plot a polygon.
#'
#' @param pts Vertices.
#' @param useShade Plot shade of the polygon.
#' @param useLines Plot lines inside the polygon.
#' @param usePoints Plot point shapes inside the polygon.
#' @param useFrame Plot a frame around the polygon.
#' @param ... Further arguments passed on the RGL plotting functions. This must be done as
#'   lists (see examples). Currently the following arguments are supported:
#'
#'   * `argsShade`: A list of arguments for [`rgl::polygon3d`] (n > 4 vertices),
#'                  [rgl::triangles3d()] (n = 3 vertices) and [rgl::quads3d()] (n = 4 vertices)
#'                  if `useShade = TRUE`.
#'   * `argsFrame`: A list of arguments for [`rgl::lines3d`] if `useFrame = TRUE`.
#'   * `argsPoints`: A list of arguments for [`rgl::shade3d`] if `usePoints = TRUE`. It is important
#'                   to give a texture using `texture`. A texture can be set using [getTexture()].
#'   * `argsLines`: A list of arguments for [rgl::persp3d()] when `useLines = TRUE`. Moreover, the list
#'                  may contain `lines`: number of lines.
#'
#' @return Object ids (invisible).
#' @importFrom stats lm predict
#' @export
#'
#' @examples
#' \donttest{
#' pts <- data.frame(x = c(1,0,0,0.4), y = c(0,1,0,0.3), z = c(0,0,1,0.3))
#' pts <- data.frame(x = c(1,0,0), y = c(0,1,0), z = c(0,0,1))
#'
#' ini3D()
#' plotPolygon3D(pts)
#' finalize3D()
#'
#' ini3D()
#' plotPolygon3D(pts, argsShade = list(color = "red", alpha = 1))
#' finalize3D()
#'
#' ini3D()
#' plotPolygon3D(pts, useFrame = TRUE, argsShade = list(color = "red", alpha = 0.5),
#'               argsFrame = list(color = "green"))
#' finalize3D()
#'
#' ini3D()
#' plotPolygon3D(pts, useFrame = TRUE, useLines = TRUE, useShade = TRUE,
#'               argsShade = list(color = "red", alpha = 0.2),
#'               argsLines = list(color = "blue"))
#' finalize3D()
#'
#' ini3D()
#' ids <- plotPolygon3D(pts, usePoints = TRUE, useFrame = TRUE,
#'               argsPoints = list(texture = getTexture(pch = 16, cex = 20)))
#' finalize3D()
#' # pop3d(id = ids) # remove object again
#'
#' # In general you have to finetune size and numbers when you use textures
#' # Different pch
#' for (i in 0:3) {
#'   fname <- getTexture(pch = 15+i, cex = 30)
#'   ini3D(TRUE)
#'   plotPolygon3D(pts, usePoints = TRUE, argsPoints = list(texture = fname))
#'   finalize3D()
#' }
#'
#' # Size of pch
#' for (i in 1:4) {
#'   fname <- getTexture(pch = 15+i, cex = 10 * i)
#'   ini3D(TRUE)
#'   plotPolygon3D(pts, usePoints = TRUE, argsPoints = list(texture = fname))
#'   finalize3D()
#' }
#'
#' # Number of pch
#' fname <- getTexture(pch = 16, cex = 20)
#' for (i in 1:4) {
#'   ini3D(TRUE)
#'   plotPolygon3D(pts, usePoints = TRUE,
#'                 argsPoints = list(texture = fname, texcoords = rbind(pts$x, pts$y, pts$z)*5*i))
#'   finalize3D()
#' }
#' }
plotPolygon3D <- function(pts, useShade = TRUE, useLines = FALSE, usePoints = FALSE,
                          useFrame = TRUE, ...) {
   args <- list(...)
   # args <- list(argsPoints=argsPoints)
   argsShade <- mergeLists(list(color = "black", col = "grey40",
                                    lwd = 2, alpha = 0.2, fill = TRUE),
                               args$argsShade)
   argsFrame <- mergeLists(list(color = "black", lwd = 1, alpha = 0.8), args$argsFrame)
   argsPoints <- mergeLists(list(color = "grey", specular = "black", alpha = 0.5,
                                 texcoords = rbind(pts[,1], pts[,2], pts[,3])*25
                                 ),
                           args$argsPoints)
   argsLines <- mergeLists(list(color = "black", lwd = 1, alpha = 0.4, lines = 50, back = 'lines',
                                front = 'lines', lit = FALSE),
                           args$argsLines)

   ids <- NULL
   pts <- .checkPts(pts, p = 3, asDF = TRUE)
   if (dimFace(pts) != 2) stop("Vertices don't define a polygon!")
   if (usePoints) useShade = FALSE

   if (useShade) {
      if (nrow(pts) > 3) ids <- c(ids, do.call(rgl::polygon3d, args = c(list(pts), argsShade)))
      if (nrow(pts) == 3) ids <- c(ids, do.call(rgl::triangles3d, args = c(list(pts), argsShade)))
   }
   if (useFrame) {
      n <- length(pts[,1])
      nas <- n+1
      prev <- 0L
      loop <- integer()
      for (i in seq_along(nas)) {
         loop <- c(loop, if (i > 1) NA, (prev + 1L):(nas[i] - 1L), prev + 1L)
         prev <- nas[i]
      }
      res <- cbind(pts[loop,1], pts[loop,2], pts[loop,3])
      ids <- c(ids, do.call(rgl::lines3d, args = c(list(res), argsFrame)))
   }
   if (usePoints) {
      if (nrow(pts) > 3)
         poly <- do.call(rgl::polygon3d, args = c(list(pts, plot = FALSE), argsPoints))
      if (nrow(pts) == 3)
         poly <- mesh3d(vertices = rbind(pts[,1], pts[,2], pts[,3], 1), triangles = 1:3)
      # poly$texcoords <- argsPoints$texcoords
      ids <- c(ids, do.call(rgl::shade3d, args = c(list(poly), argsPoints)))
   }
   if (useLines) {
      if (!rgl::cur3d()) stop("Option useLines need an open RGL window!")
      limits <- rgl::par3d()$bbox
      m <- c(limits[1], limits[3], limits[5])
      M <- c(limits[2], limits[4], limits[6])
      if (all(m == M)) stop("The RGL window's bounding box is not valid!")
      # do the mesh
      x <- seq(m[1], M[1], length.out = argsLines$lines)
      y <- seq(m[2], M[2], length.out = argsLines$lines)
      xy <- expand.grid(x, y)
      colnames(pts) <- c("x", "y", "z")
      xy[sp::point.in.polygon(xy[,1], xy[,2], pts$x, pts$y) <= 0,] <- NA
      z <- predict(lm(z ~ x + y, data = pts), newdata = data.frame(x=xy[,1], y=xy[,2]))
      ids <- c(ids, do.call(rgl::persp3d, args = c(list(x, y, z, add = TRUE), argsLines)))
   }
   print(rgl::highlevel())
   return(invisible(ids))
}


#' Save a point symbol as a temporary file.
#'
#' @param pch Point number/symbol.
#' @param cex Point size
#' @param ...  Further arguments passed to `plot`.
#'
#' @return The file name.
#' @export
#'
#' @examples
#' \donttest{
#' # Pch shapes
#' generateRPointShapes<-function(){
#'    oldPar<-par()
#'    par(font=2, mar=c(0.5,0,0,0))
#'    y=rev(c(rep(1,6),rep(2,5), rep(3,5), rep(4,5), rep(5,5)))
#'    x=c(rep(1:5,5),6)
#'    plot(x, y, pch = 0:25, cex=1.5, ylim=c(1,5.5), xlim=c(1,6.5),
#'         axes=FALSE, xlab="", ylab="", bg="blue")
#'    text(x, y, labels=0:25, pos=3)
#'    par(mar=oldPar$mar,font=oldPar$font )
#' }
#' generateRPointShapes()
#'
#' getTexture()
#' }
getTexture <- function(pch = 16, cex = 10, ...) {
   filename <- tempfile(fileext = ".png")
   grDevices::png(filename = filename)
   graphics::plot(1, 1, pch = pch, cex = cex, axes=FALSE, xlab="", ylab="", ...)
   grDevices::dev.off()
   return(filename)
}



#' Plot a cone defined by a point in 3D.
#'
#' The cones are defined as the point plus R3+.
#'
#' @param pts A matrix with a point in each row.
#' @param drawPoint Draw the points defining the cone.
#' @param drawLines Draw lines of the cone.
#' @param drawPolygons Draw polygons of the cone.
#' @param rectangle Draw the cone as a rectangle.
#' @param direction Ray direction. If i'th entry is positive, consider the i'th column of `pts`
#'   plus a value greater than on equal zero (minimize objective $i$). If negative, consider the
#'   i'th column of `pts` minus a value greater than on equal zero (maximize objective $i$).
#' @param useRGLBBox Use the RGL bounding box as ray limits for the cone.
#' @param ... Further arguments passed on the the RGL plotting functions. This must be done as
#'   lists (see examples). Currently the following arguments are supported:
#'
#'   * `argsPlot3d`: A list of arguments for [`rgl::plot3d`].
#'   * `argsSegments3d`: A list of arguments for [`rgl::segments3d`][rgl::points3d()].
#'   * `argsPolygon3d`: A list of arguments for [`rgl::polygon3d`].
#'
#' @return Object ids (invisible).
#' @export
#'
#' @examples
#' \donttest{
#' ini3D(argsPlot3d = list(xlim = c(0,6), ylim = c(0,6), zlim = c(0,6)))
#' plotCones3D(c(4,4,4), drawLines = FALSE, drawPoint = TRUE,
#'            argsPlot3d = list(col = "red", size = 10),
#'            argsPolygon3d = list(alpha = 1), rectangle = TRUE)
#' plotCones3D(c(1,1,1), rectangle = FALSE)
#' plotCones3D(matrix(c(3,3,3,2,2,2), ncol = 3, byrow = TRUE))
#' finalize3D()
#'
#' ini3D(argsPlot3d = list(xlim = c(0,6), ylim = c(0,6), zlim = c(0,6)))
#' plotCones3D(c(4,4,4), direction = 1)
#' plotCones3D(c(2,2,2), direction = -1)
#' plotCones3D(c(4,2,2), direction = c(1,-1,-1))
#' ids <- plotCones3D(c(2,2,4), direction = c(-1,-1,1))
#' finalize3D()
#' # pop3d(id = ids) # remove last cone
#' }
plotCones3D <-
   function(pts,
            drawPoint = TRUE,
            drawLines = TRUE,
            drawPolygons = TRUE,
            direction = 1,
            rectangle = FALSE,
            useRGLBBox = TRUE,
            ...) {
   args <- list(...)
   argsPlot3d <- mergeLists(list(), args$argsPlot3d)
   argsSegments3d <- mergeLists(list(lwd = 1, col = "grey40"), args$argsSegments3d)
   argsPolygon3d <- mergeLists(list(), args$argsPolygon3d)

   pts <- .checkPts(pts, p = 3)
   for (i in 1:dim(pts)[1]) {
      pt <- as.vector(pts[i, ])
      lst <- plotHull3D(pt, drawPoints = drawPoint,
                   drawLines = drawLines, drawPolygons = drawPolygons,
                   addRays = TRUE, direction = direction, drawBBoxHull = rectangle,
                   useRGLBBox = useRGLBBox,
                   argsPlot3d = argsPlot3d, argsSegments3d = argsSegments3d,
                   argsPolygon3d = argsPolygon3d
                 )
   }
   return(invisible(lst$ids))
   }



#' Plot a cone defined by a point in 2D.
#'
#' The cones are defined as the point plus/minus rays of R2.
#'
#' @param pts A matrix with a point in each row.
#' @param drawPoint Draw the points defining the cone.
#' @param drawLines Draw lines of the cone.
#' @param drawPolygons Draw polygons of the cone.
#' @param rectangle Draw the cone as a rectangle.
#' @param direction Ray direction. If i'th entry is positive, consider the i'th column of `pts`
#'   plus a value greater than on equal zero (minimize objective $i$). If negative, consider the
#'   i'th column of `pts` minus a value greater than on equal zero (maximize objective $i$).
#' @param drawPlot Draw the `ggplot`. Set to FALSE if you want to combine hulls in a single plot.
#' @param m Minimum values of the bounding box.
#' @param M Maximum values of the bounding box.
#' @param ... Further arguments passed to [plotHull2D]
#'
#' @return A `ggplot` object
#' @export
#' @import ggplot2
#' @examples
#' library(ggplot2)
#' plotCones2D(c(4,4), drawLines = FALSE, drawPoint = TRUE,
#'            argsGeom_point = list(col = "red", size = 10),
#'            argsGeom_polygon = list(alpha = 0.5), rectangle = TRUE)
#' plotCones2D(c(1,1), rectangle = FALSE)
#' plotCones2D(matrix(c(3,3,2,2), ncol = 2, byrow = TRUE))
#'
#' ## The Danish flag
#' lst <- list(argsGeom_polygon = list(alpha = 0.85, fill = "red"),
#'             drawPlot = FALSE, drawPoint = FALSE, drawLines = FALSE)
#' p1 <- do.call(plotCones2D, args = c(list(c(2,4), direction = 1), lst))
#' p2 <- do.call(plotCones2D, args = c(list(c(1,2), direction = -1), lst))
#' p3 <- do.call(plotCones2D, args = c(list(c(2,2), direction = c(1,-1)), lst))
#' p4 <- do.call(plotCones2D, args = c(list(c(1,4), direction = c(-1,1)), lst))
#' ggplot() + p1 + p2 + p3 + p4 + theme_void()
plotCones2D <-
   function(pts,
            drawPoint = TRUE,
            drawLines = TRUE,
            drawPolygons = TRUE,
            direction = 1,
            rectangle = FALSE,
            drawPlot = TRUE,
            m = apply(pts,2,min)-5,
            M = apply(pts,2,max)+5,
            ...) {
   pts <- .checkPts(pts, p = 2)
   lst <- NULL
   for (i in 1:dim(pts)[1]) {
      pt <- as.vector(pts[i, ])
      lst <- c(lst,
               plotHull2D(pt, drawPoints = drawPoint,
                 drawLines = drawLines, drawPolygons = drawPolygons,
                 addRays = TRUE, direction = direction, drawPlot = FALSE,
                 drawBBoxHull = rectangle, m = m, M = M, ...))
   }
   if (drawPlot) lst <- ggplot() + lst
   return(lst)
}


#' Plot the convex hull of a set of points in 3D.
#'
#' @param pts A matrix with a point in each row.
#' @param drawPoints Draw the points.
#' @param drawLines Draw lines of the facets.
#' @param drawPolygons Fill the facets.
#' @param addText Add text to the points. Currently `coord` (coordinates), `rownames` (rownames)
#'   and `both` supported or a vector with text.
#' @param addRays Add the ray defined by `direction`.
#' @param useRGLBBox Use the RGL bounding box when add rays.
#' @param direction Ray direction. If i'th entry is positive, consider the i'th column of `pts`
#'   plus a value greater than on equal zero (minimize objective $i$). If negative, consider the
#'   i'th column of `pts` minus a value greater than on equal zero (maximize objective $i$).
#' @param drawBBoxHull If `addRays` then draw the hull areas hitting the bounding box also.
#' @param ... Further arguments passed on the the RGL plotting functions. This must be done as
#'   lists (see examples). Currently the following arguments are supported:
#'
#'   * `argsPlot3d`: A list of arguments for [`rgl::plot3d`].
#'   * `argsSegments3d`: A list of arguments for [`rgl::segments3d`][rgl::points3d()].
#'   * `argsPolygon3d`: A list of arguments for [`rgl::polygon3d`].
#'   * `argsShade3d`: A list of arguments for [`rgl::shade3d`][rgl::mesh3d()].
#'   * `argsText3d`: A list of arguments for [`rgl::text3d`].
#'
#' @return A list with hull, `pts` classified and object ids (invisible).
#' @export
#'
#' @examples
#' \donttest{
#' ini3D()
#' pts<-matrix(c(0,0,0), ncol = 3, byrow = TRUE)
#' plotHull3D(pts) # a point
#' pts<-matrix(c(1,1,1,2,2,2,3,3,3), ncol = 3, byrow = TRUE)
#' plotHull3D(pts, drawPoints = TRUE) # a line
#' pts<-matrix(c(1,0,0,1,1,1,1,2,2,3,1,1,3,3,3), ncol = 3, byrow = TRUE)
#' plotHull3D(pts, drawLines = FALSE, argsPolygon3d = list(alpha=0.6)) # a polygon
#' pts<-matrix(c(5,5,5,10,10,5,10,5,5,5,5,10), ncol = 3, byrow = TRUE)
#' lst <- plotHull3D(pts, argsPolygon3d = list(alpha=0.9), argsSegments3d = list(color="red"))
#' finalize3D()
#' # pop3d(id = lst$ids) # remove last hull
#'
#' ## Using addRays
#' pts <- data.frame(x = c(1,3), y = c(1,3), z = c(1,3))
#' ini3D(argsPlot3d = list(xlim = c(0,max(pts$x)+10),
#'   ylim = c(0,max(pts$y)+10),
#'   zlim = c(0,max(pts$z)+10)))
#' plotHull3D(pts, drawPoints = TRUE, addRays = TRUE, , drawBBoxHull = FALSE)
#' plotHull3D(c(4,4,4), drawPoints = TRUE, addRays = TRUE)
#' finalize3D()
#'
#' pts <- data.frame(x = c(4,2.5,1), y = c(1,2.5,4), z = c(1,2.5,4))
#' ini3D(argsPlot3d = list(xlim = c(0,max(pts$x)+10),
#'   ylim = c(0,max(pts$y)+10),
#'   zlim = c(0,max(pts$z)+10)))
#' plotHull3D(pts, drawPoints = TRUE, addRays = TRUE)
#' finalize3D()
#'
#' pts <- matrix(c(
#'   0, 4, 8,
#'   0, 8, 4,
#'   8, 4, 0,
#'   4, 8, 0,
#'   4, 0, 8,
#'   8, 0, 4,
#'   4, 4, 4,
#'   6, 6, 6
#'   ), ncol = 3, byrow = TRUE)
#' ini3D(FALSE, argsPlot3d = list(xlim = c(min(pts[,1])-2,max(pts[,1])+10),
#'   ylim = c(min(pts[,2])-2,max(pts[,2])+10),
#'   zlim = c(min(pts[,3])-2,max(pts[,3])+10)))
#' plotHull3D(pts, drawPoints = TRUE, addText = "coord")
#' plotHull3D(pts, addRays = TRUE)
#' finalize3D()
#'
#' pts <- genNDSet(3, 100, dubND = FALSE)
#' pts <- as.data.frame(pts[,1:3])
#'
#' ini3D(argsPlot3d = list(
#'   xlim = c(0,max(pts[,1])+10),
#'   ylim = c(0,max(pts[,2])+10),
#'   zlim = c(0,max(pts[,3])+10)))
#' plotHull3D(pts, drawPoints = TRUE, addRays = TRUE)
#' finalize3D()
#'
#' ini3D(argsPlot3d = list(
#'   xlim = c(0,max(pts[,1])+10),
#'   ylim = c(0,max(pts[,2])+10),
#'   zlim = c(0,max(pts[,3])+10)))
#' plotHull3D(pts, drawPoints = TRUE, drawPolygons = TRUE, addText = "coord", addRays = TRUE)
#' finalize3D()
#'
#' ini3D(argsPlot3d = list(
#'   xlim = c(0,max(pts[,1])+10),
#'   ylim = c(0,max(pts[,2])+10),
#'   zlim = c(0,max(pts[,3])+10)))
#' plotHull3D(pts, drawPoints = TRUE, drawLines = FALSE,
#'   argsPolygon3d = list(alpha = 1), addRays = TRUE)
#' finalize3D()
#'
#' ini3D(argsPlot3d = list(
#'   xlim = c(0,max(pts[,1])+10),
#'   ylim = c(0,max(pts[,2])+10),
#'   zlim = c(0,max(pts[,3])+10)))
#' plotHull3D(pts, drawPoints = TRUE, argsPolygon3d = list(color = "red"), addRays = TRUE)
#' plotCones3D(pts, argsPolygon3d = list(alpha = 1), rectangle = TRUE)
#' finalize3D()
#' }
plotHull3D <- function(pts,
                       drawPoints = FALSE,
                       drawLines = TRUE,
                       drawPolygons = TRUE,
                       addText = FALSE,
                       addRays = FALSE,
                       useRGLBBox = TRUE,
                       direction = 1,
                       drawBBoxHull = TRUE,
                       ...)
{
   args <- list(...)
   argsSegments3d <- mergeLists(list(lwd = 1, col = "grey40"), args$argsSegments3d)
   argsPlot3d <- mergeLists(list(size = 5, col = "black"), args$argsPlot3d)
   argsShade3d <- mergeLists(list(col = "grey100"), args$argsShade3d)
   argsPolygon3d <- mergeLists(list(color = "grey100", lwd = 2, alpha = 0.1),
                               args$argsPolygon3d)
   argsText3d <- mergeLists(list(), args$argsText3d)

   ids <- NULL
   pts <- .checkPts(pts, p = 3)
   if (length(direction) != 3) direction = rep(direction[1],3)
   hull <- convexHull(pts, addRays = addRays, direction = direction, useRGLBBox = useRGLBBox)
   set <- hull$pts
   hull <- hull$hull
   d <- dimFace(set[,1:3, drop = FALSE])
   if (d==3) { # then poly define facets
      poly <- hull
      v <- 1:3
      names(v) <- colnames(set[,1:3])
      pN <- purrr::map_dfc(v, function(i) if (sign(direction[i]) > 0) max(set[,i]) else min(set[,i]))
      #colnames(pN) <- colnames(set[,1:3])
      for (i in 1:dim(poly)[1]) {
         tri <- poly[i,!is.na(poly[i,])]
         pt <- set[tri,1:4]
         tri <- 1:nrow(pt)
         if (drawLines) {
            if (length(tri)>3) { # then have to find the vertex sequence
               tri <- convexHull(pt[,1:3])$hull
            }
            tri1 <- NULL # use tri1 to include the addRays=TRUE case
            for (j in 2:length(tri)){
               if (!(pt[tri[j-1],4] == 0 && pt[tri[j],4] == 0)) #  && !drawBBoxHull
                  tri1 <- c(tri1,NA,tri[j-1],tri[j])
            }
            if (!(pt[tri[1], 4] == 0 && pt[tri[length(tri)], 4] == 0)) {  #  && !drawBBoxHull
               tri1 <- c(tri1, NA, tri[1], tri[length(tri)])
            }
            ids <- c(ids, do.call(rgl::polygon3d, args = c(list(pt[tri1, 1], pt[tri1, 2], pt[tri1, 3], fill = FALSE),
                    argsSegments3d)))
         }
         if (drawPolygons) {
            if (drawBBoxHull | (!drawBBoxHull & nrow(dplyr::intersect(pt[,1:3],pN)) == 0)) {
               if (length(tri)==3) {
                  ids <- c(ids, do.call(rgl::triangles3d, args = c(list(x = pt[,1:3]), argsPolygon3d) ))
               } else if (length(tri)==4){
                  tri <- convexHull(pt[,1:3])$hull # then have to find the vertex sequence
                  obj <- rgl::qmesh3d(t(pt[,1:3]),tri, homogeneous = FALSE)
                  ids <- c(ids, do.call(rgl::shade3d, args = c(list(obj), argsPolygon3d)))
               } else {
                  # idx <- apply(pt[,1:3], 2, function(x) {return(length(unique(x))==1)})
                  # idx<-which(!idx)
                  # if (length(idx)==3) coords <- 1:2 else coords <- idx
                  #plotHull3D(pt[tri,1:3])
                  tri <- convexHull(pt[,1:3])$hull
                  comb <- t(utils::combn(3,2))
                  for (j in 1:3) {
                     res <- try(
                        ids <- c(ids, do.call(rgl::polygon3d, args = c(list(
                        pt[tri, 1], pt[tri, 2], pt[tri, 3], fill = TRUE, coords = comb[j,]
                     ), argsPolygon3d))), silent = TRUE)
                     if (!inherits(res, "try-error")) break #else {cat(i, " ", j, " ")}
                  }
               }
            }
         }
      }
   }
   if (d==0) {
      ids <- c(ids, do.call(rgl::plot3d, args = c(list(set[,1], set[,2], set[,3], add=TRUE), argsPlot3d)))
   }
   if (d==1) { #segments3d(set, col = "black", lwd=10, smooth= TRUE)
      cyl <- rgl::cylinder3d(set, radius = 0.01)
      ids <- c(ids, do.call(rgl::shade3d, args = c(list(cyl), argsShade3d)))
   }
   if (d==2) {
      idx <- apply(set[,1:3], 2, function(x) {return(length(unique(x))==1)})
      tri <- hull
      if (length(tri)==3) {
         if (drawPolygons) ids <- c(ids, do.call(rgl::triangles3d,
              args = c(list(set[tri,1], set[tri,2], set[tri,3]), argsPolygon3d) ))
         if (drawLines) ids <- c(ids, do.call(rgl::polygon3d,
              args = c(list(set[tri,1], set[tri,2], set[tri,3], coords = which(!idx), fill= FALSE),
                       argsSegments3d) ))
      } else {
         if (drawPolygons) ids <- c(ids, do.call(rgl::polygon3d,
              args = c(list(set[tri,1], set[tri,2], set[tri,3], coords = which(!idx), fill= TRUE),
                       argsPolygon3d) ))
         if (drawLines) ids <- c(ids, do.call(rgl::polygon3d,
              args = c(list(set[tri,1], set[tri,2], set[tri,3], coords = which(!idx), fill= FALSE),
                       argsSegments3d) ))
      }
   }
   if (drawPoints) ids <- c(ids, do.call(plotPoints3D, args = c(list(pts), list(argsPlot3d = argsPlot3d))))


   if (length(addText) > 1) {
      if (length(addText) == nrow(pts)) {
         do.call(rgl::text3d, args = c(list(x = pts, texts = addText), argsText3d) )
      }
   } else if (length(addText) == 1) {
      if (addText == TRUE | addText == "coord") {
         text <- paste0("(",pts[,1],",",pts[,2],",",pts[,3],")")
         ids <- c(ids, do.call(rgl::text3d, args = c(list(x = pts, texts = text), argsText3d) ))
      } else if (addText == "rownames") {
         text <- rownames(as.data.frame(pts))
         ids <- c(ids, do.call(rgl::text3d, args = c(list(x = pts, texts = text), argsText3d) ))
      } else if (addText == "both") {
         text <- paste0("(",pts[,1],",",pts[,2],",",pts[,3],")")
         text <- paste(text,rownames(as.data.frame(pts)))
         ids <- c(ids, do.call(rgl::text3d, args = c(list(x = pts, texts = text), argsText3d) ))
      } else if (addText != FALSE & length(addText) == nrow(pts)) {
         ids <- c(ids, do.call(rgl::text3d, args = c(list(x = pts, texts = addText), argsText3d) ))
      }
   }
   return(invisible(list(hull = hull, pts = set, ids = ids)))
   stop("Error in plotHull3D")
}


#' Plot points in 3D.
#'
#' @param pts A vector or matrix with the points.
#' @param addText Add text to the points. Currently `coord` (coordinates), `rownames` (rownames)
#'   and `both` supported or a vector with the text.
#' @param ... Further arguments passed on the the RGL plotting functions. This must be done as
#'   lists (see examples). Currently the following arguments are supported:
#'
#'   * `argsPlot3d`: A list of arguments for [`rgl::plot3d`].
#'   * `argsPch3d`: A list of arguments for [`rgl::pch3d`].
#'   * `argsText3d`: A list of arguments for [`rgl::text3d`].
#'
#' @return Object ids (invisible).
#' @export
#'
#' @examples
#' \donttest{
#' ini3D()
#' pts<-matrix(c(1,1,1,5,5,5), ncol = 3, byrow = TRUE)
#' plotPoints3D(pts)
#' plotPoints3D(c(2,3,3), argsPlot3d = list(col = "red", size = 10))
#' plotPoints3D(c(3,2,3), argsPlot3d = list(col = "blue", size = 10, type="p"))
#' plotPoints3D(c(1.5,1.5,1.5), argsPlot3d = list(col = "blue", size = 10, type="p"))
#' plotPoints3D(c(2,2,2, 1,1,1), addText = "coord")
#' ids <- plotPoints3D(c(3,3,3, 4,4,4), addText = "rownames")
#' finalize3D()
#' rgl::rglwidget()
#' # pop3d(ids) # remove the last again
#' }
plotPoints3D <- function(pts, addText = FALSE, ...) {
   args <- list(...)
   argsPlot3d <- mergeLists(list(size = 5, col = "black", type="p"), args$argsPlot3d)
   argsPch3d <- mergeLists(list(litt = TRUE), args$argsPch3d)
   argsText3d <- mergeLists(list(), args$argsText3d)

   ids <- NULL
   p <- pts
   if (is.vector(pts)) p <- matrix(pts, ncol = 3, byrow = TRUE)
   p <- as.matrix(p)

   if (is.null(args$argsPch3d)) ids <- c(ids, do.call(rgl::plot3d, args = c(list(p, add= TRUE), argsPlot3d) ))
   if (!is.null(args$argsPch3d)) ids <- c(ids, do.call(rgl::pch3d, args = c(list(p), argsPch3d) ))

   if (length(addText) > 1) {
      if (length(addText) == nrow(p)) {
         ids <- c(ids, do.call(rgl::text3d, args = c(list(x = p, texts = addText), argsText3d) ))
      }
   } else if (length(addText) == 1) {
      if (addText == TRUE | addText == "coord") {
         text <- paste0("(",p[,1],",",p[,2],",",p[,3],")")
         ids <- c(ids, do.call(rgl::text3d, args = c(list(x = p, texts = text), argsText3d) ))
      } else if (addText == "rownames") {
         text <- rownames(as.data.frame(p))
         ids <- c(ids, do.call(rgl::text3d, args = c(list(x = p, texts = text), argsText3d) ))
      } else if (addText == "both") {
         text <- paste0("(",p[,1],",",p[,2],",",p[,3],")")
         text <- paste(text,rownames(as.data.frame(p)))
         ids <- c(ids, do.call(rgl::text3d, args = c(list(x = p, texts = text), argsText3d) ))
      } else if (addText != FALSE & length(addText) == nrow(p)) {
         ids <- c(ids, do.call(rgl::text3d, args = c(list(x = p, texts = addText), argsText3d) ))
      }
   }
   return(invisible(ids))
}


#' Plot a plane in 3D.
#'
#' @param normal Normal to the plane.
#' @param point A point on the plane.
#' @param offset The offset of the plane (only used if `point = NULL`).
#' @param useShade Plot shade of the plane.
#' @param useLines Plot lines inside the plane.
#' @param usePoints Plot point shapes inside the plane.
#' @param ... Further arguments passed on the the RGL plotting functions. This must be done as
#'   lists (see examples). Currently the following arguments are supported:
#'
#'   * `argsPlanes3d`: A list of arguments for [rgl::planes3d()] used when `useShade = TRUE`.
#'   * `argsLines`: A list of arguments for [rgl::persp3d()] when `useLines = TRUE`. Moreover, the list
#'                  may contain `lines`: number of lines.
#'
#' @return NULL (invisible)
#' @export
#'
#' @examples
#' \donttest{
#' ini3D(argsPlot3d = list(xlim = c(-1,10), ylim = c(-1,10), zlim = c(-1,10)) )
#' plotPlane3D(c(1,1,1), point = c(1,1,1))
#' plotPoints3D(c(1,1,1))
#' plotPlane3D(c(1,2,1), point = c(2,2,2), argsPlanes3d = list(color="red"))
#' plotPoints3D(c(2,2,2))
#' plotPlane3D(c(2,1,1), offset = -6, argsPlanes3d = list(color="blue"))
#' plotPlane3D(c(2,1,1), argsPlanes3d = list(color="green"))
#' finalize3D()
#'
#' ini3D(argsPlot3d = list(xlim = c(-1,10), ylim = c(-1,10), zlim = c(-1,10)) )
#' plotPlane3D(c(1,1,1), point = c(1,1,1), useLines = TRUE, useShade = TRUE)
#' ids <- plotPlane3D(c(1,2,1), point = c(2,2,2), argsLines = list(col="blue", lines = 100),
#'             useLines = TRUE)
#' finalize3D()
#' # pop3d(id = ids) # remove last plane
#' }
plotPlane3D <- function(normal, point = NULL, offset = 0, useShade = TRUE, useLines = FALSE,
                        usePoints = FALSE, ...) {
   args <- list(...)
   argsPlanes3d <- mergeLists(list(col = "grey", alpha = 0.2), args$argsPlanes3d)
   argsLines <- mergeLists(list(back = 'lines', front = 'lines', add = TRUE, lines = 50),
                           args$argsLines)

   ids <- NULL
   if (!is.null(point)) offset <- -sum(normal * point)
   if (useShade) {
      ids <- c(ids, do.call(rgl::planes3d, args = c(list(a = normal, d = offset), argsPlanes3d) ))
   }
   # else use points or lines
   if (!rgl::cur3d()) stop("Option useLines or usePoints need an open rgl window!")
   limits <- rgl::par3d()$bbox
   m <- c(limits[1], limits[3], limits[5])
   M <- c(limits[2], limits[4], limits[6])
   # do the mesh
   x <- seq(m[1], M[1], length.out = argsLines$lines)
   y <- seq(m[2], M[2], length.out = argsLines$lines)
   f <- function(x,y){-(normal[1]/normal[3])*x - (normal[2]/normal[3])*y - offset/normal[3]}
   z <- outer(x,y,f)
   z[z < m[3] | z > M[3]] <- NA
   if (useLines) {
      ids <- c(ids, do.call(rgl::persp3d, args = c(list(x, y, z), argsLines) ))
   }
   if (usePoints) {
      return(invisible(ids))
      # do.call(rgl::plot3d, args = c(list(x, y, z, add = TRUE), argsPoints) )
      # persp3d(x, y, z, back = 'lines', front = 'lines', add = TRUE)
      #
      #
      # if (nrow(pts) > 3) poly <- do.call(rgl::polygon3d, args = c(list(pts, plot = FALSE),
      #                                                             argsPolygon3d))
      # if (nrow(pts) == 3) {
      #    # poly <- do.call(rgl::triangles3d, args = c(list(pts, plot = FALSE), argsPolygon3d))
      #    poly <- tmesh3d(rbind(pts[,1], pts[,2], pts[,3], 1), indices = 1:3)
      # }
      # poly$texcoords <- argsPolygon3d$texcoords
      # do.call(rgl::shade3d, args = c(list(poly, col = "white", specular = "black"), argsPolygon3d))
   }
   return(invisible(ids))
}


#' The `ggplot` theme for the package
#'
#' @param ... Further arguments parsed to [ggplot2::theme()].
#'
#' @return The theme object.
#' @export
#'
#' @examples
#' pts <- matrix(c(1,1), ncol = 2, byrow = TRUE)
#' plotHull2D(pts)
#' pts1 <- matrix(c(2,2, 3,3), ncol = 2, byrow = TRUE)
#' pts2 <- matrix(c(1,1, 2,2, 0,1), ncol = 2, byrow = TRUE)
#' ggplot2::ggplot() +
#'   plotHull2D(pts2, drawPoints = TRUE, addText = "coord", drawPlot = FALSE) +
#'   plotHull2D(pts1, drawPoints = TRUE, drawPlot = FALSE) +
#'   gMOIPTheme() +
#'   ggplot2::xlab(expression(x[1])) +
#'   ggplot2::ylab(expression(x[2]))
gMOIPTheme <- function(...) {
   return(
      theme_bw() +
      theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
            panel.border = element_blank(),
            #axis.line = element_blank(),
            axis.line = element_line(colour = "black", size = 0.5,
                                     arrow = arrow(length = unit(0.3,"cm")) ),
            #axis.ticks = element_blank()
            #axis.text.x = element_text(margin = margin(r = 30))
            # axis.ticks.length = unit(0.5,"mm"),
            #aspect.ratio=4/3,
            legend.position="none",
            ...
      )
   )
}


#' Initialize the RGL window.
#'
#' @param new A new window is opened (otherwise the current is cleared).
#' @param clear Clear the current RGL window.
#' @param ... Further arguments passed on the the RGL plotting functions. This must be done as
#'   lists. Currently the following arguments are supported:
#'
#'   * `argsPlot3d`: A list of arguments for [`rgl::plot3d`].
#'   * `argsAspect3d`: A list of arguments for [`rgl::aspect3d`].
#'
#' @return NULL (invisible).
#' @export
#'
#' @examples
#' \donttest{
#' ini3D()
#' pts<-matrix(c(1,1,1,5,5,5), ncol = 3, byrow = TRUE)
#' plotPoints3D(pts)
#' finalize3D()
#'
#' lim <- c(-1, 7)
#' ini3D(argsPlot3d = list(xlim = lim, ylim = lim, zlim = lim))
#' plotPoints3D(pts)
#' finalize3D()
#' }
ini3D <- function(new = TRUE, clear = FALSE, ...){
   args <- list(...)
   # args <- list(argsPlot3d = list(xlim = c(0,6), ylim = c(0,6), zlim = c(0,6)))
   argsPlot3d <-
      mergeLists(list(
         xlab = '',
         ylab = '',
         zlab = '',
         box = FALSE,
         axes = FALSE
      ), args$argsPlot3d)
   argsAspect3d <- mergeLists(list(x = "iso"), args$argsAspect3d)

   if (new) rgl::open3d()
   if (clear) rgl::clear3d()
   do.call(rgl::plot3d, args = c(list(x = c(1,1), y = c(1,1), z = c(1,1), type = 'n'), argsPlot3d))
   do.call(rgl::aspect3d, args = argsAspect3d)
   rgl::highlevel()
   return(invisible(NULL))
}

#' Finalize the RGL window.
#'
#' @param ... Further arguments passed on the the RGL plotting functions. This must be done as
#'   lists. Currently the following arguments are supported:
#'
#'   * `argsAxes3d`: A list of arguments for [`rgl::axes3d`].
#'   * `argsTitle3d`: A list of arguments for [`rgl::title3d`][rgl::axes3d].
#'
#' @return The RGL object (using [rgl::highlevel()]).
#' @export
#'
#' @examples
#' \donttest{
#' ini3D()
#' pts<-matrix(c(1,1,1,5,5,5), ncol = 3, byrow = TRUE)
#' plotPoints3D(pts)
#' finalize3D()
#'
#' ini3D()
#' pts<-matrix(c(1,1,1,5,5,5), ncol = 3, byrow = TRUE)
#' plotPoints3D(pts)
#' finalize3D(argsAxes3d = list(edges = "bbox"))
#' }
finalize3D <- function(...){
   args <- list(...)
   argsAxes3d <- mergeLists(list(edges = c('x', 'y', 'z')), args$argsAxes3d)
   argsTitle3d <- mergeLists(list(xlab = expression(italic(z)[1]), ylab = expression(italic(z)[2]),
                                  zlab = expression(italic(z)[3])), args$argsTitle3d)

   do.call(rgl::axes3d, args = argsAxes3d)
   do.call(rgl::title3d, args = argsTitle3d)
   rgl::highlevel()
}


#' Convert LaTeX to a png file
#'
#' @param tex TeX string. Remember to escape backslash with \\.
#' @param viewPng View the result in the plots window.
#' @param width Width of the png.
#' @param height Height of the png (`width` are ignored).
#' @param dpi Dpi of the png. Not used if `width` or `height` are specified.
#' @param fontsize Front size used in the LaTeX document.
#' @param calcM Estimate 1 em in pixels in the resulting png.
#' @param crop Call command line program `pdfcrop` (must be installed).
#'
#' @return The filename of the png or a list if `calcM = TRUE`.
#' @export
#'
#' @examples
#' \dontrun{
#' tex <- "$\\mathbb{R}_{\\geqq}$"
#' texToPng(tex, viewPng = TRUE)
#' texToPng(tex, fontsize = 20, viewPng = TRUE)
#' texToPng(tex, height = 50, fontsize = 10, viewPng = TRUE)
#' texToPng(tex, height = 50, fontsize = 50, viewPng = TRUE)
#' tex <- "MMM"
#' texToPng(tex, dpi=72, calcM = TRUE)
#' texToPng(tex, width = 100, calcM = TRUE)
#' f <- texToPng(tex, dpi=300)
#' pngSize(f)
#' }
texToPng <- function(tex, width = NULL, height = NULL, dpi = 72, viewPng = FALSE, fontsize = 12,
                  calcM = FALSE, crop = FALSE) {
   texFile <- tempfile(fileext=".tex")
   pdfFile <- gsub("[.]tex", ".pdf", texFile)
   pngFile <- gsub("[.]tex", ".png", texFile)
   writeLines(c(paste0("\\documentclass[class=scrreprt, fontsize = ", fontsize, "pt]{standalone}"),
                "\\usepackage{amsfonts}",
                "\\usepackage{amssymb}",
                "\\begin{document}",
                tex,
                "\\end{document}"),
              texFile)
   system(paste0("pdflatex -shell-escape -output-directory=", tempdir(), " ", texFile))
   if (crop) system(paste("pdfcrop", pdfFile, pdfFile))
   tmp <- pdftools::pdf_pagesize(pdfFile)
   z <- 1
   if (!is.null(width)) z <- dpi/72 * tmp$width/width
   if (!is.null(height)) z <- dpi/72 * tmp$height/height
   dpi <- dpi * 1/z
   if (calcM) em1 <- .sizeM(dpi = dpi, fontsize = fontsize)$w
   pdftools::pdf_convert(pdfFile, format = "png", dpi = dpi, antialias = TRUE, filenames = pngFile, verbose = FALSE)
   if (viewPng) {
      img <- png::readPNG(pngFile)
      grid::grid.newpage()
      grid::grid.raster(img, just = "center", height = 0.2)
   }
   if (calcM) {
      s <- pngSize(pngFile)
      emPng <- s$w/em1
      return(list(png=pngFile, emPx = em1, emLength = emPng, w = s$w, h = s$h))
   }
   return(pngFile)
}


#' Estimate 1 em in pixels in the resulting png.
#'
#' @param ... Arguments parsed on to `texToPng`.
#'
#' @return The width and size of the png.
#' @keywords internal
.sizeM <- function(...) {
   f <- texToPng("M", ...)
   return(pngSize(f))
}


#' To size of the png file.
#'
#' @param png Png file name.
#'
#' @return A list with width and height.
pngSize <- function(png) {
   img <- png::readPNG(png)
   w <- dim(img)[2]
   h <- dim(img)[1]
   return(list(w=w,h=h))
}


#' Draw boxes, axes and other text outside the data using TeX strings.
#'
#' @param main The main title for the plot.
#' @param sub The subtitle for the plot.
#' @param xlab The axis labels for the plot .
#' @param ylab The axis labels for the plot .
#' @param zlab The axis labels for the plot .
#' @param line The ``line'' of the plot margin to draw the label on.
#' @param ... Additional parameters which are passed to \code{\link{plotMTeX3D}}.
#'
#' @details The rectangular prism holding the 3D plot has 12 edges.  They are identified
#' using 3 character strings.  The first character (`x', `y', or `z') selects
#' the direction of the axis.  The next two characters are each `-' or `+',
#' selecting the lower or upper end of one of the other coordinates.  If only
#' one or two characters are given, the remaining characters default to `-'.
#' For example \code{edge = 'x+'} draws an x-axis at the high level of y and the
#' low level of z.
#'
#' By default, \code{axes3d} uses the \code{\link{bbox3d}} function to draw the axes.
#' The labels will move so that they do not obscure the data.  Alternatively,
#' a vector of arguments as described above may be used, in which case
#' fixed axes are drawn using \code{axis3d}.
#'
#' If \code{pos} is a numeric vector of length 3, \code{edge} determines
#' the direction of the axis and the tick marks, and the values of the
#' other two coordinates in \code{pos} determine the position.  See the
#' examples.
#' @return The object IDs of objects added to the scene.
#' @export
#'
#' @examples
#' \dontrun{
#' ini3D(argsPlot3d = list(xlim = c(0, 2), ylim = c(0, 2), zlim = c(0, 2)))
#' plotTitleTeX3D(main = "\\LaTeX", sub = "subtitle $\\alpha$",
#'                xlab = "$x^1_2$", ylab = "$\\beta$", zlab = "$x\\cdot y$")
#' finalize3D()
#' }
plotTitleTeX3D <- function (main = NULL, sub = NULL, xlab = NULL, ylab = NULL,
                            zlab = NULL, line = NA, ...) {
   save <- rgl::par3d(skipRedraw = TRUE, ignoreExtent = TRUE)
   on.exit(rgl::par3d(save))
   result <- numeric(0)
   if (!is.null(main)) {
      aline <- ifelse(is.na(line), 2, line)
      result <- c(result, main = plotMTeX3D(main, "x++",
                                            line = aline, ...))
   }
   if (!is.null(sub)) {
      aline <- ifelse(is.na(line), 3, line)
      result <- c(result, sub = plotMTeX3D(sub, "x", line = aline,
                                           ...))
   }
   if (!is.null(xlab)) {
      aline <- ifelse(is.na(line), 2, line)
      result <- c(result, xlab = plotMTeX3D(xlab, "x", line = aline,
                                            ...))
   }
   if (!is.null(ylab)) {
      aline <- ifelse(is.na(line), 2, line)
      result <- c(result, ylab = plotMTeX3D(ylab, "y", line = aline,
                                            ...))
   }
   if (!is.null(zlab)) {
      aline <- ifelse(is.na(line), 2, line)
      result <- c(result, zlab = plotMTeX3D(zlab, "z", line = aline,
                                            ...))
   }
   lowlevel(result)
}


#' Get ranges of the bounding box margins
#'
#' @param expand Expand margins.
#' @param ranges The bounding box.
#'
#' @return A list with ranges.
#' @keywords internal
.getRanges <- function (expand = 1.03, ranges = par3d("bbox"))
{
   ranges <- list(xlim = ranges[1:2], ylim = ranges[3:4], zlim = ranges[5:6])
   strut <- FALSE
   ranges <- lapply(ranges, function(r) {
      d <- diff(r)
      if (d > 0)
         return(r)
      strut <<- TRUE
      if (d < 0)
         return(c(0, 1))
      else if (r[1] == 0)
         return(c(-1, 1))
      else return(r[1] + 0.4 * abs(r[1]) * c(-1, 1))
   })
   ranges$strut <- strut
   ranges$x <- (ranges$xlim - mean(ranges$xlim)) * expand +
      mean(ranges$xlim)
   ranges$y <- (ranges$ylim - mean(ranges$ylim)) * expand +
      mean(ranges$ylim)
   ranges$z <- (ranges$zlim - mean(ranges$zlim)) * expand +
      mean(ranges$zlim)
   ranges
}


#' Plot TeX in the margin
#'
#' @param tex TeX string
#' @param edge The position at which to draw the axis or text.
#' @param line The ``line'' of the plot margin to draw the label on.
#' @param at The value of a coordinate at which to draw the axis.
#' @param pos  The position at which to draw the axis or text.
#' @param ... Further arguments passed to [plotTeX3D()].
#'
#' @return The object IDs of objects added to the scene.
#' @export
plotMTeX3D <- function (tex, edge, line = 0, at = NULL, pos = NA, ...) {
   save <- rgl::par3d(ignoreExtent = TRUE)
   on.exit(rgl::par3d(save))
   ranges <- .getRanges()
   edge <- c(strsplit(edge, "")[[1]], "-", "-")[1:3]
   coord <- match(toupper(edge[1]), c("X", "Y", "Z"))
   if (coord == 2)
      edge[1] <- edge[2]
   else if (coord == 3)
      edge[1:2] <- edge[2:3]
   range <- ranges[[coord]]
   if (is.null(at))
      at <- mean(range)
   newlen <- max(length(tex), length(line), length(at))
   tex <- rep(tex, len = newlen)
   line <- rep(line, len = newlen)
   at <- rep(at, len = newlen)
   if (all(is.na(pos))) {
      pos <- matrix(NA, 3, length(at))
      if (edge[1] == "+")
         pos[1, ] <- ranges$x[2]
      else pos[1, ] <- ranges$x[1]
      if (edge[2] == "+")
         pos[2, ] <- ranges$y[2]
      else pos[2, ] <- ranges$y[1]
      if (edge[3] == "+")
         pos[3, ] <- ranges$z[2]
      else pos[3, ] <- ranges$z[1]
   }
   else pos <- matrix(pos, 3, length(at))
   pos[coord, ] <- at
   ticksize <- 0.05 * (pos[, 1] - c(mean(ranges$x), mean(ranges$y),
                                    mean(ranges$z)))
   ticksize[coord] <- 0
   plotTeX3D(pos[1, ] + 3 * ticksize[1] * line,
             pos[2, ] + 3 * ticksize[2] * line,
             pos[3, ] + 3 * ticksize[3] * line, tex,
             ...)
}



#' Plot TeX at a position.
#'
#' @param x Coordinate.
#' @param y Coordinate.
#' @param z Coordinate.
#' @param cex Expansion factor (you properly have to fine tune it).
#' @param fixedSize Fix the size of the object (no scaling when zoom).
#' @param tex TeX string.
#' @param size Size of the generated png.
#' @param ... Arguments passed on to [rgl::sprites3d()] and [texToPng()].
#'
#' @return The shape ID of the displayed object is returned.
#' @export
#'
#' @examples
#' \dontrun{
#' tex0 <- "$\\mathbb{R}_{\\geqq}$"
#' tex1 <- "\\LaTeX"
#' tex2 <- "This is a title"
#' ini3D(argsPlot3d = list(xlim = c(0, 2), ylim = c(0, 2), zlim = c(0, 2)))
#' plotTeX3D(0.75,0.75,0.75, tex0)
#' plotTeX3D(0.5,0.5,0.5, tex0, cex = 2)
#' plotTeX3D(1,1,1, tex2)
#' finalize3D()
#' ini3D(new = TRUE, argsPlot3d = list(xlim = c(0, 200), ylim = c(0, 200), zlim = c(0, 200)))
#' plotTeX3D(75,75,75, tex0)
#' plotTeX3D(50,50,50, tex1)
#' plotTeX3D(100,100,100, tex2)
#' finalize3D()
#' }
plotTeX3D <- function (x, y, z, tex, cex = graphics::par("cex"), fixedSize = FALSE, size = 480, ...) {
   f <- texToPng(tex, width = size, calcM = TRUE, ...)
   # expand png so same width and height
   system(paste0("convert ", f$png,
                 " -background none -gravity center -extent ", max(f$w,f$h), "x", max(f$w,f$h), " ", f$png))
   tmp <- rgl::par3d()$bbox
   radius <- 1/2 * cex * f$emLength * 20/f$emPx * max(c(tmp[2]-tmp[1], tmp[4]-tmp[3],tmp[6]-tmp[5]) * rgl::par3d()$scale)
   rgl::sprites3d(x, y, z, texture = f$png,
             textype = "rgba", col = "white", lit = FALSE,
             radius = radius, fixedSize = fixedSize, ...)
}
Any scripts or data that you put into this service are public.
gMOIP documentation built on May 31, 2023, 8:45 p.m.
rdrr.io home R language documentation Run R code online
CRAN packages Bioconductor packages R-Forge packages GitHub packages
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
gMOIP
Tools for 2D and 3D Plots of Single and Multi-Objective Linear/Integer Programming Models

R/plot.R
In gMOIP: Tools for 2D and 3D Plots of Single and Multi-Objective Linear/Integer Programming Models

Defines functions pngSize .sizeM texToPng finalize3D ini3D gMOIPTheme plotPlane3D plotPoints3D plotHull3D plotCones2D plotCones3D getTexture plotPolygon3D plotRectangle3D plotNDSet2D loadView saveView plotCriterion2D plotLines2D plotPolytope3D plotPolytope2D plotPolytope plotHull2D

Documented in finalize3D getTexture gMOIPTheme ini3D loadView plotCones2D plotCones3D plotCriterion2D plotHull2D plotHull3D plotLines2D plotNDSet2D plotPlane3D plotPoints3D plotPolygon3D plotPolytope plotPolytope2D plotPolytope3D plotRectangle3D pngSize saveView .sizeM texToPng

Try the gMOIP package in your browser

R Package Documentation

Browse R Packages

We want your feedback!

gMOIP Tools for 2D and 3D Plots of Single and Multi-Objective Linear/Integer Programming Models

R/plot.R In gMOIP: Tools for 2D and 3D Plots of Single and Multi-Objective Linear/Integer Programming Models

Defines functions pngSize .sizeM texToPng finalize3D ini3D gMOIPTheme plotPlane3D plotPoints3D plotHull3D plotCones2D plotCones3D getTexture plotPolygon3D plotRectangle3D plotNDSet2D loadView saveView plotCriterion2D plotLines2D plotPolytope3D plotPolytope2D plotPolytope plotHull2D

Documented in finalize3D getTexture gMOIPTheme ini3D loadView plotCones2D plotCones3D plotCriterion2D plotHull2D plotHull3D plotLines2D plotNDSet2D plotPlane3D plotPoints3D plotPolygon3D plotPolytope plotPolytope2D plotPolytope3D plotRectangle3D pngSize saveView .sizeM texToPng

Try the gMOIP package in your browser

R Package Documentation

Browse R Packages

We want your feedback!

gMOIP
Tools for 2D and 3D Plots of Single and Multi-Objective Linear/Integer Programming Models

R/plot.R
In gMOIP: Tools for 2D and 3D Plots of Single and Multi-Objective Linear/Integer Programming Models
