%%%Ejercicio examen 4 junio 2018
x0=1
xf=2
n=3
%Aplicación del espaciado de Chebyshev
X=[x0 0.5*(x0+xf)-0.5*(xf-x0)*cos(pi*(2*[1 2 3]-1)./(2*n)) xf]
%Obtendión de los valores de Y
Y=1./(X.^2)

%Angulos iniciales e intervalos de giro de las barras de entrada y salida
teta20=10
deltateta2=60
teta40=260
deltateta4=90

%transformación al dominio del ángulo
%(deltateta2/deltaX)=(teta2-teta20)/(X-Xinicial)
deltaX=X(end)-X(1)
teta2=((X-X(1))*deltateta2/deltaX)+teta20

deltaY=Y(end)-Y(1)
teta4=((Y-Y(1))*deltateta4/deltaY)+teta40


%Obtención de las Constantes de Freudenstein K1, K2 y K3

A=[cos(teta4(2:4)*pi/180)' -cos(teta2(2:4)*pi/180)'  [1 1 1]']
b=cos((teta4(2:4)-teta2(2:4))*pi/180)'
 
K=A\b;
r1=3;
r2=r1/K(1)
r4=r1/K(2)
r3=(r1^2+r2^2+r4^2-(K(3)*2*r2*r4))^0.5
[r1 r2 r3 r4]

O=[0;0]
OA=r2.*[cos(teta2*pi/180); sin(teta2*pi/180)]
B=[r1;0]
C=r4*[cos(teta4*pi/180); sin(teta4*pi/180)]

figure
plot(O(1),O(2),'o')
hold on
plot(B(1),B(2),'o')
k=1
data=[O(1) OA(1,1+k) B(1)+C(1,1+k) B(1); O(2) OA(2,1+k) B(2)+C(2,1+k) B(2)];
plot(data(1,:),data(2,:))
plot(OA(1,1+k),OA(2,1+k),'r*')
plot(data(1,3),data(2,3),'r*')

k=2
data=[O(1) OA(1,1+k) B(1)+C(1,1+k) B(1); O(2) OA(2,1+k) B(2)+C(2,1+k) B(2)];
plot(data(1,:),data(2,:),'c')
plot(OA(1,1+k),OA(2,1+k),'r*')
plot(data(1,3),data(2,3),'r*')


k=3
data=[O(1) OA(1,1+k) B(1)+C(1,1+k) B(1); O(2) OA(2,1+k) B(2)+C(2,1+k) B(2)];
plot(data(1,:),data(2,:),'g')
plot(OA(1,1+k),OA(2,1+k),'r*')
plot(data(1,3),data(2,3),'r*')