Debbuging the load and save for ALTA format in the case of RF

sources/core/plugins_manager.cpp

@@ -416,25 +416,33 @@ fitter* plugins_manager::get_fitter(const std::string& n)
 void plugins_manager::check_compatibility(data*& d, function*& f,
                                           const arguments& args)
 {
-	if(d->input_parametrization() == params::UNKNOWN_INPUT)
+	if(d->input_parametrization() == params::UNKNOWN_INPUT &&
+		f->input_parametrization() == params::UNKNOWN_INPUT)
 	{
-		std::cout << "<<WARNING>> unknown parametrization for data" << std::endl;
-	}
-
-	if(f->input_parametrization() == params::UNKNOWN_INPUT)
-	{
-		std::cout << "<<DEBUG>> function will take the parametrization of the data" << std::endl;
-		f->setParametrization(d->input_parametrization());
-	}
-	else if(d->input_parametrization() != f->input_parametrization())
-	{
-		std::cout << "<<INFO>> has to change the parametrization of the input data" << std::endl;
-		data_params* dd = new data_params(d, f->input_parametrization());
-		d = dd ;
+		std::cout << "<<WARNING>> both function and data objects have no parametrization" << std::endl;
 	}
 	else
 	{
-		std::cout << "<<DEBUG>> no change was made to the parametrization" << std::endl;
+		if(d->input_parametrization() == params::UNKNOWN_INPUT)
+		{
+			std::cout << "<<WARNING>> unknown parametrization for data" << std::endl;
+		}
+
+		if(f->input_parametrization() == params::UNKNOWN_INPUT)
+		{
+			std::cout << "<<DEBUG>> function will take the parametrization of the data" << std::endl;
+			f->setParametrization(d->input_parametrization());
+		}
+		else if(d->input_parametrization() != f->input_parametrization())
+		{
+			std::cout << "<<INFO>> has to change the parametrization of the input data" << std::endl;
+			data_params* dd = new data_params(d, f->input_parametrization());
+			d = dd ;
+		}
+		else
+		{
+			std::cout << "<<DEBUG>> no change was made to the parametrization" << std::endl;
+		}
 	}
 
 	if(f->dimY() != d->dimY())

sources/core/rational_function.cpp

@@ -366,105 +366,89 @@ vec rational_function::value(const vec& x) const
 }
 
 // IO function to text files
-void rational_function::load(std::istream& file)
+void rational_function::load(std::istream& in)
 {
-	int nX, nY ;
-	vec xmin, xmax ;
-	vec a, b;
-	int i = 0, j = 0;
-	while(file.good())
-	{
-		std::string line ;
-		std::getline(file, line) ;
-		std::stringstream linestream(line) ;
-
-		// Discard incorrect lines
-		if(linestream.peek() == '#')
-		{
-			linestream.ignore(1) ;
 
-			std::string comment ;
-			linestream >> comment ;
+	// Parse line until the next comment
+	while(in.peek() != '#')
+	{
+		char line[256];
+		in.getline(line, 256);
+	}
 
-			if(comment == std::string("DIM"))
-			{
-				linestream >> nX >> nY ;
-				setDimX(nX) ;
-				setDimY(nY) ;
+    // Checking for the comment line #FUNC nonlinear_function_lafortune
+	std::string token;
+	in >> token;
+	if(token.compare("#FUNC") != 0) 
+	{ 
+		std::cerr << "<<ERROR>> parsing the stream. The #FUNC is not the next line defined." << std::endl; 
+	}
 
-				xmin.resize(nX) ;
-				xmax.resize(nX) ;
-				for(int k=0; k<nX; ++k)
-					xmax[k] = 1.0;
+	in >> token;
+   if(token.compare("rational_function") != 0) 
+	{
+		std::cerr << "<<ERROR>> parsing the stream. function name is not the next token." << std::endl; 
+	}
 
-				setMin(xmin) ;
-				setMax(xmax) ;
-			}
-			else if(comment == std::string("NP"))
-			{
-				linestream >> np ;
-				a.resize(np);
-				std::cout << "<<DEBUG>> loading RF with np = " << np << std::endl;
-			}
-			else if(comment == std::string("NQ"))
-			{
-				linestream >> nq ;
-				b.resize(nq);
-				std::cout << "<<DEBUG>> loading RF with nq = " << nq << std::endl;
-			}
-			else if(comment == std::string("INPUT_PARAM"))
-			{
-				std::string param;
-				linestream >> param ;
+	int _np, _nq;
+	// Shoudl have the #NP [int]
+	in >> token >> _np;
+	
+	// Shoudl have the #NQ [int]
+	in >> token >> _nq;
+	setSize(_np, _nq);
 
-				setParametrization(params::parse_input(param));
-			}
-			continue ;
-		}
-		else if(line.empty())
-		{
-			continue ;
-		}
-		else if(j < nY)
-		{
-			int index ; double val ;
+	// Check for the MIN and MAX vector
+	vec min(dimX()), max(dimX());
+	in >> token;
+   if(token.compare("#MIN") != 0) 
+	{
+		std::cerr << "<<ERROR>> the min value for the input space is not defined." << std::endl; 
+	}
+	for(int k=0; k<dimX(); ++k) {in >> min[k];}
+	setMin(min);
 
-			// Accessing the index
-			for(int k=0; k<nX; ++k) {
-				linestream >> index ;
-			}
+	in >> token;
+   if(token.compare("#MAX") != 0) 
+	{
+		std::cerr << "<<ERROR>> the max value for the input space is not defined." << std::endl; 
+	}
+	for(int k=0; k<dimX(); ++k) {in >> max[k]; }
+	setMax(max);
 
-			// Accessing the value
-			linestream >> val ;
-			if(i < np)
-			{
-				a[i] = val ;
-			}
-			else
-			{
-				b[i-np] = val ;
-			}
 
+	// Check for the polynomial basis type
+	in >> token;
+   if(token.compare("#BASIS") != 0) 
+	{
+		std::cerr << "<<ERROR>> the file is not specifying the polynomial basis." << std::endl; 
+	}
+	in >> token;
+   if(token.compare("LEGENDRE") != 0) 
+	{
+		std::cerr << "<<ERROR>> the basis is different than LEGENDRE." << std::endl; 
+	}
 
-			// Update the output dimension number
-			if(i == np+nq-1) 
-			{
-				i = 0 ;
-				get(j)->update(a, b);
+	vec a(_np), b(_nq);
+	for(int i=0; i<_nY; ++i)
+	{
+		// Parse the p_i coefficients
+		for(int j=0; j<_np; ++j)
+		{
+			in >> token >> a[j];
+		}
+		
+		// Parse the q_i coefficients
+		for(int j=0; j<_nq; ++j)
+		{
+			in >> token >> b[j];
+		}
 
-#ifdef DEBUG
-				std::cout << "<<DEBUG>> loading channel " << j << " with: " << std::endl;
-				std::cout << "          " << a << std::endl;
-				std::cout << "          " << b << std::endl;
-#endif
+		std::cout << a << std::endl;
+		std::cout << b << std::endl;
 
-				++j ;
-			} 
-			else 
-			{
-				++i ;
-			}
-		}
+		// Update the i_th color channel
+		get(i)->update(a, b);
 	}
 }
 
@@ -698,6 +682,42 @@ void rational_function::save_gnuplot(const std::string& filename, const data* d,
 	file.close();
 }
 
+void rational_function::save_call(std::ostream& out, const arguments& args) const
+{
+	out << "#FUNC rational_function" << std::endl;
+	out << "#NP " << np << std::endl ;
+	out << "#NQ " << nq << std::endl ;
+	out << "#MIN "; for(int k=0; k<_nX; ++k) { out << _min[k] << " "; } out << std::endl;
+	out << "#MAX "; for(int k=0; k<_nX; ++k) { out << _max[k] << " "; } out << std::endl; 
+	out << "#BASIS LEGENDRE" << std::endl ;
+
+	for(int k=0; k<_nY; ++k)
+	{
+		rational_function_1d* rf = get(k);
+		vec a = rf->getP();
+		vec b = rf->getQ();
+
+		for(unsigned int i=0; i<np; ++i)
+		{
+			std::vector<int> index = rf->index2degree(i) ;
+			for(unsigned int j=0; j<index.size(); ++j)
+			{
+				out << index[j] << "\t" ;
+			}
+			out << a[i] << std::endl ;
+		}
+
+		for(unsigned int i=0; i<nq; ++i)
+		{
+			std::vector<int> index = rf->index2degree(i) ;
+			for(unsigned int j=0; j<index.size(); ++j)
+			{
+				out << index[j] << "\t" ;
+			}
+			out << b[i] << std::endl ;
+		}
+	}
+}
 
 void rational_function::save(const std::string& filename) const
 {

sources/core/rational_function.h

@@ -130,6 +130,27 @@ class rational_function : public function
 			rs.resize(_nY);
 		}
 
+		virtual void setMin(const vec& min)
+		{
+			function::setMin(min);
+			for(int i=0; i<dimY(); ++i)
+			{
+				get(i)->setMin(min);
+			}
+		}
+		
+		virtual void setMax(const vec& max)
+		{
+			function::setMax(max);
+			for(int i=0; i<dimY(); ++i)
+			{
+				get(i)->setMax(max);
+			}
+		}
+
+		//! \brief Save the rational function to the rational format (see \ref formating).
+		virtual void save_call(std::ostream& out, const arguments& args) const ;
+
 	protected: // functions
 
 		//! \brief Save the rational function to the rational format (see \ref formating).

sources/plugins/rational_fitter_matlab/rational_fitter.cpp

@@ -57,6 +57,8 @@ bool rational_fitter_matlab::fit_data(const data* dat, function* fit, const argu
 	{
 		QTime time ;
 		time.start() ;
+
+		r->setSize(temp_np, temp_nq);
 		
 		if(fit_data(d, temp_np, temp_nq, r))
 		{

sources/plugins/rational_fitter_quadprog/rational_fitter.cpp

@@ -57,6 +57,7 @@ bool rational_fitter_quadprog::fit_data(const data* dat, function* fit, const ar
 		QTime time ;
 		time.start() ;
 		
+		r->setSize(temp_np, temp_nq);
 		if(fit_data(d, temp_np, temp_nq, r))
 		{
 			int msec = time.elapsed() ;

sources/softs/data2moments/main.cpp

@@ -15,7 +15,8 @@
 #include <cstdlib>
 #include <cmath>
 
-/*
+#define  EPSILON 1.0E-5
+
 vec coord(vec V, vec L, vec X, vec Y, vec N)
 {
 	vec pV = V-dot(V,N)*N;
@@ -25,19 +26,24 @@ vec coord(vec V, vec L, vec X, vec Y, vec N)
 	vCoord /= (1.0+dot(V,N));
 
 	vec pL = L-dot(L,N)*N;
-	vec lCoord = vec2(dot(pL,X),dot(pL,Y));
+	vec lCoord(2);
+	lCoord[0] = dot(pL,X);
+	lCoord[1] = dot(pL,Y);
 	lCoord /= (1.0+dot(L,N));
 
-	if (length(lCoord)>EPS)	
+	if (norm(lCoord)>EPSILON)	
 	{	
-		vec2 lDir = normalize(lCoord);
-		mat2 lRot = mat2(lDir.x, lDir.y, -lDir.y, lDir.x);
-		vCoord *= lRot;
+		vec lDir = normalize(lCoord);
+
+		vec temp(2);
+		temp[0] = lDir[0]*vCoord[0] + lDir[1]*vCoord[1];
+		temp[1] = lDir[0]*vCoord[1] - lDir[1]*vCoord[0];
+		vCoord = temp;
 	}
 
 	return vCoord;
 }
-*/
+
 int main(int argc, char** argv)
 {
     QApplication app(argc, argv, false);
@@ -104,18 +110,28 @@ int main(int argc, char** argv)
                 {
                     in_angle[3] = phi_out * 0.5*M_PI / 180.0;
 
-                    vec in(d_size);
+                    vec in(d_size), cart(6), L(3), V(3);
                     params::convert(in_angle, params::SPHERICAL_TL_PL_TV_PV, data_param, &in[0]);
+                    params::convert(in_angle, params::SPHERICAL_TL_PL_TV_PV, params::CARTESIAN, &cart[0]);
+						  L[0] = cart[0];
+						  L[1] = cart[1];
+						  L[2] = cart[2];
+						  V[0] = cart[3];
+						  V[1] = cart[4];
+						  V[2] = cart[5];
 
                     // Copy the input vector
                     vec x = d->value(in);
 
+						  // Get the projected 2D coordinate
+						  vec xy = coord(V, L, X, Y, N);
+
                     for(int i=0; i<d->dimY(); ++i)
                     {
                         double val = x[i] * cos(in_angle[2]);
 
                         rawm0[i] += val;
-                        rawm1[i] += theta_out * val;
+                        rawm1[i] += val * xy[0];
                     }
                 }
             }