Add a new monotone interpolation plot handler

This handler converts each plot stream command into a curveto command, except for the first, which is converted to the previously specified action. The curveto control points are computed and adjusted so that the resulting curve respects the monotonicity of the points. If a point is a local extremum of the point set, it will also be a local extremum of the generated curve. See https://en.wikipedia.org/wiki/Monotone_cubic_interpolation for an explanation of the algorithm. The last part uses \pgfpathcurveto instead of using Hermitte base functions. The plot handler is also made available to datavisualization through the `smooth monotone line` visualizer option. Signed-off-by: Julien '_FrnchFrgg_' RIVAUD <[email protected]>
pgf-tikz · Sep 6, 2024 · fdde2c8 · fdde2c8
1 parent d43fb10
commit fdde2c8
Show file tree

Hide file tree

Showing 5 changed files with 164 additions and 0 deletions.
diff --git a/doc/generic/pgf/CHANGELOG.md b/doc/generic/pgf/CHANGELOG.md
@@ -11,6 +11,7 @@ The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/).
 ### Added
 
 - Add `RGB` and `gray` color model support for ConTeXt #1130
+- Add a new monotonic interpolation plot handler #1358
 
 ### Fixed
 

diff --git a/doc/generic/pgf/pgfmanual-en-library-plot-handlers.tex b/doc/generic/pgf/pgfmanual-en-library-plot-handlers.tex
@@ -97,6 +97,38 @@ \subsection{Curve Plot Handlers}
 \end{command}
 
 
+\begin{command}{\pgfplothandlermonotone}
+    This handler will issue a |\pgfpathcurveto| command for each point of the
+    plot, \emph{except} possibly for the first. As for the line-to handler,
+    what happens with the first point can be specified using
+    |\pgfsetmovetofirstplotpoint| or |\pgfsetlinetofirstplotpoint|.
+
+    Obviously, the |\pgfpathcurveto| command needs, in addition to the points
+    on the path, some control points. These are computed using The
+    Fritsch-Carlson monotone cubic interpolation algorithm which roughly does
+    the following:
+
+    \begin{enumerate}
+        \item Initialize tangent slopes at each input point as the
+            average of the slopes of the neighboring secant lines.
+        \item If the secants have opposite signes (because this is a local
+            extremum), set the tangent slope to~$0$ to get an horizontal
+            tangent and prevent overshoot. Same if one of the neighboring
+            secants is horizontal.
+        \item Adjust the tangent slopes to ensure strict monotonicity.
+        \item Place the control points on the computed tangent lines,
+            at one third of the $x$~distance of a secant line: if $a$, $b$
+            and~$c$ are three consecutive input points and $m$~is the
+            slope computed at~$B$, then the control points around~$b$ will
+            be $u$~and~$v$ where $x_u = x_b - \frac{1}{3}(x_b-x_a)$,
+            $y_u = y_b - \frac{1}{3} m (x_b-x_a)$,
+            $x_v = x_b + \frac{1}{3}(x_c-x_b)$ and
+            $y_v = y_b + \frac{1}{3} m (x_c-x_b)$.
+    \end{enumerate}
+\end{command}
+
+
+
 \subsection{Constant Plot Handlers}
 
 There are several plot handlers which produce piecewise constant interpolations

diff --git a/.../pgf/frontendlayer/tikz/libraries/datavisualization/tikzlibrarydatavisualization.code.tex b/.../pgf/frontendlayer/tikz/libraries/datavisualization/tikzlibrarydatavisualization.code.tex
@@ -1538,6 +1538,7 @@
   smooth cycle/.style={@set={\pgfplothandlerclosedcurve}{default label in legend closed path}},
   straight line/.style={@set={\pgfplothandlerlineto}{default label in legend path}},
   straight cycle/.style={@set={\pgfplothandlerpolygon}{default label in legend closed path}},
+  smooth monotone line/.style={@set={\pgfplothandlermonotone}{default label in legend path}},
   polygon/.style={straight cycle},% alias
   gap line/.style={@set={\pgfplothandlergaplineto}{default label in legend path}},
   gap cycle/.style={@set={\pgfplothandlergapcycle}{gap circular label in legend line}}

diff --git a/tex/generic/pgf/frontendlayer/tikz/tikz.code.tex b/tex/generic/pgf/frontendlayer/tikz/tikz.code.tex
@@ -1232,6 +1232,7 @@
 % Plot options
 \tikzoption{smooth}[]{\let\tikz@plot@handler=\pgfplothandlercurveto}%
 \tikzoption{smooth cycle}[]{\let\tikz@plot@handler=\pgfplothandlerclosedcurve}%
+\tikzoption{smooth monotone}[]{\let\tikz@plot@handler=\pgfplothandlermonotone}%
 \tikzoption{sharp plot}[]{\let\tikz@plot@handler\pgfplothandlerlineto}%
 \tikzoption{sharp cycle}[]{\let\tikz@plot@handler\pgfplothandlerpolygon}%
 

diff --git a/tex/generic/pgf/libraries/pgflibraryplothandlers.code.tex b/tex/generic/pgf/libraries/pgflibraryplothandlers.code.tex
@@ -269,6 +269,135 @@
 }%
 
 
+% This handler converts each plot stream command into a curveto
+% command, except for the first, which is converted to the previously
+% specified action. The curveto control points are computed and adjusted
+% so that the resulting curve respects the monotonicity of the points.
+% If a point is a local extremum of the point set, it will also be
+% a local extremum of the generated curve.
+% This handler treats the x and y direction differently, and expects
+% x values to be monotone.
+% See https://en.wikipedia.org/wiki/Monotone_cubic_interpolation for an
+% explanation of the algorithm. The last part uses \pgfpathcurveto
+% instead of using Hermitte base functions.
+\pgfdeclareplothandler{\pgfplothandlermonotone}{}{%
+    point macro=\pgf@plot@monotone@dopointone,
+    jump macro=\pgf@plot@monotone@dojump,
+    end macro=\pgf@plot@monotone@doend
+}%
+
+\def\pgf@plot@monotone@dojump{%
+    \pgf@plot@monotone@doend
+    \global\pgf@plot@startedfalse
+    \global\let\pgf@plotstreampoint\pgf@plot@monotone@dopointone
+}%
+
+\def\pgf@plot@monotone@dopointone#1{%
+    \pgf@process{#1}%
+    \pgf@plot@first@action{\pgfqpoint{\pgf@x}{\pgf@y}}%
+    \xdef\pgf@plot@monotone@pointone{\noexpand\pgfqpoint{\the\pgf@x}{\the\pgf@y}}%
+    \global\let\pgf@plotstreampoint=\pgf@plot@monotone@dopointtwo
+}%
+
+\def\pgf@plot@monotone@dopointtwo#1{%
+    \pgf@process{#1}%
+    \xdef\pgf@plot@monotone@pointtwo{\noexpand\pgfqpoint{\the\pgf@x}{\the\pgf@y}}%
+    % compute vector:
+    \pgf@xa=\pgf@x%
+    \pgf@ya=\pgf@y%
+    \pgf@process{\pgf@plot@monotone@pointone}%
+    \advance\pgf@xa by-\pgf@x
+    \advance\pgf@ya by-\pgf@y
+    % save left delta
+    \xdef\pgf@plot@monotone@deltaleft{\pgf@sys@tonumber\pgf@xa}%
+    % compute slope
+    \pgfmathdivide@{\the\pgf@ya}{\the\pgf@xa}%
+    \global\let\pgf@plot@monotone@rateleft\pgfmathresult
+    % init tangent one
+    \global\let\pgf@plot@monotone@slopeone\pgf@plot@monotone@rateleft
+    % prepare for next step
+    \global\let\pgf@plotstreampoint=\pgf@plot@monotone@dopointthree
+    \global\pgf@plot@startedtrue%
+}%
+
+\def\pgf@plot@monotone@dopointthree#1{%
+    \pgf@process{#1}%
+    \xdef\pgf@plot@monotone@pointthree{\noexpand\pgfqpoint{\the\pgf@x}{\the\pgf@y}}%
+    % compute vector:
+    \pgf@xa=\pgf@x%
+    \pgf@ya=\pgf@y%
+    \pgf@process{\pgf@plot@monotone@pointtwo}%
+    \advance\pgf@xa by-\pgf@x
+    \advance\pgf@ya by-\pgf@y
+    % save right delta
+    \xdef\pgf@plot@monotone@deltaright{\pgf@sys@tonumber\pgf@xa}%
+    % compute slope
+    \pgfmathdivide@{\the\pgf@ya}{\the\pgf@xa}%
+    \global\let\pgf@plot@monotone@rateright\pgfmathresult
+    % compute tangent slope
+    \pgfmathadd@{\pgf@plot@monotone@rateleft}{\pgf@plot@monotone@rateright}%
+    \pgf@xa=\pgfmathresult pt\relax
+    \pgf@xa=0.5\pgf@xa
+    % Fix slopes & draw curve
+    \pgf@plot@monotonecurveto@fixdraw
+    % Prepare next:
+    \global\let\pgf@plot@monotone@slopeone\pgf@plot@monotone@slopetwo
+    \global\let\pgf@plot@monotone@pointone\pgf@plot@monotone@pointtwo
+    \global\let\pgf@plot@monotone@pointtwo\pgf@plot@monotone@pointthree
+    \global\let\pgf@plot@monotone@deltaleft\pgf@plot@monotone@deltaright
+    \global\let\pgf@plot@monotone@rateleft\pgf@plot@monotone@rateright
+}%
+
+\def\pgf@plot@monotonecurveto@fixdraw{%
+    % fix extremums
+    % \pgf@xa is slope two
+    \pgf@ya=\pgf@plot@monotone@rateleft pt\relax
+    \pgf@yb=\pgf@plot@monotone@rateright\pgf@ya
+    \ifdim\pgf@yb<0.0001pt\relax
+        \pgf@xa=0pt\relax
+    \else
+    % fix motonicity
+        \pgf@ya=3\pgf@ya
+        \ifdim\pgf@xa>0pt\relax
+            \ifdim\pgf@xa>\pgf@ya \pgf@xa=\pgf@ya\fi
+        \else
+            \ifdim\pgf@xa<\pgf@ya \pgf@xa=\pgf@ya\fi
+        \fi
+        \pgf@ya=\pgf@plot@monotone@rateright pt\relax
+        \pgf@ya=3\pgf@ya
+        \ifdim\pgf@xa>0pt\relax
+            \ifdim\pgf@xa>\pgf@ya \pgf@xa=\pgf@ya\fi
+        \else
+            \ifdim\pgf@xa<\pgf@ya \pgf@xa=\pgf@ya\fi
+        \fi
+    \fi
+    \xdef\pgf@plot@monotone@slopetwo{\pgf@sys@tonumber\pgf@xa}%
+    % compute points for left bezier
+    % (\pgf@x,\pgf@y) is point two
+    \pgf@ya=\pgf@plot@monotone@deltaleft pt\relax
+    \pgf@ya=0.333333\pgf@ya
+    \advance \pgf@x by -\pgf@ya
+    \advance \pgf@y by -\pgf@plot@monotone@slopetwo\pgf@ya
+    \xdef\pgf@plot@monotone@supporttwo{\noexpand\pgfqpoint{\the\pgf@x}{\the\pgf@y}}%
+    \pgf@process{\pgf@plot@monotone@pointone}%
+    \advance \pgf@x by \pgf@ya
+    \advance \pgf@y by \pgf@plot@monotone@slopeone\pgf@ya
+    \xdef\pgf@plot@monotone@supportone{\noexpand\pgfqpoint{\the\pgf@x}{\the\pgf@y}}%
+    % draw left curve:
+    \pgfpathcurveto{\pgf@plot@monotone@supportone}{\pgf@plot@monotone@supporttwo}{\pgf@plot@monotone@pointtwo}%
+}%
+
+\def\pgf@plot@monotone@doend{%
+    \ifpgf@plot@started
+        % fixdraw needs (\pgf@x,\pgf@y) as point two
+        % and \pgf@xa as slope two
+        \pgf@process{\pgf@plot@monotone@pointtwo}%
+        \pgf@xa=\pgf@plot@monotone@rateleft pt\relax
+        \pgf@plot@monotonecurveto@fixdraw
+    \fi
+}%
+
+