Metapost 绘制 B-Spline 递归公式给出奇怪的结果

Metapost 绘制 B-Spline 递归公式给出奇怪的结果

我目前正在尝试使用元帖子制作通用 B 样条基绘图仪。

基函数的定义是维基百科在属性下。

在此处输入图片描述

ab 样条曲线的节点应具有模式 0,其次数等于样条曲线的顺序,一段时间内值不断增加,最后一个值重复以匹配顺序

例如对于 3 阶样条曲线:[0,0,0,1,2,3,4,5,6,6,6] 会生成一个节点向量。

为此,我尝试编写该精确公式并打印我制作的结矢量,以确保我做的事情正确:

在此处输入图片描述

上面是 2 阶样条函数之一,它应该看起来像一个线性上升然后下降的函数(金字塔),但它看起来像 1 阶基函数。我尝试更改顺序,但我的所有函数都返回相同的模式。

例如 3 阶函数: 在此处输入图片描述

我的代码如下;

\documentclass[border=10cm]{standalone}

\usepackage{luamplib}
\mplibnumbersystem{double}
\usepackage[margin=0.5cm]{geometry}

\begin{document}
{\centering
\begin{mplibcode}
u:=1cm;

numeric knots[];

vardef make_knots(expr order, num)= 

    p_num = num - 1;
    v = 0;
    for i=0 upto order-1:
        knots[i] := v;
    endfor;
    for i=0 upto p_num - order:
        v := v+1;
        knots[i + order] := v;
    endfor;
    for i=0 upto order-1:
        knots[i + p_num + 1] := v;
    endfor;

    for i=0 upto order + p_num - 1:
        tmp := round(knots[i]);
        label.top(textext("\huge$"& decimal(tmp) &"$"), (2*i*u,6*u));
    endfor;
enddef;

vardef calculate_basis(expr t, i, order)= 
    numeric ret;
    if order=1:
        if (t >= knots[i]) and (t < knots[i+1]):
            ret = 1;
        else:
            ret = 0;
        fi;
    else:
        ret = 
        ((t-knots[i]) / (knots[i] + order - knots[i])) * calculate_basis(t, i, order-1) 
        + 
        ((knots[i + order + 1] - t) / (knots[i + order + 1] - knots[i+1])) * calculate_basis(t, i+1, order-1);
    fi;
    ret
enddef;

color darkred, darkyellow, darkgreen, lightblue, brown, pink, orange; 
darkred := (0.8,0.0,0.0);
darkyellow := (1.0,0.8,0.0);
darkgreen := (0.0,0.6,0.0);
lightblue := (0.0,0.8,1.0);
brown := (0.5, 0.1, 0.1);
pink := (1, 0.0, 0.8);
orange := (1, 0.4, 0.0);
%ignore first parameter while denbugging
vardef plot_basis(expr j, order, color)= 
    res := 100;

    for i=0 upto res-1:
        fraction := (i/res) * 6; %multiply by biggest knot value
        valS := calculate_basis(fraction, 1, 3);

        save pointS;
        pair pointS;
        pointS := ((i) * 2 / res, valS) * u * 5; %scale plot to make it visible
        draw pointS withpen pencircle scaled 5bp withcolor color;

    endfor;
enddef;

% Start figure
beginfig(0);

color colors[];
colors[0] = darkred;
colors[1] = darkyellow;
colors[2] = darkgreen;
colors[3] = lightblue;
colors[4] = pink;
colors[5] = brown;
colors[6] = orange;

make_knots(2, 8);
plot_basis(2, 2, colors[0]);

endfig;
\end{mplibcode}
\par}
\end{document}

答案1

我无法完全解释原因,但我确实有一个修复程序(我认为),并且可能还修复了几个拼写错误。在您的递归中(参见注释):

vardef calculate_basis(expr t, i, order)=
    numeric ret;
    if order=0: % I think this should be zero to match your screenshot
        if (t >= knots[i]) and (t < knots[i+1]):
            ret = 1;
        else:
            ret = 0;
        fi;
    else:
        ret = 
        ((t-knots[i]) / (knots[i + order] - knots[i])) * calculate_basis(t, i, order-1) % fixed typo in denominator of first fraction
        + 
        ((knots[i + order + 1] - t) / (knots[i + order + 1] - knots[i+1])) * calculate_basis(t, i+1, order-1);
    fi;
    ret
enddef;

在 else 情况下,我认为这ret是以某种方式获取递归中评估的最终阶数为 0 的情况的值。我对 B 样条函数一无所知,也不太了解这里的分组、递归调用和方程求解如何相互作用。递归的工作原理让我感到困惑,因为它看起来扩展为类似

ret=stuff1*(ret=0)+stuff2*(ret=1)

这对我来说毫无意义。希望有人能解释原因并这样做。

无论如何,我都会这样编写递归:

vardef calculate_basis(expr t, i, order)=
    if order=0:
        if (t >= knots[i]) and (t < knots[i+1]):
                    1
        else:
                    0
        fi
    else:
        ((t-knots[i]) / (knots[i + order] - knots[i])) * calculate_basis(t, i, order-1) % typo
        + 
        ((knots[i + order + 1] - t) / (knots[i + order + 1] - knots[i+1])) * calculate_basis(t, i+1, order-1)
    fi
enddef;

并做出最后的改变后

    valS := calculate_basis(fraction, 3, 1); % changed to first order

以下似乎产生了更接近您想要的结果(我不确定绘图基础循环中剩余的固定值来自哪里)。无论哪种方式,希望这能给你一个开始:

\documentclass[border=10cm]{standalone}

\usepackage{luamplib}
\mplibnumbersystem{double}
\mplibtextextlabel{enable}
\usepackage[margin=0.5cm]{geometry}

\begin{document}
{\centering
\begin{mplibcode}
u:=1cm;

numeric knots[];

vardef make_knots(expr order, num)= 
    p_num = num - 1;
    v = 0;
    for i=0 upto order-1:
        knots[i] := v;
    endfor;
    for i=0 upto p_num - order:
        v := v+1;
        knots[i + order] := v;
    endfor;
    for i=0 upto order-1:
        knots[i + p_num + 1] := v;
    endfor;

    for i=0 upto order + p_num - 1:
        label.top("\huge$"& decimal(round(knots[i])) &"$", (2*i*u,6*u));
    endfor;
enddef;

vardef calculate_basis(expr t, i, order)=
    if order=0:
        if (t >= knots[i]) and (t < knots[i+1]):
                    1
        else:
                    0
        fi
    else:
        ((t-knots[i]) / (knots[i + order] - knots[i])) * calculate_basis(t, i, order-1) % typo
        + 
        ((knots[i + order + 1] - t) / (knots[i + order + 1] - knots[i+1])) * calculate_basis(t, i+1, order-1)
    fi
enddef;

color darkred, darkyellow, darkgreen, lightblue, brown, pink, orange; 
darkred := (0.8,0.0,0.0);
darkyellow := (1.0,0.8,0.0);
darkgreen := (0.0,0.6,0.0);
lightblue := (0.0,0.8,1.0);
brown := (0.5, 0.1, 0.1);
pink := (1, 0.0, 0.8);
orange := (1, 0.4, 0.0);

vardef plot_basis(expr order, color)= 
    res := 100;
        pair pointS;
    for i=0 upto res-1:
        fraction := (i/res) * 6; 
        valS := calculate_basis(fraction, 3, 1); % changed to first order
        pointS := ((i) * 2 / res, valS) * u * 5; 
        draw pointS withpen pencircle scaled 5bp withcolor color;
    endfor;
enddef;

% Start figure
beginfig(0);

color colors[];
colors[0] = darkred;
colors[1] = darkyellow;
colors[2] = darkgreen;
colors[3] = lightblue;
colors[4] = pink;
colors[5] = brown;
colors[6] = orange;
make_knots(2, 8);
plot_basis(1, colors[0]);
endfig;
\end{mplibcode}
\par}
\end{document}

在此处输入图片描述

相关内容