//catenary curve
rmassr=1; //mass force per unit length
numxpts=100;
x=linspace(0,1,numxpts); //half of curve since it is symmetric
z_rise=linspace(0.1,max(x),numxpts); //0 is the lowest point
a_par=acosh(1+z_rise); //parametrize the catenary curve
arclength=[];
m=[];
theta=[];
tension=[];

for(i=1:length(z_rise))
    cc=cosh(a_par(i)*x)-1; //curve from 0 on left to z_rise(i) on right
    dx=diff(x);
    dy=diff(cc);
    hp=sqrt(dx.^2+dy.^2);
    arclength(i)=sum(hp);
    m(i)=rmassr*arclength(i);
    theta(i)=atan(dy($),dx($));
    tension(i)=m(i)/sin(theta(i));
end

[a b]=min(tension);

subplot(3,1,1)
plot(z_rise,tension');
e=gce();p=e.children(1);//get the handle on the polyline
datatipCreate(p,b);
xlabel("Droop <distance>");
ylabel("Tension <force>");

subplot(3,1,2)
plot(z_rise,arclength','r');
e=gce();p=e.children(1);//get the handle on the polyline
datatipCreate(p,b);
xlabel("Droop <distance>");
ylabel("Rope length <distance>");

subplot(3,1,3)
plot(z_rise,-theta'/%pi*180,'g');
e=gce();p=e.children(1);//get the handle on the polyline
datatipCreate(p,b);
xlabel("Droop <distance>");
ylabel("Anchor angle [deg]");
