VolumeRendering(一) 光线步进 RayMarching

这边文章中的代码主要参考：
Volumetric Raymarching Sample
这个Demo中的例子比较复杂，使用SDF生成动态烟雾效果，添加了光照。但在医学领域或其他学科中，我们只是用体绘制来对三维标量场进行可视化，并不需要这么复杂的光照。其中光源可以省略，也不需要复杂的SDF生成烟雾的模型，我们只需要一个立方体就可以了。所以，我写了下面这个demo，便于理解和学习。
Volume Rendering By Ray Marching
以下是简化后的代码，主要去掉了一些动态变化

#define LARGE_NUMBER 1e20
#define MAX_VOLUME_MARCH_STEPS 160//在Volume中光线步进的最大次数
#define MAX_SDF_DETECT_STEPS 15//最大探测次数，用于确定物体边界
#define MARCH_STRIDE 0.4//固定的光线步长（在Volume中）

struct Box{
    vec3 Position;
    vec3 EdgeLength;
};

Box mBox = Box(
    vec3(0,0,0),
    vec3(16.0)
);

struct Camera
{
    vec3 Position;//
    vec3 LookAt;    

    float ImageHeight;//成像高度//根据屏幕长宽比可求宽度
    float FocalDistance;//焦距
};
Camera mCamera = Camera(
    vec3(120, 20, -165),
    vec3(0, 0, 0),
    2.0,
    7.0
);

//https://iquilezles.org/www/articles/distfunctions/distfunctions.htm
float sdBox( vec3 p /*到中心的距离*/, vec3 b/*边长*/ )
{
  vec3 q = abs(p) - b;
  return length(max(q,0.0)) + min(max(q.x,max(q.y,q.z)),0.0);
}

float QueryVolumetricDistanceField( in vec3 pos)
{
    float sdfValue= sdBox( pos-  mBox.Position, mBox.EdgeLength);
    return sdfValue;
}

float IntersectRayMarch(in vec3 rayOrigin, in vec3 rayDirection, float maxD)//用于确定物体边界
{
	float precis = MARCH_STRIDE; //这个值太大，会导致边界判断不精确，导致Volume表面有很多曲线形的分界线
    float D = 0.0f;
    for(int i=0; i<MAX_SDF_DETECT_STEPS; i++ )
    {
	    float dis = QueryVolumetricDistanceField( rayOrigin + rayDirection * D);
        if( dis < precis || D > maxD ) 
            break;
        D += dis;
    }
    return  D >= maxD  ? -1.0 : D;//没有碰到物体则返回-1，否则返回深度（据摄像机距离）
}

void SetCamera(in vec2 _uv, in float _aspectRatio, out vec3 _rayOrigin, out vec3 _rayDirection)
{
    float ImageWidth = mCamera.ImageHeight * _aspectRatio;
    
    vec3 ImagePosition = mCamera.Position;
    //vec3 ImagePosition = CameraOrbit(0.3);
    vec3 CameraView = mCamera.LookAt - ImagePosition;
    float ViewLength = length(CameraView);
    vec3 CameraViewDir = CameraView / ViewLength;

    vec3 CameraRight = cross(CameraViewDir, vec3(0, 1, 0));    
    vec3 CameraUp = cross(CameraRight, CameraViewDir);

    vec3 focalPoint = ImagePosition - mCamera.FocalDistance * CameraViewDir;//焦点位置
    
    vec3 ImagePoint = ImagePosition;//用Image的中心点初始化成像点
    //根据uv坐标偏移成像点
    ImagePoint += CameraRight * (_uv.x * 2.0 - 1.0) * ImageWidth *.5;
    ImagePoint += CameraUp * (_uv.y * 2.0 - 1.0) * mCamera.ImageHeight *.5;
    
    _rayOrigin = focalPoint;
    _rayDirection = normalize(ImagePoint - focalPoint);
}
vec3 GetAmbientLight()
{
	return vec3(0.03);
}

vec3 Render( in vec3 rayOrigin, in vec3 rayDirection)
{
    //找到volume的边界(若有外部传入的模型，可以直接得到深度。这里的模型是在shader中用代码生成的，所以使用RayMarch的方式确定边界)
    float volumeDepth = IntersectRayMarch(rayOrigin, rayDirection, LARGE_NUMBER);
    
    vec3 volumetricColor = vec3(0.0f);

    //从Volume的边界开始RayMarch
    if(volumeDepth > 0.0)//若可以触碰到Volume(没有被其他物体遮挡或Volume不在这条路径上)
    {
        float signedDistance = .0;
        for(int i = 0; i < MAX_VOLUME_MARCH_STEPS; i++)
        {
            volumeDepth += max(MARCH_STRIDE, signedDistance);//若还没有到达Volume边界，则先步进到边界处
           // volumeDepth += MARCH_STRIDE;
            vec3 position = rayOrigin + volumeDepth * rayDirection;

            signedDistance = QueryVolumetricDistanceField(position);

			if(signedDistance < 0.0f)//在Volume中
            {
                vec3 conner = mBox.Position-mBox.EdgeLength/2.0;
                float value = texture(iChannel0, (position-conner)/mBox.EdgeLength.x).x;//选取3D纹理的一部分进行采样
                float target =  0.5;//+0.2 * sin(iTime);
               //只显示0.5到0.6之间的值，相当于一个超简单的Transfer Function
                if(value<target||value>target+0.1)
                    value =0.0;
                
                volumetricColor += value * GetAmbientLight();
            }
        }
    }
    return min(volumetricColor, 1.0f);
}


void mainImage( out vec4 fragColor, in vec2 fragCoord )
{
    vec2 uv = fragCoord.xy / iResolution.xy;
    
    float aspectRatio = iResolution.x /  iResolution.y;
    
    vec3 rayOrigin,rayDirection;
    
    SetCamera(
        uv, aspectRatio, //输入
        rayOrigin, rayDirection//输出
        );
        
    vec3 color = Render(rayOrigin, rayDirection);
    
    fragColor=vec4( color, 1.0 );
}

这个直接copy到https://www.shadertoy.com上，在iChannel0中选择Volumes->GrayNoise3D就可以运行了。

这个demo中只有一个pixel shader，所以数据区域是靠SDF和RayMarching生成的，如果有顶点着色器，可以向顶点着色器传入模型数据，那样的话深度的计算就不用这么复杂了。使用SDF的好处是不需要任何顶点数据，就可以在PixelShader中利用RayMarching得到光线与复杂物体的交点，从而进行光线的计算。

接下来分析一下这个Demo中实现的体绘制算法
1.确定光线发射点和发射方向

SetCamera(uv, aspectRatio, //输入
        rayOrigin, rayDirection//输出
        );

在这个简单的场景中，我没有添加其他光源，由数据Volume自己反射环境光，由相机受光体（屏幕）接受光线。根据光线可逆原则，我们可以从屏幕发射光线，经过焦点，射到Volume数据体上。而屏幕可以按焦点对称到焦点和Volume之间。如下图所示

像素着色器执行一次就是计算一条光线，也就是计算一个像素的颜色。

2.找到光线与volume的交点，也就是光线在Volume上的入射点

float volumeDepth = IntersectRayMarch(rayOrigin, rayDirection, LARGE_NUMBER);

因为没有顶点着色器，没有顶点数据，所以使用立方体的有向距离场（signed distance field）代表数据范围。使用光线步进的方法，从RayOrigin出发，沿着RayDirection不停的向前移动，每次移动的距离就是该点到数据体表面最近的距离，这个距离由SDF得到。Ray Marching Demo For Beginner
这个例子中将算法流程演示的很清楚。
我使用立方体代表数据体。当光线步进的步长小于一个阈值precis的时候，就不再步进了。这个precis很小，此时距离数据体表面已经非常接近了，而这个点就是下一步中RayMarching的起点。

3.从Volume上的入射点开始RayMarching，直到离开数据体

if(volumeDepth > 0.0)//若可以触碰到Volume(没有被其他物体遮挡或Volume不在这条路径上)
    {
        float signedDistance = .0;
        for(int i = 0; i < MAX_VOLUME_MARCH_STEPS; i++)
        {
            volumeDepth += max(MARCH_STRIDE, signedDistance);//若还没有到达Volume边界，则先步进到边界处
           // volumeDepth += MARCH_STRIDE;
            vec3 position = rayOrigin + volumeDepth * rayDirection;

            signedDistance = QueryVolumetricDistanceField(position);

			if(signedDistance < 0.0f)//在Volume中
            {
                vec3 conner = mBox.Position-mBox.EdgeLength/2.0;
                float value = texture(iChannel0, (position-conner)/mBox.EdgeLength.x).x;//选取3D纹理的一部分进行采样
                float target =  0.5;//+0.2 * sin(iTime);
               //只显示0.5到0.6之间的值，相当于一个超简单的Transfer Function
                if(value<target||value>target+0.1)
                    value =0.0;
                
                volumetricColor += value * GetAmbientLight();
            }
        }
    }

在数据体中不断沿着光线方向步进，每到达一个位置，根据该位置在Volume立方体中的相对坐标采样3D纹理（三维数据一般保存在三维纹理中），得到一个value，也就是数据值。根据数据值的大小进行删选，上面的代码中只保留0.5~0.6的值。将这个值乘以环境光得到该位置的亮度，将这条光线上的所有亮度累加，得到该像素的值，也就是着色器的输出值。