-->// r = Beta + e(r*cos(theta))
-->// given data
-->theta = [0.88; 1.10; 1.42; 1.77; 2.14];
-->rdata = [3.00; 2.30; 1.65; 1.25; 1.01];
-->clf; // clear figure
-->plot(rdata .* cos(theta), rdata .* sin(theta),'rx')
-->// The model is r vs r*cos(theta)
-->rcostheta = rdata .* cos(theta);
-->// X is [ones rcostheta], y is rdata
-->X = [ones(theta) rcostheta], y = rdata;
 X  =
    1.    1.9114534  
    1.    1.0432711  
    1.    0.2478720  
    1.  - 0.2473610  
    1.  - 0.5443511  
-->// Normal equation X'*X * [Beta, e] = X'*y
-->A = [X'*X X'*y]
 A  =
    5.           2.4108845    9.21       
    2.4108845    5.1610149    7.6838768  
-->B=rref(A)
 B  =
    1.    0.    1.4509292  
    0.    1.    0.8110525  
-->Beta = B(1,3); e=B(2,3)
 e  =
    0.8110525  
-->// Since e < 1 it is an ellipse
-->// Solving the equation for r we get r = Beta/(1-e cos(theta))
-->rAtTheta4point6 = Beta ./ (1 - e*cos(4.6))
 rAtTheta4point6  =
    1.3299545  
-->// plot the curve
-->t2 = 0:%pi/40:2*%pi;
-->r2= Beta ./ (1 - e * cos(t2));
-->polarplot(t2,r2)
-->xs2eps(1,'6711.eps')
-->diary(0)