Add surfaceplot, isosurface (support 3d plots) (#217)

Co-authored-by: Johnny Chen <johnnychen94@hotmail.com>

.gitignore

@@ -9,3 +9,4 @@
 .DS_Store
 
 docs/imgs/
+.old/

Project.toml

@@ -9,12 +9,19 @@ version = "2.7.0"
 Contour = "d38c429a-6771-53c6-b99e-75d170b6e991"
 Crayons = "a8cc5b0e-0ffa-5ad4-8c14-923d3ee1735f"
 Dates = "ade2ca70-3891-5945-98fb-dc099432e06a"
+LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
+MarchingCubes = "299715c1-40a9-479a-aaf9-4a633d36f717"
+NaNMath = "77ba4419-2d1f-58cd-9bb1-8ffee604a2e3"
 SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
+StaticArrays = "90137ffa-7385-5640-81b9-e52037218182"
 StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91"
 
 [compat]
 Contour = "0.5"
 Crayons = "4.1"
+MarchingCubes = "0.1"
+NaNMath = "0.3"
+StaticArrays = "0.12, 1"
 StatsBase = "0.32, 0.33"
 julia = "1.6"
 

docs/generate_docs.jl

@@ -55,6 +55,20 @@ function main()
     contourplot1 = ("Contourplot", "contourplot(-3:.01:3, -7:.01:3, (x, y) -> exp(-(x / 2)^2 - ((y + 2) / 4)^2), border=:dotted)"),
     heatmap1 = ("Heatmap", """heatmap(repeat(collect(0:10)', outer=(11, 1)), zlabel="z")"""),
     heatmap2 = ("Heatmap", "heatmap(collect(0:30) * collect(0:30)', xfact=.1, yfact=.1, xoffset=-1.5, colormap=:inferno)"),
+    surfaceplot1 = ("Surfaceplot", """
+      sombrero(x, y) = 15sinc(√(x^2 + y^2) / π)
+      surfaceplot(-8:.5:8, -8:.5:8, sombrero, colormap=:jet, border=:dotted)
+      """),
+    surfaceplot2 = ("Surfaceplot", """
+      surfaceplot(
+        -8:.5:8, -8:.5:8, (x, y) -> 15sinc(√(x^2 + y^2) / π),
+        zscale=z -> 0, lines=true, azimuth=-90, elevation=90, colormap=:jet, border=:dotted
+      )
+      """),
+    isosurface = ("Isosurface", """
+      torus(x, y, z, r=0.2, R=0.5) = (√(x^2 + y^2) - R)^2 + z^2 - r^2
+      isosurface(-1:.1:1, -1:.1:1, -1:.1:1, torus, cull=true, zoom=2, elevation=50, border=:dotted)
+      """),
     width = ("Width", "lineplot(sin, 1:.5:20, width=60, border=:dotted)"),
     height = ("Height", "lineplot(sin, 1:.5:20, height=18, border=:dotted)"),
     labels = ("Labels", "lineplot(sin, 1:.5:20, labels=false, border=:dotted)"),
@@ -99,7 +113,7 @@ function main()
 
   examples = NamedTuple{keys(exs)}((
       MD(Paragraph(
-        "```julia\n$(rstrip(e[2]))\n```\n![$(e[1])]($docs_url/doc/imgs/$ver/$k.png)"
+        "```julia\n$(rstrip(e[2]))\n```\n![$(e[1])]($docs_url/$ver/$k.png)"
       )) for (k, e) in pairs(exs)
     )
   )
@@ -158,6 +172,8 @@ Here is a list of the main high-level functions for common scenarios:
   - [`densityplot`](https://github.com/JuliaPlots/UnicodePlots.jl#density-plot) (Density Plot)
   - [`contourplot`](https://github.com/JuliaPlots/UnicodePlots.jl#contour-plot) (Contour Plot)
   - [`heatmap`](https://github.com/JuliaPlots/UnicodePlots.jl#heatmap-plot) (Heatmap Plot)
+  - [`surfaceplot`](https://github.com/JuliaPlots/UnicodePlots.jl#surface-plot) (Surface Plot - 3D)
+  - [`isosurface`](https://github.com/JuliaPlots/UnicodePlots.jl#isosurface-plot) (Isosurface Plot - 3D)
 
 Here is a quick hello world example of a typical use-case:
 
@@ -258,6 +274,26 @@ The `zlabel` option and `zlabel!` method may be used to set the `z` axis (colorb
 
 $(examples.heatmap2)
 
+#### Surface Plot
+
+Plot a colored surface using height values `z` above a `x-y` plane, in three dimensions (masking values using `NaN`s is supported).
+
+$(examples.surfaceplot1)
+
+Use `lines=true` to increase the density (underlying call to `lineplot` instead of `scatterplot`).
+To plot a slice in 3D, use an anonymous function which maps to a constant value: `zscale=z -> a_constant`:
+
+$(examples.surfaceplot2)
+
+#### Isosurface Plot
+
+Uses [`MarchingCubes.jl`](https://github.com/t-bltg/MarchingCubes.jl) to extract an isosurface, where `isovalue` controls the surface isovalue.
+Using `centroid` enables plotting the triangulation centroids instead of the triangle vertices (better for small plots).
+Back face culling (hide not visible facets) can be activated using `cull=true`.
+One can use the legacy 'Marching Cubes' algorithm using `legacy=true`.
+
+$(examples.isosurface)
+
 ### Options
 
 All plots support the set (or a subset) of the following named parameters:
@@ -266,6 +302,14 @@ All plots support the set (or a subset) of the following named parameters:
 
 _Note_: If you want to print the plot into a file but have monospace issues with your font, you should probably try setting `border=:ascii` and `canvas=AsciiCanvas` (or `canvas=DotCanvas` for scatterplots).
 
+### 3D plots
+
+3d plots use a so-called "Matrix-View-Projection" transformation matrix `MVP` on input data to project 3D plots to a 2D screen.
+Use keywords`elevation`, `azimuth`, `up` or `zoom` to control the "View" matrix, a.k.a., camera.
+The `projection` type for `MVP` can be set to either `:perspective` or `orthographic`.
+Displaying the `x-`, `y-`, and `z-` axes can be controlled using the `axes3d` keyword.
+For better resolution, use wider and taller `Plot` size, which can be also be achieved using the unexported `UnicodePlots.default_size!(width=60)` for all future plots.
+
 ### Methods
 
   - `title!(plot::Plot, title::String)`
@@ -373,6 +417,7 @@ Inspired by [TextPlots.jl](https://github.com/sunetos/TextPlots.jl), which in tu
   mkpath("imgs/$ver")
   open("imgs/gen_imgs.jl", "w") do io
     println(io, """
+      # WARNING: this file has been automatically generated, please update UnicodePlots/docs/generate_docs.jl instead
       using UnicodePlots, StableRNGs, SparseArrays
       include(joinpath(dirname(pathof(UnicodePlots)), "..", "test", "fixes.jl"))
 
@@ -393,11 +438,13 @@ Inspired by [TextPlots.jl](https://github.com/sunetos/TextPlots.jl), which in tu
   open("imgs/gen_imgs.sh", "w") do io
     write(io, """
       #!/usr/bin/env bash
+      # WARNING: this file has been automatically generated, please update UnicodePlots/docs/generate_docs.jl instead
       \${JULIA-julia} gen_imgs.jl
 
       for f in $ver/*.txt; do
         html=\${f%.txt}.html
         cat \$f | \${ANSI2HTML-ansi2html} --input-encoding=utf-8 --output-encoding=utf-8 >\$html
+        sed -i "s,background-color: #000000,background-color: #1b1b1b," \$html
         \${WKHTMLTOIMAGE-wkhtmltoimage} --quiet --crop-w 800 --quality 85 \$html \${html%.html}.png
       done
 
@@ -410,7 +457,7 @@ Inspired by [TextPlots.jl](https://github.com/sunetos/TextPlots.jl), which in tu
   plain_readme = plain(readme)
   write(stdout, plain_readme)
   open("../README.md", "w") do io
-    write(io, "<!-- WARNING: this file has been automatically generated, please update Unicodeplots/docs/generate_docs.jl instead, and run \$ julia generate_docs.jl to render README.md !! -->\n")
+    write(io, "<!-- WARNING: this file has been automatically generated, please update UnicodePlots/docs/generate_docs.jl instead, and run \$ julia generate_docs.jl to render README.md !! -->\n")
     write(io, plain_readme)
   end
   return

src/UnicodePlots.jl

@@ -4,6 +4,11 @@ using Dates
 using Crayons
 using StatsBase: Histogram, fit, percentile
 using SparseArrays: AbstractSparseMatrix, findnz
+using LinearAlgebra
+using StaticArrays
+
+import MarchingCubes
+import NaNMath
 import Contour
 
 export GraphicsArea,
@@ -50,6 +55,10 @@ export GraphicsArea,
     scatterplot!,
     contourplot,
     contourplot!,
+    surfaceplot,
+    surfaceplot!,
+    isosurface,
+    isosurface!,
     stairs,
     stairs!,
     histogram,
@@ -59,6 +68,7 @@ export GraphicsArea,
     spy,
     boxplot,
     boxplot!,
+    MVP,
     savefig
 
 include("common.jl")
@@ -77,13 +87,16 @@ include("canvas/dotcanvas.jl")
 include("canvas/heatmapcanvas.jl")
 
 include("description.jl")
+include("volume.jl")
 
 include("plot.jl")
 include("colormaps.jl")
 include("interface/barplot.jl")
 include("interface/histogram.jl")
 include("interface/scatterplot.jl")
 include("interface/contourplot.jl")
+include("interface/surfaceplot.jl")
+include("interface/isosurface.jl")
 include("interface/lineplot.jl")
 include("interface/stairs.jl")
 include("interface/densityplot.jl")

src/canvas/lookupcanvas.jl

@@ -64,10 +64,10 @@ function CreateLookupCanvas(
         visible,
         pixel_width,
         pixel_height,
-        Float64(origin_x),
-        Float64(origin_y),
-        Float64(width),
-        Float64(height),
+        float(origin_x),
+        float(origin_y),
+        float(width),
+        float(height),
         xscale,
         yscale,
     )

src/colormaps.jl

@@ -1420,8 +1420,10 @@ function rgb2ansi(rgb)
 end
 
 function cmapcolor(z, minz, maxz, cmap)
-    i = if minz == maxz
+    i = if minz == maxz || z < minz
         1
+    elseif z > maxz
+        length(cmap)
     else
         1 + round(Int, ((z - minz) / (maxz - minz)) * (length(cmap) - 1))
     end

src/common.jl

@@ -186,7 +186,8 @@ function char_marker(marker::MarkerType)::Char
     if marker isa Symbol
         get(MARKERS, marker, MARKERS[:circle])
     else
-        length(marker) == 1 || throw(error("`marker` keyword has a non unit length"))
+        length(marker) == 1 ||
+            throw(ArgumentError("`marker` keyword has a non unit length"))
         marker[1]
     end
 end
@@ -208,6 +209,8 @@ function transform_name(tr, basename = "")
     string(basename, " [", name, "]")
 end
 
+as_float(x) = eltype(x) <: AbstractFloat ? x : float.(x)
+
 roundable(num::Number) = isinteger(num) & (typemin(Int) <= num < typemax(Int))
 compact_repr(num::Number) = repr(num, context = :compact => true)
 
@@ -257,25 +260,24 @@ function plotting_range(xmin, xmax)
     diff = xmax - xmin
     xmax = round_up_tick(xmax, diff)
     xmin = round_down_tick(xmin, diff)
-    Float64(xmin), Float64(xmax)
+    float(xmin), float(xmax)
 end
 
 function plotting_range_narrow(xmin, xmax)
     diff = xmax - xmin
     xmax = round_up_subtick(xmax, diff)
     xmin = round_down_subtick(xmin, diff)
-    Float64(xmin), Float64(xmax)
+    float(xmin), float(xmax)
 end
 
 extend_limits(vec, limits) = extend_limits(vec, limits, :identity)
 
 function extend_limits(vec, limits, scale::Union{Symbol,Function})
-    mi, ma = map(Float64, extrema(limits))
+    mi, ma = as_float(extrema(limits))
     if mi == 0 && ma == 0
-        mi, ma = map(Float64, extrema(vec))
+        mi, ma = as_float(extrema(vec))
     end
-    diff = ma - mi
-    if diff == 0
+    if ma - mi == 0
         ma = mi + 1
         mi = mi - 1
     end
@@ -313,6 +315,8 @@ function crayon_256_color(color::UserColorType)::ColorType
     Crayons.val(ansicolor)
 end
 
+complement(color) = (col = crayon_256_color(color)) === nothing ? nothing : ~col
+
 julia_color(color::Integer)::JuliaColorType = Int(color)
 julia_color(color::Nothing)::JuliaColorType = :normal
 julia_color(color::Symbol)::JuliaColorType = color

src/description.jl

@@ -21,6 +21,14 @@ const KEYWORDS = (
     colorbar_lim = (0, 1),
     colorbar_border = :solid,
     colormap = :viridis,
+    projection = :orthographic,
+    elevation = round(atand(1 / √2); digits = 6),
+    azimuth = 45.0,
+    axes3d = true,
+    near = 1.0,
+    far = 100.0,
+    zoom = 1.0,
+    up = :z,
     colorbar = false,
     unicode_exponent = true,
     compact = false,
@@ -35,6 +43,7 @@ const DESCRIPTION = (
     # NOTE: this named tuple has to stay ordered
     x = "horizontal position for each point",
     y = "vertical position for each point",
+    z = "depth position for each point",
     symbols = "characters used to render the bars",
     title = "text displayed on top of the plot",
     name = "current drawing annotation displayed on the right",
@@ -61,6 +70,14 @@ const DESCRIPTION = (
     grid = "draws grid-lines at the origin",
     compact = "compact plot labels",
     unicode_exponent = "use `Unicode` symbols for exponents: e.g. `10²⸱¹` instead of `10^2.1`",
+    projection = "projection for 3D plots (`:orthographic`, `:perspective`, or `Matrix-View-Projection` (MVP) matrix)",
+    axes3d = "draw 3d axes (x -> red, y -> green, z -> blue)",
+    elevation = "elevation angle above or below the `floor` plane (`-90 ≤ θ ≤ 90`)",
+    azimuth = "azimutal angle around the `up` vector (`-180° ≤ φ ≤ 180°`)",
+    zoom = "zooming factor in 3D",
+    up = "up vector (`:x`, `:y` or `:z`), prefix with `m -> -` or `p -> +` to change the sign e.g. `:mz` for `-z` axis pointing upwards",
+    near = "distance to the near clipping plane (`:perspective` projection only)",
+    far = "distance to the far clipping plane (`:perspective` projection only)",
     blend = "blend colors on the underlying canvas",
     fix_ar = "fix terminal aspect ratio (experimental)",
     visible = "visible canvas",
@@ -69,6 +86,8 @@ const DESCRIPTION = (
 const Z_DESCRIPTION =
     (:zlabel, :zlim, :colorbar, :colormap, :colorbar_lim, :colorbar_border)
 
+const PROJ_DESCRIPTION = (:projection, :azimuth, :elevation, :up, :zoom, :axes3d)
+
 const DEFAULT_KW = (
     # does not have to stay ordered
     :name,
@@ -97,6 +116,7 @@ const DEFAULT_EXCLUDED = (
     :visible,  # internals
     :fix_ar,  # experimental
     Z_DESCRIPTION...,  # by default for 2D data
+    PROJ_DESCRIPTION...,  # 3D plots
 )
 
 base_type(x) = replace(string(typeof(x).name.name), "64" => "")
@@ -129,7 +149,7 @@ function keywords(
     remove::Tuple = (),
 )
     all_kw = (; KEYWORDS..., extra...)
-    candidates = keys(extra) ∪ filter(x -> x ∈ add ∪ default, DEFAULT_KW)
+    candidates = keys(extra) ∪ filter(x -> x ∈ add ∪ default, keys(KEYWORDS))  # extra goes first !
     kw = filter(x -> x ∉ setdiff(exclude ∪ remove, add), candidates)
     @assert allunique(kw)  # extra check
     join((k isa Symbol ? "$k = $(all_kw[k] |> repr)" : k for k in kw), ", ")

src/interface/boxplot.jl

@@ -82,13 +82,13 @@ function boxplot(
 
     mean_x = (min_x + max_x) / 2
     min_x_str = compact_repr(
-        roundable(min_x) ? round(Int, Float64(min_x), RoundNearestTiesUp) : min_x,
+        roundable(min_x) ? round(Int, float(min_x), RoundNearestTiesUp) : min_x,
     )
     mean_x_str = compact_repr(
-        roundable(mean_x) ? round(Int, Float64(mean_x), RoundNearestTiesUp) : mean_x,
+        roundable(mean_x) ? round(Int, float(mean_x), RoundNearestTiesUp) : mean_x,
     )
     max_x_str = compact_repr(
-        roundable(max_x) ? round(Int, Float64(max_x), RoundNearestTiesUp) : max_x,
+        roundable(max_x) ? round(Int, float(max_x), RoundNearestTiesUp) : max_x,
     )
     label!(plot, :bl, min_x_str, color = :light_black)
     label!(plot, :b, mean_x_str, color = :light_black)
@@ -129,13 +129,13 @@ function boxplot!(
     max_x = plot.graphics.max_x
     mean_x = (min_x + max_x) / 2
     min_x_str = compact_repr(
-        roundable(min_x) ? round(Int, Float64(min_x), RoundNearestTiesUp) : min_x,
+        roundable(min_x) ? round(Int, float(min_x), RoundNearestTiesUp) : min_x,
     )
     mean_x_str = compact_repr(
-        roundable(mean_x) ? round(Int, Float64(mean_x), RoundNearestTiesUp) : mean_x,
+        roundable(mean_x) ? round(Int, float(mean_x), RoundNearestTiesUp) : mean_x,
     )
     max_x_str = compact_repr(
-        roundable(max_x) ? round(Int, Float64(max_x), RoundNearestTiesUp) : max_x,
+        roundable(max_x) ? round(Int, float(max_x), RoundNearestTiesUp) : max_x,
     )
     label!(plot, :bl, min_x_str)
     label!(plot, :b, mean_x_str)