视锥体裁剪射线的算法

射线Ray(直线情况)需要满足的条件：

1. 在视野中显示的粗细均匀，需要分段绘制，每段的粗细根据到视野的距离计算
2. 射线model的顶点尽量少以节省性能损耗

要满足条件2的话需要对射线进行裁剪，只绘制射线在视锥体内的部分，因此需要计算射线被视锥体裁剪后新的起点和终点

1. 计算三角形面积

float triangleArea(const float3& pa,const float3& pb,const float3& pc)
{
    float a = length(pb-pc);
    float b = length(pa-pc);
    float c = length(pa-pb);
    float p=(a+b+c)/2.0;
    return sqrt(p*(p-a)*(p-b)*(p-c));
}

2. 法向量法判断点是否在三角形内

代码：

bool isInTriangle2(const float3& pa,const float3& pb,const float3& pc, const float3&tP)
{
    float3 PA = pa - tP;
    float3 PB = pb - tP;
    float3 PC = pc - tP;
    float3 t1 = cross(PA,PB);
    float3 t2 = cross(PB,PC);
    float3 t3 = cross(PC,PA);
    return dot(t1,t2) >= 0 && dot(t1,t3) >= 0 && dot(t2,t3)>=0;
};

2. 计算直线穿过某平面的交点

/*
*param [in] o_orign 射线起点
*param [in] o_dir  射线方向
*param [in] fn 平面方程，点法式
*param [out] inter_pnt交点
*param [return] 是否有交点
*/
bool LineRayToPlanePnt(const float3& o_orign,const float3& o_dir,const float4& fn, float3& inter_pnt)
{
float3 N = float3(fn[0], fn[1], fn[2]);
float D = fn[3];
if (std::abs(dot(o_dir,N)) < 1e-8)
{
return false;
}
float t = -(dot(o_orign,N) + D) / (dot(o_dir,N));
inter_pnt = o_orign + t*o_dir;
    return true;
}

2.根据三点获取点法式平面方程

float4 getPlane(float3 a,float3 b,float3 c){
    float3 normal = normalize(cross(b-a,c-a));
    return float4(normal, -dot(normal,a));
}

3.计算射线和单个四边形的交点

先计算直线和平面的交点，然后缩小范围判断该交点是否在三角形内

/*
*param [in] o_orign 射线起点
*param [in] o_dir  射线方向
*param [in] pa,pb,pc,pd 四边形四个顶点
*param [out] inter_pnt 交点
*param [return] 是否有交点
*/
bool isInQuadrangle(const float3& o_orign,const float3& o_dir,const float3& pa,const float3& pb,const float3& pc ,const float3& pd, float3& inter_pnt){
    float4 plane = getPlane(pa,pb,pc);
    if(LineRayToPlanePnt(o_orign,o_dir,plane,inter_pnt)){
        if(dot(inter_pnt-o_orign,o_dir)>=0&&isInTriangle2(pa,pb,pc,inter_pnt)){
            return true;
        }
        if(dot(inter_pnt-o_orign,o_dir)>=0&&isInTriangle2(pb,pc,pd,inter_pnt)){
            return true;
        }
    }
    return false;
}

上述代码中，dot(inter_pnt-o_orign,o_dir)>=0为判断该交点是否在射线的相反方向上，是的话忽略

3.计算射线新的起点和终点

void getNewEndPoint(const std::shared_ptr<scene::Camera> camera,const float3& p1,const float3& p2, float3& new_p1,float3& new_p2){
    if(isnan(p1.x)||isnan(p1.y)||isnan(p1.z)||isnan(p2.x)||isnan(p2.y)||isnan(p2.z)){
        return;
    }

	//计算视锥体的八个顶点，这种计算比矩阵比较繁琐，计算量大
    // float3 cameraForward = normalize(camera->getForward());
    // float3 camera_pos = camera->getWorldPosition();
    // float fov_rad = camera->getVerticalFovDegrees() * filament::math::d::DEG_TO_RAD;
    // float3 view_right = -normalize(camera->getLeft());
    // float3 view_up = normalize(camera->getUp());
    // float4 bound_n = sphere::calcViewBound(camera->getVerticalFovDegrees(), camera->getAspect(), camera->getNearClipPlane());
    // float4 bound_f = sphere::calcViewBound(camera->getVerticalFovDegrees(), camera->getAspect(), camera->getFarClipPlane());
    // float3 n_left_down = camera_pos + cameraForward * camera->getNearClipPlane() + bound_n.x * view_right + view_up * bound_n.z;
    // float3 f_right_up = camera_pos + cameraForward * camera->getFarClipPlane() + bound_f.y * view_right + view_up * bound_f.w;
    // float cos_fov = cos(fov_rad);
    // float sin_fov = sin(fov_rad);
    // float tan_fov = sin(fov_rad);
    // float3 near_p = camera_pos + cameraForward * camera->getNearClipPlane();
    // float tmp = camera->getNearClipPlane()*tan_fov*cos_fov;
    // float3 normal_up = near_p + view_up * tmp - cameraForward* tmp * sin_fov;

    float3 ray_dir = p2 - p1;  //射线方向
    //视锥体八个顶点在最初NDC空间下的坐标
    float3 points[8] ;
    points[0] = float3(1, -1, 1);  //右下
    points[1] = float3(-1, -1, 1); //左下
    points[2] = float3(-1, 1, 1);  //左上
    points[3] = float3(1, 1, 1);   //右上
    points[4] = float3(1, -1, -1); //右下
    points[5] = float3(-1, -1, -1);//左下
    points[6] = float3(-1, 1, -1); //左上
    points[7] = float3(1, 1, -1);  //右上
    mat4 pro = camera->getProjectionMatrix();
    mat4 model = camera->getMainModelMatrix();
    auto inverseMat = inverse(pro* inverse(model));
    for(int i=0;i<8;i++){ //NDC空间转到世界空间
        float4 tmp = inverseMat * points[i];
        points[i] = tmp.xyz/tmp.w; 
        LOGD("ray_point view_plane_new:[%f,%f,%f]",points[i].x,points[i].y,points[i].z);
    }
    
    int findOne=0;
    auto setNewPoint = [&](float3 new_p){
        if(findOne==0){
            new_p1 = new_p;
            findOne ++;
        }else{
            new_p2 = new_p;
            findOne ++;
        }
    };

    struct quadrangle{
        float3 a,b,c,d;
    };
    quadrangle quadrangles[6]={ //视锥体八个面
        {points[0],points[1],points[3],points[2]}, 	//远平面
        {points[4],points[5],points[7],points[6]},	//近平面
        {points[1],points[2],points[5],points[6]}, //左
        {points[2],points[3],points[6],points[7]}, //上
        {points[0],points[3],points[4],points[7]}, //右
        {points[0],points[1],points[4],points[5]}, //下
    };
    float3 intersect_point;
    for(auto q:quadrangles){	//分别计算射线和八个面的焦点
        if(isInQuadrangle(p1,ray_dir,q.a,q.b,q.c,q.d,intersect_point)){
            setNewPoint(intersect_point);
        }
    }
    if(findOne==0){  //射线和视锥体无碰撞
        new_p1 = p1;
        new_p2 = p2;
        LOGD("ray_point findOne==0:start_new:[%f,%f,%f],end_new:[%f,%f,%f]",
                new_p1.x,new_p1.y,new_p1.z,new_p2.x,new_p2.y,new_p2.z);
    }else if(findOne==1){	//射线和视锥体有一个交点(射线起点在视锥体中)
        new_p2 = new_p1;
        new_p1 = p1;
        LOGD("ray_point findOne==1:start_new:[%f,%f,%f],end_new:[%f,%f,%f]",
                new_p1.x,new_p1.y,new_p1.z,new_p2.x,new_p2.y,new_p2.z);
    }else{  //射线和视锥体两个交点，判断和原始起点近的作为新起点
        if(length(new_p1-p1)>length(new_p2-p1)){
            intersect_point = new_p1;
            new_p1 = new_p2;
            new_p2 = intersect_point;
        }
        if(length(new_p2-p1)>length(p2-p1)){
            new_p2 = p2;
        }
        LOGD("ray_point findOne==%d:start_new:[%f,%f,%f],end_new:[%f,%f,%f]",findOne,
                new_p1.x,new_p1.y,new_p1.z,new_p2.x,new_p2.y,new_p2.z);
    }
}