#! /usr/bin/env python

import math

def __main__():
  nx = 50
  ny = 50
  print_metadata(nx, ny)
  for i in range(nx):
    for j in range(ny):
      z = surface_func(i, j)
      print i, j, z
    print

## a function to generate some sample densities for mayavi to visualize
def surface_func(x, y):
  off1 = (5, 20)
  off2 = (30, 10)
  f1 = math.exp(-((x-off1[0])*(x-off1[0]) + (y-off1[1])*(y-off1[0]))/25.0)
  f2 = math.exp(-((x-off2[0])*(x-off2[0]) + (y-off2[1])*(y-off2[0]))/40.0)
  return f1 + f2

def print_metadata(nx, ny):
  print '## surface data for viewing with gnuplot'
  print '##nx:', nx
  print '##ny:', ny

__main__()