60.1 概述
回旋体，大概是长这个样子：
回旋体是指曲线（称为“基本曲线”）围绕y轴转一圈得到的图形。
（基本曲线是由多段b-spline曲线段连接而成）
这里先强调一下：
上图中有两个坐标系：熟悉的空间三维xyz坐标系和曲线的局部二维uv坐标系。
对于回旋体，有如下特点：
上面这两个坐标系有这种关系：
对于回旋体上任意一点P(x,y,z)都对应一条曲线，则对应一个曲线的局部二维uv坐标系，则对应这个局部坐标系中一个坐标(u,v)。
60.2 求根
60.3 求法向量
60.4 看C++代码实现
----------------------------------------------sweeping_rotational.h ------------------------------------------
sweeping_rotational.h
/*这个文件和sweeping_translational.h差别不大，可以对比参照来理解*/
#ifndef SWEEPING_ROTATIONAL_H
#define SWEEPING_ROTATIONAL_H
#include <hitable.h>
#include <log.h>
class sweeping_rotational : public hitable
{
public:
sweeping_rotational() {}
sweeping_rotational(vec3 cen, vec3 *cp, material *m){
center = cen;
ma = m;
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;
}
}
float points_u[6], points_v[6], b_points_u[4][1], b_points_v[4][1], b_result_u[4][1], b_result_v[4][1];
int b_num = 0;
float max_u, min_v, max_v;
for (int i=0; i<6; i++) {
points_u[i] = cp[i].x();
points_v[i] = cp[i].y();
}
max_u = points_u[0];
min_v = points_v[0];
max_v = points_v[0];
for (int i=0; i<6; i++) {
if (max_u < points_u[i]) {
max_u = points_u[i];
}
if (min_v > points_v[i]) {
min_v = points_v[i];
}
if (max_v < points_v[i]) {
max_v = points_v[i];
}
}
sweeping_bl = vec3(-max_u+center.x(), min_v+center.y(), -max_u+center.z());
sweeping_bh = vec3( max_u+center.x(), max_v+center.y(), max_u+center.z());
for (int i=0; i<3; i=i+1) {
b_points_u[0][0] = points_u[i];
b_points_u[1][0] = points_u[i+1];
b_points_u[2][0] = points_u[i+2];
b_points_u[3][0] = points_u[i+3];
b_points_v[0][0] = points_v[i];
b_points_v[1][0] = points_v[i+1];
b_points_v[2][0] = points_v[i+2];
b_points_v[3][0] = points_v[i+3];
matrix_4_4_multiply_4_1(matrix_t, b_points_u, b_result_u);
matrix_4_4_multiply_4_1(matrix_t, b_points_v, b_result_v);
for (int j=0; j<4; j++) {
matrix_c_u[j][b_num] = ((fabs(b_result_u[j][0])<1e-6)? (0.0):(b_result_u[j][0]));
matrix_c_v[j][b_num] = ((fabs(b_result_v[j][0])<1e-6)? (0.0):(b_result_v[j][0]));
}
b_num ++;
}
}
virtual bool hit(const ray& r, float tmin, float tmax, hit_record& rec) const;
float matrix_c_u[4][3], matrix_c_v[4][3];
vec3 sweeping_bl, sweeping_bh, center;
material *ma;
};
#endif // SWEEPING_ROTATIONAL_H
----------------------------------------------sweeping_rotational.cpp ------------------------------------------
sweeping_rotational.cpp
#include "sweeping_rotational.h"
#include <iostream>
#include <limits>
#include "float.h"
bool ray_hit_box_sr(const ray& r, const vec3& vertex_l, const vec3& vertex_h, float& t_near, float& t_far) {
/*罪恶般地又抄了一遍这个函数，光线撞击包围回旋体外面的最小盒子*/
t_near = (numeric_limits<float>::min)();
t_far = (numeric_limits<float>::max)();
vec3 direction = r.direction();
vec3 origin = r.origin();
vec3 bl = vertex_l;
vec3 bh = vertex_h;
float array1[6];
if(direction.x() == 0) {
if((origin.x() < bl.x()) || (origin.x() > bh.x())) {
return false;
}
array1[0] = (numeric_limits<float>::min)();
array1[1] = (numeric_limits<float>::max)();
}
if(direction.y() == 0) {
if((origin.y() < bl.y()) || (origin.y() > bh.y())) {
return false;
}
array1[2] = (numeric_limits<float>::min)();
array1[3] = (numeric_limits<float>::max)();
}
if(direction.z() == 0) {
if((origin.z() < bl.z()) || (origin.z() > bh.z())) {
return false;
}
array1[4] = (numeric_limits<float>::min)();
array1[5] = (numeric_limits<float>::max)();
}
if((direction.x() != 0) && (direction.y() != 0) && (direction.z() != 0)) {
array1[0] = (bl.x()-origin.x())/direction.x();
array1[1] = (bh.x()-origin.x())/direction.x();
array1[2] = (bl.y()-origin.y())/direction.y();
array1[3] = (bh.y()-origin.y())/direction.y();
array1[4] = (bl.z()-origin.z())/direction.z();
array1[5] = (bh.z()-origin.z())/direction.z();
}
for (int i=0; i<6; i=i+2){
if(array1[i] > array1[i+1]) {
float t = array1[i];
array1[i] = array1[i+1];
array1[i+1] = t;
}
if(array1[i] >= t_near) {t_near = array1[i];}
if(array1[i+1] <= t_far) {t_far = array1[i+1];}
if(t_near > t_far) {
return false;
}
if(t_far < 0) {
return false;
}
}
if (t_near != t_near) {
t_near = t_near * 1;
}
return true;
}
bool sweeping_rotational::hit(const ray& r, float t_min, float t_max, hit_record& rec) const {
#if SWEEPING_ROTATIONAL_LOG == 1
std::cout << "-------------sweeping_rotational::hit---1-------------" << endl;
#endif // SWEEPING_ROTATIONAL_LOG
float t_near, t_far;
if (ray_hit_box_sr(r, sweeping_bl, sweeping_bh, t_near, t_far)) {
float xx0 = r.origin().x() - center.x();
float yy0 = r.origin().y() - center.y();
float zz0 = r.origin().z() - center.z();
float xxd = r.direction().x();
float yyd = r.direction().y();
float zzd = r.direction().z();
if (fabs(xxd) < 1e-6) {
xxd = 1e-6;
}
if (fabs(yyd) < 1e-6) {
yyd = 1e-6;
}
if (fabs(zzd) < 1e-6) {
zzd = 1e-6;
}
/*对应“式子60.3”*/
float aaa = xxd*xxd + zzd*zzd;
float bbb = 2 * ((xx0*xxd+zz0*zzd)*yyd - aaa*yy0);
float ccc = (xx0*yyd-xxd*yy0)*(xx0*yyd-xxd*yy0) + (zz0*yyd-zzd*yy0)*(zz0*yyd-zzd*yy0);
float ddd = yyd*yyd;
float ss6[7], roots[7], roots_t[19][3], yyv, temp0, temp1, temp2;
float tol = 1e-6;
int num_roots_t = 0;
for (int i=0; i<3; i++) {
/*“基本曲线”是有三段b-spline曲线段组成，分别对每一段求交。对应“式子60.5”*/
ss6[0] = aaa*matrix_c_v[3][i]*matrix_c_v[3][i] - ddd*matrix_c_u[3][i]*matrix_c_u[3][i];
ss6[1] = aaa*2*matrix_c_v[2][i]*matrix_c_v[3][i] - ddd*2*matrix_c_u[2][i]*matrix_c_u[3][i];
ss6[2] = aaa*(matrix_c_v[2][i]*matrix_c_v[2][i]+2*matrix_c_v[1][i]*matrix_c_v[3][i]) -
ddd*(matrix_c_u[2][i]*matrix_c_u[2][i]+2*matrix_c_u[1][i]*matrix_c_u[3][i]);
ss6[3] = aaa*2*(matrix_c_v[0][i]*matrix_c_v[3][i]+matrix_c_v[1][i]*matrix_c_v[2][i]) + bbb*matrix_c_v[3][i] -
ddd*2*(matrix_c_u[0][i]*matrix_c_u[3][i]+matrix_c_u[1][i]*matrix_c_u[2][i]);
ss6[4] = aaa*(matrix_c_v[1][i]*matrix_c_v[1][i]+2*matrix_c_v[0][i]*matrix_c_v[2][i]) + bbb*matrix_c_v[2][i] -
ddd*(matrix_c_u[1][i]*matrix_c_u[1][i]+2*matrix_c_u[0][i]*matrix_c_u[2][i]);
ss6[5] = aaa*2*matrix_c_v[0][i]*matrix_c_v[1][i] + bbb*matrix_c_v[1][i] - ddd*2*matrix_c_u[0][i]*matrix_c_u[1][i];
ss6[6] = aaa*matrix_c_v[0][i]*matrix_c_v[0][i] + bbb*matrix_c_v[0][i] - ddd*matrix_c_u[0][i]*matrix_c_u[0][i] + ccc;
roots_equation_6th(ss6, 0, 1, tol, roots);
for (int j=1; j<(int(roots[0])+1); j++) {
yyv = matrix_c_v[0][i]+matrix_c_v[1][i]*roots[j]+matrix_c_v[2][i]*roots[j]*roots[j]+matrix_c_v[3][i]*roots[j]*roots[j]*roots[j];
roots_t[num_roots_t+1][0] = (yyv-yy0)/yyd;
roots_t[num_roots_t+1][1] = i;
roots_t[num_roots_t+1][2] = roots[j];
num_roots_t ++;
}
}
roots_t[0][0] = float(num_roots_t);
if (roots_t[0][0] != 0.0) {
/*对所有的交点t从小到大排序，所以最小的值跑到了最前面*/
for (int i=1; i<int(roots_t[0][0]); i++) {
for (int j=i+1; j<int(roots_t[0][0])+1; j++) {
if (roots_t[i][0] > roots_t[j][0]) {
temp0 = roots_t[i][0];
roots_t[i][0] = roots_t[j][0];
roots_t[j][0] = temp0;
temp1 = roots_t[i][1];
roots_t[i][1] = roots_t[j][1];
roots_t[j][1] = temp1;
temp2 = roots_t[i][2];
roots_t[i][2] = roots_t[j][2];
roots_t[j][2] = temp2;
}
}
}
if ((roots_t[1][0] < t_max) && (roots_t[1][0] > t_min)) {
rec.t = roots_t[1][0];
rec.p = r.point_at_parameter(rec.t);
float nx = rec.p.x() - center.x();
float nz = rec.p.z() - center.z();
int num = int(roots_t[1][1]);
float sss = roots_t[1][2];
float nu = -(matrix_c_v[1][num]+2.0*matrix_c_v[2][num]*sss+3.0*matrix_c_v[3][num]*sss*sss);
float nv = (matrix_c_u[1][num]+2.0*matrix_c_u[2][num]*sss+3.0*matrix_c_u[3][num]*sss*sss);
float ny = sqrt(nx*nx+nz*nz)*nv/nu - center.y();
/*求法向量，对应“式子60.6”*/
rec.normal = unit_vector(vec3(nx, ny, nz));
if(dot(r.direction(), rec.normal) > 0) {
rec.normal = - rec.normal;
}
rec.mat_ptr = ma;
rec.u = -1.0;
rec.v = -1.0;
return true;
}
}
}
return false;
}
----------------------------------------------main.cpp ------------------------------------------
main.cpp
vec3 ctrl_points[6] = {vec3(-1.0, 5.0, 0.0), vec3( 2.0, 4.0, 0.0),
vec3( 2.0, 1.0, 0.0), vec3(-0.5, 1.0, 0.0),
vec3( 1.5, -3.0, 0.0), vec3( 3.0, 0.0, 0.0)};
/*“基本曲线的控制点，6个点对应3条曲线段”*/
hitable *list[1];
list[0] = new sweeping_rotational(vec3( 0.0, 0.0, 0.0), ctrl_points, new lambertian(vec3(1.0, 0.0, 0.0)));
hitable *world = new hitable_list(list,1);
vec3 lookfrom(0, 10, 20);
vec3 lookat(0, 1, 0);
float dist_to_focus = (lookfrom - lookat).length();
float aperture = 0.0;
camera cam(lookfrom, lookat, vec3(0,1,0), 20, float(nx)/float(ny), aperture, 0.7*dist_to_focus);
输出图形如下：
vec3 lookfrom(0, 10, 20);
vec3 lookfrom(10, 10, 10);
来张组合图：
vec3 ctrl_points[6] = {vec3(-1.0, 5.0, 0.0), vec3( 2.0, 4.0, 0.0),
vec3( 2.0, 1.0, 0.0), vec3(-0.5, 1.0, 0.0),
vec3( 1.5, -3.0, 0.0), vec3( 3.0, 0.0, 0.0)};
hitable *list[4];
list[0] = new sweeping_rotational(vec3( 0.0, 0, 0.0), ctrl_points, new dielectric(4.5));
list[1] = new sweeping_rotational(vec3(-5.0, 0, 0.0), ctrl_points, new lambertian(vec3(1.0, 0.0, 0.0)));
list[2] = new sweeping_rotational(vec3( 5.0, 0, 0.0), ctrl_points, new lambertian(vec3(1.0, 0.0, 1.0)));
list[3] = new sphere(vec3(0,-103,-1), 100, new metal(vec3(0.0, 1.0, 0.5), 0.0), 0);
hitable *world = new hitable_list(list,4);
vec3 lookfrom(0, 10, 20);
vec3 lookat(0, 1, 0);
float dist_to_focus = (lookfrom - lookat).length();
float aperture = 0.0;
camera cam(lookfrom, lookat, vec3(0,1,0), 25, float(nx)/float(ny), aperture, 0.7*dist_to_focus);
输出图形如下：
一个玻璃材质，两个漫射材质，最下面的大球是镜面材质的
三个全是红色漫射材质时的图形如下：
1万+
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
