Context Navigation

source: branches/1.4/gui/scripts/field.tcl @ 5459

Last change on this file since 5459 was 5459, checked in by ldelgass, 9 years ago
merge r5458 from trunk
File size: 61.6 KB

Line
1	# -- mode: tcl; indent-tabs-mode: nil --
2	# ----------------------------------------------------------------------
3	# COMPONENT: field - extracts data from an XML description of a field
4	#
5	# This object represents one field in an XML description of a device.
6	# It simplifies the process of extracting data vectors that represent
7	# the field.
8	# ======================================================================
9	# AUTHOR: Michael McLennan, Purdue University
10	# Copyright (c) 2004-2012 HUBzero Foundation, LLC
11	#
12	# See the file "license.terms" for information on usage and
13	# redistribution of this file, and for a DISCLAIMER OF ALL WARRANTIES.
14	# ======================================================================
15
16	# TODO:
17	# Vector field limits are wrong: need to compute magnitude limits and
18	# component-wise limits.
19
20	#
21	# Possible field dataset types:
22	#
23	# 1D Datasets
24	# 1D A field contained in a structure
25	# 2D Datasets
26	# vtk (range of z-axis is zero).
27	# unirect2d (deprecated)
28	# cloud (x,y point coordinates) (deprecated)
29	# mesh
30	# 3D Datasets
31	# vtk
32	# unirect3d (deprecated)
33	# cloud (x,y,z coordinates) (deprecated)
34	# mesh
35	# dx (FIXME: make dx-to-vtk converter work)
36	# ucd (AVS ucd format)
37	#
38	# Viewers:
39	# Format Dim Description Viewer Server
40	# 1D 1 structure field DeviceViewer1D N/A
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 unirect2d + extents > 1 flow flow nanovis
47	# unirect3d 3 unirect3d + extents > 1 flow flow nanovis
48	#
49	# With <views>, can specify which viewer for specific datasets. So it's OK
50	# for the same dataset to be viewed in more than one way.
51	# o Any 2D dataset can be viewed as a contour/heightmap.
52	# o Any 3D dataset can be viewed as a isosurface.
53	# o Any 2D dataset with vector data can be streamlines or flow.
54	# o Any 3D uniform rectilinear dataset can be viewed as a volume.
55	# o Any 3D dataset with vector data can be streamlines or flow.
56	#
57	# Need <views> to properly do things like qdot: volume with a polydata
58	# transparent shell. The view will combine the two objects <field> and
59	# <drawing> (??) into a single viewer.
60	#
61	package require Itcl
62	package require BLT
63
64	namespace eval Rappture {
65	# forward declaration
66	}
67
68	itcl::class Rappture::Field {
69	constructor {xmlobj path} {
70	# defined below
71	}
72	destructor {
73	# defined below
74	}
75	public method blob { cname }
76	public method components {args}
77	public method controls {option args}
78	public method extents {{cname -overall}}
79	public method fieldinfo { fname } {
80	lappend out $_fld2Label($fname)
81	lappend out $_fld2Units($fname)
82	lappend out $_fld2Components($fname)
83	return $out
84	}
85	public method fieldlimits {}
86	public method fieldnames { cname } {
87	if { ![info exists _comp2fldName($cname)] } {
88	return ""
89	}
90	return $_comp2fldName($cname)
91	}
92	public method flowhints { cname }
93	public method hints {{key ""}}
94	public method isunirect2d {}
95	public method isunirect3d {}
96	public method isvalid {} {
97	return $_isValid
98	}
99	public method limits {axis}
100	public method mesh {{cname -overall}}
101	public method numComponents {cname}
102	public method style { cname }
103	public method type {}
104	public method valueLimits { cname }
105	public method values { cname }
106	public method viewer {} {
107	return $_viewer
108	}
109	public method vtkdata {cname}
110	public method xErrorValues { cname } {
111	}
112	public method yErrorValues { cname } {
113	}
114
115	protected method Build {}
116	protected method _getValue {expr}
117	protected method GetAssociation { cname }
118	protected method GetTypeAndSize { cname }
119	protected method ReadVtkDataSet { cname contents }
120
121	private method AvsToVtk { cname contents }
122	private method DicomToVtk { cname contents }
123	private method BuildPointsOnMesh { cname }
124	private method InitHints {}
125	private method VerifyVtkDataSet { contents }
126	private method VectorLimits { vector vectorsize {comp -1} }
127	private method VtkDataSetToXy { dataset }
128
129	protected variable _dim 0; # Dimension of the mesh
130
131	private variable _xmlobj ""; # ref to XML obj with field data
132	private variable _path ""; # Path of this object in the XML
133	private variable _field ""; # This field element as XML obj
134
135	private variable _type ""; # Field type: e.g. file type
136	private variable _hints; # Hints array
137	private variable _limits; # maps axis name => {z0 z1} limits
138	private variable _units ""; # system of units for this field
139	private variable _viewer ""; # Hints which viewer to use
140	private variable _isValid 0; # Indicates if the field contains
141	# valid data.
142	private variable _isValidComponent; # Array of valid components found
143	private variable _zmax 0; # length of the device (1D only)
144
145	private variable _fld2Components; # field name => number of components
146	private variable _fld2Label; # field name => label
147	private variable _fld2Units; # field name => units
148
149	private variable _comp2fldName; # cname => field names.
150	private variable _comp2type; # cname => type (e.g. "vectors")
151	private variable _comp2size; # cname => # of components in element
152	private variable _comp2assoc; # cname => association (e.g. pointdata)
153	private variable _comp2dims; # cname => dimensionality
154	private variable _comp2xy; # cname => x,y vectors
155	private variable _comp2vtk; # cname => vtk file data
156	private variable _comp2dx; # cname => OpenDX data
157	private variable _comp2unirect2d; # cname => unirect2d obj
158	private variable _comp2unirect3d; # cname => unirect3d obj
159	private variable _comp2style; # cname => style settings
160	private variable _comp2cntls; # cname => x,y control points (1D only)
161	private variable _comp2extents; # cname => extents (Only for unirect)
162	private variable _comp2limits; # Array of limits per component
163	private variable _comp2flowhints
164	private variable _comp2mesh; # list: mesh obj, BLT vector of values
165
166	private variable _values ""; # Only for unirect2d - list of values
167
168	private common _alwaysConvertDX 0; # If set, convert DX and store as VTK,
169	# even if viewer is nanovis/flowvis
170	private common _counter 0; # counter for unique vector names
171	}
172
173	# ----------------------------------------------------------------------
174	# CONSTRUCTOR
175	# ----------------------------------------------------------------------
176	itcl::body Rappture::Field::constructor {xmlobj path} {
177	package require vtk
178	if {![Rappture::library isvalid $xmlobj]} {
179	error "bad value \"$xmlobj\": should be Rappture::library"
180	}
181	set _xmlobj $xmlobj
182	set _path $path
183	set _field [$xmlobj element -as object $path]
184	set _units [$_field get units]
185
186	set xunits [$xmlobj get units]
187	if {"" == $xunits \|\| "arbitrary" == $xunits} {
188	set xunits "um"
189	}
190
191	# determine the overall size of the device
192	set z0 [set z1 0]
193	foreach elem [$_xmlobj children components] {
194	switch -glob -- $elem {
195	box* {
196	if {![regexp {[0-9]$} $elem]} {
197	set elem "${elem}0"
198	}
199	set z0 [$_xmlobj get components.$elem.corner0]
200	set z0 [Rappture::Units::convert $z0 \
201	-context $xunits -to $xunits -units off]
202
203	set z1 [$_xmlobj get components.$elem.corner1]
204	set z1 [Rappture::Units::convert $z1 \
205	-context $xunits -to $xunits -units off]
206
207	set _limits($elem) [list $z0 $z1]
208	}
209	}
210	}
211	set _zmax $z1
212
213	# build up vectors for various components of the field
214	Build
215	InitHints
216	}
217
218	# ----------------------------------------------------------------------
219	# DESTRUCTOR
220	# ----------------------------------------------------------------------
221	itcl::body Rappture::Field::destructor {} {
222	itcl::delete object $_field
223	# don't destroy the _xmlobj! we don't own it!
224
225	foreach name [array names _comp2xy] {
226	eval blt::vector destroy $_comp2xy($name)
227	}
228	foreach name [array names _comp2unirect2d] {
229	itcl::delete object $_comp2unirect2d($name)
230	}
231	foreach name [array names _comp2unirect3d] {
232	itcl::delete object $_comp2unirect3d($name)
233	}
234	foreach name [array names _comp2flowhints] {
235	itcl::delete object $_comp2flowhints($name)
236	}
237	foreach name [array names _comp2mesh] {
238	# Data is in the form of a mesh and a vector.
239	foreach { mesh vector } $_comp2mesh($name) break
240	# Release the mesh (may be shared)
241	set class [$mesh info class]
242	${class}::release $mesh
243	# Destroy the vector
244	blt::vector destroy $vector
245	}
246	}
247
248	# ----------------------------------------------------------------------
249	# USAGE: components ?-name\|-dimensions\|-style? ?<pattern>?
250	#
251	# Returns a list of names or types for the various components of
252	# this field. If the optional glob-style <pattern> is specified,
253	# then it returns only the components with names matching the pattern.
254	# ----------------------------------------------------------------------
255	itcl::body Rappture::Field::components {args} {
256	Rappture::getopts args params {
257	flag what -name default
258	flag what -dimensions
259	flag what -style
260	flag what -particles
261	flag what -flow
262	flag what -box
263	}
264
265	set pattern *
266	if {[llength $args] > 0} {
267	set pattern [lindex $args 0]
268	set args [lrange $args 1 end]
269	}
270	if {[llength $args] > 0} {
271	error "wrong # args: should be \"components ?switches? ?pattern?\""
272	}
273
274	# There's only one dimension of the field. Components can't have
275	# different dimensions in the same field. They would by definition be
276	# using different meshes and viewers.
277	if { $params(what) == "-dimensions" } {
278	return "${_dim}D"
279	}
280	# BE CAREFUL: return component names in proper order
281	set rlist ""
282	set components {}
283	# First compile a list of valid components that match the pattern
284	foreach cname [$_field children -type component] {
285	if { ![info exists _isValidComponent($cname)] } {
286	continue
287	}
288	if { [string match $pattern $cname] } {
289	lappend components $cname
290	}
291	}
292	# Now handle the tests.
293	switch -- $params(what) {
294	-name {
295	set rlist $components
296	}
297	-style {
298	foreach cname $components {
299	if { [info exists _comp2style($cname)] } {
300	lappend rlist $_comp2style($cname)
301	}
302	}
303	}
304	}
305	return $rlist
306	}
307
308	# ----------------------------------------------------------------------
309	# USAGE: mesh ?<name>?
310	#
311	# For 1D data (curve), returns a BLT vector of x values for the field
312	# component <name>. Otherwise, this method is unused
313	# ----------------------------------------------------------------------
314	itcl::body Rappture::Field::mesh {{cname -overall}} {
315	if {$cname == "-overall" \|\| $cname == "component0"} {
316	set cname [lindex [components -name] 0]
317	}
318	if {[info exists _comp2xy($cname)]} {
319	return [lindex $_comp2xy($cname) 0] ;# return xv
320	}
321	if {[info exists _comp2vtk($cname)]} {
322	# FIXME: extract mesh from VTK file data.
323	error "method \"mesh\" is not implemented for VTK file data"
324	}
325	if {[info exists _comp2dx($cname)]} {
326	error "method \"mesh\" is not implemented for DX file data"
327	}
328	if {[info exists _comp2mesh($cname)]} {
329	# FIXME: This only works for cloud
330	set mesh [lindex $_comp2mesh($cname) 0]
331	return [$mesh points]
332	}
333	if {[info exists _comp2unirect2d($cname)]} {
334	# FIXME: unirect2d mesh is a list: xMin xMax xNum yMin yMax yNum
335	return [$_comp2unirect2d($cname) mesh]
336	}
337	if {[info exists _comp2unirect3d($cname)]} {
338	# This returns a list of x,y,z points
339	return [$_comp2unirect3d($cname) mesh]
340	}
341	error "can't get field mesh: Unknown component \"$cname\": should be one of [join [lsort [array names _comp2dims]] {, }]"
342	}
343
344	# ----------------------------------------------------------------------
345	# USAGE: values ?<name>?
346	#
347	# For 1D data (curve), returns a BLT vector of field values (y coords)
348	# for the field component <name>. Otherwise, this method is unused
349	# ----------------------------------------------------------------------
350	itcl::body Rappture::Field::values {cname} {
351	if {$cname == "component0"} {
352	set cname "component"
353	}
354	if {[info exists _comp2xy($cname)]} {
355	return [lindex $_comp2xy($cname) 1] ;# return yv
356	}
357	if { [info exists _comp2vtk($cname)] } {
358	# FIXME: extract the values from the VTK file data
359	error "method \"values\" is not implemented for VTK file data"
360	}
361	if {[info exists _comp2dx($cname)]} {
362	error "method \"values\" is not implemented for DX file data"
363	}
364	if { [info exists _comp2mesh($cname)] } {
365	set vector [lindex $_comp2mesh($cname) 1]
366	return [$vector range 0 end]
367	}
368	if {[info exists _comp2unirect2d($cname)]} {
369	return $_values
370	}
371	if {[info exists _comp2unirect3d($cname)]} {
372	return [$_comp2unirect3d($cname) values]
373	}
374	error "can't get field values. Unknown component \"$cname\": should be one of [join [lsort [array names _comp2dims]] {, }]"
375	}
376
377	# ----------------------------------------------------------------------
378	# USAGE: blob ?<name>?
379	#
380	# Returns a string representing the blob of data for the mesh and values.
381	# ----------------------------------------------------------------------
382	itcl::body Rappture::Field::blob {cname} {
383	if {$cname == "component0"} {
384	set cname "component"
385	}
386	if {[info exists _comp2dx($cname)]} {
387	return $_comp2dx($cname) ;# return gzipped, base64-encoded DX data
388	}
389	if {[info exists _comp2unirect2d($cname)]} {
390	set blob [$_comp2unirect2d($cname) blob]
391	lappend blob "values" $_values
392	return $blob
393	}
394	if {[info exists _comp2unirect3d($cname)]} {
395	return [$_comp2unirect3d($cname) blob]
396	}
397	if { [info exists _comp2vtk($cname)] } {
398	error "blob not implemented for VTK file data"
399	}
400	if {[info exists _comp2xy($cname)]} {
401	error "blob not implemented for XY data"
402	}
403	error "can't get field blob: Unknown component \"$cname\": should be one of [join [lsort [array names _comp2dims]] {, }]"
404	}
405
406	# ----------------------------------------------------------------------
407	# USAGE: valueLimits <cname>
408	#
409	# Returns an array for the requested component with a list {min max}
410	# representing the limits for each axis.
411	# ----------------------------------------------------------------------
412	itcl::body Rappture::Field::valueLimits { cname } {
413	if { [info exists _comp2limits($cname)] } {
414	return $_comp2limits($cname)
415	}
416	return ""
417	}
418
419	# ----------------------------------------------------------------------
420	# USAGE: limits <axis>
421	#
422	# Returns a list {min max} representing the limits for the specified
423	# axis.
424	# ----------------------------------------------------------------------
425	itcl::body Rappture::Field::limits {which} {
426	set min ""
427	set max ""
428
429	foreach cname [array names _comp2dims] {
430	switch -- $_comp2dims($cname) {
431	1D {
432	switch -- $which {
433	x - xlin {
434	set pos 0; set log 0; set axis x
435	}
436	xlog {
437	set pos 0; set log 1; set axis x
438	}
439	y - ylin - v - vlin {
440	set pos 1; set log 0; set axis y
441	}
442	ylog - vlog {
443	set pos 1; set log 1; set axis y
444	}
445	default {
446	error "bad axis \"$which\": should be x, xlin, xlog, y, ylin, ylog, v, vlin, vlog"
447	}
448	}
449
450	set vname [lindex $_comp2xy($cname) $pos]
451
452	if {$log} {
453	blt::vector tmp zero
454	# on a log scale, use abs value and ignore zeros
455	$vname dup tmp
456	$vname dup zero
457	zero expr {tmp == 0} ;# find the zeros
458	tmp expr {abs(tmp)} ;# get the abs value
459	tmp expr {tmp + zero*max(tmp)} ;# replace 0s with abs max
460	set axisMin [blt::vector expr min(tmp)]
461	set axisMax [blt::vector expr max(tmp)]
462	blt::vector destroy tmp zero
463	} else {
464	set axisMin [$vname min]
465	set axisMax [$vname max]
466	}
467
468	if {"" == $min} {
469	set min $axisMin
470	} elseif {$axisMin < $min} {
471	set min $axisMin
472	}
473	if {"" == $max} {
474	set max $axisMax
475	} elseif {$axisMax > $max} {
476	set max $axisMax
477	}
478	}
479	default {
480	if {[info exists _comp2limits($cname)]} {
481	array set limits $_comp2limits($cname)
482	switch -- $which {
483	x - xlin - xlog {
484	set axis x
485	foreach {axisMin axisMax} $limits(x) break
486	}
487	y - ylin - ylog {
488	set axis y
489	foreach {axisMin axisMax} $limits(y) break
490	}
491	z - zlin - zlog {
492	set axis z
493	foreach {axisMin axisMax} $limits(z) break
494	}
495	v - vlin - vlog {
496	set axis v
497	foreach {axisMin axisMax} $limits(v) break
498	}
499	default {
500	if { ![info exists limits($which)] } {
501	error "limits: unknown axis \"$which\""
502	}
503	set axis v
504	foreach {axisMin axisMax} $limits($which) break
505	}
506	}
507	} else {
508	set axisMin 0 ;# HACK ALERT! must be OpenDX data
509	set axisMax 1
510	set axis v
511	}
512	}
513	}
514	if { "" == $min \|\| $axisMin < $min } {
515	set min $axisMin
516	}
517	if { "" == $max \|\| $axisMax > $max } {
518	set max $axisMax
519	}
520	}
521	set val [$_field get "${axis}axis.min"]
522	if {"" != $val && "" != $min} {
523	if {$val > $min} {
524	# tool specified this min -- don't go any lower
525	set min $val
526	}
527	}
528	set val [$_field get "${axis}axis.max"]
529	if {"" != $val && "" != $max} {
530	if {$val < $max} {
531	# tool specified this max -- don't go any higher
532	set max $val
533	}
534	}
535	return [list $min $max]
536	}
537
538	# ----------------------------------------------------------------------
539	# USAGE: fieldlimits
540	#
541	# Returns a list {min max} representing the limits for the specified
542	# axis.
543	# ----------------------------------------------------------------------
544	itcl::body Rappture::Field::fieldlimits {} {
545	foreach cname [array names _comp2limits] {
546	array set limits $_comp2limits($cname)
547	foreach fname [fieldnames $cname] {
548	if { ![info exists limits($fname)] } {
549	puts stderr "ERROR: field \"$fname\" unknown in \"$cname\""
550	continue
551	}
552	foreach {min max} $limits($fname) break
553	if { ![info exists overall($fname)] } {
554	set overall($fname) $limits($fname)
555	continue
556	}
557	foreach {omin omax} $overall($fname) break
558	if { $min < $omin } {
559	set omin $min
560	}
561	if { $max > $omax } {
562	set omax $max
563	}
564	set overall($fname) [list $min $max]
565	}
566	}
567	if { [info exists overall] } {
568	return [array get overall]
569	}
570	return ""
571	}
572
573	# ----------------------------------------------------------------------
574	# USAGE: controls get ?<name>?
575	# USAGE: controls validate <path> <value>
576	# USAGE: controls put <path> <value>
577	#
578	# Returns a list {path1 x1 y1 val1 path2 x2 y2 val2 ...} representing
579	# control points for the specified field component <name>.
580	# ----------------------------------------------------------------------
581	itcl::body Rappture::Field::controls {option args} {
582	switch -- $option {
583	get {
584	set cname [lindex $args 0]
585	if {[info exists _comp2cntls($cname)]} {
586	return $_comp2cntls($cname)
587	}
588	return ""
589	}
590	validate {
591	set path [lindex $args 0]
592	set value [lindex $args 1]
593	set units [$_xmlobj get $path.units]
594
595	if {"" != $units} {
596	set nv [Rappture::Units::convert \
597	$value -context $units -to $units -units off]
598	} else {
599	set nv $value
600	}
601	if {![string is double $nv]
602	\|\| [regexp -nocase {^(inf\|nan)$} $nv]} {
603	error "Value out of range"
604	}
605
606	set rawmin [$_xmlobj get $path.min]
607	if {"" != $rawmin} {
608	set minv $rawmin
609	if {"" != $units} {
610	set minv [Rappture::Units::convert \
611	$minv -context $units -to $units -units off]
612	set nv [Rappture::Units::convert \
613	$value -context $units -to $units -units off]
614	}
615	# fix for the case when the user tries to
616	# compare values like minv=-500 nv=-0600
617	set nv [format "%g" $nv]
618	set minv [format "%g" $minv]
619
620	if {$nv < $minv} {
621	error "Minimum value allowed here is $rawmin"
622	}
623	}
624
625	set rawmax [$_xmlobj get $path.max]
626	if {"" != $rawmax} {
627	set maxv $rawmax
628	if {"" != $units} {
629	set maxv [Rappture::Units::convert \
630	$maxv -context $units -to $units -units off]
631	set nv [Rappture::Units::convert \
632	$value -context $units -to $units -units off]
633	}
634	# fix for the case when the user tries to
635	# compare values like maxv=-500 nv=-0600
636	set nv [format "%g" $nv]
637	set maxv [format "%g" $maxv]
638
639	if {$nv > $maxv} {
640	error "Maximum value allowed here is $rawmax"
641	}
642	}
643
644	return "ok"
645	}
646	put {
647	set path [lindex $args 0]
648	set value [lindex $args 1]
649	$_xmlobj put $path.current $value
650	Build
651	}
652	default {
653	error "bad field controls option \"$option\": should be get or put"
654	}
655	}
656	}
657
658	# ----------------------------------------------------------------------
659	# USAGE: hints ?<keyword>?
660	#
661	# Returns a list of key/value pairs for various hints about plotting
662	# this field. If a particular <keyword> is specified, then it returns
663	# the hint for that <keyword>, if it exists.
664	# ----------------------------------------------------------------------
665	itcl::body Rappture::Field::hints {{keyword ""}} {
666	if { $keyword != "" } {
667	if {[info exists _hints($keyword)]} {
668	return $_hints($keyword)
669	}
670	return ""
671	}
672	return [array get _hints]
673	}
674
675
676	# ----------------------------------------------------------------------
677	# USAGE: InitHints
678	#
679	# Returns a list of key/value pairs for various hints about plotting
680	# this field. If a particular <keyword> is specified, then it returns
681	# the hint for that <keyword>, if it exists.
682	# ----------------------------------------------------------------------
683	itcl::body Rappture::Field::InitHints {} {
684	foreach {key path} {
685	camera camera.position
686	color about.color
687	default about.default
688	group about.group
689	label about.label
690	scale about.scale
691	seeds about.seeds
692	style about.style
693	type about.type
694	xlabel about.xaxis.label
695	ylabel about.yaxis.label
696	zlabel about.zaxis.label
697	xunits about.xaxis.units
698	yunits about.yaxis.units
699	zunits about.zaxis.units
700	units units
701	updir updir
702	vectors about.vectors
703	} {
704	set str [$_field get $path]
705	if { "" != $str } {
706	set _hints($key) $str
707	}
708	}
709	foreach cname [components] {
710	if { ![info exists _comp2mesh($cname)] } {
711	continue
712	}
713	set mesh [lindex $_comp2mesh($cname) 0]
714	foreach axis {x y z} {
715	if { ![info exists _hints(${axis}units)] } {
716	set _hints(${axis}units) [$mesh units $axis]
717	}
718	if { ![info exists _hints(${axis}label)] } {
719	set _hints(${axis}label) [$mesh label $axis]
720	}
721	}
722	}
723	foreach {key path} {
724	toolid tool.id
725	toolname tool.name
726	toolcommand tool.execute
727	tooltitle tool.title
728	toolrevision tool.version.application.revision
729	} {
730	set str [$_xmlobj get $path]
731	if { "" != $str } {
732	set _hints($key) $str
733	}
734	}
735	# Set toolip and path hints
736	set _hints(path) $_path
737	if { [info exists _hints(group)] && [info exists _hints(label)] } {
738	# pop-up help for each curve
739	set _hints(tooltip) $_hints(label)
740	}
741	}
742
743	# ----------------------------------------------------------------------
744	# USAGE: Build
745	#
746	# Used internally to build up the vector representation for the
747	# field when the object is first constructed, or whenever the field
748	# data changes. Discards any existing vectors and builds everything
749	# from scratch.
750	# ----------------------------------------------------------------------
751	itcl::body Rappture::Field::Build {} {
752
753	# Discard any existing data
754	foreach name [array names _comp2xy] {
755	eval blt::vector destroy $_comp2xy($name)
756	}
757	array unset _comp2vtk
758	foreach name [array names _comp2unirect2d] {
759	eval itcl::delete object $_comp2unirect2d($name)
760	}
761	foreach name [array names _comp2unirect3d] {
762	eval itcl::delete object $_comp2unirect3d($name)
763	}
764	catch {unset _comp2xy}
765	catch {unset _comp2dx}
766	catch {unset _comp2dims}
767	catch {unset _comp2style}
768	array unset _comp2unirect2d
769	array unset _comp2unirect3d
770	array unset _comp2extents
771	array unset _dataobj2type
772	#
773	# Scan through the components of the field and create
774	# vectors for each part.
775	#
776	array unset _isValidComponent
777	foreach cname [$_field children -type component] {
778	set type ""
779	if { ([$_field element $cname.constant] != "" &&
780	[$_field element $cname.domain] != "") \|\|
781	[$_field element $cname.xy] != "" } {
782	set type "1D"
783	} elseif { [$_field element $cname.mesh] != "" &&
784	[$_field element $cname.values] != ""} {
785	set type "points-on-mesh"
786	} elseif { [$_field element $cname.vtk] != ""} {
787	set type "vtk"
788	set viewer [$_field get "about.view"]
789	if { $viewer != "" } {
790	set _viewer $viewer
791	}
792	} elseif {[$_field element $cname.opendx] != ""} {
793	global env
794	if { [info exists env(VTKVOLUME)] } {
795	set _viewer "vtkvolume"
796	}
797	set type "opendx"
798	} elseif {[$_field element $cname.dx] != ""} {
799	global env
800	if { [info exists env(VTKVOLUME)] } {
801	set _viewer "vtkvolume"
802	}
803	set type "dx"
804	} elseif {[$_field element $cname.ucd] != ""} {
805	set type "ucd"
806	} elseif {[$_field element $cname.dicom] != ""} {
807	set type "dicom"
808	}
809	set _comp2style($cname) ""
810	if { $type == "" } {
811	puts stderr "WARNING: Ignoring field component \"$_path.$cname\": no data found."
812	continue
813	}
814	# Save the extents of the component
815	if { [$_field element $cname.extents] != "" } {
816	set extents [$_field get $cname.extents]
817	} else {
818	set extents 1
819	}
820	set _comp2extents($cname) $extents
821	set _type $type
822
823	GetTypeAndSize $cname
824	GetAssociation $cname
825	if { $_comp2size($cname) > 1 } {
826	set viewer [$_field get "about.view"]
827	if { $viewer == "" } {
828	set _viewer "glyphs"
829	}
830	}
831	if {$type == "1D"} {
832	#
833	# 1D data can be represented as 2 BLT vectors,
834	# one for x and the other for y.
835	#
836	set xv ""
837	set yv ""
838	set _dim 1
839	set val [$_field get $cname.constant]
840	if {$val != ""} {
841	set domain [$_field get $cname.domain]
842	if {$domain == "" \|\| ![info exists _limits($domain)]} {
843	set z0 0
844	set z1 $_zmax
845	} else {
846	foreach {z0 z1} $_limits($domain) { break }
847	}
848	set xv [blt::vector create x$_counter]
849	$xv append $z0 $z1
850
851	foreach {val pcomp} [_getValue $val] break
852	set yv [blt::vector create y$_counter]
853	$yv append $val $val
854
855	if {$pcomp != ""} {
856	set zm [expr {0.5*($z0+$z1)}]
857	set _comp2cntls($cname) \
858	[list $pcomp $zm $val "$val$_units"]
859	}
860	} else {
861	set xydata [$_field get $cname.xy]
862	if {"" != $xydata} {
863	set xv [blt::vector create x$_counter]
864	set yv [blt::vector create y$_counter]
865	set tmp [blt::vector create \#auto]
866	$tmp set $xydata
867	$tmp split $xv $yv
868	blt::vector destroy $tmp
869	}
870	}
871
872	if {$xv != "" && $yv != ""} {
873	# sort x-coords in increasing order
874	$xv sort $yv
875	set _dim 1
876	set _comp2dims($cname) "1D"
877	set _comp2xy($cname) [list $xv $yv]
878	incr _counter
879	}
880	} elseif {$type == "points-on-mesh"} {
881	if { ![BuildPointsOnMesh $cname] } {
882	continue; # Ignore this component
883	}
884	} elseif {$type == "vtk"} {
885	set contents [$_field get $cname.vtk]
886	if { $contents == "" } {
887	puts stderr "WARNING: No data for \"$_path.$cname.vtk\""
888	continue; # Ignore this component
889	}
890	ReadVtkDataSet $cname $contents
891	set _comp2vtk($cname) $contents
892	set _comp2style($cname) [$_field get $cname.style]
893	incr _counter
894	} elseif {$type == "dx" \|\| $type == "opendx" } {
895	#
896	# Extract gzipped, base64-encoded OpenDX data
897	#
898	set viewer [$_field get "about.view"]
899	if { $viewer != "" } {
900	set _viewer $viewer
901	}
902	if { $_viewer == "" } {
903	if {[$_field element $cname.flow] != ""} {
904	set _viewer "flowvis"
905	} else {
906	set _viewer "nanovis"
907	}
908	}
909	set data [$_field get -decode no $cname.$type]
910	set contents [Rappture::encoding::decode -as zb64 $data]
911	if { $contents == "" } {
912	puts stderr "WARNING: No data for \"$_path.$cname.$type\""
913	continue; # Ignore this component
914	}
915	if 0 {
916	set f [open /tmp/$_path.$cname.dx "w"]
917	puts -nonewline $f $contents
918	close $f
919	}
920	# Convert to VTK
921	if { [catch { Rappture::DxToVtk $contents } vtkdata] == 0 } {
922	unset contents
923	# Read back VTK: this will set the field limits and the mesh
924	# dimensions based on the bounds (sets _dim). We rely on this
925	# conversion for limits even if we send DX data to nanovis.
926	ReadVtkDataSet $cname $vtkdata
927	if 0 {
928	set f [open /tmp/$_path.$cname.vtk "w"]
929	puts -nonewline $f $vtkdata
930	close $f
931	}
932	} else {
933	unset contents
934	# vtkdata variable holds error message
935	puts stderr "Can't parse DX data\n$vtkdata"
936	continue; # Ignore this component
937	}
938	if { $_alwaysConvertDX \|\|
939	($_viewer != "nanovis" && $_viewer != "flowvis") } {
940	set _type "vtk"
941	set _comp2vtk($cname) $vtkdata
942	} else {
943	set _type "dx"
944	set _comp2dx($cname) $data
945	}
946	unset data
947	unset vtkdata
948	set _comp2style($cname) [$_field get $cname.style]
949	if {[$_field element $cname.flow] != ""} {
950	set _comp2flowhints($cname) \
951	[Rappture::FlowHints ::\#auto $_field $cname $_units]
952	}
953	incr _counter
954	} elseif { $type == "dicom"} {
955	set contents [$_field get $cname.dicom]
956	if { $contents == "" } {
957	continue; # Ignore this component
958	}
959	set viewer [$_field get "about.view"]
960	if { $viewer != "" } {
961	set _viewer $viewer
962	}
963	set vtkdata [DicomToVtk $cname $contents]
964	if { $_viewer == "" } {
965	set _viewer [expr {($_dim == 3) ? "vtkvolume" : "vtkimage"}]
966	}
967	set _comp2vtk($cname) $vtkdata
968	set _comp2style($cname) [$_field get $cname.style]
969	incr _counter
970	} elseif { $type == "ucd"} {
971	set contents [$_field get $cname.ucd]
972	if { $contents == "" } {
973	continue; # Ignore this component
974	}
975	set vtkdata [AvsToVtk $cname $contents]
976	ReadVtkDataSet $cname $vtkdata
977	set _comp2vtk($cname) $vtkdata
978	set _comp2style($cname) [$_field get $cname.style]
979	incr _counter
980	}
981	set _isValidComponent($cname) 1
982	#puts stderr "Field $cname type is: $_type"
983	}
984	if { [array size _isValidComponent] == 0 } {
985	puts stderr "ERROR: All components of field \"$_path\" are invalid."
986	return 0
987	}
988	# Sanity check. Verify that all components of the field have the same
989	# dimension.
990	set dim ""
991	foreach cname [array names _comp2dims] {
992	if { $dim == "" } {
993	set dim $_comp2dims($cname)
994	continue
995	}
996	if { $dim != $_comp2dims($cname) } {
997	puts stderr "WARNING: A field can't have components of different dimensions: [join [array get _comp2dims] ,]"
998	return 0
999	}
1000	}
1001
1002	# FIXME: about.scalars and about.vectors are temporary. With views
1003	# the label and units for each field will be specified there.
1004	#
1005	# FIXME: Test that every <field><component> has the same field names,
1006	# units, components.
1007	#
1008	# Override what we found in the VTK file with names that the user
1009	# selected. We override the field label and units.
1010	foreach { fname label units } [$_field get about.scalars] {
1011	if { ![info exists _fld2Name($fname)] } {
1012	set _fld2Name($fname) $fname
1013	set _fld2Components($fname) 1
1014	}
1015	set _fld2Label($fname) $label
1016	set _fld2Units($fname) $units
1017	}
1018	foreach { fname label units } [$_field get about.vectors] {
1019	if { ![info exists _fld2Name($fname)] } {
1020	set _fld2Name($fname) $fname
1021	# We're just marking the field as vector (> 1) for now.
1022	set _fld2Components($fname) 3
1023	}
1024	set _fld2Label($fname) $label
1025	set _fld2Units($fname) $units
1026	}
1027	set _isValid 1
1028	return 1
1029	}
1030
1031	# ----------------------------------------------------------------------
1032	# USAGE: _getValue <expr>
1033	#
1034	# Used internally to get the value for an expression <expr>. Returns
1035	# a list of the form {val parameterPath}, where val is the numeric
1036	# value of the expression, and parameterPath is the XML path to the
1037	# parameter representing the value, or "" if the <expr> does not
1038	# depend on any parameters.
1039	# ----------------------------------------------------------------------
1040	itcl::body Rappture::Field::_getValue {expr} {
1041	#
1042	# First, look for the expression among the <parameter>'s
1043	# associated with the device.
1044	#
1045	set found 0
1046	foreach pcomp [$_xmlobj children parameters] {
1047	set id [$_xmlobj element -as id parameters.$pcomp]
1048	if {[string equal $id $expr]} {
1049	set val [$_xmlobj get parameters.$pcomp.current]
1050	if {"" == $val} {
1051	set val [$_xmlobj get parameters.$pcomp.default]
1052	}
1053	if {"" != $val} {
1054	set expr $val
1055	set found 1
1056	break
1057	}
1058	}
1059	}
1060	if {$found} {
1061	set pcomp "parameters.$pcomp"
1062	} else {
1063	set pcomp ""
1064	}
1065
1066	if {$_units != ""} {
1067	set expr [Rappture::Units::convert $expr \
1068	-context $_units -to $_units -units off]
1069	}
1070
1071	return [list $expr $pcomp]
1072	}
1073
1074	#
1075	# isunirect2d --
1076	#
1077	# Returns if the field is a unirect2d object.
1078	#
1079	itcl::body Rappture::Field::isunirect2d { } {
1080	return [expr [array size _comp2unirect2d] > 0]
1081	}
1082
1083	#
1084	# isunirect3d --
1085	#
1086	# Returns if the field is a unirect3d object.
1087	#
1088	itcl::body Rappture::Field::isunirect3d { } {
1089	return [expr [array size _comp2unirect3d] > 0]
1090	}
1091
1092	#
1093	# flowhints --
1094	#
1095	# Returns the hints associated with a flow vector field.
1096	#
1097	itcl::body Rappture::Field::flowhints { cname } {
1098	if { [info exists _comp2flowhints($cname)] } {
1099	return $_comp2flowhints($cname)
1100	}
1101	return ""
1102	}
1103
1104	#
1105	# style --
1106	#
1107	# Returns the style associated with a component of the field.
1108	#
1109	itcl::body Rappture::Field::style { cname } {
1110	if { [info exists _comp2style($cname)] } {
1111	return $_comp2style($cname)
1112	}
1113	return ""
1114	}
1115
1116	#
1117	# type --
1118	#
1119	# Returns the data storage type of the field.
1120	#
1121	# FIXME: What are the valid types?
1122	#
1123	itcl::body Rappture::Field::type {} {
1124	return $_type
1125	}
1126
1127	#
1128	# numComponents --
1129	#
1130	# Returns the number of components in the field component.
1131	#
1132	itcl::body Rappture::Field::numComponents {cname} {
1133	set name $cname
1134	regsub -all { } $name {_} name
1135	if {[info exists _fld2Components($name)] } {
1136	return $_fld2Components($name)
1137	}
1138	return 1; # Default to scalar.
1139	}
1140
1141	#
1142	# extents --
1143	#
1144	# Returns if the field is a unirect2d object.
1145	#
1146	itcl::body Rappture::Field::extents {{cname -overall}} {
1147	if {$cname == "-overall" } {
1148	set max 0
1149	foreach cname [$_field children -type component] {
1150	if { ![info exists _comp2extents($cname)] } {
1151	continue
1152	}
1153	set value $_comp2extents($cname)
1154	if { $max < $value } {
1155	set max $value
1156	}
1157	}
1158	return $max
1159	}
1160	if { $cname == "component0"} {
1161	set cname [lindex [components -name] 0]
1162	}
1163	return $_comp2extents($cname)
1164	}
1165
1166	itcl::body Rappture::Field::VerifyVtkDataSet { contents } {
1167	package require vtk
1168
1169	set reader $this-datasetreader
1170	vtkDataSetReader $reader
1171
1172	# Write the contents to a file just in case it's binary.
1173	set tmpfile file[pid].vtk
1174	set f [open "$tmpfile" "w"]
1175	fconfigure $f -translation binary -encoding binary
1176	puts $f $contents
1177	close $f
1178
1179	$reader SetFileName $tmpfile
1180	$reader ReadAllNormalsOn
1181	$reader ReadAllTCoordsOn
1182	$reader ReadAllScalarsOn
1183	$reader ReadAllColorScalarsOn
1184	$reader ReadAllVectorsOn
1185	$reader ReadAllTensorsOn
1186	$reader ReadAllFieldsOn
1187	$reader Update
1188	file delete $tmpfile
1189
1190	set dataset [$reader GetOutput]
1191	set dataAttrs [$dataset GetPointData]
1192	if { $dataAttrs == ""} {
1193	puts stderr "WARNING: No point data found in \"$_path\""
1194	rename $reader ""
1195	return 0
1196	}
1197	rename $reader ""
1198	}
1199
1200	itcl::body Rappture::Field::ReadVtkDataSet { cname contents } {
1201	package require vtk
1202
1203	set reader $this-datasetreader
1204	vtkDataSetReader $reader
1205
1206	# Write the contents to a file just in case it's binary.
1207	set tmpfile file[pid].vtk
1208	set f [open "$tmpfile" "w"]
1209	fconfigure $f -translation binary -encoding binary
1210	puts $f $contents
1211	close $f
1212
1213	$reader SetFileName $tmpfile
1214	$reader ReadAllNormalsOn
1215	$reader ReadAllTCoordsOn
1216	$reader ReadAllScalarsOn
1217	$reader ReadAllColorScalarsOn
1218	$reader ReadAllVectorsOn
1219	$reader ReadAllTensorsOn
1220	$reader ReadAllFieldsOn
1221	$reader Update
1222	file delete $tmpfile
1223
1224	set dataset [$reader GetOutput]
1225	set limits {}
1226	foreach {xmin xmax ymin ymax zmin zmax} [$dataset GetBounds] break
1227	# Figure out the dimension of the mesh from the bounds.
1228	set _dim 0
1229	if { $xmax > $xmin } {
1230	incr _dim
1231	}
1232	if { $ymax > $ymin } {
1233	incr _dim
1234	}
1235	if { $zmax > $zmin } {
1236	incr _dim
1237	}
1238	if { $_viewer == "" } {
1239	if { $_dim == 2 } {
1240	set _viewer contour
1241	} else {
1242	set _viewer isosurface
1243	}
1244	}
1245	set _comp2dims($cname) ${_dim}D
1246	lappend limits x [list $xmin $xmax]
1247	lappend limits y [list $ymin $ymax]
1248	lappend limits z [list $zmin $zmax]
1249	set dataAttrs [$dataset GetPointData]
1250	if { $dataAttrs == ""} {
1251	puts stderr "WARNING: No point data found in \"$_path\""
1252	rename $reader ""
1253	return 0
1254	}
1255	set vmin 0
1256	set vmax 1
1257	set numArrays [$dataAttrs GetNumberOfArrays]
1258	if { $numArrays > 0 } {
1259	for {set i 0} {$i < [$dataAttrs GetNumberOfArrays] } {incr i} {
1260	set array [$dataAttrs GetArray $i]
1261	set fname [$dataAttrs GetArrayName $i]
1262	foreach {min max} [$array GetRange -1] break
1263	if {$i == 0} {
1264	set vmin $min
1265	set vmax $max
1266	}
1267	lappend limits $fname [list $min $max]
1268	set _fld2Units($fname) ""
1269	set _fld2Label($fname) $fname
1270	# Let the VTK file override the <type> designated.
1271	set _fld2Components($fname) [$array GetNumberOfComponents]
1272	lappend _comp2fldName($cname) $fname
1273	}
1274	}
1275
1276	lappend limits v [list $vmin $vmax]
1277	set _comp2limits($cname) $limits
1278
1279	rename $reader ""
1280	}
1281
1282	#
1283	# VtkDataSetToXy --
1284	#
1285	# Attempt to convert 0 or 1 dimensional VTK DataSet to XY data (curve).
1286	#
1287	itcl::body Rappture::Field::VtkDataSetToXy { dataset } {
1288	foreach {xmin xmax ymin ymax zmin zmax} [$dataset GetBounds] break
1289	# Only X can have non-zero range. X can have zero range if there is
1290	# only one point
1291	if { $ymax > $ymin } {
1292	return 0
1293	}
1294	if { $zmax > $zmin } {
1295	return 0
1296	}
1297	set numPoints [$dataset GetNumberOfPoints]
1298	set dataAttrs [$dataset GetPointData]
1299	if { $dataAttrs == ""} {
1300	puts stderr "WARNING: No point data found"
1301	return 0
1302	}
1303	set array [$dataAttrs GetScalars]
1304	if { $array == ""} {
1305	puts stderr "WARNING: No scalar point data found"
1306	return 0
1307	}
1308	# Multi-component scalars (e.g. color scalars) are not supported
1309	if { [$array GetNumberOfComponents] != 1 } {
1310	return 0
1311	}
1312	set numTuples [$array GetNumberOfTuples]
1313	set xv [blt::vector create \#auto]
1314	for { set i 0 } { $i < $numPoints } { incr i } {
1315	set point [$dataset GetPoint $i]
1316	$xv append [lindex $point 0]
1317	}
1318	set yv [blt::vector create \#auto]
1319	for { set i 0 } { $i < $numTuples } { incr i } {
1320	$yv append [$array GetComponent $i 0]
1321	}
1322	$xv sort $yv
1323	set _comp2xy($cname) [list $xv $yv]
1324	return 1
1325	}
1326
1327	#
1328	# vtkdata --
1329	#
1330	# Returns a string representing the mesh and field data for a specific
1331	# component in the legacy VTK file format.
1332	#
1333	itcl::body Rappture::Field::vtkdata {cname} {
1334	if {$cname == "component0"} {
1335	set cname "component"
1336	}
1337	# DX: Convert DX to VTK
1338	if {[info exists _comp2dx($cname)]} {
1339	set data $_comp2dx($cname)
1340	set data [Rappture::encoding::decode $data]
1341	return [Rappture::DxToVtk $data]
1342	}
1343	# VTK file data:
1344	if { [info exists _comp2vtk($cname)] } {
1345	return $_comp2vtk($cname)
1346	}
1347	# Points on mesh: Construct VTK file output.
1348	if { [info exists _comp2mesh($cname)] } {
1349	# Data is in the form mesh and vector
1350	foreach {mesh vector} $_comp2mesh($cname) break
1351	set label $cname
1352	regsub -all { } $label {_} label
1353	append out "# vtk DataFile Version 3.0\n"
1354	append out "[hints label]\n"
1355	append out "ASCII\n"
1356	append out [$mesh vtkdata]
1357
1358	if { $_comp2assoc($cname) == "pointdata" } {
1359	set vtkassoc "POINT_DATA"
1360	} elseif { $_comp2assoc($cname) == "celldata" } {
1361	set vtkassoc "CELL_DATA"
1362	} elseif { $_comp2assoc($cname) == "fielddata" } {
1363	set vtkassoc "FIELD"
1364	} else {
1365	error "unknown association \"$_comp2assoc($cname)\""
1366	}
1367	set elemSize [numComponents $cname]
1368	set numValues [expr [$vector length] / $elemSize]
1369	if { $_comp2assoc($cname) == "fielddata" } {
1370	append out "$vtkassoc FieldData 1\n"
1371	append out "$label $elemSize $numValues double\n"
1372	} else {
1373	append out "$vtkassoc $numValues\n"
1374	if { $_comp2type($cname) == "colorscalars" } {
1375	# Must be float for ASCII, unsigned char for BINARY
1376	append out "COLOR_SCALARS $label $elemSize\n"
1377	} elseif { $_comp2type($cname) == "normals" } {
1378	# elemSize must equal 3
1379	append out "NORMALS $label double\n"
1380	} elseif { $_comp2type($cname) == "scalars" } {
1381	# elemSize can be 1, 2, 3 or 4
1382	append out "SCALARS $label double $elemSize\n"
1383	append out "LOOKUP_TABLE default\n"
1384	} elseif { $_comp2type($cname) == "tcoords" } {
1385	# elemSize must be 1, 2, or 3
1386	append out "TEXTURE_COORDINATES $label $elemSize double\n"
1387	} elseif { $_comp2type($cname) == "tensors" } {
1388	# elemSize must equal 9
1389	append out "TENSORS $label double\n"
1390	} elseif { $_comp2type($cname) == "vectors" } {
1391	# elemSize must equal 3
1392	append out "VECTORS $label double\n"
1393	} else {
1394	error "unknown element type \"$_comp2type($cname)\""
1395	}
1396	}
1397	append out [$vector range 0 end]
1398	append out "\n"
1399	if 0 {
1400	VerifyVtkDataSet $out
1401	}
1402	return $out
1403	}
1404	error "can't find vtkdata for $cname"
1405	}
1406
1407	#
1408	# BuildPointsOnMesh --
1409	#
1410	# Parses the field XML description to build a mesh and values vector
1411	# representing the field. Right now we handle the deprecated types
1412	# of "cloud", "unirect2d", and "unirect3d" (mostly for flows).
1413	#
1414	itcl::body Rappture::Field::BuildPointsOnMesh {cname} {
1415	#
1416	# More complex 2D/3D data is represented by a mesh
1417	# object and an associated vector for field values.
1418	#
1419	set path [$_field get $cname.mesh]
1420	if {[$_xmlobj element $path] == ""} {
1421	# Unknown mesh designated.
1422	puts stderr "ERROR: Unknown mesh \"$path\""
1423	return 0
1424	}
1425	set viewer [$_field get "about.view"]
1426	if { $viewer != "" } {
1427	set _viewer $viewer
1428	}
1429	set element [$_xmlobj element -as type $path]
1430	set name $cname
1431	regsub -all { } $name {_} name
1432	set _fld2Label($name) $name
1433	set label [hints zlabel]
1434	if { $label != "" } {
1435	set _fld2Label($name) $label
1436	}
1437	set _fld2Units($name) [hints zunits]
1438	set _fld2Components($name) $_comp2size($cname)
1439	lappend _comp2fldName($cname) $name
1440
1441	# Handle bizarre cases that hopefully will be deprecated.
1442	if { $element == "unirect3d" } {
1443	# Special case: unirect3d (should be deprecated) + flow.
1444	if { [$_field element $cname.extents] != "" } {
1445	set vectorsize [$_field get $cname.extents]
1446	} else {
1447	set vectorsize 1
1448	}
1449	set _type unirect3d
1450	set _dim 3
1451	if { $_viewer == "" } {
1452	set _viewer flowvis
1453	}
1454	set _comp2dims($cname) "3D"
1455	set _comp2unirect3d($cname) \
1456	[Rappture::Unirect3d \#auto $_xmlobj $_field $cname $vectorsize]
1457	set _comp2style($cname) [$_field get $cname.style]
1458	set limits {}
1459	foreach axis { x y z } {
1460	lappend limits $axis [$_comp2unirect3d($cname) limits $axis]
1461	}
1462	# Get the data limits
1463	set vector [$_comp2unirect3d($cname) values]
1464	set minmax [VectorLimits $vector $vectorsize]
1465	lappend limits $cname $minmax
1466	lappend limits v $minmax
1467	set _comp2limits($cname) $limits
1468	if {[$_field element $cname.flow] != ""} {
1469	set _comp2flowhints($cname) \
1470	[Rappture::FlowHints ::\#auto $_field $cname $_units]
1471	}
1472	incr _counter
1473	return 1
1474	}
1475	if { $element == "unirect2d" && [$_field element $cname.flow] != "" } {
1476	# Special case: unirect2d (normally deprecated) + flow.
1477	if { [$_field element $cname.extents] != "" } {
1478	set vectorsize [$_field get $cname.extents]
1479	} else {
1480	set vectorsize 1
1481	}
1482	set _type unirect2d
1483	set _dim 2
1484	if { $_viewer == "" } {
1485	set _viewer "flowvis"
1486	}
1487	set _comp2dims($cname) "2D"
1488	set _comp2unirect2d($cname) \
1489	[Rappture::Unirect2d \#auto $_xmlobj $path]
1490	set _comp2style($cname) [$_field get $cname.style]
1491	set _comp2flowhints($cname) \
1492	[Rappture::FlowHints ::\#auto $_field $cname $_units]
1493	set _values [$_field get $cname.values]
1494	set limits {}
1495	foreach axis { x y z } {
1496	lappend limits $axis [$_comp2unirect2d($cname) limits $axis]
1497	}
1498	set xv [blt::vector create \#auto]
1499	$xv set $_values
1500	set minmax [VectorLimits $xv $vectorsize]
1501	lappend limits $cname $minmax
1502	lappend limits v $minmax
1503	blt::vector destroy $xv
1504	set _comp2limits($cname) $limits
1505	incr _counter
1506	return 1
1507	}
1508	switch -- $element {
1509	"cloud" {
1510	set mesh [Rappture::Cloud::fetch $_xmlobj $path]
1511	set _type cloud
1512	}
1513	"mesh" {
1514	set mesh [Rappture::Mesh::fetch $_xmlobj $path]
1515	set _type mesh
1516	}
1517	"unirect2d" {
1518	if { $_viewer == "" } {
1519	set _viewer "heightmap"
1520	}
1521	set mesh [Rappture::Unirect2d::fetch $_xmlobj $path]
1522	set _type unirect2d
1523	}
1524	}
1525	if { ![$mesh isvalid] } {
1526	puts stderr "Mesh is invalid"
1527	return 0
1528	}
1529	set _dim [$mesh dimensions]
1530	if { $_dim == 3 } {
1531	set dim 0
1532	foreach axis {x y z} {
1533	foreach {min max} [$mesh limits $axis] {
1534	if { $min < $max } {
1535	incr dim
1536	}
1537	}
1538	}
1539	if { $dim != 3 } {
1540	set _dim $dim
1541	}
1542	}
1543
1544	if {$_dim < 2} {
1545	puts stderr "ERROR: Can't convert 1D cloud/mesh to curve. Please use curve output for 1D meshes."
1546	return 0
1547
1548	# 1D data: Create vectors for graph widget.
1549	# The prophet tool currently outputs 1D clouds with fields
1550	# Band Structure Lab used to (see isosurface1 test in rappture-bat)
1551	#
1552	# Is there a natural growth path in generating output from 1D to
1553	# higher dimensions? If there isn't, let's kill this in favor
1554	# or explicitly using a <curve> instead. Otherwise, the features
1555	# (methods such as xmarkers) of the <curve> need to be added
1556	# to the <field>.
1557	#
1558	#set xv [blt::vector create x$_counter]
1559	#set yv [blt::vector create y$_counter]
1560
1561	# This only works with a Cloud mesh type, since the points method
1562	# is not implemented for the Mesh object
1563	#$xv set [$mesh points]
1564	# TODO: Put field values in yv
1565	#set _comp2dims($cname) "1D"
1566	#set _comp2xy($cname) [list $xv $yv]
1567	#incr _counter
1568	#return 1
1569	}
1570	if {$_dim == 2} {
1571	# 2D data: By default surface or contour plot using heightmap widget.
1572	set v [blt::vector create \#auto]
1573	$v set [$_field get $cname.values]
1574	if { [$v length] == 0 } {
1575	return 0
1576	}
1577	if { $_viewer == "" } {
1578	set _viewer "contour"
1579	}
1580	set numFieldValues [$v length]
1581	set numComponentsPerTuple [numComponents $cname]
1582	if { [expr $numFieldValues % $numComponentsPerTuple] != 0 } {
1583	puts stderr "ERROR: Number of field values ($numFieldValues) not divisble by elemsize ($numComponentsPerTuple)"
1584	return 0
1585	}
1586	set numFieldTuples [expr $numFieldValues / $numComponentsPerTuple]
1587	if { $_comp2assoc($cname) == "pointdata" } {
1588	set numPoints [$mesh numpoints]
1589	if { $numPoints != $numFieldTuples } {
1590	puts stderr "ERROR: Number of points in mesh ($numPoints) and number of field tuples ($numFieldTuples) don't agree"
1591	return 0
1592	}
1593	} elseif { $_comp2assoc($cname) == "celldata" } {
1594	set numCells [$mesh numcells]
1595	if { $numCells != $numFieldTuples } {
1596	puts stderr "ERROR: Number of cells in mesh ($numCells) and number of field tuples ($numFieldTuples) don't agree"
1597	return 0
1598	}
1599	}
1600	set _comp2dims($cname) "[$mesh dimensions]D"
1601	set _comp2mesh($cname) [list $mesh $v]
1602	set _comp2style($cname) [$_field get $cname.style]
1603	if {[$_field element $cname.flow] != ""} {
1604	set _comp2flowhints($cname) \
1605	[Rappture::FlowHints ::\#auto $_field $cname $_units]
1606	}
1607	incr _counter
1608	array unset _comp2limits $cname
1609	foreach axis { x y z } {
1610	lappend _comp2limits($cname) $axis [$mesh limits $axis]
1611	}
1612	set minmax [VectorLimits $v $_comp2size($cname)]
1613	lappend _comp2limits($cname) $cname $minmax
1614	lappend _comp2limits($cname) v $minmax
1615	return 1
1616	}
1617	if {$_dim == 3} {
1618	# 3D data: By default isosurfaces plot using isosurface widget.
1619	if { $_viewer == "" } {
1620	set _viewer "isosurface"
1621	}
1622	set v [blt::vector create \#auto]
1623	$v set [$_field get $cname.values]
1624	if { [$v length] == 0 } {
1625	puts stderr "ERROR: No field values"
1626	return 0
1627	}
1628	set numFieldValues [$v length]
1629	set numComponentsPerTuple [numComponents $cname]
1630	if { [expr $numFieldValues % $numComponentsPerTuple] != 0 } {
1631	puts stderr "ERROR: Number of field values ($numFieldValues) not divisble by elemsize ($numComponentsPerTuple)"
1632	return 0
1633	}
1634	set numFieldTuples [expr $numFieldValues / $numComponentsPerTuple]
1635	if { $_comp2assoc($cname) == "pointdata" } {
1636	set numPoints [$mesh numpoints]
1637	if { $numPoints != $numFieldTuples } {
1638	puts stderr "ERROR: Number of points in mesh ($numPoints) and number of field tuples ($numFieldTuples) don't agree"
1639	return 0
1640	}
1641	} elseif { $_comp2assoc($cname) == "celldata" } {
1642	set numCells [$mesh numcells]
1643	if { $numCells != $numFieldTuples } {
1644	puts stderr "ERROR: Number of cells in mesh ($numCells) and number of field tuples ($numFieldTuples) don't agree"
1645	return 0
1646	}
1647	}
1648	set _comp2dims($cname) "[$mesh dimensions]D"
1649	set _comp2mesh($cname) [list $mesh $v]
1650	set _comp2style($cname) [$_field get $cname.style]
1651	if {[$_field element $cname.flow] != ""} {
1652	set _comp2flowhints($cname) \
1653	[Rappture::FlowHints ::\#auto $_field $cname $_units]
1654	}
1655	incr _counter
1656	foreach axis { x y z } {
1657	lappend _comp2limits($cname) $axis [$mesh limits $axis]
1658	}
1659	set minmax [VectorLimits $v $_comp2size($cname)]
1660	lappend _comp2limits($cname) $cname $minmax
1661	lappend _comp2limits($cname) v $minmax
1662	return 1
1663	}
1664	error "unhandled case in field dim=$_dim element=$element"
1665	}
1666
1667	itcl::body Rappture::Field::AvsToVtk { cname contents } {
1668	package require vtk
1669
1670	set reader $this-datasetreader
1671	vtkAVSucdReader $reader
1672
1673	# Write the contents to a file just in case it's binary.
1674	set tmpfile $cname[pid].ucd
1675	set f [open "$tmpfile" "w"]
1676	fconfigure $f -translation binary -encoding binary
1677	puts $f $contents
1678	close $f
1679	$reader SetFileName $tmpfile
1680	$reader Update
1681	file delete $tmpfile
1682
1683	set tmpfile $cname[pid].vtk
1684	set writer $this-datasetwriter
1685	vtkDataSetWriter $writer
1686	$writer SetInputConnection [$reader GetOutputPort]
1687	$writer SetFileName $tmpfile
1688	$writer SetFileTypeToBinary
1689	$writer Write
1690	rename $reader ""
1691	rename $writer ""
1692
1693	set f [open "$tmpfile" "r"]
1694	fconfigure $f -translation binary -encoding binary
1695	set vtkdata [read $f]
1696	close $f
1697	file delete $tmpfile
1698	return $vtkdata
1699	}
1700
1701	itcl::body Rappture::Field::DicomToVtk { cname path } {
1702	package require vtk
1703
1704	if { ![file exists $path] } {
1705	puts stderr "path \"$path\" doesn't exist."
1706	return 0
1707	}
1708
1709	if { [file isdir $path] } {
1710	set files [glob -nocomplain $path/*.dcm]
1711	if { [llength $files] == 0 } {
1712	set files [glob -nocomplain $path/*]
1713	if { [llength $files] == 0 } {
1714	puts stderr "no dicom files found in \"$path\""
1715	return 0
1716	}
1717	}
1718
1719	#array set data [Rappture::DicomToVtk files $files]
1720	array set data [Rappture::DicomToVtk dir $path]
1721	} else {
1722	array set data [Rappture::DicomToVtk files [list $path]]
1723	}
1724
1725	if 0 {
1726	foreach key [array names data] {
1727	if {$key == "vtkdata"} {
1728	if 0 {
1729	set f [open /tmp/$cname.vtk "w"]
1730	fconfigure $f -translation binary -encoding binary
1731	puts -nonewline $f $data(vtkdata)
1732	close $f
1733	}
1734	} else {
1735	puts stderr "$key = \"$data($key)\""
1736	}
1737	}
1738	}
1739
1740	ReadVtkDataSet $cname $data(vtkdata)
1741	return $data(vtkdata)
1742	}
1743
1744	itcl::body Rappture::Field::GetTypeAndSize { cname } {
1745	array set type2components {
1746	"colorscalars" 4
1747	"normals" 3
1748	"scalars" 1
1749	"tcoords" 2
1750	"tensors" 9
1751	"vectors" 3
1752	}
1753	set type [$_field get $cname.elemtype]
1754	if { $type == "" } {
1755	set type "scalars"
1756	}
1757	if { ![info exists type2components($type)] } {
1758	error "unknown <elemtype> \"$type\" in field"
1759	}
1760	set size [$_field get $cname.elemsize]
1761	if { $size == "" } {
1762	set size $type2components($type)
1763	}
1764	set _comp2type($cname) $type
1765	set _comp2size($cname) $size
1766	}
1767
1768	itcl::body Rappture::Field::GetAssociation { cname } {
1769	set assoc [$_field get $cname.association]
1770	if { $assoc == "" } {
1771	set _comp2assoc($cname) "pointdata"
1772	return
1773	}
1774	switch -- $assoc {
1775	"pointdata" - "celldata" - "fielddata" {
1776	set _comp2assoc($cname) $assoc
1777	return
1778	}
1779	default {
1780	error "unknown field association \"$assoc\""
1781	}
1782	}
1783	}
1784
1785	#
1786	# Compute the per-component limits or limits of vector magnitudes
1787	#
1788	itcl::body Rappture::Field::VectorLimits {vector vectorsize {comp -1}} {
1789	if {$vectorsize == 1} {
1790	set minmax [$vector limits]
1791	} else {
1792	set len [$vector length]
1793	if {[expr $len % $vectorsize] != 0} {
1794	error "Invalid vectorsize: $vectorsize"
1795	}
1796	if {$comp > $vectorsize-1} {
1797	error "Invalid vector component: $comp"
1798	}
1799	set numTuples [expr ($len/$vectorsize)]
1800	for {set i 0} {$i < $numTuples} {incr i} {
1801	if {$comp >= 0} {
1802	set idx [expr ($i * $vectorsize + $comp)]
1803	set val [$vector index $idx]
1804	} else {
1805	set idx [expr ($i * $vectorsize)]
1806	set mag 0
1807	for {set j 0} {$j < $vectorsize} {incr j} {
1808	set val [$vector index $idx]
1809	set mag [expr ($mag + $val * $val)]
1810	incr idx
1811	}
1812	set val [expr (sqrt($mag))]
1813	}
1814	if (![info exists minmax]) {
1815	set minmax [list $val $val]
1816	} else {
1817	if {$val < [lindex $minmax 0]} {
1818	lset minmax 0 $val
1819	}
1820	if {$val > [lindex $minmax 1]} {
1821	lset minmax 1 $val
1822	}
1823	}
1824	}
1825	}
1826	return $minmax
1827	}

Note: See TracBrowser for help on using the repository browser.

Download in other formats: