贝塞尔曲面在犇嘶茶壶中是如何工作的?

5
我过早地发布了一个代码高尔夫挑战,用this dataset只有茶壶)绘制犹他州茶壶。(修订和发布的茶壶挑战) 但是当我深入研究数据以便举出一个小例子时,我意识到我不知道这些数据在做什么。我很好地理解了2D中的Bezier曲线,并实现了deCasteljau。但是对于3D,它是否也是同样的工作方式?

是的!是的!

数据包含每个包含16个顶点的补丁。它们的布局是否有标准顺序?如果它们对应于2D曲线,则四个角点实际上会接触表面,其余12个是控件,对吗?

是的!

我的“原始计划”是将形状简化为矩形,将它们投影到画布上,并用灰色填充形状,其灰度值由补丁法向量与光向量的点积大小计算得出。如果我将它简化到那个程度,它还会像茶壶吗?必须使用光线追踪才能实现可识别的图像吗?

这是主观的。 :-(

虽然看起来像几个问题,但它们都是这一个问题的方面:“请教教练,如何绘制贝塞尔曲面?我需要知道什么?”


这是我目前编写的代码。(使用了这个矩阵库:mat.ps

%!
%%BoundingBox: 115 243 493 487
%-115 -243 translate

(mat.ps)run  %include matrix library

/tok{ token pop exch pop }def
/s{(,){search{ tok 3 1 roll }{ tok exit }ifelse }loop }def
/f(teapot)(r)file def
/patch[ f token pop { [ f 100 string readline pop s ] } repeat ]def
/vert[ f token pop { [ f 100 string readline pop s ] } repeat ]def
%vert == patch == %test data input

/I3 3 ident def      % 3D identity matrix
/Cam [ 0 0 10 ] def  % world coords of camera center viewpoint
/Theta [ 0 0 0 ] def % y-rotation x-rotation z-rotation
/Eye [ 0 0 15 ] def  % eye relative to camera vp
/Rot I3 def          % initial rotation seq

/makerot {
    Theta 0 get roty        % pan
    Theta 1 get rotx matmul % tilt
    Theta 2 get rotz matmul % twist
} def

/proj {
    Cam {sub} vop  % translate to camera coords
    Rot matmul  % perform camera rotation
    0 get aload pop Eye aload pop  % extract dot x,y,z and eye xyz
    4 3 roll div exch neg          % perform perspective projection
    4 3 roll add 1 index mul
    4 1 roll 3 1 roll sub mul exch % (ez/dz)(dx-ex) (ez/dz)(dy-ey)
} def


/R 20 def
/H -3 def
/ang 0 def
{
    300 700 translate
    1 70 dup dup scale div setlinewidth

    /Cam [ ang sin R mul  H  ang cos R mul ] def % camera revolves around Y axis at height H, dist R
    /Theta [ ang  H R atan  0 ] def % rotate camera back to origin
    /Rot makerot def % squash rotation sequence into a matrix

    patch {
        % Four corners
        %[ exch dup 0 get exch dup 3 get exch dup 12 get exch 15 get ]

        % Boundary curves
        [ exch
            dup 8 get exch dup 4 get exch dup 0 get exch %curveto4

            dup 14 get exch dup 13 get exch dup 12 get exch %curveto3

            dup 7 get exch dup 11 get exch dup 15 get exch %curveto2

            dup 1 get exch dup 2 get exch dup 3 get exch %curveto1

            dup 0 get exch %moveto
        pop ]

        { 1 sub vert exch get proj } forall

        moveto
            curveto curveto curveto curveto

        stroke
        %flushpage flush (%lineedit)(r)file pop
    } forall
    pstack
    showpage

    %exit
    /ang ang 10 add def
} loop

这里是原始Newell茶壶数据集

这是我的糟糕得惊人的图像:
不知道顶点的组织方式。

更新:修复了错误。也许它们实际上是“正常”布局。至少选择正确的角落可以给出对称的形状:
找到了角落。

更新:边界曲线看起来更好了。
边界曲线


基本上,你是对的:它们有4个控制点(4个角落)。同样的原则也适用于3D贝塞尔曲面片:如果你需要它们“平滑”地连接,你必须对齐控制点。一旦你“平滑”地连接了足够多的曲面片,你就可以用合理的精度表示任何想要的表面。我很久以前编写了一个贝塞尔曲面片编辑器:随意查看https://github.com/ldematte/beppe。我们用它创建了可口可乐罐、船只、飞机;因此,一个优秀的犹他茶壶的近似应该不难。 - Lorenzo Dematté
1
哦,当我们编写它时,我们真的很感激这些Gamasutra教程:http://www.gamasutra.com/view/feature/131755/curved_surfaces_using_bzier_.php?print=1 - Lorenzo Dematté
@dema80 我想使用实际的原始数据。 - luser droog
我有几个想法 :) 你在原始数据集中有什么?三角形还是体积数据?在以前的一个项目中,我有体积数据;从那里,我使用Marching Cubes提取了三角形,并且效果相当不错..确实,它被识别为一个茶壶 :) - Lorenzo Dematté
1
Foley的Intro相当简洁,但Foley&vanDam的Fundamentals则包含了相当多的内容,其中包括生成网格的伪代码。 - luser droog
显示剩余5条评论
2个回答

6
一个双三次贝塞尔曲面路径是一个4x4的3D点数组。是的,四个角落触碰到了曲面;行是贝塞尔曲线,列也是贝塞尔曲线。但deCasteljau算法是基于计算两个点之间的中值,并且在3D和2D中同样有意义。
完成上面代码的下一步是将补丁细分以覆盖更小的部分。然后上述简单的边界曲线提取变成适合的多边形网格。
首先,将补丁平展开来,直接插入顶点数据,而不是使用单独的缓存。该代码通过迭代补丁,在顶点数组中查找点并构建新的补丁数组,然后重新定义相同名称。
 /patch[ patch{ [exch { 1 sub vert exch get }forall ] }forall ]def

我们需要使用deCasteljau算法来拆分贝塞尔曲线。 vop来自于矩阵库,并且对向量的相应元素应用二元操作,产生一个新的向量作为结果。

/median { % [x0 y0 z0] [x1 y1 z1]
    {add 2 div} vop % [ (x0+x1)/2 (y0+y1)/2 (z0+z1)/2 ]
} def   

/decasteljau { % [P0] P1 P2 P3  .  P0 P1' P2' P3'  P3' P4' P5' P3
    {p3 p2 p1 p0}{exch def}forall
    /p01 p0 p1 median def
    /p12 p1 p2 median def
    /p23 p2 p3 median def
    /p012 p01 p12 median def
    /p123 p12 p23 median def
    /p0123 p012 p123 median def
    p0 p01 p012 p0123  % first half-curve
    p0123 p123 p23 p3  % second half-curve
} def   

接着需要对每个补丁的行进行堆栈操作,并将结果组装成两个新的补丁。

/splitrows { % [b0 .. b15]  .  [c0 .. c15] [d0 .. d15] 
    aload pop % b0 .. b15
    4 {                          % on each of 4 rows
        16 12 roll decasteljau   % roll the first 4 to the top
        8 4 roll                 % exch left and right halves (probably unnecessary)
        20 4 roll                % roll new curve to below the patch (pushing earlier ones lower)
    } repeat
    16 array astore              % pack the left patch
    17 1 roll 16 array astore    % roll, pack the right patch
} def

一种丑陋的实用程序让我们可以重复使用行代码用于列。堆栈注释在编写此过程时是必要的,因此它们可能阅读它时也是必要的。n j roll滚动n个元素(向左),j次; ==从1开始计数的第n个元素上面的前j个元素。因此,n不断减少,选择放置元素的位置,而j则选择将哪个元素放在那里(并将其他所有元素一起拖动)。如果应用了bind,则此过程比基于字典的过程要快得多。

% [ 0  1  2  3
%   4  5  6  7
%   8  9 10 11
%  12 13 14 15 ]
/xpose {
    aload pop % 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15
    15 12 roll % 0 4 5 6 7 8 9 10 11 12 13 14 15 1 2 3
    14 11 roll % 0 4 8 9 10 11 12 13 14 15 1 2 3 5 6 7
    13 10 roll % 0 4 8 12 13 14 15 1 2 3 5 6 7 9 10 11

    12 9 roll % 0 4 8 12 1 2 3 5 6 7 9 10 11 13 14 15
    11 9 roll % 0 4 8 12 1 5 6 7 9 10 11 13 14 15 2 3
    10 8 roll % 0 4 8 12 1 5 9 10 11 13 14 15 2 3 6 7
    9 7 roll % 0 4 8 12 1 5 9 13 14 15 2 3 6 7 10 11

    8 6 roll % 0 4 8 12 1 5 9 13 2 3 6 7 10 11 14 15
    7 6 roll % 0 4 8 12 1 5 9 13 2 6 7 10 11 14 15 3
    6 5 roll % 0 4 8 12 1 5 9 13 2 6 10 11 14 15 3 7
    5 4 roll % 0 4 8 12 1 5 9 13 2 6 10 14 15 3 7 11
    4 3 roll % 0 4 8 12 1 5 9 13 2 6 10 14 3 7 11 15
    16 array astore
} def
% [ 0 4  8 12
%   1 5  9 13
%   2 6 10 14
%   3 7 11 15 ]

/splitcols {
    xpose
    splitrows
    xpose
} def

然后将这些函数应用于补丁数据。每次重新定义补丁。
/patch[ patch{ splitrows }forall ]def
/patch[ patch{ splitrows }forall ]def
/patch[ patch{ splitcols }forall ]def
/patch[ patch{ splitcols }forall ]def

这提供了处理较小片段的能力。
细分线框图

添加一个可见性测试。

/visible { % patch  .  patch boolean 
    dup % p p
    dup 3 get exch dup 0 get exch 12 get % p p3 p0 p12
    1 index {sub} vop % p p3 p0 v0->12
    3 1 roll {sub} vop % p v0->12 v0->3
    cross /normal exch def
    dup
    [ exch dup 0 get exch dup 3 get exch dup 12 get exch 15 get ]
    { Cam {sub} vop normal dot 0 ge } forall
    %add add add 4 div 0 lt
    or or or
} def

制作中
应用可见性测试

更新:测试倒置。
已修正的可见性测试

更新:测试无效!从图像可以看出,底部未朝外,当然,背面消隐并不能阻止手柄穿过壶体。这需要进行隐藏表面移除。由于Postscript不支持Z缓冲区,我想只好使用二叉空间分割了。所以我要重新学习一下。

更新:添加模型->世界变换以将物体竖立。

/Model -90 rotx def % model->world transform

/proj {
    Model matmul 0 get             % perform model->world transform
    Cam {sub} vop                  % translate to camera coords
    ...

制作这个。

正面朝上

目前为止的完整程序。 (使用矩阵库:mat.ps.) 在ghostscript中,按住[enter]可以查看旋转动画。

    %!
    %%BoundingBox: 109 246 492 487
    %-109 -246 translate

    (mat.ps)run  %include matrix library
    (det.ps)run  %supplementary determinant function

    /tok{ token pop exch pop }def
    /s{(,){search{ tok 3 1 roll }{ tok exit }ifelse }loop }def
    /f(teapot)(r)file def
    /patch[ f token pop { [ f 100 string readline pop s ] } repeat ]def
    /vert[ f token pop { [ f 100 string readline pop s ] } repeat ]def
    /patch[ patch{ [exch { 1 sub vert exch get }forall ] }forall ]def
    %vert == patch == %test data input flush quit

    /I3 3 ident def      % 3D identity matrix
    /Cam [ 0 0 10 ] def  % world coords of camera center viewpoint
    /Theta [ 0 0 0 ] def % y-rotation x-rotation z-rotation
    /Eye [ 0 0 15 ] def  % eye relative to camera vp
    /Rot I3 def          % initial rotation seq

    /Model -90 rotx def % model->world transform

    /makerot {
            Theta 0 get roty        % pan
            Theta 1 get rotx matmul % tilt
            Theta 2 get rotz matmul % twist
    } def

    /proj {
            Model matmul 0 get             % perform model->world transform
            Cam {sub} vop  % translate to camera coords
            Rot matmul  % perform camera rotation
            0 get aload pop Eye aload pop  % extract dot x,y,z and eye xyz
            4 3 roll div exch neg          % perform perspective projection
            4 3 roll add 1 index mul
            4 1 roll 3 1 roll sub mul exch % (ez/dz)(dx-ex) (ez/dz)(dy-ey)
    } def

    /median { % [x0 y0 z0] [x1 y1 z1]
            {add 2 div} vop % [ (x0+x1)/2 (y0+y1)/2 (z0+z1)/2 ]
    } def

    /decasteljau { % [P0] P1 P2 P3  .  P0 P1' P2' P3'  P3' P4' P5' P3
            {p3 p2 p1 p0}{exch def}forall
            /p01 p0 p1 median def
            /p12 p1 p2 median def
            /p23 p2 p3 median def
            /p012 p01 p12 median def
            /p123 p12 p23 median def
            /p0123 p012 p123 median def
            p0 p01 p012 p0123
            p0123 p123 p23 p3
    } def

    /splitrows { % [b0 .. b15]  .  [c0 .. c15] [d0 .. d15]
            aload pop % b0 .. b15
            4 {
                    16 12 roll decasteljau
                    %8 4 roll
                    20 4 roll
            } repeat
            16 array astore
            17 1 roll 16 array astore
    } def

    /xpose {
            aload pop % 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15
            15 12 roll % 0 4 5 6 7 8 9 10 11 12 13 14 15 1 2 3
            14 11 roll % 0 4 8 9 10 11 12 13 14 15 1 2 3 5 6 7
            13 10 roll % 0 4 8 12 13 14 15 1 2 3 5 6 7 9 10 11

            12 9 roll % 0 4 8 12 1 2 3 5 6 7 9 10 11 13 14 15
            11 9 roll % 0 4 8 12 1 5 6 7 9 10 11 13 14 15 2 3
            10 8 roll % 0 4 8 12 1 5 9 10 11 13 14 15 2 3 6 7
            9 7 roll % 0 4 8 12 1 5 9 13 14 15 2 3 6 7 10 11

            8 6 roll % 0 4 8 12 1 5 9 13 2 3 6 7 10 11 14 15
            7 6 roll % 0 4 8 12 1 5 9 13 2 6 7 10 11 14 15 3
            6 5 roll % 0 4 8 12 1 5 9 13 2 6 10 11 14 15 3 7
            5 4 roll % 0 4 8 12 1 5 9 13 2 6 10 14 15 3 7 11
            4 3 roll % 0 4 8 12 1 5 9 13 2 6 10 14 3 7 11 14
            16 array astore
    } def

    /splitcols {
            xpose
            splitrows
            xpose exch xpose
    } def

    /patch[ patch{ splitrows }forall ]def
    /patch[ patch{ splitrows }forall ]def
    /patch[ patch{ splitrows }forall ]def
    /patch[ patch{ splitrows }forall ]def
    /patch[ patch{ splitcols }forall ]def
    /patch[ patch{ splitcols }forall ]def

    /color {normal light dot 1 add 4 div
%1 exch sub
setgray} def

    /visible { % patch  .  patch boolean
            dup % p p
            dup 3 get exch dup 0 get exch 12 get % p p3 p0 p12
            1 index {sub} vop % p p3 p0 v0->12
            3 1 roll {sub} vop % p v0->12 v0->3
            cross /normal exch def
            dup
            [ exch dup 0 get exch dup 3 get exch dup 12 get exch 15 get ]
            { Cam {sub} vop normal dot 0 ge } forall
            %add add add 4 div 0 lt
            or or or
    } def

    /drawpatch {

            % Four corners
            %[ exch dup 0 get exch dup 3 get exch dup 12 get exch 15 get ]

            visible {

                    [ exch
                            % control rows
                            %dup 4 get exch dup 5 get exch dup 6 get exch dup 7 get exch
                            %dup 11 get exch dup 10 get exch dup 9 get exch dup 8 get exch

                            % control columns
                            %dup 1 get exch dup 5 get exch dup 9 get exch dup 13 get exch
                            %dup 14 get exch dup 10 get exch dup 6 get exch dup 2 get exch

                            % Boundary curves
                            dup 8 get exch dup 4 get exch dup 0 get exch %curveto4
                            dup 14 get exch dup 13 get exch dup 12 get exch %curveto3
                            dup 7 get exch dup 11 get exch dup 15 get exch %curveto2
                            dup 1 get exch dup 2 get exch dup 3 get exch %curveto1
                            dup 0 get exch %moveto
                    pop ]

                    { proj } forall

                    moveto curveto curveto curveto curveto
                    %moveto lineto lineto lineto lineto lineto lineto lineto closepath
                    %moveto lineto lineto lineto lineto lineto lineto lineto closepath

                    stroke
                    %flushpage flush (%lineedit)(r)file pop

            }{
                    pop
            }ifelse

    } def

    /R 20 def
    /H -3 def
    /ang 10 def
    {
            300 700 translate
            1 70 dup dup scale div setlinewidth

            % camera revolves around Y axis at height H, dist R
            /Cam [ ang sin R mul  H  ang cos R mul ] def
            /Theta [ ang  H R atan  0 ] def   % rotate camera back to origin
            /Rot makerot def        % squash rotation sequence into a matrix

            patch {
                    drawpatch
            } forall
            pstack
            showpage

            %exit
            /ang ang 10 add def
    } loop

一些BSP链接:http://web.archive.org/web/20120709060316/http://maven.smith.edu/~mcharley/bsp/createbsptree.html http://web.archive.org/web/20111115163809/http://www.devmaster.net/articles/bsp-trees http://web.archive.org/web/20040625025943/http://occs.cs.oberlin.edu/faculty/jdonalds/357/lecture25.html - luser droog
math.SE 上请求有关BSP生成的帮助。 - luser droog
注意:/splitcols 中的第二个 xpose 应该是 xpose exch xpose,以便将转置应用于两个新补丁中的每一个(如果保留顺序很重要,则需要另一个 exch)。 - luser droog
上述错误可能是最后一张图中出现“鸡笼”效果的原因。通过使用斜角连接而不是斜接连接来解决问题。 - luser droog

0

根据math.StackExchange上的帮助,我被引导到了一个辅助目标,即为矩阵库添加一个计算行列式的函数。

所以,这段代码通过了一些笨拙的初始测试,但我必须承认它非常丑陋:

GS>[[1 0][0 1]] det
GS<1>=
1
GS>[[0 1][1 0]] det =
-1
GS>(mat.ps) run
GS>3 ident
GS<1>det =
1
GS>[[1 2 3][4 5 6][7 8 9]] det =
0
GS>

更新。稍微可读一点。

更新。使用点乘和叉乘更加可读。再次感谢MvG。

(mat.ps) run % use dot and cross from matrix library

/elem { % M i j
    3 1 roll get exch get % M_i_j
} def

/det {
    dup length 1 index 0 get length ne { /det cvx /typecheck signalerror } if
    1 dict begin /M exch def
        M length 2 eq {
            M 0 0 elem
            M 1 1 elem mul
            M 0 1 elem
            M 1 0 elem mul sub
        }{
            M length 3 eq {
                M aload pop cross dot
            }{ /det cvx /rangecheck signalerror } ifelse
        } ifelse
    end
} def

1
存储矩阵的名称以在该函数期间使用,并编写一个提取矩阵成员的函数如何?然后,您可以编写0 0 memb而不是dup 0 get 0 get。您可能会更系统地记录下来,例如硬编码Sarrus法则。或者,您可以将det(A,B,C)写为(A×B)∙C,即叉积后跟点积。由于det是对称的,因此您可以将列或行视为向量A、B、C。 - MvG

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接