function output ( na, ea, sa, nr, er, sr, nt, et, st, ntot )

%*****************************************************************************80
%
%% output() prints the results of the reactor shielding simulation.
%
%  Licensing:
%
%    This code is distributed under the MIT license.
%
%  Modified:
%
%    02 September 2012
%
%  Author:
%
%    Original FORTRAN77 version by Kahaner, Moler, Nash.
%    MATLAB version by John Burkardt.
%
%  Reference:
%
%    David Kahaner, Cleve Moler, Steven Nash,
%    Numerical Methods and Software,
%    Prentice Hall, 1989,
%    ISBN: 0-13-627258-4,
%    LC: TA345.K34.
%
%  Input:
%
%    integer NA, number of particles absorbed by the slab.
%
%    real EA, energy absorbed by the slab.
%
%    real SA, the sum of the squares of the 
%    absorbed energies.
%
%    integer NR, number of particles reflected by the slab.
%
%    real ER, energy reflected by the slab.
%
%    real SR, the sum of the squares of the 
%    reflected energies.
%
%    integer NT, number of particles transmitted by the slab.
%
%    real ET, energy transmitted through the slab.
%
%    real ST, the sum of the squares of the 
%    transmitted energies.
%
%    integer NTOT, the total number of particles.
%
  fprintf ( 1, '\n' );
  fprintf ( 1, '  The Reactor Shielding Problem:\n' );
  fprintf ( 1, '\n' );
  fprintf ( 1, '                           Total                   Average\n' );
  fprintf ( 1, '                    #      Energy      ' );
  fprintf ( 1, 'Percent     Energy         StDev\n' );
  fprintf ( 1, '\n' );

  etot = ea + er + et;

  if ( 0 < na )
    ea_ave = ea / na;
    sa = sqrt ( sa / na - ea_ave * ea_ave );
  else
    ea_ave = 0.0;
  end

  pa = na * 100 / ntot;

  fprintf ( 1, 'Absorbed     %8d  %14g  %6.2f  %14g  %14g\n', ...
    na, ea, pa, ea_ave, sa );

  if ( 0 < nr )
    er_ave = er / nr;
    sr = sqrt ( sr / nr - er_ave * er_ave );
  else
    er_ave = 0.0;
  end

  pr = nr * 100 / ntot;

  fprintf ( 1, 'Reflected    %8d  %14g  %6.2f  %14g  %14g\n', ...
    nr, er, pr, er_ave, sr );

  if ( 0 < nt )
    et_ave = et / nt;
    st = sqrt ( st / nt - et_ave * et_ave );
  else
    et_ave = 0.0;
  end

  pt = nt * 100 / ntot;

  fprintf ( 1, 'Transmitted  %8d  %14g  %6.2f  %14g  %14g\n', ...
    nt, et, pt, et_ave, st );

  ptot = 100.0;

  fprintf ( 1, '\n' );
  fprintf ( 1, 'Total        %8d  %14g  %6.2f\n', ...
    ntot, etot, ptot );

  return
end
