function [ids, axes] = mfs_axes_b2(elem)

# usage: [ids, axes] = mfs_axes_b2(elem)
#
# Input  elem(:)     Stucture array with element data
#
# Output ids(:)      Element identifiers
#        axes(:, 9)  Element axis vectors
#        rc          Return code: 0 = no errors
#                                 1 = errors
#
# This function computes the vectors of the beam axes.
#
# ------------------------------------------------------------------------

# Initialize

  nel  = length(elem);
  ids  = zeros(nel, 1);
  axes = zeros(nel, 9);

# Loop over elements

  for n = 1 : nel

      ids(n) = elem(n).id;

      d = elem(n).coor(2, :) - elem(n).coor(1, :);
      L = norm(d);
      ex = d' / L;

      if (isfield(elem(n).geom, "Q"))
         v = elem(n).geom.Q - elem(n).coor(1, :);
      else
         v = elem(n).geom.v;
      end

      ez = v' - (v * ex) * ex;
      Lz = norm(ez);
      ez = ez / Lz;

      ey = cross(ez, ex);

      axes(n, :) = [ex', ey', ez'];

  end

end
