source: trunk/gui/scripts/field.tcl @ 3365

Last change on this file since 3365 was 3365, checked in by gah, 12 years ago

fix for old unirect2d blob method, field limits for unirect2d
File size: 42.6 KB
Line 
1# -*- mode: tcl; indent-tabs-mode: nil -*-
2
3# ----------------------------------------------------------------------
4#  COMPONENT: field - extracts data from an XML description of a field
5#
6#  This object represents one field in an XML description of a device.
7#  It simplifies the process of extracting data vectors that represent
8#  the field.
9# ======================================================================
10#  AUTHOR:  Michael McLennan, Purdue University
11#  Copyright (c) 2004-2012  HUBzero Foundation, LLC
12#
13#  See the file "license.terms" for information on usage and
14#  redistribution of this file, and for a DISCLAIMER OF ALL WARRANTIES.
15# ======================================================================
16package require Itcl
17package require BLT
18
19namespace eval Rappture {
20    # forward declaration
21}
22
23#
24# Possible field dataset types:
25#
26# 2D Datasets
27#       vtk             (range of z-axis is zero).
28#       unirect2d       (deprecated except where extents > 1)
29#       cloud           (x,y point coordinates) (deprecated)
30#       mesh
31# 3D Datasets
32#       vtk
33#       unirect3d
34#       cloud           (x,y,z coordinates) (deprecated)
35#       mesh
36#       dx              (FIXME: make dx-to-vtk converter work)
37#       ucd avs
38#
39# Viewers:
40#       Format     Dim  Description                     Viewer          Server
41#       vtk         2   vtk file data.                  contour         vtkvis
42#       vtk         3   vtk file data.                  isosurface      vtkvis
43#       mesh        2   points-on-mesh                  heightmap       vtkvis
44#       mesh        3   points-on-mesh                  isosurface      vtkvis
45#       dx          3   DX                              volume          nanovis
46#       unirect2d   2   unirect3d + extents > 1 flow    flow            nanovis
47#       unirect3d   3   unirect2d + extents > 1 flow    flow            nanovis
48#       
49# The goal should be any 3D dataset can view a isosurface, volume,
50# streamlines, or flow (if extents > 1).  Any 2D dataset can view a
51# contour, heightmap, streamlines, or flow (if extents > 1). 
52#
53#
54itcl::class Rappture::Field {
55    private variable _dim       0;      # Dimension of the mesh
56    private variable _xmlobj ""  ;      # ref to XML obj with field data
57    private variable _limits     ;      # maps box name => {z0 z1} limits
58    private variable _field ""
59    private variable _fieldNames ""   ; # list of field names.
60    private variable _fieldUnits ""   ; # list of units of each field name.
61    private variable _fieldLabels ""  ; # list of labels of each field name.
62    private variable _viewer            ""
63    private variable _hints
64
65    constructor {xmlobj path} {
66        # defined below
67    }
68    destructor {
69        # defined below
70    }
71    public method components {args}
72    public method mesh {{cname -overall}}
73    public method values {{cname -overall}}
74    public method blob {{cname -overall}}
75    public method limits {axis}
76    public method controls {option args}
77    public method hints {{key ""}}
78    public method style { cname }
79    public method isunirect2d {}
80    public method isunirect3d {}
81    public method extents {{cname -overall}}
82    public method flowhints { cname }
83    public method type {}
84    public method viewer {} {
85        return $_viewer
86    }
87    public method vtkdata {cname}
88   
89    protected method Build {}
90    protected method _getValue {expr}
91
92    private variable _path "";          # Path of this object in the XML
93    private variable _units ""   ;      # system of units for this field
94    private variable _zmax 0     ;# length of the device
95
96    private variable _comp2dims  ;# maps component name => dimensionality
97    private variable _comp2xy    ;# maps component name => x,y vectors
98    private variable _comp2vtk   ;# maps component name => vtk file data
99    private variable _comp2dx    ;# maps component name => OpenDX data
100    private variable _comp2unirect2d ;# maps component name => unirect2d obj
101    private variable _comp2unirect3d ;# maps component name => unirect3d obj
102    private variable _comp2style ;# maps component name => style settings
103    private variable _comp2cntls ;# maps component name => x,y control points
104    private variable _comp2extents
105    private variable _comp2limits;      #  Array of limits per component
106    private variable _type ""
107    private variable _comp2flowhints
108    private variable _comp2mesh
109    private common _counter 0    ;# counter for unique vector names
110
111    private method BuildPointsOnMesh { cname }
112    private method ConvertToVtkData { cname }
113    private method ReadVtkDataSet { cname contents }
114    private method AvsToVtk { cname contents }
115    private variable _values ""
116}
117
118# ----------------------------------------------------------------------
119# CONSTRUCTOR
120# ----------------------------------------------------------------------
121itcl::body Rappture::Field::constructor {xmlobj path} {
122    if {![Rappture::library isvalid $xmlobj]} {
123        error "bad value \"$xmlobj\": should be Rappture::library"
124    }
125    set _xmlobj $xmlobj
126    set _path $path
127    set _field [$xmlobj element -as object $path]
128    set _units [$_field get units]
129
130    set xunits [$xmlobj get units]
131    if {"" == $xunits || "arbitrary" == $xunits} {
132        set xunits "um"
133    }
134
135    # determine the overall size of the device
136    set z0 [set z1 0]
137    foreach elem [$_xmlobj children components] {
138        switch -glob -- $elem {
139            box* {
140                if {![regexp {[0-9]$} $elem]} {
141                    set elem "${elem}0"
142                }
143                set z0 [$_xmlobj get components.$elem.corner0]
144                set z0 [Rappture::Units::convert $z0 \
145                    -context $xunits -to $xunits -units off]
146
147                set z1 [$_xmlobj get components.$elem.corner1]
148                set z1 [Rappture::Units::convert $z1 \
149                    -context $xunits -to $xunits -units off]
150
151                set _limits($elem) [list $z0 $z1]
152            }
153        }
154    }
155    set _zmax $z1
156
157    # build up vectors for various components of the field
158    Build
159}
160
161# ----------------------------------------------------------------------
162# DESTRUCTOR
163# ----------------------------------------------------------------------
164itcl::body Rappture::Field::destructor {} {
165    itcl::delete object $_field
166    # don't destroy the _xmlobj! we don't own it!
167
168    foreach name [array names _comp2xy] {
169        eval blt::vector destroy $_comp2xy($name)
170    }
171    foreach name [array names _comp2unirect2d] {
172        itcl::delete object $_comp2unirect2d($name)
173    }
174    foreach name [array names _comp2unirect3d] {
175        itcl::delete object $_comp2unirect3d($name)
176    }
177    foreach name [array names _comp2flowhints] {
178        itcl::delete object $_comp2flowhints($name)
179    }
180    foreach name [array names _comp2mesh] {
181        # Data is in the form of a mesh and a vector.
182        foreach { mesh vector } $_comp2mesh($name) break
183        # Release the mesh (may be shared)
184        set class [$mesh info class]
185        ${class}::release $mesh
186        # Destroy the vector
187        blt::vector destroy $vector
188    }
189}
190
191# ----------------------------------------------------------------------
192# USAGE: components ?-name|-dimensions|-style? ?<pattern>?
193#
194# Returns a list of names or types for the various components of
195# this field.  If the optional glob-style <pattern> is specified,
196# then it returns only the components with names matching the pattern.
197# ----------------------------------------------------------------------
198itcl::body Rappture::Field::components {args} {
199    Rappture::getopts args params {
200        flag what -name default
201        flag what -dimensions
202        flag what -style
203        flag what -particles
204        flag what -flow
205        flag what -box
206    }
207
208    set pattern *
209    if {[llength $args] > 0} {
210        set pattern [lindex $args 0]
211        set args [lrange $args 1 end]
212    }
213    if {[llength $args] > 0} {
214        error "wrong # args: should be \"components ?switches? ?pattern?\""
215    }
216
217    set rlist ""
218
219    # BE CAREFUL: return component names in proper order
220    foreach cname [$_field children -type component] {
221        if {[info exists _comp2dims($cname)]
222              && [string match $pattern $cname]} {
223
224            switch -- $params(what) {
225                -name { lappend rlist $cname }
226                -dimensions { lappend rlist $_comp2dims($cname) }
227                -style { lappend rlist $_comp2style($cname) }
228            }
229        }
230    }
231    return $rlist
232}
233
234# ----------------------------------------------------------------------
235# USAGE: mesh ?<name>?
236#
237# Returns a list {xvec yvec} for the specified field component <name>.
238# If the name is not specified, then it returns the vectors for the
239# overall field (sum of all components).
240# ----------------------------------------------------------------------
241itcl::body Rappture::Field::mesh {{cname -overall}} {
242    if {$cname == "-overall" || $cname == "component0"} {
243        set cname [lindex [components -name] 0]
244    }
245    if {[info exists _comp2xy($cname)]} {
246        return [lindex $_comp2xy($cname) 0]  ;# return xv
247    }
248    if { [info exists _comp2vtk($cname)] } {
249        # FIXME: extract mesh from VTK file data.
250        error "method \"mesh\" is not implemented for VTK file data"
251    }
252    if {[info exists _comp2dx($cname)]} {
253        return ""  ;# no mesh -- it's embedded in the value data
254    }
255    if {[info exists _comp2mesh($cname)]} {
256        return ""  ;# no mesh -- it's embedded in the value data
257    }
258    if {[info exists _comp2unirect2d($cname)]} {
259        set mobj [lindex $_comp2unirect2d($cname) 0]
260        return [$mobj mesh]
261    }
262    if {[info exists _comp2unirect3d($cname)]} {
263        set mobj [lindex $_comp2unirect3d($cname) 0]
264        return [$mobj mesh]
265    }
266    error "bad option \"$cname\": should be [join [lsort [array names _comp2dims]] {, }]"
267}
268
269# ----------------------------------------------------------------------
270# USAGE: values ?<name>?
271#
272# Returns a list {xvec yvec} for the specified field component <name>.
273# If the name is not specified, then it returns the vectors for the
274# overall field (sum of all components).
275# ----------------------------------------------------------------------
276itcl::body Rappture::Field::values {{cname -overall}} {
277    if {$cname == "component0"} {
278        set cname "component"
279    }
280    if {[info exists _comp2xy($cname)]} {
281        return [lindex $_comp2xy($cname) 1]  ;# return yv
282    }
283    # VTK file data
284    if { [info exists _comp2vtk($cname)] } {
285        # FIXME: extract the values from the VTK file data
286        error "method \"values\" is not implemented for vtk file data"
287    }
288    # Points-on-mesh
289    if { [info exists _comp2mesh($cname)] } {
290        set vector [lindex $_comp2mesh($cname) 1]
291        return [$vector range 0 end]
292    }
293    if {[info exists _comp2dx($cname)]} {
294        return $_comp2dx($cname)  ;# return gzipped, base64-encoded DX data
295    }
296    if {[info exists _comp2unirect2d($cname)]} {
297        return [$_comp2unirect2d($cname) values]
298    }
299    if {[info exists _comp2unirect3d($cname)]} {
300        return [$_comp2unirect3d($cname) blob]
301    }
302    error "bad option \"$cname\": should be [join [lsort [array names _comp2dims]] {, }]"
303}
304
305# ----------------------------------------------------------------------
306# USAGE: blob ?<name>?
307#
308# Returns a string representing the blob of data for the mesh and values.
309# ----------------------------------------------------------------------
310itcl::body Rappture::Field::blob {{cname -overall}} {
311    if {$cname == "component0"} {
312        set cname "component"
313    }
314    if {[info exists _comp2xy($cname)]} {
315        return ""
316    }
317    if { [info exists _comp2vtk($cname)] } {
318        error "blob not implemented for VTK file data"
319    }
320    if {[info exists _comp2dx($cname)]} {
321        return $_comp2dx($cname)  ;# return gzipped, base64-encoded DX data
322    }
323    if {[info exists _comp2unirect2d($cname)]} {
324        set blob [$_comp2unirect2d($cname) blob]
325        lappend blob "values" $_values
326        return $blob
327    }
328    if {[info exists _comp2unirect3d($cname)]} {
329        return [$_comp2unirect3d($cname) blob]
330    }
331    error "bad option \"$cname\": should be [join [lsort [array names _comp2dims]] {, }]"
332}
333
334# ----------------------------------------------------------------------
335# USAGE: limits <axis>
336#
337# Returns a list {min max} representing the limits for the specified
338# axis.
339# ----------------------------------------------------------------------
340itcl::body Rappture::Field::limits {which} {
341    set min ""
342    set max ""
343    blt::vector tmp zero
344
345    foreach cname [array names _comp2dims] {
346        switch -- $_comp2dims($cname) {
347            1D {
348                switch -- $which {
349                    x - xlin {
350                        set pos 0; set log 0; set axis x
351                    }
352                    xlog {
353                        set pos 0; set log 1; set axis x
354                    }
355                    y - ylin - v - vlin {
356                        set pos 1; set log 0; set axis y
357                    }
358                    ylog - vlog {
359                        set pos 1; set log 1; set axis y
360                    }
361                    default {
362                        error "bad option \"$which\": should be x, xlin, xlog, y, ylin, ylog, v, vlin, vlog"
363                    }
364                }
365
366                set vname [lindex $_comp2xy($cname) $pos]
367                $vname variable vec
368
369                if {$log} {
370                    # on a log scale, use abs value and ignore 0's
371                    $vname dup tmp
372                    $vname dup zero
373                    zero expr {tmp == 0}            ;# find the 0's
374                    tmp expr {abs(tmp)}             ;# get the abs value
375                    tmp expr {tmp + zero*max(tmp)}  ;# replace 0's with abs max
376                    set axisMin [blt::vector expr min(tmp)]
377                    set axisMax [blt::vector expr max(tmp)]
378                } else {
379                    set axisMin $vec(min)
380                    set axisMax $vec(max)
381                }
382
383                if {"" == $min} {
384                    set min $axisMin
385                } elseif {$axisMin < $min} {
386                    set min $axisMin
387                }
388                if {"" == $max} {
389                    set max $axisMax
390                } elseif {$axisMax > $max} {
391                    set max $axisMax
392                }
393            }
394            2D - 3D {
395                if {[info exists _comp2unirect3d($cname)]} {
396                    set limits [$_comp2unirect3d($cname) limits $which]
397                    foreach {axisMin axisMax} $limits break
398                    set axis v
399                } elseif {[info exists _comp2limits($cname)]} {
400                    array set limits $_comp2limits($cname)
401                    switch -- $which {
402                        x - xlin - xlog {
403                            set axis x
404                            foreach {axisMin axisMax} $limits(x) break
405                        }
406                        y - ylin - ylog {
407                            set axis y
408                            foreach {axisMin axisMax} $limits(y) break
409                        }
410                        z - zlin - zlog {
411                            set axis y
412                            foreach {axisMin axisMax} $limits(z) break
413                        }
414                        v - vlin - vlog {
415                            set axis v
416                            foreach {axisMin axisMax} $limits(v) break
417                        }
418                        default {
419                            if { ![info exists limits($which)] } {
420                                error "limits: unknown axis \"$which\""
421                            }
422                            set axis v
423                            foreach {axisMin axisMax} $limits($which) break
424                        }
425                    }
426                } else {
427                    set axisMin 0  ;# HACK ALERT! must be OpenDX data
428                    set axisMax 1
429                    set axis v
430                }
431            }
432        }
433        if { "" == $min || $axisMin < $min } {
434            set min $axisMin
435        }
436        if { "" == $max || $axisMax > $max } {
437            set max $axisMax
438        }
439    }
440    blt::vector destroy tmp zero
441    set val [$_field get "${axis}axis.min"]
442    if {"" != $val && "" != $min} {
443        if {$val > $min} {
444            # tool specified this min -- don't go any lower
445            set min $val
446        }
447    }
448    set val [$_field get "${axis}axis.max"]
449    if {"" != $val && "" != $max} {
450        if {$val < $max} {
451            # tool specified this max -- don't go any higher
452            set max $val
453        }
454    }
455    return [list $min $max]
456}
457
458# ----------------------------------------------------------------------
459# USAGE: controls get ?<name>?
460# USAGE: controls validate <path> <value>
461# USAGE: controls put <path> <value>
462#
463# Returns a list {path1 x1 y1 val1  path2 x2 y2 val2 ...} representing
464# control points for the specified field component <name>.
465# ----------------------------------------------------------------------
466itcl::body Rappture::Field::controls {option args} {
467    switch -- $option {
468        get {
469            set cname [lindex $args 0]
470            if {[info exists _comp2cntls($cname)]} {
471                return $_comp2cntls($cname)
472            }
473            return ""
474        }
475        validate {
476            set path [lindex $args 0]
477            set value [lindex $args 1]
478            set units [$_xmlobj get $path.units]
479
480            if {"" != $units} {
481                set nv [Rappture::Units::convert \
482                    $value -context $units -to $units -units off]
483            } else {
484                set nv $value
485            }
486            if {![string is double $nv]
487                  || [regexp -nocase {^(inf|nan)$} $nv]} {
488                error "Value out of range"
489            }
490
491            set rawmin [$_xmlobj get $path.min]
492            if {"" != $rawmin} {
493                set minv $rawmin
494                if {"" != $units} {
495                    set minv [Rappture::Units::convert \
496                        $minv -context $units -to $units -units off]
497                    set nv [Rappture::Units::convert \
498                        $value -context $units -to $units -units off]
499                }
500                # fix for the case when the user tries to
501                # compare values like minv=-500 nv=-0600
502                set nv [format "%g" $nv]
503                set minv [format "%g" $minv]
504
505                if {$nv < $minv} {
506                    error "Minimum value allowed here is $rawmin"
507                }
508            }
509
510            set rawmax [$_xmlobj get $path.max]
511            if {"" != $rawmax} {
512                set maxv $rawmax
513                if {"" != $units} {
514                    set maxv [Rappture::Units::convert \
515                        $maxv -context $units -to $units -units off]
516                    set nv [Rappture::Units::convert \
517                        $value -context $units -to $units -units off]
518                }
519                # fix for the case when the user tries to
520                # compare values like maxv=-500 nv=-0600
521                set nv [format "%g" $nv]
522                set maxv [format "%g" $maxv]
523
524                if {$nv > $maxv} {
525                    error "Maximum value allowed here is $rawmax"
526                }
527            }
528
529            return "ok"
530        }
531        put {
532            set path [lindex $args 0]
533            set value [lindex $args 1]
534            $_xmlobj put $path.current $value
535            Build
536        }
537        default {
538            error "bad option \"$option\": should be get or put"
539        }
540    }
541}
542
543# ----------------------------------------------------------------------
544# USAGE: hints ?<keyword>?
545#
546# Returns a list of key/value pairs for various hints about plotting
547# this field.  If a particular <keyword> is specified, then it returns
548# the hint for that <keyword>, if it exists.
549# ----------------------------------------------------------------------
550itcl::body Rappture::Field::hints {{keyword ""}} {
551    if { ![info exists _hints] } {
552        set _hints(fieldnames)  $_fieldNames
553        set _hints(fieldunits)  $_fieldUnits
554        set _hints(fieldlabels) $_fieldLabels
555        foreach {key path} {
556            camera          camera.position
557            color           about.color
558            default         about.default
559            group           about.group
560            label           about.label
561            fieldnames      about.fieldnames
562            fieldunits      about.fieldunits
563            fieldlabels     about.fieldlabels
564            scale           about.scale
565            seeds           about.seeds
566            style           about.style
567            type            about.type
568            xlabel          about.xaxis.label
569            ylabel          about.yaxis.label
570            zlabel          about.zaxis.label
571            xunits          about.xaxis.units
572            yunits          about.yaxis.units
573            zunits          about.zaxis.units
574            units           units
575            updir           updir
576            vectors         about.vectors
577        } {
578            set str [$_field get $path]
579            if { "" != $str } {
580                set _hints($key) $str
581            }
582        }
583        foreach {key path} {
584            toolid          tool.id
585            toolname        tool.name
586            toolcommand     tool.execute
587            tooltitle       tool.title
588            toolrevision    tool.version.application.revision
589        } {
590            set str [$_xmlobj get $path]
591            if { "" != $str } {
592                set _hints($key) $str
593            }
594        }
595        # Set toolip and path hints
596        set _hints(path) $_path
597        if { [info exists _hints(group)] && [info exists _hints(label)] } {
598            # pop-up help for each curve
599            set _hints(tooltip) $_hints(label)
600        }
601    }
602    if { $keyword != "" } {
603        if {[info exists _hints($keyword)]} {
604            return $_hints($keyword)
605        }
606        return ""
607    }
608    return [array get _hints]
609}
610
611# ----------------------------------------------------------------------
612# USAGE: Build
613#
614# Used internally to build up the vector representation for the
615# field when the object is first constructed, or whenever the field
616# data changes.  Discards any existing vectors and builds everything
617# from scratch.
618# ----------------------------------------------------------------------
619itcl::body Rappture::Field::Build {} {
620
621    # Discard any existing data
622    foreach name [array names _comp2xy] {
623        eval blt::vector destroy $_comp2xy($name)
624    }
625    array unset _comp2vtk
626    foreach name [array names _comp2unirect2d] {
627        eval itcl::delete object $_comp2unirect2d($name)
628    }
629    foreach name [array names _comp2unirect3d] {
630        eval itcl::delete object $_comp2unirect3d($name)
631    }
632    catch {unset _comp2xy}
633    catch {unset _comp2dx}
634    catch {unset _comp2dims}
635    catch {unset _comp2style}
636    array unset _comp2unirect2d
637    array unset _comp2unirect3d
638    array unset _comp2extents
639    array unset _dataobj2type
640    #
641    # Scan through the components of the field and create
642    # vectors for each part.
643    #
644    foreach cname [$_field children -type component] {
645        set type ""
646        if { ([$_field element $cname.constant] != "" &&
647              [$_field element $cname.domain] != "") ||
648              [$_field element $cname.xy] != "" } {
649            set type "1D"
650        } elseif { [$_field element $cname.mesh] != "" &&
651                   [$_field element $cname.values] != ""} {
652            set type "points-on-mesh"
653        } elseif { [$_field element $cname.vtk] != ""} {
654            set viewer [$_field get "about.view"]
655            set type "vtk"
656            if { $viewer != "" } {
657                set _viewer $viewer
658            }
659        } elseif {[$_field element $cname.opendx] != ""} {
660            global env
661            if { [info exists env(VTKVOLUME)] } {
662                set type "vtkvolume"
663            } else {
664                set type "dx"
665            }
666        } elseif {[$_field element $cname.dx] != ""} {
667            global env
668            if { [info exists env(VTKVOLUME)] } {
669                set type "vtkvolume"
670            } else {
671                set type "dx"
672            }
673        }
674        set _comp2style($cname) ""
675       
676        # Save the extents of the component
677        if { [$_field element $cname.extents] != "" } {
678            set extents [$_field get $cname.extents]
679        } else {
680            set extents 1
681        }
682        set _comp2extents($cname) $extents
683        set _type $type
684        if {$type == "1D"} {
685            #
686            # 1D data can be represented as 2 BLT vectors,
687            # one for x and the other for y.
688            #
689            set xv ""
690            set yv ""
691
692            set val [$_field get $cname.constant]
693            if {$val != ""} {
694                set domain [$_field get $cname.domain]
695                if {$domain == "" || ![info exists _limits($domain)]} {
696                    set z0 0
697                    set z1 $_zmax
698                } else {
699                    foreach {z0 z1} $_limits($domain) { break }
700                }
701                set xv [blt::vector create x$_counter]
702                $xv append $z0 $z1
703
704                foreach {val pcomp} [_getValue $val] break
705                set yv [blt::vector create y$_counter]
706                $yv append $val $val
707
708                if {$pcomp != ""} {
709                    set zm [expr {0.5*($z0+$z1)}]
710                    set _comp2cntls($cname) \
711                        [list $pcomp $zm $val "$val$_units"]
712                }
713            } else {
714                set xydata [$_field get $cname.xy]
715                if {"" != $xydata} {
716                    set xv [blt::vector create x$_counter]
717                    set yv [blt::vector create y$_counter]
718                    set tmp [blt::vector create \#auto]
719                    $tmp set $xydata
720                    $tmp split $xv $yv
721                    blt::vector destroy $tmp
722                }
723            }
724
725            if {$xv != "" && $yv != ""} {
726                # sort x-coords in increasing order
727                $xv sort $yv
728                set _comp2dims($cname) "1D"
729                set _comp2xy($cname) [list $xv $yv]
730                incr _counter
731            }
732        } elseif {$type == "points-on-mesh"} {
733            BuildPointsOnMesh $cname
734        } elseif {$type == "vtk"} {
735            set vtkdata [$_field get $cname.vtk]
736            ReadVtkDataSet $cname $vtkdata
737            set _comp2vtk($cname) $vtkdata
738            set _comp2style($cname) [$_field get $cname.style]
739            incr _counter
740        } elseif {$type == "dx" } {
741            #
742            # HACK ALERT!  Extract gzipped, base64-encoded OpenDX
743            # data.  Assume that it's 3D.  Pass it straight
744            # off to the NanoVis visualizer.
745            #
746            set _viewer "nanovis"
747            set _comp2dims($cname) "3D"
748            set _comp2dx($cname)  [$_field get -decode no $cname.dx]
749            if 1 {
750            set data  [$_field get -decode yes $cname.dx]
751            set file "/tmp/junk.dx"
752            set f [open $file "w"]
753            puts $f $data
754            close $f
755            if { [string match "<ODX>*" $data] } {
756                set data [string range $data 5 end]
757                set _comp2dx($cname) \
758                        [Rappture::encoding::encode -as zb64 $data]
759            }
760            }
761            set _comp2style($cname) [$_field get $cname.style]
762            if {[$_field element $cname.flow] != ""} {
763                set _comp2flowhints($cname) \
764                    [Rappture::FlowHints ::\#auto $_field $cname $_units]
765            }
766            incr _counter
767        } elseif {$type == "opendx"} {
768            #
769            # HACK ALERT!  Extract gzipped, base64-encoded OpenDX
770            # data.  Assume that it's 3D.  Pass it straight
771            # off to the NanoVis visualizer.
772            #
773            set _viewer "nanovis"
774            set _comp2dims($cname) "3D"
775            set data [$_field get -decode yes $cname.opendx]
776            set data "<ODX>$data"
777            set data [Rappture::encoding::encode -as zb64 $data]
778            set _comp2dx($cname) $data
779            set _comp2style($cname) [$_field get $cname.style]
780            if {[$_field element $cname.flow] != ""} {
781                set _comp2flowhints($cname) \
782                    [Rapture::FlowHints ::\#auto $_field $cname $_units]
783            }
784            incr _counter
785        } elseif {[$_field element $cname.ucd] != ""} {
786            set _viewer "isosurface"
787            set _comp2dims($cname) "3D"
788            set contents [$_field get $cname.ucd]
789            set vtkdata [AvsToVtk $cname $contents]
790            ReadVtkDataSet $cname $vtkdata
791            set _comp2vtk($cname) $vtkdata
792            set _comp2style($cname) [$_field get $cname.style]
793            incr _counter
794        }
795    }
796    # Sanity check.  Verify that all components of the field have the same
797    # dimension.
798    set dim ""
799    foreach cname [array names _comp2dims] {
800        if { $dim == "" } {
801            set dim $_comp2dims($cname)
802            continue
803        }
804        if { $dim != $_comp2dims($cname) } {
805            error "field can't have components of different dimensions: [join [array get _comp2dims] ,]"
806        }
807    }
808}
809
810# ----------------------------------------------------------------------
811# USAGE: _getValue <expr>
812#
813# Used internally to get the value for an expression <expr>.  Returns
814# a list of the form {val parameterPath}, where val is the numeric
815# value of the expression, and parameterPath is the XML path to the
816# parameter representing the value, or "" if the <expr> does not
817# depend on any parameters.
818# ----------------------------------------------------------------------
819itcl::body Rappture::Field::_getValue {expr} {
820    #
821    # First, look for the expression among the <parameter>'s
822    # associated with the device.
823    #
824    set found 0
825    foreach pcomp [$_xmlobj children parameters] {
826        set id [$_xmlobj element -as id parameters.$pcomp]
827        if {[string equal $id $expr]} {
828            set val [$_xmlobj get parameters.$pcomp.current]
829            if {"" == $val} {
830                set val [$_xmlobj get parameters.$pcomp.default]
831            }
832            if {"" != $val} {
833                set expr $val
834                set found 1
835                break
836            }
837        }
838    }
839    if {$found} {
840        set pcomp "parameters.$pcomp"
841    } else {
842        set pcomp ""
843    }
844
845    if {$_units != ""} {
846        set expr [Rappture::Units::convert $expr \
847            -context $_units -to $_units -units off]
848    }
849
850    return [list $expr $pcomp]
851}
852
853#
854# isunirect2d  --
855#
856# Returns if the field is a unirect2d object. 
857#
858itcl::body Rappture::Field::isunirect2d { } {
859    return [expr [array size _comp2unirect2d] > 0]
860}
861
862#
863# isunirect3d  --
864#
865# Returns if the field is a unirect3d object. 
866#
867itcl::body Rappture::Field::isunirect3d { } {
868    return [expr [array size _comp2unirect3d] > 0]
869}
870
871#
872# flowhints  --
873#
874# Returns the hints associated with a flow vector field. 
875#
876itcl::body Rappture::Field::flowhints { cname } {
877    if { [info exists _comp2flowhints($cname)] } {
878        return $_comp2flowhints($cname)
879    }
880    return ""
881}
882
883#
884# style  --
885#
886# Returns the style associated with a component of the field. 
887#
888itcl::body Rappture::Field::style { cname } {
889    if { [info exists _comp2style($cname)] } {
890        return $_comp2style($cname)
891    }
892    return ""
893}
894
895#
896# type  --
897#
898# Returns the style associated with a component of the field. 
899#
900itcl::body Rappture::Field::type {} {
901    return $_type
902}
903
904#
905# extents --
906#
907# Returns if the field is a unirect2d object. 
908#
909itcl::body Rappture::Field::extents {{cname -overall}} {
910    if {$cname == "-overall" } {
911        set max 0
912        foreach cname [$_field children -type component] {
913            if { ![info exists _comp2unirect3d($cname)] &&
914                 ![info exists _comp2extents($cname)] } {
915                continue
916            }
917            set value $_comp2extents($cname)
918            if { $max < $value } {
919                set max $value
920            }
921        }
922        return $max
923    }
924    if { $cname == "component0"} {
925        set cname [lindex [components -name] 0]
926    }
927    return $_comp2extents($cname)
928}
929
930itcl::body Rappture::Field::ConvertToVtkData { cname } {
931    set ds ""
932    switch -- [typeof $cname] {
933        "unirect2d" {
934            foreach { x1 x2 xN y1 y2 yN } [$dataobj mesh $cname] break
935            set spacingX [expr {double($x2 - $x1)/double($xN - 1)}]
936            set spacingY [expr {double($y2 - $y1)/double($yN - 1)}]
937           
938            set ds [vtkImageData $this-grdataTemp]
939            $ds SetDimensions $xN $yN 1
940            $ds SetOrigin $x1 $y1 0
941            $ds SetSpacing $spacingX $spacingY 0
942            set arr [vtkDoubleArray $this-arrTemp]
943            foreach {val} [$dataobj values $cname] {
944                $arr InsertNextValue $val
945            }
946            [$ds GetPointData] SetScalars $arr
947        }
948        "unirect3d" {
949            foreach { x1 x2 xN y1 y2 yN z1 z2 zN } [$dataobj mesh $cname] break
950            set spacingX [expr {double($x2 - $x1)/double($xN - 1)}]
951            set spacingY [expr {double($y2 - $y1)/double($yN - 1)}]
952            set spacingZ [expr {double($z2 - $z1)/double($zN - 1)}]
953           
954            set ds [vtkImageData $this-grdataTemp]
955            $ds SetDimensions $xN $yN $zN
956            $ds SetOrigin $x1 $y1 $z1
957            $ds SetSpacing $spacingX $spacingY $spacingZ
958            set arr [vtkDoubleArray $this-arrTemp]
959            foreach {val} [$dataobj values $cname] {
960                $arr InsertNextValue $val
961            }
962            [$ds GetPointData] SetScalars $val
963        }
964        "contour" {
965            return [$dataobj blob $cname]
966        }
967        "dx" {
968            return [Rappture::ConvertDxToVtk $_comp2dx($cname)]
969        }
970        default {
971            set mesh [$dataobj mesh $cname]
972            switch -- [$mesh GetClassName] {
973                vtkPoints {
974                    # handle cloud of points
975                    set ds [vtkPolyData $this-polydataTemp]
976                    $ds SetPoints $mesh
977                    [$ds GetPointData] SetScalars [$dataobj values $cname]
978                }
979                vtkPolyData {
980                    set ds [vtkPolyData $this-polydataTemp]
981                    $ds ShallowCopy $mesh
982                    [$ds GetPointData] SetScalars [$dataobj values $cname]
983                }
984                vtkUnstructuredGrid {
985                    # handle 3D grid with connectivity
986                    set ds [vtkUnstructuredGrid $this-grdataTemp]
987                    $ds ShallowCopy $mesh
988                    [$ds GetPointData] SetScalars [$dataobj values $cname]
989                }
990                vtkRectilinearGrid {
991                    # handle 3D grid with connectivity
992                    set ds [vtkRectilinearGrid $this-grdataTemp]
993                    $ds ShallowCopy $mesh
994                    [$ds GetPointData] SetScalars [$dataobj values $cname]
995                }
996                default {
997                    error "don't know how to handle [$mesh GetClassName] data"
998                }
999            }
1000        }
1001    }
1002
1003    if {"" != $ds} {
1004        set writer [vtkDataSetWriter $this-dsWriterTmp]
1005        $writer SetInput $ds
1006        $writer SetFileTypeToASCII
1007        $writer WriteToOutputStringOn
1008        $writer Write
1009        set out [$writer GetOutputString]
1010        $ds Delete
1011        $writer Delete
1012    } else {
1013        set out ""
1014        error "No DataSet to write"
1015    }
1016
1017    append out "\n"
1018    return $out
1019}
1020
1021itcl::body Rappture::Field::ReadVtkDataSet { cname contents } {
1022    package require vtk
1023
1024    set reader $this-datasetreader
1025    vtkDataSetReader $reader
1026
1027    # Write the contents to a file just in case it's binary.
1028    set tmpfile file[pid].vtk
1029    set f [open "$tmpfile" "w"]
1030    fconfigure $f -translation binary -encoding binary
1031    puts $f $contents
1032    close $f
1033    $reader SetFileName $tmpfile
1034    $reader ReadAllScalarsOn
1035    $reader ReadAllVectorsOn
1036    $reader ReadAllFieldsOn
1037    $reader Update
1038    set dataset [$reader GetOutput]
1039    set limits {}
1040    foreach {xmin xmax ymin ymax zmin zmax} [$dataset GetBounds] break
1041    # Figure out the dimension of the mesh from the bounds.
1042    set _dim 0
1043    if { $xmax > $xmin } {
1044        incr _dim
1045    }
1046    if { $ymax > $ymin } {
1047        incr _dim
1048    }
1049    if { $zmax > $zmin } {
1050        incr _dim
1051    }
1052    if { $_viewer == "" } {
1053        if { $_dim == 2 } {
1054            set _viewer contour
1055        } else {
1056            set _viewer isosurface
1057        }
1058    }
1059    set _comp2dims($cname) ${_dim}D
1060    lappend limits x [list $xmin $xmax]
1061    lappend limits y [list $ymin $ymax]
1062    lappend limits z [list $zmin $zmax]
1063    set dataAttrs [$dataset GetPointData]
1064    if { $dataAttrs == ""} {
1065        puts stderr "No point data"
1066    }
1067    set vmin 0
1068    set vmax 1
1069    set numArrays [$dataAttrs GetNumberOfArrays]
1070    if { $numArrays > 0 } {
1071        set array [$dataAttrs GetArray 0]
1072        foreach {vmin vmax} [$array GetRange] break
1073
1074        for {set i 0} {$i < [$dataAttrs GetNumberOfArrays] } {incr i} {
1075            set array [$dataAttrs GetArray $i]
1076            set name  [$dataAttrs GetArrayName $i]
1077            foreach {min max} [$array GetRange] break
1078            lappend limits $name [list $min $max]
1079            lappend _fieldNames $name
1080            lappend _fieldUnits ""
1081            lappend _fieldLabels $name
1082        }
1083    }
1084    lappend limits v [list $vmin $vmax]
1085    set _comp2limits($cname) $limits
1086    $reader Delete
1087    file delete $tmpfile
1088}
1089
1090#
1091# vtkdata --
1092#
1093#       Returns a string representing the mesh and field data for a specific
1094#       component in the legacy VTK file format.
1095#
1096itcl::body Rappture::Field::vtkdata {cname} {
1097    if {$cname == "component0"} {
1098        set cname "component"
1099    }
1100    # DX: Convert DX to VTK
1101    if {[info exists _comp2dx($cname)]} {
1102        return [Rappture::ConvertDxToVtk $_comp2dx($cname)]
1103    }
1104    # Unirect3d: isosurface
1105    if {[info exists _comp2unirect3d($cname)]} {
1106        return [$_comp2unirect3d($cname) vtkdata]
1107    }
1108    # VTK file data:
1109    if { [info exists _comp2vtk($cname)] } {
1110        return $_comp2vtk($cname)
1111    }
1112    # Points on mesh:  Construct VTK file output.
1113    if { [info exists _comp2mesh($cname)] } {
1114        # Data is in the form mesh and vector
1115        foreach {mesh vector} $_comp2mesh($cname) break
1116        set label [hints zlabel]
1117        if { $label == "" } {
1118            set label $cname
1119        } else {
1120            regsub -all { } $label {_} label
1121        }
1122        append out "# vtk DataFile Version 3.0\n"
1123        append out "[hints label]\n"
1124        append out "ASCII\n"
1125        append out [$mesh vtkdata]
1126        append out "POINT_DATA [$vector length]\n"
1127        append out "SCALARS $label float\n"
1128        append out "LOOKUP_TABLE default\n"
1129        append out "[$vector range 0 end]\n"
1130        return $out
1131    }
1132    error "can't find vtkdata for $cname. This method should only be called by the vtkheightmap widget"
1133}
1134
1135#
1136# BuildPointsOnMesh --
1137#
1138#       Parses the field XML description to build a mesh and values vector
1139#       representing the field.  Right now we handle the deprecated types
1140#       of "cloud", "unirect2d", and "unirect3d" (mostly for flows).
1141#
1142itcl::body Rappture::Field::BuildPointsOnMesh {cname} {
1143    #
1144    # More complex 2D/3D data is represented by a mesh
1145    # object and an associated vector for field values.
1146    #
1147    set path [$_field get $cname.mesh]
1148    if {[$_xmlobj element $path] == ""} {
1149        # Unknown mesh designated.
1150        return
1151    }
1152    set element [$_xmlobj element -as type $path]
1153    lappend _fieldNames $cname
1154    lappend _fieldLabels $cname
1155    lappend _fieldUnits ""
1156
1157    # Handle bizarre cases that hopefully will be deprecated.
1158    if { $element == "unirect3d" } {
1159        # Special case: unirect3d (should be deprecated) + flow.
1160        if { [$_field element $cname.extents] != "" } {
1161            set extents [$_field get $cname.extents]
1162        } else {
1163            set extents 1
1164        }
1165        set _dim 3
1166        set _viewer flowvis
1167        set _comp2dims($cname) "3D"
1168        set _comp2unirect3d($cname) \
1169            [Rappture::Unirect3d \#auto $_xmlobj $_field $cname $extents]
1170        set _comp2style($cname) [$_field get $cname.style]
1171        if {[$_field element $cname.flow] != ""} {
1172            set _comp2flowhints($cname) \
1173                [Rappture::FlowHints ::\#auto $_field $cname $_units]
1174        }
1175        incr _counter
1176        return
1177    }
1178    if { $element == "unirect2d" && [$_field element $cname.flow] != "" } {
1179        # Special case: unirect2d (normally deprecated) + flow.
1180        if { [$_field element $cname.extents] != "" } {
1181            set extents [$_field get $cname.extents]
1182        } else {
1183            set extents 1
1184        }
1185        set _dim 2
1186        set _viewer "flowvis"
1187        set _comp2dims($cname) "2D"
1188        set _comp2unirect2d($cname) \
1189            [Rappture::Unirect2d \#auto $_xmlobj $path]
1190        set _comp2style($cname) [$_field get $cname.style]
1191        set _comp2flowhints($cname) \
1192            [Rappture::FlowHints ::\#auto $_field $cname $_units]
1193        set _values [$_field get $cname.values]
1194        set limits {}
1195        foreach axis { x y } {
1196            lappend limits $axis [$_comp2unirect2d($cname) limits $axis]
1197        }
1198        set xv [blt::vector create \#auto]
1199        $xv set $_values
1200        lappend limits "v" [$xv limits]
1201        blt::vector destroy $xv
1202        set _comp2limits($cname) $limits
1203        incr _counter
1204        return
1205    }
1206    set _viewer "contour"
1207    switch -- $element {
1208        "cloud" {
1209            set mesh [Rappture::Cloud::fetch $_xmlobj $path]
1210        }
1211        "mesh" {
1212            set mesh [Rappture::Mesh::fetch $_xmlobj $path]
1213        }           
1214        "unirect2d" {
1215            set mesh [Rappture::Unirect2d::fetch $_xmlobj $path]
1216            set _viewer "heightmap"
1217        }
1218    }
1219    set _dim [$mesh dimensions]
1220    if {$_dim == 1} {
1221        # Is this used anywhere?
1222        #
1223        # OOPS!  This is 1D data
1224        # Forget the cloud/field -- store BLT vectors
1225        #
1226        # Is there a natural growth path in generating output from 1D to
1227        # higher dimensions?  If there isn't, let's kill this in favor
1228        # or explicitly using a <curve> instead.  Otherwise, the features
1229        # (methods such as xmarkers) or the <curve> need to be added
1230        # to the <field>.
1231        #
1232        set xv [blt::vector create x$_counter]
1233        set yv [blt::vector create y$_counter]
1234       
1235        $yv set [$mesh points]
1236        $xv seq 0 1 [$yv length]
1237        # sort x-coords in increasing order
1238        $xv sort $yv
1239       
1240        set _comp2dims($cname) "1D"
1241        set _comp2xy($cname) [list $xv $yv]
1242        incr _counter
1243        return
1244    } elseif {$_dim == 2} {
1245        set _type "heightmap"
1246        set vector [blt::vector create \#auto]
1247        $vector set [$_field get $cname.values]
1248        set _comp2dims($cname) "[$mesh dimensions]D"
1249        set _comp2mesh($cname) [list $mesh $vector]
1250        set _comp2style($cname) [$_field get $cname.style]
1251        incr _counter
1252        set label [hints zlabel]
1253        if { $label != "" } {
1254            set _fieldLabels $label
1255            set _fieldNames $label
1256            regsub -all { } $_fieldNames {_} _fieldNames
1257        }
1258        set units [hints zunits]
1259        if { $units != "" } {
1260            set _fieldUnits $units
1261        }
1262        array unset _comp2limits $cname
1263        lappend _comp2limits($cname) x [$mesh limits x]
1264        lappend _comp2limits($cname) y [$mesh limits y]
1265        lappend _comp2limits($cname) $label [$vector limits]
1266        lappend _comp2limits($cname) v [$vector limits]
1267        return
1268    } elseif {$_dim == 3} {
1269        #
1270        # 3D data: Store cloud/field as components
1271        #
1272        set values [$_field get $cname.values]
1273        set farray [vtkFloatArray ::vals$_counter]
1274        foreach v $values {
1275            if {"" != $_units} {
1276                set v [Rappture::Units::convert $v \
1277                           -context $_units -to $_units -units off]
1278            }
1279            $farray InsertNextValue $v
1280        }
1281        set _viewer "isosurface"
1282        set _type "isosurface"
1283        set vector [blt::vector create \#auto]
1284        $vector set [$_field get $cname.values]
1285        set _comp2dims($cname) "[$mesh dimensions]D"
1286        set _comp2mesh($cname) [list $mesh $vector]
1287        set _comp2style($cname) [$_field get $cname.style]
1288        incr _counter
1289        set label [hints zlabel]
1290        if { $label != "" } {
1291            set _fieldNames $label
1292            regsub -all { } $_fieldNames {_} _fieldNames
1293            set _fieldLabels $label
1294        }
1295        set units [hints zunits]
1296        if { $units != "" } {
1297            set _fieldUnits $units
1298        }
1299        lappend _comp2limits($cname) x [$mesh limits x]
1300        lappend _comp2limits($cname) y [$mesh limits y]
1301        lappend _comp2limits($cname) z [$mesh limits z]
1302        lappend _comp2limits($cname) $label [$vector limits]
1303        lappend _comp2limits($cname) v [$vector limits]
1304        return
1305    }
1306    error "unhandled case in field"
1307}
1308
1309itcl::body Rappture::Field::AvsToVtk { comp contents } {
1310    package require vtk
1311
1312    set reader $this-datasetreader
1313    vtkAVSucdReader $reader
1314
1315    # Write the contents to a file just in case it's binary.
1316    set tmpfile file[pid].vtk
1317    set f [open "$tmpfile" "w"]
1318    fconfigure $f -translation binary -encoding binary
1319    puts $f $contents
1320    close $f
1321    $reader SetFileName $tmpfile
1322    $reader Update
1323    file delete $tmpfile
1324
1325    set output [$reader GetOutput]
1326    set pointData [$output GetPointData]
1327    set _scalars {}
1328    for { set i 0 } { $i < [$pointData GetNumberOfArrays] } { incr i } {
1329        set name [$pointData GetArrayName $i]
1330        lappend _scalars $name $name "???"
1331    }
1332    set writer $this-datasetwriter
1333    vtkDataSetWriter $writer
1334    $writer SetInputConnection [$reader GetOutputPort]
1335    $writer SetFileName $tmpfile
1336    $writer Write
1337    rename $reader ""
1338    rename $writer ""
1339    set f [open "$tmpfile" "r"]
1340    fconfigure $f -translation binary -encoding binary
1341    set vtkdata [read $f]
1342    close $f
1343    file delete $tmpfile
1344    return $vtkdata
1345}
Note: See TracBrowser for help on using the repository browser.