Q100：怎么用三角形网格细分回旋体（rotational sweeping / revolution）

来源：互联网发布：软件开发思路编辑：程序博客网时间：2024/05/02 04:28

0，引入

我们在“问题六十：怎么用ray tracing画回旋体（rotational sweeping / revolution）”中已经学习了这类曲面的画法：
http://blog.csdn.net/libing_zeng/article/details/54609610

当时，是各种联立方程、微分等等，然后解一元六次方程。不得不说“三角形网格”确实是个好东西，接下来，我们看看怎么用三角形网格来细分这类曲面。

1，理论介绍

“回旋体”是指某“基本曲线”围绕某轴线旋转一周得到的曲面。我们这里考虑的“基本曲线”是由多段b-spline曲线拼接而成的曲线段；“轴线”是Y轴。（关于b-spline曲线，可以参考：
http://blog.csdn.net/libing_zeng/article/details/54565614）

接下来，我们看看用三角形网格是怎么细分回旋曲面的。（若对下方内容有疑问，请参考前面两个链接对应的文章）

这里写图片描述

图1中的曲线段是由3段“三次b-spline曲线段”组成，每一段“三次b-spline曲线段”是由四个控制点确定的。
另外，注意图一中的坐标系：曲线有自己的坐标系（一般称为局部坐标系）uv，曲线上任意点P的坐标为（u，v）。图中所示状态是曲线所在平面uov刚好旋转到xoy平面，所以此时：u坐标系和x坐标系重合；v坐标系和y坐标系重合。
对于，每一段“三次b-spline曲线段”上的任意点P（u，v），根据“三次b-spline曲线”的参数方程：
这里写图片描述

“三次b-spline曲线”的参数方程中引入了参数s，我们还需引入另一个参数phi（如图2所示）。参数phi确定曲线段在旋转过程中的位置。

曲线任意点P（u，v）在空间中的坐标（x，y，z）怎么确定呢？
根据图1，我们知道y=v；
根据图2，我们知道x=u*cos(phi)、z=u*sin(phi)。
即：
这里写图片描述

计算过程：
1，已知的控制点坐标、参数s可以确定点P在uv坐标系中的坐标。（即“式子100.1”）
2，根据“旋转”的性质（如图2）得到空间坐标（x，y，z）和局部坐标（u，v）之间的关系，从而确定空间坐标（x，y，z）

对应三角形网格细分过程：
1，细分参数s、phi（即给定这两个参数）；
2，根据参数确定细分点的空间坐标。

2，C++代码实现

像之前一样，细分函数在Grid.cpp中实现。

// ------------------------------------  tessellate_flat_rotational_sweeping  ---------------------------------------// tesselate rotational sweeping surface into flat triangles that are stored directly in the gridvoidGrid::tessellate_flat_rotational_sweeping(const int horizontal_steps, const int vertical_steps,                                          Point2D* ctrl_points, const int ctrl_points_num, const bool reverse_normal) {    int curve_num = ctrl_points_num - 3;    float   points_u[4][1],   points_v[4][1];    float matrix_c_u[4][1], matrix_c_v[4][1];    float       sss0[1][4],             phi0;    float       sss1[1][4],             phi1;    float            uuus0,            vvvs0;    float            uuus1,            vvvs1;    float             xxx0,             yyy0,             zzz0;    float             xxx1,             yyy1,             zzz1;    float             xxx2,             yyy2,             zzz2;    float             xxx3,             yyy3,             zzz3;    float matrix_t_6[4][4] = {{ 1,  4,  1, 0},                              {-3,  0,  3, 0},                              { 3, -6,  3, 0},                              {-1,  3, -3, 1}};    float matrix_t[4][4];    for (int i=0; i<4; i++) {        for (int j=0; j<4; j++) {            matrix_t[i][j] = matrix_t_6[i][j] / 6.0;        }    }    for (int i=0; i<curve_num; i++) {        int k = 0;        for (int j=i; j< (i + 4); j++) {            points_u[k][0] = ctrl_points[j].x;            points_v[k][0] = ctrl_points[j].y;            k++;        }        matrix_4_4_multiply_4_1(matrix_t, points_u, matrix_c_u);        matrix_4_4_multiply_4_1(matrix_t, points_v, matrix_c_v);        for (int k = 0; k <= vertical_steps - 1; k++) {            for (int m = 0; m <= horizontal_steps - 1; m++) {                sss0[0][0] = 1.0;                sss0[0][1] = float(k) / float(vertical_steps);                sss0[0][2] = sss0[0][1] * sss0[0][1];                sss0[0][3] = sss0[0][1] * sss0[0][2];                sss1[0][0] = 1.0;                sss1[0][1] = float(k+1) / float(vertical_steps);                sss1[0][2] = sss1[0][1] * sss1[0][1];                sss1[0][3] = sss1[0][1] * sss1[0][2];                phi0 = float(m) * TWO_PI / float(horizontal_steps);                phi1 = float(m+1) * TWO_PI / float(horizontal_steps);                matrix_1_4_multiply_4_1(sss0, matrix_c_u, uuus0);                matrix_1_4_multiply_4_1(sss0, matrix_c_v, vvvs0);                matrix_1_4_multiply_4_1(sss1, matrix_c_u, uuus1);                matrix_1_4_multiply_4_1(sss1, matrix_c_v, vvvs1);                xxx0 = uuus0 * cos(phi0);                yyy0 = vvvs0;                zzz0 = uuus0 * sin(phi0);                xxx1 = uuus0 * cos(phi1);                yyy1 = vvvs0;                zzz1 = uuus0 * sin(phi1);                xxx2 = uuus1 * cos(phi0);                yyy2 = vvvs1;                zzz2 = uuus1 * sin(phi0);                xxx3 = uuus1 * cos(phi1);                yyy3 = vvvs1;                zzz3 = uuus1 * sin(phi1);                // define the first triangle                Point3D v0(xxx2, yyy2, zzz2);                Point3D v1(xxx3, yyy3, zzz3);                Point3D v2(xxx0, yyy0, zzz0);                Triangle* triangle_ptr1 = new Triangle(v0, v2, v1, reverse_normal);                objects.push_back(triangle_ptr1);                // define the second triangle                v0 = Point3D(xxx1, yyy1, zzz1);                v1 = Point3D(xxx0, yyy0, zzz0);                v2 = Point3D(xxx3, yyy3, zzz3);                Triangle* triangle_ptr2 = new Triangle(v0, v2, v1, reverse_normal);                objects.push_back(triangle_ptr2);            }        }    }}