//*****************************************************************************************
// Truevision - a 3d modeler for gnome and povray
//
// function.cc
//
// Vincent LE PRINCE <vincentleprince@users.sourceforge.net>
// Copyright (C) 2000-2005 Vincent LE PRINCE
// This file is part of the TRUEVISION Package

//   This program is free software; you can redistribute it and/or modify
//   it under the terms of the GNU General Public License as published by
//   the Free Software Foundation; either version 2 of the License, or
//   (at your option) any later version.
//
//   This program is distributed in the hope that it will be useful,
//   but WITHOUT ANY WARRANTY; without even the implied warranty of
//   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
//   GNU General Public License for more details.
//
//   You should have received a copy of the GNU General Public License
//   along with this program; if not, write to the Free Software
//   Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.  */ 
//*******************************************************************************************
#include "include/function.h"
#include <gnome.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include <iostream>
#include "include/fnintern.h"

 /* Some <math.h> files do not define M_PI... */
#ifndef M_PI
#define M_PI 3.14159265358979323846
#endif

const char *TvFuncParserErrors[] = {
	N_("Expression expected, operator found instead"),
	N_("Unknown function : "),
	N_("Unknown keyword : "),
	N_("Numeric expression expected at end, nothing found"),
	N_("Wrong number of arguments for function : "),
};


	
const int TvFunc_Constants_num = 7;
const char *TvFunc_Constants_str[ TvFunc_Constants_num ] = {
	"false", "no", "off", "on", "pi", "true", "yes" };
const double TvFunc_Constants_vals[ TvFunc_Constants_num ] = {
	0, 0, 0, 1, M_PI, 1, 1 };

const int TvFunc_Var_num = 6;
const char *TvFunc_Var_str[ TvFunc_Var_num ] = {
	"x", "y", "z", "t", "u","v" };
const int TvFunc_Var_vals[ TvFunc_Var_num ] = { 0,1,2,2,0,1 };
	
const int TvFunc_Func_num = 104;
const int First_Internal = 28;
const char *TvFunc_Func_str[ TvFunc_Func_num ] = {
	"abs", "acos", "acosh", "asin", "asinh",
	"atan", "atanh", "ceil", "cos", "cosh",
	"degrees",  "exp", "floor", "int", "ln", 
	"log", "radians", 
	"sin", "sinh", "sqrt", "tan", "tanh",
	"atan2", "div", "mod", "pow",
	"max", "min",
	// Pov internal functions
	"f_algbr_cyl1", "f_algbr_cyl2", "f_algbr_cyl3", 	"f_algbr_cyl4", "f_bicorn",
	"f_bifolia", "f_blob", "f_blob2",
	"f_boy_surface", "f_comma", "f_cross_ellipsoids", "f_crossed_through", 
	"f_cubic_saddle", "f_cushion", "f_devils_curve", "f_devils_curve_2d",
	"f_dupin_cyclid", "f_ellipsoid", "f_enneper", "f_flange_cover", "f_folium_surface",
	"f_folium_surface_2d",  "f_glob", "f_heart",  "f_helical_torus", "f_helix1", "f_helix2",
	"f_hex_x", "f_hex_y", "f_hetero_mf", "f_hunt_surface", "f_hyperbolic_torus",
	"f_isect_ellipsoid", "f_kampyle_of_eudoxus", "f_kampyle_of_eudoxus_2d",
	"f_klein_bottle", "f_kummer_surface_v1", "f_kummer_surface_v2", "f_lemniscate_of_gerono",
	"f_lemniscate_of_gerono_2d", "f_mesh1", "f_mitre", " f_nodal_cubic", "f_odd",
	"f_ovals_of_cassini", "f_paraboloid", "f_parabolic_torus", "f_ph", "f_pillow",
	"f_piriform", "f_piriform_2d", "f_poly4", "f_polytubes", "f_quantum", "f_quartic_paraboloid",
	"f_quartic_saddle", "f_quartic_cylinder", "f_r", "f_ridge", "f_ridged_mf", "f_rounded_box",
	"f_sphere", "f_spikes", "f_spikes_2d", "f_spiral", "f_steiners_roman", "f_strophoid",
	"f_strophoid_2d", "f_superellipsoid", "f_th", "f_torus", "f_torus2", "f_torus_gumdrop",
	"f_umbrella", "f_witch_of_agnesi", "f_witch_of_agnesi_2d"
	};
const int TvFunc_Func_args[ TvFunc_Func_num ] = { 
	1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,  1, 1, 1, 1, 1, 1,
	2, 2, 2, 2,
	0, 0,
	// Pov internal functions
	8, 8, 8, 8, 5, 5, 8, 7,  5, 4, 7, 4, 4, 4, 4, 9, 9, 6, 4, 7, 6, 9, 4, 4, 13, 10, 10, 4, 4, 9, 4, 6, 7, 6, 9,
	4, 4, 7, 4, 9, 8, 4, 4, 4, 7, 4, 6, 3, 4, 4, 10, 8, 9, 4, 4, 4, 6, 3, 9, 9, 7, 4, 8, 7, 9, 4,
	7, 10, 5, 3, 5, 6, 4, 4, 5, 9
	 };
	
//*****************************************************************************************
// Function operator object
//*****************************************************************************************
TvFuncOp::TvFuncOp( char car )
{
folded = false;
if ( car == '+' ) { type = TV_FUNC_OP_PLUS; return; }
if ( car == '-' ) { type = TV_FUNC_OP_MINUS; return; }
if ( car == '*' ) { type = TV_FUNC_OP_PROD; return; }
if ( car == '/' ) { type = TV_FUNC_OP_DIV; return; }
}

TvFuncOp::~TvFuncOp()
{
Exprs.clear();
args.clear();
}


double TvFuncOp::evaluate( float x, float y, float z )
{
//cout << "\neval type =" << type;
if ( type == TV_FUNC_EQUAL ) {
	//cout << "\neval = " << Exprs[0]->evaluate( x, y, z ); cout.flush();
	return Exprs[0]->evaluate( x, y, z );
}
double arg1 = args[0]->evaluate( x, y, z );
double arg2 = args[1]->evaluate( x, y, z );
if ( type == TV_FUNC_OP_PLUS ) return arg1+arg2;
if ( type == TV_FUNC_OP_MINUS ) return arg1-arg2;
if ( type == TV_FUNC_OP_PROD ) return arg1*arg2;
if ( type == TV_FUNC_OP_DIV ) 
	{
	if ( arg2 == 0 )  return 0;
	else return arg1/arg2;
	}
return 0;
}


//*****************************************************************************************
// Function expression object
//*****************************************************************************************
TvFuncExpr::TvFuncExpr( char *str, int len )
{
string = str;
string_len = len;
status = true;
error_code = -1;
error_str_len = -1;
error_str_tmp = NULL;
invert = false;
negate = false;
	
//cout << "\ngot expr ->";
//for ( int i = 0 ; i < len ; i ++ ) cout << str[i];
//cout << "<-";
	
// Count operators
bool op_status = false;
vector<int> ops;
int brace_count = 0;
bool brace_status = false;
int par_count = 0;

	
for ( int i = 0 ; i < string_len ; i++ )
	{
	if ( string[i] == '{' ) par_count++;
	if ( string[i] == '}' ) par_count--;		
	if ( par_count != 0 ) continue;
	if ( string[i] == '(' && brace_status ) brace_count++;
	if ( string[i] == ')' ) brace_count--;
	if ( brace_count != 0 ) continue;
	
	if ( string[i] == ' ' || string[i] ==  '\t' || string[i] ==  '\n' ) continue;
	if ( string[i] == '*' )
		if (op_status == false ) { status = false; error_code = 0; break; }
		else { ops.push_back( i ); op_status = false; continue; }
	if ( string[i]== '/' )
		if (op_status == false ) { status = false; error_code = 0; break; }
		else { ops.push_back( i ); op_status = false; continue; }
	if ( string[i]== '+' )
		if (op_status == true ) { ops.push_back( i ); op_status = false;  continue; }
	if ( string[i]== '-' )
		if (op_status == true ) { ops.push_back( i ); op_status = false;  continue; }
	op_status = true;
	brace_status = true;
	}
if ( status == false ) return;

// Parsing modes
if ( ops.size() == 0 )
	{
	// mode 1 : function, constants, variables...
	mode = 1;
	parse_expression();
	}
else
	{
	// mode 2 : arythmetic expressions
	mode = 2;
	int offset = 0;
	for ( unsigned int i = 0; i < ops.size() ; i++ )
		{
		TvFuncOp * eq = new TvFuncOp( TV_FUNC_EQUAL );
		TvFuncExpr *expr = new TvFuncExpr( string+offset, ops[i] - offset );
		if ( expr->get_status() == false ) {  propagate_error( expr ); return;	}
		eq->add_expr( expr );
		AryTree.push_back( eq );
		offset = ops[i]+1;
		TvFuncOp * op = new TvFuncOp( string[ ops[i] ] );
		AryTree.push_back( op );	
		}
	TvFuncOp * eq = new TvFuncOp( TV_FUNC_EQUAL );
	TvFuncExpr *expr = new TvFuncExpr( string+offset, string_len-offset );
	if ( expr->get_status() == false ) {   propagate_error( expr ); return; }
	eq->add_expr( expr );
	AryTree.push_back( eq );		
	compile_ary_tree();
	}
ops.clear();
}

TvFuncExpr::~TvFuncExpr()
{
if ( error_str_tmp != NULL ) delete error_str_tmp;
function_args.clear();
AryTree.clear();
}


void TvFuncExpr::parse_expression()
{
int i = 0;
// skip spaces, braces and tabs, setting invert mode
while ( i < string_len )
	{
	if ( string[i] == ' ' || string[i] == '\t' || string[i] == '\n' || string[i] == '(' ) { i++; continue; }
	if ( string[i] == '-' ) { i++; invert = !invert; continue; }
	if ( string[i] == '!' ) { i++; negate = !negate; continue; }
	if ( string[i] == '+' ) { i++; continue; }
	break;
	}

// Null expression
if ( i == string_len ) 
	{
	status = false;
	error_code = 3;
	return;
	}
	
// Check for functions
bool is_func = false;
int brace_offset = 0;
for ( brace_offset = i ; brace_offset < string_len ; brace_offset++ )
	if ( string[ brace_offset ] == '(' ) { is_func = true; break; }

	
// Functions
if ( is_func )
	{
	eval_mode = 1;

	for ( int f = 0 ; f < TvFunc_Func_num ; f++ )
		{
		int len = strlen( TvFunc_Func_str[f] );
		if ( ! strncmp( string+i, TvFunc_Func_str[f], len ) )
			{
			// for extra chars after keyword
			bool wrong = false;
			for ( int ump = i + len; ump < brace_offset ; ump++ )
				{
				if ( string[ump] != '(' && string[ump] != ' ' ) { wrong = true; break; }
				}
			if ( wrong == true ) continue;
				
			function = f;
			//cout << "\nWe got a func !! -> " << TvFunc_Func_str[function] ;
				
			// Get arguments
			int offset = brace_offset;
			bool terminate = false;
			int brace_count = 0;
			for ( int u = offset+1 ; u < string_len ; u++ )
				{
				if ( string[u] == '(') brace_count++;
				if ( string[u] == ')') 
					{
					if ( brace_count == 0 ) terminate = true;
					else brace_count--;
					}
				if ( ( string[u] == ',' && brace_count == 0 ) || terminate == true )
					{
					TvFuncExpr *expr = new TvFuncExpr( string+offset+1, u-offset-1 );
					if ( expr->get_status() == false ) { propagate_error( expr ); return; }
					function_args.push_back( expr );
					offset = u;
					}
				if ( terminate ) break;
				}
			if ( function_args.size() != (unsigned int )TvFunc_Func_args[function] && TvFunc_Func_args[function] != 0 )
				{
				status = false;
				error_code  = 4;
				error_str = string+i;
				error_str_len = len;
				}
			return;	
			}
		}		
	// unknown function !!
	status = false;
	error_code = 1;
	error_str = string+i; 
	error_str_len = brace_offset - i;
	return;
	}

// Pigment block
if ( ! strncmp(	string+i, "pigment", 7 )  )
	{
	eval_mode = 4;
	//cout << "\nWe got a pig !! -> "; cout.flush();
	return;
	}
	
// Numeric constant
if ( string[i] >= '0' && string[i] <= '9' )
	{
	eval_mode = 2;
	float temp_const;
	sscanf( string+i, "%f ", &temp_const );
	constant = temp_const;
	if ( invert ) constant = - constant;
	if ( negate ) constant = ! constant;
	//cout << "\nWe got a num const !! -> " << constant;
	return;
	}
	
// Symbolic constant
for ( int u = 0 ; u < TvFunc_Constants_num ; u++ )
	{
	int len = strlen( TvFunc_Constants_str[u] );
	if ( len + i < string_len-1 & string[i+len] != ' ' ) continue;
	if ( ! strncmp( string+i, TvFunc_Constants_str[u], len ) )
		{
		eval_mode = 2;
		constant = TvFunc_Constants_vals[u];
		if ( invert ) constant = - constant;
		if ( negate ) constant = ! constant;
		//cout << "\nWe got a symb const !! -> " << constant;
		return;	
		}
	}
	
// Variable
for ( int u = 0 ; u < TvFunc_Var_num ; u++ )
	{
	int len = strlen( TvFunc_Var_str[u] );
	if ( len + i < string_len-1 & string[i+len] != ' ' ) continue;
	if ( ! strncmp( string+i, TvFunc_Var_str[u], len ) )
		{
		eval_mode = 3;
		variable = TvFunc_Var_vals[u];
		//cout << "\nWe got a var !! -> " << variable; cout.flush();
		return;	
		}
	}


	
// Unknown
status = false;
error_code = 2;
error_str = string+i;
error_str_len = string_len-i;
}


void TvFuncExpr::compile_ary_tree()
{
//cout << "\ncompiling tree"; cout.flush();
//cout << "\n\ttree size = " << AryTree.size(); cout.flush();	
bool passed = false;
while ( !passed  )
	{
	passed = true;
	for ( unsigned int i = 0 ; i < AryTree.size(); i++ )
		{
		if ( ( AryTree[i]->get_type() == TV_FUNC_OP_PROD || AryTree[i]->get_type() == TV_FUNC_OP_DIV )  && ! AryTree[i]->get_folded() ) 
			{ 
			//cout << "\none iteration";
			AryTree[i]->add_arg( AryTree[i-1] );
			AryTree[i]->add_arg( AryTree[i+1] );
			AryTree[i]->set_folded();
				//cout << "\ntree size = " << AryTree.size(); cout.flush();	
			AryTree.erase( AryTree.begin() + i +1 );
			AryTree.erase( AryTree.begin() + i -1 );
			passed = false;
			break;
			}
		}
	}
//cout << "\nfolded"; cout.flush();
//cout << "\ntree size = " << AryTree.size(); cout.flush();	
while ( AryTree.size() != 1 )
	{
	//cout << "\ntree size = " << AryTree.size(); cout.flush();
	for ( unsigned int i = 0 ; i < AryTree.size(); i++ )
		{
		if ( ( AryTree[i]->get_type() == TV_FUNC_OP_PLUS ||  AryTree[i]->get_type() == TV_FUNC_OP_MINUS ) && ! AryTree[i]->get_folded() ) 
			{ 
			//cout << "\none iteration at i=" << i ;
			AryTree[i]->add_arg( AryTree[i-1] );
			AryTree[i]->add_arg( AryTree[i+1] );
			AryTree[i]->set_folded();
			AryTree.erase( AryTree.begin() + i +1 );
			AryTree.erase( AryTree.begin() + i -1 );
			break;
			}
		}
	}	
//cout << "\ncompiled"; cout.flush();
//cout << "\ntree size = " << AryTree.size(); cout.flush();
}


double TvFuncExpr::evaluate( float x, float y, float z )
{
if ( mode == 1 )
	{
	// Eval single expression
//cout << "\neval exp"; cout.flush();
	switch ( eval_mode )
		{
		// Function
		case 1:
			{
			//cout << "\neval function"; cout.flush();
			int anum = function_args.size();
			double *fargs = new double[ anum ];
			for ( int a = 0 ; a < anum ; a++ )
				{
				fargs[a] = function_args[a]->evaluate( x, y, z );
				}
			double res = 0;

			if ( function == 0 ) res = fabs( fargs[0] );			
			if ( function == 1 ) res = acos( fargs[0] );			
			if ( function == 2 ) res = acosh( fargs[0] );			
			if ( function == 3 ) res = asin( fargs[0] );			
			if ( function == 4 ) res = asinh( fargs[0] );			
			if ( function == 5 ) res = atan( fargs[0] );			
			if ( function == 6 ) res = atanh( fargs[0] );			
			if ( function == 7 ) res = ceil( fargs[0] );			
			if ( function == 8 ) res = cos( fargs[0] );
			if ( function == 9 ) res = cosh( fargs[0] );
			if ( function == 10 ) res = fargs[0] *180.0 / M_PI;
			if ( function == 11 ) res = exp( fargs[0] );
			if ( function == 12 ) res = floor( fargs[0] );
			if ( function == 13 ) res = (int)( fargs[0] );
			if ( function == 14 ) { if ( fargs[0] != 0 ) res = log( fargs[0] ); else res = 0; }
			if ( function == 15 ) { if ( fargs[0] != 0 ) res = log10( fargs[0] ); else res = 0;  }
			if ( function == 16 ) res = fargs[0] * M_PI / 180.0;
			if ( function == 17 ) res = sin( fargs[0] );			
			if ( function == 18 ) res = sinh( fargs[0] );			
			if ( function == 19 ) res = sqrt( fargs[0] );			
			if ( function == 20 ) res = tan( fargs[0] );			
			if ( function == 21 ) res = tanh( fargs[0] );			
			if ( function == 22 ) res = atan2( fargs[0], fargs[1] );			
			if ( function == 23 ) res = (int)fargs[0] / (int)fargs[1];		
			if ( function == 24 ) res = fmod( fargs[0], fargs[1] );			
			if ( function == 25 ) res = pow( fargs[0], fargs[1] );			
			if ( function == 26 ) { res = fargs[0]; for ( int it = 0 ; it < anum ; it++ ) res = res < fargs[it] ? fargs[it] : res; }		
			if ( function == 27 ) { res = fargs[0]; for ( int it = 0 ; it < anum ; it++ ) res = res > fargs[it] ? fargs[it] : res; }		

			if ( function == 28 ) { res = f_algbr_cyl1( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6], fargs[7]  );}
			if ( function == 29 ) { res = f_algbr_cyl2( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6], fargs[7]  );}
			if ( function == 30 ) { res = f_algbr_cyl3( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6], fargs[7]  );}
			if ( function == 31 ) { res = f_algbr_cyl4( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6], fargs[7]  );}
			if ( function == 32 ) { res = f_bicorn( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4]  );}
			if ( function == 33 ) { res = f_bifolia( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4]  );}
			if ( function == 34 ) { res = f_blob( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6], fargs[7]  );}
			if ( function == 35 ) { res = f_blob2( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6]  );}
			if ( function == 36 ) { res = f_boy_surface( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4] ); }
			if ( function == 37 ) { res = f_comma( fargs[0], fargs[1], fargs[2], fargs[3] ); }
			if ( function == 38 ) { res = f_cross_ellipsoids( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6]  );}
			if ( function == 39 ) { res = f_crossed_through( fargs[0], fargs[1], fargs[2], fargs[3] ); }
			if ( function == 40 ) { res = f_cubic_saddle( fargs[0], fargs[1], fargs[2], fargs[3] ); }
			if ( function == 41 ) { res = f_cushion( fargs[0], fargs[1], fargs[2], fargs[3] ); }
			if ( function == 42 ) { res = f_devils_curve( fargs[0], fargs[1], fargs[2], fargs[3] ); }
			if ( function == 43 ) { res = f_devils_curve_2d( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6], fargs[7], fargs[8]  );}
			if ( function == 44 ) { res = f_dupin_cyclid( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6], fargs[7], fargs[8]  );}
			if ( function == 45 ) { res = f_ellipsoid( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5]  );}
			if ( function == 46 ) { res = f_enneper( fargs[0], fargs[1], fargs[2], fargs[3]  );}
			if ( function == 47 ) { res = f_flange_cover( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6]  );}
			if ( function == 48 ) { res = f_folium_surface( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5]  );}
			if ( function == 49 ) { res = f_folium_surface_2d( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6], fargs[7], fargs[8]  );}
			if ( function == 50 ) { res = f_glob( fargs[0], fargs[1], fargs[2], fargs[3]  );}
			if ( function == 51 ) { res = f_heart( fargs[0], fargs[1], fargs[2], fargs[3]  );}
			if ( function == 52 ) { res = f_helical_torus( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6], fargs[7], fargs[8] , fargs[9], fargs[10], fargs[11], fargs[12]  );}
			if ( function == 53 ) { res = f_helix1( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6], fargs[7], fargs[8] , fargs[9]  );}
			if ( function == 54 ) { res = f_helix2( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6], fargs[7], fargs[8] , fargs[9]  );}		
			if ( function == 55 ) { res = f_hex_x( fargs[0], fargs[1], fargs[2], fargs[3]  );} 				
			if ( function == 56 ) { res = f_hex_y( fargs[0], fargs[1], fargs[2], fargs[3]  );} 				
			if ( function == 57 ) { res = f_hetero_mf( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6], fargs[7], fargs[8]  );} 				
			if ( function == 58 ) { res = f_hunt_surface( fargs[0], fargs[1], fargs[2], fargs[3]  );} 					
			if ( function == 59 ) { res = f_hyperbolic_torus( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5]  );} 				
			if ( function == 60 ) { res = f_isect_ellipsoids( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6]  );}  				
			if ( function == 61 ) { res = f_kampyle_of_eudoxus( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5]  );}  				
			if ( function == 62 ) { res = f_kampyle_of_eudoxus_2d( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6], fargs[7], fargs[8] );} 				
			if ( function == 63 ) { res = f_klein_bottle(  fargs[0], fargs[1], fargs[2], fargs[3]  );}  				
			if ( function == 64 ) { res = f_kummer_surface_v1( fargs[0], fargs[1], fargs[2], fargs[3]  );}				
			if ( function == 65 ) { res = f_kummer_surface_v2(	fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6] ); }			
			if ( function == 66 ) { res = f_lemniscate_of_gerono( fargs[0], fargs[1], fargs[2], fargs[3] );} 	
			if ( function == 67 ) { res = f_lemniscate_of_gerono_2d( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6], fargs[7],  fargs[8] );} 	
			if ( function == 68 ) { res = f_mesh1( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6], fargs[7]  );} 				
			if ( function == 69 ) { res = f_mitre( fargs[0], fargs[1], fargs[2], fargs[3]  );} 				
			if ( function == 70 ) { res = f_nodal_cubic(  fargs[0], fargs[1], fargs[2], fargs[3]  );} 								
			if ( function == 71 ) { res = f_odd( fargs[0], fargs[1], fargs[2], fargs[3]  );}					
			if ( function == 72 ) { res = f_ovals_of_cassini( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6]  );}  	
			if ( function == 73 ) { res = f_paraboloid( fargs[0], fargs[1], fargs[2], fargs[3] );}  	
			if ( function == 74 ) { res = f_parabolic_torus( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5]  );}
			if ( function == 75 ) { res = f_ph( fargs[0], fargs[1], fargs[2] );}
			if ( function == 76 ) { res = f_pillow( fargs[0], fargs[1], fargs[2], fargs[3] );}
			if ( function == 77 ) { res = f_piriform( fargs[0], fargs[1], fargs[2], fargs[3] );} 
			if ( function == 78 ) { res = f_piriform_2d( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6], fargs[7], fargs[8] , fargs[9] );} 
			if ( function == 79 ) { res = f_poly4(  fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6], fargs[7]  );}
			if ( function == 80 ) { res = f_polytubes(  fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6], fargs[7], fargs[8] );}
			if ( function == 81 ) { res = f_quantum(  fargs[0], fargs[1], fargs[2], fargs[3] );}
			if ( function == 82 ) { res = f_quartic_paraboloid( fargs[0], fargs[1], fargs[2], fargs[3] );}
			if ( function == 83 ) { res = f_quartic_saddle(  fargs[0], fargs[1], fargs[2], fargs[3] );}
			if ( function == 84 ) { res = f_quartic_cylinder( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5]  );}
			if ( function == 85 ) { res = f_r( fargs[0], fargs[1], fargs[2] );}
			if ( function == 86 ) { res = f_ridge( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6], fargs[7], fargs[8] );}
			if ( function == 87 ) { res = f_ridged_mf( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6], fargs[7], fargs[8] );}
			if ( function == 88 ) { res = f_rounded_box( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6]  );} 
			if ( function == 89 ) { res = f_sphere( fargs[0], fargs[1], fargs[2], fargs[3] );}
			if ( function == 90 ) { res = f_spikes( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6], fargs[7]  );}
			if ( function == 91 ) { res = f_spikes_2d( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6]  );} 
			if ( function == 92 ) { res = f_spiral( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6], fargs[7], fargs[8] );}
			if ( function == 93 ) { res = f_steiners_roman( fargs[0], fargs[1], fargs[2], fargs[3] );}
			if ( function == 94 ) { res = f_strophoid( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6]  );} 
			if ( function == 95 ) { res = f_strophoid_2d( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6], fargs[7], fargs[8] , fargs[9] );}
			if ( function == 96 ) { res = f_superellipsoid( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4]  );} 
			if ( function == 97 ) { res = f_th( fargs[0], fargs[1], fargs[2] );}
			if ( function == 98 ) { res = f_torus( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4]  );} 
			if ( function == 99 ) { res = f_torus2( fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5]  );}
			if ( function == 100 ) { res = f_torus_gumdrop( fargs[0], fargs[1], fargs[2], fargs[3] );}
			if ( function == 101 ) { res = f_umbrella(  fargs[0], fargs[1], fargs[2], fargs[3] );}
			if ( function == 102 ) { res = f_witch_of_agnesi(  fargs[0], fargs[1], fargs[2], fargs[3], fargs[4]  );}
			if ( function == 103 ) { res = f_witch_of_agnesi_2d(  fargs[0], fargs[1], fargs[2], fargs[3], fargs[4], fargs[5], fargs[6], fargs[7], fargs[8] );}
				
			delete fargs;
			if ( invert ) res = -res;
			if ( negate ) res = ! res;
			return res;
			}
			break;
		
		// Constant
		case 2:
			//cout << "\neval const"; cout.flush();
			return constant;
			break;
		
		// Variable
		case 3:
			{
			//cout << "\neval var"; cout.flush();
			double vars[3] = { x, y, z };
			double res  = vars[ TvFunc_Var_vals[variable] ];
			if ( invert ) res = -res;
			if ( negate ) res = ! res;
			return res;			
			}
			break;
			
		// Pigment block
		default:
		case 4:
			return 0;
			break;
		}
	}
else 
	{
	// eval arythmetic tree
	//cout << "\nEvaluate tree"; cout.flush();
	//cout << "\nsize tree = " << AryTree.size(); cout.flush();
		
	return AryTree[0]->evaluate( x, y, z );
	//cout << "\nEvaluate tree"; cout.flush();
	}
return 0;
}


char *TvFuncExpr::get_error()
{
if ( error_code == -1 ) return NULL;
if ( error_str_len == -1 ) return (char*)(TvFuncParserErrors[error_code]);
int len  = strlen( TvFuncParserErrors[error_code] ) + error_str_len + 1;
error_str_tmp = new char[ len ];
strcpy( error_str_tmp, TvFuncParserErrors[error_code] );
strncat( error_str_tmp, error_str, error_str_len );
return error_str_tmp;
}


//*****************************************************************************************
// Function object
//*****************************************************************************************
TvFunction::TvFunction()
{
string = NULL;
expression = NULL;
}

TvFunction::~TvFunction()
{
if ( string != NULL ) delete string;	
if ( expression != NULL ) delete expression;
}

void TvFunction::set_expression( char *str )
{
if ( string != NULL ) delete string;
int len = strlen( str ) + 1;
string = new char[ len ];
strcpy( string, str );
if ( expression != NULL ) delete expression;
expression = new TvFuncExpr( string, len-1 );
}

double TvFunction::evaluate( float x, float y, float z )
{
return expression->evaluate( x, y, z );
}

void TvFunction::define_internals( ofstream & file )
{
for ( int i = First_Internal ; i < TvFunc_Func_num ; i++ )
	if ( TvFunc_Func_str[i] != NULL && string != NULL )
		if ( strstr( string, TvFunc_Func_str[i] ) ) 
			file  << "\n#declare " << TvFunc_Func_str[i] << " = function { internal(" << (i-First_Internal) << ") }";
}