0，引入
前续：
Q79：怎么用三角形网格（Triangle Mesh）细分曲面
http://blog.csdn.net/libing_zeng/article/details/60600404
Q80：平坦着色（Flat Shading）和平滑着色（Smooth Shading）——“Q79：怎么用三角形网格（Triangle Mesh）细分曲面”（补充）
http://blog.csdn.net/libing_zeng/article/details/60760296
【修正】问题五十五：怎么用ray tracing画Utah teapot (bicubic bezier patches)
http://blog.csdn.net/libing_zeng/article/details/69258616
这一章节，我们以Utah Teapot（Bicubic Bezier Patches）为例来进行说明。
1，理论分析
我们知道Utah Teapot是由若干个Bicubic Bezier Patch拼接而成（我们这里用到的数据是32-patches、306-vertices版本）。
关于Bezier曲面的介绍，参考：“问题五十四：怎么用ray tracing画参数方程表示的曲面（2）—— bezier surface”
http://blog.csdn.net/libing_zeng/article/details/54565264
Bezier曲面的本质是：参数方程表示的曲面。这个“参数方程”就是两条Bezier曲线的的张量积。
在求得球面的参数方程之后，用三角形网格细分该曲面的方式可以参考：
“Q79：怎么用三角形网格（Triangle Mesh）细分曲面”
http://blog.csdn.net/libing_zeng/article/details/60600404
2，C++代码实现
2.1 读取teapot数据文件
本人是将文件读取函数移植到Vector3D.cpp中
// ---------------------------------------------------------- get_teapot_data
// read utah teapot with 32 patches and 306 vertices from external file
bool get_teapot_data(int (&patches)[32][16], float (&vertices)[306][3]) {
ifstream infile( ".\\teaset\\teapot");
char str[100];
char *token;
int item_num;
int patch_num = 0;
int vertex_num = 0;
int flag = 0;
int tokens_i[16];
float tokens_f[3];
while (infile >> str) {
item_num = 0;
token = strtok(str, ",");
if ((flag == 0) || (flag == 1)) {
sscanf(token, "%d", &(tokens_i[item_num]));
if (tokens_i[0] == 32) {
flag = 1;
}
else if (tokens_i[0] == 306) {
flag = 2;
}
else {
patches[patch_num][item_num] = tokens_i[item_num];
item_num ++;
}
token = strtok(NULL, ",");
}
while (token != NULL) {
if (flag == 1) {
sscanf(token, "%d", &(tokens_i[item_num]));
patches[patch_num][item_num] = tokens_i[item_num];
}
if (flag == 2) {
sscanf(token, "%f", &(tokens_f[item_num]));
vertices[vertex_num][item_num] = tokens_f[item_num];
}
item_num ++;
token = strtok(NULL, ",");
}
if ((flag == 1) && (tokens_i[0] != 32)) {patch_num ++;}
if (flag == 2) {
if ((tokens_i[0] == 306)) {tokens_i[0] = 0;}
else {vertex_num ++;}
}
}
infile.close();
return true;
}
// ---------------------------------------------------------- matrix_4_4_multiply_4_4
void matrix_4_4_multiply_4_4(const float matrix1[4][4], const float matrix2[4][4], float (&result)[4][4]) {
for (int k=0; k<4; k++) {
for (int i=0; i<4; i++) {
result[i][k] = 0.0;
for (int j=0; j<4; j++) {
result[i][k] = result[i][k] + matrix1[i][j]*matrix2[j][k];
}
}
}
}
2.2 矩阵运算
计算过程会用到“矩阵运算”（没有使用“矩阵类”，仍旧使用二维数组表示矩阵）。矩阵运算的相关函数也是在Vector3D.cpp中实现：
// ---------------------------------------------------------- matrix_4_4_multiply_4_4
void matrix_4_4_multiply_4_4(const float matrix1[4][4], const float matrix2[4][4], float (&result)[4][4]) {
for (int k=0; k<4; k++) {
for (int i=0; i<4; i++) {
result[i][k] = 0.0;
for (int j=0; j<4; j++) {
result[i][k] = result[i][k] + matrix1[i][j]*matrix2[j][k];
}
}
}
}
// ---------------------------------------------------------- matrix_1_4_multiply_4_4
void matrix_1_4_multiply_4_4(const float matrix1[1][4], const float matrix2[4][4], float (&result)[1][4]) {
for (int k=0; k<4; k++) {
result[0][k] = 0.0;
for (int j=0; j<4; j++) {
result[0][k] = result[0][k] + matrix1[0][j]*matrix2[j][k];
}
}
}
// ---------------------------------------------------------- matrix_1_4_multiply_4_1
void matrix_1_4_multiply_4_1(const float matrix1[1][4], const float matrix2[4][1], float &result) {
result = 0.0;
for (int j=0; j<4; j++) {
result = result + matrix1[0][j]*matrix2[j][0];
}
}
2.3 细分Bezier曲面
像之前一样，细分函数在Grid.cpp中实现。
// ------------------------------------ tessellate_flat_bezier_patches ---------------------------------------
// tesselate a cubic bezier patch into flat triangles that are stored directly in the grid
void
Grid::tessellate_flat_bezier_patches(const int horizontal_steps, const int vertical_steps,
float vertices[306][3], int patches[32][16], const int patches_num) {
Vector3D patches_vertices[patches_num][16];
float matrix_c_x[4][4], matrix_c_y[4][4], matrix_c_z[4][4];
float points_x[4][4], points_y[4][4], points_z[4][4];
float points_x_t[4][4], points_y_t[4][4], points_z_t[4][4];
float uuu[1][4], vvv[4][1];
float uuu1[1][4], vvv1[4][1];
float xxx_t[1][4], yyy_t[1][4], zzz_t[1][4];
float xxx, yyy, zzz;
float xxx1, yyy1, zzz1;
float xxx2, yyy2, zzz2;
float xxx3, yyy3, zzz3;
float matrix_t[4][4] = {{ 1, 0, 0, 0},
{-3, 3, 0, 0},
{ 3, -6, 3, 0},
{-1, 3, -3, 1}};
float matrix[4][4] = {{1, -3, 3, -1},
{0, 3, -6, 3},
{0, 0, 3, -3},
{0, 0, 0, 1}};
int ip1, ip2, ipv;
for (int i=0; i<patches_num; i++) {
for (int j=0; j<16; j++) {
ip1 = int(j/4);
ip2 = int(j%4);
ipv = patches[i][j] - 1;
points_x[ip1][ip2] = vertices[ipv][0];
points_y[ip1][ip2] = vertices[ipv][1];
points_z[ip1][ip2] = vertices[ipv][2];
}
matrix_4_4_multiply_4_4(matrix_t, points_x, points_x_t);
matrix_4_4_multiply_4_4(points_x_t, matrix, matrix_c_x);
matrix_4_4_multiply_4_4(matrix_t, points_y, points_y_t);
matrix_4_4_multiply_4_4(points_y_t, matrix, matrix_c_y);
matrix_4_4_multiply_4_4(matrix_t, points_z, points_z_t);
matrix_4_4_multiply_4_4(points_z_t, matrix, matrix_c_z);
for (int k = 0; k <= vertical_steps - 1; k++) {
for (int m = 0; m <= horizontal_steps - 1; m++) {
uuu[0][0] = 1.0;
uuu[0][1] = float(m) / float(horizontal_steps);
uuu[0][2] = uuu[0][1] * uuu[0][1];
uuu[0][3] = uuu[0][1] * uuu[0][2];
vvv[0][0] = 1.0;
vvv[1][0] = float(k) / float(vertical_steps);
vvv[2][0] = vvv[0][1] * vvv[0][1];
vvv[3][0] = vvv[0][1] * vvv[0][2];
uuu1[0][0] = 1.0;
uuu1[0][1] = float(m+1) / float(horizontal_steps);
uuu1[0][2] = uuu1[0][1] * uuu1[0][1];
uuu1[0][3] = uuu1[0][1] * uuu1[0][2];
vvv1[0][0] = 1.0;
vvv1[1][0] = float(k+1) / float(vertical_steps);
vvv1[2][0] = vvv1[0][1] * vvv1[0][1];
vvv1[3][0] = vvv1[0][1] * vvv1[0][2];
matrix_1_4_multiply_4_4(uuu, matrix_c_x, xxx_t);
matrix_1_4_multiply_4_1(xxx_t, vvv, xxx);
matrix_1_4_multiply_4_4(uuu, matrix_c_y, yyy_t);
matrix_1_4_multiply_4_1(yyy_t, vvv, yyy);
matrix_1_4_multiply_4_4(uuu, matrix_c_z, zzz_t);
matrix_1_4_multiply_4_1(zzz_t, vvv, zzz);
matrix_1_4_multiply_4_4(uuu, matrix_c_x, xxx_t);
matrix_1_4_multiply_4_1(xxx_t, vvv1, xxx1);
matrix_1_4_multiply_4_4(uuu, matrix_c_y, yyy_t);
matrix_1_4_multiply_4_1(yyy_t, vvv1, yyy1);
matrix_1_4_multiply_4_4(uuu, matrix_c_z, zzz_t);
matrix_1_4_multiply_4_1(zzz_t, vvv1, zzz1);
matrix_1_4_multiply_4_4(uuu1, matrix_c_x, xxx_t);
matrix_1_4_multiply_4_1(xxx_t, vvv, xxx2);
matrix_1_4_multiply_4_4(uuu1, matrix_c_y, yyy_t);
matrix_1_4_multiply_4_1(yyy_t, vvv, yyy2);
matrix_1_4_multiply_4_4(uuu1, matrix_c_z, zzz_t);
matrix_1_4_multiply_4_1(zzz_t, vvv, zzz2);
matrix_1_4_multiply_4_4(uuu1, matrix_c_x, xxx_t);
matrix_1_4_multiply_4_1(xxx_t, vvv1, xxx3);
matrix_1_4_multiply_4_4(uuu1, matrix_c_y, yyy_t);
matrix_1_4_multiply_4_1(yyy_t, vvv1, yyy3);
matrix_1_4_multiply_4_4(uuu1, matrix_c_z, zzz_t);
matrix_1_4_multiply_4_1(zzz_t, vvv1, zzz3);
// define the first triangle
Point3D v0(xxx2, yyy2, zzz2);
Point3D v1(xxx3, yyy3, zzz3);
Point3D v2( xxx, yyy, zzz);
Triangle* triangle_ptr1 = new Triangle(v0, v1, v2);
objects.push_back(triangle_ptr1);
// define the second triangle
v0 = Point3D(xxx1, yyy1, zzz1);
v1 = Point3D( xxx, yyy, zzz);
v2 = Point3D(xxx3, yyy3, zzz3);
Triangle* triangle_ptr2 = new Triangle(v0, v1, v2);
objects.push_back(triangle_ptr2);
}
}
}
}
如上代码的截图分析如下：
3，测试图形
3.1 测试代码
3.2 输出图形
teapot_32_0,0,100_rx-90,ts(0,-2,0)_light,directional(10,-20,-20)_20,20
teapot_32_0,0,100_rx-90,ts(0,-2,0)_light,directional(10,-20,-20)_40,40_190s_ns1
teapot_32_0,0,100_rx-90,ts(0,-2,0)_light,directional(10,-20,-20)_40,40_207s_ns16
teapot_32_0,0,100_rx-90,ts(0,-2,0)_light,directional(10,-20,-20)_40,40_176s_ns16_marble,p0.1
teapot_32_0,0,100_rx-90,ts(0,-2,0)_light,directional(10,-20,-20)_40,40_183s_ns16_marble,p4
teapot_32_0,0,100_rx-90,ts(0,-2,0)_light,directional(10,-20,-20)_40,40_226s_ns16_sandstone,p0.1
teapot_32_0,0,100_rx-90,ts(0,-2,0)_light,directional(10,-20,-20)_40,40_234s_ns16_wood
143
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
