/*	sample_surface.C

		last modified on Jan/24/2000
*/

#include	<iostream.h>
#include	<stdio.h>
#include	<nurbsS.hh>

int main(int argc,char *argv[])
{
	using namespace PLib ; 

	char	chx[100];
	double	pos[3];
	int	u_deg;
	int	v_deg;
	int	u_num_poles;
	int	v_num_poles;
	int	u_num_knots;
	int	v_num_knots;

	int	nn_x;
	int	nn_y;

	int	it;
	int	jt;

	cerr << "Testing NURBS interpolation and approximation\n" ;

	cin.getline	(chx,90);
	cin.getline	(chx,90);
	sscanf	(chx,"%d%d",&nn_x,&nn_y);

	cerr << "nn_x = " << nn_x << "\tnn_y = " << nn_y << '\n';

	u_deg = 3;
	v_deg = 3;

	u_num_poles = nn_x / 3;

	if (u_num_poles < (u_deg +1))
		{
		u_num_poles = u_deg + 1;
		}

	v_num_poles = nn_y / 3;

	if (v_num_poles < (v_deg +1))
		{
		v_num_poles = v_deg + 1;
		}

	cerr << "u_num_poles = " << u_num_poles << '\t';
	cerr << "v_num_poles = " << v_num_poles << '\n';

	u_num_knots = u_deg + 1 + u_num_poles;
	v_num_knots = v_deg + 1 + v_num_poles;

	Matrix_Point3Df points(nn_x,nn_y) ;

	for (jt=0; jt< nn_y; jt++)
		for (it=0; it<nn_x; it++)
		{
		cin.getline	(chx,90);

		sscanf	(chx,"%lf%lf%lf",&pos[0],&pos[1],&pos[2]);
		points(it,jt).x() = pos[0];
		points(it,jt).y() = pos[1];
		points(it,jt).z() = pos[2];

		}

	PlNurbsSurfacef surfA;

	surfA.leastSquares (points,u_deg,v_deg,u_num_poles,v_num_poles) ;

	cerr << "D1 (degreeU) = " << surfA.degreeU() << endl ;
	cerr << "D1 (degreeV) = " << surfA.degreeV() << endl ;

	cout << "sur	20\n";

	cout << surfA.degreeU() +1 << '\t';
	cout	<< u_num_poles << '\t';
	cout	<< u_num_knots << '\n';

	cout << surfA.degreeV() +1 << '\t';
	cout	<< v_num_poles << '\t';
	cout	<< v_num_knots << '\n';

	for (it=0; it< u_num_knots; it++)
		{
		cout << surfA.knotU(it) << endl;
		}

	for (it=0; it< v_num_knots; it++)
		{
		cout << surfA.knotV(it) << endl;
		}

	for (jt=0; jt< v_num_poles; jt++)
		for (it=0; it< u_num_poles; it++)
			{
		 	cout << surfA.ctrlPnts()(it,jt) << endl ;
			}

	cerr << "Done with the interpolation and approximation test.\n" ;
}