% ---------------------------------------------------
% calculate polynome for reverse function
% 

global pan
global pbn

% GDC polynome parameters
% not normalised
a0 = 1.7839314937591553;
a1 = -4.9145695811603218e-005;
a2 = 1.4475514120704247e-009;
a3 = -2.6342226604235755e-014;
a4 = 2.2914257475971922e-019;
a5 = -3.6812716762488225e-025;
a6 = -4.023005095406736e-030;
% pa : is the polynome representation
pa = [a6;a5;a4;a3;a2;a1;a0];

% normalise the GDC circle such that we have the gain = 1
% at the outer region
% e.g. 215 is the outer reagion of interrest to us.
x_max = 215;

y_xmax = polyval(pa, x_max.^2);
% normalised version of pa
pan = pa ./ y_xmax;

printf("normalised version of address@hidden");
printf("x_max: %d\n", x_max);
printf("normalised -> pan: %e\n", van);
printf("-----------------------------------------\n");

% calculate the inverse gdc function g(x)
% define the function F(x)
function y = FF(x)
	global pan;
	y = x .* polyval(pan, x.^2);
endfunction

%define the function G(y)
function x = GG(y)
	global pbn;
	x = y .* polyval(pbn, y.^2);
endfunction

% -----------------------------------------
% calculate the inverse polynome of F(x)

% define the x range of interrest
vx_ = [0:1:215];

Fx = FF(vx_);
%plot(vx_, Fx, "-1;F(x);");
%sleep(5);

% x   = y * poly(y^2)
% x/y = poly(y^2)
rxy = vx_ ./ Fx;
% eleminate the NaN entry
rxy(1) = 1/pan(7); 

pbn = polyfit(Fx, rxy, 7);
vy_ = polyval(pbn, vx_);
plot(Fx, rxy, "-1;f(x);", vx_, vy_, "-3;g(x);" );