问题三十九：怎么用ray tracing画圆环

39.1 数学推导

书上这么说：

39.2 看C++代码实现

----------------------------------------------tori.h------------------------------------------

tori.h

 

#ifndef TORI_H
#define TORI_H

#include <hitable.h>


class tori : public hitable
{
    public:
        tori() {}
        tori(vec3 cen, float ra, float rb, material *m) : center(cen), radius_a(ra), radius_b(rb), ma(m) {}
        virtual bool hit(const ray& r, float tmin, float tmax, hit_record& rec) const;
        vec3 center;
        float radius_a;
        float radius_b;
        material *ma;
};

#endif // TORI_H

 

----------------------------------------------tori.cpp------------------------------------------

tori.cpp

#include "tori.h"
#include <iostream>
using namespace std;

bool tori::hit(const ray& r, float t_min, float t_max, hit_record& rec) const {
        vec3 oc = r.origin() - center;
        float a4 = dot(r.direction(), r.direction()) * dot(r.direction(), r.direction());
        float a3 = 4*(dot(oc, r.direction())) * dot(r.direction(), r.direction());
        float a2 = 2*dot(r.direction(), r.direction()) * (dot(oc, oc) - (radius_a*radius_a+radius_b*radius_b))
                    + 2*(dot(oc, r.direction()))
                    - 4*radius_a*radius_a*r.direction().z()*r.direction().z();
        float a1 = 4*dot(oc, r.direction())*(dot(oc, oc) - (radius_a*radius_a+radius_b*radius_b))
                    - 8*radius_a*radius_a*oc.z()*r.direction().z();
        float a0 = (dot(oc, oc) - (radius_a*radius_a+radius_b*radius_b))                    - 4*radius_a*radius_a*(radius_b*radius_b-oc.z()*oc.z());

        float *roots = roots_quartic_equation(a4, a3, a2, a1, a0);
        float temp;
        if (roots[0] > 0.0001) {
            for (int i=1; i<int(roots[0]); i++) {
                for (int j=i+1; j<int(roots[0])+1; j++) {
                    if (roots[i] > roots[j]) {
                        temp = roots[i];
                        roots[i] = roots[j];
                        roots[j] = temp;
                    }
                }
            }
            for (int k=1; k<int(roots[0])+1; k++) {
                if (roots[k] < t_max && roots[k] > t_min) {
                    rec.t = roots[k];
                    rec.p = r.point_at_parameter(rec.t);
                    vec3 pc = rec.p - center;
                    float nx = 4*pc.x()*pc.x()*pc.x() + 4*pc.x()*(pc.y()*pc.y()+pc.z()*pc.z()-(radius_a*radius_a+radius_b*radius_b));
                    float ny = 4*pc.y()*pc.y()*pc.y() + 4*pc.y()*(pc.x()*pc.x()+pc.z()*pc.z()-(radius_a*radius_a+radius_b*radius_b));
                    float nz = 4*pc.z()*pc.z()*pc.z() + 4*pc.z()*(pc.x()*pc.x()+pc.y()*pc.y()+radius_a*radius_a-radius_b*radius_b);
                    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;
                    delete [] roots;
                    return true;
                }
            }
        }
        delete [] roots;
        return false;
}

 

----------------------------------------------main.cpp------------------------------------------

main.cpp

 

        hitable *list[2];
        list[0] = new sphere(vec3(0.0,-100,0), 100, new lambertian(vec3(0.8, 0.8, 0.0)));
        list[1] = new tori(vec3(0, 3.0, 0), 2.5, 0.5, new lambertian(vec3(1.0, 0.5, 0.5)));

        hitable *world = new hitable_list(list,2);

        vec3 lookfrom(0, 3, 5);
        vec3 lookat(0, 3, 0);
        float dist_to_focus = (lookfrom - lookat).length();
        float aperture = 0.0;
        camera cam(lookfrom, lookat, vec3(0,1,0), 80, float(nx)/float(ny), aperture, 0.7*dist_to_focus);

 

输出图片如下：

我只是想问问：这是什么鬼？？？

 

39.3 修改书上错误

查来查去，算来算去。发现是书上对一元四次方程的系数算错了。坑爹啊！！！

书上给的是这个：

实际应该是这个：

 

对应代码修改如下：

       float a4 = dot(r.direction(), r.direction()) * dot(r.direction(),r.direction());

       float a3 = 4*(dot(oc, r.direction())) *dot(r.direction(), r.direction());

       float a2 = 2*dot(r.direction(), r.direction()) * (dot(oc, oc) -(radius_a*radius_a+radius_b*radius_b))

                    + 4*(dot(oc,r.direction()))*(dot(oc, r.direction()))

                    +4*radius_a*radius_a*r.direction().z()*r.direction().z();

       float a1 = 4*dot(oc, r.direction())*(dot(oc,oc) - (radius_a*radius_a+radius_b*radius_b))

                    +8*radius_a*radius_a*oc.z()*r.direction().z();

       float a0 = (dot(oc, oc) - (radius_a*radius_a+radius_b*radius_b))*(dot(oc, oc) - (radius_a*radius_a+radius_b*radius_b))

                    -4*radius_a*radius_a*(radius_b*radius_b-oc.z()*oc.z());

 

修改后输出图片如下：

 

哦啦~第一张圆环图片应该正常了……

 

再来一张：

        hitable *list[2];
        list[0] = new tori(vec3(0, 3.4, 0), 3.2, 0.2, new lambertian(vec3(0.0, 1.0, 0.0)));
        list[1] = new tori(vec3(0, 4.9, 0), 1.3, 0.2, new lambertian(vec3(1.0, 0.0, 0.0)));
        hitable *world = new hitable_list(list,2);

        vec3 lookfrom(0, 3, 5);
        vec3 lookat(0, 3, 0);
        float dist_to_focus = (lookfrom - lookat).length();
        float aperture = 0.0;
        camera cam(lookfrom, lookat, vec3(0,1,0), 80, float(nx)/float(ny), aperture, 0.7*dist_to_focus);

输出图片如下：

等等……等等……中间那个白色夹杂着红色的圈是什么鬼？？？