function haxes = ptdambar(Nodes,Bars,DImat,haxes,inodal,ngrad,DIlocs)
%
%  Name: ptdambar
% 
%  Usage: 
%
% 
%  Description:
%     Plot damage indicator DI for bar geometry
%     For nodal DI - use color gradient with ngrad increments
%     For elemental DI - use solid color for each element subdivision
%
%  Inputs:
%
%	Nodes   = List of nodes
%  Bars    = List of bar connectivity
%  DImat   = Matrix of DI values -- corresponding to rows of either Nodes or Bars
%            (If DImat has more than one column, these values apply to subdivisions
%             of each bar element)
%  haxes   = Handle of axes for plot
%  inodal  = 0, DI are elemental values
%          = 1, DI are nodal values
%  ngrad   = Number of gradient steps for each bar element (or each subdivision, if
%            applicable)
%  naxaz, naxel = 3-D viewing parameters

%  Version SWD970712
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% This matlab source code was originally     %
% developed as part of "DIAMOND" at          %
% Los Alamos National Laboratory. It may     %
% be copied, modified, and distributed in    %
% any form, provided:                        %
%  a) This notice accompanies the files and  %
%     appears near the top of all source     %
%     code files.                            %
%  b) No payment or commercial services are  %
%     received in exchange for the code.     %
%                                            %
% Original copyright is reserved by the      %
% Regents of the University of California,   %
% in addition to Scott W. Doebling, Phillip  %
% J. Cornwell, Erik G. Straser, and Charles  %
% R. Farrar.                                 %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

if ~exist('haxes')
	haxes = axes;
elseif ~ishandle(haxes),
	haxes = axes;
end

axes(haxes)

cla

%  Draw the undeformed geometry

updtgeom(haxes,Nodes)

%
%  Plot the damage indicator
%

%  Define the colors from the colormap

map = colormap;

%  Compute full range of DI

dimin = min(DImat(:));
dimax = max(DImat(:));

distep = (dimax - dimin) / (size(map,1)-1);

%  Map the colors to the DI values

%% ===================================
%%
%%    NODAL DAMAGE INDICATOR VALUES
%%
%% ===================================

if inodal == 1,
   
%  Loop over each bar element

for ibar = 1:length(Bars),
   
%  Select the end Nodes associated with that element

   n1 = Bars(ibar,2);
   n2 = Bars(ibar,3);

%  Find locations of end-nodes in Node list

   in1 = find(Nodes(:,1) == n1);
   in2 = find(Nodes(:,1) == n2);
   
%  Determine the coordinates and DI values at the nodes

   loc1 = Nodes(in1,2:4);
   loc2 = Nodes(in2,2:4);
   
   DI1  = DImat(in1);
   DI2  = DImat(in2);
   
%  Compute the geometric and color increments
   
   geoinc = (loc2 - loc1) / ngrad;
   
   diinc = (DI2 - DI1) / ngrad;

%  Loop over the number of gradient increments
   
   for igrad = 1:ngrad,
      
%    Compute the average DI value for the segment
      
      di1 = DI1 + (igrad-1) * diinc;
      di2 = DI1 + igrad * diinc;
      
      diavg = (di1 + di2) / 2;
      
%    Determine the color of the segment
      
      icolor = floor((diavg - dimin) / distep) + 1;

      segcolor = map(icolor,:);
      
%    Draw the line segment with the appropriate color

      p1 = loc1 + (igrad-1) * geoinc;
      p2 = loc1 + igrad * geoinc;
   
      h = line([p1(1) p2(1)],[p1(2) p2(2)],[p1(3) p2(3)]);
      
      set(h,'linewidth',2,'color',segcolor)

   end

end

%% ===================================
%%
%%  ELEMENTAL DAMAGE INDICATOR VALUES
%%
%% ===================================

elseif inodal == 0,
   
%  Loop over each beam element

for ibar = 1:size(DImat,1),
   
   %  Loop over each segment within the beam
   
   for iseg = 1:size(DImat,2),
      
      %  Determine the coordinates at the ends of this segment

      loc1 = [DIlocs{ibar}.xloc(iseg) DIlocs{ibar}.yloc(iseg) DIlocs{ibar}.zloc(iseg)];
      loc2 = [DIlocs{ibar}.xloc(iseg+1) DIlocs{ibar}.yloc(iseg+1) DIlocs{ibar}.zloc(iseg+1)];
   
      %  Determine the geometric increment for each subdivision
   
      geoinc = (loc2 - loc1) / size(DImat,2);
   
      %    Determine the color of the segment
      
      icolor = floor((DImat(ibar,iseg) - dimin) / distep) + 1;
         
      segcolor = map(icolor,:);
      
     %    Draw the line segment with the appropriate color
   
      h = line([loc1(1) loc2(1)],[loc1(2) loc2(2)],[loc1(3) loc2(3)]);
      
      set(h,'linewidth',3,'color',segcolor)

   end

end

end

set(haxes,'dataaspectratio',[1 1 1])

%axis equal
rotate3d on

%  Define axis labels

xlabel('X'),ylabel('Y'),zlabel('Z')

return