// Gene.cpp: implementation of the CGene class.
//
//////////////////////////////////////////////////////////////////////

#include "stdafx.h"
#include "RegCluster.h"
#include "Gene.h"

#include <math.h>

#ifdef _DEBUG
#undef THIS_FILE
static char THIS_FILE[]=__FILE__;
#define new DEBUG_NEW
#endif

//////////////////////////////////////////////////////////////////////
// Construction/Destruction
//////////////////////////////////////////////////////////////////////

CGene::CGene()
{
	SetUp();
}

CGene::~CGene()
{
	Destruct();
}

void CGene::SortCondVal()
{
	QuickSort(0, m_nCond - 1);
}

void CGene::ChainCond()
{
/*	
	for(int i = 0, j = i + 1; j < m_nCond; ++j)
	{
		if(m_pCondVal[j] == m_pCondVal[i])
		{
			m_pChain[j] = m_pChain[i];
			i = j;
		}
		else if((m_pCondVal[j] - m_pCondVal[i]) > m_dThreshold)
		{
			m_pChain[j] = i;
			i = j;
		}
	}
*/
	for(int i = 0, j = i + 1; j < m_nCond; ++j)
	{
		if((m_pCondVal[j] - m_pCondVal[i]) > m_dThreshold)
		{
			m_pChain[j] = i;
			i = j;
		}
		else
		{
			m_pChain[j] = m_pChain[i];
		}
	}
/*	
	for(int i = 1, j = i + 1; j < m_nCond; ++j)
	{
		if(m_pCondVal[j] == m_pCondVal[i])
		{
			m_pChain[j] = m_pChain[i];
			i = j;
		}
		else if((m_pCondVal[j] - m_pCondVal[i]) > m_dThreshold)
		{
			m_pChain[j] = i;
			i = j;
		}
	}
*/
}

void CGene::QuickSort(int nStart, int nEnd)
{
	if(nStart < nEnd)
	{
		int q = Partition(nStart, nEnd);
		QuickSort(nStart, q - 1);
		QuickSort(q + 1, nEnd);
	}
}

int CGene::Partition(int nStart, int nEnd)
{
	double val = m_pCondVal[nEnd];
	int i = nStart - 1;
	for(int j = nStart; j < nEnd ; ++j)	// here use '-2' as the index interval is 0 ~ n-1
	{
		if(m_pCondVal[j] < val)
		{
			i += 1;
			// exchange gene expression value
			double tempVal = m_pCondVal[i];
			m_pCondVal[i] = m_pCondVal[j];
			m_pCondVal[j] = tempVal;
			// exchange gene ID accordingly
			int tempID = m_pID[i];
			m_pID[i] = m_pID[j];
			m_pID[j] = tempID;
		}
	}

	double tempVal = m_pCondVal[i+1];
	m_pCondVal[i+1] = m_pCondVal[nEnd];
	m_pCondVal[nEnd] = tempVal;
	int tempID = m_pID[i+1];
	m_pID[i+1] = m_pID[nEnd];
	m_pID[nEnd] = tempID;

	return i+1;
}

int CGene::GetCondID(int nIndex)
{
	return m_pID[nIndex];
}

int CGene::GetCondIndex(int nCond)
{
	for(int i = 0; i < m_nCond; ++i)
	{
		if(m_pID[i] == nCond)
		{
			return i;
		}
	}
	return -1;
}

void CGene::SetCondSize(int nCond)
{
	m_nCond = nCond;
}

void CGene::InitializeGene(double *pValue)
{
	Destruct();

	m_pCondVal = new double[m_nCond];
	m_pID = new int[m_nCond];
	m_pChain = new int[m_nCond];

	for(int i = 0; i < m_nCond; ++i)
	{
		m_pCondVal[i] = pValue[i];
		m_pID[i] = i;
		m_pChain[i] = -1;	// initialize the chain value with '-1'
	}
	CalculateThreshold();
}

void CGene::CalculateThreshold()
{
	double max = m_pCondVal[0];		//INT_MIN;
	double min = m_pCondVal[0];		//INT_MAX;
	for(int i = 1; i < m_nCond; ++i)
	{
		if(m_pCondVal[i] > max)
		{
			max = m_pCondVal[i];
		}
		else if(m_pCondVal[i] < min)
		{
			min = m_pCondVal[i];
		}
	}

	m_dThreshold = m_dGama * (max - min);	
}

void CGene::SetGama(double dGama)
{
	m_dGama = dGama;
}

void CGene::SetUp()
{
	m_dGama = 0;
	m_dThreshold = 0;
	m_nCond = 0;
	m_dCohert = 0;

	m_pChain = NULL;
	m_pCondVal = NULL;
	m_pID = NULL;
}

void CGene::Destruct()
{
	if(m_pChain != NULL)
	{
		delete [] m_pChain;
		m_pChain = NULL;
	}
	if(m_pID != NULL)
	{
		delete [] m_pID;
		m_pID = NULL;
	}
	if(m_pChain != NULL)
	{
		delete [] m_pChain;
		m_pChain = NULL;
	}
}

int CGene::GetPreCond(int nCond, int *&pPreCond)
{
	int nIndex = GetCondIndex(nCond);
	int nPreIndex = m_pChain[nIndex];
	if(nPreIndex == -1)
	{
		return -1;
	}
	int nCount = 0;
	for(int i = nPreIndex - 1; i >= 0; --i)
	{
		if(m_pChain[i] == m_pChain[i+1])
		{
			++nCount;
		}
		else
		{
			break;
		}
	}
	for(i = 0; i < nCount; ++i)
	{
		pPreCond[i] = GetCondID(nPreIndex-i);
	}

	return 0;
}

double CGene::CalculateCohert(int c1, int c2, int cm, int ci)
{
	int nIndexC1 = GetCondIndex(c1);
	int nIndexC2 = GetCondIndex(c2);
	int nIndexCm = GetCondIndex(cm);
	int nIndexCi = GetCondIndex(ci);
	double dCohert = 0;

	double dNumerator = m_pCondVal[nIndexCi] - m_pCondVal[nIndexCm];
	double dDenominator = m_pCondVal[nIndexC2] - m_pCondVal[nIndexC1];

	if(nIndexC1 < nIndexCi)
	{
		dCohert = dNumerator / dDenominator;
	}
	else
	{
		dCohert = dDenominator / dNumerator;
	}
	return dCohert;
}