source: trunk/gui/scripts/mesh.tcl @ 5635

Last change on this file since 5635 was 5635, checked in by ldelgass, 9 years ago

Structured grid fixes for Rappture::Mesh
File size: 50.1 KB
Line 
1# -*- mode: tcl; indent-tabs-mode: nil -*-
2
3# ----------------------------------------------------------------------
4#  COMPONENT: mesh - represents a structured mesh for a device
5#
6#  This object represents a mesh in an XML description of simulator
7#  output.  A mesh is a structured arrangement of points, as elements
8#  composed of nodes representing coordinates.  This is a little
9#  different from a cloud, which is an unstructured arrangement
10#  (shotgun blast) of points.
11# ======================================================================
12#  AUTHOR:  Michael McLennan, Purdue University
13#  Copyright (c) 2004-2012  HUBzero Foundation, LLC
14#
15#  See the file "license.terms" for information on usage and
16#  redistribution of this file, and for a DISCLAIMER OF ALL WARRANTIES.
17# ======================================================================
18package require Itcl
19
20namespace eval Rappture {
21    # forward declaration
22}
23
24itcl::class Rappture::Mesh {
25    constructor {xmlobj path} {
26        # defined below
27    }
28    destructor {
29        # defined below
30    }
31    public method points {}
32    public method elements {}
33    public method mesh {{-type "vtk"}}
34    public method dimensions {}
35    public method limits {which}
36    public method units { axis }
37    public method label { axis }
38    public method hints {{key ""}}
39    public method isvalid {} {
40        return $_isValid
41    }
42    public proc fetch {xmlobj path}
43    public proc release {obj}
44    public method vtkdata {{what -partial}}
45    public method type {} {
46        return $_type
47    }
48    public method numpoints {} {
49        return $_numPoints
50    }
51    public method numcells {} {
52        return $_numCells
53    }
54
55    private method ReadNodesElements {path}
56    private method GetDimension { path }
57    private method GetDouble { path }
58    private method GetInt { path }
59    private method InitHints {}
60    private method ReadGrid { path }
61    private method ReadUnstructuredGrid { path }
62    private method ReadVtk { path }
63    private method WriteTriangles { path xv yv zv triangles }
64    private method WriteQuads { path xv yv zv quads }
65    private method WriteVertices { path xv yv zv vertices }
66    private method WriteLines { path xv yv zv lines }
67    private method WritePolygons { path xv yv zv polygons }
68    private method WriteTriangleStrips { path xv yv zv trianglestrips }
69    private method WriteTetrahedrons { path xv yv zv tetrahedrons }
70    private method WriteHexahedrons { path xv yv zv hexhedrons }
71    private method WriteWedges { path xv yv zv wedges }
72    private method WritePyramids { path xv yv zv pyramids }
73    private method WriteHybridCells { path xv yv zv cells celltypes }
74    private method WritePointCloud { path xv yv zv }
75    private method GetCellType { name }
76    private method GetNumIndices { type }
77
78    private variable _xmlobj "";        # Ref to XML obj with device data
79    private variable _mesh "";          # Lib obj representing this mesh
80    private variable _dim 0;            # Dimension of mesh (1, 2, or 3)
81    private variable _type "";          # Indicates the type of mesh.
82    private variable _axis2units;       # System of units for x, y, z
83    private variable _axis2labels;      #
84    private variable _hints
85    private variable _limits;           # Array of mesh limits. Keys are
86                                        # xmin, xmax, ymin, ymax, ...
87    private variable _numPoints 0;      # Number of points in mesh
88    private variable _numCells 0;       # Number of cells in mesh
89    private variable _vtkdata "";       # Mesh in vtk file format.
90    private variable _isValid 0;        # Indicates if the mesh is valid.
91
92    private common _xp2obj;             # used for fetch/release ref counting
93    private common _obj2ref;            # used for fetch/release ref counting
94}
95
96# ----------------------------------------------------------------------
97# USAGE: Rappture::Mesh::fetch <xmlobj> <path>
98#
99# Clients use this instead of a constructor to fetch the Mesh for
100# a particular <path> in the <xmlobj>.  When the client is done with
101# the mesh, he calls "release" to decrement the reference count.
102# When the mesh is no longer needed, it is cleaned up automatically.
103# ----------------------------------------------------------------------
104itcl::body Rappture::Mesh::fetch {xmlobj path} {
105    set handle "$xmlobj|$path"
106    if {[info exists _xp2obj($handle)]} {
107        set obj $_xp2obj($handle)
108        incr _obj2ref($obj)
109        return $obj
110    }
111    set obj [Rappture::Mesh ::#auto $xmlobj $path]
112    set _xp2obj($handle) $obj
113    set _obj2ref($obj) 1
114    return $obj
115}
116
117# ----------------------------------------------------------------------
118# USAGE: Rappture::Mesh::release <obj>
119#
120# Clients call this when they're no longer using a Mesh fetched
121# previously by the "fetch" proc.  This decrements the reference
122# count for the mesh and destroys the object when it is no longer
123# in use.
124# ----------------------------------------------------------------------
125itcl::body Rappture::Mesh::release {obj} {
126    if {[info exists _obj2ref($obj)]} {
127        incr _obj2ref($obj) -1
128        if {$_obj2ref($obj) <= 0} {
129            unset _obj2ref($obj)
130            foreach handle [array names _xp2obj] {
131                if {$_xp2obj($handle) == $obj} {
132                    unset _xp2obj($handle)
133                }
134            }
135            itcl::delete object $obj
136        }
137    } else {
138        error "can't find reference count for $obj"
139    }
140}
141
142# ----------------------------------------------------------------------
143# CONSTRUCTOR
144# ----------------------------------------------------------------------
145itcl::body Rappture::Mesh::constructor {xmlobj path} {
146    package require vtk
147    if {![Rappture::library isvalid $xmlobj]} {
148        error "bad value \"$xmlobj\": should be Rappture::library"
149    }
150    set _xmlobj $xmlobj
151    set _mesh [$xmlobj element -as object $path]
152
153    # Initialize mesh bounds to empty
154    foreach axis {x y z} {
155        set _limits($axis) ""
156    }
157    set units [$_mesh get units]
158    set first [lindex $units 0]
159    foreach u $units axis { x y z } {
160        if { $u != "" } {
161            set _axis2units($axis) $u
162        } else {
163            set _axis2units($axis) $first
164        }
165    }
166    foreach label [$_mesh get labels] axis { x y z } {
167        if { $label != "" } {
168            set _axis2labels($axis) $label
169        } else {
170            set _axis2labels($axis) [string toupper $axis]
171        }
172    }
173
174    # Meshes comes in a variety of flavors
175    #
176    # Dimensionality is determined from the <dimension> tag.
177    #
178    # <vtk> described mesh
179    # <element> +  <node> definitions
180    # <grid>            rectangular mesh
181    # <unstructured>    homogeneous cell type mesh.
182
183    # Check that only one mesh type was defined.
184    set subcount 0
185    foreach cname [$_mesh children] {
186        foreach type { vtk grid unstructured } {
187            if { $cname == $type } {
188                incr subcount
189                break
190            }
191        }
192    }
193    if {[$_mesh element "node"] != "" ||
194        [$_mesh element "element"] != ""} {
195        incr subcount
196    }
197
198    if { $subcount == 0 } {
199        puts stderr "WARNING: no mesh specified for \"$path\"."
200        return
201    }
202    if { $subcount > 1 } {
203        puts stderr "WARNING: too many mesh types specified for \"$path\"."
204        return
205    }
206    set result 0
207    if { [$_mesh element "vtk"] != ""} {
208        set result [ReadVtk $path]
209    } elseif {[$_mesh element "grid"] != "" } {
210        set result [ReadGrid $path]
211    } elseif {[$_mesh element "unstructured"] != "" } {
212        set result [ReadUnstructuredGrid $path]
213    } elseif {[$_mesh element "node"] != "" && [$_mesh element "element"] != ""} {
214        set result [ReadNodesElements $path]
215    }
216    set _isValid $result
217    InitHints
218}
219
220# ----------------------------------------------------------------------
221# DESTRUCTOR
222# ----------------------------------------------------------------------
223itcl::body Rappture::Mesh::destructor {} {
224    # don't destroy the _xmlobj! we don't own it!
225    itcl::delete object $_mesh
226}
227
228#
229# vtkdata --
230#
231#        This is called by the field object to generate a VTK file to send to
232#        the remote render server.  Returns the vtkDataSet object containing
233#        (at this point) just the mesh.  The field object doesn't know (or
234#        care) what type of mesh is used.  The field object will add field
235#        arrays before generating output to send to the remote render server.
236#
237itcl::body Rappture::Mesh::vtkdata {{what -partial}} {
238    if {$what == "-full"} {
239        append out "# vtk DataFile Version 3.0\n"
240        append out "[hints label]\n"
241        append out "ASCII\n"
242        append out $_vtkdata
243        return $out
244    } else {
245        return $_vtkdata
246    }
247}
248
249# ----------------------------------------------------------------------
250# USAGE: points
251#
252# Returns the vtk object containing the points for this mesh.
253# ----------------------------------------------------------------------
254itcl::body Rappture::Mesh::points {} {
255    return ""
256}
257
258#
259# units --
260#
261#       Returns the units of the given axis.
262#
263itcl::body Rappture::Mesh::units { axis } {
264    if { ![info exists _axis2units($axis)] } {
265        return ""
266    }
267    return $_axis2units($axis)
268}
269
270#
271# label --
272#
273#       Returns the label of the given axis.
274#
275itcl::body Rappture::Mesh::label { axis } {
276    if { ![info exists _axis2labels($axis)] } {
277        return ""
278    }
279    return $_axis2labels($axis)
280}
281
282# ----------------------------------------------------------------------
283# USAGE: elements
284#
285# Returns a list of the form {plist r plist r ...} for each element
286# in this mesh.  Each plist is a list of {x y x y ...} coordinates
287# for the mesh.
288# ----------------------------------------------------------------------
289itcl::body Rappture::Mesh::elements {} {
290    # build a map for region number => region type
291    foreach comp [$_mesh children -type region] {
292        set id [$_mesh element -as id $comp]
293        set regions($id) [$_mesh get $comp]
294    }
295    set regions() "unknown"
296
297    set rlist ""
298    foreach comp [$_mesh children -type element] {
299        set rid [$_mesh get $comp.region]
300
301        #
302        # HACK ALERT!
303        #
304        # Prophet puts out nodes in a funny "Z" shaped order,
305        # not in proper clockwise fashion.  Switch the last
306        # two nodes for now to make them correct.
307        #
308        set nlist [$_mesh get $comp.nodes]
309        set n2 [lindex $nlist 2]
310        set n3 [lindex $nlist 3]
311        set nlist [lreplace $nlist 2 3 $n3 $n2]
312        lappend nlist [lindex $nlist 0]
313
314        set plist ""
315        foreach nid $nlist {
316            eval lappend plist [$_mesh get node($nid)]
317        }
318        lappend rlist $plist $regions($rid)
319    }
320    return $rlist
321}
322
323# ----------------------------------------------------------------------
324# USAGE: mesh
325#
326# Returns the vtk object representing the mesh.
327# ----------------------------------------------------------------------
328itcl::body Rappture::Mesh::mesh { {type "vtk"} } {
329    switch $type {
330        "vtk" {
331            return ""
332        }
333        default {
334            error "Requested mesh type \"$type\" is unknown."
335        }
336    }
337}
338
339# ----------------------------------------------------------------------
340# USAGE: dimensions
341#
342# Returns the number of dimensions for this object: 1, 2, or 3.
343# ----------------------------------------------------------------------
344itcl::body Rappture::Mesh::dimensions {} {
345    return $_dim
346}
347
348# ----------------------------------------------------------------------
349# USAGE: limits x|y|z
350#
351# Returns the {min max} coords for the limits of the specified axis.
352# ----------------------------------------------------------------------
353itcl::body Rappture::Mesh::limits {axis} {
354    if {![info exists _limits($axis)]} {
355        error "bad axis \"$which\": should be x, y, z"
356    }
357    return $_limits($axis)
358}
359
360# ----------------------------------------------------------------------
361# USAGE: hints ?<keyword>?
362#
363# Returns a list of key/value pairs for various hints about plotting
364# this field.  If a particular <keyword> is specified, then it returns
365# the hint for that <keyword>, if it exists.
366# ----------------------------------------------------------------------
367itcl::body Rappture::Mesh::hints {{keyword ""}} {
368    if {$keyword != ""} {
369        if {[info exists _hints($keyword)]} {
370            return $_hints($keyword)
371        }
372        return ""
373    }
374    return [array get _hints]
375}
376
377# ----------------------------------------------------------------------
378# USAGE: InitHints
379#
380# Returns a list of key/value pairs for various hints about plotting
381# this mesh.  If a particular <keyword> is specified, then it returns
382# the hint for that <keyword>, if it exists.
383# ----------------------------------------------------------------------
384itcl::body Rappture::Mesh::InitHints {} {
385    foreach {key path} {
386        camera       camera.position
387        color        about.color
388        label        about.label
389        style        about.style
390        units        units
391    } {
392        set str [$_mesh get $path]
393        if {"" != $str} {
394            set _hints($key) $str
395        }
396    }
397    foreach {key path} {
398        toolid          tool.id
399        toolname        tool.name
400        toolcommand     tool.execute
401        tooltitle       tool.title
402        toolrevision    tool.version.application.revision
403    } {
404        set str [$_xmlobj get $path]
405        if { "" != $str } {
406            set _hints($key) $str
407        }
408    }
409}
410
411itcl::body Rappture::Mesh::GetDimension { path } {
412    set string [$_xmlobj get $path.dim]
413    if { $string == "" } {
414        puts stderr "WARNING: no tag <dim> found in mesh \"$path\"."
415        return 0
416    }
417    if { [scan $string "%d" _dim] == 1 } {
418        if { $_dim == 1 || $_dim == 2 || $_dim == 3 } {
419            return 1
420        }
421    }
422    puts stderr "WARNING: bad <dim> tag value \"$string\": should be 1, 2 or 3."
423    return 0
424}
425
426itcl::body Rappture::Mesh::GetDouble { path } {
427    set string [$_xmlobj get $path]
428    if { [scan $string "%g" value] != 1 } {
429        puts stderr "ERROR: can't get double value \"$string\" of \"$path\""
430        return 0.0
431    }
432    return $value
433}
434
435itcl::body Rappture::Mesh::GetInt { path } {
436    set string [$_xmlobj get $path]
437    if { [scan $string "%d" value] != 1 } {
438        puts stderr "ERROR: can't get integer value \"$string\" of \"$path\""
439        return 0.0
440    }
441    return $value
442}
443
444itcl::body Rappture::Mesh::ReadVtk { path } {
445    set _type "vtk"
446
447    if { ![GetDimension $path] } {
448        return 0
449    }
450    # Create a VTK file with the mesh in it.
451    set _vtkdata [$_xmlobj get $path.vtk]
452    append out "# vtk DataFile Version 3.0\n"
453    append out "mesh\n"
454    append out "ASCII\n"
455    append out "$_vtkdata\n"
456
457     # Write the contents to a file just in case it's binary.
458    set tmpfile file[pid].vtk
459    set f [open "$tmpfile" "w"]
460    fconfigure $f -translation binary -encoding binary
461    puts $f $out
462    close $f
463
464    # Read the data back into a vtk dataset and query the bounds.
465    set reader $this-datasetreader
466    vtkDataSetReader $reader
467    $reader SetFileName $tmpfile
468    $reader Update
469    set output [$reader GetOutput]
470    set _numPoints [$output GetNumberOfPoints]
471    set _numCells [$output GetNumberOfCells]
472    foreach { xmin xmax ymin ymax zmin zmax } [$output GetBounds] break
473    set _limits(x) [list $xmin $xmax]
474    set _limits(y) [list $ymin $ymax]
475    set _limits(z) [list $zmin $zmax]
476    file delete $tmpfile
477    rename $output ""
478    rename $reader ""
479    return 1
480}
481
482itcl::body Rappture::Mesh::ReadGrid { path } {
483    set _type "grid"
484
485    if { ![GetDimension $path] } {
486        return 0
487    }
488    set numUniform 0
489    set numRectilinear 0
490    set numCurvilinear 0
491    foreach axis { x y z } {
492        set min    [$_xmlobj get "$path.grid.${axis}axis.min"]
493        set max    [$_xmlobj get "$path.grid.${axis}axis.max"]
494        set num    [$_xmlobj get "$path.grid.${axis}axis.numpoints"]
495        set coords [$_xmlobj get "$path.grid.${axis}coords"]
496        set dim    [$_xmlobj get "$path.grid.${axis}dim"]
497        if { $min != "" && $max != "" && $num != "" && $num > 0 } {
498            set ${axis}Min $min
499            set ${axis}Max $max
500            set ${axis}Num $num
501            if {$min > $max} {
502                puts stderr "ERROR: grid $axis axis minimum larger than maximum"
503                return 0
504            }
505            incr numUniform
506        } elseif { $coords != "" } {
507            incr numRectilinear
508            set ${axis}Coords $coords
509        } elseif { $dim != "" } {
510            set ${axis}Num $dim
511            incr numCurvilinear
512        }
513    }
514    set _dim [expr $numRectilinear + $numUniform + $numCurvilinear]
515    if { $_dim == 0 } {
516        # No data found.
517        puts stderr "WARNING: bad grid \"$path\": no data found"
518        return 0
519    }
520    if { $numCurvilinear > 0 } {
521        # This is the 2D/3D curilinear case. We found <xdim>, <ydim>, or <zdim>
522        if { $numRectilinear > 0 || $numUniform > 0 } {
523            puts stderr "WARNING: bad grid \"$path\": can't mix curvilinear and rectilinear grids."
524            return 0
525        }
526        set points [$_xmlobj get $path.grid.points]
527        if { $points == "" } {
528            puts stderr "WARNING: bad grid \"$path\": no <points> found."
529            return 0
530        }
531        if { ![info exists xNum] } {
532            puts stderr "WARNING: bad grid \"$path\": invalid dimensions for curvilinear grid: missing <xdim> from grid description."
533            return 0
534        }
535        set all [blt::vector create \#auto]
536        set xv [blt::vector create \#auto]
537        set yv [blt::vector create \#auto]
538        set zv [blt::vector create \#auto]
539        $all set $points
540        set numCoords [$all length]
541        if { [info exists zNum] } {
542            if { ![info exists yNum] || ![info exists xNum] } {
543                puts stderr "WARNING: bad grid \"$path\": missing grid dimension"
544                blt::vector destroy $all $xv $yv $zv
545                return 0
546            }
547            set _dim 3
548            set _numPoints [expr $xNum * $yNum * $zNum]
549            set _numCells [expr ($xNum > 1 ? ($xNum - 1) : 1) * ($yNum > 1 ? ($yNum - 1) : 1) * ($zNum > 1 ? ($zNum - 1) : 1)]
550            if { ($_numPoints*3) != $numCoords } {
551                puts stderr "WARNING: bad grid \"$path\": invalid grid: \# of points does not match dimensions <xdim> * <ydim> * <zdim>"
552                blt::vector destroy $all $xv $yv $zv
553                return 0
554            }
555            if { ($numCoords % 3) != 0 } {
556                puts stderr "WARNING: bad grid \"$path\": wrong \# of coordinates for 3D grid"
557                blt::vector destroy $all $xv $yv $zv
558                return 0
559            }
560            $all split $xv $yv $zv
561            foreach axis {x y z} {
562                set vector [set ${axis}v]
563                set _limits($axis) [$vector limits]
564            }
565            append out "DATASET STRUCTURED_GRID\n"
566            append out "DIMENSIONS $xNum $yNum $zNum\n"
567            append out "POINTS $_numPoints double\n"
568            foreach x [$xv range 0 end] y [$yv range 0 end] z [$zv range 0 end] {
569                append out "$x $y $z\n"
570            }
571            set _vtkdata $out
572        } elseif { [info exists yNum] } {
573            if { ![info exists xNum] } {
574                puts stderr "WARNING: bad grid \"$path\": missing grid dimension"
575                blt::vector destroy $all $xv $yv $zv
576                return 0
577            }
578            set _dim 2
579            set _numPoints [expr $xNum * $yNum]
580            set _numCells [expr ($xNum > 1 ? ($xNum - 1) : 1) * ($yNum > 1 ? ($yNum - 1) : 1)]
581            if { ($_numPoints*2) != $numCoords } {
582                puts stderr "WARNING: bad grid \"$path\": \# of points does not match dimensions <xdim> * <ydim>"
583                blt::vector destroy $all $xv $yv $zv
584                return 0
585            }
586            if { ($numCoords % 2) != 0 } {
587                puts stderr "WARNING: bad grid \"$path\": wrong \# of coordinates for 2D grid"
588                blt::vector destroy $all $xv $yv $zv
589                return 0
590            }
591            $all split $xv $yv
592            foreach axis {x y} {
593                set vector [set ${axis}v]
594                set _limits($axis) [$vector limits]
595            }
596            set _limits(z) [list 0 0]
597            $zv seq 0 0 [$xv length]
598            $all merge $xv $yv $zv
599            append out "DATASET STRUCTURED_GRID\n"
600            append out "DIMENSIONS $xNum $yNum 1\n"
601            append out "POINTS $_numPoints double\n"
602            foreach x [$xv range 0 end] y [$yv range 0 end] z [$zv range 0 end] {
603                append out "$x $y $z\n"
604            }
605            set _vtkdata $out
606        } else {
607            set _dim 1
608            set _numPoints $xNum
609            set _numCells [expr $xNum - 1]
610            if { $_numPoints != $numCoords } {
611                puts stderr "WARNING: bad grid \"$path\": \# of points does not match <xdim>"
612                blt::vector destroy $all $xv $yv $zv
613                return 0
614            }
615            $all dup $xv
616            set _limits(x) [$xv limits]
617            set _limits(y) [list 0 0]
618            set _limits(z) [list 0 0]
619            $yv seq 0 0 [$xv length]
620            $zv seq 0 0 [$xv length]
621            $all merge $xv $yv $zv
622            append out "DATASET STRUCTURED_GRID\n"
623            append out "DIMENSIONS $xNum 1 1\n"
624            append out "POINTS $_numPoints double\n"
625            foreach x [$xv range 0 end] y [$yv range 0 end] z [$zv range 0 end] {
626                append out "$x $y $z\n"
627            }
628            set _vtkdata $out
629        }
630        blt::vector destroy $all $xv $yv $zv
631        return 1
632    }
633    if { $numRectilinear == 0 && $numUniform > 0} {
634        # This is the special case where all axes 2D/3D are uniform.
635        # This results in a STRUCTURED_POINTS
636        if { $_dim == 1 } {
637            if {$xNum == 1} {
638                set xSpace 0
639            } else {
640                set xSpace [expr ($xMax - $xMin) / double($xNum - 1)]
641            }
642            set _numPoints $xNum
643            set _numCells [expr $xNum - 1]
644            append out "DATASET STRUCTURED_POINTS\n"
645            append out "DIMENSIONS $xNum 1 1\n"
646            append out "ORIGIN $xMin 0 0\n"
647            append out "SPACING $xSpace 0 0\n"
648            set _vtkdata $out
649            set _limits(x) [list $xMin $xMax]
650            set _limits(y) [list 0 0]
651            set _limits(z) [list 0 0]
652        } elseif { $_dim == 2 } {
653            if {$xNum == 1} {
654                set xSpace 0
655            } else {
656                set xSpace [expr ($xMax - $xMin) / double($xNum - 1)]
657            }
658            if {$yNum == 1} {
659                set ySpace 0
660            } else {
661                set ySpace [expr ($yMax - $yMin) / double($yNum - 1)]
662            }
663            set _numPoints [expr $xNum * $yNum]
664            set _numCells [expr ($xNum > 1 ? ($xNum - 1) : 1) * ($yNum > 1 ? ($yNum - 1) : 1)]
665            append out "DATASET STRUCTURED_POINTS\n"
666            append out "DIMENSIONS $xNum $yNum 1\n"
667            append out "ORIGIN $xMin $yMin 0\n"
668            append out "SPACING $xSpace $ySpace 0\n"
669            set _vtkdata $out
670            foreach axis {x y} {
671                set _limits($axis) [list [set ${axis}Min] [set ${axis}Max]]
672            }
673            set _limits(z) [list 0 0]
674        } elseif { $_dim == 3 } {
675            if {$xNum == 1} {
676                set xSpace 0
677            } else {
678                set xSpace [expr ($xMax - $xMin) / double($xNum - 1)]
679            }
680            if {$yNum == 1} {
681                set ySpace 0
682            } else {
683                set ySpace [expr ($yMax - $yMin) / double($yNum - 1)]
684            }
685            if {$zNum == 1} {
686                set zSpace 0
687            } else {
688                set zSpace [expr ($zMax - $zMin) / double($zNum - 1)]
689            }
690            set _numPoints [expr $xNum * $yNum * $zNum]
691            set _numCells [expr ($xNum > 1 ? ($xNum - 1) : 1) * ($yNum > 1 ? ($yNum - 1) : 1) * ($zNum > 1 ? ($zNum - 1) : 1)]
692            append out "DATASET STRUCTURED_POINTS\n"
693            append out "DIMENSIONS $xNum $yNum $zNum\n"
694            append out "ORIGIN $xMin $yMin $zMin\n"
695            append out "SPACING $xSpace $ySpace $zSpace\n"
696            set _vtkdata $out
697            foreach axis {x y z} {
698                set _limits($axis) [list [set ${axis}Min] [set ${axis}Max]]
699            }
700        } else {
701            puts stderr "WARNING: bad grid \"$path\": bad dimension \"$_dim\""
702            return 0
703        }
704        return 1
705    }
706    # This is the hybrid case.  Some axes are uniform, others are nonuniform.
707    set xv [blt::vector create \#auto]
708    if { [info exists xMin] } {
709        $xv seq $xMin $xMax $xNum
710    } else {
711        $xv set [$_xmlobj get $path.grid.xcoords]
712        set xMin [$xv min]
713        set xMax [$xv max]
714        set xNum [$xv length]
715    }
716    set yv [blt::vector create \#auto]
717    if { $_dim > 1 } {
718        if { [info exists yMin] } {
719            $yv seq $yMin $yMax $yNum
720        } else {
721            $yv set [$_xmlobj get $path.grid.ycoords]
722            set yMin [$yv min]
723            set yMax [$yv max]
724            set yNum [$yv length]
725        }
726    } else {
727        set yNum 1
728    }
729    set zv [blt::vector create \#auto]
730    if { $_dim == 3 } {
731        if { [info exists zMin] } {
732            $zv seq $zMin $zMax $zNum
733        }  else {
734            $zv set [$_xmlobj get $path.grid.zcoords]
735            set zMin [$zv min]
736            set zMax [$zv max]
737            set zNum [$zv length]
738        }
739    } else {
740        set zNum 1
741    }
742    if { $_dim == 3 } {
743        set _numPoints [expr $xNum * $yNum * $zNum]
744        set _numCells [expr ($xNum > 1 ? ($xNum - 1) : 1) * ($yNum > 1 ? ($yNum - 1) : 1) * ($zNum > 1 ? ($zNum - 1) : 1)]
745        append out "DATASET RECTILINEAR_GRID\n"
746        append out "DIMENSIONS $xNum $yNum $zNum\n"
747        append out "X_COORDINATES $xNum double\n"
748        append out [$xv range 0 end]
749        append out "\n"
750        append out "Y_COORDINATES $yNum double\n"
751        append out [$yv range 0 end]
752        append out "\n"
753        append out "Z_COORDINATES $zNum double\n"
754        append out [$zv range 0 end]
755        append out "\n"
756        set _vtkdata $out
757        foreach axis {x y z} {
758            if { [info exists ${axis}Min] } {
759                set _limits($axis) [list [set ${axis}Min] [set ${axis}Max]]
760            }
761        }
762    } elseif { $_dim == 2 } {
763        set _numPoints [expr $xNum * $yNum]
764        set _numCells [expr ($xNum > 1 ? ($xNum - 1) : 1) * ($yNum > 1 ? ($yNum - 1) : 1)]
765        append out "DATASET RECTILINEAR_GRID\n"
766        append out "DIMENSIONS $xNum $yNum 1\n"
767        append out "X_COORDINATES $xNum double\n"
768        append out [$xv range 0 end]
769        append out "\n"
770        append out "Y_COORDINATES $yNum double\n"
771        append out [$yv range 0 end]
772        append out "\n"
773        append out "Z_COORDINATES 1 double\n"
774        append out "0\n"
775        foreach axis {x y} {
776            if { [info exists ${axis}Min] } {
777                set _limits($axis) [list [set ${axis}Min] [set ${axis}Max]]
778            }
779        }
780        set _limits(z) [list 0 0]
781        set _vtkdata $out
782    } elseif { $_dim == 1 } {
783        set _numPoints $xNum
784        set _numCells [expr $xNum - 1]
785        append out "DATASET RECTILINEAR_GRID\n"
786        append out "DIMENSIONS $xNum 1 1\n"
787        append out "X_COORDINATES $xNum double\n"
788        append out [$xv range 0 end]
789        append out "\n"
790        append out "Y_COORDINATES 1 double\n"
791        append out "0\n"
792        append out "Z_COORDINATES 1 double\n"
793        append out "0\n"
794        if { [info exists xMin] } {
795            set _limits(x) [list $xMin $xMax]
796        }
797        set _limits(y) [list 0 0]
798        set _limits(z) [list 0 0]
799        set _vtkdata $out
800    } else {
801        puts stderr "WARNING: bad grid \"$path\": invalid dimension \"$_dim\""
802        return 0
803    }
804    blt::vector destroy $xv $yv $zv
805    return 1
806}
807
808itcl::body Rappture::Mesh::WritePointCloud { path xv yv zv } {
809    set _type "cloud"
810    set _numPoints [$xv length]
811    append out "DATASET POLYDATA\n"
812    append out "POINTS $_numPoints double\n"
813    foreach x [$xv range 0 end] y [$yv range 0 end] z [$zv range 0 end] {
814        append out "$x $y $z\n"
815    }
816    set _vtkdata $out
817    set _limits(x) [$xv limits]
818    if { $_dim > 1 } {
819        set _limits(y) [$yv limits]
820    } else {
821        set _limits(y) [list 0 0]
822    }
823    if { $_dim == 3 } {
824        set _limits(z) [$zv limits]
825    } else {
826        set _limits(z) [list 0 0]
827    }
828    return 1
829}
830
831itcl::body Rappture::Mesh::WriteTriangles { path xv yv zv triangles } {
832    set _type "triangles"
833    set _numPoints [$xv length]
834    set _numCells 0
835    set data {}
836    set celltypes {}
837    foreach { a b c } $triangles {
838        append data "3 $a $b $c\n"
839        append celltypes "5\n"
840        incr _numCells
841    }
842    append out "DATASET UNSTRUCTURED_GRID\n"
843    append out "POINTS $_numPoints double\n"
844    foreach x [$xv range 0 end] y [$yv range 0 end] z [$zv range 0 end] {
845        append out "$x $y $z\n"
846    }
847    set count [expr $_numCells * 4]
848    append out "CELLS $_numCells $count\n"
849    append out $data
850    append out "CELL_TYPES $_numCells\n"
851    append out $celltypes
852    set _limits(x) [$xv limits]
853    set _limits(y) [$yv limits]
854    if { $_dim == 3 } {
855        set _limits(z) [$zv limits]
856    } else {
857        set _limits(z) [list 0 0]
858    }
859    set _vtkdata $out
860    return 1
861}
862
863itcl::body Rappture::Mesh::WriteQuads { path xv yv zv quads } {
864    set _type "quads"
865    set _numPoints [$xv length]
866    set _numCells 0
867    set data {}
868    set celltypes {}
869    foreach { a b c d } $quads {
870        append data "4 $a $b $c $d\n"
871        append celltypes "9\n"
872        incr _numCells
873    }
874    append out "DATASET UNSTRUCTURED_GRID\n"
875    append out "POINTS $_numPoints double\n"
876    foreach x [$xv range 0 end] y [$yv range 0 end] z [$zv range 0 end] {
877        append out "$x $y $z\n"
878    }
879    set count [expr $_numCells * 5]
880    append out "CELLS $_numCells $count\n"
881    append out $data
882    append out "CELL_TYPES $_numCells\n"
883    append out $celltypes
884    set _limits(x) [$xv limits]
885    set _limits(y) [$yv limits]
886    if { $_dim == 3 } {
887        set _limits(z) [$zv limits]
888    } else {
889        set _limits(z) [list 0 0]
890    }
891    set _vtkdata $out
892    return 1
893}
894
895itcl::body Rappture::Mesh::WriteVertices { path xv yv zv vertices } {
896    set _type "vertices"
897    set _numPoints [$xv length]
898    set _numCells 0
899    set data {}
900    set lines [split $vertices \n]
901    set count 0
902    foreach { line } $lines {
903        set numIndices [llength $line]
904        if { $numIndices == 0 } {
905            continue
906        }
907        append data "$numIndices $line\n"
908        incr _numCells
909        set count [expr $count + $numIndices + 1]
910    }
911    append out "DATASET POLYDATA\n"
912    append out "POINTS $_numPoints double\n"
913    foreach x [$xv range 0 end] y [$yv range 0 end] z [$zv range 0 end] {
914        append out "$x $y $z\n"
915    }
916    append out "VERTICES $_numCells $count\n"
917    append out $data
918    set _limits(x) [$xv limits]
919    set _limits(y) [$yv limits]
920    if { $_dim == 3 } {
921        set _limits(z) [$zv limits]
922    } else {
923        set _limits(z) [list 0 0]
924    }
925    set _vtkdata $out
926    return 1
927}
928
929itcl::body Rappture::Mesh::WriteLines { path xv yv zv polylines } {
930    set _type "lines"
931    set _numPoints [$xv length]
932    set _numCells 0
933    set data {}
934    set lines [split $polylines \n]
935    set count 0
936    foreach { line } $lines {
937        set numIndices [llength $line]
938        if { $numIndices == 0 } {
939            continue
940        }
941        append data "$numIndices $line\n"
942        incr _numCells
943        set count [expr $count + $numIndices + 1]
944    }
945    append out "DATASET POLYDATA\n"
946    append out "POINTS $_numPoints double\n"
947    foreach x [$xv range 0 end] y [$yv range 0 end] z [$zv range 0 end] {
948        append out "$x $y $z\n"
949    }
950    append out "LINES $_numCells $count\n"
951    append out $data
952    set _limits(x) [$xv limits]
953    set _limits(y) [$yv limits]
954    if { $_dim == 3 } {
955        set _limits(z) [$zv limits]
956    } else {
957        set _limits(z) [list 0 0]
958    }
959    set _vtkdata $out
960    return 1
961}
962
963itcl::body Rappture::Mesh::WritePolygons { path xv yv zv polygons } {
964    set _type "polygons"
965    set _numPoints [$xv length]
966    set _numCells 0
967    set data {}
968    set lines [split $polygons \n]
969    set count 0
970    foreach { line } $lines {
971        set numIndices [llength $line]
972        if { $numIndices == 0 } {
973            continue
974        }
975        append data "$numIndices $line\n"
976        incr _numCells
977        set count [expr $count + $numIndices + 1]
978    }
979    append out "DATASET POLYDATA\n"
980    append out "POINTS $_numPoints double\n"
981    foreach x [$xv range 0 end] y [$yv range 0 end] z [$zv range 0 end] {
982        append out "$x $y $z\n"
983    }
984    append out "POLYGONS $_numCells $count\n"
985    append out $data
986    set _limits(x) [$xv limits]
987    set _limits(y) [$yv limits]
988    if { $_dim == 3 } {
989        set _limits(z) [$zv limits]
990    } else {
991        set _limits(z) [list 0 0]
992    }
993    set _vtkdata $out
994    return 1
995}
996
997itcl::body Rappture::Mesh::WriteTriangleStrips { path xv yv zv trianglestrips } {
998    set _type "trianglestrips"
999    set _numPoints [$xv length]
1000    set _numCells 0
1001    set data {}
1002    set lines [split $trianglestrips \n]
1003    set count 0
1004    foreach { line } $lines {
1005        set numIndices [llength $line]
1006        if { $numIndices == 0 } {
1007            continue
1008        }
1009        append data "$numIndices $line\n"
1010        incr _numCells
1011        set count [expr $count + $numIndices + 1]
1012    }
1013    append out "DATASET POLYDATA\n"
1014    append out "POINTS $_numPoints double\n"
1015    foreach x [$xv range 0 end] y [$yv range 0 end] z [$zv range 0 end] {
1016        append out "$x $y $z\n"
1017    }
1018    append out "TRIANGLE_STRIPS $_numCells $count\n"
1019    append out $data
1020    set _limits(x) [$xv limits]
1021    set _limits(y) [$yv limits]
1022    if { $_dim == 3 } {
1023        set _limits(z) [$zv limits]
1024    } else {
1025        set _limits(z) [list 0 0]
1026    }
1027    set _vtkdata $out
1028    return 1
1029}
1030
1031itcl::body Rappture::Mesh::WriteTetrahedrons { path xv yv zv tetras } {
1032    set _type "tetrahedrons"
1033    set _numPoints [$xv length]
1034    set _numCells 0
1035    set data {}
1036    set celltypes {}
1037    foreach { a b c d } $tetras {
1038        append data "4 $a $b $c $d\n"
1039        append celltypes "10\n"
1040        incr _numCells
1041    }
1042    append out "DATASET UNSTRUCTURED_GRID\n"
1043    append out "POINTS $_numPoints double\n"
1044    foreach x [$xv range 0 end] y [$yv range 0 end] z [$zv range 0 end] {
1045        append out "$x $y $z\n"
1046    }
1047    set count [expr $_numCells * 5]
1048    append out "CELLS $_numCells $count\n"
1049    append out $data
1050    append out "CELL_TYPES $_numCells\n"
1051    append out $celltypes
1052    set _limits(x) [$xv limits]
1053    set _limits(y) [$yv limits]
1054    set _limits(z) [$zv limits]
1055
1056    set _vtkdata $out
1057    return 1
1058}
1059
1060itcl::body Rappture::Mesh::WriteHexahedrons { path xv yv zv hexas } {
1061    set _type "hexahedrons"
1062    set _numPoints [$xv length]
1063    set _numCells 0
1064    set data {}
1065    set celltypes {}
1066    foreach { a b c d e f g h } $hexas {
1067        append data "8 $a $b $c $d $e $f $g $h\n"
1068        append celltypes "12\n"
1069        incr _numCells
1070    }
1071    append out "DATASET UNSTRUCTURED_GRID\n"
1072    append out "POINTS $_numPoints double\n"
1073    foreach x [$xv range 0 end] y [$yv range 0 end] z [$zv range 0 end] {
1074        append out "$x $y $z\n"
1075    }
1076    set count [expr $_numCells * 9]
1077    append out "CELLS $_numCells $count\n"
1078    append out $data
1079    append out "CELL_TYPES $_numCells\n"
1080    append out $celltypes
1081    set _limits(x) [$xv limits]
1082    set _limits(y) [$yv limits]
1083    set _limits(z) [$zv limits]
1084
1085    set _vtkdata $out
1086    return 1
1087}
1088
1089itcl::body Rappture::Mesh::WriteWedges { path xv yv zv wedges } {
1090    set _type "wedges"
1091    set _numPoints [$xv length]
1092    set _numCells 0
1093    set data {}
1094    set celltypes {}
1095    foreach { a b c d e f } $wedges {
1096        append data "6 $a $b $c $d $e $f\n"
1097        append celltypes "13\n"
1098        incr _numCells
1099    }
1100    append out "DATASET UNSTRUCTURED_GRID\n"
1101    append out "POINTS $_numPoints double\n"
1102    foreach x [$xv range 0 end] y [$yv range 0 end] z [$zv range 0 end] {
1103        append out "$x $y $z\n"
1104    }
1105    set count [expr $_numCells * 7]
1106    append out "CELLS $_numCells $count\n"
1107    append out $data
1108    append out "CELL_TYPES $_numCells\n"
1109    append out $celltypes
1110    set _limits(x) [$xv limits]
1111    set _limits(y) [$yv limits]
1112    set _limits(z) [$zv limits]
1113
1114    set _vtkdata $out
1115    return 1
1116}
1117
1118itcl::body Rappture::Mesh::WritePyramids { path xv yv zv pyramids } {
1119    set _type "pyramids"
1120    set _numPoints [$xv length]
1121    set _numCells 0
1122    set data {}
1123    set celltypes {}
1124    foreach { a b c d e } $pyramids {
1125        append data "5 $a $b $c $d $e\n"
1126        append celltypes "14\n"
1127        incr _numCells
1128    }
1129    append out "DATASET UNSTRUCTURED_GRID\n"
1130    append out "POINTS $_numPoints double\n"
1131    foreach x [$xv range 0 end] y [$yv range 0 end] z [$zv range 0 end] {
1132        append out "$x $y $z\n"
1133    }
1134    set count [expr $_numCells * 6]
1135    append out "CELLS $_numCells $count\n"
1136    append out $data
1137    append out "CELL_TYPES $_numCells\n"
1138    append out $celltypes
1139    set _limits(x) [$xv limits]
1140    set _limits(y) [$yv limits]
1141    set _limits(z) [$zv limits]
1142
1143    set _vtkdata $out
1144    return 1
1145}
1146
1147itcl::body Rappture::Mesh::WriteHybridCells { path xv yv zv cells celltypes } {
1148    set _type "unstructured"
1149    set lines [split $cells \n]
1150    set numCellTypes [llength $celltypes]
1151    if { $numCellTypes == 1 } {
1152        set celltype [GetCellType $celltypes]
1153    }
1154
1155    set _numPoints [$xv length]
1156    set data {}
1157    set count 0
1158    set _numCells 0
1159    set celltypesout {}
1160    foreach line $lines {
1161        set length [llength $line]
1162        if { $length == 0 } {
1163            continue
1164        }
1165        if { $numCellTypes > 1 } {
1166            set celltype [GetCellType [lindex $celltypes $_numCells]]
1167        }
1168        set numIndices [GetNumIndices $celltype]
1169        if { $numIndices > 0 && $numIndices != $length } {
1170            puts stderr "WARNING: bad unstructured grid \"$path\": wrong \# of indices specified for celltype $celltype on line \"$line\""
1171            return 0
1172        } else {
1173            set numIndices $length
1174        }
1175        append data "$numIndices $line\n"
1176        append celltypesout "$celltype\n"
1177        incr count $length;         # Include the indices
1178        incr count;                 # and the number of indices
1179        incr _numCells
1180    }
1181    append out "DATASET UNSTRUCTURED_GRID\n"
1182    append out "POINTS $_numPoints double\n"
1183    foreach x [$xv range 0 end] y [$yv range 0 end] z [$zv range 0 end] {
1184        append out "$x $y $z\n"
1185    }
1186    append out "CELLS $_numCells $count\n"
1187    append out $data
1188    append out "CELL_TYPES $_numCells\n"
1189    append out $celltypesout
1190    set _limits(x) [$xv limits]
1191    set _limits(y) [$yv limits]
1192    if { $_dim < 3 } {
1193        set _limits(z) [list 0 0]
1194    } else {
1195        set _limits(z) [$zv limits]
1196    }
1197
1198    set _vtkdata $out
1199    return 1
1200}
1201
1202itcl::body Rappture::Mesh::ReadUnstructuredGrid { path } {
1203    set _type "unstructured"
1204
1205    if { ![GetDimension $path] } {
1206        return 0
1207    }
1208    # Step 1: Verify that there's only one cell tag of any kind.
1209    set numCells 0
1210    foreach type {
1211        cells
1212        hexahedrons
1213        lines
1214        polygons
1215        pyramids
1216        quads
1217        tetrahedrons
1218        triangles
1219        trianglestrips
1220        vertices
1221        wedges
1222    } {
1223        set data [$_xmlobj get $path.unstructured.$type]
1224        if { $data != "" } {
1225            incr numCells
1226        }
1227    }
1228    # The generic <cells> tag requires there be a <celltypes> tag too.
1229    set celltypes [$_xmlobj get $path.unstructured.celltypes]
1230    if { $numCells == 0 && $celltypes != "" } {
1231        puts stderr "WARNING: bad unstuctured grid \"$path\": no <cells> description found."
1232        return 0
1233    }
1234    if { $numCells > 1 } {
1235        puts stderr "WARNING: bad unstructured grid \"$path\": too many <cells>, <triangles>, <quads>... descriptions found: should be only one."
1236        return 0
1237    }
1238    foreach type {
1239        cells
1240        hexahedrons
1241        lines
1242        polygons
1243        pyramids
1244        quads
1245        tetrahedrons
1246        triangles
1247        trianglestrips
1248        vertices
1249        wedges
1250    } {
1251        set data [$_xmlobj get $path.unstructured.$type]
1252        if { $data != "" } {
1253            break
1254        }
1255    }
1256    # Step 2: Allow points to be specified as <points> or
1257    #         <xcoords>, <ycoords>, <zcoords>.  Split and convert into
1258    #         3 vectors, one for each coordinate.
1259    if { $_dim == 1 } {
1260        set xcoords [$_xmlobj get $path.unstructured.xcoords]
1261        set ycoords [$_xmlobj get $path.unstructured.ycoords]
1262        set zcoords [$_xmlobj get $path.unstructured.zcoords]
1263        set data    [$_xmlobj get $path.unstructured.points]
1264        if { $ycoords != "" } {
1265            put stderr "can't specify <ycoords> with a 1D mesh"
1266            return 0
1267        }
1268        if { $zcoords != "" } {
1269            put stderr "can't specify <zcoords> with a 1D mesh"
1270            return 0
1271        }
1272        if { $xcoords != "" } {
1273            set xv [blt::vector create \#auto]
1274            $xv set $xcoords
1275        } elseif { $data != "" } {
1276            Rappture::ReadPoints $data dim points
1277            if { $points == "" } {
1278                puts stderr "WARNING: bad unstructured grid \"$path\": no <points> found."
1279                return 0
1280            }
1281            if { $dim != 1 } {
1282                puts stderr "WARNING: bad unstructured grid \"$path\": \# of coordinates per point is \"$dim\": does not agree with dimension specified for mesh \"$_dim\""
1283                return 0
1284            }
1285            set xv [blt::vector create \#auto]
1286            $xv set $points
1287        } else {
1288            puts stderr "WARNING: bad unstructured grid \"$path\": no points specified."
1289            return 0
1290        }
1291        set yv [blt::vector create \#auto]
1292        set zv [blt::vector create \#auto]
1293        $yv seq 0 0 [$xv length];       # Make an all zeroes vector.
1294        $zv seq 0 0 [$xv length];       # Make an all zeroes vector.
1295    } elseif { $_dim == 2 } {
1296        set xcoords [$_xmlobj get $path.unstructured.xcoords]
1297        set ycoords [$_xmlobj get $path.unstructured.ycoords]
1298        set zcoords [$_xmlobj get $path.unstructured.zcoords]
1299        set data    [$_xmlobj get $path.unstructured.points]
1300        if { $zcoords != "" } {
1301            put stderr "can't specify <zcoords> with a 2D mesh"
1302            return 0
1303        }
1304        if { $xcoords != "" && $ycoords != "" } {
1305            set xv [blt::vector create \#auto]
1306            set yv [blt::vector create \#auto]
1307            $xv set $xcoords
1308            $yv set $ycoords
1309        } elseif { $data != "" } {
1310            Rappture::ReadPoints $data dim points
1311            if { $points == "" } {
1312                puts stderr "WARNING: bad unstructured grid \"$path\": no <points> found."
1313                return 0
1314            }
1315            if { $dim != 2 } {
1316                puts stderr "WARNING: bad unstructured grid \"$path\": \# of coordinates per point is \"$dim\": does not agree with dimension specified for mesh \"$_dim\""
1317                return 0
1318            }
1319            set all [blt::vector create \#auto]
1320            set xv [blt::vector create \#auto]
1321            set yv [blt::vector create \#auto]
1322            $all set $points
1323            $all split $xv $yv
1324            blt::vector destroy $all
1325        } else {
1326            puts stderr "WARNING: bad unstructured grid \"$path\": no points specified."
1327            return 0
1328        }
1329        set zv [blt::vector create \#auto]
1330        $zv seq 0 0 [$xv length];       # Make an all zeroes vector.
1331    } elseif { $_dim == 3 } {
1332        set xcoords [$_xmlobj get $path.unstructured.xcoords]
1333        set ycoords [$_xmlobj get $path.unstructured.ycoords]
1334        set zcoords [$_xmlobj get $path.unstructured.zcoords]
1335        set data    [$_xmlobj get $path.unstructured.points]
1336        if { $xcoords != "" && $ycoords != "" && $zcoords != "" } {
1337            set xv [blt::vector create \#auto]
1338            set yv [blt::vector create \#auto]
1339            set zv [blt::vector create \#auto]
1340            $xv set $xcoords
1341            $yv set $ycoords
1342            $zv set $zcoords
1343        } elseif { $data != "" }  {
1344            Rappture::ReadPoints $data dim points
1345            if { $points == "" } {
1346                puts stderr "WARNING: bad unstructured grid \"$path\": no <points> found."
1347                return 0
1348            }
1349            if { $dim != 3 } {
1350                puts stderr "WARNING: bad unstructured grid \"$path\": \# of coordinates per point is \"$dim\": does not agree with dimension specified for mesh \"$_dim\""
1351                return 0
1352            }
1353            set all [blt::vector create \#auto]
1354            set xv [blt::vector create \#auto]
1355            set yv [blt::vector create \#auto]
1356            set zv [blt::vector create \#auto]
1357            $all set $points
1358            $all split $xv $yv $zv
1359            blt::vector destroy $all
1360        } else {
1361            puts stderr "WARNING: bad unstructured grid \"$path\": no points specified."
1362            return 0
1363        }
1364    }
1365    set _numPoints [$xv length]
1366
1367    # Step 3: Write the points and cells as vtk data.
1368    if { $numCells == 0 } {
1369        set result [WritePointCloud $path $xv $yv $zv]
1370    } elseif { $type == "cells" } {
1371        set cells [$_xmlobj get $path.unstructured.cells]
1372        if { $cells == "" } {
1373            puts stderr "WARNING: bad unstructured grid \"$path\": no cells found."
1374            return 0
1375        }
1376        set result [WriteHybridCells $path $xv $yv $zv $cells $celltypes]
1377    } else {
1378        set cmd "Write[string totitle $type]"
1379        set cells [$_xmlobj get $path.unstructured.$type]
1380        if { $cells == "" } {
1381            puts stderr "WARNING: bad unstructured grid \"$path\": no cells found."
1382            return 0
1383        }
1384        set result [$cmd $path $xv $yv $zv $cells]
1385    }
1386    # Clean up.
1387    blt::vector destroy $xv $yv $zv
1388    return $result
1389}
1390
1391# ----------------------------------------------------------------------
1392# USAGE: ReadNodesElements <path>
1393#
1394# Used internally to build a mesh representation based on nodes and
1395# elements stored in the XML.
1396# ----------------------------------------------------------------------
1397itcl::body Rappture::Mesh::ReadNodesElements {path} {
1398    set _type "nodeselements"
1399    set count 0
1400
1401    # Scan for nodes.  Each node represents a vertex.
1402    set data {}
1403    foreach cname [$_xmlobj children -type node $path] {
1404        append data "[$_xmlobj get $path.$cname]\n"
1405    }
1406    Rappture::ReadPoints $data _dim points
1407    if { $_dim == 2 } {
1408        set all [blt::vector create \#auto]
1409        set xv [blt::vector create \#auto]
1410        set yv [blt::vector create \#auto]
1411        set zv [blt::vector create \#auto]
1412        $all set $points
1413        $all split $xv $yv
1414        set _numPoints [$xv length]
1415        set _limits(x) [$xv limits]
1416        set _limits(y) [$yv limits]
1417        set _limits(z) [list 0 0]
1418        # 2D Dataset. All Z coordinates are 0
1419        $zv seq 0.0 0.0 $_numPoints
1420        $all merge $xv $yv $zv
1421        set points [$all range 0 end]
1422        blt::vector destroy $all $xv $yv $zv
1423    } elseif { $_dim == 3 } {
1424        set all [blt::vector create \#auto]
1425        set xv [blt::vector create \#auto]
1426        set yv [blt::vector create \#auto]
1427        set zv [blt::vector create \#auto]
1428        $all set $points
1429        $all split $xv $yv $zv
1430        set _numPoints [$xv length]
1431        set _limits(x) [$xv limits]
1432        set _limits(y) [$yv limits]
1433        set _limits(z) [$zv limits]
1434        set points [$all range 0 end]
1435        blt::vector destroy $all $xv $yv $zv
1436    } else {
1437        error "bad dimension \"$_dim\" for nodes mesh"
1438    }
1439    array set node2celltype {
1440        3 5
1441        4 10
1442        8 12
1443        6 13
1444        5 14
1445    }
1446    set count 0
1447    set _numCells 0
1448    set celltypes {}
1449    set data {}
1450    # Next scan for elements.  Each element represents a cell.
1451    foreach cname [$_xmlobj children -type element $path] {
1452        set nodeList [$_mesh get $cname.nodes]
1453        set numNodes [llength $nodeList]
1454        if { ![info exists node2celltype($numNodes)] } {
1455            puts stderr "WARNING: bad nodes/elements mesh \$path\": unknown number of indices \"$_numNodes\": should be 3, 4, 5, 6, or 8"
1456            return 0
1457        }
1458        set celltype $node2celltype($numNodes)
1459        append celltypes "$celltype\n"
1460        if { $celltype == 12 } {
1461            # Formerly used voxels instead of hexahedrons. We're converting
1462            # it here to be backward compatible and still fault-tolerant of
1463            # non-axis aligned cells.
1464            set newList {}
1465            foreach i { 0 1 3 2 4 5 7 6 } {
1466                lappend newList [lindex $nodeList $i]
1467            }
1468            set nodeList $newList
1469        }
1470        append data "$numNodes $nodeList\n"
1471        incr _numCells
1472        incr count $numNodes
1473        incr count;                        # One extra for the VTK celltype id.
1474    }
1475
1476    append out "DATASET UNSTRUCTURED_GRID\n"
1477    append out "POINTS $_numPoints double\n"
1478    append out $points
1479    append out "\n"
1480    append out "CELLS $_numCells $count\n"
1481    append out $data
1482    append out "CELL_TYPES $_numCells\n"
1483    append out $celltypes
1484    set _vtkdata $out
1485    set _isValid 1
1486}
1487
1488itcl::body Rappture::Mesh::GetCellType { name } {
1489    array set name2type {
1490        "vertex"          1
1491        "polyvertex"      2
1492        "line"            3
1493        "polyline"        4
1494        "triangle"        5
1495        "trianglestrip"   6
1496        "polygon"         7
1497        "pixel"           8
1498        "quad"            9
1499        "tetrahedron"     10
1500        "voxel"           11
1501        "hexahedron"      12
1502        "wedge"           13
1503        "pyramid"         14
1504        "pentagonalprism" 15
1505        "hexagonalprism"  16
1506    }
1507    if { [info exists name2type($name)] } {
1508        return $name2type($name)
1509    }
1510    if { [string is int $name] == 1 && $name >= 1 && $name <= 16 } {
1511        return $name
1512    }
1513    error "invalid celltype \"$name\""
1514}
1515
1516itcl::body Rappture::Mesh::GetNumIndices { type } {
1517    array set type2indices {
1518        1       1
1519        2       0
1520        3       2
1521        4       0
1522        5       3
1523        6       0
1524        7       0
1525        8       4
1526        9       4
1527        10      4
1528        11      8
1529        12      8
1530        13      6
1531        14      5
1532        15      10
1533        16      12
1534    }
1535    if { [info exists type2indices($type)] } {
1536        return $type2indices($type)
1537    }
1538    error "invalid celltype \"$type\""
1539}
Note: See TracBrowser for help on using the repository browser.