Source code for pygal.interpolate

# -*- coding: utf-8 -*-
# This file is part of pygal
#
# A python svg graph plotting library
# Copyright © 2012-2016 Kozea
#
# This library is free software: you can redistribute it and/or modify it under
# the terms of the GNU Lesser General Public License as published by the Free
# Software Foundation, either version 3 of the License, or (at your option) any
# later version.
#
# This library is distributed in the hope that it will be useful, but WITHOUT
# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
# FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public License for more
# details.
#
# You should have received a copy of the GNU Lesser General Public License
# along with pygal. If not, see <http://www.gnu.org/licenses/>.
"""
Interpolation functions

These functions takes two lists of points x and y and
returns an iterator over the interpolation between all these points
with `precision` interpolated points between each of them

"""
from __future__ import division

from math import sin


[docs]def quadratic_interpolate(x, y, precision=250, **kwargs): """ Interpolate x, y using a quadratic algorithm https://en.wikipedia.org/wiki/Spline_(mathematics) """ n = len(x) - 1 delta_x = [x2 - x1 for x1, x2 in zip(x, x[1:])] delta_y = [y2 - y1 for y1, y2 in zip(y, y[1:])] slope = [delta_y[i] / delta_x[i] if delta_x[i] else 1 for i in range(n)] # Quadratic spline: a + bx + cx² a = y b = [0] * (n + 1) c = [0] * (n + 1) for i in range(1, n): b[i] = 2 * slope[i - 1] - b[i - 1] c = [(slope[i] - b[i]) / delta_x[i] if delta_x[i] else 0 for i in range(n)] for i in range(n + 1): yield x[i], a[i] if i == n or delta_x[i] == 0: continue for s in range(1, precision): X = s * delta_x[i] / precision X2 = X * X yield x[i] + X, a[i] + b[i] * X + c[i] * X2
[docs]def cubic_interpolate(x, y, precision=250, **kwargs): """ Interpolate x, y using a cubic algorithm https://en.wikipedia.org/wiki/Spline_interpolation """ n = len(x) - 1 # Spline equation is a + bx + cx² + dx³ # ie: Spline part i equation is a[i] + b[i]x + c[i]x² + d[i]x³ a = y b = [0] * (n + 1) c = [0] * (n + 1) d = [0] * (n + 1) m = [0] * (n + 1) z = [0] * (n + 1) h = [x2 - x1 for x1, x2 in zip(x, x[1:])] k = [a2 - a1 for a1, a2 in zip(a, a[1:])] g = [k[i] / h[i] if h[i] else 1 for i in range(n)] for i in range(1, n): j = i - 1 l = 1 / (2 * (x[i + 1] - x[j]) - h[j] * m[j]) if x[i + 1] - x[j] else 0 m[i] = h[i] * l z[i] = (3 * (g[i] - g[j]) - h[j] * z[j]) * l for j in reversed(range(n)): if h[j] == 0: continue c[j] = z[j] - (m[j] * c[j + 1]) b[j] = g[j] - (h[j] * (c[j + 1] + 2 * c[j])) / 3 d[j] = (c[j + 1] - c[j]) / (3 * h[j]) for i in range(n + 1): yield x[i], a[i] if i == n or h[i] == 0: continue for s in range(1, precision): X = s * h[i] / precision X2 = X * X X3 = X2 * X yield x[i] + X, a[i] + b[i] * X + c[i] * X2 + d[i] * X3
[docs]def hermite_interpolate( x, y, precision=250, type='cardinal', c=None, b=None, t=None ): """ Interpolate x, y using the hermite method. See https://en.wikipedia.org/wiki/Cubic_Hermite_spline This interpolation is configurable and contain 4 subtypes: * Catmull Rom * Finite Difference * Cardinal * Kochanek Bartels The cardinal subtype is customizable with a parameter: * c: tension (0, 1) This last type is also customizable using 3 parameters: * c: continuity (-1, 1) * b: bias (-1, 1) * t: tension (-1, 1) """ n = len(x) - 1 m = [1] * (n + 1) w = [1] * (n + 1) delta_x = [x2 - x1 for x1, x2 in zip(x, x[1:])] if type == 'catmull_rom': type = 'cardinal' c = 0 if type == 'finite_difference': for i in range(1, n): m[i] = w[i] = .5 * ((y[i + 1] - y[i]) / (x[i + 1] - x[i]) + (y[i] - y[i - 1]) / (x[i] - x[i - 1]) ) if x[i + 1] - x[i] and x[i] - x[i - 1] else 0 elif type == 'kochanek_bartels': c = c or 0 b = b or 0 t = t or 0 for i in range(1, n): m[i] = .5 * ((1 - t) * (1 + b) * (1 + c) * (y[i] - y[i - 1]) + (1 - t) * (1 - b) * (1 - c) * (y[i + 1] - y[i])) w[i] = .5 * ((1 - t) * (1 + b) * (1 - c) * (y[i] - y[i - 1]) + (1 - t) * (1 - b) * (1 + c) * (y[i + 1] - y[i])) if type == 'cardinal': c = c or 0 for i in range(1, n): m[i] = w[i] = (1 - c) * (y[i + 1] - y[i - 1]) / ( x[i + 1] - x[i - 1] ) if x[i + 1] - x[i - 1] else 0 def p(i, x_): t = (x_ - x[i]) / delta_x[i] t2 = t * t t3 = t2 * t h00 = 2 * t3 - 3 * t2 + 1 h10 = t3 - 2 * t2 + t h01 = -2 * t3 + 3 * t2 h11 = t3 - t2 return ( h00 * y[i] + h10 * m[i] * delta_x[i] + h01 * y[i + 1] + h11 * w[i + 1] * delta_x[i] ) for i in range(n + 1): yield x[i], y[i] if i == n or delta_x[i] == 0: continue for s in range(1, precision): X = x[i] + s * delta_x[i] / precision yield X, p(i, X)
[docs]def lagrange_interpolate(x, y, precision=250, **kwargs): """ Interpolate x, y using Lagrange polynomials https://en.wikipedia.org/wiki/Lagrange_polynomial """ n = len(x) - 1 delta_x = [x2 - x1 for x1, x2 in zip(x, x[1:])] for i in range(n + 1): yield x[i], y[i] if i == n or delta_x[i] == 0: continue for s in range(1, precision): X = x[i] + s * delta_x[i] / precision s = 0 for k in range(n + 1): p = 1 for m in range(n + 1): if m == k: continue if x[k] - x[m]: p *= (X - x[m]) / (x[k] - x[m]) s += y[k] * p yield X, s
[docs]def trigonometric_interpolate(x, y, precision=250, **kwargs): """ Interpolate x, y using trigonometric As per http://en.wikipedia.org/wiki/Trigonometric_interpolation """ n = len(x) - 1 delta_x = [x2 - x1 for x1, x2 in zip(x, x[1:])] for i in range(n + 1): yield x[i], y[i] if i == n or delta_x[i] == 0: continue for s in range(1, precision): X = x[i] + s * delta_x[i] / precision s = 0 for k in range(n + 1): p = 1 for m in range(n + 1): if m == k: continue if sin(0.5 * (x[k] - x[m])): p *= sin(0.5 * (X - x[m])) / sin(0.5 * (x[k] - x[m])) s += y[k] * p yield X, s
INTERPOLATIONS = { 'quadratic': quadratic_interpolate, 'cubic': cubic_interpolate, 'hermite': hermite_interpolate, 'lagrange': lagrange_interpolate, 'trigonometric': trigonometric_interpolate } if __name__ == '__main__': from pygal import XY points = [(.1, 7), (.3, -4), (.6, 10), (.9, 8), (1.4, 3), (1.7, 1)] xy = XY(show_dots=False) xy.add('normal', points) xy.add('quadratic', quadratic_interpolate(*zip(*points))) xy.add('cubic', cubic_interpolate(*zip(*points))) xy.add('lagrange', lagrange_interpolate(*zip(*points))) xy.add('trigonometric', trigonometric_interpolate(*zip(*points))) xy.add( 'hermite catmul_rom', hermite_interpolate(*zip(*points), type='catmul_rom') ) xy.add( 'hermite finite_difference', hermite_interpolate(*zip(*points), type='finite_difference') ) xy.add( 'hermite cardinal -.5', hermite_interpolate(*zip(*points), type='cardinal', c=-.5) ) xy.add( 'hermite cardinal .5', hermite_interpolate(*zip(*points), type='cardinal', c=.5) ) xy.add( 'hermite kochanek_bartels .5 .75 -.25', hermite_interpolate( *zip(*points), type='kochanek_bartels', c=.5, b=.75, t=-.25 ) ) xy.add( 'hermite kochanek_bartels .25 -.75 .5', hermite_interpolate( *zip(*points), type='kochanek_bartels', c=.25, b=-.75, t=.5 ) ) xy.render_in_browser()