--- a
+++ b/maxima-pre59/bin/omplotdata
@@ -0,0 +1,4649 @@
+#!/bin/sh
+# comment \
+exec wish8.0 "$0" "$@"
+#############################################
+##### Copyright William Schelter 1997 #######
+#############################################
+set ws_openMath(date) 11/08/1998
+
+###### plotting.tcl ######
+## source plotconf.tcl
+
+###### plotconf.tcl ######
+
+## source private.tcl
+
+###### private.tcl ######
+
+# a private way of storing variables on a window by window
+# basis
+
+proc makeLocal { win args } {
+  foreach v $args {
+     uplevel 1 set  $v \[oget $win $v\]
+ }
+}
+
+proc linkLocal { win args } {
+  foreach v $args {
+      uplevel 1 upvar #0 _WinInfo${win}\($v) $v
+ }
+}
+
+proc clearLocal { win } {
+    global _WinInfo$win
+       # puts "clearing info for $win in [info level 1]"
+
+    catch { unset _WinInfo$win }
+}
+ 
+
+proc oset { win var val } {
+  global _WinInfo$win
+  set _WinInfo[set win]($var) $val
+}
+
+proc oarraySet { win vals } {
+  global _WinInfo$win
+  array set  _WinInfo$win $vals
+}
+
+proc oloc { win var } {
+  return _WinInfo[set win]($var)
+}
+
+proc oarray { win  } {
+  return _WinInfo[set win]
+}
+
+proc oget { win var } {
+  global _WinInfo$win
+  return [set _WinInfo[set win]($var)]
+}
+
+
+## endsource private.tcl
+## source parse.tcl
+
+###### parse.tcl ######
+
+## source getopt.tcl
+
+###### getopt.tcl ######
+## source macros.tcl
+
+###### macros.tcl ######
+proc desetq {lis lis2} {
+    set i 0
+    foreach v $lis {
+	uplevel 1 set $v [list [lindex $lis2 $i]]
+	set i [expr {$i + 1}]
+    }
+}
+
+proc assoc { key lis args } {
+    foreach { k val } $lis {
+	if { "$k" == "$key" } {
+	    return $val }
+    }
+    return [lindex $args 0]
+}
+
+proc delassoc { key lis } {
+    foreach { k val } $lis {
+	if { "$k" != "$key" } {
+	lappend new $k $val
+	}
+    }
+    return $new
+}
+
+
+proc putassoc {key lis value } {
+    set done 0
+    foreach { k val } $lis {
+	if { "$k" == "$key" } {
+	    set done 1
+	    set val $value
+	}
+	lappend new $k $val
+    }
+    if { !$done } {
+	lappend new $key $value
+    }
+    return $new
+}
+
+proc intersect { lis1 lis2 } {
+    set new ""
+    foreach v $lis1 { set there($v) 1 }
+    foreach v $lis2 { if { [info exists there($v)] } { lappend new $v }}
+    return $new
+}
+
+
+
+#
+ #-----------------------------------------------------------------
+ #
+ # ldelete --  remove all copies of ITEM from LIST
+ #
+ #  Results: new list without item
+ #
+ #  Side Effects: 
+ #
+ #----------------------------------------------------------------
+#
+proc ldelete { item list } {
+    while { [set ind [lsearch $list $item]] >= 0  } {
+	set list [concat [lrange $list 0 [expr {$ind -1}]] [lrange $list [expr {$ind +1}] end]]
+    }
+    return $list
+}
+
+# apply f a1 a2 a3 [list  u1 u2 ..un]   , should call
+# f with n+3 arguments.
+proc apply {f args } {
+    set lis1 [lrange $args 0 [expr {[llength $args] -2}]]
+    foreach v [lindex $args end] { lappend lis1 $v}
+    set lis1 [linsert $lis1  0 $f]
+    uplevel 1 $lis1
+}
+
+
+
+
+
+
+## endsource macros.tcl
+
+#####sample option list.   Error will be signalled if "Required" option
+##### not given.
+#set dfplotOptions {
+#    {xdot Required {specifies dx/dt = xdot.  eg -xdot "x+y+sin(x)^2"} }
+#    {ydot Required {specifies dy/dt = ydot.  eg -ydot "x-y^2+exp(x)"} }
+#    {xradius 10 "Width in x direction of the x values" }
+#    {yradius 10 "Height in y direction of the y values"}
+#}
+
+
+#
+ #-----------------------------------------------------------------
+ #
+ # optLoc --  if $usearray is not 0, then the OPTION is stored 
+ # in a hashtable, otherwise in the variable whose name is the
+ # same as OPTION.  
+ #  Results: a form which when 'set' will allow storing value.
+ #
+ #  Side Effects: none
+ #
+ #----------------------------------------------------------------
+#
+proc optLoc { op ar }  {
+  #  puts "$ar,[lindex $op 0]"
+  # puts "return=$ar\([lindex $op 0]\)"
+   if { "$ar" == 0 } {
+       return [lindex $op 0]
+   } else {
+     #puts "$ar\([lindex $op 0]\)"
+    return "$ar\([lindex $op 0]\)" }
+
+}
+
+
+#
+ #-----------------------------------------------------------------
+ #
+ # getOptions --  given OPTLIST a specification for the options taken,
+ # parse the alternating keyword1 value1 keyword2 value2 options_supplied
+ # to make sure they are allowed, and not just typos, and to supply defaults
+ # for ones not given.   Give an error message listing options.
+ # a specification is  { varname default_value "doc string" }
+ # and optlist, is a list of these.   the key should be -varname
+ # 
+ #  -debug 1 "means print the values on standard out"
+ #  -allowOtherKeys 1 "dont signal an error if -option is supplied but not in
+ #                     the list"
+ #  -usearray "should give a NAME, so that options are stored in NAME(OPTION)
+ #  -setdefaults "if not 0 (default is 1) do `set OPTION dflt' for all options"
+ # If a key is specified twice eg.  -key1 val1 -key1 val2, then the first
+ # value val1 will be used 
+ #  Results:
+ #
+ #  Side Effects: set the values in the callers environment
+ #
+ #----------------------------------------------------------------
+#
+
+proc getOptions { optlist options_supplied args } {
+   # global  getOptionSpecs
+    set ar [assoc -usearray $args 0]
+    set help [assoc -help $args ""]
+    if { "$ar" != "0" } { global $ar }
+    set debug [assoc -debug $args 0]
+    set allowOtherKeys [assoc -allowOtherKeys $args 0]
+    set setdefaults [assoc -setdefaults $args 1]
+    set supplied ""
+
+    foreach {key val } $options_supplied {
+	if { [info exists already($key)] } { continue }
+	set already($key) 1
+	set found 0
+	foreach op $optlist {
+	    if { "$key" == "-[lindex $op 0]" } {
+		uplevel 1 set [optLoc $op $ar] [list $val]
+		
+		append supplied " [lindex $op 0]"
+		set found 1
+		break
+	    }
+	}
+	set caller global
+
+	if { $found == 0  && !$allowOtherKeys } {
+	    catch {set caller [lindex [info level -1] 0]}
+             error "`$caller' does not take the key `$key':\n[optionHelpMessage $optlist]\n$help"
+
+	    }
+	}
+    foreach op $optlist {
+ 	if { [lsearch $supplied [lindex $op 0]] < 0 } {
+           if { "[lindex $op 1]" == "Required" } {
+	       catch {set caller [lindex [info level -1] 0]}	       
+              error "`-[lindex $op 0]' is required option for `$caller':\n[optionHelpMessage $optlist]"
+            }
+       	   if { $setdefaults }  {
+            
+	    uplevel 1 set [optLoc $op $ar] [list [lindex $op 1]]
+          }
+	}
+	# for debugging see them.
+	# if { $debug } {   uplevel 1 puts "[optLoc $op $ar]=\$[optLoc $op $ar]"}
+	 if { $debug } {   puts "[optLoc $op $ar]=[safeValue [optLoc $op $ar] 2]"}
+	
+    }
+}
+
+proc getOptionDefault { key optionList } {
+ foreach v $optionList {
+  if { "[lindex $v 0]" == "$key" } { return [lindex $v 1]}
+   }
+ return ""
+}
+
+proc assq {key list {dflt ""}} {
+    foreach v $list { if { "[lindex $v 0]" == "$key" } { return $v }}
+    return $dflt
+}
+
+proc safeValue { loc level} {
+  if { ![catch { set me [uplevel $level set $loc] } ] } {
+     return $me  }  else {return "`unset'" }
+}
+     
+  
+
+proc optionFirstItems { lis } {
+    set ans ""
+    foreach v $lis { append ans " [list [lindex $v 0]]" }
+    return $ans
+}
+
+proc optionHelpMessage { optlist } {
+    set msg ""
+    foreach op $optlist  {	
+	append msg \
+		" -[lindex $op 0] \[ [lindex $op 1] \] --[lindex $op 2]\n"
+	    }
+    return $msg
+	}
+	    
+
+
+## endsource getopt.tcl
+
+catch { unset Parser }
+
+foreach v  { { ( 120 } { \[ 120 } { ) 120 } { \] 120 }  { ^ 110}
+         {* 100} { / 100} {% 100}  {- 90 } { + 90 }
+           { << 80} { >> 80 } { < 70 } { > 70 } { <= 70 } {>= 70}
+	   { == 60 } { & 50} { | 40 } { , 40 } {= 40}
+	   { && 30 } { || 20 } { ? 10 } { : 10 }  { ; 5 }}  {
+	       set parse_table([lindex $v 0]) [lindex $v 1]
+	       set getOp([lindex $v 0]) doBinary
+	       
+	   }
+
+proc binding_power {s} {
+    global parse_table billy
+    set billy $s
+    if { [catch { set tem $parse_table($s) }] } { return 0 } else { return $tem }
+}
+
+proc getOneMatch { s inds } {
+    return [string range $s [lindex $inds 0] [lindex $inds 1]]
+}
+proc parseTokenize { str } {
+    regsub  -all {[*][*]} $str "^" str
+    set ans ""
+    while { [string length $str ] > 0 } {
+#    puts "ans=$ans,str=$str"	
+    set str [string trimleft $str " \t\n" ]
+    set s [string range $str 0 1]
+    set bp [binding_power $s]
+    if { $bp > 0 } { append ans " $s"
+	   set str [string range $str 2 end]
+	continue
+    } else {
+	set s [string range $s 0 0]
+        set bp [binding_power $s]
+        if { $bp > 0 } { append ans " $s"
+	   set str [string range $str 1 end]
+	continue
+   }
+  }
+  if { "$s" == "" } {
+      return $ans
+  }
+  if { [regexp -indices {^[0-9.]+[eE]*[0-9]*} $str all] } {
+      append ans " { number [getOneMatch $str $all] }"
+     # append ans " [getOneMatch $str $all]"
+      set str [string range $str [expr {1+ [lindex $all 1]}] end]
+  }  elseif { [regexp -indices {^[$a-zA-Z][a-zA-Z0-9]*} $str all] } {
+       append ans " { id [getOneMatch $str $all] } "
+      # append ans " [getOneMatch $str $all]"
+      set str [string range $str [expr {1+ [lindex $all 1]}] end]
+  }  else { error "parser unrecognized: $str"
+  }
+  }
+  return $ans
+}
+
+set Parser(reserved) " acos cos hypo sinh asin cosh log sqrt atan exp log10 tan atan2 floor pow tanh ceil fmod sin abs double int round"
+
+set Parser(help) [join [list {
+The syntax is like C except that it is permitted to write x^n
+instead of pow(x,n).
+} "\nFunctions: $Parser(reserved)\n\nOperators: == % & || ( << <= ) : * >=  + && , | < >> - > ^ ? /" ] ""]
+
+    
+
+proc nexttok { } {
+    global Parser
+    set x [lindex $Parser(tokenlist) [incr Parser(tokenind) ]]
+    # puts "nexttok=$x"
+    if {[llength $x ] > 1 } {
+	set Parser(tokenval) [lindex $x 1]
+	return [lindex $x 0]
+    } else { return $x }
+}
+
+
+#
+ #-----------------------------------------------------------------
+ #
+ # parseToSuffixLists -- Convert EXPR1; EXPR2; ..
+ # to a list of suffix lists.  Each suffix list is suitable for
+ # evaluating on a stack machine (like postscript) or for converting
+ # further into another form.  see parseFromSuffixList.
+ #  "1+2-3^4;" ==>
+ #   {number 1} {number 2} + {number 3} {number 4} ^ -
+ #  Results: suffix list form of the original EXPR
+ #
+ #  Side Effects: none
+ #
+ #----------------------------------------------------------------
+#
+proc parseToSuffixLists { a }  {
+    global    Parser 
+    set Parser(result) ""
+    set Parser(tokenlist) [parseTokenize $a]
+    set Parser(tokenind) -1
+    set Parser(lookahead)  [nexttok]
+    set ans ""
+    while { "$Parser(lookahead)" != ""  } {
+      getExpr  ; parseMatch ";"
+      #puts "here: $Parser(result) "	
+      append ans "[list	$Parser(result)] "
+      set Parser(result) "" 	
+    }
+    return $ans
+}
+
+proc parseMatch { t } {
+    global Parser
+    if { "$t" == "$Parser(lookahead)" } {
+	set Parser(lookahead)  [nexttok]
+    } else { error "syntax error: wanted $t"}
+}
+
+proc emit { s args } {
+    global Parser
+    if { "$args" == "" }   {
+	append Parser(result) " $s"
+	# puts " $s "
+    } else {
+	append Parser(result) " {[lindex $args 0 ] $s}"
+	#puts " {[lindex $args 0 ] $s} "
+    }
+}
+
+proc getExpr { } { getExprn 0 }
+
+proc getExprn { n } {
+    global Parser
+    if { $n == 110 } {
+      getExpr120
+      return
+     }
+
+    incr n 10
+    if  { $n == 110 } {
+       if { "$Parser(lookahead)" == "-" || "$Parser(lookahead)" == "+"  } {
+            if { "$Parser(lookahead)" == "-" } { set this PRE_MINUS } else {
+       	    set this PRE_PLUS }
+	    parseMatch $Parser(lookahead)
+            getExprn $n
+	    #puts "l=$Parser(lookahead),pl=$Parser(result)"
+            emit $this
+            return
+           } 
+	       
+    }
+
+    getExprn $n
+    while { 1 } {
+	if { [binding_power $Parser(lookahead)] == $n } {
+	    set this $Parser(lookahead)
+	    parseMatch $Parser(lookahead)
+            getExprn $n
+	    if { $n == 110 } {
+		set toemit ""
+		while { "$this" == "^" &&  "$Parser(lookahead)" == "^" } {
+		    # puts "p=$Parser(result),$
+		    set this $Parser(lookahead)
+		    append toemit " $this"
+		    parseMatch $Parser(lookahead)
+		    getExprn $n
+		}
+		foreach v $toemit { emit $v }
+	    }
+	    emit $this
+           
+	} else { return }
+    }
+}
+
+proc getExpr120 { } {
+    global Parser
+    while { 1 } {
+	if { "$Parser(lookahead)" == "(" } {
+	    parseMatch $Parser(lookahead)
+	    getExpr
+	    parseMatch ")"
+	    break;
+	} elseif { $Parser(lookahead) == "id" } {
+	    emit $Parser(tokenval) id
+
+	    parseMatch $Parser(lookahead)
+	    if { "$Parser(lookahead)" == "(" } {
+		getExpr120
+		emit funcall
+	    }
+	    break;
+	} elseif { $Parser(lookahead) == "number" } {
+	    emit $Parser(tokenval) number
+	    parseMatch $Parser(lookahead)
+	    break;
+	} else { error "syntax error" }
+    }
+}
+
+set getOp(PRE_PLUS) doPrefix
+set getOp(PRE_MINUS) doPrefix
+set getOp(funcall) doFuncall
+set getOp(^) doPower
+set getOp(:) doConditional
+set getOp(?) doConditional
+
+proc doBinary { } {
+    uplevel 1 {set s $nargs; incr nargs -1 ;
+    if { "$x" == "," } {    set a($nargs) "$a($nargs) $x $a($s)"} else { 
+
+	set a($nargs) "($a($nargs) $x $a($s))"} }
+}
+
+proc doPower { } {
+    uplevel 1 {set s $nargs; incr nargs -1 ; set a($nargs) "pow($a($nargs),$a($s))" }
+}
+proc doFuncall {} {
+    uplevel 1 {
+	#puts nargs=$nargs
+	set s $nargs; incr nargs -1 ; set a($nargs) "$a($nargs)($a($s))"}
+}
+
+proc doPrefix {} {
+  uplevel 1  { if { "$x" == "PRE_MINUS" } { set a($nargs) "-$a($nargs)" } }
+}
+
+proc doConditional { } {
+    set x [uplevel 1 set x]
+    if { "$x" == "?" } { return }
+    # must be :
+    uplevel 1 { 
+	set s $nargs ; incr nargs -2 ;
+	set a($nargs) "($a($nargs) ? $a([expr {$nargs + 1}]) : $a($s))"
+ }
+}
+
+
+#
+ #-----------------------------------------------------------------
+ #
+ # parseFromSuffixList --  takes a token list, and turns
+ # it into a suffix form.  eg: 1 + 2 - 3 ^ 4 --> 1 2 + 3 4 ^ -
+ #  Results:
+ #
+ #  Side Effects: 
+ #
+ #----------------------------------------------------------------
+#
+proc parseFromSuffixList { list } {
+    global getOp
+  set stack ""
+  set lim [llength $list]
+  set i 0
+  set nargs 0
+  while { $i < $lim } {
+    set x [lindex $list $i ]
+    set bp [binding_power $x]
+    incr i
+   # all binary
+    if { [llength $x] > 1 } {
+	
+	set a([incr nargs]) [lindex $x 1]
+
+     } else {
+	 $getOp($x) }
+   }
+
+  return $a(1)
+}
+	
+
+#
+ #-----------------------------------------------------------------
+ #
+ # parseConvert --  given an EXPRESSION, parse it and find out
+ # what are the variables, and convert a^b to pow(a,b).  If
+ # -variables "x y" is given, then x and y will be replaced by $x $y
+ #  doall 1 is giv 
+ #  Results:
+ #
+ #  Side Effects: 
+ #
+ #----------------------------------------------------------------
+#
+set Parser(convertOptions) {
+    { doall 0 "convert all variables x to \$x" }
+    { variables "" "list of variables to change from x to \$x" }
+}
+proc parseConvert { expr args } {
+    global   Parser 
+    getOptions $Parser(convertOptions) $args
+    if { "$expr" == "" } { return [list {} {}] }
+    set parselist [parseToSuffixLists "$expr;"]
+    #puts "parselist=$parselist"
+    catch { unset allvars }
+	set new ""
+    set answers ""
+    foreach lis $parselist {
+      
+     foreach v $lis {
+
+	if { ("[lindex $v 0]" == "id")
+	&& ([llength $v] == 2)
+	&& ([lsearch  $Parser(reserved) [set w [lindex $v 1]]] < 0)
+    } {  
+	if { ($doall != 0)  || ([lsearch  $variables $w] >= 0) } {
+	    append new " {id \$$w}"
+	    set allvars(\$$w) 1
+	} else {
+	    set allvars($w) 1
+	    append new " {$v}"
+    }   }  else {
+	if { [llength $v] > 1 } { 
+	    append new " {$v}" } else {
+		append new " $v" }
+	    }
+    }
+    #puts "new=$new"
+    append answers "[list [parseFromSuffixList $new]] "
+    set new ""
+ }
+    return [list $answers [array names allvars]]
+}
+
+proc test { s } {
+    set me [parseFromSuffixList [lindex [parseToSuffixLists "$s;"] 0]]
+    puts $me
+    return "[eval expr $s] [eval expr $me]"
+}
+     
+# 
+# Local Variables:
+# mode: tcl
+# version-control: t
+# End:
+
+
+## endsource parse.tcl
+## source textinsert.tcl
+
+###### textinsert.tcl ######
+
+proc mkTextItem { c x y args  } {
+    set font [assoc -font $args {Helvetica 14}]
+    set tags [assoc -tags $args {}]
+    set item [$c create text $x $y -text " " -width 440 -anchor n -font $font -justify left]
+    append tags text
+    foreach v $tags { $c addtag $v withtag $item}
+    $c bind text <1> "textB1Press $c %x %y"
+    $c bind text <B1-Motion> "textB1Move $c %x %y"
+    $c bind text <Shift-1> "$c select adjust current @%x,%y"
+    $c bind text <Shift-B1-Motion> "textB1Move $c %x %y"
+    $c bind text <KeyPress> "textInsert $c %A"
+    $c bind text <Return> "textInsert $c \\n"
+    $c bind text <Control-h> "textBs $c"
+    $c bind text <BackSpace> "textBs $c"
+    $c bind text <Delete> "textDel $c"
+    $c bind text <2> "textPaste $c @%x,%y" 
+}
+
+
+## endsource textinsert.tcl
+## source printops.tcl
+
+###### printops.tcl ######
+
+### fix a4 size !
+set paperSizes {{letter 8.5 11} { A4 8.5 11} {legal 8.5 13}} 
+
+set printOptions { 
+    { landscape  1 "Non zero means use landscape mode in printing" }
+    { tofile 1 "Non zero means print to file" }
+    { pagewidth "" "Figure width" }
+    { pageheight "" "Figure height" }
+    { papersize letter "letter, legal or A4"}
+    { hoffset .5 "Left margin for printing"}
+    { voffset .5 "Right margin for printing"}
+    { xticks 20 "Rough number of ticks on x axis"}
+    { yticks 20 "Rough number of ticks on y axis"}
+    { domargin 1 "Print the frame and the margin ticks"}
+    { printer "" "Printer to print to, eg lw8b " }
+    { title "" "Title" }
+    { psfilename "~/sdfplot.ps" "Postscript filename" }
+    { gsview "gsview32" "postscript viewer, used for printing under Windows" }
+    { centeronpage 1 ""} 
+}
+
+# proc getPageOffsets { widthbyheight} {
+#     global printOption paperSizes
+#     puts "wbh=$widthbyheight"
+#     set pwid 8.5
+#     set phei 11.0
+
+#     foreach v $paperSizes {
+# 	if { "[lindex $v 0]" == "$printOption(papersize)" } {
+# 	    set pwid [lindex $v 1]
+# 	    set phei [lindex $v 2]
+# 	}
+#     }
+#     set wid [expr {$pwid - 2* $printOption(hoffset)}]
+#     set hei [expr {$phei - 2* $printOption(voffset)}]
+# #    if { $printOption(landscape) } {set widthbyheight [expr  {1.0 /$widthbyheight}]}
+# #    set w $wid ; set hei $wid ; set wid $w
+
+#     puts "pw=$wid,ph=$hei,w/h=$widthbyheight,hh=[expr {$hei * $widthbyheight}], ww=[expr {$wid / $widthbyheight}]"
+
+#     set fac   $widthbyheight
+#     puts "fac=$fac"
+#     if { $fac * $hei < $wid } {
+    # 	set iwid [expr {$fac *$hei}]
+# 	set ihei $hei
+
+#     } else {
+    # 	set ihei [expr {$wid / $fac}]
+
+# 	set iwid $wid
+
+#     }
+
+#     if { $printOption(landscape) } { set fac1 [expr {1/$fac}] }
+#     if { $wid/$hei > $fac } {
+# 	set ihei $hei
+    # 	set iwid   [expr {$hei / $fac }]
+
+#     } else {
+# 	 set iwid $wid
+    # 	 set ihei [expr {$wid * $fac }]
+#     }
+
+#     #-pagex = left margin (whether landscape or not)
+#     #-pagey = right margin (whether landscape or not)
+#     #-pagewidth becomes vertical height if landscape
+#     #-pageheight becomes horiz width if landscape
+    
+#     set xoff [expr {($pwid-$iwid)/2.0}]
+#     set yoff [expr  {($phei-$ihei)/2.0}]
+
+#     if { $printOption(landscape) } {
+# 	set h $ihei
+# 	set ihei $iwid
+# 	set iwid $h
+#     }
+
+#     puts "phei=$phei,ihei=$ihei,yoff=$yoff,voff=$printOption(voffset)"
+#     set ans "-pagex [set xoff]i -pagey [set yoff]i \
+# 	    -pagewidth [set iwid]i -pageheight [set ihei]i"
+#     set ans "-pagex [set xoff]i -pagey [set yoff]i \
+# 	    -pagewidth [set iwid]i -pageheight [set ihei]i"    
+#     return $ans
+# }
+
+proc swap { a b } {
+    set me [uplevel 1 set $b]
+    uplevel 1 set $b \[set $a\]
+    uplevel 1 set $a [list $me]
+}
+
+proc getPageOffsets { widthbyheight} {
+    global printOption paperSizes
+    #puts "wbh=$widthbyheight"
+    set pwid 8.5
+    set phei 11.0
+
+    foreach v $paperSizes {
+	if { "[lindex $v 0]" == "$printOption(papersize)" } {
+	    set pwid [lindex $v 1]
+	    set phei [lindex $v 2]
+	}
+    }
+    set wid [expr {$pwid - 2* $printOption(hoffset)}]
+    set hei [expr {$phei - 2* $printOption(voffset)}]
+    if { $printOption(landscape) } {
+	swap wid hei
+#	swap pwid phei
+    }
+    if { $wid / $hei  < $widthbyheight  } {
+	# width dominates
+	set iwid $wid
+	set ihei [expr {$wid / $widthbyheight }]
+	append opts " -pagewidth [set wid]i"
+    } else {
+	set ihei $hei
+	set iwid [expr {$hei * $widthbyheight }]
+	append opts " -pageheight [set hei]i"
+    }
+
+    #-pagex = left margin (whether landscape or not)
+    #-pagey = right margin (whether landscape or not)
+    #-pagewidth becomes vertical height if landscape
+    #-pageheight becomes horiz width if landscape
+    
+    append opts " -pagex [expr {$pwid / 2.0}]i -pagey [expr {$phei / 2.0}]i "
+
+	if { $printOption(landscape) } {
+	    append opts " -rotate $printOption(landscape)" 
+	}
+    return $opts
+}
+
+set printOption(setupDone) 0
+
+proc getEnv { name } {
+  global env
+ if { [catch { set tem $env($name) } ] } { return "" }
+ return $tem
+}
+proc setPrintOptions { lis } {
+    global browser_version
+   global printOptions printOption printSetUpDone 
+    if { !$printOption(setupDone) } {
+	set printOption(setupDone) 1
+	getOptions $printOptions $lis -allowOtherKeys 1 \
+		-setdefaults [catch { source [getEnv HOME]/.printOptions }] -usearray printOption
+        if { "$printOption(printer)" == "" } {set printOption(printer) [getEnv PRINTER] } else { set printOption(printer) lw8b }
+	
+    }
+    if { [info exists browser_version] } { set printOption(tofile) 2 }
+}
+
+proc mkentryPr { w var text buttonFont }  {
+     set fr $w ; frame $fr
+    uplevel 1 append topack [list " $fr"]
+    label $fr.lab -text "$text" -font $buttonFont
+    entry $fr.e -width 20 -textvariable $var -font $buttonFont
+    pack $fr.lab $fr.e -side left -expand 1 -padx 3 -fill x
+}
+
+
+proc mkPrintDialog { name args } {
+    global printSet argv env printOptions printOption printSetUpDone paperSizes buttonfont
+
+    set canv [assoc -canvas $args ]
+    set buttonFont [assoc -buttonfont $args $buttonfont]
+    catch { destroy $name }
+    set dismiss "destroy $name"
+    if { "$canv" == "" } {
+     catch {destroy $name}
+    toplevel $name
+    wm geometry $name -0+20
+   
+    } else {
+        $canv delete printoptions
+        set name [winfo parent $canv].printoptions
+	# set name $canv.fr1
+        catch {destroy $name}
+	frame $name -borderwidth 2 -relief raised
+	
+	set item [$canv create window [$canv canvasx 10] [$canv canvasy  10] -window $name -anchor nw -tags printoptions]
+        $canv raise printoptions
+	set dismiss "$canv delete $item; destroy $name "
+    }
+	
+    frame $name.fr
+
+    set w $name.fr
+    label $w.msg  -wraplength 600 -justify left -text "Printer Setup"
+    pack $w
+    pack $w.msg
+    set wb $w.buttons
+    frame $wb 
+    pack $wb -side left -fill x -pady 2m
+    set topack ""
+    catch { set printOption(psfilename) \
+	    [file nativename $printOption(psfilename)]}
+    button $wb.ok -text "ok" -font $buttonFont  -command "destroy $name ; $canv delete printoptions"
+    radiobutton $wb.b0 -text "Save via ftp" -variable printOption(tofile) -relief flat -value 2 -command {set writefile "Save"} -font $buttonFont  -highlightthickness 0 
+    radiobutton $wb.b1 -text "Save as Postscript File" -variable printOption(tofile) -relief flat -value 1 -command {set writefile "Save"} -font $buttonFont  -highlightthickness 0 
+    radiobutton $wb.b2 -text "Print To Printer" -variable printOption(tofile) -relief flat -value 0 -command {set writefile "Print"} -font $buttonFont -highlightthickness 0 
+    checkbutton $wb.b3 -text "Center on Page" -variable printOption(centeronpage) -relief flat -font $buttonFont -highlightthickness 0 
+    checkbutton $wb.b4 -text "Landscape Mode" -variable printOption(landscape) -relief flat -font $buttonFont -highlightthickness 0 
+
+    mkentryPr  $wb.pagewidth printOption(pagewidth) "Figure width" $buttonFont
+    mkentryPr  $wb.pageheight printOption(pageheight) "Figure height" $buttonFont
+    mkentryPr  $wb.hoffset printOption(hoffset) "Left margin for printing" $buttonFont
+    mkentryPr  $wb.voffset printOption(voffset) "bottom margin for printing" $buttonFont
+    mkentryPr  $wb.psfilename printOption(psfilename) "postscript filename" $buttonFont
+    mkentryPr  $wb.printer printOption(printer) "Printer to print to" $buttonFont
+    mkentryPr  $wb.gsview printOption(gsview) "postscript viewer, used for printing under Windows" $buttonFont
+   mkentryPr  $wb.xticks printOption(xticks) "Rough number of xticks" $buttonFont
+   mkentryPr  $wb.yticks printOption(yticks) "Rough number of yticks" $buttonFont
+    eval pack $wb.ok $wb.b0 $wb.b1 $wb.b2 $wb.b3 $wb.b4
+    eval pack $topack -expand 1
+
+    foreach v  $paperSizes {
+	set papersize [lindex $v 0]
+        set lower [string tolower $papersize]
+        radiobutton $wb.$lower -text [lindex $v 0] -variable printOption(papersize) \
+	   -value [lindex $v 0] -font $buttonFont -highlightthickness 0 
+    pack $wb.$lower -pady 2 -anchor w -fill x
+    }
+    checkbutton $wb.domargin -variable printOption(domargin) -text "do margin" 
+    pack $wb.domargin -pady 2 -anchor w -fill x
+
+    frame $w.grid
+    pack $w.grid -expand yes -fill both -padx 1 -pady 1
+    grid rowconfig    $w.grid 0 -weight 1 -minsize 0
+    grid columnconfig $w.grid 0 -weight 1 -minsize 0
+}
+
+proc markToPrint { win tag title } {
+    # puts "$win $tag"
+   # bind $win <1> "bindBeginDrag $win %x %y $tag [list $title]"
+    pushBind $win <1> "$win delete printrectangle ; popBind $win <1>"
+    pushBind $win <1> "bindBeginDrag $win %x %y $tag [list $title]; popBind $win <1>"    
+}
+
+proc bindBeginDrag { win x y tag title } {
+    $win delete $tag printrectangle
+    set beginRect "[$win canvasx $x] [$win canvasy $y]"
+    set it1 [eval $win create rectangle $beginRect $beginRect -tags $tag -width 3]
+    set old [bind $win <B1-Motion>]
+    set new "eval $win coords $it1 \
+	    $beginRect \[$win canvasx %x\] \[$win canvasy %y\]; \
+	    "
+    if { "$old" == "$new" } {set old ""}
+    bind $win <B1-Motion> $new
+    bind $win <ButtonRelease-1> "bind $win <B1-Motion> [list $old];\
+	    bind $win <ButtonRelease-1> {} ; unbindAdjustWidth $win $tag [list $title];"
+}
+
+proc unbindAdjustWidth { canv tag title } {
+    set win [winfo parent $canv]
+    global printOption
+
+    set it [$canv find withtag $tag]
+    set co1 [$canv coords $tag]
+    set co [$canv coords $it]
+   # if { "$co" != "$co1" } {puts differ,$co1,$co}
+    desetq "x1 y1 x2 y2" $co
+    set center [expr { ($x1+$x2 )/2}]
+   set h [expr {$y2 - $y1}]
+    set it [$canv find withtag $tag]
+   set new [$canv create rectangle $x1 $y1 $x2 $y2 -outline white -width [expr {$h* .04}] -tags [concat $tag bigger] ]
+
+    # puts "<marginTicks $canv $x1 $y1 $x2 $y2 printrectangle>"
+    marginTicks $canv [storx$win $x1] [story$win $y2] [storx$win $x2] [story$win $y1] "printrectangle marginticks"
+    desetq "a1 b1 a2 b2" [$canv bbox $new]
+   set textit [$canv create text $center [expr {$y1 - $h *.03}] \
+	    -font [font create -family Courier -size 14 -weight bold] -text "$title" \
+	    -anchor s -tags [concat $tag bigger title]]
+
+    set bb [$canv bbox $textit]
+   $canv create rectangle $a1 [lindex $bb 1]  $a2 [expr {$y1 - 0.02 * $h}]  -tags $tag -fill white -outline {}
+   $canv itemconfig $it -width [expr {$h *.002}]
+    $canv raise $it
+    $canv raise $textit
+    $canv raise marginticks
+    if { $printOption(domargin) == 0 } {
+	$canv delete marginticks
+    }
+
+    $canv create text [expr {($a1 + $a2)/2.0}] [expr {$y2 + .01*$h  }] -anchor nw -text "For [getEnv USER] [clock format [clock seconds]]" -font [font create -family Courier -size 10 -weight normal] -tag $tag
+    # puts h=$h
+
+}
+    
+
+proc getPSBbox  { } {
+  set fi [open /home/wfs/sdfplot.ps r]
+  set me [read $fi 500]
+  regexp {BoundingBox: (-*[0-9]+) (-*[0-9]+) (-*[0-9]+) (-*[0-9]+)} $me junk x1 y1 x2 y2
+    set w [expr {72 * 8.5}]
+    set h [expr {72 * 11}]
+    # puts "hei=[expr {$y2-$y1}],tm=[expr {$h - $y2}],bm=$y1"
+    # puts "wid=[expr {$x2-$x1}],lm=$x1,rm=[expr {$w - $x2}]"
+    # puts "hei=[expr {($y2-$y1)/72.0}],tm=[expr {($h - $y2)/72.0}],bm=([expr {$y1/72.0}])"
+    #puts "wid=[expr {($x2-$x1)/72.0}],lm=([expr {$x1/72.0}]),rm=[expr {($w - $x2)/72.0}]"    
+  close $fi
+}
+
+
+## endsource printops.tcl
+# set font {Courier 8}
+set fontCourier8 "-*-Courier-Medium-R-Normal--*-120-*-*-*-*-*-*"
+
+if { "[winfo screenvisual .]" == "staticgray" } { set axisGray black
+}     else  { set axisGray gray60}
+
+set writefile  "Save"
+# make printing be by ftp'ing a file..
+
+if {[catch { set doExit }] } { set doExit ""}
+set width_ [winfo screenwidth .]
+if { $width_ >= 1280 } { set fontSize 12
+  } elseif { $width_ <= 640} { set fontSize 8 } else {
+    set fontSize 10}
+unset width_    
+
+proc makeFrame { w type } {
+    global   writefile doExit fontSize buttonfont ws_openMath   
+    set win $w
+    if { "$w" == "." } {
+        set w "" } else {
+	    catch { destroy $w}
+	    
+	    frame $w
+	    # toplevel $w
+	    # set w $w.new
+            # frame $w
+           # puts "making $w"	
+	    
+    }
+
+    set dismiss "destroy $win"
+    catch { set  parent [winfo parent $win] 
+    if { "$parent" == "." } {
+	set dismiss "destroy ."
+    }
+    }
+    
+    if { "$doExit" != "" } {set dismiss $doExit } 	
+    oset $w type $type
+
+    frame $w.grid
+   #positionWindow $w
+    set c $w.c
+    oset $win c $c
+    if { [catch { set buttonfont} ] } {
+	set buttonfont [font create -family Helvetica -size $fontSize]
+    }
+    set buttonFont $buttonfont    
+    oset $win buttonFont $buttonfont
+
+#    puts "children wb=[winfo children $w]"
+    set wb $w.buttons
+    frame $wb
+    set dismiss [concat $dismiss "; clearLocal $win "]
+
+    button $wb.dismiss -text Dismiss -command $dismiss -font $buttonFont
+    setBalloonhelp $win $wb.dismiss {Close this plot window}
+    button $wb.zoom -text "Zoom" -command "showZoom $w" -font $buttonFont
+    setBalloonhelp $win $wb.zoom {Magnify the plot.  Causes clicking with the left mouse button on the plot, to magnify (zoom in) the plot where you click.  Also causes Shift+Click to  it to unmagnify (zoom out) at that point}
+    oset $w position "" 
+#    button $w.position -textvariable [oloc $w position] -font $buttonFont -width 10
+    label $w.position  -textvariable [oloc $w position] -font $buttonFont -width 10
+    setBalloonhelp $win $w.position {Position of the pointer in real x y coordinates.  For 3d it is the position of the nearest vertex of the polygon the pointer is over.}
+
+    button $wb.help -text "Help" -command "doHelp$type $win" -font $buttonFont
+    setBalloonhelp $win $wb.help {Give more help about this plot window}
+    button $wb.postscript -textvariable writefile -command "writePostscript $w" -font $buttonFont
+    setBalloonhelp $win $wb.postscript {Prints or Saves the plot in postscript format.  The region to be printed is marked using Mark.   Other print options can be obtained by using "Print Options" in the Config menu }
+    
+    button $wb.markrect -text "Mark" -command "markToPrint $c printrectangle \[eval \[oget $win maintitle\]\]" -font $buttonFont
+    setBalloonhelp $win $wb.markrect {Mark the region to be printed.  Causes the left mouse button to allow marking of a rectangle by clicking at the upper left corner, and dragging the mouse to the lower right corner.  The title can be set under "Print Options" under Config}
+    button $wb.replot -text "Replot" -command "replot$type $win" -font $buttonFont
+    setBalloonhelp $win $wb.replot {Use the current settings and recompute the plot.  The settings may be altered in Config}
+    
+	
+    
+    button $wb.config -text "Config" -command "doConfig$type $win" -font $buttonFont
+    setBalloonhelp $win $wb.config {Configure various options about the plot window.  After doing this one may do replot.  Hint: you may leave the config menu on the screen and certain actions take place immediately, such as rotating or computing a trajectory at a point.  To make room for the window you might slide the graph to the right, and possibly shrink it using the unzoom feature}    
+    
+
+
+    bind $win.position <Enter> "+place $win.buttons -in $win.position -x 0 -rely 1.0 ;  after cancel lower $win.position ; raise $win.buttons "
+    bind $win.buttons <Leave> "place forget $win.buttons"
+
+    # pack $wb
+    scrollbar $w.hscroll -orient horiz -command "$c xview"
+    scrollbar $w.vscroll -command "$c yview"
+    # -relief sunken
+    canvas $c  -borderwidth 2 \
+	    -scrollregion {-1200 -1200 1200 1200} \
+	-xscrollcommand "$w.hscroll set" \
+	-yscrollcommand "$w.vscroll set" -cursor arrow -background white
+    # puts "$c config  -height [oget $win height] -width [oget $win width] "
+    set buttonsLeft 1
+    set wid [oget $win width]
+    catch {$c config  -height [oget $win height] -width  $wid
+           oset $win oldCheight [oget $win height]
+           oset $win oldCwidth $wid
+     }
+    # puts "$c height =[$c cget   -height],$c width =[$c cget   -width]"
+    # bind $c <2> "$c scan mark %x %y"
+    bind $c <B3-Motion> "$c scan dragto %x %y"
+    bind $c <3> "$c scan mark %x %y"
+    bind $c <B3-Motion> "$c scan dragto %x %y"    
+    bind $c <Motion> "showPosition $w %x %y"
+    bind $c <Configure> "reConfigure $c %w %h"
+    bind $c <Enter> "raise $win.position"
+    bind $c <Leave> "after 200 lower $win.position"
+    $w.position config -background [$c cget -background]
+     
+     
+    pack  $wb.dismiss $wb.help $wb.zoom   \
+	    $wb.postscript $wb.markrect $wb.replot $wb.config -side top -expand 1 -fill x
+    if { 0 } {
+	pack $w.hscroll -side bottom -expand 1 -fill x
+	pack $w.vscroll -side right -expand 1 -fill y
+    }
+	pack $w.c -side right -expand 1 -fill both
+    
+    pack $w
+    place $w.position -in $w -x 2 -y 2 -anchor nw
+    oset $w position "Menu Here"
+    if { ![info exists ws_openMath(showedplothelp)] ||
+    [llength $ws_openMath(showedplothelp)] < 2 } {
+	lappend ws_openMath(showedplothelp) 1
+	
+	after 100 balloonhelp $w $w.position [list \
+		"Initial help: Moving the mouse over the position \
+		window (top left corner), will bring up a menu.  Holding down \
+		right mouse button and dragging will translate the plot"]
+	after 2000 $w.c delete balloon
+	
+
+    }
+    
+    raise $w.position
+    
+    pack [winfo parent $wb]
+   # update
+#    set wid [ winfo width $win]
+#    if { $wid > [      $c cget -width ] } {
+#    $c config -width $wid
+#	    oset $win width $wid
+#    }
+
+   bind $w <Configure> "resizePlotWindow $w %w %h"
+    return $w    
+}
+
+proc mkentry { newframe textvar text buttonFont } {
+    frame $newframe
+    set parent $newframe
+    set found 0
+    while { !$found } {
+	set parent [winfo parent $parent]
+	if { "$parent" == "" } { break }
+	if { ![catch {  set type [oget $parent type] } ] } {
+	    global plot[set type]Options
+	    foreach v [set plot[set type]Options] {
+		if { "[oloc $parent [lindex $v 0]]" == "$textvar" } {
+		     setBalloonhelp $parent $newframe [lindex $v 2]
+		    set found 1
+		     break
+
+		}
+	    }
+    }
+    }
+    label $newframe.lab1 
+    label $newframe.lab -text "$text:" -font $buttonFont -width 0
+    entry $newframe.e -width 20 -textvariable $textvar -font $buttonFont
+    pack $newframe.lab1 -side left -expand 1 -fill x 
+    pack $newframe.lab -side left
+    pack $newframe.e -side right -padx 3 -fill x
+   # pack $newframe.lab $newframe.e -side left -padx 3 -expand 1 -fill x
+}
+    
+
+proc doHelp { win msg } {
+    makeLocal $win c
+    set atx [$c canvasx 0]
+    set aty [$c canvasy 0]
+    $c create rectangle [expr {$atx -1000}] [expr  {$aty -1000}] 10000 10000 -fill white -tag help
+
+    $c create text [expr {$atx +10}] [expr {$aty + 10.0}] -tag help  -anchor nw  -width 400 -text $msg 
+
+    pushBind $c <1> "$c delete help; popBind $c <1>"
+}
+
+## source push.tcl
+
+###### push.tcl ######
+
+
+
+#
+ #-----------------------------------------------------------------
+ #
+ # pushl --  push VALUE onto a stack stored under KEY
+ #
+ #  Results:
+ #
+ #  Side Effects: 
+ #
+ #----------------------------------------------------------------
+#
+
+global __pushl_ar
+proc pushl { val key  } {
+    global __pushl_ar
+  append __pushl_ar($key) " [list $val]"
+}
+
+
+#
+ #-----------------------------------------------------------------
+ #
+ # peekl --  if a value has been pushl'd under KEY return the 
+ # last value otherwise return DEFAULT.   If M is supplied, get the
+ # M'th one pushed... M == 1 is the last one pushed.
+ #  Results:  a previously pushed value or DEFAULT
+ #
+ #  Side Effects: none
+ #
+ #----------------------------------------------------------------
+#
+proc peekl {key default {m 1}} {
+    global __pushl_ar
+    if { [catch { set val [set __pushl_ar($key) ] } ] } {
+	return $default } else {
+	    set n [llength $val]
+	    if { $m > 0 && $m <= $n } {
+		return [lindex $val [incr n -$m]]
+	    } else { return $default }
+	}
+    }
+    
+    
+
+#
+ #-----------------------------------------------------------------
+ #
+ # popl --  pop off  last value stored under KEY, or else return DFLT
+ #
+ #  Results: last VALUE stored or DEFAULT
+ #
+ #  Side Effects: List stored under KEY becomes one shorter
+ #
+ #----------------------------------------------------------------
+#
+proc popl { key  dflt} {
+    global __pushl_ar
+    
+    if { [catch { set val [set __pushl_ar($key) ] } ] } {
+	return $dflt } else {
+	    set n [llength $val]
+   	    set result [lindex $val [incr n -1]]
+
+	    if { $n > 0 } {
+		set __pushl_ar($key) [lrange $val 0 [expr {$n -1}]]
+	    } else {unset __pushl_ar($key) }
+	    return $result
+	}
+    }
+
+
+#
+ #-----------------------------------------------------------------
+ #
+ # clearl --  clear the list stored under KEY
+ # 
+ #  Result: none
+ #
+ #  Side Effects:  clear the list stored under KEY
+ #
+ #----------------------------------------------------------------
+#
+proc clearl { key } {
+    global __pushl_ar
+    catch { unset __pushl_ar($key) }
+}
+    
+
+
+## endsource push.tcl
+proc pushBind { win key action } {
+    pushl [bind $win $key] [list $win $key ] 
+    bind $win $key $action
+}
+
+proc popBind { win key  } {
+    set binding [popl [list $win $key] {}]
+   
+    bind $win $key $binding
+}
+
+# exit if not part of openmath browser
+proc maybeExit { n } {
+    if { "[info proc OpenMathOpenUrl]" != "" } {
+	uplevel 1 return
+    } else { exit 0 }
+}
+
+proc showPosition { win x y } {
+   # global position c
+    makeLocal $win c
+    # we catch so that in case have no functions or data..
+    catch {
+    oset $win position \
+      "[format {(%.2f,%.2f)}  [storx$win [$c canvasx $x]] [story$win [$c canvasy $y]]]"
+}   }
+
+proc showZoom  { win } {
+  #  global c position
+    makeLocal $win c
+    oset $win position "Click to Zoom\nShift+Click Unzoom"
+     
+    bind $c <1> "doZoom $win %x %y 1"
+    bind $c  <Shift-1> "doZoom $win %x %y -1"
+}
+
+proc doZoom { win x y direction } {
+    set zf [oget $win zoomfactor]
+    if { $direction < 0 } {
+	set zf 	"[expr {1/[lindex $zf 0]}] [expr {1/[lindex $zf 1]}]"
+    }
+    eval doZoomXY $win $x $y $zf
+}
+    
+
+
+#
+ #-----------------------------------------------------------------
+ #
+ # doZoomXY --  given screen coordinates (x,y) and factors (f1,f2)
+ #  perform a scaling on the canvas, centered at (x,y) so that
+ #  the distance in the x direction from this origin is multiplied by f1
+ #  and similarly in the y direction
+ #  Results:
+ #
+ #  Side Effects: scale the canvas, and set new transforms for translation
+ #   from real to canvas coordinates.
+ #----------------------------------------------------------------
+#
+
+proc doZoomXY { win x y facx facy } {
+    makeLocal $win c transform
+    
+    set x [$c canvasx $x]
+    set y [$c canvasy $y]
+
+    $c scale all $x $y $facx $facy
+
+    set ntransform [composeTransform \
+	    "$facx 0 0 $facy [expr {(1-$facx)* $x}] [expr {(1-$facy)* $y}]" \
+	    $transform  ]
+    oset $win transform $ntransform
+    getXtransYtrans $ntransform rtosx$win rtosy$win
+    getXtransYtrans [inverseTransform $ntransform] storx$win story$win
+    axisTicks $win $c
+}
+
+
+#
+ #-----------------------------------------------------------------
+ #
+ # scrollPointTo --  attempt to scroll the canvas so that point
+ #  x,y on the canvas appears at screen (sx,sy)
+ #
+ #  Results: none
+ #
+ #  Side Effects: changes x and y view of canvas 
+ #
+ #----------------------------------------------------------------
+#
+proc scrollPointTo { c x y sx sy } {
+    desetq "x0 y0 x1 y1" [$c cget -scrollregion]
+    $c xview moveto [expr { 1.0*($x-$x0-$sx)/($x1-$x0)} ]
+    $c yview moveto [expr { 1.0*($y-$y0-$sy)/($y1-$y0)} ]
+}
+
+
+
+#
+ #-----------------------------------------------------------------
+ #
+ # reConfigure --  
+ #
+ #  Results:
+ #
+ #  Side Effects: 
+ #
+ #----------------------------------------------------------------
+#
+
+proc reConfigure { c width height  } {
+    set w [winfo parent $c]
+    if { [catch { makeLocal $w oldCwidth oldCheight } ] } {
+	oset $w oldCwidth $width
+	oset $w oldCheight $height
+	return
+    }
+    set oldx [$c canvasx [expr {$oldCwidth/2.0}]]
+    set oldy [$c canvasy [expr {$oldCheight/2.0}]]
+    doZoomXY $w [expr {$oldCwidth/2.0}] [expr {$oldCheight/2.0}] \
+	    [expr {1.0*$width/$oldCwidth}] [expr {1.0*$height/$oldCheight}]
+    scrollPointTo $c $oldx $oldy [expr {$width/2.0}] [expr {$height/2.0}]
+   # update
+    oset $w oldCwidth $width
+    oset $w oldCheight $height
+}
+
+proc writePostscript { win } {
+    global  printOption argv
+    makeLocal $win c transform transform0 xmin ymin xmax ymax
+    set rtosx rtosx$win ; set rtosy rtosy$win
+    drawPointsForPrint $c
+    if { "[$c find withtag printrectangle]" == "" } {
+	# $c create rectangle [$rtosx $xmin] [$rtosy $ymin] [$rtosx $xmax] [$rtosy $ymax] -tags printrectangle -width .5
+	$c create rectangle [$c canvasx 0] [$c canvasy 0] [$c canvasx [$c cget -width ]] [$c canvasy [$c cget -height ]]   -tags printrectangle -width .5	
+	unbindAdjustWidth $c printrectangle [eval [oget $win maintitle]]
+    }
+    $c delete balloon
+	
+	
+    set bbox [eval $c bbox [$c find withtag printrectangle]]
+    desetq "x1 y1 x2 y2" $bbox
+#     set title "unknown plot"
+#     catch { set title [eval $printOption(maintitle)] }
+
+#     $c create text [expr {($x1 + $x2)/2}]  [expr {$y1 + .04 * ($y2 - $y1)}] \
+# 	    -anchor center -text $title -tag title
+
+    update
+set diag [vectorlength [expr {$y1-$x1}] [expr {$y2-$x2}]]
+#  get rid of little arrows that creep onto the outside, ie let
+#  the blank rectangle cover them.
+set x1 [expr {$x1+.01 * $diag}]
+set x2 [expr {$x2-.01 * $diag}]
+set y1 [expr {$y1+.01 * $diag}]
+set y2 [expr {$y2-.01 * $diag}]
+
+    set com "$c postscript  \
+      	    -x  $x1  -y $y1 \
+	    -width [expr {($x2 - $x1)}] \
+            -height [expr {($y2 - $y1)}] \
+	    [getPageOffsets [expr {($x2 - $x1)/(1.0*($y2 - $y1))}] ] "
+
+    #puts com=$com
+    set output [eval $com]
+    switch $printOption(tofile) {
+	0 { global tcl_platform
+	    set usegsview 0  
+	    if { "$tcl_platform(platform)" == "windows" } {
+		set usegsview 1
+	    }
+	    if { $usegsview } {
+		set fi [open $printOption(psfilename) w]
+		puts $fi $output
+		close $fi
+		exec "$printOption(gsview) /S $printOption(psfilename)"
+	    } else {
+	    set fi [open "|lpr -P[set printOption(printer)]" w]
+	    puts $fi $output
+	    close $fi
+	    }
+	}
+	1 { set fi [open $printOption(psfilename) w]
+	puts $fi $output
+	close $fi }
+	2 { global ftpInfo
+	    set ftpInfo(data) $output
+	    ftpDialog $win
+	}
+    }
+#    if { $printOption(tofile) } {
+#	set fi [open $printOption(psfilename) w]
+#    } else { set fi [open "|lpr -P[set printOption(printer)]" w] }
+ #   puts $fi $output
+#    close $fi
+}
+
+
+#
+ #-----------------------------------------------------------------
+ #
+ # ftpDialog --  open up a dialog to send ftpInfo(data) to a file
+ # via http and ftp.   The http server can be specified.
+ #
+ #  Results:
+ #
+ #  Side Effects: 
+ #
+ #----------------------------------------------------------------
+#
+
+set ftpInfo(host) genie1.ma.utexas.edu
+set ftpInfo(viahost) genie1.ma.utexas.edu
+
+proc ftpDialog { win args } {
+    global ftpInfo buttonFont fontSize
+    set fr ${win}plot
+    set usefilename [assoc -filename $args  0]
+    if { "$usefilename" != "0"} {
+	set ftpInfo(filename) $usefilename
+	set usefilename 1
+    }
+    catch { destroy $fr }
+    set ftpInfo(percent) 0
+    set buttonFont [font create -family Courier -size $fontSize]
+    frame $fr -borderwidth 2 -relief raised
+    if { [catch { set ftpInfo(directory) } ] } { set ftpInfo(directory) homework }
+    label $fr.title -text "Ftp Dialog Box" -font [font create -family Helvetica -size [expr {2+ $fontSize}]]
+    mkentry $fr.host ftpInfo(host) "host to write file on" $buttonFont
+    mkentry $fr.viahost ftpInfo(viahost) "host to write to via" $buttonFont
+    mkentry $fr.username ftpInfo(username) "Your User ID on host" $buttonFont
+    mkentry $fr.password ftpInfo(password) "Your password on host" $buttonFont
+    $fr.password.e config -show *
+    mkentry $fr.directory ftpInfo(directory) "remote subdirectory for output" $buttonFont
+
+    if { $usefilename } {
+	mkentry $fr.filename ftpInfo(filename) "filename " $buttonFont
+    } else {
+    mkentry $fr.chapter ftpInfo(chapter) "chapter " $buttonFont
+    mkentry $fr.section ftpInfo(section) "section" $buttonFont
+    mkentry $fr.problemnumber ftpInfo(number) "Problem number" $buttonFont
+    }
+    scale   $fr.scale -orient horizontal -variable ftpInfo(percent) -length 100 
+    button $fr.doit -text "Send it" -command "doFtpSend $fr" -font $buttonFont
+    button $fr.cancel -text "Cancel" -command "destroy $fr" -font $buttonFont
+    set ftpInfo(message) ""
+    label $fr.message  -width 30 -height 3 -textvariable ftpInfo(message) -font $buttonFont
+    eval pack  [winfo  children $fr] -side top 
+    raise $fr
+    place $fr -in $win -relx .5 -rely .5 -anchor center
+   }
+
+proc doFtpSend { fr } {
+    global ftpInfo om_ftp
+    set error ""
+    if { [winfo exists $fr.filename] } {
+	set filename $ftpInfo(filename)
+	set check "host username directory filename"
+    } else {
+	set check "host username directory chapter section number"
+    }
+    foreach v $check {
+	if { $ftpInfo($v) == "" } {
+	    if  { "$error" == "" } { set error "Failed to specify $v " } else {
+		append error ", $v"}
+	}   
+    }
+    if { "$error" != "" } {
+	set ftpInfo(message) $error
+	return -1
+    }
+    if { [winfo exists $fr.chapter] } {
+	set filename "$ftpInfo(chapter).$ftpInfo(section)-$ftpInfo(number).ps"
+    }
+    
+    
+    set res [submitFtp $ftpInfo(viahost) $ftpInfo(host) $ftpInfo(username) $ftpInfo(password) $ftpInfo(directory) $filename]
+    if { "$res" == 1 }  {
+	   after 1000 "destroy $fr"
+    }
+    return $res
+    
+#    set counter [ ftp $ftpInfo(host) $ftpInfo(username) $ftpInfo(password)]
+#    if { $counter < 0 } {
+#	set ftpInfo(message) [concat "Failed:" $om_ftp($counter,log)]
+#	return -1
+#    }
+
+#     if { [ftpDoCd $counter $ftpInfo(directory)] < 0 &&
+#          [ftpDoMkdir $counter $ftpInfo(directory)] > -10 &&
+#        [ftpDoCd $counter $ftpInfo(directory)] < 0 } {
+# 	set ftpInfo(message) [concat "Failed:" $om_ftp($counter,log)]
+# 	return -1
+#     }
+
+    
+#     set res [ftpDoStore $counter $ftpInfo(chapter).$ftpInfo(section)-$ftpInfo(number).ps $ftpInfo(data)]
+#     if { $res < 0 } {
+# 	set ftpInfo(message) "Failed: $om_ftp($counter,log)"
+# 	return -1
+#     } else {
+# 	set ftpInfo(message) "Wrote $ftpInfo(directory)/$ftpInfo(chapter).$ftpInfo(section)-$ftpInfo(number).ps"
+# 	after 1000 destroy $fr
+#     }
+#     ftpClose $counter
+}
+
+proc vectorlength { a b } {
+    return [expr {sqrt($a*$a + $b * $b)} ]
+}
+
+proc setupCanvas { win } {
+  makeLocal $win   xcenter xradius ycenter yradius
+  oset $win xmin [expr {$xcenter - $xradius}]
+  oset $win xmax [expr { $xcenter + $xradius}]
+  oset $win ymin [expr { $ycenter - $yradius}]
+  oset $win ymax [expr { $ycenter + $yradius} ]
+
+}
+
+
+#
+ #-----------------------------------------------------------------
+ #
+ # compose --  A and B are transformations of the form "origin scalefac"
+ # and composing them means applying first b then a, as in a.b.x
+ #  "o s" . x ==> (x-o)*s + o
+ #  Results: the "origin scalefac" which corresponds to the composition.
+ #
+ #  Side Effects: 
+ #
+ #----------------------------------------------------------------
+#
+proc compose { a b } {
+  return  "[expr {-[lindex $a 1]*[lindex $b 0]*[lindex $b 1] \
+      +[lindex $a 1]*[lindex $b 0]-[lindex $a 0]*[lindex $a 1] \
+      +[lindex $a 0]}] [expr {[lindex $a 1]*[lindex $b 1]}]"
+}
+
+proc sparseList { s } {
+    if  { [catch {
+	set val [parseConvert "$s" -variables "x y t"] } err ] } {
+	    tkerror "Syntax error with `$s'\n $err"
+	}
+	return [lindex $val 0]
+    }
+
+proc sparse { s } {
+    set val [sparseList $s]
+    set first $val
+    if { [llength $first] != 1 } {
+	tkerror "only one function wanted" }
+	
+	return [lindex $first 0]
+    }
+
+	
+
+
+#
+ #-----------------------------------------------------------------
+ #
+ # setUpTransforms --  set up transformations for the canvas of WINDOW
+ # so that the image is on FACTOR fractionof the window
+ # these transforms are used for real to screen and vice versa.
+ #  Results: 
+ #
+ #  Side Effects: transform functions rtosx$win rtosy$win storx$win story$win
+ #  are defined.
+ #
+ #----------------------------------------------------------------
+#    
+proc setUpTransforms { win fac } {
+    makeLocal $win xcenter ycenter xradius yradius c
+
+    set delx [$c cget -width]
+    set dely [$c cget -height]
+    set f1 [expr {(1 - $fac)/2.0}]
+    
+    set x1 [expr {$f1 *$delx}]
+    set y1 [expr {$f1 *$dely}]
+    set x2 [expr {$x1 + $fac*$delx}]
+    set y2 [expr {$x1 + $fac*$dely}]
+
+    
+    
+    set xmin [expr {$xcenter - $xradius}]
+    set xmax [expr {$xcenter + $xradius}]
+    set ymin [expr {$ycenter - $yradius}]
+    set ymax [expr {$ycenter + $yradius}]
+    
+    oset $win xmin $xmin
+    oset $win xmax $xmax
+    oset  $win ymin $ymin
+    oset $win ymax $ymax
+    
+    oset $win transform [makeTransform "$xmin $ymin $x1 $y2" "$xmin $ymax $x1 $y1 " "$xmax $ymin $x2 $y2"]
+    set transform [makeTransform "$xmin $ymin $x1 $y2" "$xmin $ymax $x1 $y1 " "$xmax $ymin $x2 $y2"]
+    oset $win transform $transform
+    oset $win transform0 $transform
+    
+    getXtransYtrans $transform rtosx$win rtosy$win
+    getXtransYtrans [inverseTransform $transform] storx$win story$win 
+    
+}
+
+proc inputParse { in } {
+  if { [regexp -indices \
+       {D\[([a-zA-Z][0-9a-zA-Z]*[ ]*),([a-zA-Z][0-9a-zA-Z]*[ ]*)\] *=} \
+	  $in all1 i1 i2] } {
+   set v1 [getOneMatch $in $i1]
+   set v2 [getOneMatch $in $i2]
+   set s1 [string range $in [lindex $all1 1] end]
+
+     if { [regexp -indices {,[ \n]*D\[([a-zA-Z][0-9a-zA-Z]*[ ]*),([a-zA-Z][0-9a-zA-Z]*[ ]*)\] *=} \
+	  $s1 all2 i1 i2] } {
+   set v3  [getOneMatch $s1 $i1]
+   set v4 [getOneMatch $s1 $i2]
+   set end [string first \} $s1 ]
+      set form2 [string range $s1 [expr {1 + [lindex $all2 1]}] [expr {$end -1}]]
+    if { "$v4" != "$v2" } {error "different variable $v2 and $v4"}
+
+    set form1 [string range $in [expr {1 + [lindex $all1 1]}] [expr {[lindex $all2 0] + -1 + [lindex $all1 1]}]]
+    return [list  $v2 $v1 $v3 $form1 $form2]
+    # puts "v1=$v1,form1=$form1,form2=$form2"  
+  } 
+ }
+}
+
+proc composeTransform { t1 t2  } {
+    desetq "a11 a12 a21 a22 e1 e2" $t1
+    desetq "b11 b12 b21 b22 f1 f2" $t2
+   return  [list \
+	   [expr {$a11*$b11+$a12*$b21}] \
+	   [expr {$a11*$b12+$a12*$b22}] \
+	   [expr {$a21*$b11+$a22*$b21}] \
+	   [expr {$a22*$b22+$a21*$b12}] \
+	   [expr {$a11*$f1+$a12*$f2+$e1}] \
+	   [expr {$a21*$f1+$a22*$f2+$e2}] ]
+}
+    
+
+
+#
+ #-----------------------------------------------------------------
+ #
+ # makeTransform --  Given three points mapped to three other points
+ # write down the affine transformation (A.X+B) which performs this.
+ # the arguments are of the form "x1 y1 u1 v1" "x2 y2 u2 v2" "x3 y3 u3 v3"
+ # where (x1,y1) --> (u1,v1)  etc.
+ #  Results: an affine transformation "a b c d e f" which is
+ #     [ a  b ]  [ x1 ] + [ e ]     
+ #     [ c  d ]  [ y1 ]   [ f ]
+ #  Side Effects: none
+ #
+ #----------------------------------------------------------------
+#
+proc makeTransform { P1 P2 P3 } {
+    desetq  "X1 Y1 U1 V1" $P1
+    desetq  "X2 Y2 U2 V2" $P2
+    desetq  "X3 Y3 U3 V3" $P3
+    set tem [expr {double((($X2-$X1)*$Y3+($X1-$X3)*$Y2+($X3-$X2)*$Y1))}]
+    set A [expr {(($U2-$U1)*$Y3+($U1-$U3)*$Y2+($U3-$U2)*$Y1) \
+	    /$tem}]
+    set B [expr {-(($U2-$U1)*$X3+($U1-$U3)*$X2+($U3-$U2)*$X1) \
+	    /$tem}]
+    set E [expr {(($U1*$X2-$U2*$X1)*$Y3+($U3*$X1-$U1*$X3)*$Y2+($U2*$X3-$U3*$X2)*$Y1) \
+	    /$tem}]
+    set C [expr {(($V2-$V1)*$Y3+($V1-$V3)*$Y2+($V3-$V2)*$Y1) \
+	    /$tem}]
+    set D [expr {-(($V2-$V1)*$X3+($V1-$V3)*$X2+($V3-$V2)*$X1) \
+	    /$tem}]
+    set F [expr {(($V1*$X2-$V2*$X1)*$Y3+($V3*$X1-$V1*$X3)*$Y2+($V2*$X3-$V3*$X2)*$Y1) \
+	    /$tem}]
+    set xf ""
+    set yf ""
+    if { $B == 0  && $C == 0 } {
+	set xf "$A*\$X+$E"
+	set yf "$D*\$Y+$F"
+    }
+    return [list $A $B $C $D $E $F]
+}
+
+
+#
+ #-----------------------------------------------------------------
+ #
+ # getXtransYtrans --   If the x coordinate transforms independently
+ #  of the y and vice versa, give expressions suitable for building a
+ # proc. 
+ #  Results:
+ #
+ #  Side Effects: 
+ #
+ #----------------------------------------------------------------
+#
+proc getXtransYtrans { transform p1 p2 } {
+    desetq "a b c d e f"  $transform
+    if { $b == 0  && $c == 0 } {
+	proc $p1 { x } "return \[expr {$a*\$x+$e}\]" 
+	proc $p2 { y } "return \[expr {$d*\$y+$f} \]"
+	return 1
+    }
+    return 0
+}
+
+
+#
+ #-----------------------------------------------------------------
+ #
+ # inverseTransform --   Find the inverse of an affine transformation.
+ #
+ #  Results:
+ #
+ #  Side Effects: 
+ #
+ #----------------------------------------------------------------
+#
+proc inverseTransform { transform } {
+    desetq "a b c d e f" $transform
+    set det [expr {double($a*$d - $b*$c)}]
+    return [list [expr {$d/$det}] [expr {- $b / $det }] [expr {- $c / $det}] [expr {$a / $det}]  [expr {($b*$f-$d*$e)/ $det }] [expr {-($a*$f-$c*$e)/ $det}]]
+
+}
+
+
+#
+ #-----------------------------------------------------------------
+ #
+ # getTicks --  given an interval (a,b) subdivide it and 
+ # calculate where to put the ticks and what to print there.
+ # we want DESIRED number of ticks, but we also want the ticks
+ # to be at points in the real coords of the form .2*10^i or .5*10^j
+ #  Results: the ticks
+ #
+ #  Side Effects: 
+ #
+ #----------------------------------------------------------------
+#
+
+proc getTicks { a b n } {
+    set len [expr {(($b - $a))}]
+    if { $len < [expr {pow(10,-40)}] } { return ""}
+    set best 0
+    foreach v { .1 .2 .5 } {
+	# want $len/(.1*10^i) == $n
+	set val($v)  [expr {ceil(log10($len/(double($n)*$v)))}]
+	set use [expr {$v*pow(10,$val($v))}]
+	set fac [expr {1/$use}]
+	set aa [expr {$a * $fac + .03}]
+	set bb [expr {$b * $fac -.03}]
+	set j [expr {round(ceil($aa)) }]
+	set upto [expr {floor($bb) }]
+	set ticks ""
+	while { $j <= $upto } {
+	    set tt [expr {$j / $fac}]
+	    if { $j%5 == 0 } {
+		append ticks " { $tt $tt }"
+	    } else  {
+		append ticks " $tt"
+	    }
+	    incr j   
+	}
+	set answer($v) $ticks
+	set this [llength $ticks]
+	if { $this  > $best } {
+	    set best $this
+	    set at $v
+	}
+	#puts "for $v [llength $ticks] ticks"
+    }
+    #puts "using $at [llength $answer($at)]"
+     
+    return $answer($at)
+}
+     
+proc axisTicks { win c }  {
+    $c delete axisTicks
+    if { ![catch {oget $win noaxisticks}] } { return }
+    set swid [$c cget -width]
+    set shei [$c cget -height]
+    set x1 [storx$win [$c canvasx 0]]
+    set y1 [story$win [$c canvasy 0]]
+    set x2 [storx$win [$c canvasx $swid]]
+    set y2 [story$win [$c canvasy $shei]]
+    #puts "x1=$x1,y1=$y1,y2=$y2,x2=$x2"
+    if { $y1 > 0  &&  $y2 < 0 } {
+	set ticks [getTicks $x1 $x2 [expr {$swid/50}] ]
+	#puts "ticks=$ticks"
+	set eps [expr {.005 * abs($y1 - $y2)}]
+	set neps [expr {-.005 * abs($y1 - $y2)}]
+	set donext 0
+	foreach v $ticks {
+	    set x [lindex $v 0]
+	    set text [lindex $v 1]
+	    if { $donext } {set text [lindex $v 0] ; set donext 0 }
+	    if { [lindex $v 0] == 0 } { set text "" ; set donext 1 }
+	    #puts " drawTick $c $x 0 0 $neps 0 $eps  $text axisTicks"
+	    drawTick $c $x 0 0 $neps 0 $eps  $text axisTicks
+	    }
+	}
+    if { 0 < $x2 && 0 > $x1 } {
+	set ticks [getTicks $y2 $y1 [expr {$shei/50}]]
+	set eps [expr {.005 * ($x2 - $x1)}]
+	set neps [expr {-.005 * ($x2 - $x1)}]
+	set donext 0
+	foreach v $ticks {
+	    set y [lindex $v 0]
+	    set text [lindex $v 1]
+	    if { $donext } {set text [lindex $v 0] ; set donext 0}
+	    if { [lindex $v 0] == 0 } { set text "" ; set donext 1}
+
+	    drawTick $c 0 $y $neps 0 $eps 0  $text axisTicks
+	    }
+	}
+
+    }
+
+
+#
+ #-----------------------------------------------------------------
+ #
+ # marginTicks --  draw ticks around the border of window
+ #  x1,y1  top left x2,y2 bottom right.
+ #
+ #  Results:
+ #
+ #  Side Effects: 
+ #
+ #----------------------------------------------------------------
+#    
+proc marginTicks { c x1 y1 x2 y2 tag }  {
+    global printOption
+    set win [winfo parent $c]
+
+    if { ![catch {oget $win noaxisticks}] } { return }
+    $c delete marginTicks
+    set ticks [getTicks $x1 $x2 $printOption(xticks)]
+    # puts "x=$x1 $x2"
+    set eps [expr {.008 * ($y1 - $y2)}]
+    set neps [expr {-.008 * ($y1 - $y2)}]
+    foreach v $ticks {
+	set x [lindex $v 0]
+	set text [lindex $v 1]
+	drawTick $c $x $y1 0 0 0 $neps  $text $tag
+	drawTick $c $x $y2 0 0 0 $eps  $text $tag
+	
+    }
+    #puts "y=$y2,$y1"
+    set ticks [getTicks $y1 $y2 $printOption(yticks)]
+    set eps [expr {.005 * ($x2 - $x1)}]
+    set neps [expr {-.005 * ($x2 - $x1)}]
+    set donext 0
+    foreach v $ticks {
+	set y [lindex $v 0]
+	set text [lindex $v 1]
+	drawTick $c $x1 $y 0 0 $eps 0  $text $tag
+	drawTick $c $x2 $y 0 0 $neps 0  $text $tag
+	    }
+	}
+
+proc drawTick {c x y dx dy ex ey n tags} {
+    global axisGray     fontCourier8
+    set win [winfo parent $c]
+    set rtosx rtosx$win ; set rtosy rtosy$win
+    set it [$c create line [$rtosx [expr {$x +$dx}]] [$rtosy [expr {$y +$dy}]] [$rtosx [expr {$x +$ex}]] [$rtosy [expr {$y +$ey}]] -fill $axisGray -tags $tags]
+    $c lower $it
+    
+    if { "$n" != "" } {
+   if { $ey > 0 } { set anch s
+    } elseif { $ex > 0 } {set anch w 
+    } elseif { $ex < 0 } {set anch e
+    } elseif { $ey < 0 } {set anch n}
+    
+    $c create text  [$rtosx [expr {$x +1.5*$ex}]] [$rtosy [expr {$y +1.5*$ey}]] \
+		-text [format "%.8g" $n] -font $fontCourier8 -tags $tags \
+		-anchor $anch
+}   }
+
+proc doConfig { win }  {
+    makeLocal $win c buttonFont
+    $c delete configoptions
+    set canv $c
+   # set w $c.config
+     set w $win.config
+    catch {destroy $w}
+    frame $w -borderwidth 2 -relief raised
+
+    label $w.msg  -wraplength 600 -justify left -text "Plot Setup" -font $buttonFont
+    pack $w
+    pack $w.msg -side top
+    set wb1 $w.choose1
+    frame $wb1
+    set wb2 $w.choose2
+    frame $wb2
+    pack $wb1 $wb2 -side left -fill x -pady 2m
+    set item [$canv create window [$canv canvasx 10] [$canv canvasy  10] -window $w -anchor nw -tags configoptions]
+    button $wb1.dismiss -command  "$canv delete $item; destroy $w " -text "ok" -font $buttonFont
+    button $wb1.printoptions -text "Print Options" -command "mkPrintDialog .dial -canvas $c -buttonfont $buttonFont " -font $buttonFont
+
+    pack $wb1.dismiss  $wb1.printoptions -side top 
+    return "$wb1 $wb2"
+}
+# mkentry { newframe textvar text } 
+
+set show_balloons 1
+
+proc balloonhelp { win subwin msg } {
+    global show_balloons
+    if { $show_balloons == 0 } return;
+    makeLocal $win c buttonFont
+    set x0 [winfo rootx $win]
+    set y0 [winfo rooty $win]
+    
+    
+    set atx [expr {[winfo rootx $subwin] + [winfo width $subwin] - $x0} ]
+    set aty [expr {[winfo rooty $subwin] + [winfo height $subwin] - $y0} ]
+
+    set wid [$c cget -width]
+    set wid2 [expr {round ($wid /2.0)}]
+    set wid10 [expr {round ($wid /10.0)}]
+
+    if { $aty <=1 } { set aty 30 } 
+    incr aty 10
+    incr atx 10
+    set atx [$c canvasx $atx]
+    set aty [$c canvasy $aty]
+    #puts "$atx $aty"
+    $c delete balloon
+    $c create text $atx $aty -anchor nw -text $msg -font $buttonFont -width $wid2 -fill white -fill black -tags "balloon btext"
+    desetq "x1 y1 x2 y2" [$c bbox btext]
+
+    set x1 [expr {$x1 - .3*($x2-$x1)}]
+    set x2 [expr {$x2 + .3*($x2-$x1)}]
+    
+    set y1 [expr {$y1 - .3*($y2-$y1)}]
+    set y2 [expr {$y2 + .3*($y2-$y1)}]
+
+    eval $c create polygon $x1 $y1  $x2 $y1 $x2 $y2 $x1 $y2  -fill beige -tags balloon -smooth 1
+    $c raise btext
+    
+}
+
+proc setBalloonhelp { win subwin msg } {
+    makeLocal $win c
+    bind $subwin <Enter> "balloonhelp $win $subwin [list $msg]"
+    bind $subwin <Leave> "$c delete balloon"
+}
+    
+
+
+#
+ #-----------------------------------------------------------------
+ #
+ # minMax --  Compute the max and min of the arguments, which may
+ # be vectors or numbers
+ #
+ #  Results: list of MIN and MAX
+ #
+ #  Side Effects: none
+ #
+ #----------------------------------------------------------------
+#
+proc minMax { args } {
+    set max [lindex [lindex $args 0] 0] ; set min $max ;
+    foreach vec $args {
+	foreach v $vec {
+	    if { $v > $max } {set max $v }
+	    if { $v < $min} {set min $v }
+	}
+    }
+    return [list $min $max]
+}
+
+proc matrixMinMax { list } {
+# compute the min max of the list    
+    set min +10e300
+    set max -10e300
+    foreach mat $list {
+	foreach row $mat {
+	    foreach v [ldelete nam $row] {
+		if { $v > $max } {catch  { set max [expr {$v + 0}] }}
+		if { $v < $min} {catch  { set min [expr {$v + 0}] }}
+		}
+	    }
+	}
+    list $min $max
+}
+	
+proc omPlotAny { data args } {
+    # puts "data=<[lindex $data 0]>"
+    set command [list [lindex [lindex $data 0] 0]  -data [lindex $data 0] ]
+    if { "[lindex $command 0]" == "plot2d" } {
+	lappend command -xfun {}
+    }
+    foreach v $args { [lappend command $v] }
+    eval $command
+    #eval [lindex [lindex $data 0] 0] -xfun [list {}] -data [list [lindex $data 0]] $args
+}
+
+
+proc resizeSubPlotWindows { win wid height } {
+    set at [$win yview "@0,0"]
+    foreach w [winfo children $win] {
+	if { [string match plot* [lindex [split $w .] end]] } {
+	    resizePlotWindow $w [winfo width $w] $height
+    }
+  }
+  if { "$at" != "" } { $win yview $at} 
+ }
+
+
+	    
+proc resizePlotWindow  { w width height } {
+    if { [winfo width $w.c] <= 1 } {
+	after 100 update ;
+	return }
+	if { ![catch { set tem [oget $w lastResize] } ] && [expr {[clock seconds] - $tem }] < 2 } { return
+} else { oset $w lastResize [clock seconds ]
+    }
+    #puts "resizePlotWindow $w $width $height"
+	
+   # return
+   set par [winfo parent $w]
+   set facx 1.0
+   set facy 1.0    
+    set wid [winfo width $par]
+    set hei [winfo height $par]
+    
+   if { "[winfo class $par]" == "Text" } {
+   set dif 10
+
+   set wid1 $wid ; set hei1 $hei
+   #puts "now w=$w"
+   #set wid1 [getPercentDim [oget $w widthDesired] width $par]        
+   catch {set wid1 [getPercentDim [oget $w widthDesired] width $par] }
+   catch {set hei1 [getPercentDim [oget $w heightDesired] height $par] }
+   set wid [expr {($wid1 > $wid - 30 ? $wid - 30 : $wid1 )}]
+   set hei [expr {($hei1 > $hei - 30 ? $hei - 30 : $hei1 )}]
+   } else {
+       set dif 10
+
+   }
+   
+
+    #puts "width arg=$width,width $w=[winfo width $w],wid of $par=$wid,height=$height,hei=$hei,\[winfo width \$w.c\]=[winfo width $w.c]"
+#     if { $width > $wid -20 || $wid > $width -20 }
+    if { (abs($width-$wid) > $dif ||  abs($height-$hei) > $dif)
+&&  [winfo width $w.c] > 1 } {
+    set eps [expr {2 * [$w.c cget -insertborderwidth] + [$w.c cget -borderwidth] }]
+    set epsx $eps
+    set epsy $eps
+	#puts "reconfiguring: w=$w,par=$par,dif=$dif,widths=$wid, \
+    $width,[winfo width $par],[winfo width $w],[winfo width $w.c]\
+	heights=$hei,$height,[winfo height $par],[winfo height $w],\
+    [winfo height $w.c]"
+    
+    set extrawidth [expr {([winfo width $w] - [winfo width  $w.c]) +$epsx}]
+    set extraheight [expr {([winfo height $w] - [winfo height  $w.c]) +$epsy}]
+    set nwidth [expr {$wid - ($extrawidth > 0  ? $extrawidth : 0)}]
+    set nheight [expr {$hei - ($extraheight > 0  ? $extraheight : 0)}]
+    
+    #puts "$w.c config -width $nwidth  -height $nheight, extraheight=$extraheight,epsy=$epsy"
+    $w.c config -width $nwidth  -height $nheight
+
+		}
+    		
+ }
+
+
+
+
+
+
+## endsource plotconf.tcl
+## source plotdf.tcl
+
+###### plotdf.tcl ######
+#######################################################################
+#######  Copyright William F. Schelter.  All rights reserved.  ########
+#######################################################################
+
+set plotdfOptions {
+    {dxdt "x-y^2+sin(x)*.3" {specifies dx/dt = dxdt.  eg -dxdt "x+y+sin(x)^2"} }
+    {dydt "x+y" {specifies dy/dt = dydt.  eg -dydt "x-y^2+exp(x)"} }
+    {dydx "" { may specify dy/dx = x^2+y,instead of dy/dt = x^2+y and dx/dt=1 }}
+    {xradius 10 "Width in x direction of the x values" }
+    {yradius 10 "Height in y direction of the y values"}
+    {width 500 "Width of canvas in pixels"}
+    {height 500 "Height of canvas in pixels" }
+    {scrollregion {} "Area to show if canvas is larger" }
+    {xcenter 0.0 {(xcenter,ycenter) is the origin of the window}}
+    {ycenter 0.0 "see xcenter"}
+    {tinitial 0.0 "The initial value of variable t"}
+    {nsteps 100 "Number of steps to do in one pass"}
+    {xfun "" "A semi colon separated list of functions to plot as well"}
+    {tstep "" "t step size"}
+    {direction "both" "May be both, forward or backward" }
+    {windowname ".dfplot" "window name"}
+    {linecolors {blue green red brown gray black} "colors to use for lines in data plots"}
+    {doTrajectoryAt "1.0 1.0" "Place to calculate trajectory"}
+    {linewidth "0.6" "Width of integral lines" }
+    {nolines 0 "If not 0, plot points and nolines"}
+    {bargraph 0 "If not 0 this is the width of the bars on a bar graph" }
+    {plotpoints 0 "if not 0 plot the points at pointsize" }
+    {pointsize 2 "radius in pixels of points" }
+    {autoscale "x y" "Set {x,y}center and {x,y}range depending on data and function. "}
+    {zoomfactor "1.6 1.6" "Factor to zoom the x and y axis when zooming.  Zoom out will be reciprocal" }
+    {labelposition "10 35" "Position for the curve labels nw corner"}
+}
+
+if { "[info proc makeFrame]" == "" } { source "plotconf.tcl" }
+
+    
+
+ proc makeFrameDf { win } {
+   set w [makeFrame $win df]
+    makeLocal $win c dydx
+
+    set top $win
+   # puts "w=$w,win=$win"
+    catch { set top [winfo parent $win]}
+    catch {
+
+    wm title $top "Direction Fields"
+    wm iconname $top "DF plot"
+#    wm geometry $top 750x700-0+20
+   }
+    set wb $w.buttons
+   makeLocal $win buttonFont 
+   label $w.msg  -wraplength 600 -justify left -text "A direction field plotter by William Schelter" -font $buttonFont
+   
+  button $wb.integrate -text "Integrate" -command "setForIntegrate $w" -font $buttonFont
+   setBalloonhelp $win $wb.integrate {Causes clicking on the  plot with the left mouse button at a point, to draw a trajectory passing through that point.   Under Config there is an entry box which allows entering exact x,y coordinates, and which also records the place of the last trajectory computed.}
+   
+  setForIntegrate $w
+  pack $wb.integrate -side right -expand 1 -fill x
+ # pack $w.msg -side top
+  pack $w
+  return $win
+}
+
+proc swapChoose {win msg winchoose } {
+   # global dydx dxdt dydt
+    
+    if { "$msg" == "dydt" } {
+	pack $winchoose.dxdt -before $winchoose.dydt -side bottom
+	oset $win dydx ""
+	$winchoose.dydt.lab config -text "dy/dt"
+    } else {
+	pack forget $winchoose.dxdt
+	oset $win dxdt 1
+	oset $win dydx " "
+	$winchoose.dydt.lab config -text "dy/dx"
+    }
+}
+	
+
+proc doHelpdf { win } {
+    global Parser
+ doHelp $win [join [list \
+{
+William Schelter's solver/plotter for ode systems.
+
+To QUIT this HELP click here.
+
+Clicking at a point computes the trajectory
+(x(t),y(t)) starting at that point, and satisfying
+the differential equation
+	
+  dx/dt = dxdt
+  dy/dt = dydt
+
+By clicking on Zoom, the mouse now allows you to zoom
+in on a region of the plot.  Each click near a point
+magnifies the plot, keeping the center at the point
+you clicked.  Depressing the SHIFT key while clicking
+zooms in the opposite direction.
+
+To resume computing trajectories click on Integrate.
+
+To change the differential equation, click on Config and
+enter new values in the entry windows, and then click on
+Replot in the main menu bar.
+
+Holding the right mouse button down allows you to drag
+(translate) the plot sideways or up and down.
+
+Additional parameters such as the number of steps (nsteps),
+the initial t value (tinitial), and the x and y centers
+and radii, may be set under the  Config menu.
+
+You may print to a postscript printer, or save the plot \
+as a postscript file, by clicking on save.   To change \
+between printing and saving see the Print Options under Config.
+	
+} $Parser(help)]]
+}
+
+proc setForIntegrate { win} {
+    makeLocal $win c
+    $c delete printrectangle
+    bind $c  <1> "doIntegrateScreen $win %x %y "
+}
+
+## source rk.tcl
+
+###### rk.tcl ######
+#######################################################################
+#######  Copyright William F. Schelter.  All rights reserved.  ########
+#######################################################################
+
+#proc try { } {
+    #  proc ff { a b c } { return [expr {$b + $c}] }
+    #  proc gg { a b c } { return [expr {$b - $c}] }
+#  rungeKutta ff gg 0.2 0.2 0 .01 10
+#}
+
+proc rungeKutta { f g t0 x0 y0  h nsteps } {
+  set n $nsteps
+  set ans "$x0 $y0"
+  set xn $x0
+  set yn $y0
+  set tn $t0
+    set h2 [expr {$h / 2 }]
+    set h6 [expr {$h / 6.0 }]
+ catch {
+  while { [incr nsteps -1] >= 0 } {
+  
+
+  set kn1 [$f $tn $xn $yn]
+  set ln1 [$g $tn $xn $yn]
+
+      set arg "[expr {$tn + $h2}] [expr {$xn + $h2 * $kn1}] [expr {$yn + $h2*$ln1}]"
+  set kn2 [eval $f $arg]
+  set ln2 [eval $g $arg]
+
+      set arg "[expr {$tn + $h2}] [expr {$xn + $h2 * $kn2}] [expr {$yn +$h2*$ln2}]"
+  set kn3 [eval $f $arg]
+  set ln3 [eval $g $arg]
+
+      set arg "[expr {$tn + $h}] [expr {$xn + $h * $kn3}] [expr {$yn + $h*$ln3}]"
+  set kn4 [eval $f $arg]
+  set ln4 [eval $g $arg]
+
+      set xn [expr {$xn + $h6 * ($kn1+2*$kn2+2*$kn3+$kn4)}]
+      set yn [expr {$yn + $h6 * ($ln1+2*$ln2+2*$ln3+$ln4)}]
+      set tn [expr {$tn+ $h}]
+
+  append ans " $xn $yn"
+  }
+ }
+
+ # puts "ans=$ans"
+ return $ans 
+}
+
+proc pathLength { list } {
+  set sum 0
+  foreach { x y } $list {
+      set sum [expr {$sum + sqrt($x*$x+$y*$y)}]
+ }
+  return $sum
+}
+proc rungeKuttaA { f g t0 x0 y0  h nsteps } {
+  set ans [rungeKutta $f $g $t0 $x0 $y0 $h $nsteps]
+  set count 0
+  # puts "retrying([llength $ans]) .."
+  while { [llength $ans] < $nsteps * .5  && $count < 7 } {
+       incr count
+       #set leng [pathLength $ans]
+       #if { $leng == 0 } {set leng .001}
+       set th [expr {$h / 3.0}]
+       if { $th  < $h }  { set h $th }
+       set ans  [rungeKutta $f $g $t0 $x0 $y0 $h $nsteps]
+      # puts -nonewline "..(h=[format "%.5f" $h],pts=[llength $ans])"
+       # flush stdout
+  }
+  return $ans
+}
+
+  
+
+## endsource rk.tcl
+
+# sample procedures
+# proc xff { t x y } { return [expr {$x + $y }] }
+# proc yff { t x y } { return [expr {$x - $y }] }
+
+proc doIntegrateScreen { win sx sy  } {
+    makeLocal $win c
+    doIntegrate $win [storx$win [$c canvasx $sx]] [story$win [$c canvasy $sy]]
+}
+
+proc doIntegrate { win x y } {
+    # global xradius yradius c tstep  nsteps
+    makeLocal $win xradius yradius c tstep  nsteps direction linewidth
+    set rtosx rtosx$win ; set rtosy rtosy$win      
+    oset $win doTrajectoryAt [format "%.10g  %.10g" $x $y]
+    # puts "doing at $doTrajectoryAt"
+    set steps $nsteps
+    if { "$tstep" == "" } {  
+	set h [expr {[vectorlength $xradius $yradius] / 200.0}]
+	set tstep $h
+    } else {set h $tstep }
+    
+    # puts h=$h
+    set todo $h
+    switch $direction {
+	forward { set todo "$h" }
+	backward { set todo "[expr {- $h}]" }
+	both { set todo "$h [expr {- $h}]" }
+    }
+
+    foreach h $todo {
+	
+	set ans [rungeKuttaA xff yff 0 $x $y $h $steps]
+	set i -1
+	set xn1 [$rtosx [lindex $ans [incr i]]]
+	set yn1 [$rtosy [lindex $ans [incr i]]]
+	set lim [expr {$steps * 2}]
+	set mee [expr {pow(10.0,9)}]
+	catch  { 
+	    while { $i <= $lim } {
+		set xn2  [$rtosx [lindex $ans [incr i]]]
+		set yn2  [$rtosy [lindex $ans [incr i]]]
+		# puts "$xn1 $yn1"
+		# xxxxxxxx following is for a bug in win95 version
+		if { abs($xn1) + abs($yn1) +abs($xn2)+abs($yn2) < $mee    } {
+		    $c create line $xn1 $yn1 $xn2 $yn2 -tags path -width $linewidth
+		}
+		# puts "$xn1 $yn1"
+		set xn1 $xn2
+		set yn1 $yn2
+	    }
+    }   }
+}
+
+
+#
+ #-----------------------------------------------------------------
+ #
+ # $rtosx,$rtosy --  convert Real coordinate to screen coordinate
+ #
+ #  Results: a window coordinate
+ #
+ #  Side Effects: 
+ #
+ #----------------------------------------------------------------
+
+
+#
+ #-----------------------------------------------------------------
+ #
+ # $storx,$story --  Convert a screen coordinate to a Real coordinate.
+ #
+ #  Results:
+ #
+ #  Side Effects: 
+ #
+ #----------------------------------------------------------------
+#
+
+proc drawArrowScreen { c atx aty dfx dfy } {
+
+    set x1 [expr {$atx + $dfx}]
+    set y1 [expr {$aty + $dfy}]
+    #   set x2 [expr {$atx + .8*$dfx +.1* $dfy}]
+    #   set y2 [expr {$aty + .8*$dfy - .1* $dfx}]
+    #   set x3 [expr {$atx + .8*$dfx -.1* $dfy}]
+    #   set y3 [expr {$aty + .8*$dfy + .1* $dfx}]
+  $c create line $atx $aty $x1 $y1 -tags arrow -fill blue -arrow last -arrowshape {3 5 2}
+#  $c create line $x2 $y2  $x1 $y1 -tags arrow -fill red
+#  $c create line $x3 $y3 $x1 $y1 -tags arrow -fill red
+}
+
+proc drawDF { win tinitial } {
+    global  axisGray
+    makeLocal  $win xmin xmax   xcenter ycenter c ymin ymax transform
+
+    # flush stdout
+    set rtosx rtosx$win ; set rtosy rtosy$win
+    set storx   storx$win  ;   set story   story$win  
+    set stepsize 30
+    set min 100000000000.0
+    set max 0.0
+    set t0 $tinitial
+    set xfactor [lindex $transform 0]
+    set yfactor [lindex $transform 3]
+    set extra $stepsize  
+    set uptox [expr {[$rtosx $xmax] + $extra}]
+    set uptoy [expr {[$rtosy $ymin] + $extra}]
+    # draw the axes:
+    #puts "draw [$rtosx $xmin] to $uptox"
+    for { set x [expr {[$rtosx $xmin] - $extra}] } { $x < $uptox } { set x [expr {$x +$stepsize}] } {
+	for { set y [expr {[$rtosy $ymax] - $extra}] } { $y < $uptoy } { set y [expr {$y + $stepsize}] } {
+	    set args "$t0 [$storx $x] [$story $y]"
+	    set dfx [expr {$xfactor * [eval xff $args]}]
+	    # screen y is negative of other y
+	    set dfy [expr  {$yfactor * [eval yff $args]}]
+	    #     puts "$dfx $dfy"
+	    set len  [vectorlength $dfx $dfy]
+	    append all " $len $dfx $dfy "
+	    if { $min > $len } { set min $len }
+	    if { $max < $len } {set  max $len}
+    }   }
+    set fac [expr {($stepsize -5 -8)/($max - $min)}]
+    set arrowmin 8
+    set arrowrange [expr {$stepsize -4 - $arrowmin}]
+    set s1 [expr {($arrowrange*$min+$arrowmin*$min-$arrowmin*$max)/($min-$max)}]
+    set s2 [expr {$arrowrange/($max-$min) }]
+    # we calculate fac for each length, so that
+    # when we multiply the vector times fac, its length
+    # will fall somewhere in [arrowmin,arrowmin+arrowrange].
+    # vectors of length min and max resp. should get mapped
+    # to the two end points.
+    # To do this we set fac [expr {$s1/$len + $s2}]
+    # puts "now to draw,s1=$s1 s2=$s2,max=$max,min=$min"
+    # puts "xfactor=$xfactor,yfactor=$yfactor"
+
+    
+    set i -1
+    for { set x [expr {[$rtosx $xmin] - $stepsize}] } { $x < $uptox } { set x [expr {$x +$stepsize}] } {
+	for { set y [expr {[$rtosy $ymax] - $stepsize}] } { $y < $uptoy } { set y [expr {$y + $stepsize}] } {
+	    
+	    
+	    set len [lindex $all [incr i]]
+	    
+	    set fac [expr {$s1/$len + $s2}]
+	    set dfx [lindex $all [incr i]]
+	    set dfy [lindex $all [incr i]]
+	    #puts "[$storx $x] [$story $y] x=$x y=$y dfx=$dfx dfy=$dfy fac=$fac"
+	    # puts "$len $dfx $dfy" 
+	    drawArrowScreen $c $x $y [expr {$fac * $dfx}] [expr {$fac * $dfy}]
+        }
+    }
+    
+    $c create line [$rtosx 0 ] [$rtosy -1000] [$rtosx 0] [$rtosy 1000] \
+	    -fill $axisGray
+    $c create line [$rtosx -1000] [$rtosy 0] [$rtosx 1000] [$rtosy 0] \
+	    -fill $axisGray
+    axisTicks $win $c
+}
+
+proc parseOdeArg {  s } {
+    set w "\[ ]*"
+    set exp "\[dD]$w\\($w\(\[xyz])$w,$w\(\[xyt])$w\\)$w=(\[^;]+)"
+    while { [regexp $exp $s junk x t expr ] } {
+	lappend ans  -d${x}d$t
+	lappend ans $expr
+	regexp -indices $exp $s junk x t expr
+	set s [string range $s [lindex $junk 1] end]
+    }
+    if { ![info exists ans] } {
+	error "bad -ode argument: $s\nwant d(y,x)=f(x,y) \n   OR d(x,t)=f(x,y,t) d(y,t)=g(x,y,t) "
+    }
+    return $ans
+}
+
+proc plotdf { args } {
+    global plotdfOptions   printOption printOptions plot2dOptions
+    # puts "args=$args"
+    # to see options add: -debug 1
+    set win [assoc -windowname $args]
+    if { "$win" == "" } {set win [getOptionDefault windowname $plotdfOptions] }
+    if { "[lindex $args 0]" == "-ode" } {
+	set tem [parseOdeArg [lindex $args 1]]
+	set args [lrange $args 2 end]
+	set args [concat $tem  $args]
+    }
+    global [oarray $win]
+    getOptions $plotdfOptions $args -usearray [oarray $win] 
+    makeLocal $win dydx
+    if { "$dydx" !="" } { oset $win dxdt 1 ; oset $win dydt $dydx }
+    setPrintOptions $args
+    
+    makeFrameDf $win
+    oset $win maintitle [concat "makeLocal $win  dxdt dydt dydx ;"  \
+	    {if { "$dydx" == "" } { concat "dx/dt = $dxdt , dy/dt = $dydt"}  else {
+	concat "dy/dx = $dydt" } } ]
+	replotdf $win
+    }
+
+proc replotdf { win } {
+    global plotdfOptions 
+    makeLocal $win c dxdt dydt tinitial nsteps xfun
+    
+    setUpTransforms $win 1.0
+    setXffYff $dxdt $dydt
+    $c delete all
+    setForIntegrate $win
+    oset $win curveNumber -1
+    drawDF $win $tinitial
+    foreach v [sparseList $xfun] {
+	proc _xf {  x  } "return \[expr { $v } \]"
+	set listPts [calculatePlot $win _xf $nsteps]
+	drawPlot $win [RealtoScreen $win $listPts]  -fill  [nextColor $win] \
+		-tags line[oget $win curveNumber]
+    }
+    
+}    
+
+proc setXffYff { dxdt dydt } {
+   proc xff { t x y } "return \[expr { [sparse $dxdt] } \]"
+   proc yff { t x y } "return \[expr { [sparse $dydt] } \]"
+}
+
+proc doConfigdf { win } {
+    desetq "wb1 wb2" [doConfig $win]
+    makeLocal $win buttonFont 
+    frame $wb1.choose1
+    set frdydx $wb1.choose1
+    button $frdydx.dydxbut -command "swapChoose $win dydx $frdydx " \
+	    -text "dy/dx" -font $buttonFont
+    button $frdydx.dydtbut -command "swapChoose $win dydt $frdydx" \
+	    -text "dy/dt,dx/dt" -font $buttonFont
+    mkentry $frdydx.dxdt [oloc $win dxdt] "dx/dt" $buttonFont
+    mkentry $frdydx.dydt [oloc $win dydt] "dy/dt" $buttonFont
+    pack $frdydx.dxdt  $frdydx.dydt -side bottom  -fill x -expand 1
+    pack $frdydx.dydxbut $frdydx.dydtbut -side left -fill x -expand 1
+    
+    foreach w {linewidth xradius yradius xcenter ycenter tinitial nsteps tstep direction xfun linecolors } {
+	mkentry $wb1.$w [oloc $win $w] $w $buttonFont
+	pack $wb1.$w -side bottom -expand 1 -fill x
+    }
+    mkentry $wb1.doTrajectoryAt [oloc $win doTrajectoryAt] \
+	    "Trajectory at" $buttonFont
+    bind $wb1.doTrajectoryAt.e <KeyPress-Return> \
+	    "eval doIntegrate $win \[oget $win doTrajectoryAt\] "
+    pack  $wb1.doTrajectoryAt   $frdydx    -side bottom -expand 1 -fill x
+    if { "[oget $win dydx]" != "" } { swapChoose $win dydx $frdydx }
+    setForIntegrate $win
+    
+}
+
+
+## endsource plotdf.tcl
+## source plot2d.tcl
+
+###### plot2d.tcl ######
+
+
+set p .plot
+catch { destroy $p }
+
+set plot2dOptions { 
+    {xradius 10 "Width in x direction of the x values" }
+    {yradius 10 "Height in y direction of the y values"}
+    {width 500 "Width of canvas in pixels"}
+    {height 500 "Height of canvas in pixels" }
+    {xcenter 0.0 {(xcenter,ycenter) is the origin of the window}}
+    {xfun "" {function of x to plot eg: sin(x) or "sin(x);x^2+3" }}
+    {nsteps "100" "mininmum number of steps in x direction"}
+    {ycenter 0.0 "see xcenter"}
+    {screenwindow "20 20 700 700" "Part of canvas on screen"}
+    {windowname ".plot2d" "window name"}
+    {nolines 0 "If not 0, plot points and nolines"}
+    {bargraph 0 "If not 0 this is the width of the bars on a bar graph" }
+    {plotpoints 0 "if not 0 plot the points at pointsize" }
+    {pointsize 2 "radius in pixels of points" }
+    {linecolors {blue green red brown gray black} "colors to use for lines in data plots"}
+    {labelposition "10 35" "Position for the curve labels nw corner"}
+    {xaxislabel "" "Label for the x axis"}
+    {yaxislabel "" "Label for the y axis"}
+    {autoscale "x y" "Set {x,y}center and {x,y}range depending on data and function. "}
+    {zoomfactor "1.6 1.6" "Factor to zoom the x and y axis when zooming.  Zoom out will be reciprocal" }
+    {errorbar 0 "If not 0 width in pixels of errorbar.  Two y values supplied for each x: {y1low y1high y2low y2high  .. }"} 
+    {data "" "List of data sets to be plotted.  Has form { {xversusy {x1 x2 ... xn} {y1 .. yn ... ym}} .. {againstIndex {y1 y2 .. yn}}  .. }"}
+}
+    
+proc mkPlot2d { args } {
+    global plot2dOptions  printOption axisGray
+    #puts "args=<$args>"
+    # global  screenwindow c xmax xmin ymin ymax 
+    # eval global [optionFirstItems $plot2dOptions]
+    set win [assoc -windowname $args]
+    if { "$win" == "" } {
+	set win [getOptionDefault windowname $plot2dOptions] }
+    global  [oarray $win]
+    set data [assoc -data $args ]
+    # puts ranges=[plot2dGetDataRange $data]
+   if {( "[assoc -yradius $args nothere]" != "nothere" ||
+	"[assoc -ycenter $args nothere]" != "nothere" ) &&
+	"[assoc -autoscale $args nothere]" == "nothere" } {
+	    lappend args -autoscale "x"
+	}
+    getOptions $plot2dOptions $args -usearray [oarray $win]
+    oset $win curveNumber -1
+    setPrintOptions $args
+    oset $win maintitle ""	
+    setupCanvas $win 
+    catch { destroy $windowname }
+    
+    makeFrame2d $win
+    makeLocal $win c
+    return $win
+    
+}
+
+proc  makeFrame2d  { win } {
+    set w [makeFrame $win 2d]
+    set top $w
+    catch { set top [winfo parent $w]}
+    catch {
+    wm title $top "Schelter's 2d Plot Window"
+    wm iconname $top "2d plot"
+   # wm geometry $top 750x700-0+20
+   }
+    pack $w
+   return $w
+
+}
+
+proc doConfig2d { win } {
+    desetq "wb1 wb2" [doConfig $win]
+    makeLocal $win buttonFont 
+    mkentry $wb1.nsteps [oloc $win nsteps]  "Number of mesh grids"  $buttonFont
+    mkentry $wb1.xfun [oloc $win xfun]  "y=f(x)"  $buttonFont
+    # button .jim.buttons.rot "rotate" -command "bindForRotation"
+    # pack .jim.buttons.rot
+    pack $wb1.xfun  $wb1.nsteps -expand 1 -fill x
+    foreach w {xradius yradius xcenter ycenter linecolors autoscale } {
+	mkentry $wb1.$w [oloc $win $w] $w $buttonFont
+	pack $wb1.$w -side bottom -expand 1 -fill x
+    }
+}
+
+proc doHelp2d {win } {
+ global Parser
+
+doHelp $win [join [list \
+{
+
+William Schelter's plotter for two dimensional graphics.
+
+To QUIT this HELP click here.
+
+By clicking on Zoom, the mouse now allows you to zoom \
+in on a region of the plot.  Each click near a point \
+magnifies the plot, keeping the center at the point \
+you clicked.  Depressing the SHIFT key while clicking \
+zooms in the opposite direction. 
+
+To change the functions plotted, click on Config and \
+enter new values in the entry windows, and then click on \
+Replot in the main menu bar.
+
+Holding the right mouse button down allows you to drag
+(translate) the plot sideways or up and down.
+
+Additional parameters such as the number of steps (nsteps), \
+and the x and y centers and radii, may be set under the \
+Config menu.
+
+You may print to a postscript printer, or save the plot \
+as a postscript file, by clicking on save.   To change \
+between printing and saving see the Print Options under Config.
+	
+
+
+
+
+} $Parser(help)]]
+}
+
+
+set   plot(numberPlots) 4
+proc mkExtraInfo { name args } {
+    # global plot 	
+    catch { destroy $name }
+
+    toplevel $name
+    wm geometry $name -10+10
+ # pack $name
+    set canv [assoc -canvas $args ]
+    set i 0
+    set w $name
+    frame $w.grid
+    pack $w.grid -expand yes -fill both -padx 1 -pady 1
+    grid $w.grid
+    grid rowconfig    $w.grid 0 -weight 1 -minsize 0
+    grid columnconfig $w.grid 0 -weight 2 -minsize 0
+    
+    set i 0
+    label $w.title -text "Extra Plotting Information" -width 50
+    grid $w.title -in $w.grid -columnspan 2 -row 0 -column 0
+    incr i
+    label $w.labppl -text "Plot Function f(x)"
+    label $w.labcol -text "plot color"
+    grid $w.labppl -padx 1 -in $w.grid  -pady 1 -row $i -column 0 -sticky news
+    grid $w.labcol -padx 1 -in $w.grid  -pady 1 -row $i -column 1 -sticky news
+    incr i
+    set k 1
+    proc mkPlotEntry { w k i } {
+      entry $w.plot$k -textvariable plot(fun$k)
+      entry $w.color$k -textvariable plot(col$k)
+      grid $w.plot$k -padx 10 -in $w.grid  -pady 1 -row $i -column 0 -sticky news
+      grid $w.color$k -padx 4 -in $w.grid  -pady 1 -row $i -column 1 -sticky news
+    }
+    while { $k <= $plot(numberPlots) } { mkPlotEntry $w $i $k ; incr i ; incr k}
+   }
+
+proc calculatePlot { win fun  nsteps } {
+  #  global xmin xmax  ymax ymin
+    makeLocal $win xmin xmax  ymax ymin
+    set h0 [expr {($xmax - $xmin)/double($nsteps )}]
+    set x0 $xmin
+    set res ""
+    set limit [expr {100 * (abs($ymax)> abs($ymin) ? abs($ymax) : abs($ymin))}]
+    while { $x0 < $xmax } {
+	set lastx0 $x0
+	append res " " [calculatePlot1 $win $x0 $h0 $fun $limit]
+	if { $x0 <= $lastx0 }	{
+	    # puts "x0=$x0,($lastx0)"
+	    set x0 [expr {$x0 + $h0/4}]
+	    #error "how is this?"
+	}
+    }
+    # puts "plength=[llength $res]"
+    return $res
+}
+
+
+#
+ #-----------------------------------------------------------------
+ #
+ # calculatePlot1 --   must advance x0 in its caller
+ #
+ #  Results: one connected line segment as "x0 y0 x1 y1 x2 y2 .."
+ #
+ #  Side Effects: must advance x0 in its caller
+ #
+ #----------------------------------------------------------------
+#
+proc  calculatePlot1 { win x0 h0 fun  limit } {
+     #puts "calc:$win $x0 $h0 $limit $fun"
+    makeLocal $win xmax 
+    set ansx ""
+    set ansy ""
+   while { [catch { set y0 [$fun $x0] } ] && $x0 <= $xmax }  {
+       set x0 [expr {$x0 + $h0}] }
+    if { $x0 > $xmax } {
+	# puts "catching {$fun $x0}"
+	uplevel 1 set x0 $x0
+	return ""
+    }
+    set ans "$x0 $y0"
+    set delta 0
+    set littleLimit [expr {$limit/50.0 }]
+    set veryLittleLimit [expr {$littleLimit * 10}]
+    # now have one point..
+    # this is really set below for subsequent iterations.
+    set count 10
+    set heps [expr {$h0/pow(2,6)}]
+    set h2 [expr {$h0 *2 }]
+    set ii 0
+    set x1 [expr {$x0 + $h0}]
+    while { $x1 <= $xmax  && $ii < 5000 } {
+	    # puts $x1
+	    incr ii
+	if { [catch { set y1 [$fun $x1] } ] } {
+	    	#puts "catching1 {$fun $x1}"
+	    if { $count > 0 } {
+		# try a shorter step.
+		set x1 [expr {($x1 -$x0)/2 + $x0}]
+		incr count -1
+		continue
+	    } else {
+		uplevel 1 set x0 [expr {$x0 + $heps}]
+		return [list $ansx $ansy]
+	    }
+	}
+	# ok have x1,y1
+	# do this on change in slope!! not change in limit..
+
+	set nslope [expr {($y1-$y0)/($x1-$x0)}]
+	catch { set delta [expr {($slope * $nslope < 0 ? abs($slope-$nslope) : .1*abs($slope-$nslope))}]} 
+	# catch { set delta [expr {abs($slope - ($y1-$y0)/($x1-$x0))}] }
+	
+	if { $count > 0 && (abs($y1 - $y0) > $h2 || $delta > $h2)  && (0 || abs($y1) < $littleLimit)
+	        } {
+		    #puts "too  big $y1 [expr {abs($y1-$y0)}] at $x1"
+		    set x1 [expr {($x1 -$x0)/2 + $x0}]
+		incr count -1
+	         continue 
+             } elseif { abs($y1) > $limit || abs($y1-$y0) > $limit
+		    || $delta > $littleLimit } {
+		incr ii
+		if { $count == 0 } {
+		    uplevel 1 set x0 [expr {$x0 + $heps}]
+		    return [list $ansx $ansy]
+		} else {
+		    
+		    set x1 [expr {($x1 -$x0)/2 + $x0}]
+		    incr count -1
+		    continue
+		}
+	    } else {
+    	 if {   abs($y1-$y0) > $limit/4} {
+	     
+	    # puts "x0=$x0,x1=$x1,y0=$y0,y1=$y1"
+	    uplevel 1 set x0 $x1
+	    return [list $ansx $ansy]
+	}
+
+	
+		# hopefully common case!!
+		# puts "got it: $x1,$y1,"
+	        lappend ansx $x1
+	        lappend ansy $y1
+		#append ans " $x1 $y1"
+              	set slope [expr {($y1-$y0)/($x1-$x0)} ]
+		set x0 $x1
+		set y0 $y1
+	set x1 [expr {$x0 + $h0}]
+		set count 4
+	    }
+	}
+	uplevel 1 set x0 $x1
+	return [list $ansx $ansy]
+    }
+
+    
+
+	
+
+proc setup_xf { vars form } {
+    set s [sparse $form ] 
+    proc _xf [list $vars ]  "return \[ expr { $s } \]"
+}
+
+
+#
+ #-----------------------------------------------------------------
+ #
+ # nextColor --  get next COLOR and advance the curveNumber
+ #
+ #  Results: a color
+ #
+ #  Side Effects: the local variable for WIN called curveNumber is incremented
+ #
+ #----------------------------------------------------------------
+#
+proc nextColor { win } {
+    makeLocal $win linecolors 
+    if { [catch { set i [oget $win curveNumber] } ] } { set i -1 }
+    set color [lindex $linecolors [expr {[incr i]%[llength $linecolors]}]]
+    oset $win curveNumber $i
+    return $color
+}
+    
+
+proc plot2d {args } {
+    #puts "args=$args"
+    set win [apply mkPlot2d $args]
+    replot2d $win
+    return $win
+}
+
+proc replot2d {win } {
+    global printOption axisGray plot2dOptions
+    linkLocal $win xfundata data
+    foreach v $data {
+	if { "[assq [lindex $v 0] $plot2dOptions notthere]" != "notthere" } {
+	    oset $win [lindex $v 0] [lindex $v 1]
+	}
+    }
+    makeLocal $win xfun nsteps c linecolors xaxislabel yaxislabel  autoscale
+    set xfundata ""
+    foreach v [sparseList $xfun] {
+	proc _xf {  x  } "return \[expr { $v } \]"
+	regsub "\\$" $v "" label
+	lappend xfundata [list label $label] \
+	  [linsert [calculatePlot $win _xf $nsteps]  \
+		0 xversusy]
+    }
+
+    # in case only functions and no y autoscale dont bother.
+    if { "$data" != "" || [lsearch $autoscale y]>=0  } {
+	set ranges [plot2dGetDataRange [concat $data $xfundata]]	
+	foreach {v k} [eval plot2dRangesToRadius $ranges] {
+	    if { [lsearch $autoscale [string index $v 1] ] >= 0 } {
+		oset $win [string range $v 1 end] $k
+	    }
+	}
+    }
+    
+    setUpTransforms $win 1.0
+    set rtosx rtosx$win ; set rtosy rtosy$win
+    $c del axes
+    $c create line [$rtosx 0 ] [$rtosy -1000] [$rtosx 0] [$rtosy 1000] -fill $axisGray -tags axes
+    $c create line [$rtosx -1000] [$rtosy 0] [$rtosx 1000] [$rtosy 0] -fill $axisGray -tags axes
+    axisTicks $win $c
+
+    if { "$xfun" != "" } {
+	 oset $win maintitle [concat list "Plot of y = \[oget $win xfun\]" ]
+    }
+    $c del path
+    $c del label
+    oset  $win curveNumber -1
+    redraw2dData $win -tags path
+    $c create text    [expr {[$rtosx 0] + 10}] [expr {[$rtosy [oget $win ymax]] +20}] -text [oget $win yaxislabel] -anchor nw
+    $c create text     [expr {[$rtosx [oget $win xmax]] -20}] [expr {[$rtosy 0] - 10}] -text [oget $win xaxislabel] -anchor se
+
+    
+}
+
+
+
+#
+ #-----------------------------------------------------------------
+ #  Should change name to plotData since works for 3d to now..
+ # plot2dData --  create WIN and plot 2d OR 3d DATA which is a list of 
+ #  data sets.  Each data set must begin with xversusy or againstIndex
+ #  In the first case the data set looks like:
+ #       { xversusy {x1 x2 ...xn} {y1 ... yn yn+1 ... ym} }
+ #  and will be plotted as m/n curves : (x1,y1) (x2,y2) .. (xn,yn)
+ #  and (x1,yn+1) (x2,yn+2) ..
+ #  In the againstIndex case the x values are replace by the indices
+ #  0,1,2,... [length $yvalues]-1 
+ #  Results: none
+ #
+ #  Side Effects: curves draw
+ #
+ #----------------------------------------------------------------
+#
+proc plot2dData { win data args } {
+   clearLocal $win
+    #puts "data=$data, [regexp plot2d $data junk ]"
+    if { [regexp plot2d $data junk] } {
+  # eval plot2d $args -windowname $win  [plot2dGetRanges $data] -xfun [list {}] -data [list $data]
+   eval plot2d $args -windowname $win   -xfun [list {}] -data [list $data]	
+    } else {
+	# puts data=$data
+	set com [concat \
+		plot3d $args -windowname $win -zfun {{}} -data [lrange $data 1 end]]
+	# puts com=$com
+	eval $com
+    }
+ }
+
+
+
+proc plot2dGetDataRange { data } {
+    set rangex ""
+    set rangey ""
+     #puts "data=$data"
+    set extra ""
+    foreach d $data {
+	#puts first=[lindex $d 0]
+      if { [catch { 	
+	switch -exact -- [lindex $d 0] {
+	   xversusy {
+	       foreach { xx yy } [lrange $d 1 end] {
+		  # puts "hi xx=[llength $xx],yy=[llength $yy]"
+		   if { [llength $xx] > 0 } {
+		       set rangex [minMax $xx $rangex]
+		       set rangey [minMax $yy $rangey]
+		   }
+	       }
+	       #puts "rangex=$rangex,rangey=$rangey"
+	   }
+	   againstIndex {
+	       set rangex [minMax [list 0 [llength [lindex $d 1]]] $rangex]
+	       set rangey [minMax [lindex $d 1] $rangey]
+	   }
+	   default {
+	       set vv [lindex $d 0]
+	       if { [lsearch {xrange yrange   } $vv] >= 0 } {
+		   set radius [expr {([lindex $d 2] -[lindex $d 1])/2.0 }]
+		   set center [expr {([lindex $d 2] +[lindex $d 1])/2.0 }]
+		   set var [string range $vv 0 0]
+		   lappend extra -${var}radius $radius -${var}center $center
+	       }
+	       if { [lsearch bargraph $vv] >= 0 } {
+		       set rangey [minMax 0 $rangey]
+		   }
+
+
+	       if { [lsearch {xradius yradius xcenter ycenter } $vv] >= 0 } {
+		   lappend extra -$vv [list [lindex $d 1]]
+	       }
+
+	    }
+       }
+   } errmsg ] } {
+       set com [list error "bad data: [string range $d 0 200].." $errmsg]
+       after 1 $com
+   }
+ }
+   list $rangex $rangey $extra
+}
+
+
+
+proc plot2dRangesToRadius  { rangex rangey extra } {
+   set ranges ""
+  # puts "extra=$extra"
+   foreach u { x y } {
+       if { "[assoc -[set u]radius $extra]" == "" } {
+	   desetq "min max" [set range$u]
+	   if { "$min" == "$max" } {
+	       set min [expr {$min - .5}]
+	       set max [expr {$max + .5}]
+	   }
+	   #puts "$u has $min,$max"
+	   # use 1.7 to get a bit bigger radius than really necessary.
+	   if { "$max" != "" } {
+	   
+	       lappend extra -[set u]radius [expr {($max-$min)/1.7}] \
+		       -[set u]center [expr {($max+$min)/2.0}]
+	   }
+   }
+ }
+ # puts "extra=$extra"
+ return $extra
+}
+   
+
+proc redraw2dData { win  args } {
+   makeLocal $win c linecolors data xfundata errorbar
+   set tags [assoc -tags $args {} ]
+   set rtosx rtosx$win ; set rtosy rtosy$win  
+   set i -1
+   set label _default
+   append data " " $xfundata
+    
+   #puts "data=$data"
+   foreach d $data {
+       set type [lindex $d 0]
+       switch  $type {
+	   xversusy {
+	        #puts "starting .. [oget $win curveNumber]"
+	       set curvenumber [oget $win curveNumber]
+	       # the data can be multiple lists and each list
+	       # will not be line connected to previous
+	       foreach {xvalues yvalues} [lrange $d 1 end] {
+		   # puts "xvalues=$xvalues"
+		   #puts "here:$curvenumber,[oget $win curveNumber]"
+		   oset $win curveNumber $curvenumber
+		   set n [expr {[llength $xvalues] -1}]
+		   while { [llength $yvalues] > 0 } {
+		       set ans ""
+		       set color [nextColor $win]
+		       catch { set color [oget $win color] }
+		       
+		       if { [info exists didLabel([oget $win curveNumber])] } {
+			   set label "" } else { set didLabel([oget $win curveNumber]) 1
+		       }
+		       set errorbar [oget $win errorbar]
+		       # puts "errorbar=$errorbar"
+		       if { $errorbar != 0 } {
+			   set j 0
+			  # puts "xvalues=$xvalues,yvalues=$yvalues"
+     		       for { set i 0 } { $i <= $n } {incr i} {
+			       set x [lindex $xvalues $i]
+			       set y1 [lindex $yvalues [expr {$i * 2}]]
+			       set y2 [lindex $yvalues [expr { $i * 2 +1}]]
+			   if { 1 } {
+			      # puts "x=$x,y1=$y1,y2=$y2"
+			       set xx [$rtosx $x]
+			       set y1 [$rtosy $y1]
+			       set y2 [$rtosy $y2]
+			       $c create line [expr {$xx - $errorbar}] $y1 [expr {$xx +$errorbar}] $y1 $xx $y1 $xx $y2 [expr {$xx -$errorbar}] $y2 [expr {$xx + $errorbar}] $y2  -tags [list [concat $tags line[oget $win curveNumber]]]  -fill $color
+			   }
+		       }
+		       
+			       
+		       set yvalues [lrange $yvalues [llength $xvalues] end]
+		       } else { 
+
+		       foreach x $xvalues y [lrange $yvalues 0 $n] {
+			   append ans "[$rtosx $x] [$rtosy $y] "
+		       
+		       }
+
+		       drawPlot $win [list $ans] -tags [list [concat $tags line[oget $win curveNumber]]]  -fill $color -label $label
+		       }
+		       set label _default
+
+		       set yvalues [lrange $yvalues [llength $xvalues] end]
+		   }
+
+	     } }
+	   againstIndex {
+               set color [nextColor $win]
+	       set ind 0
+	       set ans ""
+	       foreach y [lindex $d 1] {
+		   append ans "[$rtosx $ind] [$rtosy $y] "
+		   incr ind
+	       }
+	       
+	       drawPlot $win [list $ans] -tags \
+		       [list [concat $tags line[oget $win curveNumber]]] \
+		       -fill $color -width .2 -label $label
+	       set label _default
+
+#	       eval $c create line $ans -tags \
+#		        [list [concat $tags line[oget $win curveNumber]]] \
+#		       -fill $color -width .2
+	   }
+	   label {
+	       set label [lindex $d 1]
+	   }
+	   default {
+
+	       # puts "$type,[lindex $d 1]"
+	       if { [lsearch { xfun color plotpoints linecolors pointsize nolines bargraph errorbar maintitle
+	       labelposition
+	       xaxislabel yaxislabel } $type] >= 0 } {
+		   # puts "setting oset $win $type [lindex $d 1]"
+		   oset $win $type [lindex $d 1]
+	       } elseif { "$type" == "text" } {
+		   desetq "x y text" [lrange $d 1 end]
+		   $c create text [$rtosx $x] [$rtosy $y] -anchor nw -text $text -tags "text all" -font times-roman
+	       }
+
+	   }
+
+       }
+   }
+
+}
+
+proc plot2dDrawLabel { win label color } {
+    makeLocal $win c labelposition
+    #puts "$win $label $color"
+    if { "$label" == ""} {return }
+    set bb [$c bbox label]
+    desetq "a0 b0" $labelposition
+    if { "$bb" == "" } { set bb "$a0 $b0 $a0 $b0" }
+    desetq "x0 y0 x1 y1" $bb
+    set leng  15
+    set last [$c create text [expr {$a0 +$leng +4}] \
+	    [expr {2 + $y1}] \
+	    -anchor nw       -text "$label" -tags label]
+    desetq "ux0 uy0 ux1 uy1" [$c bbox $last]
+    $c create line $a0 [expr {($uy0+$uy1) /2}] [expr {$a0 +$leng}] [expr {($uy0+$uy1) /2}]   -tags "label" -fill $color
+}
+	
+
+proc RealtoScreen { win listPts } {
+    set rtosx rtosx$win ; set rtosy rtosy$win  
+    set ans ""
+    if { [llength [lindex $listPts  0]] != 1 } {
+	foreach v $listPts {
+	    append ans " {"
+	    append ans [RealtoScreen $win $v]
+	    append ans "}"
+	}
+    }    else {
+	foreach {x y } $listPts {
+	    append ans " [$rtosx $x] [$rtosy $y]"
+	}
+    }
+    return $ans
+}
+
+proc drawPlot {win listpts args } {
+    makeLocal $win  c nolines plotpoints  pointsize bargraph
+    
+    # puts ll:[llength $listpts]
+    set tags [assoc -tags $args ""]
+    if { [lsearch $tags path] < 0 } {lappend tags path}
+    set fill [assoc -fill $args black]
+    set label [assoc -label $args ""]
+    if { "$label" == "_default" } {
+	set label line[oget $win curveNumber]
+    }
+
+    catch { set fill [oget $win color] }
+    
+    if { $nolines == 1 && $plotpoints == 0 && $bargraph == 0} {
+	set plotpoints 1
+    }
+
+    catch { 
+    foreach pts $listpts {
+	if { $bargraph } {
+	    set rtosy rtosy$win
+	    set rtosx rtosx$win
+	    set width [expr {abs([$rtosx $bargraph] - [$rtosx 0])}]
+	    set w2 [expr {$width/2.0}]
+	    # puts "width=$width,w2=$w2"
+	    set ry0 [$rtosy 0]
+	    foreach { x y } $pts {
+		$c create rectangle [expr {$x-$w2}] $y  [expr {$x+$w2}] \
+			$ry0 -tags $tags -fill $fill }
+	    } else {
+		if { $plotpoints } {
+		    set im [getPoint $pointsize $fill]
+		    
+		    # there is no eval, so we need this.
+		    if { "$im" != "" } {
+			foreach { x y } $pts {
+			$c create image $x $y -image $im -anchor center \
+				-tags "$tags point"
+			}
+		    } else {
+		    foreach { x y } $pts {
+			$c create oval [expr {$x -$pointsize}] \
+				[expr {$y -$pointsize}] [expr {$x +$pointsize}] \
+				[expr {$y +$pointsize}] -tags $tags \
+				-fill $fill -outline {}
+			
+		    }
+		}
+		}
+		
+		if { $nolines == 0 } {
+		    set n [llength $pts]
+		    set i 0
+		    set res "$win create line "
+		    #puts npts:[llength $pts]
+		    if { $n >= 6 } {
+			eval $c create line $pts  	-tags [list $tags] -width .2 -fill $fill
+		    }
+		}
+	    }
+	    
+	}
+    }
+    plot2dDrawLabel $win $label $fill
+}
+
+    
+
+proc drawPointsForPrint { c } {
+    global ws_openMath
+    foreach v [$c find withtag point] {
+	set tags [ldelete point [$c gettags $v]]
+	desetq "x y" [$c coords $v]
+	
+	
+	desetq "pointsize fill" $ws_openMath(pointimage,[$c itemcget $v -image])
+	catch { 
+	    $c create oval [expr {$x -$pointsize}] \
+		    [expr {$y -$pointsize}] [expr {$x +$pointsize}] \
+		    [expr {$y +$pointsize}] -tags $tags \
+				-fill $fill -outline {}
+         $c delete $v			
+	}
+
+
+    }
+
+}
+
+array set ws_openMath { bitmap,disc4 {#define disc4_width 4
+#define disc4_height 4
+static unsigned char disc4_bits[] = {
+    0x06, 0x0f, 0x0f, 0x06};}
+    bitmap,disc6 {#define disc_width 6
+#define disc_height 6
+static unsigned char disc_bits[] = {
+    0xde, 0xff, 0xff, 0xff, 0xff, 0xde};}
+}
+ 
+
+proc getPoint { size color } {
+    global ws_openMath
+    set im ""
+    if { ![catch { set im $ws_openMath(pointimage,$size,$color) }] } {
+	return $im
+    }
+    catch { set data $ws_openMath(bitmap,disc[expr {$size * 2}]) 
+    set im [image create bitmap -data $data -foreground $color]
+    set ws_openMath(pointimage,$size,$color) $im
+    set ws_openMath(pointimage,$im) "$size $color"
+   }
+    return $im
+}
+    
+
+
+
+## endsource plot2d.tcl
+## source plot3d.tcl
+
+###### plot3d.tcl ######
+
+
+set ws_openMath(speed) [expr {(9700.0 / (1 + [lindex [time {set i 0 ; while { [incr i] < 1000} {}} 1] 0]))}]
+
+
+set plot3dOptions { 
+    {xradius 1 "Width in x direction of the x values" }
+    {yradius 1 "Height in y direction of the y values"}
+
+    {width 500 "Width of canvas in pixels"}
+    {height 500 "Height of canvas in pixels" }
+    {xcenter 0.0 {(xcenter,ycenter) is the origin of the window}}
+    {ycenter 0.0 "see xcenter"}
+    {zcenter 0.0 "see xcenter"}
+    {zradius auto " Height in z direction of the z values"}
+    {az 60 "azimuth angle" }
+    {el 30 "elevantion angle" }
+    
+    {thetax 10.0 "ignored is obsolete: use az and el"}
+    {thetay 20.0 "ignored is obsolete: use az and el"}
+    {thetaz 30.0 "ignored is obsolete: use az and el"}
+    {flatten 0 "Flatten surface when zradius exceeded" }
+    {zfun "" "a function of z to plot eg: x^2-y^2"}
+    {data  "" "a data set of type { variable_grid xvec yvec zmatrix}
+    or {matrix_mesh xmat ymat zmat} or {grid {xmin xmax} {ymin ymax} zmatrix}"}
+    {nsteps "10 10" "steps in x and y direction"}
+    {rotationcenter "" "Origin about which rotation will be done"}
+    {zoomfactor "1.6 1.6" "Factor to zoom the x and y axis when zooming.  Zoom out will be reciprocal" }
+    {screenwindow "20 20 700 700" "Part of canvas on screen"}
+    {windowname ".plot3d" "window name"}
+}
+    
+
+## source matrix.tcl
+
+###### matrix.tcl ######
+
+# In this file a matrix is represented by a list of M*N entries together
+# with an integer N giving the number of columns: {1 0 0 1} 2  would give
+# the two by two identity
+
+proc comment {args } { }
+  set mee " } \] \[ expr { "  
+
+proc mkMultLeftExpr { mat n prefix { constant "" } } {
+#create a function body that does MAT (prefix1,prefix2,..) + constant 
+    global mee
+    set all ""
+    
+    set vars ""
+    for { set i 0} { $i < $n} {incr i} { append vars " $prefix$i" }
+    set j 0
+    set k 0
+    
+    foreach v $mat {
+	if { $j == 0 } {
+	    set ro ""
+	    # append ans ""
+	    set op ""
+	}
+        append ro " $op $v*\$$prefix$j"
+	set op "+"
+	if { $j == [expr {$n -1}] } {
+	     append ans " "
+	    if { "[lindex $constant $k]" != "" } {
+		append ro " + [lindex $constant $k] "
+	    }
+	    incr k
+	    append ans [concat \[ expr [list $ro] \]]
+	    set j -1
+	}
+	incr j
+    }
+     return [list $vars $ans]
+}
+
+proc mkMultLeftFun { mat n name { constant ""} } {
+    set expr [mkMultLeftExpr $mat $n _a $constant]
+    set bod1 [string trim [lindex $expr 1] " "]
+    set bod "return \"$bod1\""
+    proc $name [lindex $expr 0] $bod
+}
+
+proc rotationMatrix { th ph } {
+   return [list \
+	   [expr {cos($ph)*cos($th)}] [expr {- cos($ph)*sin($th)}] [expr {sin($ph)}] \
+	   [expr {sin($th)}] [expr {cos($th)}] 0.0 \
+	   [expr {- sin($ph)*cos($th)}] [expr {sin($ph)*sin($th)}] [expr {cos($ph)}]]
+}
+
+# proc rotationMatrix { thx thy thz } {
+#   return [list  \
+#  [expr { cos($thy)*cos($thz)} ]  \
+#  [expr { cos($thy)*sin($thz)} ]  \
+#  [expr { sin($thy)} ]  \
+#  [expr { sin($thx)*sin($thy)*cos($thz)-cos($thx)*sin($thz)} ]  \
+#  [expr { sin($thx)*sin($thy)*sin($thz)+cos($thx)*cos($thz)} ]  \
+#  [expr { -sin($thx)*cos($thy)} ]  \
+#  [expr { -sin($thx)*sin($thz)-cos($thx)*sin($thy)*cos($thz)} ]  \
+#  [expr { -cos($thx)*sin($thy)*sin($thz)+sin($thx)*cos($thz)} ]  \
+#  [expr { cos($thx)*cos($thy)} ] ]
+# }
+
+proc rotationMatrix { thx thy thz } {
+    return \
+ [list  \
+ [expr { cos($thy)*cos($thz) } ] \
+ [expr { cos($thy)*sin($thz) } ] \
+ [expr { sin($thy) } ] \
+ [expr { sin($thx)*sin($thy)*cos($thz)-cos($thx)*sin($thz) } ] \
+ [expr { sin($thx)*sin($thy)*sin($thz)+cos($thx)*cos($thz) } ] \
+ [expr { -sin($thx)*cos($thy) } ] \
+ [expr { -sin($thx)*sin($thz)-cos($thx)*sin($thy)*cos($thz) } ] \
+ [expr { sin($thx)*cos($thz)-cos($thx)*sin($thy)*sin($thz) } ] \
+ [expr { cos($thx)*cos($thy) } ] ]
+}
+
+# cross [a,b,c] [d,e,f] == [B*F-C*E,C*D-A*F,A*E-B*D]
+# cross_product([a,b,c],[d,e,f]):=[B*F-C*E,C*D-A*F,A*E-B*D]
+# cross_product(u,v):=sublis([a=u[1],b=u[2],c=u[3],d=v[1],e=v[2],f=v[3]],[B*F-C*E,C*D-A*F,A*E-B*D]);
+# the rotation by azimuth th, and elevation ph
+# MATRIX([COS(TH),SIN(TH),0],[-COS(PH)*SIN(TH),COS(PH)*COS(TH),SIN(PH)],
+#	    [SIN(PH)*SIN(TH),-SIN(PH)*COS(TH),COS(PH)]);
+
+proc rotationMatrix { th ph {ignore {} } } {
+    return \
+[list \
+[	    expr {cos($th)   } ]\
+[expr {sin($th)   } ]\
+0 \
+[expr {-cos($ph)*sin($th)   } ]\
+[expr {cos($ph)*cos($th)   } ]\
+[expr {sin($ph)   } ]\
+[expr {sin($ph)*sin($th)   } ]\
+[expr {-sin($ph)*cos($th)   } ]\
+[expr {cos($ph)   } ]]
+}
+
+proc setMatFromList {name lis n} {
+    set i 1
+    set j 1
+    foreach v $lis {
+	uplevel 1 set [set name]($i,$j) $v
+	if { $j == $n } {set j 1; incr i} else { incr j}
+}   }
+
+proc matRef { mat cols i j } { [lindex $mat [expr {$i*$cols + $j}]] }
+proc matTranspose { mat cols } {
+    set j 0
+    set m [expr {[llength $mat ] / $cols}]
+    while { $j < $cols} {
+	set i 0
+	while { $i < $m } {
+	    append ans " [lindex $mat [expr {$i*$cols + $j}]]"
+	    incr i
+	}
+	incr j
+    }
+    return $ans
+}
+
+
+proc matMul { mat1 cols1 mat2 cols2 } {
+    mkMultLeftFun $mat1 $cols1 __tem
+    set tr [matTranspose $mat2 $cols2]
+    set rows1 [expr {[llength $mat1] / $cols1}]
+    #puts "tr=$tr"
+    set upto [llength $tr]
+    set j 0
+    set ans ""
+    set i 0
+    while { $j < $cols2  } {
+	append ans " [eval __tem [lrange $tr $i [expr {$i+$cols1 -1}]]]"
+	incr i $cols1
+	incr j
+    }
+ #   return $ans
+   # puts "matTranspose $ans $rows1"
+    return [matTranspose $ans $rows1]
+}
+
+
+
+proc invMat3 { mat } {
+    setMatFromList xx $mat 3
+    set det [expr { double($xx(1,1))*($xx(2,2)*$xx(3,3)-$xx(2,3)*$xx(3,2))-$xx(1,2)* \
+	    ($xx(2,1)*$xx(3,3)-$xx(2,3)*$xx(3,1))+$xx(1,3)*($xx(2,1)*$xx(3,2)\
+	    -$xx(2,2)*$xx(3,1)) }]
+    
+    return [list   [expr { ($xx(2,2)*$xx(3,3)-$xx(2,3)*$xx(3,2))/$det}] \
+	    [expr { ($xx(1,3)*$xx(3,2)-$xx(1,2)*$xx(3,3))/$det}] \
+	    [expr { ($xx(1,2)*$xx(2,3)-$xx(1,3)*$xx(2,2))/$det}] \
+	    \
+	    [expr { ($xx(2,3)*$xx(3,1)-$xx(2,1)*$xx(3,3))/$det}] \
+	    [expr { ($xx(1,1)*$xx(3,3)-$xx(1,3)*$xx(3,1))/$det}] \
+	    [expr { ($xx(1,3)*$xx(2,1)-$xx(1,1)*$xx(2,3))/$det}] \
+	    \
+	    [expr { ($xx(2,1)*$xx(3,2)-$xx(2,2)*$xx(3,1))/$det}] \
+	    [expr { ($xx(1,2)*$xx(3,1)-$xx(1,1)*$xx(3,2))/$det}] \
+	    [expr { ($xx(1,1)*$xx(2,2)-$xx(1,2)*$xx(2,1))/$det}]]
+}
+
+
+proc vectorOp { a op b} {
+    set i [llength $a]
+    set k 0
+    set ans [expr [list [lindex $a 0]  $op [lindex $b 0]]]
+    while { [incr k] < $i } {
+	lappend ans  [expr  [list [lindex $a $k] $op [lindex $b $k]]]
+    }
+    return $ans
+}
+## endsource matrix.tcl
+
+proc transformPoints { pts fun } {
+  set ans ""
+  foreach { x y z } $pts {
+      append ans " "
+      append ans [$fun $x $y $z]
+  }
+  return $ans
+}
+
+proc calculatePlot3d {win fun  nx ny } {
+    global plot3dMeshes$win
+    set meshes  plot3dMeshes$win
+    makeLocal $win xradius xmin yradius ymin zradius zcenter flatten
+    
+    set stepx [expr { 2*$xradius / double($nx)}]
+    set stepy [expr { 2*$yradius / double($ny)} ]
+    set i 0
+    set j 0
+    set zmax -1000000000
+    set zmin 1000000000
+    # check if zradius is a number
+    set dotruncate [expr ![catch {expr {$zradius + 1} }]]
+    if { $dotruncate } {
+	if { $flatten } { set dotruncate 0 }
+	set zzmax [expr {$zcenter + $zradius}]
+	set zzmin [expr {$zcenter - $zradius}]
+	#puts "zzmax=$zzmax,$zzmin"
+   } else { set flatten 0 }
+
+    catch { unset  $meshes }
+    set k 0
+    for {set i 0} { $i <= $nx } { incr i} {
+	set x [expr { $xmin + $i * $stepx }]
+	for {set j 0} { $j <= $ny } { incr j} {
+	    set y [expr { $ymin + $j *$stepy }]
+	   if { [catch {  set z [$fun $x $y] }] } {
+	       set z nam
+	   } elseif { $dotruncate  &&  ($z > $zzmax || $z < $zzmin) } {
+	       set z nam
+
+	    } else {
+		if { $flatten } {
+		    if { $z > $zzmax } { set z $zzmax } elseif {
+			$z < $zzmin } { set z $zzmin }}
+			
+	       if { $z < $zmin }  { set zmin $z } elseif {
+		$z > $zmax } { set zmax $z }
+		if { $j != $ny && $i != $nx } {
+		    set [set meshes]($k) \
+                      "$k [expr { $k+3 }] [expr { $k+3+($ny+1)*3 }] \
+		      [expr { $k+($ny+1)*3 }]"} else {
+		  # set plot3dMeshes($k) ""
+	      }
+	  }
+	      incr k 3
+	    append ans " $x $y $z"
+	}
+    }
+    oset $win zmin $zmin
+    oset $win zmax $zmax
+    oset $win points $ans
+    oset $win nx $nx
+    oset $win ny $ny
+    oset $win colorfun plot3dcolorFun
+    addAxes $win
+    setupPlot3dColors $win
+}
+
+proc calculatePlot3data {win fun  nx ny } {
+# calculate the 3d data from function:
+    makeLocal $win xradius xmin xmax ymax yradius ymin zradius zcenter flatten
+    
+    set rowx [linspace $xmin $xmax $nx]
+    set rowy [linspace $ymin $ymax $ny]
+    foreach  y $rowy {
+	set row ""
+	foreach x $rowx {
+	    if { [catch {  set z [$fun $x $y] }] } {
+		set z nam
+	    }
+	    lappend row $z
+	}
+	lappend matrix $row
+    }
+    return [list variable_grid $rowx $rowy $matrix ]
+    
+}
+
+
+proc addAxes { win } {
+    #global plot3dPoints plot3dMeshes xradius yradius xcenter ycenter
+    global [oarray $win] plot3dMeshes$win
+    linkLocal $win lmesh
+    makeLocal $win   xradius yradius xcenter ycenter  points zmax zcenter zmin
+    set meshes plot3dMeshes$win
+    set ll [llength $points]
+
+    # puts "oset $win  axisstart  $ll"
+    oset $win  axisstart  $ll
+    set nx2 5
+    set ny2 5
+    set xstep [expr { 1.2 * $xradius/double($nx2) }]
+    set ystep [expr { 1.2 * $yradius/double($ny2) }]
+    set nz2 $ny2
+
+    set ans " "
+    set x0 $xcenter
+    set y0 $ycenter
+    set z0 $zcenter
+
+    set k $ll
+    for { set i 0 } { $i < $nx2 } { incr i } {
+	append ans "[expr {$x0 +$i * $xstep}] $y0 $z0 "
+	lappend lmesh [list $k [incr k 3]]
+	#set [set meshes]($k) "$k [incr k 3]"
+    }
+    append ans "[expr {$x0 +$nx2 * $xstep}] $y0 $z0 "
+    incr k 3
+    # set plot3dMeshes($k) ""
+
+    for { set i 0 } { $i < $ny2 } { incr i } {
+	append ans "$x0 [expr {$y0 +$i * $ystep}] $z0 "
+	lappend lmesh [list $k [incr k 3]]
+	#set [set meshes]($k) "$k [incr k 3]"
+    }
+    append ans "$x0 [expr {$y0 +$ny2 * $ystep}] $z0 "
+    incr k 3
+    # set $meshes($k) ""
+
+    set zstep [expr {1.2 * $zmax/double($nz2)}]
+    if { $zstep < $ystep } { set zstep $ystep }
+    
+    for { set i 0 } { $i < $ny2 } { incr i } {
+	append ans "$x0 $y0 [expr {$z0 +$i * $zstep}] "
+	# puts "set [set meshes]($k) \"$k [incr k 3]\""
+	lappend lmesh [list $k [incr k 3]]
+	# set [set meshes]($k) "$k [incr k 3]"
+    }
+    append ans "$x0 $y0 [expr {$z0 +$nz2 * $zstep}] "
+    incr k 3
+    # puts "ans=$ans"
+    append [oloc $win points] $ans
+    
+   # set $meshes($k) ""
+
+}
+
+proc addBbox { win } {
+    global plot3dMeshes$win
+    makeLocal $win xmin xmax ymin ymax zmin zmax  cmap
+    linkLocal $win points lmesh 
+    set ll [llength $points]
+     append points " $xmin $ymin $zmin \
+	    $xmax $ymin $zmin \
+            $xmin $ymax $zmin \
+            $xmax $ymax $zmin \
+            $xmin $ymin $zmax \
+	    $xmax $ymin $zmax \
+            $xmin $ymax $zmax \
+            $xmax $ymax $zmax "
+    foreach  { a b } { 0 1 0 2 2 3 3 1
+                       4 5 4 6 6 7 7 5
+    0 4 1 5 2 6 3 7  }  {
+	set k [expr {$a*3 + $ll}]
+	set l [expr {$b*3 + $ll}]
+	# set plot3dMeshes${win}($k) [list $k $l]
+	lappend lmesh [list $k $l]
+    }
+    lappend lmesh [list $ll]
+    oset $win $cmap,[list $ll [expr {$ll + 3}]] red
+    oset $win $cmap,[list $ll [expr {$ll + 6}]] blue
+    oset $win $cmap,[list $ll [expr {$ll + 12}]] green
+    
+    oset $win special($ll) "drawOval [oget $win c] 3 -fill red -tags axis"
+}
+
+proc drawOval { c radius args } {
+    set ll [llength $args]
+    set x [lindex $args [expr {$ll -2}]]
+    set y [lindex $args [expr {$ll -1}]]
+    set rest [lrange $args 0 [expr {$ll -3}]]
+    set com [concat $c create oval [expr {$x - $radius}]  [expr {$y - $radius}] [expr {$x + $radius}]  [expr {$y + $radius}] $rest]
+    eval $com
+}
+    
+
+proc plot3dcolorFun {win z } {
+    makeLocal $win zmin zmax
+    set ncolors 180
+    set tem [expr {(180/$ncolors)*round(($z - $zmin)*$ncolors/($zmax - $zmin+.001))}]
+    #puts "tem=$tem,z=[format %3g $z],[format "#%.2x%.2x%.2x" 50 50 $tem]"
+    return [format "#%.2x%.2x%.2x" [expr {180 -$tem}] [expr {240 - $tem}] $tem]
+}
+
+proc setupPlot3dColors { win } {
+    upvar #0 [oarray $win] wvar
+    # the default prefix for cmap
+    set wvar(cmap) c1
+    set k 0
+    makeLocal $win colorfun points
+    foreach { x y z } $points {
+	catch { set wvar(c1,$k) [$colorfun $win $z] }
+	incr k 3 
+    }
+}
+
+proc calculateRotated { win } {
+    set pideg [expr {3.14159/180.0}]
+    linkLocal $win scale
+    makeLocal $win az el rotationcenter xradius zradius yradius
+    set rotmatrix [rotationMatrix [expr {$az * $pideg }] \
+	    [expr {$el * $pideg }] \
+      ]
+
+    # shrink by .2 on z axis
+    # set fac [expr  {[vectorlength $xradius $yradius] / (sqrt(2) * $zradius)}]
+
+    set rotmatrix [ matMul  $rotmatrix 3 $scale 3 ]
+    set tem [matMul $scale 3 $rotationcenter 1]
+    
+    mkMultLeftFun  $rotmatrix 3 _rot$win
+    set rot _rot$win
+    set ans ""
+    # puts "points=[oget $win points]"
+    if { "$rotationcenter" != "" } {
+	#puts "rotationcenter = $rotationcenter"
+	set constant [vectorOp $tem - [eval $rot $rotationcenter]]
+	mkMultLeftFun  $rotmatrix 3 _rot$win $constant
+    }
+    foreach { x y z } [oget $win points] {
+	if { [catch { append ans " [$rot $x $y $z]" } ] } {
+	    append ans "  nam nam nam " }
+	}
+    oset $win rotatefun $rot
+    oset $win rotated $ans 
+}
+
+proc getOrderedMeshIndices { win } {
+ #   global  plot3dMeshes$win
+ #    set meshes plot3dMeshes$win
+    linkLocal $win lmesh
+    # puts "array names $meshes =[array names $meshes ]"
+    set pts [oget $win rotated]
+    set i 0
+    foreach tem $lmesh {
+        set k [llength $tem]
+       if { [catch {
+        if {  $k == 4 } {
+	    set z [expr { ([lindex $pts [expr {[lindex $tem 0] +2}]] \
+		    +[lindex $pts [expr {[lindex $tem 1] +2}]] \
+		    + [lindex $pts [expr {[lindex $tem 2] +2}]] \
+		    + [lindex $pts [expr {[lindex $tem 3] +2}]])/4.0 }]
+	} elseif { $k == 2 } {
+            set z [expr { ([lindex $pts [expr {[lindex $tem 0] +2}]] \
+		    +[lindex $pts [expr {[lindex $tem 1] +2}]])/2.0 }]
+	} else {
+	    set z 0
+            foreach w $tem {
+		set z [expr {$z + [lindex $pts [expr {$w + 2}]] } ]
+		
+	    }	
+	    set z [expr { $z/double($k)}]
+        }
+	append pp($z) "$i "
+	incr i
+	
+         } ]} {
+	    set lmesh [lreplace $lmesh $i $i]
+	 }
+    }
+    #puts "pp=[array get pp]"
+    
+    set tem [lsort -real [array names pp]]
+    set ans ""
+    foreach v $tem {
+	append ans $pp($v)
+    }
+    oset $win meshes $ans
+    return
+}
+
+# we dont want xy plane to fill whole page!
+proc setUpTransforms3d { win } {
+    global screenwindow
+    #set scr $screenwindow
+    setUpTransforms $win .7
+    # set screenwindow $scr
+
+}
+
+proc setUpTransforms3d { win } {
+    global screenwindow
+    #set scr $screenwindow
+    # setUpTransforms $win .7
+    # set screenwindow $scr
+    linkLocal $win scale
+        makeLocal $win xcenter ycenter xradius yradius c zmin zmax xmin xmax ymin ymax zradius
+    set fac .5
+
+    set delx [$c cget -width]
+    set dely [$c cget -height]
+    set f1 [expr {(1 - $fac)/2.0}]
+    set scale [list [expr {1.5/($xradius)}] 0 0 0 [expr {1.5/($yradius)}] 0 0 0 [expr {1.5/($zradius)}] ]
+    set x1 [expr {$f1 *$delx}]
+    set y1 [expr {$f1 *$dely}]
+    set x2 [expr {$x1 + $fac*$delx}]
+    set y2 [expr {$y1 + $fac*$dely}]
+    # set xmin [expr {($xcenter - $xradius) * 1.5/ ($xradius)}]
+    # set ymin [expr {($ycenter - $yradius) * 1.5/ ($yradius)}]
+    # set xmax [expr {($xcenter + $xradius) * 1.5/ ($xradius)}]
+    # set ymax [expr {($ycenter + $yradius) * 1.5/ ($yradius)}]
+    #puts "RANGES=$xmin,$xmax $ymin,$ymax $zmin,$zmax"
+    desetq "xmin ymin" [matMul $scale 3 "$xmin $ymin 0" 1]
+    desetq "xmax ymax" [matMul $scale 3 "$xmax $ymax 0" 1]
+    #puts "RANGES=$xmin,$xmax $ymin,$ymax $zmin,$zmax"
+    # set transform [makeTransform "$xmin $ymin $x1 $y2" "$xmin $ymax $x1 $y1 " "$xmax $ymin $x2 $y2"]
+   # desetq "xmin xmax ymin ymax" "-2 2 -2 2"
+
+    set transform [makeTransform "$xmin $ymin $x1 $y2" "$xmin $ymax $x1 $y1 " "$xmax $ymin $x2 $y2"]
+    oset $win transform $transform
+    oset $win transform0 $transform
+    
+    getXtransYtrans $transform rtosx$win rtosy$win
+    getXtransYtrans [inverseTransform $transform] storx$win story$win 
+
+}
+
+# 
+
+proc plot3d { args } {
+    global  plot3dOptions
+    set win [assoc -windowname $args]
+    if { "$win" == "" } {
+	set win [getOptionDefault windowname $plot3dOptions] }
+    clearLocal $win
+    apply mkPlot3d  $win $args
+    replot3d $win
+}
+
+proc replot3d { win } {
+    global   printOption
+    makeLocal $win nsteps zfun data
+    oset $win maintitle    "concat \"Plot of z = [oget $win zfun]\""
+    if { [llength $nsteps] == 1 }    {
+	oset $win nsteps \
+		[set nsteps  [list [lindex $nsteps 0] [lindex $nsteps 0]]]
+    }
+    if { "$zfun" != "" } {
+	    proc _xf {  x  y } "return \[expr { [sparse $zfun] } \]"
+	addOnePlot3d $win [calculatePlot3data $win _xf  [lindex $nsteps 0] [lindex $nsteps 1]]
+	# calculatePlot3d $win _xf [lindex $nsteps 0] [lindex $nsteps 1]
+    }
+    if { "$data" != "" } {
+	addOnePlot3d $win $data
+    }
+    setUpTransforms3d $win
+    oset $win colorfun plot3dcolorFun
+#    addAxes $win
+    oset $win cmap c1
+    setupPlot3dColors $win
+    addBbox $win
+    # grab the bbox just as itself
+    global ws_openMath
+    linkLocal $win lmesh
+    if { [llength $lmesh] >   100 * $ws_openMath(speed)  } {
+	# if we judge that rotation would be too slow, we make a secondary list
+	# of meshes (random) including the bbox, and display those. 
+	linkLocal $win  points lmeshBbox pointsBbox
+	set n [llength $lmesh]
+	set lmeshBbox [lrange $lmesh [expr {$n -13}] end]
+	set i 0 ;
+	while { [incr i ] < ( 35*$ws_openMath(speed)) } {
+	    set j [expr {round(floor(rand()*($n-13))) }]
+	    if { ![info exists temm($j)] } {
+		lappend lmeshBbox [lindex $lmesh $j ]
+		set temm(j) 1
+	    }
+	}
+	resetPtsForLmesh $win
+    }
+    oset $win lastAnglesPlotted ""
+    setView $win ignore
+} 
+
+proc setView { win ignore } {
+    global timer
+    foreach v [after info] {
+	 #puts "$v=<[after info $v]>"
+	if { "[lindex [after info $v] 0]" == "setView1" } {
+	    after cancel $v
+	}
+    }
+    after 2 setView1 $win
+}
+
+proc setView1 { win  } {
+    linkLocal $win lastAnglesPlotted points
+    set new [list [oget  $win az] [oget  $win el] ]
+    if { "$new" != "$lastAnglesPlotted" } {
+   	makeLocal $win c
+	calculateRotated $win
+	getOrderedMeshIndices $win
+	drawMeshes $win $c
+	oset $win lastAnglesPlotted $new
+    }
+}
+
+proc setQuick { win on } {
+  linkLocal $win  lmesh  points savedData cmap 	lmeshBbox pointsBbox
+    if { $on } {
+	if { ![info exists savedData] && [info exists lmeshBbox] } {
+	    set savedData [list $lmesh $points $cmap]
+	    set lmesh $lmeshBbox
+	    set points $pointsBbox
+	    set cmap c2
+	}
+    } else {
+	if { [info exists savedData] } {
+	    desetq "lmesh points cmap" $savedData
+	    unset savedData
+	    oset $win lastAnglesPlotted ""
+}   }   }
+
+
+# reduce the set of pointsBbox to include only those needed by lmeshBbox
+proc resetPtsForLmesh { win } {
+    upvar 1 lmeshBbox lmeshBbox
+    upvar 1 pointsBbox pointsBbox
+    upvar 1 points points
+    upvar #0 [oarray $win] wvar
+    set k 0
+    foreach v $lmeshBbox {
+	if { [llength $v] == 1 } {
+	    lappend nmesh $v
+	} else {
+	    set s ""
+	    foreach w $v {
+		if { [info exists tem($w)] } {
+		    lappend s $tem($w)
+		} else {
+		    set tem($w) $k
+		    lappend s $k
+		    lappend pointsBbox \
+			    [lindex $points $w] \
+			    [lindex $points [expr {$w +1}]] \
+			    [lindex $points [expr {$w +2}]]
+		    catch {set wvar(c2,$k) $wvar(c1,$w)}
+		    incr k 3
+		
+		}
+			
+	    }
+	    lappend nmesh $s
+	    if { [info exists wvar(c1,$v)] } {
+		set wvar(c2,$s) $wvar(c1,$v)
+	    }
+	}
+    }
+    set lmeshBbox  $nmesh
+}
+
+proc drawMeshes {win canv} {
+    # $canv delete poly
+    # only delete afterwards, to avoid relinquishing the colors
+    $canv addtag oldpoly withtag poly 
+    $canv delete axis
+    foreach k [oget $win meshes] {
+	#puts "drawOneMesh $win $canv $k"
+	#puts "drawOneMesh $win $canv $k"
+	drawOneMesh $win $canv $k
+    }
+    $canv delete oldpoly
+}
+
+
+#
+ #-----------------------------------------------------------------
+ # plot3dMeshes  --  given K an index in plot3dPoints(points) 
+ # if this is the index of a lower grid corner, return the other points.
+ # k takes values 0,3,6,9,... the values returned all have a 3 factor,  
+ # and so are true lindex indices into the list of points.
+ # returns {} if this is not a mesh point.
+ #  Results:
+ #
+ #  Side Effects: none... NOTE we should maybe cash this in an array. 
+ #
+ #----------------------------------------------------------------
+#
+
+proc drawOneMesh { win  canv k } {
+    #k=i*(ny+1)+j
+    # k,k+1,k+1+nyp,k+nyp
+    set pts [oget $win rotated]
+    set coords ""
+    set n 0
+ #   set meshes plot3dMeshes$win
+    set rtosx rtosx$win ; set rtosy rtosy$win
+    makeLocal $win lmesh cmap
+
+    set mesh [lindex $lmesh $k]
+    foreach kk $mesh {
+	incr n
+	append coords "[$rtosx [lindex $pts $kk]] [$rtosy [lindex $pts [expr {$kk + 1}]]] "
+    }
+    if { $n <= 2 } {
+	#puts "drawing $k,n=$n $coords, points $mesh "
+	#desetq "a b" $mesh
+	#puts "<[lrange $points $a [expr {$a +2}]]> <[lrange $points $b [expr {$b +2}]]"
+	if { $n == 2 } {
+	    set color gray70
+	    catch { set color [oget $win $cmap,$mesh]}
+	    eval $canv create line $coords -tags "{axis mesh.$k}" \
+		    -fill $color -width 5 
+	} else {
+	   # puts "doing special $mesh, $coords"
+	    catch { set tem [oget $win special([lindex $mesh 0])]
+	    eval [concat $tem $coords]
+	}
+	}
+    } else {
+	if { [catch { set color [oget $win $cmap,[lindex $mesh 0]] }] } { set color white}
+	
+	eval $canv create polygon $coords -tags "{poly mesh.$k}" \
+		-fill $color \
+		-outline black
+    }
+}
+
+proc doHelp3d { win } {
+ global Parser
+ doHelp $win [join [list \
+{ 
+
+William Schelter's plotter for three dimensional graphics.
+
+To QUIT this HELP click here.
+
+By clicking on Zoom, the mouse now allows you \
+to zoom in on a region of the plot.  Each click \
+near a point magnifies the plot, keeping the \
+center at the point you clicked.  Depressing \
+the SHIFT key while clicking zooms in the \
+opposite direction.
+
+Clicking on Rotate, makes the left mouse button  \
+cause rotation of the image.   The current position \
+can be determined by azimuth and elevation angles \
+which are given under the Config menu.   They may also \
+be specified on the command line.
+
+To change the equations enter in the entry \
+windows, and click on replot.
+
+You may print to a postscript printer, or save the plot \
+as a postscript file, by clicking on save.   To change \
+between printing and saving see the Print Options under Config.
+	
+Clicking with the right mouse button and dragging may be used \
+instead of the scroll bars to slide the plot \
+around.
+
+
+} $Parser(help)]]
+}
+
+proc     makeFrame3d { win } {
+  global plot3dPoints
+   set w [makeFrame $win 3d]
+    set top $w
+    catch { set top [winfo parent $w]}
+    catch {
+
+    wm title $top "Schelter's 3d Plot Window"
+    wm iconname $top "DF plot"
+ #   wm geometry $top 750x700-0+20
+   }
+  
+    pack $w
+
+}
+    
+proc mkPlot3d { win  args } {
+    global plot3dOptions  printOption [oarray $win] axisGray
+
+    getOptions $plot3dOptions $args -usearray [oarray $win]
+    #puts "$win width=[oget $win width],args=$args"
+    setPrintOptions $args
+    set printOption(maintitle) ""
+    set wb $win.buttons
+    setupCanvas $win
+   # catch { destroy $win }
+    makeFrame3d $win
+   oset $win noaxisticks 1
+   
+   makeLocal $win buttonFont c
+    bind $c <Motion> "showPosition3d $win %x %y"
+    button $wb.rotate -text "Rotate" -command "setForRotate $win" -font $buttonFont
+   setBalloonhelp $win $wb.rotate {Dragging the mouse with the left button depressed will cause the object to rotate.  The rotation keeps the z axis displayed in an upright position (ie parallel to the sides of the screen), but changes the viewpoint.   Moving right and left changes the azimuth (rotation about the z axis), and up and down changes the elevation (inclination of z axis).   The red,blue and green sides of the bounding box are parallel to the X, Y and Z axes, and are on the smaller side.} 
+
+   $win.position config -width 15
+    pack $wb.rotate -expand 1 -fill x
+   setForRotate $win
+
+    
+}   
+
+proc doConfig3d { win } {
+
+    
+    desetq "wb1 wb2" [doConfig $win]
+
+    makeLocal $win buttonFont
+
+    mkentry $wb1.zfun [oloc $win zfun]  "z=f(x,y)" $buttonFont 
+    mkentry $wb1.nsteps [oloc $win nsteps]  "Number of mesh grids"  $buttonFont 
+    # button .jim.buttons.rot "rotate" -command "bindForRotation"
+    # pack .jim.buttons.rot
+    pack $wb1.zfun  $wb1.nsteps
+    pack	    $wb1.zfun  $wb1.nsteps 
+   foreach w {xradius yradius xcenter ycenter zcenter zradius  } {
+	mkentry $wb1.$w [oloc $win $w] $w $buttonFont
+	pack $wb1.$w 
+    }
+
+    scale $wb1.rotxscale -label "azimuth"  \
+	    -orient horizontal -length 150 -from -180 -to 180 -resolution 1 \
+	    -command "setView $win" -variable [oloc $win az] -tickinterval 120 -font $buttonFont 
+
+    scale $wb1.rotyscale -label "elevation"  \
+	    -orient horizontal -length 150 -from -180 -to 180 -resolution 1 \
+	    -command "setView $win" -variable [oloc $win el] -tickinterval 120 -font $buttonFont 
+
+
+#    scale $wb1.rotzscale -label "thetaz"  \
+#	    -orient horizontal -length 150 -from -180 -to 180 \
+#	    -command "setView $win" -variable [oloc $win thetaz] -tickinterval 120 -font $buttonFont 
+
+    pack   $wb1.rotxscale   $wb1.rotyscale   
+
+}
+
+
+proc showPosition3d { win x y } {
+   # global position c
+    makeLocal $win c
+    set x [$c canvasx $x]
+    set y [$c canvasy $y]
+    set it [ $c find closest $x $y]
+    set tags [$c gettags $it]
+    if { [regexp {mesh[.]([0-9]+)} $tags junk k] } {
+	set i 0
+	set min 1000000
+	set at 0
+	# find closest.
+	foreach {x1 y1} [$c coords $it] {
+	    set d [expr {($x1 - $x)*($x1 - $x)+($y1 - $y)*($y1 - $y)}]
+	    if { $d < $min} { set at $i ; set min $d }
+	    incr i
+	}
+	set mesh [lindex [oget $win lmesh] $k]
+	set ll [lindex $mesh $at]
+	set pt [lrange [oget $win points] $ll [expr {$ll + 2}]]
+	# puts pt=$pt
+	catch { oset $win position [eval [concat "format {(%.2f %.2f %.2f)}" $pt]] }
+    }
+#    oset $win position [format {(%.1f %.1f)} $x $y]
+#    oset $win position \
+#      "[format {(%.2f,%.2f)}  [storx$win [$c canvasx $x]] [story$win [$c canvasy $y]]]"
+}
+
+
+
+#
+ #-----------------------------------------------------------------
+ #
+ # rotateRelative --  do a rotation indicated by a movement
+ # of dx,dy on the screen.  
+ #
+ #  Results:
+ #
+ #  Side Effects: 
+ #
+ #----------------------------------------------------------------
+#
+
+proc rotateRelative { win x1 x2 y1 y2 } {
+    makeLocal $win c az el rotatefun
+    set x1 [$c canvasx $x1]
+    set x2 [$c canvasx $x2]
+    set y1 [$c canvasy $y1]
+    set y2 [$c canvasy $y2]
+    set xx [expr {$x2-$x1}]
+    set yy [expr {($y2-$y1)}]
+    set res [$rotatefun 0 0 1]
+    set res1 [$rotatefun 0 0 0]
+    set fac [expr {([lindex $res 1] > [lindex $res1 1] ? -1 : 1) }] ;
+   # puts "fac=$fac,[lindex $res 1],[lindex $res1 1]"
+    oset $win az [reduceMode360 [expr   {round($az + $fac *  $xx /2.0) }]]
+    oset $win el [reduceMode360 [expr   {round($el -  $yy /2.0) }]]
+    setView $win ignore
+
+
+}
+
+proc reduceMode360 { n } {
+  return [  expr fmod(($n+180+5*360),360)-180]
+
+}
+
+proc setForRotate { win} {
+    makeLocal $win c
+    $c delete printrectangle
+    bind $c  <Button-1> "setQuick $win 1 ; doRotateScreen $win %x %y "
+    bind $c  <ButtonRelease-1> "setQuick $win 0 ; setView $win ignore"
+}
+proc doRotateScreen { win x y } {
+    makeLocal $win c
+    oset $win lastx $x
+    oset $win lasty $y
+    bind $c <B1-Motion> "doRotateScreenMotion $win %x %y"
+    
+
+}
+
+proc doRotateScreenMotion {win x y } {
+    makeLocal $win lastx lasty
+    set dx [expr {$x - $lastx}]
+    set dy [expr {$y - $lasty}]
+    if { [vectorlength $dx $dy] < 4 } { return }
+    rotateRelative $win $lastx $x $lasty $y
+    oset $win lastx $x
+    oset $win lasty $y
+    
+}
+
+    
+    
+    
+## endsource plot3d.tcl
+## source nplot3d.tcl
+
+###### nplot3d.tcl ######
+# source plotting.tcl ; source nplot3d.tcl ; catch { destroy .plot3d} ;  plot3d -zfun "" -data $sample -xradius 10 -yradius 10
+# newidea:
+# { plot3d
+#  { gridequal {minx maxx} {miny maxy}
+#   {{z00 z01 z02 .. z0n } { z10 z11 z12 .. z1n} {..  } ...}
+#  { grid {x0 x1  xm} {y0 y1 yn } miny maxy}
+#   {{z00 z01 z02 .. z0n } { z10 z11 z12 .. z1n} {..  } ...}
+#  { xyzgrid {{x00 y00 z00 x01 y01 z01 .. x0  }{x0 x1  xm} {y0 y1 yn } miny maxy}
+#   {{z00 z01 z02 .. z0n } { z10 z11 z12 .. z1n} {..  } ...}
+# tclMesh(2*[0,0,0,0,0;1,1,1,1,1]-1,2*[0,1,1,0,0;0,1,1,0,0]-1,2*[0,0,1,1,0;0,0,1,1,0]-1)
+  
+#     { gridequal { 
+
+# z00 z01 .. all belong to x=minx and y = miny,.... up y=maxy in n+1 steps
+#{ grid {minx maxx} {miny maxy}
+#  {{z00 z01 z02 .. z0n } { z10 z11 z12 .. z1n} {..  } ...}
+# }
+# where a mesh(1) {z00 z01 z11 z10} above 
+
+
+
+# { mesh {{{x00 y00 z00 } { x01 y01 z01} { x02 y02 z02}  ..}{{x10 y10 z10} {x11 y11 z11} ......} ..}}
+# mesh(1) = P00 P01 P11 P10
+
+set sample { variable_grid { 0 1 2 } { 3 4 5} { {21 111 2} {3 4 5 } {6 7 8 }}}
+set sample { variable_grid { 0 1 2 } { 3 4 5} { {0 1 2} {3 4 5 } {6 7 8 }}}
+set sample { matrix_mesh {{0 1} { 2 3 } {4 5 }}  {{0 1} { 2 3 } {4 5 }}  {{0 1} { 2 3 } {4 5 }} }
+set sample { matrix_mesh {{0 1 2} {0 1 2 } {0 1 2 }} {{3 4 5} {3 4 5} {3 4 5}} { {0 1 2} {3 4 5 } {6 7 8 }}}
+set sample1 { variable_grid  { 1 2 3 4 5 6 7 8 9 10 }
+ { 1 2 3 }
+ {  { 0 0 0 0 0 0 0 0 0 0 }
+ { 0 0.68404 1.28558 1.73205 1.96962 1.96962 1.73205 1.28558 0.68404 2.44921e-16 }
+ { 0 1.36808 2.57115 3.4641 3.93923 3.93923 3.4641 2.57115 1.36808 4.89843e-16 }
+ }  }
+
+set sample { matrix_mesh  {  { 0 0 0 0 0 }
+ { 1 1 1 1 1 }
+ }  {  { 0 1 1 0 0 }
+ { 0 1 1 0 0 }
+ }  {  { 0 0 1 1 0 }
+ { 0 0 1 1 0 }
+ }  } 
+
+    
+proc  fixupZ { } {
+    uplevel 1 {
+	if { [catch { expr $z + 0 } ] } {
+	    set z nam
+	}  elseif { $dotruncate  &&  ($z > $zzmax || $z < $zzmin) } {
+	    set z nam
+	    
+	} else {
+	    if { $flatten } {
+		if { $z > $zzmax } { set z $zzmax } elseif {
+		    $z < $zzmin } { set z $zzmin }}
+		    
+		    if { $z < $zmin }  { set zmin $z } elseif {
+			$z > $zmax } { set zmax $z }
+		    }
+		}
+}
+
+
+proc vectorLength { v } {
+    expr { sqrt(1.0 * [lindex $v 0]*[lindex $v 0] + [lindex $v 1]*[lindex $v 1] + [lindex $v 2]*[lindex $v 2]) }
+}
+
+proc normalizeToLengthOne { v } {
+    set norm [expr { sqrt(1.0 * [lindex $v 0]*[lindex $v 0] + [lindex $v 1]*[lindex $v 1] + [lindex $v 2]*[lindex $v 2]) }]
+    if { $norm != 0.0 } {
+	return [list [expr { [lindex $v 0] / $norm  } ] \
+		[expr { [lindex $v 1] / $norm  } ] \
+		[expr { [lindex $v 2] / $norm  } ] ]
+	 
+    } else { return "1.0 0.0 0.0 " }
+}
+    
+    
+
+proc vectorCross { x1 x2 }  {
+     list \
+      [expr { [lindex $x1 1]*[lindex $x2 2]- [lindex $x2 1]*[lindex $x1 2]}] \
+      [expr { [lindex $x1 2]*[lindex $x2 0]- [lindex $x2 2]*[lindex $x1 0] } ] \
+      [expr { [lindex $x1 0]*[lindex $x2 1]- [lindex $x2 0]*[lindex $x1 1] }]
+}
+    
+proc linspace { a b n } {
+    if { $n < 2 } { error "from $a to $b requires at least 2 points" }
+    set del [expr {($b - $a)*1.0/($n -1)  }]
+    for { set i 0 } { $i < $n } { incr i } {
+	lappend ans [expr {$a + $del * $i}]
+    }
+    return $ans
+}
+
+
+proc addOnePlot3d { win data } {
+    upvar #0 plot3dMeshes$win meshes 
+    #puts " adding meshes = plot3dMeshes$win"
+    #puts "data=$data"
+    linkLocal $win points zmax zmin zcenter zradius rotationcenter xradius yradius xmin xmax ymin ymax lmesh
+    makeLocal $win flatten 
+    catch { unset  meshes }
+    set points ""
+
+
+    set dotruncate [expr ![catch {expr {$zradius + 1} }]]
+    set k [llength $points]
+    set type [lindex $data 0]
+    # in general the data should be a list of plots..
+    if { [lsearch {grid mesh variable_grid matrix_mesh }  $type ]>=0 } {
+	set alldata [list $data]
+    } else {set alldata $data}
+    foreach data $alldata {	
+	set type [lindex $data 0]
+    if { "$type" == "grid" } {
+	desetq "xmin xmax" [lindex $data 1]
+	desetq "ymin ymax" [lindex $data 2]
+	set pts [lindex $data 3]
+	
+	set ncols [llength $pts]
+	set nrows  [llength [lindex $pts 0]]
+	set data [list variable_grid [linspace $xmin $xmax $ncols] \
+		[linspace $ymin $ymax $nrows] \
+		$pts ]
+    }
+    if { "$type" == "variable_grid" } {
+	desetq "xrow yrow zmat" [lrange $data 1 end]
+	# puts "xrow=$xrow,yrow=$yrow,zmat=$zmat"
+	set nx [expr {[llength $xrow] -1}]
+	set ny [expr {[llength $yrow] -1}]
+	#puts "nx=$nx,ny=$ny"
+	set xmin [lindex $xrow 0]
+	set xmax [lindex $xrow $nx]
+	set ymin [lindex $yrow 0]
+	set ymax [lindex $yrow $ny]
+	desetq "zmin zmax" [matrixMinMax [list $zmat]]
+	
+	if { $dotruncate } {
+	    if { $flatten } { set dotruncate 0 }
+
+	    set zzmax [expr {$zcenter + $zradius}]
+	    set zzmin [expr {$zcenter - $zradius}]
+	    #puts "zzmax=$zzmax,$zzmin"
+	} else { set flatten 0 }
+
+
+
+	for {set j 0} { $j <= $ny } { incr j} {
+	    set y [lindex $yrow $j]
+	    set row [lindex $zmat $j]    
+	for {set i 0} { $i <= $nx } { incr i} {
+	    set x [lindex $xrow $i]
+	    set z [lindex $row $i]
+	    #puts "x=$x,y=$y,z=$z, at ($i,$j)"
+	    fixupZ 
+	    if { $j != $ny && $i != $nx } {
+		lappend lmesh [list $k [expr { $k+3 }] \
+			[expr { $k+3+($nx+1)*3 }] \
+		      [expr { $k+($nx+1)*3 }]]
+	    }
+	      incr k 3
+	  lappend points $x $y $z
+	  }
+	}
+    } elseif { "$type" == "matrix_mesh" } {
+	
+	desetq "xmat ymat zmat" [lrange $data 1 end]
+	foreach v {x y z} {
+	    
+	    
+	    desetq "${v}min ${v}max" [matrixMinMax [list [set ${v}mat]]]
+	    
+	}
+	#puts "zrange=$zmin,$zmax"
+	set nj [expr {[llength [lindex $xmat 0]] -1 }]
+	set ni [expr {[llength $xmat ] -1 }]
+	set i -1
+	set k [llength $points]
+	foreach rowx $xmat rowy $ymat rowz $zmat {
+	    set j -1
+	    incr i
+	    if { [llength $rowx] != [llength $rowy] } {
+		error "mismatch rowx:$rowx,rowy:$rowy"
+	    }
+	    if { [llength $rowx] != [llength $rowz] } {
+		error "mismatch rowx:$rowx,rowz:$rowz"
+	    }
+	    foreach x $rowx y $rowy z $rowz {
+		incr j
+		if { $j != $nj && $i != $ni } {
+		#puts "tes=($i,$j) $x, $y, $z"
+		    lappend lmesh [ list \
+			    $k [expr { $k+3 } ] [expr { $k + 3  + ($nj+1)*3}] \
+			    [expr { $k+($nj+1)*3 }] ]
+		}
+		incr k 3
+		lappend points $x $y $z
+	    }
+	}
+    } elseif { 0 && "$type" == "mesh" } {
+  # walk thru compute the xmin, xmax, ymin , ymax...
+  # and then go thru setting up the mesh array..
+  # and maybe setting up the color map for these meshes..
+  #
+    # { mesh {{{x00 y00 z00 } { x01 y01 z01} { x02 y02 z02}  ..}{{x10 y10 z10} {x11 y11 z11} ......} ..}}
+# mesh(1) = P00 P01 P11 P10
+    set mdata [lindex $data 1]
+    set nx [llength $mdata]
+    set ny [llength [lindex $mdata 0]]
+    
+    for {set i 0} { $i <= $nx } { incr i} {
+	set pts [lindex $mdata $i]
+	set j 0
+	foreach { x y z} $pts {
+	    fixupZ $z
+	    if { $j != $ny && $i != $nx } {
+		lappend lmesh [list 
+			$k [expr { $k+3 }] [expr { $k+3+($ny+1)*3 }] \
+			[expr { $k+($ny+1)*3 }] ]
+		}
+	    }
+	    incr k 3
+	    lappend points $x $y $z
+	    incr j
+	}
+    }
+   }
+    foreach v { x y z } {
+	set a [set ${v}min]
+	set b  [set ${v}max]
+	if { $a == $b } {
+	    set ${v}min [expr {$a -1}]
+	    set ${v}max [expr {$a +1}]
+	}
+	set ${v}radius [expr {($b - $a)/2.0}]
+	set ${v}center [expr {($b + $a)/2.0}]
+    }
+    if { "$rotationcenter" == "" } {
+	set rotationcenter "[expr {.5*($xmax + $xmin)}] [expr {.5*($ymax + $ymin)}]   [expr {.5*($zmax + $zmin)}] "
+    }
+    
+    #puts "meshes data=[array get meshes]"
+    #global plot3dMeshes.plot3d
+    #puts "array names plot3dMeshes.plot3d = [array names plot3dMeshes.plot3d]"
+}
+
+proc vectorDiff { x1 x2 } {
+    list [expr { [lindex $x1 0] - [lindex $x2 0] }] \
+	    [expr { [lindex $x1 1] - [lindex $x2 1] }] \
+	    [expr { [lindex $x1 2] - [lindex $x2 2] }]
+}
+
+
+proc oneCircle { old2 old1 pt radius nsides { angle 0 } } {
+    set dt  [expr {  3.14159265358979323*2.0/($nsides-1.0) + $angle }]
+    for  { set i 0 } { $i < $nsides } { incr i } {
+	set t [expr {$dt*$i }]
+	lappend ans [expr { $radius*([lindex $old2 0]*cos($t) + [lindex $old1 0] * sin($t)) + [lindex $pt 0] } ] \
+		[expr { $radius*([lindex $old2 1]*cos($t) + [lindex $old1 1] * sin($t)) + [lindex $pt 1] } ] \
+		[expr { $radius*([lindex $old2 2]*cos($t) + [lindex $old1 2] * sin($t)) + [lindex $pt 2] } ]
+    }
+    return $ans
+}
+
+proc curve3d { xfun yfun zfun trange } {
+    foreach u { x y z} {
+	set res [parseConvert [set ${u}fun] -variables t]
+	proc _${u}fun { t } [list expr [lindex [lindex $res 0] 0]]
+}   }
+
+proc tubeFromCurveData { pts nsides radius } {
+    set n [llength $pts] ;
+    set closed [ expr { [vectorLength [vectorDiff [lindex $pts 0] [lindex $pts end]]] < .02} ]
+    if { $closed } {
+	set f1 [expr {$n -2}]
+	set f2 1
+    } else { set f1 0
+	 set f2 1
+    }
+    set delta [vectorDiff [lindex $pts $f2] [lindex $pts $f1]]
+    if { [lindex $delta 0] == 0 && [lindex $delta 1] == 0 && [lindex $delta 2] == 0 } { set delta "0 0 1.0" }
+    set old ".6543654 0.0765456443 0.2965433"
+    set old1 [normalizeToLengthOne [vectorCross $delta $old]]
+    set n1 $old1
+    set n2 [normalizeToLengthOne [vectorCross $delta $old1]]
+    set first1 $n1 ; set first2 $n2
+    
+    lappend ans [oneCircle $n2   old1 [lindex $pts 0]]
+    for { set j 1 } { $j < $n -1 } { incr j } {
+	set delta [vectorDiff [lindex $pts $j] [lindex $pts [expr {$j+1}]]]
+	if { [lindex $delta 0] == 0 && [lindex $delta 1] == 0 && [lindex $delta 2] == 0 } { set delta $old
+    }
+    set old $delta
+    set old1 [normalizeToLengthOne [vectorCross $delta $n1]]
+    set old2 [normalizeToLengthOne [vectorCross $delta $n2]]
+    set n2 $old1
+    set n1 $old2
+    lappend ans [oneCircle $n2 $n1 [lindex $pts $j] $radius $nsides]
+}
+    if { $closed } {
+	set f2 1 ; set f1 [expr {$n -2}] ; set f3 0
+    } else {
+	set f1 [expr {$n -2}] ; set f2 [expr {$n-1}] ; set f3 $f2
+    }
+
+    set delta [vectorDiff [lindex $pts $f2] [lindex $pts $f1]]
+    if { [lindex $delta 0] == 0 && [lindex $delta 1] == 0 && \
+	    [lindex $delta 2] == 0 } { set delta $old }
+    set old1 [normalizeToLengthOne [vectorCross delta $n1]]
+    set old2 [normalizeToLengthOne [vectorCross $n2 $delta]]
+    set n2 $old1 ; set n1 $old2
+    if { $closed } {
+	set angle [vangle $first1 $n1]
+	set n1 $first1 ; st n2 $first2;
+    }
+    lappend ans [oneCircle $n2 $n1 [lindex $pts $f3] $radius $nsides $angle]
+   return $ans
+}
+
+
+#
+ #-----------------------------------------------------------------
+ #
+ # vangle --  angle between two unit vectors
+ #
+ #  Results: an angle
+ #
+ #  Side Effects: none.
+ #
+ #----------------------------------------------------------------
+#
+proc vangle { x1 x2 } {
+    set dot [expr { [lindex $x1 0]*[lindex $x2 0] +\
+	     [lindex $x1 1]*[lindex $x2 1] +\
+	     [lindex $x1 2]*[lindex $x2 2]} ]
+    if { $dot >= 1 } { return 0.0 }
+    if { $dot <= -1.0 } { return 3.141592653589 }
+    return [expr { acos($dot) } ]
+}
+
+## endsource nplot3d.tcl
+
+# from shell 
+# wish8.0 plotting.tcl -eval {plot2d -xfun  x^2+3}
+# or in html
+# <embed src=plotting.tcl eval="plot2d -xfun  x^2+3" >
+#
+
+
+
+omPlotAny [exec cat [lindex $argv 0]]