## Automatically adapted for scipy Oct 31, 2005 by

# Copyright (c) 1996, 1997, The Regents of the University of California.
# All rights reserved.  See Legal.htm for full text and disclaimer.

from scipy import *
from numpy.core.umath import *

GraphicsError = "GraphicsError"
import os
try:
    graphics = os.environ["PYGRAPH"]
except KeyError:
    graphics = "Gist"

if graphics [0:3] == "Nar" :
    import NarPlotter
elif graphics == "Gist" :
    import GistPlotter
else :
    raise GraphicsError , \
       graphics + " is an unknown graphics package. Check PYGRAPH " + \
          "environment variable."

if graphics [0:3] == "Nar" :
    print "This is a test of the Python interface to the Limeil Lab graphics"
    print "package, Narcisse. You  need Narcisse to be running. Fire it up by"
    print "doing setenv PORT_SERVEUR 0, typing  /dist/basis/Narcisse/bin/Narcisse &,"
    print "and then do another senetv PORT_SERVEUR to the port number which"
    print "appears in the Narcisse GUI."
else :
    print "This is a test of the Python interface to the Gist graphics package."

print "Invoke function demo () to run."

from plane import *
from surface import Surface
from graph3d import Graph3d
from mesh3d import *
from numpy.corelinalg import rand

def paws ( ) :
    i = raw_input ("Type in any string to continue; ^C to return to prompt. ")
    return

def demo () :
#  vsf = 0.
#  c = 1
    s = 1000.
    kmax = 25
    lmax = 35
    # The following computations define an interesting 3d surface.
    xr = multiply.outer (
       arange (1, kmax + 1, dtype = Float), ones (lmax, Float))
    yr = multiply.outer (
       ones (kmax, Float), arange (1, lmax + 1, dtype = Float))
    zt = 5. + xr + .2 * rand (kmax, lmax)   # ranf (xr)
    rt = 100. + yr + .2 * rand (kmax, lmax)   # ranf (yr)
    z = s * (rt + zt)
    z = z + .02 * z * rand (kmax, lmax)   # ranf (z)
    ut = rt / sqrt (rt ** 2 + zt ** 2)
    vt = zt / sqrt (rt ** 2 + zt ** 2)
    ireg =  multiply.outer ( ones (kmax), ones (lmax))
    ireg [0:1, 0:lmax] = 0
    ireg [0:kmax, 0:1] = 0
    ireg [1:15, 7:12] = 2
    ireg [1:15, 12:lmax] = 3
    ireg [3:7, 3:7] = 0
    freg = ireg + .2 * (1. - rand (kmax, lmax))  # ranf (ireg))
    freg = array (freg, Float)
    #rt [4:6, 4:6] = -1.e8
    z [3:10, 3:12] = z [3:10, 3:12] * .9
    z [5, 5] = z [5, 5] * .9
    z [17:22, 15:18] = z [17:22, 15:18] * 1.2
    z [16, 16] = z [16, 16] * 1.1

    s1 = Surface (z = z, mask = "max", opt_3d = ["wm", "i3"])
    g1 = Graph3d ( s1 , titles = "Surface with contour lines",
                   xyequal = 1., hardcopy = "talk.ps",
                   theta = 45., phi = 10., roll = 0.)
    g1.plot ( )
    paws ( )
    s1.set (opt_3d = ["wm", "f3"])
    g1.change (titles = "Flat mode")
    g1.plot ()
    paws ()
    s1.set (opt_3d = ["wm", "s3"])
    g1.change (titles = "Smooth mode")
    g1.plot ()
    paws ()
    s1.new (x = xr, y = yr, z = z - z [kmax/2, lmax/2],
            mask = "max", opt_3d = ["wm" , "i3"])
    s11 = sslice (s1, Plane (array ( [0., 0., 1.]), array ( [0., 0., 0.])))
    g1.replace (1, s11)
## g1.change ( titles = "Part of surface above xy plane", phi = 30.,
##          y_factor = 2.0,
##          axis_limits = [[0., 0.],[0., 0.], [0., 100000.]])
    g1.change ( titles = "Part of surface above xy plane",
             y_factor = 2.0,
             axis_limits = [[0., 0.],[0., 0.], [0., 100000.]])
    g1.plot ( )
    paws ( )
    s1.set (opt_3d = ["wm", "f3"])
    s11 = sslice (s1, Plane (array ( [0., 0., 1.]), array ( [0., 0., 0.])))
    g1.replace (1, s11)
    g1.change (titles = "Flat mode")
    g1.plot ()
    paws ()
    s1.set (opt_3d = ["wm", "s3"])
    s11 = sslice (s1, Plane (array ( [0., 0., 1.]), array ( [0., 0., 0.])))
    g1.replace (1, s11)
    g1.change (titles = "Smooth mode")
    g1.plot ()
    paws ()
    s1.set (z = z, c = freg, opt_3d = ["wm", "s4"])
## s1 = Surface (x = xr, y = yr, z = z, c = freg, opt_3d = ["wm", "s4"])
    g1.replace (1, s1)
    g1.change ( titles = "Surface colored by mesh values", phi = 20.,
            xyequal = 1,
            axis_limits = [[0., 0.], [0., 0.], [0., 0.], [0., 0.]])
    g1.plot ( )
    paws ( )
    xr1 = zeros ( (16, 7), Float)
    yr1 = zeros ( (16, 7), Float)
    z1 = zeros ( (16, 7), Float)
    zs1 = zeros ( (16, 7), Float)
    for i in range (16) :
        for j in range (7) :
            xr1 [i, j] = xr [i, j+6]
            yr1 [i, j] = yr [i, j+6]
            z1 [i, j] = z [i, j+6]
            zs1 [i, j] = freg [i, j+6]
    s1.set (x = xr1, y = yr1, z = z1, c = zs1)
    g1.change ( titles = "Region 2 colored by mesh values", phi = 10.)
    g1.plot ( )
    paws ( )
    s1.new (x = xr1, y = yr1, z = z1, opt_3d = ["wm", "s3"], mask = "max")
    g1.change ( titles = "Region 2 with mesh and contours", phi = 10.)
    g1.plot ( )
    paws ()
    s1.set (opt_3d = "s3")
    g1.change ( titles = "Region 2 with contours alone")
    g1.plot ()
    paws ( )
    zs1 = zeros ( (16, lmax - 7), Float)
    for i in range (16) :
        for j in range (lmax - 7) :
            zs1 [i, j] = z [i, j+6]
    s1.new (z = zs1, opt_3d = ["wm", "s3"], mask = "max")
    g1.change ( titles = "Regions 2 and 3, mesh and contours",
             theta = 70., phi = 10., roll = 0.)
    g1.plot ( )
    paws ( )
    s1.set (opt_3d = "s3")
    g1.change ( titles = "Regions 2 and 3, contours alone")
    g1.plot ( )