Context Navigation

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

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

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

Download in other formats: