function y = universal_egg ( B, L, w, D, x )

%*****************************************************************************80
%
%% universal_egg() is a universal formula for the shape of bird eggs.
%
%  Licensing:
%
%    This code is distributed under the MIT license.
%
%  Modified:
%
%    08 June 2011
%
%  Author:
%
%    John Burkardt
%
%  Reference:
%
%    Valeriy Narushin, Michael Romanov, Darren Griffin,
%    "Egg and math: introducing a universal formula for egg shape",
%    Annals of the New York Academy of Sciences,
%    Number 1505, pages 169-177, 23 August 2021.
%
%  Input:
%
%    real B, the maximum breadth of the egg (y direction);
%
%    real L, the maximum length of the egg (x direction).
%
%    real w: the x-distance between the location of the maximum
%    breadth and the midpoint of the egg's length (x=0).
%
%    real D: the diameter of the egg at a quarter of the length
%    from the pointed end (perhaps x=0.25?).
%
%    real x(*): sample points between -L/2 and +L/2.
%
%  Output:
%
%    real y(*): the height of the egg at each x value.
%    The egg is assumed to be vertically symmetric, and 
%    rotationally symmetric about the x axis.
%
  y_chicken = chicken_egg ( B, L, w, x );

  s1 = sqrt ( 5.5 * L^2 + 11.0 * L * w + 4.0 * w^2 );
  s2 = sqrt ( 3.0 ) * B * L - 2.0 * D * sqrt ( L^2 + 2.0 * w * L + 4.0 * w^2 );
  s3 = sqrt ( 5.5 * L^2 + 11.0 * L * w + 4.0 * w^2 );
  s4 = 2.0 * sqrt ( L^2 + 2.0 * w * L + 4.0 * w^2 );
% t2 = 1.0 - ( s1 * s2 ) / ( sqrt ( 3.0 ) * B * L * ( s3 - s4 ) );
  t2 = ( s1 * s2 ) ./ ( sqrt ( 3.0 ) * ( s3 - s4 ) );

  s5 = L * ( L^2 + 8.0 * w * x + 4.0 * w^2 );
  s6 = 2.0 * ( L - 2.0 * w ) * x.^2 ...
   + ( L^2 + 8.0 * L * w - 4.0 * w^2 ) * x ...
    + 2.0 * L * w^2 + L^2 * w + L^3;
  t3 = 1.0 - sqrt ( s5 ./ s6 )

  y = y_chicken .* ( 1.0 - t2 .* t3 );

  return
end

