シェーディング処理

小一時間で実装してみた。unsigned char *bufに既に色が入っていて、double *zbufにz値(手前が小さい)が入っている場合、周囲48方向の最大仰角を平均して空が見える立体角をざっくり計算する。それに比例して明るさを調節する。サンプル方向が少ないため境界が見えるが、まあ筋は悪くない。速度もそれなり。
以下コード

double maxSlope(int x0, int y0, int dx, int dy)
{
    double z0 = zbuf[x0+y0*XSIZE];

    double maxsl=-1000.0;

    int i;
    for(i=1;i<XSIZE;i++)
    {
	int x=x0+dx*i;
	int y=y0+dy*i;
	if(x<0 || x>=XSIZE || y<0 || y>=YSIZE) break;
	double z=-(zbuf[x+y*XSIZE]-z0)/i;
	if (z>maxsl)maxsl=z;
    }
    
    // t=tan(max lookup angle)
    double t= maxsl/sqrt(dx*dx+dy*dy) /( (xmax-xmin)/XSIZE);
    // back to degree
    double deg = atan(t);
    // Visible solid angle =int_deg^pi/2 dphi dtheta cos(theta) = dphi (1-sin(deg))
    // Normalize max to 1
    return (1.0-sin(deg))/2;
}

double GetShade(int x, int y, int nx, int ny)
{
    return 
	(
	 +maxSlope(x,y, nx, ny)
	 +maxSlope(x,y, -nx, ny)
	 +maxSlope(x,y, nx, -ny)
	 +maxSlope(x,y, -nx, -ny)
	 +maxSlope(x,y, ny, nx)
	 +maxSlope(x,y,-ny, nx)
	 +maxSlope(x,y, ny,-nx)
	 +maxSlope(x,y,-ny,-nx)
	)/8;
}

void AddShade()
{
    int x,y;
    for(x=0;x<XSIZE;x++){
	printf("%d/%d\n",x,XSIZE);
    for(y=0;y<YSIZE;y++){
	int ad=x+y*XSIZE;
        if(zbuf[ad]>9999999.0)continue;
	double shd= (
	    +maxSlope(x,y, -1, 0)
	    +maxSlope(x,y,  1, 0)
	    +maxSlope(x,y,  0, 1)
	    +maxSlope(x,y,  0,-1)
	    +maxSlope(x,y, -1, -1)
	    +maxSlope(x,y,  1, -1)
	    +maxSlope(x,y, -1,  1)
	    +maxSlope(x,y,  1,  1) )/8
	    +GetShade(x,y, 1,2)
	    +GetShade(x,y, 1,3)
	    +GetShade(x,y, 2,3)
	    +GetShade(x,y, 1,4)
	    +GetShade(x,y, 3,4)
	    ;
	shd = shd/6 * 2;
	shd = 0.2+0.8*shd;
	if(shd>1.0)shd=1.0;

	buf[ad*3+0] = (unsigned char)(buf[ad*3+0]*shd);
	buf[ad*3+1] = (unsigned char)(buf[ad*3+1]*shd);
	buf[ad*3+2] = (unsigned char)(buf[ad*3+2]*shd);
    }}
}