#! /usr/bin/env python3
#
def corkscrew_plot3d ( ):

#*****************************************************************************80
#
## grid_surface() draws a surface from a table Z(X,Y).
#
#  Discussion:
#
#    Here, the grid is a 41x41 array of equally spaced points in [-2,2]x[-2,2].
#
#    and Z(X,Y) is a function which can be evaluated by:
#
#      Z = exp ( - ( X^2 + Y^2 ) ) * cos ( 0.25 * X )
#        * sin ( Y ) * cos ( 2 * (X^2 + Y^2 ) )
#
#  Licensing:
#
#    This code is distributed under the MIT license.
#
#  Modified:
#
#    27 September 2026
#
#  Author:
#
#    John Burkardt
#
  import matplotlib
  import numpy as np
  import platform

  print ( '' )
  print ( 'grid_surface():' )
  print ( '  matplotlib version: ' + matplotlib.__version__ )
  print ( '  numpy version:      ' + np.version.version )
  print ( '  python version:     ' + platform.python_version ( ) )
  print ( '  Make a surface plot Z(X,Y) for a 41x41 table.' )

  header = 'grid'
#
#  Create the data file.
#
  m = 41
  n = 41

  data_filename = header + '_data.txt'
  data_unit = open ( data_filename, 'wt' )

  for i in range ( 0, m ):

    if ( 0 < i ):
      data_unit.write ( "\n" )

    x = ( ( m - i - 1 ) * ( -2.0 ) \
        + (     i     ) * (  2.0 ) ) \
        / ( m     - 1 )

    for j in range ( 0, n ):

      y = ( ( n - j - 1 ) * ( -2.0 ) \
          + (     j     ) * (  2.0 ) ) \
          / ( n     - 1 )
      r2 = x * x + y * y
      z = np.exp ( - r2 ) * np.cos ( 0.25 * x ) * np.sin ( y ) * np.cos ( 2 * r2 )
      s = ( '  %g  %g  %g\n' % ( x, y, z ) )
      data_unit.write ( s )

  data_unit.close ( )
  print ( '  Created data file "' + data_filename + '"' )
#
#  Create the command file.
#
  command_filename = header + '_commands.txt'
  command_unit = open ( command_filename, 'wt' )
  command_unit.write ( '# ' + command_filename + '\n' )
  command_unit.write ( '#\n' )
  command_unit.write ( '# Usage:\n' )
  command_unit.write ( '#  gnuplot < ' + command_filename + '\n' )
  command_unit.write ( '#\n' )
  command_unit.write ( 'set term png\n' )
  command_unit.write ( 'set output "' + header + '_surface.png"\n' )
  command_unit.write ( 'set xlabel "X"\n' )
  command_unit.write ( 'set ylabel "Y"\n' )
  command_unit.write ( 'set zlabel "Z"\n' )
  command_unit.write ( 'set title "Grid surface"\n' )
  command_unit.write ( 'set view 77, 77\n' )
  command_unit.write ( 'set xyplane -0.5\n' )
  command_unit.write ( 'set hidden3d\n' )
  command_unit.write ( 'set timestamp\n' )
  command_unit.write ( 'set grid\n' )
  command_unit.write ( 'set style data lines\n' )
  command_unit.write ( 'splot "' + data_filename + '" with lines\n' )
  command_unit.write ( 'quit\n' )
  command_unit.close ( )
  print ( '  Created command file "' + command_filename + '"' )

  return

def timestamp ( ):

#*****************************************************************************80
#
## timestamp() prints the date as a timestamp.
#
#  Licensing:
#
#    This code is distributed under the MIT license. 
#
#  Modified:
#
#    21 August 2019
#
#  Author:
#
#    John Burkardt
#
  import time

  t = time.time ( )
  print ( time.ctime ( t ) )

  return

if ( __name__ == '__main__' ):
  timestamp ( )
  corkscrew_plot3d ( )
  timestamp ( )