fix stuff

Don't run gnuplot process connected to stdout

.travis.yml

@@ -6,15 +6,15 @@ os:
 #  - osx
 
 julia:
-  - 1.2
-  - nightly
+  - 1.5
+#  - nightly
 
 notifications:
   email: false
 
-matrix:
-  allow_failures:
-    - julia: nightly
+#matrix:
+#  allow_failures:
+#    - julia: nightly
 
 addons:
   apt:

ChangeLog.md

@@ -1,3 +1,68 @@
+# Version 1.4.1 (released on: )
+- New features:
+	* Implicit recipes can now returns a `Vector{PlotElement}`;
+
+	* Allow using single quotes in output file names (#52);
+
+	* New function: `palette_levels()` can be used to modify palette levels before passing them to gnuplot;
+
+- Bugfix:
+	* Fixed `BoundsErrors` in `hist()` (#49);
+
+	* Fixed problem when generating documentation (#51);
+
+
+# Version 1.4.0 (released on: May 5, 2021)
+- New features:
+    * Missing values are accepted if the input arrays have `eltype <:
+      AbstractFloat`;
+
+    * Missing values are also accepted in calls to `hist`;
+
+	* VSCode and Pluto sessions are now properly handled (#35 and #43);
+
+- Bugfix:
+	* Multiplot were not displayed in Jupyter (#25);
+
+	* `gpvars()` fails if gnuplot character encoding is utf8
+      (#24);
+
+
+# Version 1.3.0 (released on: Apr. 29, 2020)
+
+- New features:
+    * The new `dgrid3d()` allows to interpolate scattered 2D data on a
+       2D regular grid;
+
+    * The `Options` structure features a new `mime` field containing a
+      dictionary to map a MIME type to gnuplot terminals;
+
+    * The `Options` structure features a new `gpviewer` field allowing
+      to choose the display behaviour (using either gnuplot
+      interactive terminals or anexternal viewer such as Jupyter or
+      Juno);
+
+    * The `save()` function now accepts a `MIME` argument in place of
+      the `term=` keyword.  The actual terminal is retrieved from the
+      `Options.mime` dictionary;
+
+    * The `contourlines()` function now accepts `AbstractVector` and
+      `AbstractMatrix` as arguments, rather than `Vector` and
+      `Matrix`;
+
+    * The `contourlines()` function now accepts a `fractions` input to
+      generate contours encompassing given fractions of the total
+      counts in a 2D histogram;
+
+    * The `palette()` function now accept a boolean `smooth` keyword,
+      allowing to interpolate a discrete palette into a continuous one.
+
+- Breaking changes:
+    * The `Options` structure no longer provides the `term_svg` and
+      `term_png` fields.  They have been replaced by the `mime`
+      dictionary.
+
+
 # Version 1.2.0 (released on: Apr. 20, 2020)
 
 - New features:
@@ -30,6 +95,8 @@
       returned by `hist()`), contour lines (as returned by
       `contourlines()`) and images;
 
+    * Implemented automatic display of plots in both Jupyter and Juno;
+
     * Documentation updated;
 
 

Project.toml

@@ -1,6 +1,6 @@
 name = "Gnuplot"
 uuid = "dc211083-a33a-5b79-959f-2ff34033469d"
-version = "1.2.0"
+version = "1.4.1"
 
 [deps]
 ColorSchemes = "35d6a980-a343-548e-a6ea-1d62b119f2f4"
@@ -14,11 +14,11 @@ StructC14N = "d2514e9c-36c4-5b8e-97e2-51e7675c221c"
 Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
 
 [compat]
-ColorSchemes = "^3.5"
-ColorTypes = "^0.10"
+ColorSchemes = "^3.12"
+ColorTypes = "^0.11"
 Colors = "^0.12"
-DataStructures = "^0.17"
+DataStructures = "^0.18"
 ReplMaker = "^0.2"
 StatsBase = "^0.33"
 StructC14N = "^0.3"
-julia = "^1.2"
+julia = "^1.6"

README.md

@@ -3,10 +3,10 @@
 
 [![Build Status](https://travis-ci.org/gcalderone/Gnuplot.jl.svg?branch=master)](https://travis-ci.org/gcalderone/Gnuplot.jl)
 [![License](http://img.shields.io/badge/license-MIT-brightgreen.svg?style=flat)](LICENSE.md)
-[![DocumentationStatus](https://img.shields.io/badge/docs-stable-blue.svg?style=flat)](https://gcalderone.github.io/Gnuplot.jl/v1.2.0/)
+[![DocumentationStatus](https://img.shields.io/badge/docs-stable-blue.svg?style=flat)](https://gcalderone.github.io/Gnuplot.jl/v1.4.0/index.html)
 
 
-**Gnuplot.jl** is a simple package able to send both data and commands from Julia to an underlying [gnuplot](http://gnuplot.sourceforge.net/) process.  Its main purpose it to provide a fast and powerful data visualization framework, using an extremely concise Julia syntax.
+**Gnuplot.jl** is a simple package able to send both data and commands from Julia to an underlying [gnuplot](http://gnuplot.sourceforge.net/) process.  Its main purpose it to provide a fast and powerful data visualization framework, using an extremely concise Julia syntax. It also has automatic display of plots in Jupyter, Juno and VS Code.
 
 
 ## Installation
@@ -21,9 +21,9 @@ You may check the installed **Gnuplot.jl** version with:
 ```julia
 ]st Gnuplot
 ```
-If the displayed version is not v1.2.0 you are probably having a dependency conflict. In this case try forcing installation of the latest version with:
+If the displayed version is not v1.4.0 you are probably having a dependency conflict. In this case try forcing installation of the latest version with:
 ```julia
-]add Gnuplot@1.2.0
+]add Gnuplot@1.4.0
 ```
 and check which package is causing the conflict.
 
@@ -37,7 +37,7 @@ test_terminal()
 
 
 ## Quick start
-The following examples are supposed to be self-explaining.  See [documentation](https://gcalderone.github.io/Gnuplot.jl/v1.2.0/) for further informations.
+The following examples are supposed to be self-explaining.  See [documentation](https://gcalderone.github.io/Gnuplot.jl/v1.4.0/) for further informations.
 
 ### A simple parabola
 ```julia

docs/make.jl

@@ -1,4 +1,5 @@
 using Documenter, Gnuplot
+empty!(Gnuplot.options.mime)
 
 makedocs(sitename="Gnuplot.jl",
          authors = "Giorgio Calderone",

docs/src/advanced.md

@@ -255,6 +255,70 @@ saveas("advanced014e") # hide
 ![](assets/advanced014e.png)
 
 
+The [`contourlines()`](@ref) function also allows to calculate the contour lines encompassing a given fraction of the total counts of a 2D histogram.  E.g. to plot the contours corresponding to 1, 2, and 3 $\sigma$ of a 2D Gaussian distribution:
+```@example abc
+x = randn(10^5);
+y = randn(10^5);
+h = hist(x, y, nbins1=20, nbins2=20);
+
+# Calculate probability within 0 < r < σ
+p(σ) = round(1 - exp(-(σ^2) / 2), sigdigits=3)
+
+# Draw contour lines at 1, 2 and 3 σ
+clines = contourlines(h, p.(1:3));
+@gp palette(:beach, smooth=true, rev=true) "set grid front" "set size ratio -1" h clines
+saveas("advanced014f") # hide
+```
+![](assets/advanced014f.png)
+
+
+## Interpolation of 2D scattered data
+The `dgrid3d()` function allows to interpolate 2D scattered data onto a 2D regular grid, e.g.:
+```@example abc
+x = (rand(200) .- 0.5) .* 3;
+y = (rand(200) .- 0.5) .* 3;
+z = exp.(-(x.^2 .+ y.^2));
+
+# Interpolate on a 20x30 regular grid with splines
+gx, gy, gz = dgrid3d(x, y, z, "20,30 splines")
+
+@gsp "set size ratio -1" "set xyplane at 0" xlab="X" ylab="Y" :-
+@gsp :-  x  y  z "w p t 'Scattered data' lc pal"
+@gsp :- gx gy gz "w l t 'Interpolation on a grid' lc pal"
+saveas("advanced015a") # hide
+```
+![](assets/advanced015a.png)
+
+
+!!! warn
+    The `splines` algorithm may be very slow on large datasets.  An alternative option is to use a smoothing kernel, such as `gauss`.
+
+
+The interpolated data in scarcely sampled regions are poorly constrained, i.e. they are actually *extrapolated values*.  By using the `extra=false` keyword all extrapolated values are set to `NaN`:
+
+```@example abc
+x = randn(2000) .* 0.5;
+y = randn(2000) .* 0.5;
+rsq = x.^2 + y.^2;
+z = exp.(-rsq) .* sin.(y) .* cos.(2 * rsq);
+
+@gsp "set size ratio -1" palette(:balance, smooth=true) "set view map" "set pm3d" :-
+@gsp :- "set multiplot layout 1,3" xr=[-2,2] yr=[-2,2] :-
+@gsp :- 1 tit="Scattered data"  x  y  z "w p notit lc pal"
+
+# Show extrapolated values
+gx, gy, gz = dgrid3d(x, y, z, "40,40 gauss 0.1,0.1")
+@gsp :- 2 tit="Interpolation on a grid\\n(extrapolated values are shown)"  gx gy gz "w l notit lc pal"
+
+# Hide exrapolated values
+gx, gy, gz = dgrid3d(x, y, z, "40,40 gauss 0.1,0.1", extra=false)
+@gsp :- 3 tit="Interpolation on a grid\\n(extrapolated values are hidden)" gx gy gz "w l notit lc pal"
+save(term="pngcairo size 1000,400 fontscale 1.0", output="assets/advanced015b.png")  # hide
+```
+![](assets/advanced015b.png)
+
+
+
 ## Animations
 
 The [Multiplot](@ref) capabilities can also be used to stack plots one above the other in order to create an animation, as in the following example:

docs/src/api.md

@@ -13,6 +13,7 @@ The list of **Gnuplot.jl** exported symbols is as follows:
 boxxy
 contourlines
 dataset_names
+dgrid3d
 gpexec
 gpmargins
 gpranges
@@ -20,6 +21,7 @@ gpvars
 hist
 linetypes
 palette
+palette_levels
 palette_names
 recipe
 save

docs/src/basic.md

@@ -21,7 +21,7 @@ saveas("basic000") # hide
 ```
 ![](assets/basic000.png)
 
-The plots are displayed either in an interactive window (if running in the Julia REPL), as an inline image (if running in Jupyter) or in the plot pane (if running in Juno).  See [Options](@ref) and [Jupyter and Juno](@ref) for further informations.
+The plots are displayed either in an interactive window (if running in the Julia REPL), as an inline image (if running in Jupyter) or in the plot pane (if running in Juno).  See [Display options](@ref) for further informations.
 
 
 Both the [`@gp`](@ref) and [`@gsp`](@ref) macros accept any number of arguments, whose meaning is interpreted as follows:
@@ -210,6 +210,16 @@ saveas("basic008a") # hide
 ```
 ![](assets/basic008a.png)
 
+The palette levels may be easily stretched by using the [`palette_levels()`](@ref) and modifying the numeric levels, e.g.:
+```@example abc
+x = 0:0.1:10pi
+v, l, n = palette_levels(:viridis)
+@gsp palette(v.^0.25, l, n) cbr=[-1,1].*30 :-
+@gsp :-  x  x.*sin.(x)  x.*cos.(x)  x./20  "w p pt 7 ps var lc pal"
+saveas("basic008b") # hide
+```
+![](assets/basic008b.png)
+
 The list of all available palette can be retrieved with [`palette_names()`](@ref):
 ```@repl abc
 palette_names()

docs/src/index.md

@@ -21,6 +21,8 @@ The **Gnuplot.jl** package allows easy and fast use of [gnuplot](http://gnuplot.
 
 - enhanced support for contour plots;
 
+- 2D interpolation of scattered data on a regular grid;
+
 - export to a huge number of formats such as `pdf`, `png`, `gif`, ``\LaTeX``, `svg`, etc. (actually all those supported by gnuplot);
 
 - compatibility with Jupyter and Juno;

docs/src/install.md

@@ -19,11 +19,11 @@ Check **Gnuplot.jl** version with:
 ```julia-repl
 julia> ]st Gnuplot
 Status `~/.julia/environments/v1.4/Project.toml`
-  [dc211083] Gnuplot v1.2.0
+  [dc211083] Gnuplot v1.4.1
 ```
-If the displayed version is not `v1.2.0` you are probably having a dependency conflict.  In this case try forcing installation of the latest version with:
+If the displayed version is not `v1.4.1` you are probably having a dependency conflict.  In this case try forcing installation of the latest version with:
 ```julia-repl
-julia> ]add Gnuplot@1.2.0
+julia> ]add Gnuplot@1.4.1
 ```
 and check which package is causing the conflict.
 

docs/src/options.md

@@ -9,6 +9,19 @@ push!( Gnuplot.options.init, linetypes(:Set1_5, lw=1.5, ps=1.5))
 saveas(file) = save(term="pngcairo size 550,350 fontscale 0.8", output="assets/$(file).png")
 ```
 
+# Display options
+
+The display behaviour of **Gnuplot.jl** depends on the value of the `Gnuplot.options.gpviewer` flag:
+
+- if `true` the plot is displayed in a gnuplot window, using one of the interactive terminals such as `wxt`, `qt` or `aqua`.  This is the default setting when running a Julia REPL session; The terminal options can be customized using `Gnuplot.options.term`;
+
+- if `false` the plot is displayed through the Julia [multimedia interface](https://docs.julialang.org/en/v1/base/io-network/#Multimedia-I/O-1), i.e. it is exported as either a `png`, `svg` or `html` file, and displayed in an external viewer.  This is the default setting when running a Jupyter, JupyterLab or Juno session.  The terminal options can be customized using the `Gnuplot.options.mime` dictionary.
+
+The `Gnuplot.options.gpviewer` flag is automatically set when the package is first loaded according to the runtime environment, however the user can change its value at any time to fit specific needs.  Further informations and examples for both options are available in this Jupyter [notebook](https://github.com/gcalderone/Gnuplot.jl/blob/gh-pages/v1.3.0/options/display.ipynb).
+
+
+
+
 # Package options and initialization
 
 ## Options
@@ -18,14 +31,21 @@ The package options are stored in a global structure available in Julia as `Gnup
 
 - `cmd::String`: command to start the gnuplot process, default value is `"gnuplot"`.  Use this field to specify a custom path to the gnuplot executable;
 
+- `gpviewer::Bool`: use a gnuplot terminal as main plotting device (if `true`) or an external viewer (if `false`);
+
 - `term::String`: default terminal for interactive use (default is an empty string, i.e. use gnuplot settings).  A custom terminal can be set with, e.g.:
 ```@repl abc
 Gnuplot.options.term = "wxt size 700,400";
 ```
 
-- `term_svg::String`: terminal to save png files (default: `"svg background rgb 'white' dynamic"`);
+- `mime::Dict{MIME, String}`: dictionary of MIME types and corresponding gnuplot terminals.  Used to export images with either [`save()`](@ref) or `show()` (see [Display options](@ref)).  Default values are:
+  - `MIME"application/pdf" => "pdfcairo enhanced"`
+  - `MIME"image/jpeg"      => "jpeg enhanced"`
+  - `MIME"image/png"       => "pngcairo enhanced"`
+  - `MIME"image/svg+xml"   => "svg enhanced mouse standalone dynamic background rgb 'white'"`
+  - `MIME"text/html"       => "svg enhanced mouse standalone dynamic"`
+  - `MIME"text/plain"      => "dumb enhanced ansi"`
 
-- `term_png::String`: terminal to save png files (default: `"pngcairo"`);
 
 - `init::Vector{String}`: commands to initialize the session when it is created or reset.  It can be used to, e.g., set a custom linetypes or palette:
 ```@repl abc
@@ -48,46 +68,47 @@ gpexec("set term unknown");                                # hide
 Each line reports the package name (`GNUPLOT`), the session name (`default`), the command or string being sent to gnuplot process, and the returned response (line starting with `->`).  Default value is `false`;
 
 
-## Jupyter and Juno
-
-**Gnuplot.jl** can display plots in Jupyter and Juno by exporting images in the PNG and SVG formats.  To customize the terminals used to export the images set the `term_png` or `term_svg` fields of the [`Gnuplot.Options`](@ref) structure, e.g.:
-```@repl abc
-Gnuplot.options.term_png = "pngcairo size 700,400 linewidth 2";
-Gnuplot.options.term_svg = "svg dynamic";
-```
-
-
 ## Package initialization
 
-If you use **Gnuplot.jl** frequently you may find convenient to collect all the package settings ([Options](@ref)) in a single place, to quickly recall them in a Julia session.  I suggest to put the following code in the `~/.julia/config/startup.jl` initialization file (further info [here](https://docs.julialang.org/en/v1/stdlib/REPL/)):
+If you use **Gnuplot.jl** frequently you may find convenient to automatically apply the package settings ([Options](@ref)) whenever the package is loaded.  A possibility is to use [Requires.jl](https://github.com/JuliaPackaging/Requires.jl) and put the following code in the `~/.julia/config/startup.jl` initialization file (further info [here](https://docs.julialang.org/en/v1/stdlib/REPL/)):
 ```julia
-macro gnuplotrc()
-    return :(
-        using Revise, Gnuplot;
-
-        # Uncomment following to true if you don't have the gnuplot
+using Requires
+@require Gnuplot="dc211083-a33a-5b79-959f-2ff34033469d" begin
+        @info "Custom Gnuplot initialization"
+        # Uncomment the following if you don't have the gnuplot
         # executable installed on your platform:
         #Gnuplot.options.dry = true;
 
-        # Uncomment the following and set the proper path if the
-        # gnuplot executable is not available in your $PATH
+        # Set the proper path if the gnuplot executable is not
+        # available in your $PATH
         #Gnuplot.options.cmd = "/path/to/gnuplot";
 
+        # Force a specific display behaviour (see documentation).  If
+        # not given explicit Gnuplot.jl will choose the best option
+        # according to your runtime environment.
+        #Gnuplot.options.gpviewer = true
+
         # Set the default terminal for interacitve use
         Gnuplot.options.term = "wxt size 700,400";
 
+        # Set the terminal options for the exported MIME types:
+        #Gnuplot.options.mime[MIME"image/png"] = "";
+        #Gnuplot.options.mime[MIME"image/svg+xml"] = "svg enhanced standalone dynamic";
+        #Gnuplot.options.mime[MIME"text/html"] = "svg enhanced standalone mouse dynamic";
+
+        # Set the terminal to plot in a terminal emulator:
+        # (try with `save(MIME"text/plain")`):
+        #Gnuplot.options.mime[MIME"text/plain"] = "sixelgd enhanced"; # requires vt340 emulation
+
         # Set the default linetypes
         empty!(Gnuplot.options.init);
-        push!(Gnuplot.options.init, linetypes(:Set1_5, lw=1.5, ps=1.5));
+        push!(Gnuplot.options.init, Gnuplot.linetypes(:Set1_5, lw=1.5, ps=1.5));
 
         # Initialize the gnuplot REPL using the provided `start_key`.
-        # Comment the following to disable the REPL.
+        if Gnuplot.options.gpviewer;
             Gnuplot.repl_init(start_key='>');
-    )
+        end;
 end
 ```
-At the Julia prompt you may load the package and the associated settings by typing:
-```julia
-julia> @gnuplotrc
-```
-and you're ready to go.
+
+The above code will be automatically when you first load the package with `using Gnuplot`.

docs/src/recipes.md

@@ -125,7 +125,7 @@ end
 with only one mandatory argument.  In order to exploit the optional keyword we can explicitly invoke the recipe as follows:
 ```@example abc
 img = testimage("walkbridge");
-@gp palette(:gray) recipe(img, "flipy rot=15deg")
+@gp palette(:gray1) recipe(img, "flipy rot=15deg")
 saveas("recipes007c") # hide
 ```
 ![](assets/recipes007c.png)

docs/src/style.md

@@ -1,6 +1,6 @@
 # Style Guide
 
-The **Gnuplot.jl** loose syntax allows to create a plot using very different approaches.  While this was one of the initial purposes for the package, it may lead to decreased code readability if not used judiciously.
+The **Gnuplot.jl** loose syntax allows to create a plot using very different approaches.  While this was one of the initial purposes for the package, it may lead to a reduced code readability if not used judiciously.
 
 Here I will summarize a few, non-mandatory, guidelines which allows to maintain a neat syntax and a high readability:
 

docs/src/terminals.md

@@ -4,11 +4,7 @@ Gnuplot provides dozens of terminals to display plots or export them into files
 
 To use a specific terminal for interactive use you may either set it as initialization command for all new session with (see [Options](@ref)):
 ```julia
-Gnuplot.options.term = "wxt")
-```
-or directly send the command to a specific session (see [Direct command execution](@ref))
-```julia
-gpexec("set term wxt")
+Gnuplot.options.term = "wxt"
 ```
 See official [gnuplot documentation](http://gnuplot.sourceforge.net/documentation.html) for further info on terminals and their options.
 
@@ -16,13 +12,13 @@ See official [gnuplot documentation](http://gnuplot.sourceforge.net/documentatio
 ## Interactive terminals (`wxt` and `qt`)
 The multiplatform `wxt` and `qt` terminals are among the most widely used ones for their nicely looking outputs on display and for their interactive capabilities.
 
-You may set them as terminal with:
+You may use such terminals with:
 ```
-"set term wxt size 800,600"
+Gnuplot.options.term = "wxt size 800,600"
 ```
 or
 ```
-"set term qt  size 800,600"
+Gnuplot.options.term = "qt  size 800,600"
 ```
 (the `size 800,600` is optional and can be omitted).
 
@@ -30,9 +26,35 @@ Press the `h` key on the window to display an help message with all available ke
 
 
 ## Plot in a terminal application (`dumb`, `sixel` and `sixelgd`)
-Gnuplot supports plotting in a terminal application, with no need for X11 or other GUI support, via the `dumb`, `sixel` and `sixelgd` terminals.  These are extremely useful when you run Julia on a remote shell through `ssh`, with no X11 forwarding.
+Gnuplot supports plotting in a terminal application, with no need for X11 or other GUI support, via the `dumb`, `sixel` and `sixelgd` terminals.  These are extremely useful when you run Julia on a remote shell through `ssh`, with no X11 forwarding. The `dumb` terminal uses ASCII characters to draw a plot, while `sixel` and `sixelgd` actually use bitmaps (but require Sixel support to be enabled in the terminal, e.g. `xterm -ti vt340`).
 
-The `dumb` terminal uses ASCII characters to draw a plot, while `sixel` and `sixelgd` actually use bitmaps (but require Sixel support to be enabled in the terminal, e.g. `xterm -ti vt340`).  A sixel plot on `xterm` looks as follows:
+Dumb terminal can be used as follows:
+
+```jldoctest; setup = :(using Gnuplot)
+julia> origterm = Gnuplot.options.term;
+
+julia> Gnuplot.options.term = "dumb size 60,15";
+
+julia> @gp "plot sin(x)"
+
+    1 +-------------------------------------------------+
+  0.8 |-+         *+ *         + ** **      +    *  * +-|
+  0.6 |-+        *   **          *   *   sin(x) *******-|
+  0.4 |*+        *    *         *     *         *    *+-|
+  0.2 |*+       *      *        *     *        *      *-|
+    0 |*+       *      *       *       *       *      *-|
+      | *      *       *       *       *      *        *|
+ -0.2 |-*      *        *     *        *      *       +*|
+ -0.4 |-+*    *         *     *         *    *        +*|
+ -0.6 |-+*    *          *   *          **   *        +-|
+ -0.8 |-+ *  *     +     ** ** +         *  *         +-|
+   -1 +-------------------------------------------------+
+     -10          -5           0            5           10
+
+julia> Gnuplot.options.term = origterm;
+
+```
+A sixel plot on `xterm` looks as follows:
 ![](assets/sixelgd.png)
 
 The above terminals are available if gnuplot has been compiled with the `--with-bitmap-terminals` option enabled and Libgd (only for `sixelgd`).
@@ -41,7 +63,7 @@ The above terminals are available if gnuplot has been compiled with the `--with-
 ## Export to image files
 
 Gnuplot provides dozens of terminals able to export on files.  Examples are:
-- `cairopng` to export PNG files;
+- `pngcairo` to export PNG files;
 - `pdfcairo` for PDF;
 - `jpeg` for JPG;
 - `gif` for GIF (see [Animations](@ref)).
@@ -58,7 +80,7 @@ Gnuplot is also able to export vector (i.e. non-raster) plots through the `svg`
 The `cairolatex` terminal allows to produce high quality plots by splitting the output into a PDF file (containing a rasterized image of a plot) and a `.tex` file (containing all the text as ``\LaTeX`` code).  The following example shows how to write plot tics and an equation in ``\LaTeX``:
 ```julia
 x = LinRange(-2pi, 2pi, 1000)
-@gp t="Polynomial approximation of sin(x)"  "set style fill transparent solid 0.6 noborder"
+@gp tit="Polynomial approximation of sin(x)"  "set style fill transparent solid 0.6 noborder"
 @gp :- raw"""set xtics ('$-\pi$' -pi, '$-\pi/2$' -pi/2, 0, '$\pi/2$' pi/2, '$\pi$' pi)"""
 @gp :- xr=3.8.*[-1, 1] yr=[-1.5,1.5] key="box opaque left horiz" linetypes(:Blues_3) "set grid front"
 latex = raw"""\begin{minipage}[c]{\textwidth}\begin{equation*}""" *
@@ -74,7 +96,7 @@ approx = fill(0., length(x));
 save(term="cairolatex pdf input color dashed size 5in,3.3in", output="test.tex")
 ```
 !!! warning
-    If you add a path in the `output=` keyword this will also be copied in the the `.tex` file.  I suggest to use just filenames, with no path, in order to avoid possible errors when compiling ``\LaTeX`` code.
+    If you add a path in the `output=` keyword this will also be copied in the the `.tex` file, and may generate errors when compiling ``\LaTeX`` code.  The simplest way to solve this problem is to use just filenames, with no paths.
 
 The two output files (`test.tex` and `test.pdf`) can then be included in a ``\LaTeX`` file as follows:
 ```latex

src/Gnuplot.jl

@@ -7,10 +7,10 @@ import Base.reset
 import Base.write
 import Base.show
 
-export session_names, dataset_names, palette_names, linetypes, palette,
+export session_names, dataset_names, palette_names, linetypes, palette_levels, palette,
     terminal, terminals, test_terminal,
     stats, @gp, @gsp, save, gpexec,
-    boxxy, contourlines, hist, recipe, gpvars, gpmargins, gpranges
+    boxxy, contourlines, dgrid3d, hist, recipe, gpvars, gpmargins, gpranges
 
 
 # ╭───────────────────────────────────────────────────────────────────╮
@@ -184,7 +184,6 @@ mutable struct GPSession <: Session
     plots::Vector{SinglePlot}           # commands and plot commands (one entry for each plot of the multiplot)
     curmid::Int                         # current multiplot ID
     pin::Base.Pipe;
-    pout::Base.Pipe;
     perr::Base.Pipe;
     proc::Base.Process;
     channel::Channel{String};
@@ -202,8 +201,8 @@ Structure containing the package global options, accessible through `Gnuplot.opt
 - `cmd::String`: command to start the Gnuplot process (default: `"gnuplot"`)
 - `default::Symbol`: default session name (default: `:default`)
 - `term::String`: default terminal for interactive use (default: empty string, i.e. use gnuplot settings);
-- `term_svg::String`: terminal to save png files (default `"svg background rgb 'white' dynamic"`);
-- `term_png::String`: terminal to save png files (default `"pngcairo"`);
+- `mime::Dict{DataType, String}`: dictionary of MIME types and corresponding gnuplot terminals.  Used to export images with either [`save()`](@ref) or `show()` (see [Display options](@ref));
+- `gpviewer::Bool`: use a gnuplot terminal as main plotting device (if `true`) or an external viewer (if `false`);
 - `init::Vector{String}`: commands to initialize the session when it is created or reset (e.g., to set default palette);
 - `verbose::Bool`: verbosity flag (default: `false`)
 - `preferred_format::Symbol`: preferred format to send data to gnuplot.  Value must be one of:
@@ -216,8 +215,14 @@ Base.@kwdef mutable struct Options
     cmd::String = "gnuplot"
     default::Symbol = :default
     term::String = ""
-    term_svg::String = "svg background rgb 'white' dynamic"
-    term_png::String = "pngcairo"
+    mime::Dict{DataType, String} = Dict(
+        MIME"application/pdf" => "pdfcairo enhanced",
+        MIME"image/jpeg"      => "jpeg enhanced",
+        MIME"image/png"       => "pngcairo enhanced",
+        MIME"image/svg+xml"   => "svg enhanced mouse standalone dynamic background rgb 'white'",
+        MIME"text/html"       => "svg enhanced mouse standalone dynamic",  # canvas mousing
+        MIME"text/plain"      => "dumb enhanced ansi")
+    gpviewer::Bool = false
     init::Vector{String} = Vector{String}()
     verbose::Bool = false
     preferred_format::Symbol = :auto
@@ -225,11 +230,26 @@ end
 
 
 # ╭───────────────────────────────────────────────────────────────────╮
-# │                         GLOBAL VARIABLES                          │
+# │           GLOBAL VARIABLES AND MODULE INITIALIZATION              │
 # ╰───────────────────────────────────────────────────────────────────╯
 const sessions = OrderedDict{Symbol, Session}()
 const options = Options()
 
+function __init__()
+    # Check whether we are running in an IJulia, Juno, VSCode or Pluto session.
+    # (copied from Gaston.jl).
+    options.gpviewer = !(
+        ((isdefined(Main, :IJulia)  &&  Main.IJulia.inited)  ||
+         (isdefined(Main, :Juno)    &&  Main.Juno.isactive()) ||
+         (isdefined(Main, :VSCodeServer)) ||
+         (isdefined(Main, :PlutoRunner)) )
+    )
+    if isdefined(Main, :VSCodeServer)
+        # VS Code shows "dynamic" plots with fixed and small size :-(
+        options.mime[MIME"image/svg+xml"] = replace(options.mime[MIME"image/svg+xml"], "dynamic" => "")
+    end
+end
+
 
 # ╭───────────────────────────────────────────────────────────────────╮
 # │                         LOW LEVEL FUNCTIONS                       │
@@ -263,11 +283,11 @@ function parseKeywords(; kwargs...)
     ismissing(kw.zrange ) || (push!(out, replace("set zrange  [" * join(kw.zrange , ":") * "]", "NaN"=>"*")))
     ismissing(kw.cbrange) || (push!(out, replace("set cbrange [" * join(kw.cbrange, ":") * "]", "NaN"=>"*")))
     ismissing(kw.key    ) || (push!(out, "set key " * kw.key  * ""))
-    ismissing(kw.title  ) || (push!(out, "set title   \"" * kw.title  * "\""))
-    ismissing(kw.xlabel ) || (push!(out, "set xlabel  \"" * kw.xlabel * "\""))
-    ismissing(kw.ylabel ) || (push!(out, "set ylabel  \"" * kw.ylabel * "\""))
-    ismissing(kw.zlabel ) || (push!(out, "set zlabel  \"" * kw.zlabel * "\""))
-    ismissing(kw.cblabel) || (push!(out, "set cblabel \"" * kw.cblabel * "\""))
+    ismissing(kw.title  ) || (push!(out, "set title   '" * kw.title  * "'"))
+    ismissing(kw.xlabel ) || (push!(out, "set xlabel  '" * kw.xlabel * "'"))
+    ismissing(kw.ylabel ) || (push!(out, "set ylabel  '" * kw.ylabel * "'"))
+    ismissing(kw.zlabel ) || (push!(out, "set zlabel  '" * kw.zlabel * "'"))
+    ismissing(kw.cblabel) || (push!(out, "set cblabel '" * kw.cblabel * "'"))
     ismissing(kw.xlog   ) || (push!(out, (kw.xlog  ?  ""  :  "un") * "set logscale x"))
     ismissing(kw.ylog   ) || (push!(out, (kw.ylog  ?  ""  :  "un") * "set logscale y"))
     ismissing(kw.zlog   ) || (push!(out, (kw.zlog  ?  ""  :  "un") * "set logscale z"))
@@ -284,7 +304,6 @@ function parseKeywords(; kwargs...)
     ismissing(kw.rmargin) ||  push!(out, (kw.rmargin == ""  ?  "unset rmargin"  :  "set rmargin at screen $(kw.rmargin)"))
     ismissing(kw.bmargin) ||  push!(out, (kw.bmargin == ""  ?  "unset bmargin"  :  "set bmargin at screen $(kw.bmargin)"))
     ismissing(kw.tmargin) ||  push!(out, (kw.tmargin == ""  ?  "unset tmargin"  :  "set tmargin at screen $(kw.tmargin)"))
-
     return join(out, ";\n")
 end
 
@@ -439,60 +458,68 @@ end
 
 
 # ---------------------------------------------------------------------
+
+function readTask(gp::GPSession)
     pagerTokens() = ["Press return for more:"]
 
-function GPSession(sid::Symbol)
-    function readTask(sid, stream, channel)
-        function gpreadline(stream)
+    captureID = 0
+    function gpreadline()
         line = ""
         while true
-                c = read(stream, Char)
+            c = read(gp.perr, Char)
             (c == '\r')  &&  continue
             (c == '\n')  &&  break
-                if c == Char(0x1b)  # sixel
-                    buf = Vector{UInt8}()
-                    push!(buf, UInt8(c))
-                    while true
-                        c = read(stream, Char)
-                        push!(buf, UInt8(c))
-                        (c == Char(0x1b))  &&  break
-                    end
-                    c = read(stream, Char)
-                    push!(buf, UInt8(c))
-                    write(stdout, buf)
-                    continue
-                end
             line *= c
             for token in pagerTokens()  # handle pager interaction
                 if (length(line) == length(token))  &&  (line == token)
-                        return line
+                    # GNUPLOT_CAPTURE_END maybe lost when pager is
+                    # running: send it again.
+                    captureID += 1
+                    write(gp.pin, "\nprint 'GNUPLOT_CAPTURE_END $(captureID)'\n")
+                    line = ""
                 end
             end
         end
         return line
     end
 
+    try
         saveOutput = false
-        while isopen(stream)
-            line = gpreadline(stream)
+        while isopen(gp.perr)
+            line = gpreadline()
+
             if line == "GNUPLOT_CAPTURE_BEGIN"
                 saveOutput = true
-            elseif line == "GNUPLOT_CAPTURE_END"
-                put!(channel, line)
+            elseif line == "GNUPLOT_CAPTURE_END $(captureID)"
                 saveOutput = false
+                put!(gp.channel, "GNUPLOT_CAPTURE_END")
+                captureID = 0
+            elseif !isnothing(findfirst("GNUPLOT_CAPTURE_END", line))
+                ;# old GNUPLOT_CAPTURE_END, ignore it
             else
                 if line != ""
                     if options.verbose  ||  !saveOutput
-                        printstyled(color=:cyan, "GNUPLOT ($sid) -> $line\n")
+                        printstyled(color=:cyan, "GNUPLOT ($(gp.sid)) -> $line\n")
                     end
                 end
-                (saveOutput)  &&  (put!(channel, line))
+                (saveOutput)  &&  (put!(gp.channel, line))
             end
         end
-        delete!(sessions, sid)
-        return nothing
+    catch err
+        if isopen(gp.perr)
+            @error "Error occurred in readTask for session $(gp.sid)"
+            @show(err)
+        else
+            put!(gp.channel, "GNUPLOT_CAPTURE_END")
+        end
+    end
+    if options.verbose
+        printstyled(color=:red, "GNUPLOT ($(gp.sid)) Process terminated\n")
+    end
+    delete!(sessions, gp.sid)
 end
 
+function GPSession(sid::Symbol)
     session = DrySession(sid)
     if !options.dry
         try
@@ -518,14 +545,25 @@ function GPSession(sid::Symbol)
     Base.start_reading(pout.out)
     Base.start_reading(perr.out)
 
-    # Start reading tasks
-    @async readTask(sid, pout, chan)
-    @async readTask(sid, perr, chan)
-
     out = GPSession(getfield.(Ref(session), fieldnames(DrySession))...,
-                    pin, pout, perr, proc, chan)
+                    pin, perr, proc, chan)
     sessions[sid] = out
 
+    # Start reading tasks
+    @async readTask(out)
+    @async while !eof(pout) # see PR #51
+        write(stdout, readavailable(pout))
+    end
+
+    # Read gnuplot default terminal
+    if options.term == ""
+        options.term = terminal()
+    end
+
+    # The stderr of the gnuplot process goes to Julia which can parse
+    # UTF8 characters (regardless of the terminal).
+    gpexec("set encoding utf8")
+
     return out
 end
 
@@ -611,33 +649,15 @@ function writeread(gp::GPSession, str::AbstractString)
     write(gp, str)
 
     options.verbose = false
-    write(gp, "print 'GNUPLOT_CAPTURE_END'")
+    write(gp, "print 'GNUPLOT_CAPTURE_END 0'")
     options.verbose = verbose
 
     out = Vector{String}()
     while true
         l = take!(gp.channel)
-        if l in pagerTokens()
-            # Consume all data from the pager
-            while true
-                write(gp, "")
-                sleep(0.5)
-                if isready(gp.channel)
-                    while isready(gp.channel)
-                        push!(out, take!(gp.channel))
-                    end
-                else
-                    options.verbose = false
-                    write(gp, "print 'GNUPLOT_CAPTURE_END'")
-                    options.verbose = verbose
-                    break
-                end
-            end
-        else
         l == "GNUPLOT_CAPTURE_END"  &&  break
         push!(out, l)
     end
-    end
     return out
 end
 
@@ -693,7 +713,7 @@ end
 function DatasetBin(cols::Vararg{AbstractVector, N}) where N
     source = "binary record=$(length(cols[1])) format='"
     types = Vector{DataType}()
-    (length(cols) == 1)  &&  (source *= "%int")
+    #(length(cols) == 1)  &&  (source *= "%int")
     for i in 1:length(cols)
         @assert length(cols[1]) == length(cols[i])
         if     isa(cols[i][1], Int32);   push!(types, Int32);   source *= "%int"
@@ -710,7 +730,7 @@ function DatasetBin(cols::Vararg{AbstractVector, N}) where N
     (path, io) = mktemp()
     source = " '$path' $source"
     for row in 1:length(cols[1])
-        (length(cols) == 1)  &&  (write(io, convert(Int32, row)))
+        #(length(cols) == 1)  &&  (write(io, convert(Int32, row)))
         for col in 1:length(cols)
             write(io, convert(types[col], cols[col][row]))
         end
@@ -741,14 +761,6 @@ end
 # │              PRIVATE FUNCTIONS TO MANIPULATE SESSIONS             │
 # ╰───────────────────────────────────────────────────────────────────╯
 # ---------------------------------------------------------------------
-function enableExportThroughShow()
-    # Trick to check whether we are running in a IJulia or Juno
-    # session.  Copied from Gaston.jl.
-    return ((isdefined(Main, :IJulia)  &&  Main.IJulia.inited)  ||
-            (isdefined(Main, :Juno)    &&  Main.Juno.isactive()))
-end
-
-
 function reset(gp::Session)
     delete_binaries(gp)
     gp.datas = OrderedDict{String, Dataset}()
@@ -758,11 +770,8 @@ function reset(gp::Session)
     gpexec(gp, "set output")
     gpexec(gp, "reset session")
 
-    # When the `show()` method is enabled ignore options.term and set
-    # the unknown terminal
-    if enableExportThroughShow()
-        gpexec(gp, "set term unknown")
-    else
+    if options.gpviewer
+        # Use gnuplot viewer
         (options.term != "")  &&  gpexec(gp, "set term " * options.term)
 
         # Set window title (if not already set)
@@ -773,11 +782,14 @@ function reset(gp::Session)
                 writeread(gp, "set term $term $opts title 'Gnuplot.jl: $(gp.sid)'")
             end
         end
+    else
+        # Use external viewer
+        gpexec(gp, "set term unknown")
     end
 
     # Note: the reason to keep Options.term and .init separate are:
-    # - .term can be overriden by enableExportThroughShow()
-    # - .init is dumped in scripts, while .term is not
+    # - .term can be overriden by "unknown" (if options.gpviewer is false);
+    # - .init is dumped in scripts, while .term is not;
     add_cmd.(Ref(gp), options.init)
     return nothing
 end
@@ -806,7 +818,7 @@ function useBinaryMethod(args...)
     elseif options.preferred_format == :auto
         if (length(args) == 1)  &&  isa(args[1], AbstractMatrix)
             binary = true
-        elseif all(ndims.(args) .== 1)
+        elseif all(ndims.(args) .== 1)  &&  all(Base.:<:.(eltype.(args), Real))
             s = sum(length.(args))
             if s > 1e4
                 binary = true
@@ -863,7 +875,6 @@ end
 
 function quit(gp::GPSession)
     close(gp.pin)
-    close(gp.pout)
     close(gp.perr)
     wait( gp.proc)
     exitCode = gp.proc.exitcode
@@ -915,9 +926,12 @@ function execall(gp::GPSession; term::AbstractString="", output::AbstractString=
     if term != ""
         former_term = writeread(gp, "print GPVAL_TERM")[1]
         former_opts = writeread(gp, "print GPVAL_TERMOPTIONS")[1]
+        gpexec(gp, "unset multiplot")
         gpexec(gp, "set term $term")
     end
-    (output != "")  &&  gpexec(gp, "set output '$output'")
+    (output != "")  &&  gpexec(gp, "set output '$(replace(output, "'" => "''"))'")
+
+    # printstyled("Plotting with terminal: " * terminal() * "\n", color=:blue, bold=true)
 
     for i in 1:length(gp.plots)
         d = gp.plots[i]
@@ -930,7 +944,7 @@ function execall(gp::GPSession; term::AbstractString="", output::AbstractString=
             gpexec(gp, s)
         end
     end
-    (length(gp.plots) > 1)  &&  gpexec(gp, "unset multiplot")
+    gpexec(gp, "unset multiplot")
     (output != "")  &&  gpexec(gp, "set output")
     if term != ""
         gpexec(gp, "set term $former_term $former_opts")
@@ -1039,35 +1053,49 @@ function parseArguments(_args...)
         return (isplot, is3d, string(cmd))
     end
 
-    # First pass: check for `:-` and session names
-    sid = options.default
-    doDump  = true
-    doReset = true
-    if length(_args) == 0
-        return (sid, doReset, doDump, Vector{PlotElement}())
+    # First pass: check for session names and `:-`
+    sid = nothing
+    args = Vector{Any}([_args...])
+    pos = 1
+    while pos <= length(args)
+        arg = args[pos]
+        if  (typeof(arg) == Symbol)  &&
+            (arg != :-)
+            @assert isnothing(sid) "Only one session at a time can be addressed"
+            sid = arg
+            deleteat!(args, pos)
+            continue
         end
-    for iarg in 1:length(_args)
-        arg = _args[iarg]
+        pos += 1
+    end
+    isnothing(sid)  &&  (sid = options.default)
 
-        if typeof(arg) == Symbol
-            if arg == :-
-                if iarg == 1
+    if length(args) == 0
         doReset = false
-                elseif iarg == length(_args)
+    else
+        doReset = true
+    end
+    doDump = true
+    pos = 1
+    while pos <= length(args)
+        arg = args[pos]
+        if typeof(arg) == Symbol
+            @assert arg == :-
+            if pos == 1
+                doReset = false
+            elseif pos == length(args)
                 doDump  = false
             else
-                    @warn "Symbol `:-` at position $iarg in argument list has no meaning."
-                end
-            else
-                @assert (sid == options.default) "Only one session at a time can be addressed"
-                sid = arg
+                @warn "Symbol `:-` has a meaning only if it is at first or last position in argument list."
             end
+            deleteat!(args, pos)
+            continue
         end
+        pos += 1
     end
 
     # Second pass: check data types, run implicit recipes and splat
     # Vector{PlotElement}
-    args = Vector{Any}([_args...])
     pos = 1
     while pos <= length(args)
         arg = args[pos]
@@ -1093,15 +1121,24 @@ function parseArguments(_args...)
             end
             insert!(args, pos, string(strip(arg[1])) => nothing)
         elseif isa(arg, AbstractArray) &&            # ==> a dataset column
-            ((valtype(arg) <: Real)    ||
-             (valtype(arg) <: AbstractString))  ;
+            ((nonmissingtype(eltype(arg)) <: Real)    ||
+             (nonmissingtype(eltype(arg)) <: AbstractString))  ;
         elseif isa(arg, Real)                        # ==> a dataset column with only one row
             args[pos] = [arg]
         elseif isa(arg, Dataset)                ;    # ==> a Dataset object
         elseif hasmethod(recipe, tuple(typeof(arg))) # ==> implicit recipe
             # @info which(recipe, tuple(typeof(arg)))  # debug
             deleteat!(args, pos)
-            insert!(args, pos, recipe(arg))
+            pe = recipe(arg)
+            if isa(pe, PlotElement)
+                insert!(args, pos, pe)
+            elseif isa(pe, Vector{PlotElement})
+                for i in 1:length(pe)
+                    insert!(args, pos, pe[i])
+                end
+            else
+                error("Recipe must return a PlotElement or Vector{PlotElement}")
+            end
             continue
         elseif isa(arg, Vector{PlotElement})         # ==> explicit recipe (vector)
             deleteat!(args, pos)
@@ -1116,40 +1153,54 @@ function parseArguments(_args...)
         pos += 1
     end
 
-    # Third pass: convert data into Dataset objetcs
+    # Third pass: convert data into Dataset objects
     pos = 1
+    accum = Vector{AbstractArray}()
     while pos <= length(args)
         arg = args[pos]
-        if isa(arg, AbstractArray)   &&      # ==> beginning of a dataset
-            ((valtype(arg) <: Real)  ||
-             (valtype(arg) <: AbstractString))
+        taken = false
 
-            # Collect all data
-            accum = Vector{AbstractArray}()
-            while isa(arg, AbstractArray)  &&
-                ((valtype(arg) <: Real)    ||
-                 (valtype(arg) <: AbstractString))
+        if isa(arg, AbstractArray)
+            if nonmissingtype(eltype(arg)) != eltype(arg)
+                @assert nonmissingtype(eltype(arg)) <: AbstractFloat "Missing values are supported only on arrays of floats"
+                arg = replace(arg, missing => NaN)
+            end
+            tt = eltype(arg)
+
+            # Try to convert into Int, Float64 and String
+            if (tt  <: Integer)  &&  !(tt <: Int)
+                arg = convert(Array{Int}, arg)
+            elseif (tt  <: AbstractFloat)  &&  !(tt <: Float64)
+                arg = convert(Array{Float64}, arg)
+            elseif (tt  <: AbstractString)  &&  !(tt <: String)
+                arg = convert(Array{String}, arg)
+            end
+
+            tt = eltype(arg)
+            if  (tt <: Real)  ||
+                (tt <: AbstractString)
                 push!(accum, arg)
                 deleteat!(args, pos)
-                if pos <= length(args)
-                    arg = args[pos]
-                else
-                    break
+                taken = true
             end
         end
 
+        if !taken  ||  (pos > length(args))
+            if length(accum) > 0
                 mm = extrema(length.(accum))
-            if mm[1] == 0
-                # empty Dataset
+                if mm[1] == 0   # empty Dataset
                     @assert mm[1] == mm[2] "At least one input array is empty, while other(s) are not"
                     d = DatasetEmpty()
                 else
                     d = Dataset(accum)
                 end
                 insert!(args, pos, d)
+                empty!(accum)
             end
+
             pos += 1
         end
+    end
 
     # Fourth pass: collect PlotElement objects
     mid = 0
@@ -1240,12 +1291,6 @@ function driver(_args...; is3d=false)
         return source
     end
 
-    if length(_args) == 0
-        gp = getsession()
-        execall(gp)
-        return SessionID(gp.sid, true)
-    end
-
     (sid, doReset, doDump, elems) = parseArguments(_args...)
     gp = getsession(sid)
     doReset  &&  reset(gp)
@@ -1269,6 +1314,7 @@ function driver(_args...; is3d=false)
         end
     end
 
+    # Add commands and plot specifications
     for elem in elems
         (elem.mid > 0)  &&  setmulti(gp, elem.mid)
         gp.plots[gp.curmid].is3d = (is3d | elem.is3d)
@@ -1299,7 +1345,9 @@ function driver(_args...; is3d=false)
         end
     end
 
-    (doDump)  &&  (execall(gp))
+    if options.gpviewer  &&  doDump
+        execall(gp)
+    end
     return SessionID(gp.sid, doDump)
 end
 
@@ -1312,7 +1360,7 @@ end
 
 Return the **Gnuplot.jl** package version.
 """
-version() = v"1.2.0"
+version() = v"1.4.1"
 
 # ---------------------------------------------------------------------
 """
@@ -1476,66 +1524,82 @@ end
 
 # --------------------------------------------------------------------
 """
-    save(sid::Symbol; term="", output="")
-    save(sid::Symbol, script_filename::String, ;term="", output="")
-    save(; term="", output="")
-    save(script_filename::String ;term="", output="")
+    save([sid::Symbol]; term="", output="")
+    save([sid::Symbol,] mime::Type{T}; output="") where T <: MIME
+    save([sid::Symbol,] script_filename::String, ;term="", output="")
 
-Export a (multi-)plot into the external file name provided in the `output=` keyword.  The gnuplot terminal to use is provided through the `term=` keyword.
+Export a (multi-)plot into the external file name provided in the `output=` keyword.  The gnuplot terminal to use is provided through the `term=` keyword or the `mime` argument.  In the latter case the proper terminal is set according to the `Gnuplot.options.mime` dictionary.
 
 If the `script_filename` argument is provided a *gnuplot script* will be written in place of the output image.  The latter can then be used in a pure gnuplot session (Julia is no longer needed) to generate exactly the same original plot.
 
-If the `sid` argument is provided the operation applies to the corresponding session.
+If the `sid` argument is provided the operation applies to the corresponding session, otherwise the default session is considered.
+
+Example:
+```julia
+@gp hist(randn(1000))
+save(MIME"text/plain")
+save(term="pngcairo", output="output.png")
+save("script.gp")
+```
 """
 save(           ; kw...) = execall(getsession()   ; kw...)
 save(sid::Symbol; kw...) = execall(getsession(sid); kw...)
 save(             file::AbstractString; kw...) = savescript(getsession()   , file, kw...)
 save(sid::Symbol, file::AbstractString; kw...) = savescript(getsession(sid), file, kw...)
 
+save(mime::Type{T}; kw...) where T <: MIME = save(options.default, mime; kw...)
+function save(sid::Symbol, mime::Type{T}; kw...) where T <: MIME
+    if mime in keys(options.mime)
+        term = strip(options.mime[mime])
+        if term != ""
+            return save(sid; term=term, kw...)
+        end
+    end
+    @error "No terminal is defined for $mime.  Check `Gnuplot.options.mime` dictionary."
+end
+
+
 # ╭───────────────────────────────────────────────────────────────────╮
 # │                     Interfacing Julia's show                      │
 # ╰───────────────────────────────────────────────────────────────────╯
 # --------------------------------------------------------------------
 
-#=
-# Define a display that will be used when Gnuplot.jl is used
-# in the Julia REPL (see PGFPlotsX.jl).
-struct GnuplotDisplay <: AbstractDisplay end
-function __init__()
-    pushdisplay(GnuplotDisplay())
-    atreplinit(i -> begin
-        if PlotDisplay() in Base.Multimedia.displays
-            popdisplay(GnuplotDisplay())
+function show(obj::SessionID)
+    gp = getsession(obj.sid)
+    @info "Gnuplot session" sid=gp.sid datasets=length(gp.datas) plots=length(gp.plots)
+    nothing
 end
-        pushdisplay(GnuplotDisplay())
-    end)
+show(io::IO, gp::SessionID) = nothing
+
+
+function showable(mime::Type{T}, gp::SessionID) where T <: MIME
+    if gp.dump  &&  !options.gpviewer
+        if mime in keys(options.mime)
+            term = strip(options.mime[mime])
+            if term != ""
+                return true
             end
-function Base.display(d::GnuplotDisplay, gp::Session)
-    execall(gp)
-    return
         end
-=#
-Base.show(gp::SessionID) = nothing
-Base.show(io::IO, gp::SessionID) = nothing
-function Base.show(io::IO, ::MIME"image/svg+xml", gp::SessionID)
-    if gp.dump  &&  enableExportThroughShow()
-        tmpfile = tempname()*".svg"
-        save(gp.sid; term=options.term_svg, output=tmpfile)
-        write(io, read(tmpfile))
-        rm(tmpfile; force=true)
-    end
-    nothing
-end
-function Base.show(io::IO, ::MIME"image/png", gp::SessionID)
-    if gp.dump  &&  enableExportThroughShow()
-        tmpfile = tempname()*".png"
-        save(gp.sid; term=options.term_png, output=tmpfile)
-        write(io, read(tmpfile))
-        rm(tmpfile; force=true)
+    end
+    return false
+end
+
+function internal_show(io::IO, mime::Type{T}, gp::SessionID) where T <: MIME
+    if showable(mime, gp)
+        term = strip(options.mime[mime])
+        file = tempname()
+        save(gp.sid, term=term, output=file)
+        write(io, read(file))
+        rm(file; force=true)
     end
     nothing
 end
 
+show(io::IO, mime::MIME"image/svg+xml", gp::SessionID) = internal_show(io, typeof(mime), gp)
+show(io::IO, mime::MIME"image/png"    , gp::SessionID) = internal_show(io, typeof(mime), gp)
+show(io::IO, mime::MIME"text/html"    , gp::SessionID) = internal_show(io, typeof(mime), gp)
+
+
 # ╭───────────────────────────────────────────────────────────────────╮
 # │                     HIGH LEVEL FACILITIES                         │
 # ╰───────────────────────────────────────────────────────────────────╯
@@ -1648,28 +1712,52 @@ end
 
 
 """
-    palette(cmap::ColorScheme; rev=false)
-    palette(s::Symbol; rev=false)
+    palette_levels(cmap::ColorScheme; rev=false, smooth=false)
+    palette_levels(s::Symbol; rev=false, smooth=false)
 
-Convert a `ColorScheme` object into a string containing the gnuplot commands to set up the corresponding palette.
+Convert a `ColorScheme` object into a `Tuple{Vector{Float64}, Vector{String}, Int}` containing:
+- the numeric levels (between 0 and 1 included) corresponding to colors in the palette;
+- the corresponding colors (as hex strings);
+- the total number of different colors in the palette.
 
-If the argument is a `Symbol` it is interpreted as the name of one of the predefined schemes in [ColorSchemes](https://juliagraphics.github.io/ColorSchemes.jl/stable/basics/#Pre-defined-schemes-1). If `rev=true` the palette is reversed.
+If the argument is a `Symbol` it is interpreted as the name of one of the predefined schemes in [ColorSchemes](https://juliagraphics.github.io/ColorSchemes.jl/stable/basics/#Pre-defined-schemes-1).
+
+If `rev=true` the palette is reversed.  If `smooth=true` the palette is interpolated in 256 levels.
 """
-palette(s::Symbol; rev=false) = palette(colorschemes[s], rev=rev)
-function palette(cmap::ColorScheme; rev=false)
-    levels = Vector{String}()
-    for x in LinRange(0, 1, length(cmap.colors))
+palette_levels(s::Symbol; kwargs...) = palette_levels(colorschemes[s]; kwargs...)
+function palette_levels(cmap::ColorScheme; rev=false, smooth=false)
+    levels = OrderedDict{Float64, String}()
+    for x in LinRange(0, 1, (smooth  ?  256  : length(cmap.colors)))
         if rev
             color = get(cmap, 1-x)
         else
             color = get(cmap, x)
         end
-        push!(levels, "$x '#" * Colors.hex(color) * "'")
+        levels[x] = "#" * Colors.hex(color)
     end
-    return "set palette defined (" * join(levels, ", ") * ")\nset palette maxcol $(length(cmap.colors))\n"
+    return (collect(keys(levels)), collect(values(levels)), length(cmap.colors))
 end
 
 
+"""
+    palette(cmap::ColorScheme; rev=false, smooth=false)
+    palette(s::Symbol; rev=false, smooth=false)
+
+Convert a `ColorScheme` object into a string containing the gnuplot commands to set up the corresponding palette.
+
+If the argument is a `Symbol` it is interpreted as the name of one of the predefined schemes in [ColorSchemes](https://juliagraphics.github.io/ColorSchemes.jl/stable/basics/#Pre-defined-schemes-1).
+
+If `rev=true` the palette is reversed.  If `smooth=true` the palette is interpolated in 256 levels.
+"""
+function palette(values::Vector{Float64}, levels::Vector{String}, ncolors::Int)
+    str = string.(values) .* " '" .* levels .* "'"
+    return "set palette defined (" * join(str, ", ") * ")\nset palette maxcol $(ncolors)\n"
+end
+palette(s::Symbol; kwargs...) = palette(colorschemes[s]; kwargs...)
+palette(cmap::ColorScheme; kwargs...) =
+    palette(palette_levels(cmap; kwargs...)...)
+
+
 # --------------------------------------------------------------------
 """
     terminals()
@@ -1793,16 +1881,20 @@ function hist(v::Vector{T}; range=[NaN,NaN], bs=NaN, nbins=0, pad=true) where T
         if sum(hh.weights) < length(i)
             j = findall(v[i] .== range[2])
             @assert length(j) == (length(i) - sum(hh.weights))
+            if length(hh.weights) > 0
                 hh.weights[end] += length(j)
+            else
+                push!(hh.weights, length(j))
+            end
         end
     else
         hh = fit(Histogram, v[i], closed=:left)
     end
     @assert sum(hh.weights) == length(i)
     x = collect(hh.edges[1])
-    x = (x[1:end-1] .+ x[2:end]) ./ 2
+    binsize = isfinite(bs) ? bs : x[2] - x[1]
+    length(x) > 1 && (x = (x[1:end-1] .+ x[2:end]) ./ 2)
     h = hh.weights
-    binsize = x[2] - x[1]
     if pad
         x = [x[1]-binsize, x..., x[end]+binsize]
         h = [0, h..., 0]
@@ -1866,6 +1958,35 @@ function hist(v1::Vector{T1}, v2::Vector{T2};
     return Histogram2D(x1, x2, hh.weights, binsize1, binsize2)
 end
 
+# Allow missing values in input
+function hist(v::Vector{Union{Missing,T}}; kw...) where T <: Real
+    ii = findall(.!ismissing.(v))
+    @info "Neglecting missing values ($(length(v) - length(ii)))"
+    hist(convert(Vector{T}, v[ii]); kw...)
+end
+
+function hist(v1::Vector{Union{Missing,T1}}, v2::Vector{T2}; kw...) where {T1 <: Real, T2 <: Real}
+    ii = findall(.!ismissing.(v1)  .&
+                 .!ismissing.(v2)  )
+    @info "Neglecting missing values ($(length(v1) - length(ii)))"
+    hist(convert(Vector{T1}, v1[ii]), convert(Vector{T2}, v2[ii]), kw...)
+end
+
+function hist(v1::Vector{T1}, v2::Vector{Union{Missing, T2}}; kw...) where {T1 <: Real, T2 <: Real}
+    ii = findall(.!ismissing.(v1)  .&
+                 .!ismissing.(v2)  )
+    @info "Neglecting missing values ($(length(v1) - length(ii)))"
+    hist(convert(Vector{T1}, v1[ii]), convert(Vector{T2}, v2[ii]), kw...)
+end
+
+function hist(v1::Vector{Union{Missing,T1}}, v2::Vector{Union{Missing,T2}}; kw...) where {T1 <: Real, T2 <: Real}
+    ii = findall(.!ismissing.(v1)  .&
+                 .!ismissing.(v2)  )
+    @info "Neglecting missing values ($(length(v1) - length(ii)))"
+    hist(convert(Vector{T1}, v1[ii]), convert(Vector{T2}, v2[ii]), kw...)
+end
+
+
 
 # --------------------------------------------------------------------
 """
@@ -1935,6 +2056,8 @@ struct IsoContourLines
     paths::Vector{Path2d}
     data::Dataset
     z::Float64
+    prob::Float64
+end
 function IsoContourLines(paths::Vector{Path2d}, z)
     @assert length(z) == 1
     # Prepare Dataset object
@@ -1944,14 +2067,14 @@ struct IsoContourLines
         push!(data, "")
         push!(data, "")
     end
-        return new(paths, DatasetText(data), z)
-    end
+    return IsoContourLines(paths, DatasetText(data), z, NaN)
 end
 
 
 """
-    contourlines(x::Vector{Float64}, y::Vector{Float64}, z::Matrix{Float64}, cntrparam="level auto 10")
-    contourlines(h::Histogram2D, cntrparam="level auto 10")
+    contourlines(x, y, z, cntrparam="level auto 4")
+    contourlines(x, y, z, fractions)
+    contourlines(h::Histogram2D, ...)
 
 Compute paths of contour lines for 2D data, and return a vector of [`IsoContourLines`](@ref) object.
 
@@ -1959,16 +2082,19 @@ Compute paths of contour lines for 2D data, and return a vector of [`IsoContourL
     This feature is not available in *dry* mode and will raise an error if used.
 
 # Arguments:
-- `x`, `y`: Coordinates;
-- `z`: the levels on which iso contour lines are to be calculated
-- `cntrparam`: settings to compute contour line paths (see gnuplot documentation for `cntrparam`).
+- `x`, `y` (as `AbstractVector{Float64}`): Coordinates;
+- `z::AbstractMatrix{Float64}`: the levels on which iso-contour lines are to be calculated;
+- `cntrparam::String`: settings to compute contour line paths (see gnuplot documentation for `cntrparam`);
+- `fractions::Vector{Float64}`: compute contour lines encompassing these fractions of total counts;
+- `h::Histogram2D`: use histogram bins and counts to compute contour lines.
+
 
 # Example
 ```julia
-x = randn(5000);
-y = randn(5000);
+x = randn(10^5);
+y = randn(10^5);
 h = hist(x, y, nbins1=20, nbins2=20);
-clines = contourlines(h, "levels discrete 15, 30, 45");
+clines = contourlines(h, "levels discrete 500, 1500, 2500");
 
 # Use implicit recipe
 @gp clines
@@ -1978,14 +2104,45 @@ clines = contourlines(h, "levels discrete 15, 30, 45");
 for i in 1:length(clines)
     @gp :- clines[i].data "w l t '\$(clines[i].z)' lw \$i dt \$i"
 end
+
+# Calculate probability within 0 < r < σ
+p(σ) = round(1 - exp(-(σ^2) / 2), sigdigits=3)
+
+# Draw contour lines at 1, 2 and 3 σ
+clines = contourlines(h, p.(1:3));
+@gp palette(:beach, smooth=true, rev=true) "set grid front" "set size ratio -1" h clines
 ```
 """
 contourlines(h::Histogram2D, args...) = contourlines(h.bins1, h.bins2, h.counts, args...)
-function contourlines(x::Vector{Float64}, y::Vector{Float64}, z::Matrix{Float64},
-                      cntrparam="level auto 10")
+function contourlines(x::AbstractVector{Float64}, y::AbstractVector{Float64}, z::AbstractMatrix{Float64},
+                      fraction::Vector{Float64})
+    @assert minimum(fraction) > 0
+    @assert maximum(fraction) < 1
+    @assert length(fraction) >= 1
+    sorted_fraction = sort(fraction, rev=true)
+
+    i = sortperm(z[:], rev=true)
+    topfrac = cumsum(z[i]) ./ sum(z)
+    selection = Int[]
+    for f in sorted_fraction
+        push!(selection, minimum(findall(topfrac .>= f)))
+    end
+    levels = z[i[selection]]
+    clines = contourlines(x, y, z, "levels discrete " * join(string.(levels), ", "))
+    @assert issorted(getfield.(clines, :z))
+
+    if  length(clines) == length(fraction)
+        out = [IsoContourLines(clines[i].paths, clines[i].data, clines[i].z,
+                               sorted_fraction[i]) for i in 1:length(clines)]
+        return out
+    end
+    return clines
+end
+
+function contourlines(x::AbstractVector{Float64}, y::AbstractVector{Float64}, z::AbstractMatrix{Float64},
+                      cntrparam="level auto 4")
     lines = gp_write_table("set contour base", "unset surface",
                            "set cntrparam $cntrparam", x, y, z, is3d=true)
-
     level = NaN
     path = Path2d()
     paths = Vector{Path2d}()
@@ -2031,6 +2188,105 @@ function contourlines(x::Vector{Float64}, y::Vector{Float64}, z::Matrix{Float64}
 end
 
 
+# --------------------------------------------------------------------
+"""
+    dgrid3d(x, y, z, opts=""; extra=true)
+
+Interpolate non-uniformly spaced 2D data onto a regular grid.
+
+!!! note
+    This feature is not available in *dry* mode and will raise an error if used.
+
+# Arguments:
+- `x`, `y`, `z` (as `AbstractVector{Float64}`): coordinates and values of the function to interpolate;
+- `opts`: interpolation settings (see gnuplot documentation for `dgrid3d`);
+- `extra`: if `true` (default) compute inerpolated values in all regions, even those which are poorly constrained by input data (namely, extrapolated values).  If `false` set these values to `NaN`.
+
+# Return values:
+A tuple with `x` and `y` coordinates on the regular grid (as `Vector{Float64}`), and `z` containing interpolated values (as `Matrix{Float64}`).
+
+# Example
+```julia
+x = (rand(200) .- 0.5) .* 3;
+y = (rand(200) .- 0.5) .* 3;
+z = exp.(-(x.^2 .+ y.^2));
+
+# Interpolate on a 20x30 regular grid with splines
+gx, gy, gz = dgrid3d(x, y, z, "20,30 splines")
+
+@gsp "set size ratio -1" "set xyplane at 0" xlab="X" ylab="Y" :-
+@gsp :-  x  y  z "w p t 'Scattered data' lc pal"
+@gsp :- gx gy gz "w l t 'Interpolation on a grid' lc pal"
+```
+!!! warn
+    The `splines` algorithm may be very slow on large datasets.  An alternative option is to use a smoothing kernel, such as `gauss`:
+
+```julia
+x = randn(2000) .* 0.5;
+y = randn(2000) .* 0.5;
+rsq = x.^2 + y.^2;
+z = exp.(-rsq) .* sin.(y) .* cos.(2 * rsq);
+
+@gsp "set size ratio -1" palette(:balance, smooth=true) "set view map" "set pm3d" :-
+@gsp :- "set multiplot layout 1,3" xr=[-2,2] yr=[-2,2] :-
+@gsp :- 1 tit="Scattered data"  x  y  z "w p notit lc pal"
+
+# Show extrapolated values
+gx, gy, gz = dgrid3d(x, y, z, "40,40 gauss 0.1,0.1")
+@gsp :- 2 tit="Interpolation on a grid\\\\n(extrapolated values are shown)"  gx gy gz "w l notit lc pal"
+
+# Hide exrapolated values
+gx, gy, gz = dgrid3d(x, y, z, "40,40 gauss 0.1,0.1", extra=false)
+@gsp :- 3 tit="Interpolation on a grid\\\\n(extrapolated values are hidden)" gx gy gz "w l notit lc pal"
+```
+"""
+function dgrid3d(x::AbstractVector{Float64},
+                 y::AbstractVector{Float64},
+                 z::AbstractVector{Float64},
+                 opts::String="";
+                 extra=true)
+    c = Gnuplot.gp_write_table("set dgrid3d $opts", x, y, z, is3d=true)
+    gx = Vector{Float64}()
+    gy = Vector{Float64}()
+    gz = Vector{Float64}()
+    ix = 0
+    iy = 0
+    for l in c
+        l = string(strip(l))
+        if l == "# x y z type"
+            ix += 1
+            iy = 1
+        else
+            (l == "" )     &&  continue
+            (l[1] == '#')  &&  continue
+            n = Meta.parse.(split(l)[1:3])
+            (iy == 1)  &&  push!(gx, n[1])
+            (ix == 1)  &&  push!(gy, n[2])
+            if n[3] == :NaN
+                push!(gz, NaN)
+            else
+                push!(gz, n[3])
+            end
+            iy += 1
+        end
+    end
+    gz = collect(reshape(gz, length(gy), length(gx))')
+    if !extra
+        dx = abs(gx[2]-gx[1]) / 2
+        dy = abs(gy[2]-gy[1]) / 2
+        for ix in 1:length(gx)
+            for iy in 1:length(gy)
+                n = length(findall(((gx[ix] - dx) .< x .< (gx[ix] + dx))  .&
+                                   ((gy[iy] - dy) .< y .< (gy[iy] + dy))))
+                (n == 0)  &&  (gz[ix, iy] = NaN)
+            end
+        end
+    end
+    return (gx, gy, gz)
+end
+
+
+
 # ╭───────────────────────────────────────────────────────────────────╮
 # │                        GNUPLOT REPL                               │
 # ╰───────────────────────────────────────────────────────────────────╯
@@ -2092,7 +2348,7 @@ function gpvars(sid::Symbol)
         if length(s) == 2
             key = Symbol(s[1])
             if s[2][1] == '"'
-                out[key] = s[2][2:end-1]
+                out[key] = s[2][2:prevind(s[2], end, 1)]
             else
                 try
                     out[key] = Meta.parse(s[2])

src/recipes.jl

@@ -29,7 +29,12 @@ recipe(h::Histogram2D) =
 
 Implicit recipes to visualize iso-contour lines.
 """
-recipe(c::IsoContourLines) = PlotElement(data=c.data, plot="w l t '$(c.z)'")
+function recipe(c::IsoContourLines)
+    if isnan(c.prob)
+        return PlotElement(data=c.data, plot="w l t '$(c.z)'")
+    end
+    return PlotElement(data=c.data, plot="w l t '$(round(c.prob * 100, sigdigits=6))%'")
+end
 recipe(v::Vector{IsoContourLines}) = recipe.(v)
 
 

test/runtests.jl

@@ -1,16 +1,15 @@
 using Test, Gnuplot
 try
-    @info "Gnuplot version: " * string(Gnuplot.gpversion())
+    @info "Gnuplot.jl version: " * string(Gnuplot.version())
+    @info "gnuplot    version: " * string(Gnuplot.gpversion())
 catch
     Gnuplot.options.dry = true
 end
-Gnuplot.options.term = "unknown"
+Gnuplot.options.gpviewer = true
 
 x = [1, 2, 3]
 y = [4, 5, 6]
 
-
-
 s = Gnuplot.arrays2datablock(x)
 @test all(s .== [" 1"   ,
                  " 2"   ,
@@ -90,10 +89,30 @@ s = Gnuplot.arrays2datablock(1:3, 1:3, ["One", "Two", "Three"])
 
 
 #-----------------------------------------------------------------
+dummy = palette_names()
 pal = palette(:deepsea)
 @test pal == "set palette defined (0.0 '#2B004D', 0.25 '#4E0F99', 0.5 '#3C54D4', 0.75 '#48A9F8', 1.0 '#C5ECFF')\nset palette maxcol 5\n"
-ls = linetypes(:Set1_5)
-@test ls == "unset for [i=1:256] linetype i\nset linetype 1 lc rgb '#E41A1C' lw 1 dt solid pt 1 ps 1\nset linetype 2 lc rgb '#377EB8' lw 1 dt solid pt 2 ps 1\nset linetype 3 lc rgb '#4DAF4A' lw 1 dt solid pt 3 ps 1\nset linetype 4 lc rgb '#984EA3' lw 1 dt solid pt 4 ps 1\nset linetype 5 lc rgb '#FF7F00' lw 1 dt solid pt 5 ps 1\nset linetype cycle 5\n"
+ls = linetypes(:Set1_5, lw=1.5, ps=2)
+@test ls == "unset for [i=1:256] linetype i\nset linetype 1 lc rgb '#E41A1C' lw 1.5 dt solid pt 1 ps 2\nset linetype 2 lc rgb '#377EB8' lw 1.5 dt solid pt 2 ps 2\nset linetype 3 lc rgb '#4DAF4A' lw 1.5 dt solid pt 3 ps 2\nset linetype 4 lc rgb '#984EA3' lw 1.5 dt solid pt 4 ps 2\nset linetype 5 lc rgb '#FF7F00' lw 1.5 dt solid pt 5 ps 2\nset linetype cycle 5\n"
+
+dummy = terminals()
+# if "sixelgd" in terminals()
+#     Gnuplot.options.term = "sixelgd enhanced"
+# elseif "sixel" in terminals()
+#     Gnuplot.options.term = "sixel enhanced"
+# elseif "dumb" in terminals()
+#     Gnuplot.options.term = "dumb"
+# else
+#     Gnuplot.options.term = "unknown"
+# end
+# Gnuplot.quitall()
+
+# Force unknown on Travis CI
+Gnuplot.options.term = "unknown"
+
+@gp 1:9
+@info "using terminal: " terminal()
+#test_terminal("unknown")
 
 #-----------------------------------------------------------------
 # Test wth empty dataset
@@ -247,8 +266,19 @@ Gnuplot.quitall()
 	"splot x9, v, (u<0.5) ? -1 : sinc(x9,v) notitle")
 
 
+
+x = randn(5000);
+y = randn(5000);
+h = hist(x, y, nbins1=20, nbins2=20);
+clines = contourlines(h, "levels discrete 15, 30, 45");
+@gp clines
+@gp "set size ratio -1"
+for i in 1:length(clines)
+    @gp :- clines[i].data "w l t '$(clines[i].z)' lw $i dt $i"
+end
+
+
 Gnuplot.options.verbose = true
-# Gnuplot.options.term = "sixel"  not vailable in Travis CI
 @gp randn(10^6) randn(10^6)
 @gp :- 0. 0.
 Gnuplot.quit(:default)
@@ -256,6 +286,8 @@ Gnuplot.quit(:default)
 Gnuplot.options.dry = true
 @gp hist(randn(1000))
 
-t = terminals()
+# Various hist() corner cases
+@gp hist([1,2,3], bs=2)
+@gp hist([1,1,1], bs=1)
 
 Gnuplot.quitall()