如何在犹他茶壶贝塞尔工作？(How do Bezier Patches work in the Ut

我过早地公布了代码高尔夫挑战赛吸引了犹他茶壶使用此数据集（只是茶壶）。（修订和发布茶壶挑战），但是当我在数据以掀起一个小例子看起来更深，我意识到我不知道发生了什么事情与该数据。我有贝塞尔曲线的2D有很好的理解，实现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

这里是原来的纽维尔茶壶数据集。

这是我的壮观坏形象：

更新：修正错误。也许他们毕竟是奠定了“正常的”。选择正确的角落至少给出了一个对称的形状：

更新：边界曲线看起来更好。

Answer 1:

双三次Bezier曲面补丁是3D点的一个4×4阵列。是的，四角触摸表面上; 和行是贝塞尔曲线，和列也贝塞尔曲线。但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

然后一些堆栈操作适用于补丁的每一行和组装结果存入2个新补丁。

/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

一个丑陋的工具可以让我们重用列行代码。堆栈的意见是必要的，写这个程序，所以他们很可能需要阅读。 nj roll辊n个元素（向左）中，j倍; ==第n个元素上面的顶部元件Ĵ（从1开始计数）。因此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

生产

更新：测试是倒退。

更新： 测试是没用！ 可以从该底部件不向外定向的图像看，当然，背面-剔除不防止把手通过锅表示。这就要求隐面消除。而且，由于后记不具有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

Answer 2:

基于对援助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