レイをトレースしたい　第五回　アンチエイリアス

前回までの画像は境界付近がギザギザしてるのでその部分をならすアンチエイリアシングをしていこうと思います。

アンチエイリアスの方法は、描画対象のピクセルに対して光線を複数飛ばすことによって実現します。
描画対象のピクセルには大きさがあるので、その中に収まる座標をランダムにｎ点選び、そのn点分の交差判定の結果を平均することで平滑化をかけます。こうすることで２つ以上の物体がかぶさるピクセルは２つの物体の色が混じって表示され、私達の目がそれに対して補完をかけるためなめらかに見えます。

収まる座標は

$$(中心座標)+ k\cdot(ピクセル幅)$$

としてkを-0.5から0.5の値からランダムに選んで実現します。
ここでランダムに選ぶ本数がふえるほど室は良くなりますが、計算不可につながるのでそこそこにすることが大事です。

実際の実装は以下のようになります。

```cpp
#include <bits/stdc++.h>
#include "Eigen/Core"
using namespace std;
using namespace Eigen;
/*
Vector3d
Matrix3d
(i,j)の形で要素アクセス
A.array().sin()で要素ごと演算
A.norm()
A.normalize()
a.dot(b)
*/

const double PI = 3.141592653589793238;
const double PI2 = 6.28318530718;
const double eps = 1e-6;




double hit_sphere(Vector3d center, Vector3d ray, Vector3d camera_pos, double size){
    Vector3d oc = camera_pos-center;
    double a = ray.dot(ray);
    double b = ray.dot(oc)*2.0;
    double c = (oc).dot(oc)-size*size;
    double D = b*b-4*a*c;
    if (D < 0){
        return -1e30;
    }else{
        return (-b-sqrt(D))/(2.0*a);
    }
}

Vector3d norm_sphere(Vector3d center, Vector3d hit_point){
    Vector3d temp = (hit_point-center).normalized();
    return temp;
}

double saturate(double d){
    double ans = 0;

    if(d < 0){
        ans = 0.0;
    }else if(d > 1){
        ans = 1.0;
    }else{
        ans = d;
    }
    return ans;
}

void write_ppm(vector<vector<vector<int>>> pic){
    ofstream out;
    out.open("test.pnm");
    out << "P3" << endl;
    out << pic.size() << " " << pic[0].size() << endl;
    out << 255 << endl;
    for (int i = 0 ; i < pic.size(); ++i){
        for (int j = 0; j < pic[0].size(); ++j){
            out << pic[i][j][0] << " " << pic[i][j][1] << " " << pic[i][j][2] << endl;
        }
    }

}

Matrix3d Rod_rot(Vector3d axis, double angle){
    double theta = angle*PI/180;
    double c = cos(theta);
    double s = sin(theta);
    double x = axis(0);
    double y = axis(1);
    double z = axis(2);
    
    Matrix3d ret;

    ret << c+(1-c)*x*x    , x*y*(1-c)-z*s , x*z*(1-c)+y*s,
           y*x*(1-c)+z*s  , c+(1-c)*y*y   , y*z*(1-c)-x*s,
           z*x*(1-c)-y*s  , z*y*(1-c)+x*s , c+z*z*(1-c);
    return ret;
}

class camera{
    public:
        camera();
        Vector3d direction;
        Vector3d position;
        double angle;
        double rad_ang_half;
    Vector3d ray(double x, double y){
        Vector3d temp;
        temp << sin(rad_ang_half)*x, sin(rad_ang_half)*y, cos(rad_ang_half);
        temp = temp.normalized();
        return temp;
    }
    private:
        
};

camera::camera(){
    angle = 60;
    rad_ang_half = angle*0.5*PI/180;
    direction << 0,0,1;
    position << 0,0,-3;
}

int main(void){
    random_device rd{};
    int seed = rd();
    mt19937 mt(rd());
    uniform_real_distribution<> dist(0.01,0.99);
    int height = 512;
    int width = 512;
    double xstep = 2.0/width;
    double ystep = 2.0/height;
    const int ray_per_px = 10;
    vector<vector<vector<int>>> pic;
    pic = vector<vector<vector<int>>>(height,vector<vector<int>>(width,vector<int>(3,255)));
    camera cam ;
    Vector3d sp_pos;
    sp_pos << 0,0,2;
    Vector3d sp_pos2;
    sp_pos2 << 0,-20,10;
    for (int i = 0; i<height; ++i){
        for(int j = 0; j < width; ++ j){
            vector<double> buf = vector<double>(3,0);
            double u = xstep*j-1.0;
            double v = 1.0-ystep*i;
            for (int k = 0; k < ray_per_px; ++ k){
                double dx = dist(mt)-0.5;
                double dy = dist(mt)-0.5;
                double temp_u = u+xstep*dx;
                double temp_v = v+ystep*dy;
                Vector3d ray = cam.ray(temp_u,temp_v);
                double t1 = hit_sphere(sp_pos,ray,cam.position,1);
                double t2 = hit_sphere(sp_pos2,ray,cam.position,20);
                Vector3d hit_point = cam.position + ray*t1;
                Vector3d hit_point2 = cam.position + ray*t2;
                if(t1 > 0 || t2 > 0){
                if (abs(t1) < abs(t2)){
                    Vector3d norm = norm_sphere(sp_pos,hit_point);
                    buf[0] += saturate(norm(0))*255;
                    buf[1] += saturate(norm(1))*255;
                    buf[2] += saturate(-norm(2))*255;
                }else{
                    buf[0] += 0;
                    buf[1] += 255;
                    buf[2] += 0;
                }
                }else{
                    buf[0] += 255;
                    buf[1] += 255;
                    buf[2] += 255;
                } 

            }
            pic[i][j][0] = buf[0]/ray_per_px;
            pic[i][j][1] = buf[1]/ray_per_px;
            pic[i][j][2] = buf[2]/ray_per_px;
        }
    }
    write_ppm(pic);
    return 0;
}
```
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/244678/65b45a38-24da-38fe-c149-5d125ae9273c.png)

なめらかに補完されました。