#include "MGCLStdAfx.h"
#include "mg/Tolerance.h"
#include "mg/Transf.h"
#include "topo/Edge.h"
#include "topo/Loop.h"
#include "topo/Face.h"
#include "Tl2/TL2parameter.h"
#include "Tl2/TL2Fans.h"
#include "Tl2/TL2LPline.h"
#include "Tl2/TL2Face.h"

#if defined(_DEBUG)
#define new DEBUG_NEW
#undef THIS_FILE
static char THIS_FILE[] = __FILE__;
#endif


/****************************************************************/
/*   Copyright (c) 2017 by System fugen G.K.                */
/*                       All rights reserved.                   */
/****************************************************************/

///Check the angle of all the vertices, the concavity or sharp angle.
///Function's return value is -1 if concavity found.
///                            1 if sharp angle found.
///                            0 if both concavity ro sharp angle not found.
int checkAngle(
	int nline,///<number of lines of LPline.
	const mgTL2LPline* LPline,///<The target lines in LPline[.].
	int& vertexID///<if function's return value is 1 or -1, the vertex id will be returned.
){
	const MGVector& N=MGDefault::z_unit_vector();
	double angle[4];//angle around the normal at the 4 veritces will be stored.
	const mgTL2parameter& para=LPline[0].TL2param();
	double uerr=para.get_UError(), verr=para.get_VError();
	double paraerr=para.get_UVError();
	int nlinem1=nline-1;
	int idmin, idmax;
	double angMin=mgDBLPAI, angMax=-mgDBLPAI;
	for(int i=0; i<nline; i++){
		const mgTL2LPline& currentL=LPline[i];
		//std::cout<<currentL<<std::endl;
		int preID=(i+nlinem1)%nline;
		const mgTL2LPline& preL=LPline[preID];
		MGPosition uv=currentL.uv(0);//uv of start point.
		MGVector Tpre=preL.eval(1.,1);
		MGVector Tnow=currentL.eval(0.,1);
		if(Tpre.len()<=paraerr || Tnow.len()<=paraerr)
			continue;
		double ang=Tpre.angle2pai(Tnow,N);
		if(ang>=mgPAI)
			ang-=mgDBLPAI;
		angle[i]=ang;
		if(ang<angMin){
			angMin=ang;
			idmin=i;
		}
		if(ang>angMax){
			angMax=ang;
			idmax=i;
		}
	}
	if(CONCAVEANGLE2<=angMin && angMin<=LOOSE_ZERO_ANGLE){
		vertexID=idmin;
		return -1;
	}else if(angMax>SHARPANGLE){
		vertexID=idmax;
		return 1;
	}
	return 0;
}

//Interpolate pline[idEdge] and pline[(idEdge+2)%4] for the tessellation.
//Let id0=idEdge, id1=(idEdge+1)%4, id2=(idEdge+2)%4, id3=(idEdge+3)%4;
//Then,the direction of the output is the same as pline[id2],
//that is, always is opposite to pline[id0].
//Interpolation is done from vetex id1V of pline[id1]
//to the vetex id3V of pline[id3].
std::unique_ptr<mgTL2Polyline> interpolateWithVertexNum(
	int idEdge,	//id of edge that interpolate.
	const mgTL2LPline pline[4],
	int id1V,	//id of pline[id1]'s vertex.
	int id3V	//id of pline[id3]'s vertex.
){
	const mgTL2LPline& lp0=pline[idEdge];
	const mgTL2LPline& lp1=pline[(idEdge+1)%4];
	const mgTL2LPline& lp2=pline[(idEdge+2)%4];
	const mgTL2LPline& lp3=pline[(idEdge+3)%4];

	std::unique_ptr<mgTL2Polyline> uvpolyline
		=std::unique_ptr<mgTL2Polyline>(new mgTL2Polyline(lp0.TL2param()));

	double w2=double(id1V)/double(lp1.number_of_points()-1);
	double w1=1.-w2;
	int n0=lp0.number_of_points(), n2=lp2.number_of_points();
	//int nbd=int(double(n0)*w1+double(n2)*w2+.5);
	int nbd=int(double(n0)*w1+double(n2)*w2+.05);
		//Interpolated mgTL2Polyline's vertex number.
	assert(nbd>=2);
	MGBPointSeq& uvbp=uvpolyline->line_bcoef();
	uvbp.resize(nbd,2);
	MGPosition uvS=lp1.uv(id1V),uvE=lp3.uv(id3V);
	uvbp.store_at(0,uvS);

	MGKnotVector& uvKnotv=uvpolyline->knot_vector();
	uvKnotv.size_change(2,nbd);
	uvKnotv(0)=uvKnotv(1)=0.;

	int nbdm1=nbd-1;
	double tnbdm1=double(nbdm1);
	if(nbd>2){
		mgTL2LPline lp0R(lp0);
		lp0R.reverse();
		MGPosition uv00=lp0R.uv(0),uv01=lp0R.uv(n0-1);
		MGTransf tr1(uv00,uv01,uvS,uvE);
		MGPosition uv20=lp2.uv(0),uv21=lp2.uv(n2-1);
		MGTransf tr2(uv20,uv21,uvS,uvE);

		for(int j=1; j<nbdm1; j++){
			double tj=double(j);
			double tjOtnbdm1=tj/tnbdm1;
			MGPosition uv1=lp0R.eval(tjOtnbdm1);
			uv1*=tr1;
			MGPosition uv2=lp2.eval(tjOtnbdm1);
				uv2*=tr2;
			MGPosition uv=uv1*w1+uv2*w2;
			uvbp.store_at(j,uv);
			uvKnotv(j+1)=tj;//Set the uniform knot configuration of uvpolyline.
		}
	}
	uvbp.store_at(nbdm1,uvE);
	uvKnotv(nbd+1)=uvKnotv(nbd)=tnbdm1;

	short idBP[3];
	if(lp1.get_id_from_VertexID(id1V,idBP))
		uvpolyline->set_startID(idBP);
	if(lp3.get_id_from_VertexID(id3V,idBP))
		uvpolyline->set_endID(idBP);
	return uvpolyline;
}

//When pline[i] is null and pline[j] for j>i is not null, replace [i] with [j].
void compressNull(mgTL2LPline pline[4]){
	for(int i=0; i<4; i++){
		if(pline[i].is_null()){
			int i2=i;
			for(int j=i+1; j<4; j++){
				if(pline[j].is_null())
					continue;
				else{
					pline[i2++]=pline[j];
					pline[j].setNull();
				}
			}
		}
	}
}

///Let id0=idEdge, id2=(idEdge+2)%4, then
///Subdivide rectangle pline[id0] at the vertex id0V, and
///subdivide rectangle pline[id2] at the vertex id2V.
///Then generate the two subdivided rectangles pline1 and pline2.
///Subdividing pline is input as midpline whose direction is from id0 to id2.
///pline[idEdge+1] or pline[idEdge+3] can be null mgTL2LPline for 3 or 2 edges pline.
void subdivideWithMidline(
	const std::unique_ptr<mgTL2Polyline>& midpline, //mgTL2Polyline that subdivide pline.
	int idEdge,	//id of edge that indicates which edge be subdivided.
	const mgTL2LPline pline[4],
	int id0V,
	int id2V,
	mgTL2LPline pline1[4],	//1st subdivided rectangle.
	mgTL2LPline pline2[4]	//2nd subdivided rectangle.
){
	int id0=idEdge;
	int id1=(id0+1)%4;
	int id2=(id1+1)%4;
	int id3=(id2+1)%4;
	//build Edge0.
	pline[idEdge].subdivide(id0V,pline1[id0],pline2[id0]);
	//std::cout<<pline1[id0]<<pline2[id0]<<std::endl;

	//build Edge1.
	pline1[id1]=mgTL2LPline(midpline.get());
	pline2[id1]=pline[(idEdge+1)%4];//Can be null

	//build Edge2.
	const mgTL2LPline& plineOpo=pline[(idEdge+2)%4];
	if(!plineOpo.is_null()){
		int nOpo=plineOpo.number_of_points();
		if(nOpo>2)
			plineOpo.subdivide(id2V,pline2[id2],pline1[id2]);
		else
			pline2[id2]=plineOpo;
		//std::cout<<pline1[id2]<<pline2[id2]<<std::endl;
	}

	//build Edge3.
	pline1[id3]=pline[(idEdge+3)%4];//can be null
	pline2[id3]=mgTL2LPline(midpline.get()); pline2[id3].reverse();
	compressNull(pline1);
	compressNull(pline2);
	//std::cout<<pline1[0]<<pline1[1]<<pline1[2]<<pline1[3]<<std::endl;
	//std::cout<<pline2[0]<<pline2[1]<<pline2[2]<<pline2[3]<<std::endl;
}

///Let id1=(idEdge+1)%4, id3=(idEdge+3)%4, then
///Subdivide rectangle pline[id1] at the vertex id1V, and
///subdivide rectangle pline[id3] at the vertex id3V.
///And interpolate from the vertex id1V to the vertex id3V.
///Then generate the two subdivided rectangles.
///Subdivided pline will be output into Pline1 and 2.
void subdivide(
	int idEdge,	//id of edge that indicates which edge be subdivided.
	const mgTL2LPline pline[4],
	int id1V,
	int id3V,
	mgTL2LPline pline1[4],	//1st subdivided rectangle.
	mgTL2LPline pline2[4],	//2nd subdivided rectangle.
	std::unique_ptr<mgTL2Polyline>& midpline //mgTL2Polyline that subdivide pline.
){
	midpline=interpolateWithVertexNum(idEdge,pline,id1V,id3V);
	subdivideWithMidline(midpline,(idEdge+1)%4,pline,id1V,id3V,pline1,pline2);
}

///Subdivide rectangle pline[.] at the middle of u(when idEdge%2 =1) or v(when idEdge%2 =0).
///Subdivided pline will be output into Pline1 and 2.
void subdivideMiddle(
	int idEdge,	//id of edge that indicates which edge be subdivided.
	const mgTL2LPline pline[4],
	mgTL2LPline pline1[4],	//1st subdivided rectangle.
	mgTL2LPline pline2[4],	//2nd subdivided rectangle.
	std::unique_ptr<mgTL2Polyline>& midpline //mgTL2Polyline that subdivide pline.
){
	const mgTL2LPline& lp1=pline[(idEdge+1)%4];
	const mgTL2LPline& lp3=pline[(idEdge+3)%4];

	int n1=lp1.number_of_points(), n3=lp3.number_of_points();
	int id1V=n1/2;
	int id3V=n3-1;id3V-=n3/2;
		//This is because lp1 has the opposite direction of lp3's.

	subdivide(idEdge,pline,id1V,id3V,pline1,pline2,midpline);
}

#define ALLOWANCE 1.4
///Tessellate a rectangle that has 2 non-continuous edges of only 2 vertices.
///Tessellated triangles will be appended.
///pl.number_of_points()>=pl_opposite.number_of_points() is assumed.
void mgTL2Triangles::tessellate2n2n(
	const mgTL2LPline& pl,//an edge whose number of vetices is more than 2.
	const mgTL2LPline& pl_opposite//Opposite edge of pl.
){
	int np1=pl.number_of_points();assert(np1>=3);
	int np2=pl_opposite.number_of_points();assert(np2>=3);
	int np1half=(np1+1)/2;
	int np1r=np1-np1half+1;
	mgTL2LPline pl10(pl,0,np1half);///HALF

	int np2half=(np2+1)/2;
	if(np2half+1<np2){
		MGVector P1=pl.xyz(np1half-1);
		MGVector V21=pl_opposite.xyz(np2half)-P1;//MGPosition uv1=pl_opposite.uv(np2half);
		MGVector V22=pl_opposite.xyz(np2half-1)-P1;//MGPosition uv2=pl_opposite.uv(np2half-1);	
		if(V21%V21<V22%V22)
			np2half+=1;
	}
	int np2r=np2-np2half+1;
	mgTL2LPline pl21(pl_opposite,np2half-1,np2r);///Other half

	if(np1half==2){
		if(np2r==2)
			tessellate2by2(pl10,pl21);
		else
			tessellate222n(pl21,pl10);
	}else{
		if(np2r==2)
			tessellate222n(pl10,pl21);
		else{
			if(np1half<np2r)
				tessellate2n2n(pl21,pl10);
			else
				tessellate2n2n(pl10,pl21);
		}
	}

	mgTL2LPline pl11(pl,np1half-1,np1r);
	mgTL2LPline pl20(pl_opposite,0,np2half);
	if(np1r==2){
		if(np2half==2)
			tessellate2by2(pl11,pl20);
		else
			tessellate222n(pl20,pl11);
	}else{
		if(np2half==2)
			tessellate222n(pl11,pl20);
		else
			tessellate2n2n(pl11,pl20);
	}
}

///Tessellate a rectangle in case that one of the edges has
///2 or 3 vertices,and the other 3 edges have 3, or more than 3 vertices.
///Tessellated triangles will be appended.
void mgTL2Triangles::tessellate3nnn(
	int id3,	//id of edge that has 2 or 3 vertices.
	const mgTL2LPline pline[4]
){
	mgTL2LPline pline1[4];
	mgTL2LPline pline2[4];
	std::unique_ptr<mgTL2Polyline> midpline;
	subdivideMiddle(id3,pline,pline1,pline2,midpline);
	tessellate4(pline1);
	tessellate4(pline2);
}

//Tessellate a rectangle of a General case.
//All of the four edges have 4 or more than 4 vertices.
//Tessellated triangles will be appended.
void mgTL2Triangles::tessellate_nnnn(
	int minID,
	const mgTL2LPline pline[4]
){
	mgTL2LPline pline1[4];
	mgTL2LPline pline2[4];
	std::unique_ptr<mgTL2Polyline> midpline;
	subdivideMiddle(minID,pline,pline1,pline2,midpline);
	tessellate4(pline1);
	tessellate4(pline2);
}

///Tessellate a polygon pline whose two opposing edges' number of vertices
///are the same and the other edges' numbers of vertices differ at most 1.
///Let id0=maxID, id1=(id0+1)%4, id2=(id0+2)%4, id3=(id0+3)%4.
///Then the number of vertices of id1 and id3 are the same.
///The number of id0 is greater or equal to id2, and the difference is
///at most 1.
///Tessellated triangles will be appended.
void mgTL2Triangles::tessellate_nm(
	int maxID,
	const mgTL2LPline pline[4]
){
	int n1m1=pline[(maxID+1)%4].number_of_points()-1;

	mgTL2LPline pline0[4];
	mgTL2LPline pline1[4];	//1st subdivided rectangle.
	mgTL2LPline pline2[4];	//2nd subdivided rectangle.
	for(int i=0; i<4; i++)
		pline2[i]=pline[(maxID+i)%4];
	std::unique_ptr<mgTL2Polyline> midpline1, midpline2; //mgTL2Polyline that subdivide pline.
	for(int i=0; i<n1m1; i++){
		for(int j=0; j<4; j++)
			pline0[j]=pline2[j];
		int npl01=pline0[1].number_of_points();
		if(npl01>2){
			subdivide(0,pline0,1,npl01-2,pline1,pline2,midpline1);
		}else{
			pline1[0]=pline0[0];
			pline1[2]=pline0[2];
		}
		makeStrip(pline1[0], pline1[2]);
		midpline2=std::move(midpline1);
	}
}

///Do perform the tessellation for a concave rectangle.
///The result will be appended onto triangles.
///Tessellated triangles will be appended.
void mgTL2Triangles::tessellateConcave4(
	int concaveID,	//vertiex id at which vertex concavity is detected.
	const mgTL2LPline LPline[4]
){
	const mgTL2LPline& lp0=LPline[concaveID];
	const mgTL2LPline& lp1=LPline[(concaveID+1)%4];
	const mgTL2LPline& lp2=LPline[(concaveID+2)%4];
	const mgTL2LPline& lp3=LPline[(concaveID+3)%4];
	MGVector v1=lp0.eval(0.,1);
	MGVector v2=lp3.eval(1.,1);
	MGStraight slmid(v2,v1,lp0.uv(0));
	int n2=lp2.number_of_points();
	const mgTL2LPline* lp=&lp1;
	int idisec1=lp->isectSlTl(slmid);
	if(idisec1==0)
		idisec1=1; 
	if(idisec1<0){
		lp=&lp2;
		idisec1=lp->isectSlTl(slmid);
		if(idisec1==n2-1)
			idisec1=n2-2;
		else if(idisec1<0)
			idisec1=0;
	}

	std::unique_ptr<mgTL2Polyline> midLine=lp0.polygonizeSL(*lp,0,idisec1);
	mgTL2Polyline* midLineP=midLine.get();
	mgTL2LPline LPline1[4], LPline2[4];
	int nV=midLine->number_of_points();
	if(lp==&lp1){
		LPline1[0]=lp0;
		LPline1[1]=mgTL2LPline(lp1,0,idisec1+1);
		LPline1[2]=mgTL2LPline(midLineP,nV-1,-nV);
		tessellate3(LPline1);

		LPline2[0]=mgTL2LPline(midLineP,0,nV);
		int id=1, n1=lp1.number_of_points();
		if(idisec1<n1-1)
			LPline2[id++]=mgTL2LPline(lp1,idisec1,n1-idisec1);
		LPline2[id++]=lp2;
		LPline2[id++]=lp3;
		if(id==3)
			tessellate3(LPline2);
		else
			tessellate4(LPline2);
	}else{
		LPline1[0]=lp0;
		LPline1[1]=lp1;
		int id=2;
		if(idisec1>0)
			LPline1[id++]=mgTL2LPline(lp2,0,idisec1+1);
		LPline1[id++]=mgTL2LPline(midLineP,nV-1,-nV);
		if(id==3)
			tessellate3(LPline1);
		else
			tessellate4(LPline1);

		LPline2[0]=mgTL2LPline(midLineP,0,nV);
		LPline2[1]=mgTL2LPline(lp2,idisec1,n2-idisec1);
		LPline2[2]=lp3;
		tessellate3(LPline2);
	}
}

///Do perform the tessellation for a sharp vertex rectangle.
///The result will be appended onto triangles.
///Tessellated triangles will be appended.
void mgTL2Triangles::tessellateSharp4(
	int sharpID,	//vertiex id at which sharpness is detected.
	const mgTL2LPline LPline[4]
){
	const mgTL2LPline& lp0=LPline[sharpID];
	int n0=lp0.number_of_points();
	const mgTL2LPline& lp1=LPline[(sharpID+1)%4];
	int n1=lp1.number_of_points();
	const mgTL2LPline& lp2=LPline[(sharpID+2)%4];
	int n2=lp2.number_of_points();
	const mgTL2LPline& lp3=LPline[(sharpID+3)%4];
	int n3=lp3.number_of_points();

	int idlp0=(n0-1)/2, idlp3=n3/2;
	if(n0<=2){
		idlp0=1;
		idlp3=0;
	}
	if(n3<=2 || (n1==2 && n2==2))
		idlp3=0;
	std::unique_ptr<mgTL2Polyline> midLine=lp0.polygonizeSL(lp3,idlp0,idlp3);
	int nV=midLine->number_of_points();

	mgTL2Polyline* midLineP=midLine.get();
	mgTL2LPline LPline1[4];
	LPline1[0]=mgTL2LPline(lp0,0,idlp0+1);	
	LPline1[1]=mgTL2LPline(midLineP,0,nV);
	LPline1[2]=mgTL2LPline(lp3,idlp3,n3-idlp3);
	tessellate3(LPline1);
	if(idlp3==0){
		if(n0<=2){
			LPline1[0]=mgTL2LPline(lp1);	
			LPline1[1]=mgTL2LPline(lp2);
			LPline1[2]=mgTL2LPline(midLineP,nV-1,-nV);
			tessellate3(LPline1);
		}else{
			LPline1[0]=mgTL2LPline(lp0,idlp0,n0-idlp0);	
			LPline1[1]=mgTL2LPline(lp1);
			LPline1[2]=mgTL2LPline(lp2);
			LPline1[3]=mgTL2LPline(midLineP,nV-1,-nV);
			tessellate4(LPline1);
		}
		return;
	}

	mgTL2LPline LPline2[4];
	if(n2>=n1){
		int idlp2=(n2-1)/2;
		std::unique_ptr<mgTL2Polyline> midLine2=lp0.polygonizeSL(lp2,idlp0,idlp2);
		int nV2=midLine2->number_of_points();
		mgTL2Polyline* midLine2P=midLine2.get();
		LPline2[0]=mgTL2LPline(lp3,0,idlp3+1);	
		LPline2[1]=mgTL2LPline(midLineP,nV-1,-nV);
		LPline2[2]=mgTL2LPline(midLine2P,0,nV2);
		LPline2[3]=mgTL2LPline(lp2,idlp2,n2-idlp2);
		tessellate4(LPline2);

		LPline2[0]=mgTL2LPline(lp0,idlp0,n0-idlp0);	
		LPline2[1]=mgTL2LPline(lp1);
		LPline2[2]=mgTL2LPline(lp2,0,idlp2+1);
		LPline2[3]=mgTL2LPline(midLine2P,nV2-1,-nV2);
		tessellate4(LPline2);
	}else{
		int idlp1=(n1-1)/2;
		std::unique_ptr<mgTL2Polyline> midLine2=lp3.polygonizeSL(lp1,idlp3,idlp1);
		int nV2=midLine2->number_of_points();
		mgTL2Polyline* midLine2P=midLine2.get();
		LPline2[0]=mgTL2LPline(lp1,0,idlp1+1);	
		LPline2[1]=mgTL2LPline(midLine2P,nV2-1,-nV2);
		LPline2[2]=mgTL2LPline(midLineP,nV-1,-nV);
		LPline2[3]=mgTL2LPline(lp0,idlp0,n0-idlp0);
		tessellate4(LPline2);

		LPline2[0]=mgTL2LPline(lp1,idlp1,n1-idlp1);	
		LPline2[1]=mgTL2LPline(lp2);
		LPline2[2]=mgTL2LPline(lp3,0,idlp3+1);
		LPline2[3]=mgTL2LPline(midLine2P,0,nV2);
		tessellate4(LPline2);
	}
}

///Do perform the tessellation for the mgTL2Plygon that has 4 edges as the boundary.
///The result will be appended onto triangles.
///When triangles.is_uv()=false, all of the element of the triangle position data
///has normal data as (x,y,z,xn,yn,zn). Here (x,y,z) is the position data and
///(xn,yn,zn) is the normal vector at the position (x,y,z).
///When triangles.is_uv()=true, all of the element of the triange position data are (u,v).
///Tessellated triangles will be appended.
void mgTL2Triangles::tessellate4(
	const mgTL2LPline LPline[4]
){
	//std::cout<<LPline[0]<<LPline[1]<<LPline[2]<<LPline[3]<<std::endl;//////****
	int nVert[4];//Vertex number of pline[i] will be stored.
	int num2=0;
	int nVMin=-1,minID=0;//minimum number of vertices and the edge id will be stored.
	int nvMax=0, eidMaxN=0;//Maximum vertex number and the id of pline[i] will be stored.
	for(int i=0; i<4; i++){
		const mgTL2LPline& currentL=LPline[i];
		int nVi=nVert[i]=currentL.number_of_points();
		if(nVi==2)
			num2++;
		if(nVMin<0 || nVi<nVMin){
			nVMin=nVi;
			minID=i;
		}
		if(nVi>nvMax){
			nvMax=nVi;
			eidMaxN=i;
		}
	}

	int eidNext=(eidMaxN+1)%4;
	int eidOpo=(eidMaxN+2)%4;
	int eidNextOpo=(eidMaxN+3)%4;
	int otherNum=LPline[eidNext].number_of_points();
	int otherNum2=LPline[eidNextOpo].number_of_points();
	if(otherNum<otherNum2)
		otherNum=otherNum2;
	int nvOpo=LPline[eidOpo].number_of_points();

	if(nvOpo>=otherNum*DEVIDE_RATIO){
		const mgTL2LPline& lp0=LPline[eidMaxN];
		const mgTL2LPline& lp2=LPline[eidOpo];
		//std::cout<<lp0<<lp2<<std::endl;
		int id0V, id2V;
		mgTL2LPline LPline1[4], LPline2[4];
		std::unique_ptr<mgTL2Polyline> midLine=lp0.getPolygonizedMidLine(lp2,id0V,id2V);
		subdivideWithMidline(midLine,eidMaxN,LPline,id0V,id2V,LPline1,LPline2);
		tessellate4(LPline1);
		tessellate4(LPline2);
		return;
	}

	if(num2==4){//Case that all of the edges have only 2 vertices.
		tessellate2by2(LPline[0],LPline[2]);
		return;
	}

	if(num2==3){//Case that 3 of the edges have only 2 vertices.
		int i=0;
		for(; i<4; i++){
			if(nVert[i]>=3)
				break;
		}
		tessellate222n(LPline[i],LPline[(i+2)%4]);
		return;
	}

	if(nVert[0]==nVert[2]){
		int ndif=nVert[1]-nVert[3];
		if(-1<=ndif && ndif<=1){
			if(ndif>=0)
				tessellate_nm(1,LPline);
			else
				tessellate_nm(3,LPline);
			return;
		}
	}
	if(nVert[1]==nVert[3]){
		int ndif=nVert[0]-nVert[2];
		if(-1<=ndif && ndif<=1){
			if(ndif>=0)
				tessellate_nm(0,LPline);
			else
				tessellate_nm(2,LPline);
			return;
		}
	}

	int minIDp3=(minID+3)%4;
	if(num2==2){
		if(nVert[(minID+1)%4]==2){//Case that 2 continuous edges have only 2 vertices.
			tessellate22nn(minID,LPline);
		}else if(nVert[minIDp3]==2){//Case that 2 continuous edges have only 2 vertices.
			tessellate22nn(minIDp3,LPline);
		}else{
			const mgTL2LPline& lp1=LPline[minIDp3];
			const mgTL2LPline& lp2=LPline[(minIDp3+2)%4];
			if(lp1.number_of_points()>=lp2.number_of_points())
				tessellate2n2n(lp1,lp2);
			else
				tessellate2n2n(lp2,lp1);
		}
		return;
	}

	int vertexID;
	int sharpOrConcave=checkAngle(4,LPline, vertexID);
	if(sharpOrConcave==-1){
		tessellateConcave4(vertexID,LPline);
		return;
	}else if(sharpOrConcave==1){
		tessellateSharp4(vertexID,LPline);
		return;
	}

	if(nVMin<=3){
		//Case that one of the edges has at most 3(that is 2 or 3) vertices,
		//and the other 3 edges have 2, 3, or more than 3 vertices.
		tessellate3nnn(minID,LPline);
		return;
	}

	//General case. All of the four edges have 4 or more than 4 vertices.
	tessellate_nnnn(minID,LPline);
}

///Do perform the tessellation for the mgTL2Plygon that has 2 edges as the boundary.
///The result will be appended onto triangles.
///Tessellated triangles will be appended.
void mgTL2Triangles::tessellate2(
	const mgTL2LPline LPline[2]
){
	//std::cout<<LPline[0]<<LPline[1]<<LPline[2]<<std::endl;//////****
	const mgTL2LPline& lp0=LPline[0];
	const mgTL2LPline& lp1=LPline[1];
	int n0=lp0.number_of_points();
	int n1=lp1.number_of_points();
	if(n0<=2 && n1<=2){
		makeFan(lp0,lp1);
		return;
	}
	int id0V=n0/2;
	int id1V=(n1-1)/2;
	std::unique_ptr<mgTL2Polyline> midLine=lp0.polygonizeSL(lp1,id0V,id1V);
	mgTL2LPline LPlineTemp[4];
	LPlineTemp[0]=lp0; LPlineTemp[2]=lp1;

	mgTL2LPline LPline1[4], LPline2[4];
	subdivideWithMidline(midLine,0,LPlineTemp,id0V,id1V,LPline1,LPline2);
	if(LPline1[2].is_null())
		tessellate2(LPline1);
	else
		tessellate3(LPline1);
	tessellate3(LPline2);
}

///Do perform the tessellation for the mgTL2Plygon that has 3 edges as the boundary.
///The result will be appended onto triangles. Tessellated triangles will be appended.
void mgTL2Triangles::tessellate3(
	const mgTL2LPline LPline[3]
){
	//std::cout<<LPline[0]<<LPline[1]<<LPline[2]<<std::endl;//////****
	int num2=0;
	int nVMin=LPline[0].number_of_points();//minimum number of vertices will be stored.
	if(nVMin==2)
		num2++;
	int nVMax=nVMin;
	int minID=0, maxID=0;
	for(int i=1; i<3; i++){
		int nVi=LPline[i].number_of_points();
		if(nVi==2)
			num2++;
		if(nVi<nVMin){
			nVMin=nVi;
			minID=i;
		}
		if(nVi>nVMax){
			nVMax=nVi;
			maxID=i;
		}
	}

	if(num2>=2){
		int id=(minID+1)%3;
		if(LPline[id].number_of_points()==2)
			id=(minID+2)%3;

		const mgTL2LPline& lpn=LPline[id];
		const mgTL2LPline& lpnPre=LPline[(id+2)%3];
		makeFan(lpn,lpnPre);
	}else{
		mgTL2LPline pline1[4], pline2[4];
		int vertexID;
		int Concave=checkAngle(3,LPline, vertexID);
		if(Concave==-1){
			const mgTL2LPline& lp0=LPline[vertexID];
			const mgTL2LPline& lp1=LPline[(vertexID+1)%3];
			const mgTL2LPline& lp2=LPline[(vertexID+2)%3];
			int n1=lp1.number_of_points();
			if(n1>2){
				int idlp1=n1/2;
				pline1[0]=lp0;
				std::unique_ptr<mgTL2Polyline> midpline=lp0.polygonizeSL(lp1,0,idlp1);
				pline1[2]=mgTL2LPline(midpline.get()); pline1[2].reverse();
				lp1.subdivide(idlp1,pline1[1],pline2[0]);
				pline2[1]=lp2;
				pline2[2]=mgTL2LPline(midpline.get());
				tessellate3(pline1);
				tessellate3(pline2);
				return;
			}
		}

		const mgTL2LPline& lp1=LPline[(minID+1)%3];
		const mgTL2LPline& lp2=LPline[(minID+2)%3];
		int idlp1=lp1.number_of_points()/2;
		int idlp2=(lp2.number_of_points()-1)/2;
		pline1[0]=LPline[minID];
		std::unique_ptr<mgTL2Polyline> midpline=lp1.polygonizeSL(lp2,idlp1,idlp2);
		lp1.subdivide(idlp1,pline1[1],pline2[0]);
		pline1[2]=mgTL2LPline(midpline.get());
		lp2.subdivide(idlp2,pline2[1],pline1[3]);
		pline2[2]=mgTL2LPline(midpline.get()); pline2[2].reverse();
		tessellate3(pline2);
		tessellate4(pline1);
	}
}

///Do perform the tessellation for the mgTL2Plygon that has 3 or 4 edges outer loop, and
///that has no inner loops.
///The result will be appended onto triangles.
///When triangles.is_uv()=false, all of the element of the triangle position data has normal data as
///(x,y,z,xn,yn,zn). Here (x,y,z) is the position data and (xn,yn,zn) is the normal vector
///at the position (x,y,z).
///When triangles.is_uv()=true, all of the element of the triange position data are (u,v).
///Tessellated triangles will be appended.
void mgTL2Triangles::tessellate1234edgesLoop(
	const MGLoop& oloop//Outer boundary loop of the face of mgTL2Plygon
					//to tessellate whose number of edges is 4.
){
	mgTL2LPline LPline[4];
	MGComplex::const_pcellItr ei=oloop.pcell_begin();
	int nedges=oloop.number_of_edges();
	int nvMax=0, eidMaxN=0;//Maximum vertex number and the id of pline[i] will be stored.
	for(int i=0; i<nedges; ei++, i++){
		const mgTL2Polyline* pli=TL2Polyline(ei);
		int nVi=pli->number_of_points();
		LPline[i]=mgTL2LPline(pli,0,nVi);
		if(nVi>nvMax){
			nvMax=nVi;
			eidMaxN=i;
		}
	}

	if(nedges==1){
		const mgTL2LPline& lp0=LPline[0];//std::cout<<lp0<<std::endl;
		int id0V=lp0.number_of_points()/2;
		std::unique_ptr<mgTL2Polyline> midLine=lp0.polygonizeSL(lp0,id0V,0);//std::cout<<*midLine<<std::endl;
		mgTL2LPline LPline1[4], LPline2[4];
		subdivideWithMidline(midLine,0,LPline,id0V,0,LPline1,LPline2);
		tessellate2(LPline1);
		tessellate2(LPline2);
	}else if(nedges==2){
		tessellate2(LPline);
	}else if(nedges==3){
		tessellate3(LPline);
	}else{
		int eidNext=(eidMaxN+1)%4;
		int eidOpo=(eidMaxN+2)%4;
		int eidNextOpo=(eidMaxN+3)%4;
		int otherNum=LPline[eidNext].number_of_points();
		int otherNum2=LPline[eidNextOpo].number_of_points();
		if(otherNum<otherNum2)
			otherNum=otherNum2;
		int nvOpo=LPline[eidOpo].number_of_points();

		if(nvOpo>=5){//This means nvMax>=nvOpo>=5.
			const mgTL2LPline& lp0=LPline[eidMaxN];
			const mgTL2LPline& lp2=LPline[eidOpo];
			int id0V, id2V;
			std::unique_ptr<mgTL2Polyline> midLine=lp0.getPolygonizedMidLine(lp2,id0V,id2V);
			if(midLine.get()){
				if(midLine->number_of_points()>otherNum*3){
					mgTL2LPline LPline1[4], LPline2[4];
					subdivideWithMidline(midLine,eidMaxN,LPline,id0V,id2V,LPline1,LPline2);
					tessellate4(LPline1);
					tessellate4(LPline2);
					return;
				}
			}
		}
		tessellate4(LPline);
	}
}
