Bezier Patches 如何在犹他州茶壶中工作? [英] How do Bezier Patches work in the Utah Teapot?
问题描述
我过早地发布了一个代码高尔夫挑战，使用 .)在 ghostscript 中，您可以通过按住 [enter] 来查看动画旋转.
<代码> %!%% 边界框:109 246 492 487%-109 -246 翻译(mat.ps) 运行 %include 矩阵库(det.ps) 运行 % 补充行列式函数/tok{ token pop exch pop }def/s{(,){search{ tok 3 1 roll }{ tok exit }ifelse }loop }def/f(茶壶)(r)文件定义/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 数据输入刷新退出/I3 3 ident def % 3D 单位矩阵/Cam [ 0 0 10 ] def % 相机中心视点的世界坐标/Theta [ 0 0 0 ] def % y-rotation x-rotation z-rotation/Eye [ 0 0 15 ] def % 眼睛相对于相机 vp/Rot I3 def % 初始旋转序列/Model -90 rotx def % 模型->世界变换/makerot {Theta 0 get roty % panTheta 1 获得 rotx matmul % 倾斜Theta 2 得到 rotz matmul % 扭曲定义/项目{模型 matmul 0 get % perform model->world transformCam {sub} vop % 转换为相机坐标Rot matmul % 执行相机旋转0 get aload pop Eye aload pop % 提取点 x,y,z 和 eye xyz4 3 roll div exch neg % 执行透视投影4 3 卷加 1 索引 mul4 1 roll 3 1 roll sub mul exch % (ez/dz)(dx-ex) (ez/dz)(dy-ey)定义/median { % [x0 y0 z0] [x1 y1 z1]{add 2 div} vop % [ (x0+x1)/2 (y0+y1)/2 (z0+z1)/2 ]定义/decasteljau { % [P0] P1 P2 P3 .P0 P1' P2' P3' P3' P4' P5' P3{p3 p2 p1 p0}{exch def}forall/p01 p0 p1 中值定义/p12 p1 p2 中值定义/p23 p2 p3 中值定义/p012 p01 p12 中值定义/p123 p12 p23 中值定义/p0123 p012 p123 中值定义p0 p01 p012 p0123p0123 p123 p23 p3定义/splitrows { % [b0 .. b15] .[c0 .. c15] [d0 .. d15]aload pop % b0 .. b154 {16 12 卷 decasteljau%8 4 卷20 4 卷} 重复16 阵列存储17 1 卷 16 阵列存储定义/xpose {aload pop % 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 1515 12 卷% 0 4 5 6 7 8 9 10 11 12 13 14 15 1 2 314 11 卷% 0 4 8 9 10 11 12 13 14 15 1 2 3 5 6 713 10 卷% 0 4 8 12 13 14 15 1 2 3 5 6 7 9 10 1112 9 卷% 0 4 8 12 1 2 3 5 6 7 9 10 11 13 14 1511 9 卷% 0 4 8 12 1 5 6 7 9 10 11 13 14 15 2 310 8 卷% 0 4 8 12 1 5 9 10 11 13 14 15 2 3 6 79 7 滚动 % 0 4 8 12 1 5 9 13 14 15 2 3 6 7 10 118 6 卷% 0 4 8 12 1 5 9 13 2 3 6 7 10 11 14 157 6 卷% 0 4 8 12 1 5 9 13 2 6 7 10 11 14 15 36 5 卷% 0 4 8 12 1 5 9 13 2 6 10 11 14 15 3 75 4 卷% 0 4 8 12 1 5 9 13 2 6 10 14 15 3 7 114 3 滚动 % 0 4 8 12 1 5 9 13 2 6 10 14 3 7 11 1416 阵列存储定义/splitcols {摆姿势分流xpose 交易所 xpose定义/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 {普通光点1加4格%1 交换子setgray} 定义/可见{％补丁.补丁布尔值复制 % p pdup 3 得到 exch dup 0 得到 exch 12 得到 % p p3 p0 p121 索引 {sub} vop % p p3 p0 v0->123 1 卷 {sub} vop % p v0->12 v0->3交叉/正常交换定义重复[ exch dup 0 得到 exch dup 3 得到 exch dup 12 得到 exch 15 ]{ Cam {sub} vop normal dot 0 ge } forall%add 添加添加 4 div 0 lt或 或 或定义/画布{% 四个角%[ exch dup 0 得到 exch dup 3 得到 exch dup 12 得到 exch 15 ]可见的 {[ 交易所% 控制行%dup 4 得到 exch dup 5 得到 exch dup 6 得到 exch dup 7 得到 exch%dup 11 得到 exch dup 10 得到 exch dup 9 得到 exch dup 8 得到 exch% 控制列%dup 1 得到 exch dup 5 得到 exch dup 9 得到 exch dup 13 得到 exch%dup 14 得到 exch dup 10 得到 exch dup 6 得到 exch dup 2 得到 exch% 边界曲线dup 8 得到 exch dup 4 得到 exch dup 0 得到 exch %curveto4dup 14 得到 exch dup 13 得到 exch dup 12 得到 exch %curveto3dup 7 得到 exch dup 11 得到 exch dup 15 得到 exch %curveto2dup 1 得到 exch dup 2 得到 exch dup 3 得到 exch %curveto1dup 0 get exch %moveto流行音乐 ]{ 项目} forall移动到曲线到曲线到曲线到曲线到%moveto lineto lineto lineto lineto lineto lineto lineto closepath%moveto lineto lineto lineto lineto lineto lineto lineto closepath中风%flushpage 刷新 (%lineedit)(r) 文件弹出}{流行音乐}如果别的定义/R 20 防御/H -3 定义/ang 10 定义{300 700 翻译1 70 dup dup scale div setlinewidth% 相机围绕 Y 轴在高度 H, dist R 处旋转/Cam [ang sin R mul Hang cos R mul ] def/Theta [ ang H R atan 0 ] def % 将相机旋转回原点/Rot makerot def % 将旋转序列挤压成矩阵修补 {抽签} 对所有人堆栈显示页面％出口/ang ang 10 添加定义} 环形I prematurely posted a code golf challenge to draw the Utah Teapot using this dataset (just the teapot). (Revised and Posted teapot challenge) But when I looked deeper at the data in order to whip up a little example, I realized I have no idea what's going on with that data. I have a good understanding of Bezier curves in 2D, implemented deCasteljau. But for 3D does it work the same?
Yes! It does!
The data contains patches containing 16 vertices each. Is there a standard ordering for how these are laid out? And if they correspond to the 2D curves, then the four corner points actually touch the surface and the remaining 12 are controls, right?
Yes!
My "original plan" was to simplify the shape to rectangles, project them to the canvas, and draw filled shapes in a grey computed by the magnitude of the dot product of the patch normal to a light vector. If I simplify it that far, will it even look like a teapot? Does one have to use raytracing to achieve a recognizable image?
That's subjective. :-(
While this may look like several questions, but they are all aspects of this one: "Please, kindly Guru, school me on some Bezier Patches? What do I need to know to draw the teapot?"
Here's the code I've written so far. (uses this matrix library: 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
Here's the original Newell Teapot dataset.
And here's my spectacularly bad image:
Update: bugfix. Maybe they are laid out 'normally' after all. Selecting the correct corners at least gives a symmetrical shape:
Update: boundary curves looks better.
A Bi-Cubic Bezier surface patch is a 4x4 array of 3D points. Yes, the four corners touch the surface; and the rows are Bezier curves, and columns are also Bezier curves. But the deCasteljau algorithm is based on calculating the median between two points, and is equally meaningful in 3D as in 2D.
The next step in completing the above code is subdividing the patches to cover smaller portions. Then the simple boundary curve extraction above becomes a suitable polygon mesh.
Start by flattening the patches, inserting the vertex data directly instead of using a separate cache. This code iterates through the patches, looking up points in the vertex array and constructs a new array of patches which is then redefined with the same name.
/patch[ patch{ [exch { 1 sub vert exch get }forall ] }forall ]def
Then we need the deCasteljau algorithm to split Bezier curves. vop comes from the matrix library and applies a binary operation upon corresponding elements of a vector and produces a new vector as the result.
/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
Then some stack manipulation to apply to each row of a patch and assemble the results into 2 new patches.
/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
An ugly utility lets us reuse the row code for columns. The stack comments were necessary to write this procedure, so they're probably necessary to read it. n j roll rolls n elements (to the left), j times; == the top j elements above n-th element (counting from 1). So n steady decreases, selecting where to put the element, and j selects which element to put there (dragging everything else with it). If bind were applied, this procedure would be substantially faster than a dictionary-based procedure.
% [ 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
Then apply these functions to the patch data. Again, redefining patch each time.
/patch[ patch{ splitrows }forall ]def
/patch[ patch{ splitrows }forall ]def
/patch[ patch{ splitcols }forall ]def
/patch[ patch{ splitcols }forall ]def
This gives the ability to deal with smaller fragments.
Add a visibility test.
/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
Producing
Update: test was backwards.
Update: Test is Useless! You can see from the image that the bottom piece is not oriented outward, and of course, backface-culling doesn't prevent the handle from showing through the pot. This calls for hidden surface removal. And since Postscript doesn't have support for a Z-buffer, I guess it'll have to be a Binary Space Partition. So it's back to the books for me.
Update: Add a Model->World transform to turn the thing upright.
/Model -90 rotx def % model->world transform
/proj {
Model matmul 0 get % perform model->world transform
Cam {sub} vop % translate to camera coords
...
Producing this.
Complete program so far. (uses matrix library:mat.ps.) In ghostscript, you can view an animated rotation by holding [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
这篇关于Bezier Patches 如何在犹他州茶壶中工作?的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！
