【Python】グラフから値を抽出し，関数にする【SVG】

はじめに

研究を行っていると，先行研究と比較するときなど，グラフから実数を得たい場面は意外と多いです．
以前，画像から半自動的に抽出する記事を書いたのですが↓今一つ使い勝手が良くなかったため（色分けされた実線グラフのみ対応），今回は別のアプローチを試してみます．
【Python】（線）グラフの画像から値を抽出する

方法

グラフをなぞったSVG画像からベジェ曲線を計算し，線形補間で関数化します．

SVGとはベクター画像のフォーマットのひとつで，XMLベースで記述されています．
有償ソフトでは Adobe Illustrator，無償ソフトでは Inkscape などで作成することができます．
ベクター画像で図を作成できると，論文やスライドなどを見栄え良くできるので，スキルとして習得しておくと何かと役に立つと思います．

今回は線をなぞるだけですので，難しくはありません．

グラフ画像の準備

Inkscape を例に説明します．

適当ですが，このようなグラフがあったとします．

Inkscape で開いて，ペンツールを選択します．

欲しい線をなるべくきれいになぞります（複数のグラフでもOKです，1つのグラフは1本の繋がった線でなぞってください）．
今回は緑をなぞりました．

どこでも良いので2点，なるべく確実な座標間に直線を引きます．
この２点の座標は後でスケール合わせに使用するため覚えておいてください．
今回は(0, 3)と(8, 0)の間に直線を引きました．

画像を好きな名前のSVGファイルとして保存します．

これで準備は完了です．

Pythonで抽出する

関数たち

2/18 xyの入れ替え，logスケール，点の座標のみ取得に対応しました．

import numpy as np
import re
from bs4 import BeautifulSoup

# ベジェ曲線の計算
def Bezier(p, t):
    *_, n = p.shape
    J = np.vectorize(lambda i, t: (n/np.r_[1:i+1]-1).prod() * t**i * (1-t)**(n-1-i), signature='(),(m)->(m)')
    return p @ J(np.r_[:n], t)

# スケール合わせ
def _scale(l, ax, pt1, pt2):
    e1, e3 = np.eye(3), np.eye(3)
    ax_ = _parse_svg(ax)[0]
    e1[:2, -1] = -ax_[:, 0]
    e2 = np.diag(np.r_[(np.array(pt2) - np.array(pt1)) / np.diff(ax_, 1)[:, 0], [1]])
    e3[:2, -1] = pt1
    e = e3@e2@e1    
    return [(e @ np.pad(i, ((0, 1), (0, 0)), constant_values=1))[:2] for i in l]

# SVGからpathを取得
def _get_paths(svg_path):
    with open(svg_path) as f:
        dat = f.read()
    svg = BeautifulSoup(dat, "xml")
    paths = []
    for path in svg("path"):
        v = np.array([0.0, 0.0])
        d = []
        for i in re.findall("([a-zA-Z])|(-?[\d\.]+)", path["d"]):
            c = i[0] or c
            if i[1]:
                v = float(i[1])
                d.append((c, v))
        paths.append([(i[0], np.array([i[1], -j[1]])) for i, j in zip(d[::2], d[1::2])])
    return paths

# pathから制御点の座標をセグメント毎に抽出
def _parse_svg(d):
    o = d[0][1]
    i = 1
    l = []
    while i < len(d):
        c = d[i][0]
        cc = {"c": 3, "q": 2, "l": 1, "m": 1}[c.lower()]
        l.append(np.r_[[o]+[d[i+j][1]+o*int(c.islower()) for j in range(cc)]].T)
        i += cc
        o = l[-1][:,-1]
    return l

# pathから制御点の座標を抽出
def _parse_svgp(d):
    o = d[0][1]
    i = 1
    l = [o.tolist()]
    while i < len(d):
        c = d[i][0]
        cc = {"c": 3, "q": 2, "l": 1, "m": 1}[c.lower()]
        l.append((d[i+cc-1][1]+o*int(c.islower())).tolist())
        i += cc
        o = l[-1]
    return np.array(l).T[None]

# 線形補間で関数化
def _funcfy(beziers, resolution, swapxy):
    t = np.linspace(0, 1, resolution)
    x, y = np.hstack([Bezier(b, t) for b in beziers])
    if swapxy:
        x, y = y, x
    idx = x.argsort()
    return lambda t: np.interp(t, x[idx], y[idx])

# メインの関数
def svg2curvefunc(svg_path, pt1=None, pt2=None, resolution=256, swapxy=False, logx=False, logy=False, loglog=False):
    d = _get_paths(svg_path)
    ax, *cvs = sorted(d, key=lambda i: abs(len(i)-2))
    curves = []
    if loglog:
        logx = logy = True
    for curve in cvs:
        l = _parse_svg(curve)
        if pt1 is not None and pt2 is not None:
            pt1 = [np.log(pt1[0]) if logx else pt1[0], np.log(pt1[1]) if logy else pt1[1]]
            pt2 = [np.log(pt2[0]) if logx else pt2[0], np.log(pt2[1]) if logy else pt2[1]]
            l = _scale(l, ax, pt1, pt2)
            l = [np.c_[np.exp(i[0]) if logx else i[0], np.exp(i[1]) if logy else i[1]].T for i in l]
        curves.append(_funcfy(l, resolution, swapxy))
    return curves

# メインの関数(点座標を取得)
def svg2points(svg_path, pt1=None, pt2=None, resolution=256, swapxy=False, logx=False, logy=False, loglog=False):
    d = _get_paths(svg_path)
    ax, *cvs = sorted(d, key=lambda i: abs(len(i)-2))
    curves = []
    if loglog:
        logx = logy = True
    for curve in cvs:
        l = _parse_svgp(curve)
        if pt1 is not None and pt2 is not None:
            pt1 = [np.log(pt1[0]) if logx else pt1[0], np.log(pt1[1]) if logy else pt1[1]]
            pt2 = [np.log(pt2[0]) if logx else pt2[0], np.log(pt2[1]) if logy else pt2[1]]
            l = _scale(l, ax, pt1, pt2)
            l = [np.c_[np.exp(i[0]) if logx else i[0], np.exp(i[1]) if logy else i[1]].T for i in l]
        curves.append(np.hstack(l))
    return curves

上の関数群を実行して，

# svg2curvefunc(SVGファイルのパス, 先程の座標１, 先程の座標２)
svg2curvefunc('../test_curve.svg', [0, 3], [8, 0])

とするだけで関数化されたものがリストとして得られます．

今回は１つだけなぞったので，0番目の関数を取り出してプロットしてみましょう．

import matplotlib.pyplot as plt
import numpy as np

t = np.linspace(-1, 9, 1000)
f = svg2curvefunc('../test_curve.svg', [0, 3], [8, 0])[0]
plt.plot(t, f(t))
plt.grid(alpha=.5)
plt.show()

元画像

きれいに抽出できました．

なお，範囲外の値は端点の値で固定になります．

スケールがおかしい場合は，直線を引いた座標が逆になっている可能性がありますので，座標1と座標2を入れ替えて実行してみてください．

# svg2curvefunc(SVGファイルのパス, 先程の座標２, 先程の座標１)
svg2curvefunc('../test_curve.svg', [8, 0], [0, 3])

おわりに

SVGのパース部分はあまり知識のない人間が自力で書いたものなので，場合によっては上手く座標が計算できない可能性があります．ご了承下さい．

曲線をなぞった際の各点間が大きく離れている場合は，ベジェ曲線の計算の解像度が不十分になる可能性が考えられます．
引数resolutionを用意していますので（デフォルト256点）適宜大きい値を設定してください．