rxysegment_temp.cpp

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

								bBottom_Connectivity = true;
							}
						} else bBottom_Connectivity = false;
					} else bBottom_Connectivity = false;
				}
				if(bUp_Check){
					if(pLabel[nZUpIndex+nIndexNeighbor]==0 ){
						int nRefData = pImg[nZUpIndex+nIndexNeighbor];
						if((nMin<=nRefData) && (nRefData<=nMax)){
							if(!bUp_Connectivity){
								Q.Enqueue(nZUpIndex+nIndexNeighbor);
								
								bUp_Connectivity = true;
							}
						} else bUp_Connectivity = false;
					}  else bUp_Connectivity = false;
				}
				if(bDown_Check){
					if( pLabel[nZDownIndex+nIndexNeighbor]==0 ){
						int nRefData = pImg[nZDownIndex+nIndexNeighbor];
						if((nMin<=nRefData) && (nRefData<=nMax)){
							if(!bDown_Connectivity){
								Q.Enqueue(nZDownIndex+nIndexNeighbor);								
								bDown_Connectivity = true;
							}
						} else bDown_Connectivity = false;
					} else bDown_Connectivity = false;
				}
			} 
			else	break;
			
		}
	}while(!Q.IsEmpty());// Q俊 甸绢埃 point啊 家柳瞪 锭鳖瘤

	return bRet;
}



/////////////////////////////////
// Edge detection
/////////////////////////////////
BOOL RxySegment::EdgeDetection(unsigned short *pOriginalVolume, int iVolX, int iVolY, int iSliceNum, int nMode,
							   RxProgressWnd* pWndProgress, RxPoint3D<float> ptMaxBound, RxPoint3D<float> ptMinBound)
{
	TRACE(_T("Start Edge Detection!\n"));
	RxxFilters Filter;
	LPBYTE pbyCanny;
	int i, j, k;
	int w, x, y, z;
	long iIndex;
	
	RxVolumeData* pVolume = RxGetVolumeData(nMode);
	int iSliceSize = iVolX*iVolY;

	RxMaskVolume* pMask = pVolume->GetEdgeVolume();
	pMask->ResetMask(0x00);

/*	int tempIndex = (int)(iSliceSize*iSliceNum*0.5);
	short* xArr = new short[tempIndex]; 
	short* yArr = new short[tempIndex];
	short* zArr = new short[tempIndex];
	int nNumEdge = 0;
*/
	unsigned short* pData = new unsigned short[iSliceSize];			
	for( k = 0; k < iSliceNum; k++) {
		memcpy(pData, &(pOriginalVolume[k*iSliceSize]), sizeof(unsigned short)*iSliceSize);
		
		pbyCanny = Filter.SobelProcessing(pData, iVolX, iVolY);

		for (j = 0; j< iVolY; j++){
			for (i = 0; i< iVolX; i++) {
				iIndex = j*iVolX+i;
			
				if(pbyCanny[iIndex]) {
					pMask->SetMask(i+ptMinBound.x, j+ptMinBound.y, k+ptMinBound.z);		
		/*			xArr[nNumEdge] = i+ptMinBound.x;
					yArr[nNumEdge] = j+ptMinBound.y;
					zArr[nNumEdge] = k+ptMinBound.z;
					nNumEdge++;
*/
				}
			}
		}
		delete pbyCanny;

		if(k%(iSliceNum/4)==0)
			ControlProgress(pWndProgress);
	}
/*	int dilation = 2;
	for(w = 0; w < nNumEdge; w++)
	{			
		i = xArr[w];
		j = yArr[w];
		k = zArr[w];
		for (x = -dilation; x <= dilation; x++)
			for (y = -dilation; y <= dilation; y++)
				for (z = -dilation; z <= dilation; z++)
				{
					if ((0 <= (i + x)) && ((i + x) < iVolX) &&
						(0 <= (j + y)) && ((j + y) < iVolY) &&
						(0 <= (k + z)) && ((k + z) < iSliceNum))
						pMask->SetMask(i + x, j + y, k + z);
				}
	}

	delete [] xArr;
	delete [] yArr;
	delete [] zArr;
*/
		// draw seed point
/*	
	for( j=0; j<5; j++)
		for( i=0; i<5; i++)
			pMask->SetMask( ptMinBound.x + m_ptSeed[k].x + i, ptMinBound.y + m_ptSeed[k].y + j, 
			ptMinBound.z + k );

*/
	delete[] pData;
	
	return TRUE;
}


///////////////////////////
// 2D Morphology
// 2004.02.24.
///////////////////////////

BOOL RxySegment::Dilate(unsigned char *pSrc, unsigned char *pDst, int width, int height, int iIter)
{
	ASSERT(width % 4 == 0);

	IplImage* srcImage;
	IplImage* dstImage;

	srcImage = iplCreateImageHeader(
			1, 0, IPL_DEPTH_8U, // data of byte type
			"Gray", "Gray", // color order
			//"Binary", "Binary", // for binary morphology
			IPL_DATA_ORDER_PIXEL, // channel arrangement
			IPL_ORIGIN_TL, // top left orientation
			IPL_ALIGN_QWORD, // 8 bytes align
			width, height, // image height
			NULL, NULL, NULL, NULL); // not tiled
	dstImage = iplCreateImageHeader(
			1, 0, IPL_DEPTH_8U, // data of byte type
			"Gray", "Gray", // color order
			//"Binary", "Binary", // for binary morphology
			IPL_DATA_ORDER_PIXEL, // channel arrangement
			IPL_ORIGIN_TL, // top left orientation
			IPL_ALIGN_QWORD, // 8 bytes align
			width, height, // image height
			NULL, NULL, NULL, NULL); // not tiled

	srcImage->imageData = (char *)pSrc;
	dstImage->imageData = (char *)pDst;

	iplDilate(srcImage, dstImage, iIter);
	
	iplDeallocate( srcImage, IPL_IMAGE_HEADER );
	iplDeallocate( dstImage, IPL_IMAGE_HEADER );

	return TRUE;
}

BOOL RxySegment::Dilate_CircleMask(unsigned char *pSrc, unsigned char *pDst, int width, int height, int iIter)
{
	int i, j, k;
	int row, column, temp;

	for( k=0; k<iIter; k++) {
		memset(pDst, 0x00, sizeof(unsigned char)*width*height);

		for( row=0; row<width; row++) {
			for( column=0; column<height; column++) {
				if(row >= 3 && row< width-3 && column >=3 && column < height-3) {
					temp = 1;
					for(j=0; j<7; j++) {
						for(i=0; i<7; i++) {
							if( (i-3)*(i-3) + (j-3)*(j-3) <= 10 ) {
								if(pSrc[(i + row-3) + width*(j + column-3)] == 0 ) { temp*=1; }			
								else { 	temp*=0; break; }
							}
						}
					}
					if(temp) pDst[row + width*column] = 0;			
					else pDst[row + width*column] = 255;	
										
				}
				// 版拌 贸府
				else if( ((row==2 || row==width-3) && (column>=2 && column<=height-3))
					|| ((column==2 || column==height-3) && (row>=2 && row <=width-3)) ) {
					temp = 1;
					for(j=0; j<5; j++) {
						for(i=0; i<5; i++) {
							if( (i-2)*(i-2) + (j-2)*(j-2) <= 5 ) {
								if( pSrc[(i + row-2) + width*(j + column-2)] == 0 ) { temp*=1; }
								else { temp*=0; break; }
							}
						}
					}
					if(temp)  pDst[row + width*column] = 0;						
					else pDst[row + width*column] = 255;					
				}
				else if( ((row==1 || row==width-2) && (column>=1 && column<=height-2))
					|| ((column==1 || column==height-2) && (row>=1 && row <=width-2)) ) {
					temp = 1;
					for(j=0; j<3; j++) {
						for(i=0; i<3; i++) {
							//if( (i-1)*(i-1) + (j-1)*(j-1) <= 1 ) {
								if( pSrc[(i + row-1) + width*(j + column-1)] == 0 ) { temp*=1; }
								else { temp*=0; break; }
							//}
						}
					}
					if(temp)  pDst[row + width*column] = 0;						
					else pDst[row + width*column] = 255;
				}
			}
		}
		memcpy(pSrc, pDst, width*height*sizeof(unsigned char) );
	}

	return TRUE;
}

BOOL RxySegment::Erose_CircleMask(unsigned char *pSrc, unsigned char *pDst, int width, int height, int iIter)
{
	int i, j, k;
	int row, column, temp;

	for( k=0; k<iIter; k++) {
		memset(pDst, 0x00, sizeof(unsigned char)*width*height);

		for( row=0; row<width; row++) {
			for( column=0; column<height; column++) {
				if(row >= 3 && row< width-3 && column >=3 && column < height-3) {
					temp = 1;
					for(j=0; j<7; j++) {
						for(i=0; i<7; i++) {
							if( (i-3)*(i-3) + (j-3)*(j-3) <= 10 ) {
								if( pSrc[(i + row-3) + width*(j + column-3)] == 255 ) { temp*=1; }							
								else { temp*=0; break; }								
							}
						}
					}
					if(temp)	pDst[row + width*column] = 255;
					else		pDst[row + width*column] = 0;
				}

				else if( ((row==2 || row==width-3) && (column>=2 && column<=height-3))
					|| ((column==2 || column==height-3) && (row>=2 && row <=width-3)) ) {
					temp = 1;
					for(j=0; j<5; j++) {
						for(i=0; i<5; i++) {
							if( (i-2)*(i-2) + (j-2)*(j-2) <= 5 ) {
								if( pSrc[(i + row-2) + width*(j + column-2)] == 255 ) { temp*=1; }
								else { temp*=0;  break; }								
							}
						}
					}
					if(temp) 	pDst[row + width*column] = 255;
					else		pDst[row + width*column] = 0;
				}
				else if( ((row==1 || row==width-2) && (column>=1 && column<=height-2))
					|| ((column==1 || column==height-2) && (row>=1 && row <=width-2)) ) {
					temp = 1;
					for(j=0; j<3; j++) {
						for(i=0; i<3; i++) {
							//if( (i-1)*(i-1) + (j-1)*(j-1) <= 1 ) {
								if(pSrc[(i + row-1) + width*(j + column-1)] == 255 ) { temp*=1; }
								else { temp*=0; break; }
							//}
						}
					}
					if(temp)  pDst[row + width*column] = 255;						
					else pDst[row + width*column] = 0;
				}
			}
		}

		memcpy(pSrc, pDst, width*height*sizeof(unsigned char) );
	}

	return TRUE;
}

BOOL RxySegment::Erose(unsigned char *pSrc, unsigned char *pDst, int width, int height, int iIter)
{
	ASSERT(width % 4 == 0);

	IplImage* srcImage;
	IplImage* dstImage;

	srcImage = iplCreateImageHeader(
			1, 0, IPL_DEPTH_8U, // data of byte type
			"Gray", "Gray", // color order
			//"Binary", "Binary", // for binary morphology
			IPL_DATA_ORDER_PIXEL, // channel arrangement
			IPL_ORIGIN_TL, // top left orientation
			IPL_ALIGN_QWORD, // 8 bytes align
			width, height, // image height
			NULL, NULL, NULL, NULL); // not tiled
	dstImage = iplCreateImageHeader(
			1, 0, IPL_DEPTH_8U, // data of byte type
			"Gray", "Gray", // color order
			//"Binary", "Binary", // for binary morphology
			IPL_DATA_ORDER_PIXEL, // channel arrangement
			IPL_ORIGIN_TL, // top left orientation
			IPL_ALIGN_QWORD, // 8 bytes align
			width, height, // image height
			NULL, NULL, NULL, NULL); // not tiled

	srcImage->imageData = (char *)pSrc;
	dstImage->imageData = (char *)pDst;

	iplErode(srcImage, dstImage, iIter);
	
	iplDeallocate( srcImage, IPL_IMAGE_HEADER );
	iplDeallocate( dstImage, IPL_IMAGE_HEADER );

	return TRUE;
}

BOOL RxySegment::Close(unsigned char *pSrc, unsigned char *pDst, int width, int height, int iIter)
{
	ASSERT(width % 4 == 0);

	IplImage* srcImage;
	IplImage* dstImage;

	srcImage = iplCreateImageHeader(
			1, 0, IPL_DEPTH_8U, // data of byte type
			"Gray", "Gray", // color order
		//	"Binary", "Binary",
			IPL_DATA_ORDER_PIXEL, // channel arrangement
			IPL_ORIGIN_TL, // top left orientation
			IPL_ALIGN_QWORD, // 8 bytes align
			width, height, // image height
			NULL, NULL, NULL, NULL); // not tiled
	dstImage = iplCreateImageHeader(
			1, 0, IPL_DEPTH_8U, // data of byte type
			"Gray", "Gray", // color order - for graylevel morphology
		//	"Binary", "Binary", // for binary morphology
			IPL_DATA_ORDER_PIXEL, // channel arrangement
			IPL_ORIGIN_TL, // top left orientation
			IPL_ALIGN_QWORD, // 8 bytes align
			width, height, // image height
			NULL, NULL, NULL, NULL); // not tiled

	srcImage->imageData = (char *)pSrc;
	dstImage->imageData = (char *)pDst;

	iplClose(srcImage, dstImage, iIter);

	iplDeallocate( srcImage, IPL_IMAGE_HEADER );
	iplDeallocate( dstImage, IPL_IMAGE_HEADER );

	return TRUE;
}

void RxySegment::Segmentation(int nMode)
{
	TRACE(_T("Start Segmentation\n"));

	// for calculating elapsed time
	char TimeStr[255];
	DWORD CurrentTime, ElapsedTime, SRGTime, ErosionTime, DilationTime, CCTime, CorrectionTime;

	// yeni - 2004.02.04
	int i, j, k;
	int label = 0;
	int thresholdmin, thresholdmax;
	int iVolX, iVolY, iVolZ, iSliceNum, iSliceSize;
	unsigned int iVolumeSize;
	int iIndexX, iIndexY,
		iIndexZ, iSeedIndex;

	RxPoint3D<float> ptMaxBound(-100, -100, -100);
	RxPoint3D<float> ptMinBound(5000, 5000, 5000);

	bool bSrgSuccess = false;	

	// 1. for getting original volume 
	RxVolumeData *pVolumeData = RxGetVolumeData(nMode); // 0 : Reference, 1 : float, 2 : combine
	unsigned short *pVolume = pVolumeData->GetBigVolume();
	pVolumeData->GetBigVolumeDmsn(&iVolX, &iVolY, &iSliceNum);
	iVolumeSize = iVolX*iVolY*iSliceNum;
	iSliceSize = iVolX*iVolY;

	// Progress bar
	RxProgressWnd*	pWndProgress;
	pWndProgress = new RxProgressWnd;
	int nSetRange = 18 + iSliceNum;
	pWndProgress->GoModal();
	pWndProgress->SetWindowSize(0);
	pWndProgress->SetRange(0, nSetRange);
	pWndProgress->SetText(_T("Segmenting..."));
 	ControlProgress(pWndProgress);

	CurrentTime = GetTickCount ();

	// 2. code for searching bounding box
	for( i = 0; i < iVolumeSize; i++) {
		if( pVolume[i] > THRESHOLD ) // for chest
		{	
			iIndexX = (i%(iSliceSize))%iVolX;
			iIndexY = (i%(iSliceSize))/iVolX;
			iIndexZ = i/(iSliceSize);

			if(ptMaxBound.x < iIndexX)
				ptMaxBound.x = iIndexX;
			if(ptMaxBound.y < iIndexY)
				ptMaxBound.y = iIndexY;
			if(ptMaxBound.z < iIndexZ)
				ptMaxBound.z = iIndexZ;
			if(ptMinBound.x > iIndexX)
				ptMinBound.x = iIndexX;
			if(ptMinBound.y > iIndexY)
				ptMinBound.y = iIndexY;
			if(ptMinBound.z > iIndexZ)
				ptMinBound.z = iIndexZ;
		}
	}
	ControlProgress(pWndProgress);

	int nX = ptMaxBound.x - ptMinBound.x + 1;
	int nY = ptMaxBound.y - ptMinBound.y + 1;
	iSliceNum = ptMaxBound.z - ptMinBound.z + 1;
	int iRemainderX = nX%4;
	int iRemainderY = nY%4;
	
	if( iRemainderX!=0 ) { 
		if(ptMaxBound.x + (4-iRemainderX) < iVolX) {
			iVolX = nX + (4-iRemainderX);
			ptMaxBound.x += (4-iRemainderX);
		}
		else if( ptMinBound.x - (4-iRemainderX) >= 0 ) {
			ptMinBound.x -=(4-iRemainderX);
			iVolX = nX + (4-iRemainderX);
		}
	} else iVolX = nX;

	if( iRemainderY!=0 ) { 
		if(ptMaxBound.y + (4-iRemainderY) < iVolY) {
			iVolY = nY + (4-iRemainderY);
			ptMaxBound.y += (4-iRemainderY);
		}
		else if( ptMinBound.y - (4-iRemainderY) >= 0 ) {
			ptMinBound.y -=(4-iRemainderY);
			iVolY = nY + (4-iRemainderY);
		}
	} else iVolY = nY;

	m_pLabel = new unsigned char[iVolX*iVolY*iSliceNum];
	memset(m_pLabel, 0x00, sizeof(unsigned char)*iVolX*iVolY*iSliceNum);
	unsigned short*	pTempVolume = new unsigned short[iVolX*iVolY*iSliceNum];

	int iCount = 0;
	for(i=0; i<iVolumeSize; i++) {
		iIndexX = (i%iSliceSize)%iVolX;
		iIndexY = (i%iSliceSize)/iVolX;
		iIndexZ = i/iSliceSize;
		if(iIndexX >= ptMinBound.x && iIndexX <= ptMaxBound.x && 
			iIndexY >= ptMinBound.y && iIndexY <= ptMaxBound.y &&
			iIndexZ >= ptMinBound.z && iIndexZ <= ptMaxBound.z) 
			pTempVolume[iCount++] = pVolume[i];
	}
	iVolumeSize = iVolX*iVolY*iSliceNum;
	iSliceSize = iVolX*iVolY;
	// SRG threshold for background	
	thresholdmin = 0, thresholdmax = THRESHOLD;

	// 3. First Seeded Region Growing