连载：徒手写一个DICOM阅图软件（005）之图像的插值算法

本文链接：https://blog.csdn.net/microzhe/article/details/135885486

本文主要介绍“图像的操作插值算法”.

本文为了简明体现整个操作过程，因此程序框架及软件设计并不规范，所以代码仅供参考。

上一篇介绍了图像坐标系的转换，见：连载：徒手写一个DICOM阅图软件（004）之图像的操作-CSDN博客

0. 临近插值的问题

上篇中，我们完成了对图像的操作（移动、缩放、旋转）。尤其在图像的缩放和旋转中，图像会表现出马赛克状，如下图：

原因是在图像显示时，咱们采用的是临近插值算法，因此，在图像显示时，采用双线性插值即可改善此问题。

1.双线性插值的代码

代码如下：

void ShowWork::DrawInPixel( )
{
	short* pPixelBuffer = mData.GetPixelBuffer( ); // DICOM原始像素
	if ( pPixelBuffer == NULL )
		return;

	int pixelWidth = mData.GetWidth( ); // DICOM图像的宽度
	int pixelHeight = mData.GetHeight( ); // DICOM图像的高度

	float a, b;
	mData.GetWLFactor( a, b ); // 窗宽窗位转换系数

	int w = mShowSize[0]; // 显示尺寸
	int h = mShowSize[1];

	M3DMatrix44f ScreenToDicom;
	m3dInvertMatrix44( ScreenToDicom, mDicomToScreen );

	M3DVector3f v1, v2;

	UCHAR c = 0;
	for ( int y = 0; y < h; y++ )
	{
		for ( int x = 0; x < w; x++ )
		{
			m3dLoadVector3( v1, float( x ), float( y ), 0 ); // 当前窗口位置
			m3dTransformVector3( v2, v1, ScreenToDicom ); // 计算出原始图像中位置

			// 测试是否越界
			if ( v2[0] < 1 || v2[0] >= pixelWidth-1 ||
				 v2[1] < 1 || v2[1] >= pixelHeight-1 )
			{
				mShowPixel[x + y * w ] = 0; // 如果越界，则置0
				continue;
			}

			// 未越界，则取出CT值

/*			// 旧算法：临近插值
			int posX = int( v2[0] + 0.5f );
			int posY = int( v2[1] + 0.5f );
			short pixel = pPixelBuffer[ posX + posY * pixelWidth]; 
			float v = pixel * a + b; // 将原始数据进行转换，y=x*a+b
*/

			// 以下是双线性插值的算法
			int intX = int( v2[0] );
			int intY = int( v2[1] );

			float fx = v2[0] - intX;
			float fy = v2[1] - intY;

			short* tempP = pPixelBuffer + intX + intY*pixelWidth;

			// 取出周围4个点
			short v00 = *tempP;
			short v01 = *(tempP+1);
			short v10 = *(tempP+pixelWidth);
			short v11 = *(tempP+pixelWidth+1);

			float fxInv = 1.0f - fx;
			float fyInv = 1.0f - fy;

			float v0001 = v00*fxInv + v01*fx;
			float v1011 = v10*fxInv + v11*fx;
			float v00011011 = v0001 * fyInv + v1011 * fy;

			float v = v00011011 * a + b; // 计算显示的像素值
			if ( v < 0.0f )
				c = 0;
			else if ( v > 255.0f )
				c = 255;
			else
				c = UCHAR( v );

			mShowPixel[x + y * w] = RGB( c, c, c ); // 显示像素为32bit，因为显示灰度，所以3通道颜色一样。
		}
	}
}