%
%  macm 202 -- 09 mar 03 -- djm
%
%  w10setup.m:  set-up truss & plot 
%     - run w10mytruss#.m first (define:  nodes,join,loads)

Nnodes = size(nodes,1);

%  plot truss

figure(1);  clf
plot(nodes(:,1),nodes(:,2),'*')
hold on
plot(nodes(1:2,1),nodes(1:2,2),'o')
axis([min(nodes(:,1))-1.5 max(nodes(:,1))+1 min(nodes(:,2))-1 max(nodes(:,2))+1])
axis equal
shift = 1/8;
for j = 1:Nnodes
    text(nodes(j,1)+shift,nodes(j,2)+shift,num2str(j),'color','b')
end

%  process strut geometry
%
%  join shows both connections
%  angles:  the alpha^j and alphabar^j
%  length:  of course
%  (Xstrut,Ystrut) = vector coordinates of each strut

join2  = join + transpose(join);
    
length = zeros(size(join));  
angles = length;  struts = length;  
Xstrut = length;  Ystrut = length;
Nstrut = 0;

for j = 1:Nnodes
    for k = 1:Nnodes
        if (join2(j,k)==1)
            Xstrut(j,k) = nodes(k,1)-nodes(j,1);
            Ystrut(j,k) = nodes(k,2)-nodes(j,2);

            if (j<=k)
                Nstrut = Nstrut + 1;
                struts(j,k) = Nstrut;
                plot([nodes(j,1) nodes(k,1)],[nodes(j,2) nodes(k,2)],'k')
            else
                struts(j,k) = struts(k,j);
            end
        end
    end
end

[angles,length] = cart2pol(Xstrut,Ystrut);

%  matlab doesn't plot arrowed vectors automatically

ang = 3*pi/4;
for j = 1:Nnodes
    for k = j:Nnodes
        if (join(j,k)==1)
            xmid = (nodes(k,1)+nodes(j,1))/2;
            xarr = xmid ...
                 + [shift*cos(angles(j,k)+ang) 0 shift*cos(angles(j,k)-ang)];
            ymid = (nodes(k,2)+nodes(j,2))/2;
            yarr = ymid ...
                 + [shift*sin(angles(j,k)+ang) 0 shift*sin(angles(j,k)-ang)];
            plot(xarr,yarr,'r')
            text(xmid+shift,ymid+shift,num2str(struts(j,k)),'color','r')
        end
    end
end

title(['\bf truss plot (' num2str(Nnodes) ...
        ' nodes & ' num2str(Nstrut) ' struts)'])

%  plot loads

[loadA,loadR] = cart2pol(loads(:,1),loads(:,2));

for j = 1:Nnodes
    if (loadR(j)~=0)
        plot(nodes(j,1)+[0 3*shift*cos(loadA(j))], ...
             nodes(j,2)+[0 3*shift*sin(loadA(j))],'g');
    end
end