Thursday 4 April 2019

Graphing Derivatives

A Series Of Points

Previously:
Differentiating f(x)^g(x)
The Differentiator
Note: The latest code can be found at https://github.com/dogbiscuituk/Sid/.
I've been inundated by a request from one reader, for a concrete example illustrating the plotting of graphs and their derivatives using my Differentiator class. This here quick and filthy Winforms application is my terse response. It uses no modern graphics library, just oldschool GDI+. Starting with a new Winforms app containing the FormulaBuilder class library from last time, add a new FormulaGrapher library to contain the two classes I'm about to describe.

Visual graphing controls normally provide a capability to plot multiple traces on the same canvas, so the first thing we need is a class to represent just one of these traces. Again following tradition, we'll call this the Series class. Its main responsibility will be to store the Expression used to calculate the points on its curve. This is passed into the constructor as the Formula parameter, and immediately lambda-compiled into the callable function Func.

For demo purposes, we'll also have this class store the resolution StepCount (the number of tiny line segments in the output, normally one less than the number of points computed); the PenColour used to draw the trace; and a cache containing the current plot limits, together with all of the computed point coordinates. After all, you might have any number of nested trigonometric functions to compute per point, and as long as the logical plot limits don't change, there's no need to recalculate these in response to e.g. resizing the output canvas in terms of physical pixels, or switching between screen and printer output.

using System;
using System.Collections.Generic;
using System.Drawing;
using System.Linq.Expressions;
using FormulaBuilder;

namespace FormulaGrapher
{
    public class Series
    {
        public Series(Expression formula, int stepCount)
        {
            Func = formula.AsFunction();
            StepCount = stepCount;
        }

        public Color PenColour { get; set; }

        public void Draw(Graphics g, RectangleF limits, float penWidth)
        {
            if (Limits != limits)
            {
                Limits = limits;
                ComputePoints();
            }
            using (var pen = new Pen(PenColour, penWidth))
                PointLists.ForEach(p => g.DrawLines(pen, p.ToArray()));
        }

        private Func<double, double> Func;
        private int StepCount;
        private RectangleF Limits;
        private List<List<PointF>> PointLists = new List<List<PointF>>();

        private void ComputePoints()
        {
            PointLists.Clear();
            List<PointF> points = null;
            float
                x1 = Limits.Left, y1 = Limits.Top, y2 = Limits.Bottom,
                w = Limits.Width, h = 8 * Limits.Height;
            var skip = true;
            for (var step = 0; step <= StepCount; step++)
            {
                float x = x1 + step * w / StepCount, y = (float)Func(x);
                if (float.IsInfinity(y) || float.IsNaN(y) || y < y1 - h || y > y2 + h)
                    skip = true;
                else
                {
                    if (skip)
                    {
                        skip = false;
                        points = new List<PointF>();
                        PointLists.Add(points);
                    }
                    points.Add(new PointF(x, y));
                }
            }
            // Every segment of the trace must include at least 2 points.
            PointLists.RemoveAll(p => p.Count < 2);
        }
    }
}

So our first class essentially consists of just this Draw method, which first checks whether the supplied logical limits are the same as the previous call, if any. When it decides recomputation is required, it iterates over x values between the supplied limits, splitting the interval into StepCount equal segments, and computing and storing the corresponding y values.

One interesting aspect of the structure used to store these points is its type: a List of Lists of PointFs. This particular contortion is needed because some graphs have discontinuities, or worse still, whole interior intervals on which they are entirely undefined. Hence the need to break the list into segments, whenever a y value is found to be infinite, or not a number, or just too damned big for our paltry coordinate system.

A Group Of Series

Next up, we have our supervisor class, which I've called Graph. This will take full control of its multiple individual Series, so we will only need to interact directly with the Graph class itself. Let's start by examining its constructor parameters.
  • limits is a RectangleF structure specifying the logical window for the plots. Suppose we want the x domain to go from -9 to +9, and y to range from -7 to +7. The two minimum (negative) values can be supplied as-is, but the GDI+ convention requires width and height values to specify the window size, rather than the corresponding maximum (positive) values of x and y. So in our example, we would pass the window bounds as new RectangleF(-9, -7, 18, 14).
  • stepCount has already been described; it is simply passed through to the Series objects for use when drawing.
As with the previous Series class, the interesting details are in the Draw method and its minions:

using System;
using System.Collections.Generic;
using System.Diagnostics;
using System.Drawing;
using System.Drawing.Drawing2D;
using System.Linq.Expressions;
using FormulaBuilder;

namespace FormulaGrapher
{
    public class Graph
    {
        public Graph(RectangleF limits, int stepCount)
        {
            Limits = limits;
            StepCount = stepCount;
        }

        public Series AddSeries(Expression formula)
        {
            Debug.WriteLine(formula.AsString());
            var series = new Series(formula, StepCount);
            Series.Add(series);
            return series;
        }

        public void Clear()
        {
            Series.Clear();
        }

        public void Draw(Graphics g, Rectangle r)
        {
            if (r.Width == 0 || r.Height == 0)
                return; // Nothing to draw!
            g.SmoothingMode = SmoothingMode.AntiAlias;
            g.Transform = GetMatrix(r);
            g.FillRectangle(Brushes.LightYellow, Limits);
            var penWidth = (Limits.Width / r.Width + Limits.Height / r.Height);
            DrawGrid(g, penWidth);
            Series.ForEach(s => s.Draw(g, Limits, penWidth));
        }

        public PointF ScreenToGraph(Point p, Rectangle r)
        {
            var points = new[] { new PointF(p.X, p.Y) };
            var matrix = GetMatrix(r);
            matrix.Invert();
            matrix.TransformPoints(points);
            return points[0];
        }

        private RectangleF Limits;
        private int StepCount;
        private List<Series> Series = new List<Series>();

        private void DrawGrid(Graphics g, float penWidth)
        {
            float x1 = Limits.Left, y1 = Limits.Bottom, x2 = Limits.Right, y2 = Limits.Top;
            using (Pen gridPen = new Pen(Color.LightGray, penWidth) { DashStyle = DashStyle.Dot },
                axisPen = new Pen(Color.DarkGray, penWidth))
            using (var font = new Font("Arial", 5 * penWidth))
            using (var format = new StringFormat(StringFormat.GenericTypographic) { Alignment = StringAlignment.Far })
            {
                var brush = Brushes.DarkGray;
                for (int phase = 0; phase < 4; phase++)
                {
                    var vertical = (phase & 1) != 0;
                    if (phase < 2)
                    {
                        var size = Math.Log10(Math.Abs(y2 - y1));
                        var step = Math.Floor(size);
                        size -= step;
                        step = (size < 0.3 ? 2 : size < 0.7 ? 5 : 10) * Math.Pow(10, step - 1);
                        for (var y = step; y <= Math.Max(Math.Abs(y1), Math.Abs(y2)); y += step)
                        {
                            DrawGridLine(g, gridPen, font, brush, x1, x2, (float)y, format, vertical);
                            DrawGridLine(g, gridPen, font, brush, x1, x2, -(float)y, format, vertical);
                        }
                    }
                    else
                        DrawGridLine(g, axisPen, font, brush, x1, x2, 0, format, vertical);
                    var t = x1; x1 = y1; y1 = t;
                    t = x2; x2 = y2; y2 = t;
                }
            }
        }

        private void DrawGridLine(Graphics g, Pen pen, Font font, Brush brush,
            float x1, float x2, float y, StringFormat format, bool vertical)
        {
            float x = 0, y1 = y, y2 = y, z = -y;
            if (vertical)
            {
                var t = x; x = y; y = t;
                t = x1; x1 = y1; y1 = t;
                t = x2; x2 = y2; y2 = t;
                z = -z;
            }
            g.DrawLine(pen, x1, y1, x2, y2);
            g.ScaleTransform(1, -1);
            g.DrawString(z.ToString(), font, brush, x - pen.Width, y + pen.Width, format);
            g.ScaleTransform(1, -1);
        }

        private Matrix GetMatrix(Rectangle r)
        {
            return new Matrix(Limits,
                new[] { new PointF(r.Left, r.Bottom), new PointF(r.Right, r.Bottom), r.Location });
        }
    }
}

The GDI+ Transform property of a Graphics context is a Matrix, which handles all conversions between our logical units, e.g. the 'one' that's the sine of half a pie (π/2), and physical units, e.g. the 'one' that's the size of a pixel on the display device. Here we basically set the Transform to equal the ratio between these coordinate systems, and forget it...

...until, that is, we come to specify a pen width for our plots. Here we can't just say 1, hoping for a single pixel width line. What we would get instead is a line as thick as the height of the sine graph we're trying to plot, making it look as if drawn with a giant crayon. Go ahead, try it, it's comical! The solution is just to scale the desired pixel pen width by the inverse of the coordinate system ratio, so undoing the effect of the Transform on our pen size.

DrawGrid() prepares the canvas, giving it the appearance of a musty old sheet of graph paper with x-y axes, and grid lines with numerical labels. Regardless of whether you scale your limits up or down by a factor of 100, each axis is divided into between 4 and 10 sensibly sized strips, by the magic of logarithms (initiates will quickly spot the usual four lines of code doing this essential work). The strips are then labelled in ones, twos or fives, as appropriate. A total of four phases are used to draw the x grid, y grid, x axis, and y axis in turn, ensuring that the grid is always behind the axes (and the axes behind the traces, which will be drawn next).

Note the inversions of the y axis near the end of DrawGridLine(). Since our graph's y values increase from bottom up, while device y values increase from top down, there is necessarily a flip in the transform used to map one to the other. If ignored, this flip would cause all text to be printed upside down! One of several simple solutions is to invert the y axis before printing text, and then again after, to restore the original coordinate mapping.

GetMatrix() has been factored out from the original Draw() method, to let external clients access logical coordinates via ScreenToGraph(). The main form below will use this feature to show on a ToolTip the updated logical position, i.e. graph world coordinates, of the moving mouse.
A Form Of Pencils

Finally it's time to add a main form to our Winforms app. Give it a context menu with an item for each formula you want to draw, following the pattern shown here. The code in the DrawGraphs() method will draw your selected function in black, then on the same canvas it will draw the first, second and third derivatives in red, green and blue respectively.

using System.Drawing;
using System.Linq.Expressions;
using System.Windows.Forms;
using FormulaBuilder;
using FormulaGrapher;

namespace Sid
{
    public partial class MainForm : Form
    {
        public MainForm()
        {
            InitializeComponent();
        }

        private void MainForm_Paint(object sender, PaintEventArgs e)
        {
            DrawGraph(e.Graphics);
        }

        private void MainForm_Resize(object sender, System.EventArgs e)
        {
            Invalidate();
        }

        private void MainForm_MouseMove(object sender, MouseEventArgs e)
        {
            var p = Graph.ScreenToGraph(e.Location, ClientRectangle);
            ToolTip.SetToolTip(this, string.Format("({0}, {1})", p.X, p.Y));
        }

        private void MenuSinX_Click(object sender, System.EventArgs e)
        {
            DrawGraphs(sender, x.Sin());
        }

        private void MenuSinhXover2_Click(object sender, System.EventArgs e)
        {
            DrawGraphs(sender, x.Over(2).Sinh());
        }

        private void Menu1overX_Click(object sender, System.EventArgs e)
        {
            DrawGraphs(sender, x.Power(-1));
        }

        private void MenuLnXsquaredMinus1_Click(object sender, System.EventArgs e)
        {
            DrawGraphs(sender, x.Squared().Minus(1).Log());
        }

        private Graph Graph = new Graph(new RectangleF(-9, -7, 18, 14), 16000);
        private Expression x = Differentiator.x;

        private void DrawGraph(Graphics g)
        {
            Graph.Draw(g, ClientRectangle);
        }

        private void DrawGraphs(object Sender, Expression y)
        {
            Text = ((ToolStripMenuItem)Sender).Text;
            Graph.Clear();
            Graph.AddSeries(y).PenColour = Color.Black;
            y = y.Differentiate();
            Graph.AddSeries(y).PenColour = Color.Red;
            y = y.Differentiate();
            Graph.AddSeries(y).PenColour = Color.Green;
            y = y.Differentiate();
            Graph.AddSeries(y).PenColour = Color.Blue;
            Invalidate();
        }
    }
}

The results are anisotropic, so as you resize the form, the x and y scalings vary independently. This affects the pen too - horizontal and vertical lines are generally of different thicknesses. This is only proof of concept code, which you might extend for your own purposes.

One last thing to note is the appearance of the sample graph at the head of this article. Its primary function y=ln(x²-1) is undefined for all values of x between -1 and +1, as is borne out by the shape of the two black curves, resembling a funnel, with no connection between the left and right sides. So far, so good. But if it's undefined in that region, why does it seem to have perfectly well defined first, second and third derivatives there, as indicated by the red, green and blue curves in the centre area? (Hint: look for a complex explanation!)

Next time I'll be returning to the business of simplifying the algebraic expressions returned by the differentiator.

No comments:

Post a Comment