#version 3.8;
global_settings {
	assumed_gamma 1.0
}

/*
** find rational approximation to given real number
** David Eppstein / UC Irvine / 8 Aug 1993
**
** With corrections from Arno Formella, May 2008
**
** usage: a.out r d
**   r is real number to approx
**   d is the maximum denominator allowed
**
** based on the theory of continued fractions
** if x = a1 + 1/(a2 + 1/(a3 + 1/(a4 + ...)))
** then best approximation is found by truncating this series
** (with some adjustments in the last term).
**
** Note the fraction can be recovered as the first column of the matrix
**  ( a1 1 ) ( a2 1 ) ( a3 1 ) ...
**  ( 1  0 ) ( 1  0 ) ( 1  0 )
** Instead of keeping the sequence of continued fraction terms,
** we just keep the last partial product of these matrices.
*/

#declare Num = array;
#declare Den = array;

#macro MakeRatio (r, d)
	//#local r = r/180;							// for converting angles to fractions of pi radians
	#local m = array [2][2];

	#local (startx, X) = (r, r);
	#local maxden = d;
	
	// initialize matrix
	#local (m[0][0], m[1][1]) = (1, 1);			// tuple-style assignment
	#local (m[0][1], m[1][0]) = (0, 0);
	
	#local ai = 1;								// 
	
	#local index = 0;
	
	// loop finding terms until denom gets too big
	#while (m[1][0] * (ai = X) + m[1][1] <= maxden)
		#local ai = floor (ai);
		
		// swap entries in top row
		#local T = m[0][0] * ai + m[0][1];
		#local m[0][1] = m[0][0];
		#local m[0][0] = T;
		
		// swap entries in bottom row
		#local T = m[1][0] * ai + m[1][1];
		#local m[1][1] = m[1][0];
		#local m[1][0] = T;
		
		// debugging
		#declare Num [index+0] = "|";
		#declare Num [index+1] = concat (str (m[0][0], 3, 4), " ");
		#declare Num [index+2] = str (m[0][1], 3, 4);
		#declare Num [index+3] = "|";
		
		#declare Den [index+0] = "|";
		#declare Den [index+1] = concat (str (m[1][0], 3, 4), " ");
		#declare Den [index+2] = str (m[1][1], 3, 4);
		#declare Den [index+3] = "|";
		
		
		
		#if (X = ai)
			#debug "X = ai \n"
			#break							// AF: division by zero
		#end
		
		#local X = 1/(X - ai);
		#local index = index + 4;
		
		
		
		#if (X > 8388607)
			#break							// AF: representation failure
		#end
		
	#end // end while
	
	// now remaining x is between 0 and 1/ai
	// approx as either 0 or 1/m where m is max that will fit in maxden
	// first try zero
	// printf("%ld/%ld, error = %e\n", m[0][0], m[1][0], startx - ((double) m[0][0] / (double) m[1][0]));
	
	#local Div = "/";
	//#local Div = "pi/";
	
	#local Result1 = concat (str (m[0][0], 0, 0), Div, str(m[1][0], 0, 0) )
	#local Result1_error = startx - (m[0][0] / m[1][0]);
	#debug concat ("First Fraction = ", Result1, "     error =", str (Result1_error, 0, 8), "\n")
	
	// now try other possibility
	#local ai = (maxden - m[1][1]) / m[1][0];
	#local m[0][0] = m[0][0] * ai + m[0][1];
	#local m[1][0] = m[1][0] * ai + m[1][1];
	// printf("%ld/%ld, error = %e\n", m[0][0], m[1][0], startx - ((double) m[0][0] / (double) m[1][0]));
	#local Result2 = concat (str (m[0][0], 0, 0), Div, str(m[1][0], 0, 0) )
	#local Result2_error = startx - (m[0][0] / m[1][0]);
	#debug concat ("Second Fraction = ", Result2, "     error =", str (Result2_error, 0, 8), "\n\n")
	
	#if (abs (Result1_error) < abs (Result2_error))
		Result1
	#else
		Result2
	#end
	
	
#end // end macro



#macro _MakeRatio (decNum, Accuracy)
	
	#local (A, Y, R, P, Q) = (array, array, array, array, array);
	#local prevRatio = 0;
	
	#for (k, 0, 24)
	
		#if (k = 0)
			
			#local A[0] = floor (decNum); // function Int in VB is equivalent to Floor</em> 
			#local R[0] = decNum - A[0];
			
			#if (R[0] !=0)
				#local Y[0] = 1 / R[0]; 
			#end
			
			#local P[0] = A[0]; 
			#local Q[0] = 1 ;
			#local Ratio = P[0] / Q[0]; 
			#local Error = abs (decNum - Ratio);
			
		#else 
			
			#local A[k] = floor (Y[k - 1]);
			#local R[k] = Y[k - 1] - A[k]; 
			#if (R[k]=0)
				#break
			#end
			#local Y[k] = 1 / R[k]; 
			
			#if (k = 1)
				#local P[1] = A[1] * P[0] + 1 ;
				#local Q[1] = A[1] * Q[0]; 
			
			#else
				#local P[k] = A[k] * P[k - 1] + P[k - 2]; 
				#local Q[k] = A[k] * Q[k - 1] + Q[k - 2]; 
			
			#end // end if
			
			#local Ratio = P[k] / Q[k]; 
			#local Error = abs (Ratio - prevRatio);
			
		#end // end if
		
		#if (Error < Accuracy)
			#break
		#else
		
			#local prevRatio = Ratio;
		#end // end if
	
	#end // end for k
	
	// report the results
	#local nominator = P[k - 1];
	#local denominator = Q[k - 1];
	#local remainder = R[k - 1];
	#local Decimal2Fraction = nominator / denominator ;
	
	#local Result = concat (str (nominator, 0, 0), "/", str(denominator, 0, 0) );
	
	Result
	
#end // end macro
  
  
  
  
  


#declare N = 1/32;

//#declare Fraction = MakeRatio (N, 0.000001)




#debug concat ("Best Fraction = ", MakeRatio (N, 0.000001) "\n")


/*
#declare NumSize = dimension_size (Num, 1)-1;
#declare DenSize = dimension_size (Den, 1)-1;

#debug concat (
	#for (i, 0, NumSize)
		Num[i],
	#end
	"\n"
)

#debug concat (
	#for (i, 0, DenSize)
		Den[i],
	#end
	"\n"
)
*/


#error "No objects in scene"

// Original code:
/*
main(ac, av)
int ac;
char ** av;
{
    double atof();		// casts char to float
    int atoi();		// casts char to int
    void exit();

    long m[2][2];
    double x, startx;
    long maxden;
    long ai;

    /* read command line arguments 
    if (ac != 3) {
        fprintf(stderr, "usage: %s r d\n",av[0]);  // AF: argument missing
        exit(1);
    }
    startx = x = atof(av[1]);		// command line argument 1 (r)
    maxden = atoi(av[2]);		// command line argument 2 (d)

    /* initialize matrix 
    m[0][0] = m[1][1] = 1;
    m[0][1] = m[1][0] = 0;

    /* loop finding terms until denom gets too big 
    while (m[1][0] *  ( ai = (long)x ) + m[1][1] <= maxden) {
        long t;
        t = m[0][0] * ai + m[0][1];
        m[0][1] = m[0][0];
        m[0][0] = t;
        t = m[1][0] * ai + m[1][1];
        m[1][1] = m[1][0];
        m[1][0] = t;
        if(x==(double)ai) break;     // AF: division by zero
        x = 1/(x - (double) ai);
        if(x>(double)0x7FFFFFFF) break;  // AF: representation failure
    } 

    /* now remaining x is between 0 and 1/ai 
    /* approx as either 0 or 1/m where m is max that will fit in maxden 
    /* first try zero 
    printf("%ld/%ld, error = %e\n", m[0][0], m[1][0],
           startx - ((double) m[0][0] / (double) m[1][0]));

    /* now try other possibility 
    ai = (maxden - m[1][1]) / m[1][0];
    m[0][0] = m[0][0] * ai + m[0][1];
    m[1][0] = m[1][0] * ai + m[1][1];
    printf("%ld/%ld, error = %e\n", m[0][0], m[1][0],
           startx - ((double) m[0][0] / (double) m[1][0]));
}

*/