const Error *Path::AddLine( uunit length )
{
   // Length must be a positive value.
   CML_ASSERT( length >= 0 );

   // Calculate the ending position based on the
   // current position and direction of travel
   Point<PATH_MAX_DIMENSIONS> p;
   p.setDim( GetDim() );

   p = posEnd;
   p[0] += cos(dirEnd) * length;
   if( GetDim() > 1 )
      p[1] += sin(dirEnd) * length;

   // add the line segment to my path
   LineSeg *seg = new LineSeg( posEnd, dirEnd, length );
   if( !seg ) return &PathError::Alloc;

   const Error *err = AddSegment( seg );
   if( err ) return err;

   // Update my ending position & direction
   posEnd = p;
   return 0;
}

void GTRProfileProcess::Create(int innsite, int indim)	{
	if (! rr)	{
		ProfileProcess::Create(innsite,indim);
		Nrr = GetDim() * (GetDim()-1) / 2;
		rr = new double[Nrr];
		SampleRR();
	}
}

// MPI: these two functions are responsible for broadcasting/receiving the current state of the parameter vector
// are model dependent
// should be implemented in .cpp file
void CodonMutSelSBDPPhyloProcess::SlaveUpdateParameters()	{

	// SlaveBroadcastTree();

	int i,j,L1,L2,ni,nd,nbranch = GetNbranch(),nnucrr = GetNnucrr(),nnucstat = 4;
	L1 = GetNmodeMax();
	L2 = GetDim();
	//nd = nbranch + nnucrr + nnucstat + L2 + L1*(L2+1); // check if these last terms are correct in this context...
	nd = nbranch + nnucrr + nnucstat + L1*L2 + GetDim() + 1;
	ni = 1 + ProfileProcess::GetNsite();
	int* ivector = new int[ni];
	double* dvector = new double[nd];
	MPI_Bcast(ivector,ni,MPI_INT,0,MPI_COMM_WORLD);
	MPI_Bcast(dvector,nd,MPI_DOUBLE,0,MPI_COMM_WORLD);
	int index = 0;
	for(i=0; i<nbranch; ++i) {
		blarray[i] = dvector[index];
		index++;
	}
	for(i=0; i<nnucrr; ++i) {
		nucrr[i] = dvector[index];
		index++;
	}
	for(i=0; i<nnucstat; ++i) {
		nucstat[i] = dvector[index];
		index++;
	}
	for(i=0; i<L1; ++i) {
		for(j=0; j<L2; ++j) {
			profile[i][j] = dvector[index];
			index++;
		}
	}
	for (int i=0; i<GetDim(); i++)	{
		dirweight[0][i] = dvector[index];
		index++;
	}
	kappa = dvector[index];
	index++;
	
	Ncomponent = ivector[0];
	for(i=0; i<ProfileProcess::GetNsite(); ++i) {
		SBDPProfileProcess::alloc[i] = ivector[1+i];
	}
	//GetBranchLengthsFromArray();
	delete[] dvector;
	delete[] ivector;
	// this one is really important
	// in those cases where new components have appeared, or some old ones have disappeared
	// during allocation move on the master node.
	// 
	// note that CreateMatrices() in fact creates only those that are not yet allocated
	// and also deletes those that are now obsolete
	// CreateMatrices();
	UpdateMatrices();
}

static void DrawRegistrationMarks( drawCmd_p d )
{
	long x, y, delta, divisor;
	coOrd p0, p1, qq, q0, q1;
	POS_T len;
	char msg[10];
	wFont_p fp;
	wFontSize_t fs;
	fp = wStandardFont( F_TIMES, FALSE, FALSE );
	if ( units==UNITS_METRIC ) {
		delta = 10;
		divisor = 100;
	} else {
		delta = 3;
		divisor = 12;
	}
	for ( x=delta; (POS_T)x<PutDim(mapD.size.x); x+=delta ) {
		qq.x = p0.x = p1.x = (POS_T)GetDim(x);
		p0.y = 0.0;
		p1.y = mapD.size.y;
		if (!ClipLine( &p0, &p1, d->orig, d->angle, d->size ))
			continue;
		for ( y=(long)(ceil(PutDim(p0.y)/delta)*delta); (POS_T)y<PutDim(p1.y); y+=delta ) {
			qq.y = (POS_T)GetDim(y);
			q0.x = q1.x = qq.x;
			if ( x%divisor == 0 && y%divisor == 0 ) {
				len = 0.25;
				fs = 12.0;
			} else {
				len = 0.125;
				fs = 8.0;
			}
			q0.y = qq.y-len;
			q1.y = qq.y+len;
			DrawLine( d, q0, q1, 0, wDrawColorBlack );
			q0.y = q1.y = qq.y;
			q0.x = qq.x-len;
			q1.x = qq.x+len;
			DrawLine( d, q0, q1, 0, wDrawColorBlack );
			q0.x = qq.x + len/4;;
			q0.y = qq.y + len/4;;
			if (units == UNITS_METRIC)
				sprintf( msg, "%0.1fm", (DOUBLE_T)x/100.0 );
			else
				sprintf( msg, "%ld\' %ld\"", x/12, x%12 );
			DrawString( d, q0, 0.0, msg, fp, fs, wDrawColorBlack );
			q0.y = qq.y - len*3/4;
			if (units == UNITS_METRIC)
				sprintf( msg, "%0.1fm", (DOUBLE_T)y/100.0 );
			else
				sprintf( msg, "%ld\' %ld\"", y/12, y%12 );
			DrawString( d, q0, 0.0, msg, fp, fs, wDrawColorBlack );
		}
	}
}

void CSceneNodeAABB::Render() const{
	if(!GetSingleton<CProgram>()->IsInWireFrameMode())
	{
		//glPushAttrib(GL_ALL_ATTRIB_BITS);
		////glBegin();
		//glEnable(GL_MAP2_VERTEX_3);
		//vector3 r[4] = { 
		//	m_aabb.first , 
		//	vector3(m_aabb.first.x,m_aabb.second.y,m_aabb.first.z) ,
		//	vector3(m_aabb.second.x,m_aabb.first.y,m_aabb.first.z) , 
		//	vector3(m_aabb.second.x,m_aabb.second.y,m_aabb.first.z)
		//};
		//glMap2f( GL_MAP2_VERTEX_3,0.0f,1.0f,6,2,0.0f,1.0f,3,2,(GLfloat*)r);
		//glMapGrid2f(10,0,1,20,0,1);
		//glEvalMesh2(GL_FILL,0,10,0,20);
		////glEnd();
		//glPopAttrib();
		vector3 dim(GetDim(m_aabb));
		vector3 dx(dim.x,0,0);
		vector3 dy(0,dim.y,0);
		vector3 dz(0,0,dim.z);
		unsigned int nx(dim.x*0.02+1);
		unsigned int ny(dim.y*0.02+1);
		unsigned int nz(dim.z*0.02+1);
		DrawFace(m_aabb.first,dy,dx,ny,nx);
		DrawFace(m_aabb.first,dx,dz,nx,nz);
		DrawFace(m_aabb.first,dz,dy,nz,ny);
		DrawFace(m_aabb.second,-dx,-dy,nx,ny);
		DrawFace(m_aabb.second,-dy,-dz,ny,nz);
		DrawFace(m_aabb.second,-dz,-dx,nz,nx);
		return;
	}
	glPushAttrib(GL_TEXTURE_BIT);
	glPushMatrix();
	vector3 center(GetCenter(m_aabb)),dim(GetDim(m_aabb));
	glTranslatef(center.x,center.y,center.z);
	{
		GLenum coord[2]={GL_S,GL_T};
		for(int iCoord(0); iCoord<2; ++iCoord)
		{
			float f[4];
			glGetTexGenfv(coord[iCoord],GL_OBJECT_PLANE,f);
			f[0]*=dim.x*0.001f;
			f[1]*=dim.y*0.001f;
			f[2]*=dim.z*0.001f;
			//f[2]*=1000;
			glTexGenfv(coord[iCoord],GL_OBJECT_PLANE,f);
		}
	}
	glScalef(dim.x,dim.y,dim.z);
	glutSolidCube(1.0);
	glPopMatrix();
	glPopAttrib();
}

void RASCATFiniteGammaPhyloProcess::SlaveUpdateParameters()	{
	int i,j,L1,L2,ni,nd,nbranch = GetNbranch();
	L1 = GetNmodeMax();
	L2 = GetDim();
	nd = 3 + nbranch + L1*L2 + GetDim() + 1;
	ni = 1 + GetNsite();
	int* ivector = new int[ni];
	double* dvector = new double[nd];
	MPI_Bcast(ivector,ni,MPI_INT,0,MPI_COMM_WORLD);
	MPI_Bcast(dvector,nd,MPI_DOUBLE,0,MPI_COMM_WORLD);

	int index = 0;
	SetAlpha(dvector[index]);
	index++;
	branchalpha = dvector[index];
	index++;
	branchbeta = dvector[index];
	index++;
	for(i=0; i<nbranch; ++i) {
		blarray[i] = dvector[index];
		index++;
	}

	for(i=0; i<L1; ++i) {
		for(j=0; j<L2; ++j) {
			profile[i][j] = dvector[index];
			index++;
		}
	}
	for (int i=0; i<GetDim(); i++)	{
		dirweight[i] = dvector[index];
		index++;
	}

	Ncomponent = ivector[0];
	for(i=0; i<GetNsite(); ++i) {
		FiniteProfileProcess::alloc[i] = ivector[1+i];
	}

	MPI_Bcast(weight,GetNcomponent(),MPI_DOUBLE,0,MPI_COMM_WORLD);

	UpdateZip();

	delete[] dvector;
	delete[] ivector;

	UpdateZip();
}

CComSafeArray::CComSafeArray(LPCVARIANT pSrc)
{
	_SafeArrayInit(this);
	*this			= pSrc;
	m_dwDims		= GetDim();
	m_dwElementSize = GetElemSize();
}

/* Market::GetSpace
 * Takes a Market, specifically its D_States, and constructs the array of possible ending states.
 * JSpace represents D_State[0] x D_State[1] x ... x D_State[n-1] possible outcomes
 */
Space Market::GetSpace(bool Verbose) const
{
	std::vector<uint32_t> Dim = GetDim();
	if (Verbose) {
		printf("Market Dimention(s):\n");
		for(uint32_t i=0; i < D_State.size(); i++) {
			const char *str = (i+1==D_State.size())? "": ",";
			printf("%d%s", (uint32_t)D_State[i].size(),str);
		}
		printf("\n");
	}
	Space S;
	S.D_State = D_State;
	for(uint32_t i=0; i < D_State.size(); i++) {
		std::vector<std::string> Names;
		for(uint32_t j=0; j < D_State[i].size()+1; j++) {
			char name[16];
			sprintf(name, "d%u.%s", i+1, ((j==0)? "No": "Yes"));
			Names.push_back(std::string(name));
		}
		S.Names.push_back(Names);
	}
	uint32_t MaxN = 1;
	for(uint32_t i=0; i < D_State.size(); i++)
		MaxN *= D_State[i].size() + 1;
	for(uint32_t i=0; i < MaxN; i++)
		S.data.push_back(i+1);
	return S;
}

void TextRender::SetExpanse(const Vec2& _expanse, bool _onlyScaleDown)
{
    int len = (int)m_text.length();
    int lineCount = 1;
    int charCount = 0;
    int charCountMax = 0;
    for (int i = 0; i < len; i++)
    {
        charCount++;
        if (m_text[i] == '\n')
        {
            if (charCountMax < charCount)
                charCountMax = charCount;
            charCount = 0;
            lineCount++;
        }
    }

    //if the text contains no linebreaks
    if (!charCountMax) charCountMax = charCount;

    Vec2 captionDim = GetDim();

    // in case that the text is to large in any direction scale it down
    if (!_onlyScaleDown || (captionDim[0] * charCountMax * GetDefaultSize() >= _expanse[0] || captionDim[1] * lineCount >= _expanse[1]))
    {
        SetDefaultSize(min((float)(_expanse[0] / (captionDim[0] * charCountMax)),
                           (float)(_expanse[1] / (captionDim[1] * lineCount))));
    }

}

void TMultinomial::AddFtr(const TStr& Str, TIntFltKdV& SpV, int& Offset) const {
    const int FtrId = FtrGen.GetFtr(Str);
    if (FtrId != -1) {
        SpV.Add(TIntFltKd(Offset + FtrId, 1.0));
    }
    Offset += GetDim();
}

void TDateWnd::AddFtr(const TTm& Val, TIntFltKdV& SpV, int& Offset) const {
    const int Ftr = GetFtr(Val);
    for (int FtrN = 0; FtrN < WndSize; FtrN++) {
        SpV.Add(TIntFltKd(Offset + Ftr + FtrN, Wgt));
    }
    Offset += GetDim();
}

void PoissonSubstitutionProcess::CreateZip()	{
	if (! zipstat)	{
		zipstat = new double*[GetNsite()];
		for (int i=0; i<GetNsite(); i++)	{
			zipstat[i] = new double[GetDim()];
		}
	}
}

void PartitionGTRMixtureProfileProcess::CreateMatrix(int k)	{
	if (matrixarray[k])	{
		cerr << "error matrixarray is not 0\n";
		cerr << k << '\t' << GetNmodeMax() << '\t' << matrixarray[k]->GetNstate() << '\n';
		exit(1);
	}
	matrixarray[k] = new GTRSubMatrix(GetDim(),rr[k],profile[k],false);
}

int PoissonSubstitutionProcess::GetRandomStateFromZip(int site, int zipstate)	{
	int truestate = 0;
	if ((GetZipSize(site) != GetOrbitSize(site)) && (zipstate == GetOrbitSize(site)))	{
		double v = rnd::GetRandom().Uniform();
		double u = zipstat[site][GetOrbitSize(site)] * v;
		double total = 0;
		double* pi = GetProfile(site);

		int choose = -1;
		while ((choose < GetDim()) && (total < u))	{
			choose ++;
			if (choose == GetDim())	{
				cerr << "error in getstatefromzip\n";
				cerr << choose << '\t' << GetDim() << '\n';
				cerr << total << '\n';
				cerr << v << '\t' << zipstat[site][GetOrbitSize(site)] << '\t' << u << '\n';
				cerr << total - zipstat[site][GetOrbitSize(site)] << '\n';
				double newtotal = 0;
				cerr << '\n';
				for (int k=0; k<GetDim(); k++)	{
					cerr << pi[k] << '\t';
					cerr << InOrbit(site,k) << '\n';
					if (! InOrbit(site,k))	{
						newtotal += pi[k];
					}
				}
				cerr << '\n';
				for (int k=0; k<=GetOrbitSize(site); k++)	{
					cerr << zipstat[site][k] << '\n';
				}
				cerr << "new total : " << newtotal << '\t' << newtotal - total << '\t' << total - choose << '\n';
				exit(1);
			}
			if (!InOrbit(site,choose))	{
				total += pi[choose];
			}
		}
		truestate = choose;
	}
	else	{
		truestate = GetStateFromZip(site,zipstate);
	}
	return truestate;
}

void CComSafeArray::ResizeOneDim(DWORD dwElements)
{
	ATLASSERT(GetDim() == 1);

	SAFEARRAYBOUND rgsabound;

	rgsabound.cElements = dwElements;
	rgsabound.lLbound	= 0;

	Redim(&rgsabound);
}

void TBagOfWords::AddFtr(const TStrV& TokenStrV, TFltV& FullV, int& Offset) const {
    // create sparse vector
    TIntFltKdV ValSpV; AddFtr(TokenStrV, ValSpV);
    // add to the full feature vector and increase offset count
    for (int ValSpN = 0; ValSpN < ValSpV.Len(); ValSpN++) {
        const TIntFltKd& ValSp = ValSpV[ValSpN];
        FullV[Offset + ValSp.Key] = ValSp.Dat;
    }
    // increase the offset by the dimension
    Offset += GetDim();    
}