function panres = mfs_results_vlmh(panels, pancols, cref, freqresp)

# usage: panres = mfs_results_vlmh(panels, pancols, cref, freqresp)
#
# Input  panels      Structure with panel data
#        pancols     Structure with panel column data
#        cref        Reference chord length
#        freqresp    Structure with frequency response
#
# Output panres      Structure with panel results
#
# Structure panres has the following fields:
#    nofpan              Number of panels
#    P(3, nofpan)        Coordinates of point between A and B
#    f(3, nofpan, nfreq) Panel force vectors
#    p(nofpan, nfreq)    Panel pressure
#
# --------------------------------------------------------------------

# Check arguments

  if (nargin != 4 || nargout != 1)
     print_usage();
  end

# Number of panels and reduced frequencies

  nofpan = panels.nofpan;
  nfreq  = freqresp.nfreq;

# Coordinates of points the aerodynamic force is applied to

  P = 0.5 * (panels.A + panels.B);

# Panel lengths

  dx = (panels.C(1, :) - P(1, :))';

# Pressure

  p = freqresp.G ./ dx;
  s = (4 / cref) * i * diag(freqresp.kred);

  for n = 1 : pancols.npcol
      pids = pancols.panids{n};
      Gs   = [zeros(1, nfreq); 
              cumsum(freqresp.G(pids(1 : end-1), :))];
      p(pids, :) += Gs * s;
  end

# Force

  f    = complex(zeros(3, nofpan, nfreq));
  nvec = (panels.nvec .* panels.area)';

  for k = 1 : nfreq
      f(:, :, k) = (p(:, k) .* nvec)';
  end

# Build structure with results

  panres = struct("nofpan", nofpan, "P", P, "f", f, "p", p);

end
