myrevkin1(basepoint, endpoint, lengths)

function [j_vars] = myrevkin1(basepoint, endpoint, lengths)
% Returns 3*pi if there is no fesible soution 
% else returns values between 0 and 2*pi
l1=lengths(1); l2=lengths(2);
% bx=2; by=3; ex=6; ey=6;
bx=basepoint(1);
by=basepoint(2);
ex=endpoint(1);
ey=endpoint(2);
syms th1 th2;
f(1,1)=bx+l1*cos(th1)+l2*cos(th1+th2)-ex;
f(2,1)=by+l1*sin(th1)+l2*sin(th1+th2)-ey;
var=[th1;th2];
numvar=length(var);
numeq=length(f);
for j=1:numeq
for i=1:numvar
    J(j,i)=diff(f(j),var(i));
end
end
n=0; nmax=30;
var1=[0.1;0.1];
dmain=[0.1;0.1];
tol=0.0001;
while(((abs(dmain(1))>tol)||(abs(dmain(2)>tol)))&&(n<nmax))
    f1=subs(f,var,var1);
    J1=subs(J,var,var1);
    if(abs(det(J1))>tol*tol)
    dmain=inv(J1)*f1;
    end
    var1=var1-dmain;
    n=n+1;
end
if(n~=nmax)
j_vars(1,:)=var1';
% Now find second position based on geometrical data and first position
slope=atan2(ey-by,ex-bx);
j_vars(2,:)=[2*slope-var1(1) -var1(2)];
for j=1:numvar
for i=1:2
if(j_vars(i,j)<0)
    j_vars(i,j)=j_vars(i,j)-round(j_vars(i,j)/2/3.1416)*2*3.1416+2*3.1416;
else
    j_vars(i,j)=j_vars(i,j)-round(j_vars(i,j)/2/3.1416)*2*3.1416;
end
end
end
% subs(f,var,j_vars(:,1))
% subs(f,var,j_vars(:,2))
else
    disp('no feasible solution');
    j_vars=[10 10; 10 10];
end
return