#include "Stdafx.h"
#include "GridFaces.h"
#include "SplitAllInsectionCurves.h"
//#include "atlstr.h"
#include "TColStd_Array1OfInteger.hxx"
#include "Geom_BSplineSurface.hxx"
#include "TColgp_Array1OfPnt.hxx"
#include "TColStd_Array1OfReal.hxx"



GridFaces::GridFaces( void )
{
	isDone = false;
}

GridFaces::GridFaces( map<int, vector<TopoDS_Edge>> wires )
{
	this->wireSet = wires;
	isDone = false;
}

GridFaces::GridFaces( vector<TopoDS_Edge> wires) //20180101 modified by czb~{#,~} ~{Tv<SPr:EWwN*2NJ}~}
{
	CSplitAllInsectionCurves m_tool;
	vector<std::string> strNames;
	for(int i = 0; i < wires.size(); i++)
	{
		std::string linename;
		linename = std::to_string(i+1);
		strNames.push_back(linename);
	}
	m_tool.Init(wires, strNames);
	m_tool.Build(false);
	this->wireSet = m_tool.GetResult();
	wireOldSet = m_tool.GetOldResult();
}

GridFaces::GridFaces( vector<TopoDS_Edge> wires,std::vector<std::string> seqNo) //20180101 modified by czb~{#,~} ~{Tv<SPr:EWwN*2NJ}~}
{
	CSplitAllInsectionCurves m_tool;
	//vector<std::string> strNames;
	//for(int i = 0; i < wires.size(); i++)
	//{
	//	std::string linename;
	//	linename.Format(_T("%d"),i);
	//	strNames.push_back(linename);
	//}
	//m_tool.Init(wires, strNames);
	m_tool.Init(wires, seqNo);
	m_tool.Build(false);
	this->wireSet = m_tool.GetResult();
	wireOldSet = m_tool.GetOldResult();
}

// ADDED BY XUEFENG 20180731
// FOR CODES COMBINE
//ADDED BY XUEFENG 201805 
// ADDED exp wires NO
GridFaces::GridFaces( vector<TopoDS_Edge> wires,std::vector<std::string> seqNo, vector<int> expSeqNo) 
{
	CSplitAllInsectionCurves m_tool;
	//vector<std::string> strNames;
	//for(int i = 0; i < wires.size(); i++)
	//{
	//	std::string linename;
	//	linename.Format(_T("%d"),i);
	//	strNames.push_back(linename);
	//}
	//m_tool.Init(wires, strNames);
	
	//m_tool.Init(wires, seqNo);//DELETED BY XUEFENG 201805 
	m_tool.ExpInit(wires, seqNo, expSeqNo); //ADDED BY XUEFENG 201805 
	m_tool.Build(false);
	this->wireSet = m_tool.GetResult();
	wireOldSet = m_tool.GetOldResult();
}
// END ADDED
//END ADDED

GridFaces::~GridFaces(void)
{
}


TopoDS_Shape GridFaces::GetResult()
{
	TopoDS_Shape shape;
	if(isDone)
	{
		shape = this->result;
	}
	return shape;
}

void GridFaces::Build()
{
	TopoDS_Shell shell;
	BRep_Builder b;
	b.MakeShell(shell);
	static int ii = 0;
	ii++;
	int k = 0;
	vector<TopoDS_Shape>  faceSet;
	for(map<int, vector<TopoDS_Edge>>::iterator iter = this->wireSet.begin(); iter != this->wireSet.end(); iter++)
	{
		k++;
		vector<TopoDS_Edge> edgeSet = iter->second;
		if(k==27)
		{
			k=k;
			////continue;
			//TopoDS_Compound c;
			//BRep_Builder bb;
			//bb.MakeCompound(c);
			//for(int t = 0; t < edgeSet.size(); t++)
			//{
			//	bb.Add(c,edgeSet[t]);
			//}
			//BRepTools::Write(c,"C:\\makeface1.brep");
			//throw	c;
		}
		if (edgeSet.size() < 3 || edgeSet.size()>4) //20170601 by czb~{#,2;TJPm=vA=Lu1_99TlMx8qCf#,2;TJPm4sSZ~}4~{Lu1_99TlGzCf~}
		{
			continue;
		}
		
		//TopoDS_Face face = MakeFace(edgeSet); // DELETE BY XUEFENG 201810
		TopoDS_Shape face = MakeFace(edgeSet); // ADDED BY XUEFENG 201810
		if(face.IsNull())
		{
			continue;
		}
		//b.Add(shell, face);
		faceSet.push_back(face);
	}
	shell = SewShell(faceSet);

	this->result = shell;
	isDone = true;
}

// DELETE BY XUEFENG 201803
//TopoDS_Face GridFaces::MakeFace(vector<TopoDS_Edge> edgeSet)
// END DELETE

// ADDED BY XUEFENG 201803
TopoDS_Shape GridFaces::MakeFace(vector<TopoDS_Edge> edgeSet)
// END ADDED
{
	vector<Handle(Geom_BSplineCurve)> vctBSplineCurves;

	BRep_Builder b;
	TopoDS_Compound c;
	b.MakeCompound(c);

	for(int i = 0; i < edgeSet.size(); i++)
	{
		b.Add(c, edgeSet[i]);
	}
	//BRepTools::Write(c, "c:\\f.brep");
	vector<gp_Pnt> ps;
	for(int i = 0; i < edgeSet.size(); i++)
	{
		gp_Pnt p1,p2;
		TopoDS_Edge edge = edgeSet[i];
		GetEdgeStartEndPoint(edge, p1,p2);
		//CSplitAllInsectionCurves csp;
		//edge = csp.RebuildEdge(edge, p1,p2,false);
		//edge = BRepBuilderAPI_MakeEdge(p1,p2);
		Handle(Geom_BSplineCurve) curve = GetBSplineCurve(edge);

		//TColStd_Array1OfInteger M1(1, curve->NbKnots());
		//curve->Multiplicities(M1);

		vctBSplineCurves.push_back(curve);
		ps.push_back(p1);
		ps.push_back(p2);

		std::cout<<p1.X()<<"  "<<p1.Y()<<"  "<<p1.Z()<<std::endl;
		std::cout<<p2.X()<<"  "<<p2.Y()<<"  "<<p2.Z()<<std::endl;
	}
   std::cout<<std::endl<<std::endl;
   
   // DELETE BY XUEFENG 201803
   //TopoDS_Face face ;
    TopoDS_Shape face;  // ADDED BY XUEFENG 201803
   try
   {
	   // XUEFENG DELETE 20210319
	   //_ASSERTE( _CrtCheckMemory( ) );  //~{?IRT6(N;DZ4fP9B65DPP~}        //pj
	   face = MakeFace(vctBSplineCurves);
   }
   catch (Standard_Failure e)
   {
	   static int tt = 0;
	   TopoDS_Wire c;
	   BRep_Builder b;
	   b.MakeWire(c);
	   for(int i = 0; i < edgeSet.size(); i++)
	   {
		   //if(tt != i)
		   //	continue;
		   b.Add(c, edgeSet[i]);
	   }
	   tt++;
	   //BRepTools::Write(c,"C:\\makeface.brep");
	   //throw e.GetMessageString();
	   throw c;
   }
   return face;

	//BRep_Builder b;
	//TopoDS_Wire w;
	//b.MakeWire(w);
	//for(int i = 0; i < edgeSet.size(); i++)
	//{
	//	gp_Pnt p1,p2;
	//	GetEdgeStartEndPoint(edgeSet[i], p1,p2);
	//	b.Add(w, edgeSet[i]);
	//}
	//TopoDS_Face face = BRepBuilderAPI_MakeFace(w);
	//return face;
}

// XUEFENG ADDED FOR 5 CURVES FACE 201810
Handle(Geom_BSplineCurve) GridFaces::getBSplineCurveMultiPoints(gp_Pnt points[], int len)
{
	Handle(TColgp_HArray1OfPnt) SPLPoints = new TColgp_HArray1OfPnt(1,len);
	for(int i=1;i<=len;i++)
		SPLPoints->SetValue(i, points[i-1]);

	GeomAPI_Interpolate PtB(SPLPoints, Standard_False, 0.0);
	PtB.Perform();
	Handle(Geom_BSplineCurve) SPL;
	if(PtB.IsDone())
		SPL=PtB.Curve();

	return SPL;
}

void GridFaces::formFiveCurvesSeq(Handle(Geom_BSplineCurve) SPLs[], gp_Pnt startPoints[])
{
	startPoints[0] = SPLs[0]->StartPoint();

	Handle(Geom_BSplineCurve) finalSPLs[5];
	bool sig[5]={true,true,true,true,true};
	finalSPLs[0] = SPLs[0];
	sig[0] = false;
	gp_Pnt foreEnd = SPLs[0]->EndPoint();

	int n=1;
	while(n<5)
	{
		for(int i=0;i<5;i++)
		{
			if(sig[i]==true)
			{
				if(foreEnd.IsEqual(SPLs[i]->StartPoint(), 0.001)==true)
				{
					finalSPLs[n] = SPLs[i];
					sig[i] = false;
					foreEnd = SPLs[i]->EndPoint();

					startPoints[n] = finalSPLs[n]->StartPoint();

					break;
				}

				if(foreEnd.IsEqual(SPLs[i]->EndPoint(), 0.001)==true)
				{
					SPLs[i]->Reverse();
					finalSPLs[n] = SPLs[i];
					sig[i] = false;
					foreEnd = SPLs[i]->EndPoint();

					startPoints[n] = finalSPLs[n]->StartPoint();

					break;
				}
			}

		}

		n++;
	}

	for(int i=0;i<5;i++)
		SPLs[i] = finalSPLs[i];
}

int GridFaces::chooseSmoothPoint(gp_Pnt startPoints[], Handle(Geom_BSplineCurve) SPLs[], int divNum)
{
	int startPointSeq = 0;
	double minAngleSub = 999999999.0;

	for(int i=0;i<5;i++)
	{
		gp_Pnt tmpPnt[5];
		Handle(Geom_BSplineCurve) tmpSPLs[5];

		int num = 0;
		for(int j=0;j<5;j++)
		{
			num = i+j;
			if(num>4)
				num = i+j-5;

			tmpPnt[j]=startPoints[num];
			tmpSPLs[j]=SPLs[num];
		}

		double totalAngelSum = this->getAngleSum(tmpPnt, tmpSPLs, divNum);
		if(totalAngelSum<minAngleSub)
		{
			minAngleSub = totalAngelSum;
			startPointSeq = i;
		}
	}

	return startPointSeq;
}

double GridFaces::getAngleSum(gp_Pnt startPoints[], Handle(Geom_BSplineCurve) SPLs[], int divNum)
{
	double allDiveLinesAngles = 0.0;

	// start point is i;
	// 
	std::vector<gp_Pnt> pair1_Curve1_Points;	// SPLs[0]
	std::vector<gp_Pnt> pair1_Curve2_Points;	// SPLs[4]
	{
		double start = SPLs[0]->FirstParameter();
		double end = SPLs[0]->LastParameter();
		double tmpU = (end-start)/(double)divNum;
		for(int i=1;i<divNum;i++)
		{
			gp_Pnt tmpPnt = SPLs[0]->Value(start+((double)i)*tmpU);
			pair1_Curve1_Points.push_back(tmpPnt);
		}
	}
	{
		double start = SPLs[4]->FirstParameter();
		double end = SPLs[4]->LastParameter();
		double tmpU = (end-start)/(double)divNum;
		for(int i=1;i<divNum;i++)
		{
			gp_Pnt tmpPnt = SPLs[4]->Value(end - ((double)i)*tmpU);
			pair1_Curve2_Points.push_back(tmpPnt);
		}
	}

	gp_Dir foreDir = gp_Dir(pair1_Curve2_Points[0].X()-pair1_Curve1_Points[0].X(),
		pair1_Curve2_Points[0].Y()-pair1_Curve1_Points[0].Y(),
		pair1_Curve2_Points[0].Z()-pair1_Curve1_Points[0].Z());;
	for(int i=0;i<divNum-1;i++)
	{
		gp_Dir nowDir = gp_Dir(pair1_Curve2_Points[i].X()-pair1_Curve1_Points[i].X(),
			pair1_Curve2_Points[i].Y()-pair1_Curve1_Points[i].Y(),
			pair1_Curve2_Points[i].Z()-pair1_Curve1_Points[i].Z());
		allDiveLinesAngles+=nowDir.Angle(foreDir);
	}

	std::vector<gp_Pnt> pair2_Curve1_Points;	// SPLs[1]
	std::vector<gp_Pnt> pair2_Curve2_Points;	// SPLs[3]
	{
		double start = SPLs[1]->FirstParameter();
		double end = SPLs[1]->LastParameter();
		double tmpU = (end-start)/(double)divNum;
		for(int i=1;i<divNum;i++)
		{
			gp_Pnt tmpPnt = SPLs[1]->Value(start+((double)i)*tmpU);
			pair2_Curve1_Points.push_back(tmpPnt);
		}
	}

	{
		double start = SPLs[3]->FirstParameter();
		double end = SPLs[3]->LastParameter();
		double tmpU = (end-start)/(double)divNum;
		for(int i=1;i<divNum;i++)
		{
			gp_Pnt tmpPnt = SPLs[3]->Value(end - ((double)i)*tmpU);
			pair2_Curve2_Points.push_back(tmpPnt);
		}
	}
	for(int i=0;i<divNum-1;i++)
	{
		gp_Dir nowDir = gp_Dir(pair2_Curve2_Points[i].X()-pair2_Curve1_Points[i].X(),
			pair2_Curve2_Points[i].Y()-pair2_Curve1_Points[i].Y(),
			pair2_Curve2_Points[i].Z()-pair2_Curve1_Points[i].Z());
		allDiveLinesAngles+=nowDir.Angle(foreDir);
		foreDir = nowDir;
	}

	return allDiveLinesAngles;
}

//Handle(Geom_Surface) GridFaces::getSurfaceFromFiveCurves(Handle(Geom_BSplineCurve) SPLs[], int divNum, BRep_Builder myBuilder) // divNum>=2
//TopoDS_Shape GridFaces::getSurfaceFromFiveCurves(Handle(Geom_BSplineCurve) SPLs[], int divNum)
TopoDS_Shape GridFaces::getSurfaceFromFiveCurves(Handle(Geom_BSplineCurve) SPLs[], int divNum)
{
	//Handle(Geom_Surface) retSurf = NULL;

	gp_Pnt halfLineStartPoint = SPLs[0]->StartPoint();
	gp_Pnt halfLineEndPoint = SPLs[2]->Value((SPLs[2]->FirstParameter()+SPLs[2]->LastParameter())/2.0);
	gp_Pnt halfLinePoints[] = {halfLineStartPoint, halfLineEndPoint};
	Handle(Geom_BSplineCurve) halfLine =  this->getBSplineCurveMultiPoints(halfLinePoints, 2);
	//ShowRedShape(this->getEdgeFromPnt(myBuilder, halfLineStartPoint, halfLineEndPoint, halfLine));

	std::vector<gp_Pnt> pair1_Curve1_Points;	// SPLs[0]
	std::vector<gp_Pnt> pair1_Curve2_Points;	// SPLs[4]
	{
		double start = SPLs[0]->FirstParameter();
		double end = SPLs[0]->LastParameter();
		double tmpU = (end-start)/(double)divNum;
		for(int i=1;i<divNum;i++)
		{
			gp_Pnt tmpPnt = SPLs[0]->Value(start+((double)i)*tmpU);
			pair1_Curve1_Points.push_back(tmpPnt);
		}
	}

	{
		double start = SPLs[4]->FirstParameter();
		double end = SPLs[4]->LastParameter();
		double tmpU = (end-start)/(double)divNum;
		for(int i=1;i<divNum;i++)
		{
			gp_Pnt tmpPnt = SPLs[4]->Value(end - ((double)i)*tmpU);
			//this->ShowGreenShape(getVertexFromPnt(myBuilder, tmpPnt));
			pair1_Curve2_Points.push_back(tmpPnt);
		}
	}

	std::vector<Handle(Geom_BSplineCurve)> pair1_Curves;
	{
		for(int i=0;i<divNum-1;i++)
		{
			gp_Pnt connLinePoints[] = {pair1_Curve1_Points[i], pair1_Curve2_Points[i]};
			Handle(Geom_BSplineCurve) connLine =  this->getBSplineCurveMultiPoints(connLinePoints, 2);
			pair1_Curves.push_back(connLine);

			//ShowBlueShape(this->getEdgeFromPnt(myBuilder, pair1_Curve1_Points[i],  pair1_Curve2_Points[i], connLine));	
		}
	}

	std::vector<gp_Pnt> pair2_Curve1_Points;	// SPLs[1]
	std::vector<gp_Pnt> pair2_Curve2_Points;	// SPLs[3]
	{
		double start = SPLs[1]->FirstParameter();
		double end = SPLs[1]->LastParameter();
		double tmpU = (end-start)/(double)divNum;
		for(int i=1;i<divNum;i++)
		{
			gp_Pnt tmpPnt = SPLs[1]->Value(start+((double)i)*tmpU);
			pair2_Curve1_Points.push_back(tmpPnt);
		}
	}

	{
		double start = SPLs[3]->FirstParameter();
		double end = SPLs[3]->LastParameter();
		double tmpU = (end-start)/(double)divNum;
		for(int i=1;i<divNum;i++)
		{
			gp_Pnt tmpPnt = SPLs[3]->Value(end - ((double)i)*tmpU);
			//this->ShowGreenShape(getVertexFromPnt(myBuilder, tmpPnt));
			pair2_Curve2_Points.push_back(tmpPnt);
		}
	}

	std::vector<Handle(Geom_BSplineCurve)> pair2_Curves;
	{
		for(int i=0;i<divNum-1;i++)
		{
			gp_Pnt connLinePoints[] = {pair2_Curve1_Points[i], pair2_Curve2_Points[i]};
			Handle(Geom_BSplineCurve) connLine =  this->getBSplineCurveMultiPoints(connLinePoints, 2);
			pair2_Curves.push_back(connLine);

			//ShowBlueShape(this->getEdgeFromPnt(myBuilder, pair2_Curve1_Points[i],  pair2_Curve2_Points[i], connLine));
		}
	}

	GeomAPI_ExtremaCurveCurve ecc;
	std::vector<gp_Pnt> divPair1Points;
	std::vector<gp_Pnt> divPair2Points;
	for(int i=0;i<divNum-1;i++)
	{
		gp_Pnt halfPoint, pair1Point, pair2Point;
		ecc.Init(halfLine, pair1_Curves[i]);
		ecc.NearestPoints(halfPoint, pair1Point);

		divPair1Points.push_back(pair1Point);
		//this->ShowGreenShape(getVertexFromPnt(myBuilder, pair1Point));

		ecc.Init(halfLine, pair2_Curves[i]);
		ecc.NearestPoints(halfPoint, pair2Point);
	
		divPair2Points.push_back(pair2Point);
		//this->ShowGreenShape(getVertexFromPnt(myBuilder, pair2Point));
	}

	std::vector<gp_Pnt> finalHalfCurvePoints;
	finalHalfCurvePoints.push_back(halfLineStartPoint);
	for(int i=0;i<divNum-1;i++)
		finalHalfCurvePoints.push_back(divPair1Points[i]);
	for(int i=0;i<divNum-1;i++)
		finalHalfCurvePoints.push_back(divPair2Points[i]);
	finalHalfCurvePoints.push_back(halfLineEndPoint);

	Handle(Geom_BSplineCurve) finalHalfCurve = this->getBSplineCurveFromVectorPnt(finalHalfCurvePoints, 2+2*(divNum-1));
	//ShowRedShape(this->getEdgeFromPnt(myBuilder, halfLineStartPoint, halfLineEndPoint, finalHalfCurve));

	Handle(Geom_BSplineCurve) halfSPL2ToSPL01; 
	Handle(Geom_BSplineCurve) halfSPL2ToSPL34; 
	{
		Handle(Geom_BSplineCurve) divCurveFor01 = Handle(Geom_BSplineCurve)::DownCast(SPLs[2]->Copy());
		Handle(Geom_BSplineCurve) divCurveFor34 = Handle(Geom_BSplineCurve)::DownCast(SPLs[2]->Copy());
		double startU = divCurveFor01->FirstParameter();
		double endU = divCurveFor01->LastParameter();
		double midU = (startU+endU)/2.0;
		divCurveFor01->Segment(startU, midU);
		divCurveFor34->Segment(midU, endU);
		halfSPL2ToSPL01 = divCurveFor01;
		halfSPL2ToSPL34 = divCurveFor34;
	}

	Handle(Geom_BSplineCurve) halfFace1SPL[] = {finalHalfCurve, SPLs[0], SPLs[1], halfSPL2ToSPL01};

	Handle(Geom_Surface) halfSurface1 = this->getBSplineSurface(halfFace1SPL);
	
	//BRepBuilderAPI_MakeFace MF1( halfSurface1 );      // XUEFENG DELETE 202101 
	BRepBuilderAPI_MakeFace MF1(halfSurface1, 0.0001);  // XUEFENG ADDED 202101 OCC 7.1.0

	TopoDS_Face halfBSplineFace1 = MF1.Face();
	//ShowShape(MF1.Face());

	Handle(Geom_BSplineCurve) halfFace2SPL[] = {finalHalfCurve, SPLs[4], SPLs[3], halfSPL2ToSPL34};
	Handle(Geom_Surface) halfSurface2 = this->getBSplineSurface(halfFace2SPL);
	
	//BRepBuilderAPI_MakeFace MF2( halfSurface2 );   // XUEFENG DELETE 202009 6.3.0
	BRepBuilderAPI_MakeFace MF2(halfSurface2, 0.0001);  // XUEFENG ADDED 202001 7.1.0

	TopoDS_Face halfBSplineFace2 = MF2.Face();
	//ShowShape(MF2.Face());

	BRepBuilderAPI_Sewing sewMethod(0.1);
	sewMethod.Add(halfBSplineFace1);
	sewMethod.Add(halfBSplineFace2);
	sewMethod.Perform();
	TopoDS_Shape retShape = sewMethod.SewedShape();
	//TopoDS_Shell retShell =  TopoDS::Shell( retShape );
	//TopoDS_Face retFace =  TopoDS::Face( retShape );
	//BRepOffsetAPI_Sewing aMethodSew;
	//aMethodSew.SetTolerance(0.1);
	//aMethodSew.Add(halfBSplineFace1);
	//aMethodSew.Add(halfBSplineFace2);
	//aMethodSew.Perform();
	//TopoDS_Shape retShape =  aMethodSew.SewedShape();
	//TopoDS_Face retFace =  TopoDS::Face( aMethodSew.SewedShape() );

	return retShape;
}

Handle(Geom_Surface) GridFaces::getBSplineSurface(Handle(Geom_BSplineCurve) SPL[])
{
	GeomFill_FillingStyle fillType = GeomFill_StretchStyle;
	GeomFill_BSplineCurves aGeomFill(SPL[0],SPL[1],SPL[2],SPL[3], fillType);
	Handle(Geom_BSplineSurface) aBSplineSurface = aGeomFill.Surface();
	Handle(Geom_Surface) aSurf = aBSplineSurface;

	return aSurf;
}

Handle(Geom_BSplineCurve) GridFaces::getBSplineCurveFromVectorPnt(std::vector<gp_Pnt> points, int len)
{
	TColgp_Array1OfPnt SPLPoints(1,len);
	for(int i=1;i<=len;i++)
		SPLPoints(i)=points[i-1];

	GeomAPI_PointsToBSpline aPointToSPL(SPLPoints);
	Handle(Geom_BSplineCurve) SPL = aPointToSPL.Curve();

	return SPL;
}
// END

// CZB ADD 20200602
// ADDED BY XUEFENG 20200416
#include "GeomFill_ConstrainedFilling.hxx"
//#include "GeomAdaptor_HCurve.hxx"  // NOT EXIST IN 7.6.0
#include "Adaptor3d_Curve.hxx"		// Added for OCC 7.6.0 XUEFENG 2024.01
#include "GeomFill_SimpleBound.hxx"
// END ADDED
// CZB END ADDED 20200602

//TopoDS_Face GridFaces::MakeFace(vector<Handle(Geom_BSplineCurve)> vctBSplineCurves)
// END DELETE

// ADDED BY XUEFENG 201803
TopoDS_Shape GridFaces::MakeFace(vector<Handle(Geom_BSplineCurve)> vctBSplineCurves)
//END ADDED
{
	//TopoDS_Face face;		//DELETE BY XUEFENG 201803
	TopoDS_Shape face;

	GeomFill_BSplineCurves fill; 
	bool isFill = false;



	//Standard_Integer MaxDeg = 8;                                                               
	//Standard_Integer MaxSeg = 1000;                                                               
	//GeomFill_ConstrainedFilling aConstrainedFilling(MaxDeg, MaxSeg);          
	try
	{
		static int ii = 0;
		if(vctBSplineCurves.size() == 2)
		{
			Handle(Geom_BSplineCurve)& bc1 = vctBSplineCurves[0];
			Handle(Geom_BSplineCurve)& bc2 = vctBSplineCurves[1];
			TColgp_Array1OfPnt      Poles(1,2);
			TColStd_Array1OfReal    Knots(1,2);
			TColStd_Array1OfInteger Mults(1,2);
			Poles( 1) = bc1->StartPoint();
			Poles( 2) = bc1->StartPoint();
			Knots( 1) = bc2->Knot(bc2->FirstUKnotIndex());
			Knots( 2) = bc2->Knot(bc2->LastUKnotIndex());
			Mults( 1) = Mults( 2) = 2;
			Handle(Geom_BSplineCurve)  bc3 = new Geom_BSplineCurve( Poles, Knots, Mults, 1);
			gp_Pnt p1 = bc1->StartPoint();
			gp_Pnt p2 = bc1->EndPoint();
			gp_Pnt p3 = bc2->StartPoint();
			gp_Pnt p4 = bc2->EndPoint();
			if(bc1->StartPoint().Distance(bc3->StartPoint()) < Precision::Confusion() || bc1->StartPoint().Distance(bc3->EndPoint()) < Precision::Confusion())
			{
				bc1->Reverse();
			}
			fill=GeomFill_BSplineCurves(bc1,bc3,bc2,GeomFill_StretchStyle);
			//fill=GeomFill_BSplineCurves(vctBSplineCurves[0],vctBSplineCurves[1],GeomFill_StretchStyle);
			Handle(Geom_BSplineSurface) BSS = fill.Surface();
			
			//face = BRepBuilderAPI_MakeFace(BSS);  // XUEFENG DELETE 202101
			face = BRepBuilderAPI_MakeFace(BSS, Precision::Confusion()); // XUEFENG ADDED 202101
			
			isFill = true;
		}

		else if(vctBSplineCurves.size() == 3)
		{
			// XUEFENG ADDED 20230827 
			// MODIFIED BACK TO OLD VERSION
			Handle(Geom_BSplineCurve)& c1 = vctBSplineCurves[0];
			Handle(Geom_BSplineCurve)& c2 = vctBSplineCurves[1];

			if(c1->StartPoint().Distance(c2->StartPoint()) < Precision::Confusion() || c1->StartPoint().Distance(c2->EndPoint()) < Precision::Confusion())
			{
				c1->Reverse();
			}
			//TriangleEdge(vctBSplineCurves);
			fill=GeomFill_BSplineCurves(vctBSplineCurves[0],vctBSplineCurves[1],vctBSplineCurves[2],GeomFill_StretchStyle);
			face = BRepBuilderAPI_MakeFace(fill.Surface(), Precision::Confusion()*100);
			isFill = true;
			
			// END ADDED TEST

			// XUEFENG ADDED 20230802
			// RETURN BACK TO SIMPLE WAY
			/*
			fill=GeomFill_BSplineCurves(vctBSplineCurves[0],vctBSplineCurves[1],vctBSplineCurves[2],GeomFill_StretchStyle);
			isFill = true;
			if(isFill)
			{
				Handle(Geom_BSplineSurface) BSS = fill.Surface();
				//face = BRepBuilderAPI_MakeFace(BSS);  // XUEFENG DELETE 202101
				face = BRepBuilderAPI_MakeFace(BSS, Precision::Confusion());   // XUEFENG ADDED 202101
			}
			*/
			//TriangleEdge(vctBSplineCurves);
			// XUEFENG DELETE 20200421
			/*
			Handle(Geom_BSplineCurve)& c1 = vctBSplineCurves[0];
			Handle(Geom_BSplineCurve)& c2 = vctBSplineCurves[1];

			if(c1->StartPoint().Distance(c2->StartPoint()) < Precision::Confusion() || c1->StartPoint().Distance(c2->EndPoint()) < Precision::Confusion())
			{
				c1->Reverse();
			}
			*/

			//TriangleEdge(vctBSplineCurves);
			
			// XUEFENG DELETE 20200421
			//fill=GeomFill_BSplineCurves(vctBSplineCurves[0],vctBSplineCurves[1],vctBSplineCurves[2],vctBSplineCurves[3],GeomFill_StretchStyle);
			//face = BRepBuilderAPI_MakeFace(fill.Surface());
			// END DELETE
			
			/* DELETE 20230801
			// CZB ADD 20200602
			// XUEFENG ADDED 20200421
			// 
			Standard_Integer MaxDeg = 8;                                                               
			Standard_Integer MaxSeg = 1000;                                                               
			GeomFill_ConstrainedFilling aConstrainedFilling(MaxDeg, MaxSeg);
			Handle(Geom_BSplineCurve) SPL1 = vctBSplineCurves[0];                  
			Handle(Geom_BSplineCurve) SPL2 = vctBSplineCurves[1];                  
			Handle(Geom_BSplineCurve) SPL3 = vctBSplineCurves[2];                  

			Handle(GeomAdaptor_HCurve) SPL1Adaptor = new GeomAdaptor_HCurve(SPL1);                     
			Handle(GeomFill_SimpleBound) B1 =                                                          
				new GeomFill_SimpleBound(SPL1Adaptor,Precision::Approximation(),Precision::Angular());  
			Handle(GeomAdaptor_HCurve) SPL2Adaptor = new GeomAdaptor_HCurve(SPL2);                     
			Handle(GeomFill_SimpleBound) B2 =                                                          
				new GeomFill_SimpleBound(SPL2Adaptor,Precision::Approximation(),Precision::Angular());  
			Handle(GeomAdaptor_HCurve) SPL3Adaptor = new GeomAdaptor_HCurve(SPL3);                     
			Handle(GeomFill_SimpleBound) B3 =                                                          
				new GeomFill_SimpleBound(SPL3Adaptor,Precision::Approximation(),Precision::Angular());  
			
			Standard_Boolean NoCheck= Standard_False;

			// XUEFENG ADDED TEST 20211029
			int knotsNum1 = SPL1->NbKnots();
			int polesNum1 = SPL1->NbPoles();
			int knotsNum2 = SPL2->NbKnots();
			int polesNum2 = SPL2->NbPoles();
			int knotsNum3 = SPL3->NbKnots();
			int polesNum3 = SPL3->NbPoles();
			*/
			// DELETE 20230801

			/*
			if(knotsNum1>=knotsNum2 && knotsNum1>=knotsNum3)
			{
				if(polesNum2>polesNum3)
					aConstrainedFilling.Init(B3,B1,B2,NoCheck); 
				else// if(polesNum2<=polesNum3)
					aConstrainedFilling.Init(B2,B1,B3,NoCheck); 

			}
			else if(knotsNum2>=knotsNum1 && knotsNum2>=knotsNum3)
			{
				if(polesNum1>polesNum3)
					aConstrainedFilling.Init(B3,B2,B1,NoCheck); 
				else// if(polesNum1<=polesNum3)
					aConstrainedFilling.Init(B1,B2,B3,NoCheck); 

				//aConstrainedFilling.Init(B1,B2,B3,NoCheck); 
			}
			else if(knotsNum3>=knotsNum1 && knotsNum3>=knotsNum2)
			{
				if(polesNum1>polesNum2)
					aConstrainedFilling.Init(B2,B3,B1,NoCheck); 
				else// if(polesNum1<=polesNum2)
					aConstrainedFilling.Init(B1,B3,B2,NoCheck); 

				//aConstrainedFilling.Init(B2,B3,B1,NoCheck); 
			}
			else
			*/
	 
			/* DELETE 20230801
			if(knotsNum1>=knotsNum2 && knotsNum1>=knotsNum3)
			{
				if(polesNum2>polesNum3)
					aConstrainedFilling.Init(B3,B2,B1,NoCheck); 
				else// if(polesNum2<=polesNum3)
					aConstrainedFilling.Init(B2,B3,B1,NoCheck); 

			}
			else if(knotsNum2>=knotsNum1 && knotsNum2>=knotsNum3)
			{
				if(polesNum1>polesNum3)
					aConstrainedFilling.Init(B3,B1,B2,NoCheck); 
				else// if(polesNum1<=polesNum3)
					aConstrainedFilling.Init(B1,B3,B2,NoCheck); 

				//aConstrainedFilling.Init(B1,B2,B3,NoCheck); 
			}
			else if(knotsNum3>=knotsNum1 && knotsNum3>=knotsNum2)
			{
				if(polesNum1>polesNum2)
					aConstrainedFilling.Init(B2,B1,B3,NoCheck); 
				else// if(polesNum1<=polesNum2)
					aConstrainedFilling.Init(B1,B2,B3,NoCheck); 

				//aConstrainedFilling.Init(B2,B3,B1,NoCheck); 
			}
			else
			
			// END ADDED
			aConstrainedFilling.Init(B1,B2,B3,NoCheck);                                                

			
			Handle(Geom_BSplineSurface) aBSplineSurface = aConstrainedFilling.Surface();               
			//face = BRepBuilderAPI_MakeFace(aBSplineSurface);  // XUEFENG DELETE 202009  OCC 6.3.0
			
			face = BRepBuilderAPI_MakeFace(aConstrainedFilling.Surface(), Precision::Confusion());  // XUEFENG ADDED 202009 OCC 7.1.0

			isFill = true;
			*/
			// DELETE 20230801
			// END ADDED
			// CZB END ADDED 20200602

			//END DELETE COMPLEX THE 3 LINES FACE 20230801
		}
		else if(vctBSplineCurves.size() == 4)
		{
			// ~{4K4&Hg:NIhVC;rP^8DTKKc>+6H#?~}
			//fill=GeomFill_BSplineCurves(vctBSplineCurves[0],vctBSplineCurves[1],vctBSplineCurves[2],vctBSplineCurves[3],GeomFill_StretchStyle);
			//face = BRepBuilderAPI_MakeFace(fill.Surface()); // XUEFENG DELETE 202101
			//face = BRepBuilderAPI_MakeFace(fill.Surface(), Precision::Confusion());   // XUEFENG ADDED 202101

			// IMPORTANT!!!! MODIFIED BACK AGAIN!
			// XUEFENG MODIFIED 20210420 AGAIN
			// XUEFENG TEST

			Standard_Integer MaxDeg = 8;                                                               
			Standard_Integer MaxSeg = 1000;                                                               
			GeomFill_ConstrainedFilling aConstrainedFilling(MaxDeg, MaxSeg);
			Handle(Geom_BSplineCurve) SPL1 = vctBSplineCurves[0];                  
			Handle(Geom_BSplineCurve) SPL2 = vctBSplineCurves[1];                  
			Handle(Geom_BSplineCurve) SPL3 = vctBSplineCurves[2];                  
			Handle(Geom_BSplineCurve) SPL4 = vctBSplineCurves[3];         

			//Handle(GeomAdaptor_HCurve) SPL1Adaptor = new GeomAdaptor_HCurve(SPL1);  // NOT EXIST IN 7.6.0
			Handle(GeomAdaptor_Curve) SPL1Adaptor = new GeomAdaptor_Curve(SPL1);  // Added for OCC 7.6.0

			Handle(GeomFill_SimpleBound) B1 =                                                          
				new GeomFill_SimpleBound(SPL1Adaptor,Precision::Approximation(),Precision::Angular());  

			//Handle(GeomAdaptor_HCurve) SPL2Adaptor = new GeomAdaptor_HCurve(SPL2);    // NOT EXIST IN 7.6.0
			Handle(GeomAdaptor_Curve) SPL2Adaptor = new GeomAdaptor_Curve(SPL2);  // Added for OCC 7.6.0   

			Handle(GeomFill_SimpleBound) B2 =                                                          
				new GeomFill_SimpleBound(SPL2Adaptor,Precision::Approximation(),Precision::Angular());  

			//Handle(GeomAdaptor_HCurve) SPL3Adaptor = new GeomAdaptor_HCurve(SPL3); // NOT EXIST IN 7.6.0
			Handle(GeomAdaptor_Curve) SPL3Adaptor = new GeomAdaptor_Curve(SPL3);  // Added for OCC 7.6.0                    
			Handle(GeomFill_SimpleBound) B3 =                                                          
				new GeomFill_SimpleBound(SPL3Adaptor,Precision::Approximation(),Precision::Angular());  

			//Handle(GeomAdaptor_HCurve) SPL4Adaptor = new GeomAdaptor_HCurve(SPL4);   // NOT EXIST IN 7.6.0
			Handle(GeomAdaptor_Curve) SPL4Adaptor = new GeomAdaptor_Curve(SPL4);  // Added for OCC 7.6.0                  
			Handle(GeomFill_SimpleBound) B4 =                                                          
				new GeomFill_SimpleBound(SPL4Adaptor,Precision::Approximation(),Precision::Angular());  

			Standard_Boolean NoCheck= Standard_False;     
			aConstrainedFilling.Init(B1,B2,B3,B4,NoCheck);                                                

			Handle(Geom_BSplineSurface) aBSplineSurface = aConstrainedFilling.Surface();               
			face = BRepBuilderAPI_MakeFace(aBSplineSurface, Precision::Approximation());

			isFill = true;
			// END MODIFIED!
	}
	// XUEFENG ADDED 201810
	else if(vctBSplineCurves.size() == 5)
	{
		Handle(Geom_BSplineCurve) inputBSPL[5] = { vctBSplineCurves[0],vctBSplineCurves[1],vctBSplineCurves[2],vctBSplineCurves[3],vctBSplineCurves[4]}; 
		// Form the BSplines Curves Seq
		gp_Pnt startPoints[5];
		this->formFiveCurvesSeq(inputBSPL,startPoints);

		int chosenPoint = this->chooseSmoothPoint(startPoints, inputBSPL, 8);

		gp_Pnt finalPoints[5];
		Handle(Geom_BSplineCurve) finalBSPL[5];

		int finalPos=0;
		for(int i=0;i<5;i++)
		{
			finalPos = i+chosenPoint;
			if(finalPos>4)
				finalPos = i+chosenPoint-5;
			finalPoints[i] = startPoints[finalPos];
			finalBSPL[i] = inputBSPL[finalPos];
		}
		//this->getSurfaceFromFiveCurves(finalBSPL, 8, myBuilder);
		face = this->getSurfaceFromFiveCurves(finalBSPL, 8);
		//TopoDS_Face halfBSplineFace2 = MF2.Face();

		// End Face
		isFill = true;
	}
	// END ADDED

	int i = ii;
	if(i == 202)
	{
		i = i;
	}
	}
	catch (Standard_Failure e)
	{
		string myError = e.GetMessageString();
		//throw myError; // XUEFENG DELETE 20230827
		//throw "~{GzCfF,Iz3IRl3##,GkH7HOMx8qO_7b1U#!~}";
	}
	return face;
}

map<int, vector<TopoDS_Edge>> GridFaces::GetWires()
{
	return wireSet;
	//return wireOldSet;
}