File: surftest3d.py

package info (click to toggle)
python-scipy 0.6.0-12
  • links: PTS, VCS
  • area: main
  • in suites: lenny
  • size: 32,016 kB
  • ctags: 46,675
  • sloc: cpp: 124,854; ansic: 110,614; python: 108,664; fortran: 76,260; objc: 424; makefile: 384; sh: 10
file content (156 lines) | stat: -rw-r--r-- 5,381 bytes parent folder | download | duplicates (2)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
## 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 ( )