/* This file is part of the source code for the following publication:
 * Assembling Self-Supporting Structures, Deuss et al., SIGGRAPH Asia 2014
 *
 * Copyright (C) 2014 Mario Deuss <mario.deuss@epfl.ch>
 *
 * This Source Code Form is subject to the terms of the Mozilla Public
 * License, v. 2.0. If a copy of the MPL was not distributed with this
 * file, You can obtain one at http://mozilla.org/MPL/2.0/. */

#include "PolygonPartition.h"
#include "../polypartition/polypartition.h"
#include <iostream>

namespace SelfAssembly{

int PolygonPartition::partition(int method, Eigen::Matrix3Xd &f, std::vector<std::vector<int> > &newF)
{
	//fit a plane to vertices of both faces of the interface
	Eigen::Vector3d mean = f.rowwise().mean();
	for(unsigned int i = 0; i < f.cols(); ++i) f.col(i) -= mean;
	Eigen::JacobiSVD<Eigen::Matrix3Xd > jSVD;
	jSVD.compute(f, Eigen::ComputeFullU);
	
	//project points to basis
	f = jSVD.matrixU().transpose()*f;
	f.row(2) = Eigen::VectorXd::Constant(f.cols(), 0.0);
	
	//partition
	TPPLPoly p;
	p.Init(f.cols());
	for(unsigned int i = 0; i < f.cols(); ++i){
		p[i].x= f(0,i);
		p[i].y= f(1,i);
	}
	bool inv=(p.GetOrientation() == TPPL_CW);
	if(inv){  p.Invert(); }
	
	TPPLPartition tp;
	std::list<TPPLPoly> p_;
	int ret=-1;
	switch ( method )
	{
	case 0:
		ret=tp.Triangulate_EC(&p,&p_);
		break;
	case 1:
		ret=tp.Triangulate_OPT(&p,&p_);
		break;
	case 2:
		ret=tp.Triangulate_MONO(&p,&p_);
		break;
	case 3:
		ret=tp.ConvexPartition_HM(&p,&p_);
		break;
	case 4:
		ret=tp.ConvexPartition_OPT(&p,&p_);
		break;
	default:
		std::cerr << __FUNCTION__ << ": Invalid method nr. " << method << std::endl;
		break;
	}
	//dolphinOut() << __FUNCTION__ << ": Number of new faces: " << p_.size() << std::endl;
	
	//collect vertex handles of new faces
	newF.clear();
	for( std::list<TPPLPoly>::const_iterator it=p_.begin();it!=p_.end();++it){
		TPPLPoly pi=*it;
		//pi.SetOrientation(p.GetOrientation());
		if(inv){  pi.Invert(); }
		long n=pi.GetNumPoints();
		newF.push_back(std::vector<int>());
		for( unsigned int j=0;j<n;++j){
			TPPLPoint pt = pi.GetPoint(j);
			
			//Since TPPL does not store indices, we need to find the row index of pt in f:
			int ch=-1;
			for(unsigned int c=0; c<f.cols();++c ){
				double d= std::fabs(pt.x - f(0,c)) + std::fabs(pt.y - f(1,c));
				if( 1e-10 > d ){
					ch=c;
					newF.back().push_back(c);
					break;
				}
			}
			assert(ch>=0);
		}
	}
	
	return ret;
}
}
