voxelshader.h

3D reconstruction, medical image processing from colons, using intel image processing for based clas

#ifndef _VOXELSHADER_H_
#define _VOXELSHADER_H_

#pragma once

#include "DensityTable.h"

#define NORM_N (44)				/*   N parameter for normal computation */
static const int PIXEDPOINT_ACC=16;
static const int MAX_NORMAL=7924;	/*   maximum normal index */

class RxVoxelShader : public RxDensityTable
{
	// member variables
protected:
	unsigned short		*m_pVolume;
	unsigned short		*m_pnOriginalVolume;
	
	float		m_pfIntensity[MAX_NORMAL];


	short		*m_pnNormalPy;
	short		m_anNormalPyStorage[1 + 2 * NORM_N];		/* storage for NormalPy */
	short		*m_pnNormalXLimit;							/* max abs(x) allowed for a given y */
	short		m_anNormalXLimitStorage[1 + 2 * NORM_N];	/* storage for NormalXLimit */

	BOOL		m_bShade;


	// viewing
	int			m_iBright;
	int			m_iContrast;
	RxMatrix4D	m_mtxView;

	// volume size
	int			m_iVoxX;
	int			m_iVoxY;
	int			m_iVoxZ;


	// voi
	int			m_iVOIXm, m_iVOIXM;
	int			m_iVOIYm, m_iVOIYM;
	int			m_iVOIZm, m_iVOIZM;


	// voxel size
	double		m_fVoxelSizeX;
	double		m_fVoxelSizeY;
	double		m_fVoxelSizeZ;
	int			m_iRatioZ;

	// inverse
	int			m_iInverseVoxelSizeX;
	int			m_iInverseVoxelSizeY;
	int			m_iInverseVoxelSizeZ;

	int			m_iInterpolatedZ;


	// offset
	int			m_iOffsetX;
	int			m_iOffsetY;
	int			m_iOffsetZ;


	// view
	int			m_iPrincipalAxis;
	int			m_iDistance[2][3];
	float		m_fSamplingDistance;


	int			m_iShadeBufferSize;

	unsigned char	*m_pbyShadedWindowing;	// 寇何俊辑 且寸
	unsigned char	*m_pbyUnshadedWindowing;	// 寇何俊辑 且寸
	unsigned char	*m_pLocalWindowing;			// 寇何俊辑 且寸
	unsigned int	m_aiSATUnshadedWindowing[DENSITY_RANGE+1];
	unsigned int	m_aiSATShadedWindowing[DENSITY_RANGE+1];

	__m128			*m_aByte2M128Table;

public:

	RxVoxelShader();
	~RxVoxelShader();

	void operator=(const RxVoxelShader& info);

	// member functions
protected:
	BOOL		InitNormalTables();
	int			GetNormalIndexFromPosition(int iVoxelOffset);
	int			GetNormalIndexFromPosition(int iVoxelOffset, int iInternalNum);
	int			GetGradientMagnitude(int nOffset);
	BOOL		MakeShadeTable();		// from setviewing
	
public:
	int GetNormalIndexFromNormalVector(int iX, int iY, int iZ) {
		int n, x, y, xlimit;
		int denom;

		if((iX | iY | iZ) == 0)
			return 0;

		denom = (iX < 0) ? -iX : iX;
		denom += (iY < 0) ? -iY : iY;
		denom += (iZ < 0) ? -iZ : iZ;

		__m128 vDenom = _mm_cvtsi32_ss(_mm_setzero_ps(), denom);
		__m128 vDenomInv = _mm_rcp_ss(vDenom);
		vDenomInv = _mm_mul_ss(vDenomInv, _mm_set_ss(NORM_N));

		__m128 vX = _mm_cvtsi32_ss(_mm_setzero_ps(), iX);
		vX = _mm_mul_ss(vDenomInv, vX);
		x = _mm_cvtss_si32(vX);

		__m128 vY = _mm_cvtsi32_ss(_mm_setzero_ps(), iY);
		vY = _mm_mul_ss(vDenomInv, vY);
		y = _mm_cvtss_si32(vY);


	//	denominv = (float)NORM_N / denom;
	//	x = (int)(iX * denominv + 0.5);
	//	y = (int)(iY * denominv + 0.5);

		/* clamp x */
		xlimit = m_pnNormalXLimit[y];
		if (x < 0) {
		if (-x > xlimit)
		  x = -xlimit;
		}
		else {
		if (x > xlimit)
		  x = xlimit;
		}

		n = m_pnNormalPy[y] + (x << 1);
		if (iZ < 0)
		n |= 1;

		return n;
	}


	float		GetAlpha(int iI, int iJ, int iK);
	__m128		GetShadedVector(__m128 &vColor, int iOffset);
	__m128		GetShadedVector(__m128 &vColor, int iOffset, int iInternalNum);


	// windowing
	__m128		GetCutPlaneColor(unsigned short unDensity) {
		__m128 vColor = GetColor(unDensity);
		__m128 vWindow = m_aByte2M128Table[m_pbyShadedWindowing[unDensity]];

		return _mm_mul_ps(vColor, vWindow);
	};

	BOOL				GetWindowingColor(unsigned short unDensity, float color[4]) {
		__m128 vColor = GetColor(unDensity);
		__m128 vWindow = m_aByte2M128Table[m_pbyUnshadedWindowing[unDensity]];

		__m128 vRetval = _mm_mul_ps(vColor, vWindow);
		_mm_storeu_ps(color, vRetval);

		return TRUE;
	}
	BOOL				GetPreIntegratedWindowingColor(int nPrevDensity, int nCurrentDensity, float color[4]) {
		int nMaxDensity = __max(nPrevDensity, nCurrentDensity);
		int nMinDensity = __min(nPrevDensity, nCurrentDensity);

		int nRange = nMaxDensity - nMinDensity + 1;
		__m128 vRange = _mm_set1_ps((float)nRange);	// cvt肺 函券
		vRange = _mm_rcp_ps(vRange);

		int iWindowing = m_aiSATUnshadedWindowing[nMaxDensity+1] -  m_aiSATUnshadedWindowing[nMinDensity];
		__m128 vWindowing = _mm_set1_ps((float)iWindowing / 255.0f);
		vWindowing = _mm_move_ss( _mm_mul_ps(vWindowing, vRange), _mm_set_ss(1.0f));
//		vWindowing = m_aByte2M128Table[m_pbyUnshadedWindowing[nCurrentDensity]];

		__m128 vColor = _mm_sub_ps(m_pPreIntegratedTable[nMaxDensity], m_pPreIntegratedTable[nMinDensity-1]);
		_mm_storeu_ps(color, _mm_mul_ps(vWindowing, _mm_mul_ps(vColor, vRange)));

		return TRUE;
	}

	__m128				GetWindowingColor(unsigned short unDensity) {
		__m128 vColor = GetColor(unDensity);
		__m128 vWindow = m_aByte2M128Table[m_pbyUnshadedWindowing[unDensity]];

		__m128 vRetval = _mm_mul_ps(vColor, vWindow);

		return vRetval;
	};

	__m128				GetLocalWindowing(unsigned short unDensity) {
		return m_aByte2M128Table[m_pLocalWindowing[unDensity]];
	}

	__m128				GetWindowingColor(__m128 &vColor, unsigned short unDensity) {
		return _mm_mul_ps(vColor, m_aByte2M128Table[m_pbyUnshadedWindowing[unDensity]]);
	}

	unsigned char inline	GetBYTELocalWindowing(unsigned short unDensity) {
		return m_pLocalWindowing[unDensity];
	}

	__m128		GetRawWindowing(int iI, int iJ, int iK, float fAlpha=1.0f);
	__m128		GetRawWindowing(int iI, int iJ, int iK, int iInternalNum, float fAlpha=1.0f);


	BOOL		MakeDensityTable();
	
	BOOL		SetVolumeParameters(unsigned short *p3DVolume, unsigned short *pnOriginalVolume,
						const int iVolumeSizeX, const int iVolumeSizeY, const int iVolumeSizeZ,
						const double fVoxelSizeX, const double fVoxelSizeY, const double fVoxelSizeZ,
						int iRatioZ=1);
	BOOL		SetViewing(const RxMatrix4D &mtrxView);
	BOOL		SetShadingMode(BOOL bMode);
	BOOL		GetShadingMode() { return m_bShade; };
	BOOL		SetBrightContrast(int iBright, int iContrast);
	BOOL		SetWindowing(unsigned char *pShadedWindowing, unsigned char *pUnshadedWindowing, unsigned char *pLocalWindowing);
	BOOL		SetWindowingTable(unsigned char *pShadedWindowing, unsigned char *pUnshadedWindowing);

	BOOL		Set3DVOI(const int iVOIXm, const int iVOIXM,
						const int iVOIYm, const int iVOIYM,
						const int iVOIZm, const int iVOIZM);

	unsigned short GetOriginalDensity(int iI, int iJ, int iK) {
		int nOffset1 = ((iK*m_iVoxY) + iJ)*m_iVoxX + iI;
		return m_pnOriginalVolume[nOffset1];
	}
	unsigned short GetInterpolatedOriginalDensity(int iI, int iJ, int iK, int iInternalNum) {
		int nOffset1 = ((iK*m_iVoxY) + iJ)*m_iVoxX + iI;
		int nOffset2 = nOffset1 + m_iVoxY*m_iVoxX;
		unsigned short unDensity1 = m_pnOriginalVolume[nOffset1];
		unsigned short unDensity2 = m_pnOriginalVolume[nOffset2];
		unsigned int uiDensity = iInternalNum*unDensity2 + ((1<<PIXEDPOINT_ACC)-iInternalNum)*unDensity1;

		return (unsigned short)(uiDensity>>PIXEDPOINT_ACC);
	}

	unsigned short GetInterpolatedDensity(int iI, int iJ, int iK, int iInternalNum) {
		int nOffset1 = ((iK*m_iVoxY) + iJ)*m_iVoxX + iI;
		int nOffset2 = nOffset1 + m_iVoxY*m_iVoxX;
		unsigned short unDensity1 = m_pVolume[nOffset1];
		unsigned short unDensity2 = m_pVolume[nOffset2];
		unsigned int uiDensity = iInternalNum*unDensity2 + ((1<<PIXEDPOINT_ACC)-iInternalNum)*unDensity1;

		return (unsigned short)(uiDensity>>PIXEDPOINT_ACC);
	}

//	int			IsTransparent(int nOffset) {
//		return RxDensityTable::IsTransparent(m_pVolume[nOffset]);
//	}

//	int			IsTransparent(int iI, int iJ, int iK) {
//		int nOffset = ((iK*m_iVoxY) + iJ)*m_iVoxX + iI;
//		unsigned short unDensity = m_pVolume[nOffset];
//		return RxDensityTable::IsTransparent(unDensity);
//	}

//	int			IsTransparent(int iI, int iJ, int iK, int iInternalNum) {
//		unsigned short unDensity = GetInterpolatedDensity(iI, iJ, iK, iInternalNum);
//		return RxDensityTable::IsTransparent(unDensity);
//	}
};

#endif