每条光线相当于一条射线，具有两个固定属性（起点o以及方向d，此外参数t表示光线的长度）
本节中所学习的光线类型为摄影机光线或主光线。对图像中的每一个像素，我们需要构造一条相机射线，然后将其投影到场景中。如果光线与对象相交，我们这些交点计算对象颜色，并将这些颜色分配给相应的像素。（我们需要区分主光线（投射到场景中的第一条光线，其起点是相机的原点）和副光线之间的区别。）
在几乎所有的3D应用程序中，创建摄像机时的默认位置都是世界原点，世界原点坐标为（0，0，0）。在光线追踪中，通常将平面位置放在距离相机原点正好1个单位的位置（此距离永远不变）。并且，我们还将摄像头朝着负Z轴的方向。（光线这里定义的都是从眼睛到物体，从物体到光源）
我们首先将一个6×6像素的图像投影到 2×2大小的区域中。
我们需要将栅格空间中的像素对于到世界空间中，寻找两个之间的关系。
![[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-4Yk2f8S9-1595991732684)(C:\Users\Administrator\AppData\Roaming\Typora\typora-user-images\image-20200719155407579.png)]](https://img-blog.csdnimg.cn/20200810192540454.png?x-oss-process=image/watermark,type_ZmFuZ3poZW5naGVpdGk,shadow_10,text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3doeTE4NzY3MTgzMDg2,size_16,color_FFFFFF,t_70#pic_center)
我们首先使用框架的尺寸归一化这些像素坐标，这些新的坐标在NDC（Normalized Device Coordinates）空间中定义。
我们将像素位置添加小偏移0.5，是我们希望射线穿过像素的中心。NDC空间表示的像素坐标正在[0,1]范围内（这与在光栅化时的NDC([-1,1])不同)。
再从NDC空间转化到Screen space[-1,1]中。
但是由于在Raster space和NDC space中，左上角为(0,0)，因此正在y轴向下为正方向。所以式子修改如下：
具体步骤如下图：将像素中间点的坐标转换成世界坐标
该点的坐标在Raster space中表示（像素坐标加上偏移量0.5），然后转化为NDC space（坐标重新
映射到[0,1]范围），然后转换为屏幕空间重新映射到[-1,1]，最后应用相机到世界转换4×4矩阵将屏幕空间中的坐标转换成世界空间。
到目前为止，我们一直处理图像为正方形的情况。如果有一张尺寸为7×5像素的图像。当他被缩放时，如果还是按照比例计算出Screen space的话，因为横向的像素点较多，像素就无法方正，所以无法正确投影像素。我们可以将Screen space的横坐标范围根据图像的宽高比由原来的[-1,1]变成[-1.4,1.4]。
最后我们来考虑一下视角大小对投影的影响。到目前为止，屏幕空间中定义的任何点的y坐标都在[-1,1]的范围内。我们也知道图像平面与相机原点的距离相距一个单位。如果从侧面看摄像机的设置，则可以通过将摄像机的原点连接到胶片平面的顶部和底部边缘来绘制三角形。
因为我们知道从照相机原点到胶片平面的距离（1个单位）和胶片平面的高度（2个单位，因为它从y = 1到y = -1），所以我们可以使用一些简单的三角函数来找到角度直角ABC的一半，它是垂直角度的一半α。
随着视野的增大，我们能看到的场景也会变多，**相当于更往外的物品投影到这个像素点。**因此，我们可以将屏幕像素坐标乘以该数字来进行放大与缩小。
因为相机朝着负Z轴方向，所以像素在图像平面上的最终坐标为：
如果我们希望可以从任何特定角度渲染场景图像。将摄影机从原始位置移动后，可以使用4×4矩阵表示摄影机的平移与旋转（相机到世界 矩阵）。
代码如下：
float imageAspectRatio = imageWidth / imageHeight; // assuming width > height
float Px = (2 * ((x + 0.5) / imageWidth) - 1) * tan(fov / 2 * M_PI / 180) * imageAspectRatio;
float Py = (1 - 2 * ((y + 0.5) / imageHeight) * tan(fov / 2 * M_PI / 180);
Vec3f rayOrigin = Point3(0, 0, 0);
Matrix44f cameraToWorld;
cameraToWorld.set(...); // set matrix
Vec3f rayOriginWorld, rayPWorld;
cameraToWorld.multVectMatrix(rayOrigin, rayOriginWorld);
cameraToWorld.multVectMatrix(Vec3f(Px, Py, -1), rayPWorld);
Vec3f rayDirection = rayPWorld - rayOriginWorld;
rayDirection.normalize(); // it's a direction so don't forget to normalize
参数曲面：如下公式表示一个球体
隐式曲面：如下公式为球面的隐式方程
- 射线与隐式曲面的交点通常比其他几何类型更快地进行计算。
- 隐式曲面可以用作边界体积，加速射线几何交点测试
如果光线与球体有交点，如p与p’。
代码
bool solveQuadratic(const float &a, const float &b, const float &c, float &x0, float &x1)
{
float discr = b * b - 4 * a * c;
if (discr < 0) return false;
else if (discr == 0) x0 = x1 = - 0.5 * b / a;
else {
float q = (b > 0) ?
-0.5 * (b + sqrt(discr)) :
-0.5 * (b - sqrt(discr));
x0 = q / a;
x1 = c / q;
}
if (x0 > x1) std::swap(x0, x1);
return true;
}
bool intersect(const Ray &ray) const
{
float t0, t1; // solutions for t if the ray intersects
#if 0
// geometric solution
Vec3f L = center - orig;
float tca = L.dotProduct(dir);
// if (tca < 0) return false;
float d2 = L.dotProduct(L) - tca * tca;
if (d2 > radius2) return false;
float thc = sqrt(radius2 - d2);
t0 = tca - thc;
t1 = tca + thc;
#else
// analytic solution
Vec3f L = orig - center;
float a = dir.dotProduct(dir);
float b = 2 * dir.dotProduct(L);
float c = L.dotProduct(L) - radius2;
if (!solveQuadratic(a, b, c, t0, t1)) return false;
#endif
if (t0 > t1) std::swap(t0, t1);
if (t0 < 0) {
t0 = t1; // if t0 is negative, let's use t1 instead
if (t0 < 0) return false; // both t0 and t1 are negative
}
t = t0;
return true;
}
bool intersectPlane(const Vec3f &n, const Vec3f &p0, const Vec3f &l0, const Vec3f &l, float &t)
{
// assuming vectors are all normalized
float denom = dotProduct(n, l);
if (denom > 1e-6) {
Vec3f p0l0 = p0 - l0;
t = dotProduct(p0l0, n) / denom;
return (t >= 0);
}
return false;
}
bool intersectDisk(const Vec3f &n, const Vec3f &p0, const float &radius, const Vec3f &l0, const Vec3 &l)
{
float t = 0;
if (intersectPlane(n, p0, l0, l, t)) {
Vec3f p = l0 + l * t;
Vec3f v = p - p0;
float d2 = dot(v, v);
return (sqrtf(d2) <= radius);
// or you can use the following optimisation (and precompute radius^2)
// return d2 <= radius2; // where radius2 = radius * radius
}
return false;
}
光线与盒子的交点，要判断光线何时进入盒子，何时离开。
以二维平面为例，如图：
[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-v9LuwROD-1595991732715)
class Box3
{
public:
Box3(cont Vec3f &vmin, const Vec3f &vmax)
{
bounds[0] = vmin;
bounds[1] = vmax;
}
Vec3f bounds[2];
};
bool intersect(const Ray &r)
{
//进入、离开x轴的时间
float tmin = (min.x - r.orig.x) / r.dir.x;
float tmax = (max.x - r.orig.x) / r.dir.x;
if (tmin > tmax) swap(tmin, tmax);
//进入、离开y轴的时间
float tymin = (min.y - r.orig.y) / r.dir.y;
float tymax = (max.y - r.orig.y) / r.dir.y;
if (tymin > tymax) swap(tymin, tymax);
if ((tmin > tymax) || (tymin > tmax))
return false;
if (tymin > tmin)
tmin = tymin;
if (tymax < tmax)
tmax = tymax;
//进入、离开y轴的时间
float tzmin = (min.z - r.orig.z) / r.dir.z;
float tzmax = (max.z - r.orig.z) / r.dir.z;
if (tzmin > tzmax) swap(tzmin, tzmax);
if ((tmin > tzmax) || (tzmin > tmax))
return false;
if (tzmin > tmin)
tmin = tzmin;
if (tzmax < tmax)
tmax = tzmax;
return true;
}
优化：
1.根据射线方向获取最大最小值，减少不必要的判断
2.优化射线的方向，当dir为0或-0时，tmin和tmax为无穷大，将dir = 1/dir
class Ray
{
public:
Ray(const Vec3f &orig, const Vec3f &dir) : orig(orig), dir(dir)
{
invdir = 1 / dir;
sign[0] = (invdir.x < 0);
sign[1] = (invdir.y < 0);
sign[2] = (invdir.z < 0);
}
Vec3 orig, dir; // ray orig and dir
Vec3 invdir;
int sign[3];
};
bool intersect(const Ray &r) const
{
float tmin, tmax, tymin, tymax, tzmin, tzmax;
tmin = (bounds[r.sign[0]].x - r.orig.x) * r.invdir.x;
tmax = (bounds[1-r.sign[0]].x - r.orig.x) * r.invdir.x;
tymin = (bounds[r.sign[1]].y - r.orig.y) * r.invdir.y;
tymax = (bounds[1-r.sign[1]].y - r.orig.y) * r.invdir.y;
if ((tmin > tymax) || (tymin > tmax))
return false;
if (tymin > tmin)
tmin = tymin;
if (tymax < tmax)
tmax = tymax;
tzmin = (bounds[r.sign[2]].z - r.orig.z) * r.invdir.z;
tzmax = (bounds[1-r.sign[2]].z - r.orig.z) * r.invdir.z;
if ((tmin > tzmax) || (tzmin > tmax))
return false;
if (tzmin > tmin)
tmin = tzmin;
if (tzmax < tmax)
tmax = tzmax;
return true;
}
#include <cstdio>
#include <cstdlib>
#include <memory>
#include <vector>
#include <utility>
#include <cstdint>
#include <iostream>
#include <fstream>
#include <cmath>
#include <limits>
#include <random>
#include "geometry.h"
const float kInfinity = std::numeric_limits<float>::max();
std::random_device rd;
std::mt19937 gen(rd());
std::uniform_real_distribution<> dis(0, 1);
inline
float clamp(const float &lo, const float &hi, const float &v)
{ return std::max(lo, std::min(hi, v)); }
inline
float deg2rad(const float °)
{ return deg * M_PI / 180; }
inline
Vec3f mix(const Vec3f &a, const Vec3f& b, const float &mixValue)
{ return a * (1 - mixValue) + b * mixValue; }
struct Options
{
uint32_t width;
uint32_t height;
float fov;
Matrix44f cameraToWorld;
};
//基类对象
class Object
{
public:
Object() : color(dis(gen), dis(gen), dis(gen)) {}
virtual ~Object() {}
//计算交点，return true or false
virtual bool intersect(const Vec3f &, const Vec3f &, float &) const = 0;
// Method to compute the surface data such as normal and texture coordnates at the intersection point.
// See method implementation in children class for details
virtual void getSurfaceData(const Vec3f &, Vec3f &, Vec2f &) const = 0;
Vec3f color;
};
//计算二次方程根
bool solveQuadratic(const float &a, const float &b, const float &c, float &x0, float &x1)
{
float discr = b * b - 4 * a * c;
if (discr < 0) return false;
else if (discr == 0) {
x0 = x1 = - 0.5 * b / a;
}
else {
float q = (b > 0) ?
-0.5 * (b + sqrt(discr)) :
-0.5 * (b - sqrt(discr));
x0 = q / a;
x1 = c / q;
}
return true;
}
//继承类，球类
class Sphere : public Object
{
public:
Sphere(const Vec3f &c, const float &r) : radius(r), radius2(r *r ), center(c) {}
//射线相交测试
bool intersect(const Vec3f &orig, const Vec3f &dir, float &t) const
{
float t0, t1; // solutions for t if the ray intersects
#if 0
// geometric solution
Vec3f L = center - orig;
float tca = L.dotProduct(dir);
if (tca < 0) return false;
float d2 = L.dotProduct(L) - tca * tca;
if (d2 > radius2) return false;
float thc = sqrt(radius2 - d2);
t0 = tca - thc;
t1 = tca + thc;
#else
// analytic solution
Vec3f L = orig - center;
float a = dir.dotProduct(dir);
float b = 2 * dir.dotProduct(L);
float c = L.dotProduct(L) - radius2;
if (!solveQuadratic(a, b, c, t0, t1)) return false;
#endif
if (t0 > t1) std::swap(t0, t1);
if (t0 < 0) {
t0 = t1; // if t0 is negative, let's use t1 instead
if (t0 < 0) return false; // both t0 and t1 are negative
}
t = t0;
return true;
}
/*
在表面给定点设置表面数据：法线和纹理坐标
Phit为表面的点，Nhit为法线，tex为纹理坐标
*/
void getSurfaceData(const Vec3f &Phit, Vec3f &Nhit, Vec2f &tex) const
{
Nhit = Phit - center;
Nhit.normalize();
// In this particular case, the normal is simular to a point on a unit sphere
// centred around the origin. We can thus use the normal coordinates to compute
// the spherical coordinates of Phit.
// atan2 returns a value in the range [-pi, pi] and we need to remap it to range [0, 1]
// acosf returns a value in the range [0, pi] and we also need to remap it to the range [0, 1]
tex.x = (1 + atan2(Nhit.z, Nhit.x) / M_PI) * 0.5;
tex.y = acosf(Nhit.y) / M_PI;
}
float radius, radius2;
Vec3f center;
};
/*如果光线与对象相交，则返回true。
变量tNear设置为最接近的相交距离，并且hitObject是指向相交对象的指针。
如果未找到交集，则将变量tNear设置为infinity并将hitObject设置为null。*/
bool trace(const Vec3f &orig, const Vec3f &dir, const std::vector<std::unique_ptr<Object>> &objects, float &tNear, const Object *&hitObject)
{
tNear = kInfinity;
std::vector<std::unique_ptr<Object>>::const_iterator iter = objects.begin();
for (; iter != objects.end(); ++iter) {
float t = kInfinity;
if ((*iter)->intersect(orig, dir, t) && t < tNear) {
hitObject = iter->get();
tNear = t;
}
}
return (hitObject != nullptr);
}
//计算相交点的颜色（如果有）（否则返回背景颜色）
Vec3f castRay(
const Vec3f &orig, const Vec3f &dir,
const std::vector<std::unique_ptr<Object>> &objects)
{
Vec3f hitColor = 0;
const Object *hitObject = nullptr; // this is a pointer to the hit object
float t; // 光源到交点的距离
if (trace(orig, dir, objects, t, hitObject)) {
Vec3f Phit = orig + dir * t;
Vec3f Nhit;
Vec2f tex;
hitObject->getSurfaceData(Phit, Nhit, tex);
// Use the normal and texture coordinates to shade the hit point.
// The normal is used to compute a simple facing ratio and the texture coordinate
// to compute a basic checker board pattern
float scale = 4;
float pattern = (fmodf(tex.x * scale, 1) > 0.5) ^ (fmodf(tex.y * scale, 1) > 0.5);
hitColor = std::max(0.f, Nhit.dotProduct(-dir)) * mix(hitObject->color, hitObject->color * 0.8, pattern);
}
return hitColor;
}
//主要的渲染功能。我们在此迭代图像中的所有像素，生成主光线并将这些光线投射到场景中。
//帧缓冲区的内容保存到文件中。
void render(
const Options &options,
const std::vector<std::unique_ptr<Object>> &objects)
{
Vec3f *framebuffer = new Vec3f[options.width * options.height];
Vec3f *pix = framebuffer;
float scale = tan(deg2rad(options.fov * 0.5));
float imageAspectRatio = options.width / (float)options.height;
//使用相机到世界矩阵将射线原点（也就是相机原点，通过将坐标（0,0,0）的点转换到世界空间来进行变换）。
Vec3f orig;
options.cameraToWorld.multVecMatrix(Vec3f(0), orig);
for (uint32_t j = 0; j < options.height; ++j) {
for (uint32_t i = 0; i < options.width; ++i) {
/*生成主光线方向。计算射线在屏幕空间中的x和y位置。这在z = 1的图像平面上给出了一个点。
从那里开始，我们只需通过归一化vec3f变量来计算方向。这类似于在图像平面上的点和相机原点之间获取矢量，在相机空间中该矢量为（0,0,0）：
ray.dir = normalize（Vec3f（x，y，-1）-Vec3f（0））; */
#ifdef MAYA_STYLE
float x = (2 * (i + 0.5) / (float)options.width - 1) * scale;
float y = (1 - 2 * (j + 0.5) / (float)options.height) * scale * 1 / imageAspectRatio;
#elif
float x = (2 * (i + 0.5) / (float)options.width - 1) * imageAspectRatio * scale;
float y = (1 - 2 * (j + 0.5) / (float)options.height) * scale;
#endif
//不要忘记使用相机到世界的矩阵来改变射线方向。
Vec3f dir;
options.cameraToWorld.multDirMatrix(Vec3f(x, y, -1), dir);
dir.normalize();
*(pix++) = castRay(orig, dir, objects);
}
}
// Save result to a PPM image (keep these flags if you compile under Windows)
std::ofstream ofs("./out.ppm", std::ios::out | std::ios::binary);
ofs << "P6\n" << options.width << " " << options.height << "\n255\n";
for (uint32_t i = 0; i < options.height * options.width; ++i) {
char r = (char)(255 * clamp(0, 1, framebuffer[i].x));
char g = (char)(255 * clamp(0, 1, framebuffer[i].y));
char b = (char)(255 * clamp(0, 1, framebuffer[i].z));
ofs << r << g << b;
}
ofs.close();
delete [] framebuffer;
}
//在程序的主要功能中，我们创建场景（创建对象）并设置渲染选项（图像宽度和高度等）。然后，我们调用render函数（）。
int main(int argc, char **argv)
{
// creating the scene (adding objects and lights)
std::vector<std::unique_ptr<Object>> objects;
// generate a scene made of random spheres
uint32_t numSpheres = 32;
gen.seed(0);
for (uint32_t i = 0; i < numSpheres; ++i) {
Vec3f randPos((0.5 - dis(gen)) * 10, (0.5 - dis(gen)) * 10, (0.5 + dis(gen) * 10));
float randRadius = (0.5 + dis(gen) * 0.5);
objects.push_back(std::unique_ptr<Object>(new Sphere(randPos, randRadius)));
}
// setting up options
Options options;
options.width = 640;
options.height = 480;
options.fov = 51.52;
options.cameraToWorld = Matrix44f(0.945519, 0, -0.325569, 0, -0.179534, 0.834209, -0.521403, 0, 0.271593, 0.551447, 0.78876, 0, 4.208271, 8.374532, 17.932925, 1);
// finally, render
render(options, objects);
return 0;
}
当光线与一个物体相交的时候，不必遍历物体的所有三角形面，可以利用一个包围盒包住该物体，先求光线与该物体的包围盒相交，如果与包围盒不相交说明与物体也不相交。
首先如上图最左边所示，求出光线与x平面的交点，将先进入的交点(偏小的那个)记为 tmin, 后出去的交点(偏大的那个)记为 tmax，紧接着如中间图所示计算出光线与y平面的两个交点同样记为另外一组tmin, tmax，当然计算的过程中要注意如果任意的 t < 0，那么这代表的是光线反向传播与对应平面的交点。
BVH与前几种方法最显著的区别就是，不再以空间作为划分依据，而是从对象的角度考虑，即三角形面。
