rxxpowellmethod.cpp

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

// RxxPowellMethod.cpp: implementation of the RxxPowellMethod class.
//
//////////////////////////////////////////////////////////////////////

#include "stdafx.h"
#include "RxxPowellMethod.h"
#include "Matrix4D.h"
#include "ProgressWnd.h"

#include "JJDistanceMap.h"

#include <math.h>

#define DEGREETORADIAN		3.141592/180
#define RADIANTODEGREE		180/3.141592
//////////////////////////////////////////////////////////////////////
// Construction/Destruction
//////////////////////////////////////////////////////////////////////
RxxPowellMethod::RxxPowellMethod()
{
	m_lfMaxCC = 0.0f;	

	for (int i = 0; i < 1000; i++)
		m_arrNoduleMapping[i] = i;
}

RxxPowellMethod::RxxPowellMethod(JJDistanceMap* pRefDistanceMap, JJBinaryVolume* pFloatBinaryVolume,
								 short* pListX, short* pListY, short* pListZ, unsigned long nNum, double BoundingBoxX, double BoundingBoxY, double BoundingBoxZ)
{
// translation initialization
	m_lfTransX			= 0.0;
	m_lfTransY			= 0.0;
	m_lfTransZ			= 0.0;

// rotation initialization
	m_lfRadX			= 0.0;		
	m_lfRadY			= 0.0;		
	m_lfRadZ			= 0.0;
	m_lfRotX			= 0.0;
	m_lfRotY			= 0.0;
	m_lfRotZ			= 0.0;	

	m_lfOptimizedScaleX = g_ScaleX;
	m_lfOptimizedScaleY = g_ScaleY;
	m_lfOptimizedScaleZ = g_ScaleZ;

// shear initialization
	m_lfShearX = 0.0;
	m_lfShearY = 0.0;
	m_lfShearZ = 0.0;

// bounding box movement
	m_lfBoundingBoxX	= BoundingBoxX;
	m_lfBoundingBoxY	= BoundingBoxY;
	m_lfBoundingBoxZ	= BoundingBoxZ;

// distance map & volume initialization
	m_pRefDistanceMap			= pRefDistanceMap;
	m_pFloatBinaryVolume		= pFloatBinaryVolume;	

// edge list initialization
	m_pEdgeListX = pListX;
	m_pEdgeListY = pListY;
	m_pEdgeListZ = pListZ;
	m_nNumEdge = nNum;
	
	m_lfMaxCC = 10000.0f;	

	for (int i = 0; i < 1000; i++)
		m_arrNoduleMapping[i] = i;
}

RxxPowellMethod::RxxPowellMethod(_MODAL_INFO* pFloat, _MODAL_INFO* pRefer, double lfGlobalScaleX, double lfGlobalScaleY, double lfGlobalScaleZ, double lfGlobalRotateX, double lfGlobalRotateY, double lfGlobalRotateZ, int nLabelIdx, UINT nInterpID)
{

	m_lfGlobalScaleX = lfGlobalScaleX;
	m_lfGlobalScaleY = lfGlobalScaleY;
	m_lfGlobalScaleZ = lfGlobalScaleZ;	

	m_plfJointHisto = NULL;
	m_lfMaxMI = 0.0;	

//	m_lfOptimizedScaleX = 0.0;
//	m_lfOptimizedScaleY = 0.0;
//	m_lfOptimizedScaleZ = 0.0;

	m_lfOptimizedScaleX = lfGlobalScaleX;
	m_lfOptimizedScaleY = lfGlobalScaleY;
	m_lfOptimizedScaleZ = lfGlobalScaleZ;

	m_lfTranslateXofFlt = 0.0;
	m_lfTranslateYofFlt = 0.0;
	m_lfTranslateZofFlt = 0.0;

	m_lfThetaX = 0.0;		
	m_lfThetaY = 0.0;		
	m_lfThetaZ = 0.0;

	m_lfRotX = 0.0;
	m_lfRotY = 0.0;
	m_lfRotZ = 0.0;

//	m_ptrCenter = CPoint(-1,-1);
	m_lfRegionInPointCnt = 0.0;

//	m_szFloatImg = CSize(512,512);
//	m_szRefImg = CSize(2048,1216);

	m_lfGlobalRotateX = lfGlobalRotateX*CIRCLEDEGREE;
	m_lfGlobalRotateY = lfGlobalRotateY*CIRCLEDEGREE;
	m_lfGlobalRotateZ = lfGlobalRotateZ*CIRCLEDEGREE;

	m_pFloat = pFloat;
	m_pRefer = pRefer;

	m_lfGradientInfo = 0.0;

	m_nLabelIdx = nLabelIdx;

	m_nInterpID = nInterpID;

	for (int i = 0; i < 1000; i++)
		m_arrNoduleMapping[i] = i;
}

RxxPowellMethod::~RxxPowellMethod()
{

}

void RxxPowellMethod::FindTransformVector()
{		
	DWORD CurrentTime = GetTickCount();	
	 
	double lfStartptMI = RxxMutualInfoAgent::ComputeMI(m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ, 0, 0, 0, 0.0, 0.0, 0.0);
	RxxMutualInfoAgent::KeepParameter(m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ, 0, 0, 0, 0.0, 0.0, 0.0);

	double lfIncS = 0.0;
	double lfIncT = 0.2;
	double lfIncR;
	int nIteration = 3;	
	double lfDetailNum = 1.0;	

	for(int i=0 ; i<nIteration ; i++){		
		
		lfIncR = 0.1*lfDetailNum*CIRCLEDEGREE;

		lfStartptMI = IncreaseTransXvector(lfStartptMI, lfIncT);
		lfStartptMI = DecreaseTransXvector(lfStartptMI, lfIncT);

		lfStartptMI = IncreaseTransYvector(lfStartptMI, lfIncT);
		lfStartptMI = DecreaseTransYvector(lfStartptMI, lfIncT);		

		lfStartptMI = CWRotZvector(lfStartptMI, lfIncR);
		lfStartptMI = CCWRotZvector(lfStartptMI, lfIncR);		

		lfStartptMI = CWRotXvector(lfStartptMI, lfIncR);
		lfStartptMI = CCWRotXvector(lfStartptMI, lfIncR);

		lfStartptMI = CWRotYvector(lfStartptMI, lfIncR);
		lfStartptMI = CCWRotYvector(lfStartptMI, lfIncR);

		lfStartptMI = IncreaseTransZvector(lfStartptMI, lfIncT);
		lfStartptMI = DecreaseTransZvector(lfStartptMI, lfIncT);		

		lfDetailNum *= 0.1;
	}
	
	DWORD ElapsedTime = GetTickCount () - CurrentTime;		

	double lfTransXMM = m_lfTranslateXofFlt*m_pFloat->lfVoxelX;
	double lfTransYMM = m_lfTranslateYofFlt*m_pFloat->lfVoxelY;
	double lfTransZMM = m_lfTranslateZofFlt*m_pFloat->lfVoxelZ;	
		
	CString str;
	str.Format(_T("S = (%0.3lf, %0.3lf, %0.3lf)\r\nT = (%0.6lf, %0.6lf, %0.6lf)\r\nT = (%0.3lfmm, %0.3lfmm, %0.3lfmm)\r\nR(Rad) = (%0.3lf, %0.3lf, %0.3lf)\r\nR(Degree) = (%0.3lf, %0.3lf, %0.3lf)\r\n\r\nMaxMI=%lf\r\n\r\nTotal Processing Time =%dms"), 
				m_lfOptimizedScaleX, m_lfOptimizedScaleY, m_lfOptimizedScaleZ, 
				m_lfTranslateXofFlt, m_lfTranslateYofFlt, m_lfTranslateZofFlt,
				lfTransXMM, lfTransYMM, lfTransZMM,
				m_lfThetaX, m_lfThetaY, m_lfThetaZ, 
				m_lfRotX, m_lfRotY, m_lfRotZ, lfStartptMI, ElapsedTime);
	AfxMessageBox(str);
}

double RxxPowellMethod::IncreaseTransXvector(double lfPresentMI, double lfIncT)
{	
	double MI =0.0;
	int nIsTranslateCOUNTERState = EMPTYCOUNTER;
	for(double lfTransX=m_lfTranslateXofFlt+lfIncT ; ; lfTransX +=lfIncT){
		MI = RxxMutualInfoAgent::ComputeMI(m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ, lfTransX, m_lfTranslateYofFlt, m_lfTranslateZofFlt, m_lfThetaX, m_lfThetaY, m_lfThetaZ);		
//		TRACE(_T("Scale=(%lf,%lf,%lf), Translate=(%lf,%lf,%lf), Theta=(%lf,%lf,%lf), MI=%d\n"),m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ,lfTransX,m_lfTranslateYofFlt,m_lfTranslateZofFlt, m_lfThetaX, m_lfThetaY, m_lfThetaZ, MI);
		if(MI>lfPresentMI){
			RxxMutualInfoAgent::KeepParameter(m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ, lfTransX, m_lfTranslateYofFlt, m_lfTranslateZofFlt, m_lfThetaX, m_lfThetaY, m_lfThetaZ);
			lfPresentMI = MI;
			nIsTranslateCOUNTERState = EMPTYCOUNTER;
		}
		else if(nIsTranslateCOUNTERState!=FULLCOUNTER)	nIsTranslateCOUNTERState++;		
		else	break;		
	}
	return lfPresentMI;
}

double RxxPowellMethod::IncreaseTransYvector(double lfPresentMI, double lfIncT)
{
	double MI =0.0;
	int nIsTranslateCOUNTERState = EMPTYCOUNTER;
	for(double lfTransY=m_lfTranslateYofFlt+lfIncT ; ; lfTransY +=lfIncT){
		MI = RxxMutualInfoAgent::ComputeMI(m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ, m_lfTranslateXofFlt, lfTransY, m_lfTranslateZofFlt, m_lfThetaX, m_lfThetaY, m_lfThetaZ);
//		TRACE(_T("Scale=(%lf,%lf,%lf), Translate=(%lf,%lf,%lf), Theta=(%lf,%lf,%lf), MI=%d\n"),m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ,m_lfTranslateXofFlt,lfTransY,m_lfTranslateZofFlt, m_lfThetaX, m_lfThetaY, m_lfThetaZ, MI);
		if(MI>lfPresentMI){
			RxxMutualInfoAgent::KeepParameter(m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ, m_lfTranslateXofFlt, lfTransY, m_lfTranslateZofFlt, m_lfThetaX, m_lfThetaY, m_lfThetaZ);
			lfPresentMI = MI;
			nIsTranslateCOUNTERState = EMPTYCOUNTER;
		}
		else if(nIsTranslateCOUNTERState!=FULLCOUNTER)	nIsTranslateCOUNTERState++;		
		else	break;		
	}
	return lfPresentMI;
}

double RxxPowellMethod::IncreaseTransZvector(double lfPresentMI, double lfIncT)
{
	double MI =0.0;
	int nIsTranslateCOUNTERState = EMPTYCOUNTER;
	for(double lfTransZ=m_lfTranslateZofFlt+lfIncT ; ; lfTransZ +=lfIncT){
		MI = RxxMutualInfoAgent::ComputeMI(m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ, m_lfTranslateXofFlt, m_lfTranslateYofFlt, lfTransZ, m_lfThetaX, m_lfThetaY, m_lfThetaZ);		

//		TRACE(_T("Scale=(%lf,%lf,%lf), Translate=(%lf,%lf,%lf), Theta=(%lf,%lf,%lf), MI=%lf\n"),m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ,m_lfTranslateXofFlt, m_lfTranslateYofFlt, lfTransZ, m_lfThetaX, m_lfThetaY, m_lfThetaZ, MI);
		if(MI>lfPresentMI){
			RxxMutualInfoAgent::KeepParameter(m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ, m_lfTranslateXofFlt, m_lfTranslateYofFlt, lfTransZ, m_lfThetaX, m_lfThetaY, m_lfThetaZ);
			lfPresentMI = MI;
			nIsTranslateCOUNTERState = EMPTYCOUNTER;
		}
		else if(nIsTranslateCOUNTERState!=FULLCOUNTER)	nIsTranslateCOUNTERState++;		
		else	break;		
	}
	return lfPresentMI;
}

double RxxPowellMethod::DecreaseTransXvector(double lfPresentMI, double lfIncT)
{
	double MI =0.0;
	int nIsTranslateCOUNTERState = EMPTYCOUNTER;
	for(double lfTransX=m_lfTranslateXofFlt-lfIncT ; ; lfTransX -=lfIncT){
		MI = RxxMutualInfoAgent::ComputeMI(m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ, lfTransX, m_lfTranslateYofFlt, m_lfTranslateZofFlt, m_lfThetaX, m_lfThetaY, m_lfThetaZ);
//		TRACE(_T("Scale=(%lf,%lf,%lf), Translate=(%lf,%lf,%lf), Theta=(%lf,%lf,%lf), MI=%lf\n"),m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ,lfTransX,m_lfTranslateYofFlt,m_lfTranslateZofFlt, m_lfThetaX, m_lfThetaY, m_lfThetaZ, MI);
		if(MI>lfPresentMI){
			RxxMutualInfoAgent::KeepParameter(m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ, lfTransX, m_lfTranslateYofFlt, m_lfTranslateZofFlt, m_lfThetaX, m_lfThetaY, m_lfThetaZ);
			lfPresentMI = MI;
			nIsTranslateCOUNTERState = EMPTYCOUNTER;
		}
		else if(nIsTranslateCOUNTERState!=FULLCOUNTER)	nIsTranslateCOUNTERState++;		
		else	break;		
	}
	return lfPresentMI;
}

double RxxPowellMethod::DecreaseTransYvector(double lfPresentMI, double lfIncT)
{
	double MI =0.0;
	int nIsTranslateCOUNTERState = EMPTYCOUNTER;
	for(double lfTransY=m_lfTranslateYofFlt-lfIncT ; ; lfTransY -=lfIncT){
		MI = RxxMutualInfoAgent::ComputeMI(m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ, m_lfTranslateXofFlt, lfTransY, m_lfTranslateZofFlt, m_lfThetaX, m_lfThetaY, m_lfThetaZ);
//		TRACE(_T("Scale=(%lf,%lf,%lf), Translate=(%lf,%lf,%lf), Theta=(%lf,%lf,%lf), MI=%lf\n"),m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ,m_lfTranslateXofFlt,lfTransY,m_lfTranslateZofFlt, m_lfThetaX, m_lfThetaY, m_lfThetaZ, MI);
		if(MI>lfPresentMI){
			RxxMutualInfoAgent::KeepParameter(m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ, m_lfTranslateXofFlt, lfTransY, m_lfTranslateZofFlt, m_lfThetaX, m_lfThetaY, m_lfThetaZ);
			lfPresentMI = MI;
			nIsTranslateCOUNTERState = EMPTYCOUNTER;
		}
		else if(nIsTranslateCOUNTERState!=FULLCOUNTER)	nIsTranslateCOUNTERState++;		
		else	break;		
	}
	return lfPresentMI;

}

double RxxPowellMethod::DecreaseTransZvector(double lfPresentMI, double lfIncT)
{
	double MI =0.0;
	int nIsTranslateCOUNTERState = EMPTYCOUNTER;
	for(double lfTransZ=m_lfTranslateZofFlt-lfIncT ; ; lfTransZ -=lfIncT){
		MI = RxxMutualInfoAgent::ComputeMI(m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ, m_lfTranslateXofFlt, m_lfTranslateYofFlt, lfTransZ, m_lfThetaX, m_lfThetaY, m_lfThetaZ);
//		TRACE(_T("Scale=(%lf,%lf,%lf), Translate=(%lf,%lf,%lf), Theta=(%lf,%lf,%lf), MI=%lf\n"),m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ,m_lfTranslateXofFlt, m_lfTranslateYofFlt, lfTransZ, m_lfThetaX, m_lfThetaY, m_lfThetaZ, MI);
		if(MI>lfPresentMI){
			RxxMutualInfoAgent::KeepParameter(m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ, m_lfTranslateXofFlt, m_lfTranslateYofFlt, lfTransZ, m_lfThetaX, m_lfThetaY, m_lfThetaZ);
			lfPresentMI = MI;
			nIsTranslateCOUNTERState = EMPTYCOUNTER;
		}
		else if(nIsTranslateCOUNTERState!=FULLCOUNTER)	nIsTranslateCOUNTERState++;		
		else	break;		
	}
	return lfPresentMI;
}

double RxxPowellMethod::CWRotXvector(double lfPresentMI, double lfIncR)
{
	double MI =0.0;
	int nIsThetaCOUNTERState = EMPTYCOUNTER;
	for(double lfThetaX=m_lfThetaX+lfIncR ; ; lfThetaX +=lfIncR){
		MI = RxxMutualInfoAgent::ComputeMI(m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ, m_lfTranslateXofFlt, m_lfTranslateYofFlt, m_lfTranslateZofFlt, lfThetaX, m_lfThetaY, m_lfThetaZ);
//		TRACE(_T("Scale=(%lf,%lf,%lf), Translate=(%lf,%lf,%lf), Theta=(%lf,%lf,%lf), MI=%lf\n"),m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ,m_lfTranslateXofFlt, m_lfTranslateYofFlt, m_lfTranslateZofFlt, lfThetaX, m_lfThetaY, m_lfThetaZ, MI);
		if(MI>lfPresentMI){
			RxxMutualInfoAgent::KeepParameter(m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ, m_lfTranslateXofFlt, m_lfTranslateYofFlt, m_lfTranslateZofFlt, lfThetaX, m_lfThetaY, m_lfThetaZ);
			lfPresentMI = MI;
			nIsThetaCOUNTERState = EMPTYCOUNTER;
		}
		else if(nIsThetaCOUNTERState!=FULLCOUNTER)	nIsThetaCOUNTERState++;		
		else	break;		
	}
	return lfPresentMI;
}

double RxxPowellMethod::CWRotYvector(double lfPresentMI, double lfIncR)
{
	double MI =0.0;
	int nIsThetaCOUNTERState = EMPTYCOUNTER;
	for(double lfThetaY=m_lfThetaY+lfIncR ; ; lfThetaY +=lfIncR){
		MI = RxxMutualInfoAgent::ComputeMI(m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ, m_lfTranslateXofFlt, m_lfTranslateYofFlt, m_lfTranslateZofFlt, m_lfThetaX, lfThetaY, m_lfThetaZ);
//		TRACE(_T("Scale=(%lf,%lf,%lf), Translate=(%lf,%lf,%lf), Theta=(%lf,%lf,%lf), MI=%lf\n"),m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ,m_lfTranslateXofFlt, m_lfTranslateYofFlt, m_lfTranslateZofFlt, m_lfThetaX, lfThetaY, m_lfThetaZ, MI);
		if(MI>lfPresentMI){
			RxxMutualInfoAgent::KeepParameter(m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ, m_lfTranslateXofFlt, m_lfTranslateYofFlt, m_lfTranslateZofFlt, m_lfThetaX, lfThetaY, m_lfThetaZ);
			lfPresentMI = MI;
			nIsThetaCOUNTERState = EMPTYCOUNTER;
		}
		else if(nIsThetaCOUNTERState!=FULLCOUNTER)	nIsThetaCOUNTERState++;		
		else	break;		
	}
	return lfPresentMI;
}

double RxxPowellMethod::CWRotZvector(double lfPresentMI, double lfIncR)
{
	double MI =0.0;
	int nIsThetaCOUNTERState = EMPTYCOUNTER;
	for(double lfThetaZ=m_lfThetaZ+lfIncR ; ; lfThetaZ +=lfIncR){
		MI = RxxMutualInfoAgent::ComputeMI(m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ, m_lfTranslateXofFlt, m_lfTranslateYofFlt, m_lfTranslateZofFlt, m_lfThetaX, m_lfThetaY, lfThetaZ);
//		TRACE(_T("Scale=(%lf,%lf,%lf), Translate=(%lf,%lf,%lf), Theta=(%lf,%lf,%lf), MI=%lf\n"),m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ,m_lfTranslateXofFlt, m_lfTranslateYofFlt, m_lfTranslateZofFlt, m_lfThetaX, m_lfThetaY, lfThetaZ, MI);
		if(MI>lfPresentMI){
			RxxMutualInfoAgent::KeepParameter(m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ, m_lfTranslateXofFlt, m_lfTranslateYofFlt, m_lfTranslateZofFlt, m_lfThetaX, m_lfThetaY, lfThetaZ);
			lfPresentMI = MI;
			nIsThetaCOUNTERState = EMPTYCOUNTER;
		}
		else if(nIsThetaCOUNTERState!=FULLCOUNTER)	nIsThetaCOUNTERState++;		
		else	break;		
	}
	return lfPresentMI;
}

double RxxPowellMethod::CCWRotXvector(double lfPresentMI, double lfIncR)
{
	double MI =0.0;
	int nIsThetaCOUNTERState = EMPTYCOUNTER;
	for(double lfThetaX=m_lfThetaX-lfIncR ; ; lfThetaX -=lfIncR){
		MI = RxxMutualInfoAgent::ComputeMI(m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ, m_lfTranslateXofFlt, m_lfTranslateYofFlt, m_lfTranslateZofFlt, lfThetaX, m_lfThetaY, m_lfThetaZ);
//		TRACE(_T("Scale=(%lf,%lf,%lf), Translate=(%lf,%lf,%lf), Theta=(%lf,%lf,%lf), MI=%lf\n"),m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ,m_lfTranslateXofFlt, m_lfTranslateYofFlt, m_lfTranslateZofFlt, lfThetaX, m_lfThetaY, m_lfThetaZ, MI);
		if(MI>lfPresentMI){
			RxxMutualInfoAgent::KeepParameter(m_lfGlobalScaleX, m_lfGlobalScaleY, m_lfGlobalScaleZ, m_lfTranslateXofFlt, m_lfTranslateYofFlt, m_lfTranslateZofFlt, lfThetaX, m_lfThetaY, m_lfThetaZ);
			lfPresentMI = MI;
			nIsThetaCOUNTERState = EMPTYCOUNTER;
		}
		else if(nIsThetaCOUNTERState!=FULLCOUNTER)	nIsThetaCOUNTERState++;		
		else	break;		
	}
	return lfPresentMI;
}

double RxxPowellMethod::CCWRotYvector(double lfPresentMI, double lfIncR)
{
	double MI =0.0;
	int nIsThetaCOUNTERState = EMPTYCOUNTER;