home *** CD-ROM | disk | FTP | other *** search
/ MacFormat 1995 January / macformat-020.iso / Shareware City / Education / RLaB / rlib / plplot.r < prev    next >
Encoding:
Text File  |  1994-07-31  |  41.4 KB  |  1,854 lines  |  [TEXT/ttxt]
  1. #
  2. # New plot.r for use with PLPLOT library.
  3. # The help files for these functions are in
  4. # misc/plhelp
  5. #
  6.  
  7. # plplot.r
  8.  
  9. #  This file is a part of RLaB ("Our"-LaB)
  10. #  Copyright (C) 1994  Ian R. Searle
  11.  
  12. #  This program is free software; you can redistribute it and/or modify
  13. #  it under the terms of the GNU General Public License as published by
  14. #  the Free Software Foundation; either version 2 of the License, or
  15. #  (at your option) any later version.
  16.  
  17. #  This program is distributed in the hope that it will be useful,
  18. #  but WITHOUT ANY WARRANTY; without even the implied warranty of
  19. #  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  20. #  GNU General Public License for more details.
  21.  
  22. #  You should have received a copy of the GNU General Public License
  23. #  along with this program; if not, write to the Free Software
  24. #  Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
  25.  
  26. #  See the file ../COPYING
  27.  
  28.  
  29. static (WIN)            # The static plot window structure
  30. static (P)            # The active/current plot window
  31.  
  32. static (create_plot_object)
  33. static (check_plot_object)
  34. static (xy_scales)
  35. static (xyz_scales)
  36. static (list_scales)
  37. static (hist_scales)
  38. static (plot_matrix)
  39. static (plot_list)
  40. static (check_3d_list)
  41. static (find_char)
  42. static (get_style)
  43. static (make_legend)
  44. static (plhold_first)
  45.  
  46. #
  47. # Defaults
  48. #
  49.  
  50. static (grid_x_default, grid_y_default)
  51. static (grid_3x_default, grid_3y_default, grid_3z_default)
  52.  
  53. grid_x_default = "bcnst";
  54. grid_y_default = "bcnstv";
  55. grid_3x_default = "bnstu";
  56. grid_3y_default = "bnstu";
  57. grid_3z_default = "bcdmnstuv";
  58.  
  59. #
  60. # Create the default plot-object.
  61. # Initialize to all the default values
  62. #
  63.  
  64. WIN = <<>>;    # Create the plot-object list
  65.  
  66. create_plot_object = function ( N, nx, ny )
  67. {
  68.   local (i, j, pobj)
  69.   
  70.   if (!exist (N)) { N = 0; }
  71.   
  72.   pobj = <<>>;
  73.   pobj.subplot = 0;        # The current subplot no.
  74.   pobj.nplot = nx*ny;        # Total no. of plots on window
  75.   
  76.   pobj.fontld = 0;        # Loaded extended fonts?
  77.   
  78.   for (i in 1:(nx*ny))
  79.   {
  80.     pobj.style.[i] = "line";    # The type/style of plot to draw
  81.     pobj.pstyle[i] = 1;            # The point-style
  82.     pobj.nbin[i] = inf();    # The number of bins for a histogram
  83.     pobj.width[i] = 1;        # The pen width for current plot
  84.     pobj.font[i] = 1;        # The current font
  85.     pobj.xlabel[i] = "";
  86.     pobj.ylabel[i] = "";
  87.     pobj.zlabel[i] = "";
  88.     pobj.title[i] = "";
  89.     pobj.orientation[i] = "portrait";
  90.     pobj.desc.[i] = "default";        # The legend description
  91.     pobj.gridx[i] =  grid_x_default;    # Plot axes style, 2D-X
  92.     pobj.gridy[i] =  grid_y_default;    # Plot axes style, 2D-Y
  93.     pobj.grid3x[i] = grid_3x_default;    # Plot axes style, 3D-X
  94.     pobj.grid3y[i] = grid_3y_default;    # Plot axes style, 3D-Y
  95.     pobj.grid3z[i] = grid_3z_default;    # Plot axes style, 3D-Z
  96.     pobj.aspect[i] = 0;        # Plot aspect style
  97.     pobj.alt[i] = 60;
  98.     pobj.az[i] = 45;
  99.     
  100.     pobj.xmin[i] = inf();
  101.     pobj.xmax[i] = inf();
  102.     pobj.ymin[i] = inf();
  103.     pobj.ymax[i] = inf();
  104.     pobj.zmin[i] = inf();
  105.     pobj.zmax[i] = inf();
  106.     
  107.     pobj.page.xp = 0;
  108.     pobj.page.yp = 0;
  109.     pobj.xleng = 400;
  110.     pobj.yleng = 300;
  111.     pobj.xoff = 200;
  112.     pobj.yoff = 200;
  113.  
  114.     for (j in 1:14) { pobj.color[i;j] = j; }       
  115.     for (j in 1:8)  { pobj.lstyle[i;j] = j; }
  116.   }
  117.   
  118.   #
  119.   # Save the newly generated plot-object
  120.   # in a list of plot-objects.
  121.   
  122.   WIN.[N] = pobj;
  123. };
  124.  
  125. #
  126. # Check to make sure a plot-object exists. If one
  127. # does not exist, create it.
  128. #
  129.  
  130. check_plot_object = function ()
  131. {
  132.   if (length (WIN) == 0)
  133.   {
  134.     pstart();
  135.     return 0;
  136.   }
  137.   return 1;
  138. };
  139.  
  140. #
  141. # Set the current plot window
  142. # Default value = 0
  143. #
  144.  
  145. pwin = function ( N )
  146. {
  147.   check_plot_object ();
  148.   if (!exist (N)) { N = 0; }
  149.   
  150.   # Check to make sure N is valid
  151.   
  152.   for (i in members (WIN))
  153.   {
  154.     if (N == strtod (i))
  155.     {
  156.       _plsstrm (N);
  157.       return P = N;
  158.     }
  159.   }
  160.   printf ("pwin: invalid argument, N = %i\n", N);
  161.   printf ("      valid values are:\n");
  162.   WIN?
  163. };
  164.  
  165. showpwin = function ()
  166. {
  167.   check_plot_object ();
  168.   
  169.   printf ("Curren plot-window is: %i\n", P);
  170.   printf ("Available plot windows are:\n");
  171.   WIN?
  172. };
  173.  
  174. #
  175. # Set/start/select the plot device
  176. #
  177.  
  178. pstart = function ( nx, ny, dev )
  179. {
  180.   if (!exist (nx)) { nx = 1; }
  181.   if (!exist (ny)) { ny = 1; }
  182.   if (!exist (dev)) { dev = "?"; }
  183.   
  184.   # Create the plot-object
  185.   # First, figure out the index
  186.   if (!exist (P))
  187.   {
  188.     P = 0;
  189.   else
  190.     P = P + 1;
  191.   }
  192.   
  193.   create_plot_object (P, nx, ny);
  194.   _plsstrm (P);
  195.   
  196.   #_plscolbg (255,255,255); # white
  197.   # Default window size for X
  198.   _plspage (0, 0, 400, 300, 200, 200);
  199.   
  200.   # Start up the plot-window
  201.   _plstart (dev, nx, ny);
  202.   
  203.   _plwid (8);
  204.   
  205.   # Turn between plot pause off
  206.   _plspause (0);
  207.   _pltext ();
  208.   
  209.   return P;
  210. };
  211.  
  212. #
  213. # Close a plot device. We must destroy the current plot-object
  214. # And switch the output stream back to the default.
  215. #
  216.  
  217. pclose = function ()
  218. {
  219.   local (n)
  220.   if (size (WIN) > 1)
  221.   {   
  222.     clear (WIN.[P]);
  223.     _plend1 ();
  224.     _plsstrm (strtod (members (WIN)[1]));
  225.     P = strtod (members (WIN)[1]);
  226.   else
  227.     if (exist (WIN.[P])) { clear (WIN.[P]); }
  228.     _plend1 ();
  229.   }
  230. };
  231.  
  232. #
  233. # Close ALL the plot-windows
  234. #
  235.  
  236. pend = function ()
  237. {
  238.   _plend ();
  239.   if (exist (WIN)) { clear (WIN); }
  240.   if (exist (P)) { clear (P); }
  241.   WIN = <<>>;
  242. };
  243.  
  244. ##############################################################################
  245. #
  246. # Plot the columns of a matrix (X-Y plot).
  247. #
  248. ##############################################################################
  249.  
  250. plot = function ( data )
  251. {
  252.   local (hscale, xmin, xmax, ymin, ymax, i, p)
  253.   
  254.   check_plot_object ();
  255.   
  256.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  257.   
  258.   #
  259.   # Draw the graph
  260.   # Step through the matrix plotting
  261.   # each column versus the 1st
  262.   #
  263.   
  264.   if (class (data) == "num")
  265.   {
  266.     _plgra ();
  267.     _plcol (1);
  268.     _pllsty (1);
  269.     _plfont (WIN.[P].font[p]);
  270.     _plwid (WIN.[P].width[p]);
  271.     xy_scales ( real(data), p, xmin, xmax, ymin, ymax );
  272.     
  273.     _pladv (0);
  274.     _plvsta ();
  275.     _plwind (xmin, xmax, ymin, ymax);
  276.     _plbox (WIN.[P].gridx[p], 0, 0, WIN.[P].gridy[p], 0, 0);
  277.     if (plot_matrix ( data, p, 0, xmin, xmax, ymin, ymax, ymax-ymin ) < 0) 
  278.     {
  279.       return -1;
  280.     }
  281.     
  282.     else if (class (data) == "list") {
  283.       
  284.       _plgra ();
  285.       _plcol (1);
  286.       _pllsty (1);
  287.       _plfont (WIN.[P].font[p]);
  288.       _plwid (WIN.[P].width[p]);
  289.       list_scales ( data , p, xmin, xmax, ymin, ymax );
  290.       _pladv (0);
  291.       _plvsta ();
  292.       _plwind (xmin, xmax, ymin, ymax);
  293.       _plbox (WIN.[P].gridx[p], 0, 0, WIN.[P].gridy[p], 0, 0);
  294.       if (plot_list ( data, p, xmin, xmax, ymin, ymax ) < 0) 
  295.       { 
  296.     return -1;
  297.       }
  298.  
  299.     else
  300.       error ("plot: un-acceptable argument");
  301.   } }
  302.  
  303.   _pllsty (1);
  304.   _plcol (1);
  305.   _pllab (WIN.[P].xlabel[p], WIN.[P].ylabel[p], WIN.[P].title[p]);  
  306.   _plflush ();
  307.   _pltext ();
  308.   
  309.   #
  310.   # Increment the plot no. so that next time
  311.   # we use the correct settings.
  312.   #
  313.   
  314.   WIN.[P].subplot = WIN.[P].subplot + 1;
  315.   return P;
  316. };
  317.  
  318. #
  319. # plhold:
  320. # Plot some data, and "hold" on for more.
  321. # Plot the data, setting up the plot as usual the first time.
  322. # On subsequent plots do not do any setup, just plot some
  323. # more. plhold_off must be called to finish up.
  324. #
  325.  
  326. plhold_first = 1;    # True (1) if plhold() has NOT been used.
  327.                         # Or if plhold_off() has been used.
  328.             # False (0) if plhold is in use
  329.  
  330. static (hxmin, hxmax, hymin, hymax)
  331.  
  332. plhold = function ( data )
  333. {
  334.   local (hscale, xmin, xmax, ymin, ymax, i, p)
  335.   
  336.   check_plot_object ();
  337.   
  338.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  339.   
  340.   if (plhold_first)
  341.   {
  342.     if (class (data) == "num")
  343.     {
  344.       #
  345.       # Do the setup ONCE
  346.       #
  347.       hxmin = hxmax = hymin = hymax = 0;
  348.       _plgra ();
  349.       _plcol (1);
  350.       _pllsty (1);
  351.       _plfont (WIN.[P].font[p]);
  352.       _plwid (WIN.[P].width[p]);
  353.       xy_scales ( real(data), p, hxmin, hxmax, hymin, hymax );
  354.       
  355.       _pladv (0);
  356.       _plvsta ();
  357.       _plwind (hxmin, hxmax, hymin, hymax);
  358.       _plbox (WIN.[P].gridx[p], 0, 0, WIN.[P].gridy[p], 0, 0);
  359.       _pllab (WIN.[P].xlabel[p], WIN.[P].ylabel[p], WIN.[P].title[p]);
  360.     else
  361.       error ("plot: un-acceptable argument");
  362.     }
  363.     plhold_first = 0;
  364.   }
  365.   
  366.   if (plot_matrix ( data, p, 0, hxmin, hxmax, hymin, hymax, hymax-hymin ) < 0) 
  367.   { 
  368.     return -1; 
  369.   }
  370.  
  371.   _plcol (1);
  372.   _plflush ();
  373.   _pltext ();
  374.  
  375.   return P;
  376. };
  377.  
  378. #
  379. # Clean up the plotting environment and get ready
  380. # for normal interactive usage.
  381. #
  382.  
  383. plhold_off = function ( )
  384. {
  385.   local (p)
  386.   
  387.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  388.   plhold_first = 1;  
  389.   _plcol (1);
  390.   _plflush ();
  391.   _pltext ();
  392.   WIN.[P].subplot = WIN.[P].subplot + 1;
  393.   return P;
  394. };
  395.  
  396. ##############################################################################
  397. #
  398. # Plot a 3-Dimensional graph. The data is composed in a list, with
  399. # elements `x', `y', and `z'. x and y are single-dimension arrays
  400. # (row or column matrices), and z is a two-dimensional array. The
  401. # array z, is a function of x and y: z = f(x,y). Thus, the values in
  402. # the array x can be thought of a "row-labels", and the values of y
  403. # can be thought of as "column-lables" for the 2-dimensioal array z.
  404. #
  405. # At present plot3 can plot 3 distinct lists. Each list may have
  406. # different X, Y, and Z scales.
  407. #
  408. ##############################################################################
  409.  
  410. plot3 = function ( L31, L32, L33 )
  411. {
  412.   local (p, xmin, xmax, ymin, ymax, zmin, zmax, basex, basey, height, ...
  413.          xmin2d, xmax2d, ymin2d, ymax2d, alt, az, ...
  414.          Xmin, Xmax, Ymin, Ymax, Zmin, Zmax)
  415.  
  416.   check_plot_object ();
  417.  
  418.   #
  419.   # 1st check list contents
  420.   #
  421.   
  422.   if (exist (L31)) { check_3d_list (L31); }
  423.   if (exist (L32)) { check_3d_list (L32); }
  424.   if (exist (L33)) { check_3d_list (L33); }
  425.   
  426.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  427.  
  428.   #
  429.   # Figure out the scale limits. 
  430.   # Needs improvement!
  431.   #
  432.   
  433.   xmin = xmax = ymin = ymax = zmin = zmax = 0;
  434.   if (exist (L31)) 
  435.   {
  436.     xyz_scales (L31, p, Xmin, Xmax, Ymin, Ymax, Zmin, Zmax);
  437.     if (Xmin < xmin) { xmin = Xmin; } if (Xmax > xmax) { xmax = Xmax; }
  438.     if (Ymin < ymin) { ymin = Ymin; } if (Ymax > ymax) { ymax = Ymax; }
  439.     if (Zmin < zmin) { zmin = Zmin; } if (Zmax > zmax) { zmax = Zmax; }
  440.   }
  441.   if (exist (L32)) 
  442.   {
  443.     xyz_scales (L32, p, Xmin, Xmax, Ymin, Ymax, Zmin, Zmax);
  444.     if (Xmin < xmin) { xmin = Xmin; } if (Xmax > xmax) { xmax = Xmax; }
  445.     if (Ymin < ymin) { ymin = Ymin; } if (Ymax > ymax) { ymax = Ymax; }
  446.     if (Zmin < zmin) { zmin = Zmin; } if (Zmax > zmax) { zmax = Zmax; }
  447.   }   
  448.   if (exist (L33)) 
  449.   {
  450.     xyz_scales (L33, p, Xmin, Xmax, Ymin, Ymax, Zmin, Zmax);
  451.     if (Xmin < xmin) { xmin = Xmin; } if (Xmax > xmax) { xmax = Xmax; }
  452.     if (Ymin < ymin) { ymin = Ymin; } if (Ymax > ymax) { ymax = Ymax; }
  453.     if (Zmin < zmin) { zmin = Zmin; } if (Zmax > zmax) { zmax = Zmax; }
  454.   }
  455.   
  456.   _plgra ();
  457.   _plcol (1);
  458.   _pllsty (1);
  459.   _plfont (WIN.[P].font[p]);
  460.   _plwid (WIN.[P].width[p]);
  461.   
  462.   basex = 2; basey = 2; height = 4;
  463.   xmin2d = -2.0; xmax2d = 2.0;
  464.   ymin2d = -2.5; ymax2d = 5.0;
  465.   
  466.   _plenv (xmin2d, xmax2d, ymin2d, ymax2d, 0, -2);
  467.   _plw3d (basex, basey, height, xmin, xmax, ymin, ymax, ...
  468.           zmin, zmax, WIN.[P].alt[p], WIN.[P].az[p]);
  469.   _plbox3 (WIN.[P].grid3x[p], WIN.[P].xlabel[p], 0, 0, ...
  470.            WIN.[P].grid3y[p], WIN.[P].ylabel[p], 0, 0, ...
  471.            WIN.[P].grid3z[p], WIN.[P].zlabel[p], 0, 0);
  472.   _plmtex ("t", 1.0, 0.5, 0.5, WIN.[P].title[p]);
  473.  
  474.   if (exist (L31))
  475.   {
  476.     _plcol (2);
  477.     _plot3d (real(L31.x), real(L31.y), real(L31.z), ...
  478.              L31.x.n, L31.y.n, 3, 0);
  479.   }
  480.   if (exist (L32))
  481.   {
  482.     _plcol (3);
  483.     _pllsty (2);
  484.     _plot3d (real(L32.x), real(L32.y), real(L32.z), ...
  485.              L32.x.n, L32.y.n, 3, 0);
  486.   }
  487.   if (exist (L33)) 
  488.   {
  489.     _plcol (4);
  490.     _pllsty (3);
  491.     _plot3d (real(L33.x), real(L33.y), real(L33.z), ...
  492.              L33.x.n, L33.y.n, 3, 0);
  493.   }
  494.   
  495.   _plflush ();
  496.   _pltext ();
  497.   
  498.   #
  499.   # Increment the plot no. so that next time
  500.   # we use the correct settings.
  501.   #
  502.   
  503.   WIN.[P].subplot = WIN.[P].subplot + 1;
  504.   
  505.   return P;
  506. };
  507.  
  508. ##############################################################################
  509. #
  510. # Plot a 3-Dimensional graph. The data is composed in a list, with
  511. # elements `x', `y', and `z'. x and y are single-dimension arrays
  512. # (row or column matrices), and z is a two-dimensional array. The
  513. # array z, is a function of x and y: z = f(x,y). Thus, the values in
  514. # the array x can be thought of a "row-labels", and the values of y
  515. # can be thought of as "column-lables" for the 2-dimensioal array z.
  516. #
  517. # At present plmesh can plot 3 distinct lists. Each list may have
  518. # different X, Y, and Z scales.
  519. #
  520. ##############################################################################
  521.  
  522. plmesh = function ( L31, L32, L33 )
  523. {
  524.   local (p, xmin, xmax, ymin, ymax, zmin, zmax, basex, basey, height, ...
  525.          xmin2d, xmax2d, ymin2d, ymax2d, alt, az, ...
  526.          Xmin, Xmax, Ymin, Ymax, Zmin, Zmax)
  527.  
  528.   check_plot_object ();
  529.   
  530.   #
  531.   # 1st check list contents
  532.   #
  533.   
  534.   if (exist (L31)) { check_3d_list (L31); }
  535.   if (exist (L32)) { check_3d_list (L32); }
  536.   if (exist (L33)) { check_3d_list (L33); }
  537.   
  538.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  539.  
  540.   #
  541.   # Figure out the scale limits. 
  542.   # Needs improvement!
  543.   #
  544.   
  545.   xmin = xmax = ymin = ymax = zmin = zmax = 0;
  546.   if (exist (L31)) 
  547.   {
  548.     xyz_scales (L31, p, Xmin, Xmax, Ymin, Ymax, Zmin, Zmax);
  549.     if (Xmin < xmin) { xmin = Xmin; } if (Xmax > xmax) { xmax = Xmax; }
  550.     if (Ymin < ymin) { ymin = Ymin; } if (Ymax > ymax) { ymax = Ymax; }
  551.     if (Zmin < zmin) { zmin = Zmin; } if (Zmax > zmax) { zmax = Zmax; }
  552.   }
  553.   if (exist (L32)) 
  554.   {
  555.     xyz_scales (L32, p, Xmin, Xmax, Ymin, Ymax, Zmin, Zmax);
  556.     if (Xmin < xmin) { xmin = Xmin; } if (Xmax > xmax) { xmax = Xmax; }
  557.     if (Ymin < ymin) { ymin = Ymin; } if (Ymax > ymax) { ymax = Ymax; }
  558.     if (Zmin < zmin) { zmin = Zmin; } if (Zmax > zmax) { zmax = Zmax; }
  559.   }   
  560.   if (exist (L33)) 
  561.   {
  562.     xyz_scales (L33, p, Xmin, Xmax, Ymin, Ymax, Zmin, Zmax);
  563.     if (Xmin < xmin) { xmin = Xmin; } if (Xmax > xmax) { xmax = Xmax; }
  564.     if (Ymin < ymin) { ymin = Ymin; } if (Ymax > ymax) { ymax = Ymax; }
  565.     if (Zmin < zmin) { zmin = Zmin; } if (Zmax > zmax) { zmax = Zmax; }
  566.   }
  567.  
  568.   _plgra ();
  569.   _plcol (1);
  570.   _pllsty (1);
  571.   _plfont (WIN.[P].font[p]);
  572.   _plwid (WIN.[P].width[p]);
  573.   
  574.   basex = 2; basey = 2; height = 4;
  575.   xmin2d = -2.0; xmax2d = 2.0;
  576.   ymin2d = -2.5; ymax2d = 5.0;
  577.   
  578.   _plenv (xmin2d, xmax2d, ymin2d, ymax2d, 0, -2);
  579.   _plw3d (basex, basey, height, xmin, xmax, ymin, ymax, ...
  580.           zmin, zmax, WIN.[P].alt[p], WIN.[P].az[p]);
  581.   _plbox3 (WIN.[P].grid3x[p], WIN.[P].xlabel[p], 0, 0, ...
  582.            WIN.[P].grid3y[p], WIN.[P].ylabel[p], 0, 0, ...
  583.            WIN.[P].grid3z[p], WIN.[P].zlabel[p], 0, 0);
  584.   _plmtex ("t", 1.0, 0.5, 0.5, WIN.[P].title[p]);
  585.  
  586.   if (exist (L31))
  587.   {
  588.     _plcol (2);
  589.     _plmesh (real(L31.x), real(L31.y), real(L31.z), ...
  590.              L31.x.n, L31.y.n, 3);
  591.   }
  592.   if (exist (L32))
  593.   {
  594.     _plcol (3);
  595.     _pllsty (2);
  596.     _plmesh (real(L32.x), real(L32.y), real(L32.z), ...
  597.              L32.x.n, L32.y.n, 3);
  598.   }
  599.   if (exist (L33)) 
  600.   {
  601.     _plcol (4);
  602.     _pllsty (3);
  603.     _plmesh (real(L33.x), real(L33.y), real(L33.z), ...
  604.              L33.x.n, L33.y.n, 3);
  605.   }
  606.  
  607.   _plflush ();
  608.   _pltext ();
  609.   
  610.   #
  611.   # Increment the plot no. so that next time
  612.   # we use the correct settings.
  613.   #
  614.   
  615.   WIN.[P].subplot = WIN.[P].subplot + 1;
  616.   
  617.   return P;
  618. };
  619.  
  620. ##############################################################################
  621. #
  622. # Plot contours. The data is composed in a list, with
  623. # elements `x', `y', and `z'. x and y are single-dimension arrays
  624. # (row or column matrices), and z is a two-dimensional array. The
  625. # array z, is a function of x and y: z = f(x,y). Thus, the values in
  626. # the array x can be thought of a "row-labels", and the values of y
  627. # can be thought of as "column-lables" for the 2-dimensioal array z.
  628. #
  629. ##############################################################################
  630.  
  631. plcont = function ( CL )
  632. {
  633.   local (p, xmin, xmax, ymin, ymax, zmin, zmax, ...
  634.          Xmin, Xmax, Ymin, Ymax, Zmin, Zmax, i, v, ...
  635.          xl, yl, clevel, llevel, k, j, l)
  636.      
  637.   check_plot_object ();
  638.  
  639.   #
  640.   # 1st check list contents
  641.   #
  642.   
  643.   if (exist (CL)) { check_3d_list (CL); }
  644.  
  645.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  646.   
  647.   #
  648.   # Figure out the scale limits. 
  649.   # Needs improvement!
  650.   #
  651.  
  652.   xmin = xmax = ymin = ymax = zmin = zmax = 0;
  653.   if (exist (CL)) 
  654.   {
  655.     xyz_scales (CL, p, Xmin, Xmax, Ymin, Ymax, Zmin, Zmax);
  656.     if (Xmin < xmin) { xmin = Xmin; } if (Xmax > xmax) { xmax = Xmax; }
  657.     if (Ymin < ymin) { ymin = Ymin; } if (Ymax > ymax) { ymax = Ymax; }
  658.     if (Zmin < zmin) { zmin = Zmin; } if (Zmax > zmax) { zmax = Zmax; }
  659.   }
  660.  
  661.   _plgra ();
  662.   _plcol (1);
  663.   _pllsty (1);
  664.   _plfont (WIN.[P].font[p]);
  665.   _plwid (WIN.[P].width[p]);
  666.  
  667.   #
  668.   # Set up the 1st viewport for drawing the plot.
  669.   #
  670.  
  671.   _pladv (0);
  672.   _plvpor (0.15, 0.75, 0.15, 0.85);
  673.   _plwind (xmin, xmax, ymin, ymax);
  674.   _plbox (WIN.[P].gridx[p], 0, 0, WIN.[P].gridy[p], 0, 0);
  675.  
  676.   if (exist (CL.clevel))
  677.   {
  678.     clevel = CL.clevel;
  679.   else
  680.     clevel = linspace(zmin, zmax, 10);
  681.   }
  682.  
  683.   #
  684.   # Draw the contours
  685.   #
  686.  
  687.   l = 1;
  688.   for (i in 1:clevel.n)
  689.   {
  690.     k = mod (i-1, 14) + 1;
  691.     j = mod (i-1, 8) + 1;
  692.     _pllsty(j);
  693.     _plcol (1+k);
  694.     if (_plcont (CL.x, CL.y, CL.z, 1, CL.x.n, 1, CL.y.n, clevel[i]))
  695.     {
  696.       llevel[l] = clevel[i];
  697.       l = l + 1;
  698.     }
  699.   }
  700.  
  701.   #
  702.   # Reset color and draw the labels.
  703.   #
  704.  
  705.   _plcol (1);
  706.   _pllab (WIN.[P].xlabel[p], WIN.[P].ylabel[p], WIN.[P].title[p]);  
  707.  
  708.   #
  709.   # Draw the contour legend. Use a new viewport to the right
  710.   # of the contour plot.
  711.   #
  712.  
  713.   _plvpor (0.75, 1.0, 0.15, 0.85);
  714.   _plwind (0, 1, 0, 1);
  715.  
  716.   v = 1 - 1/(2*llevel.n);
  717.  
  718.   for (i in 1:llevel.n)
  719.   {
  720.     xl = [0.1, 0.2, 0.3]';
  721.     yl = [v, v, v]';
  722.     v = v - 1/llevel.n;
  723.  
  724.     k = mod (i-1, 14) + 1;
  725.     j = mod (i-1, 8) + 1;
  726.  
  727.     _plcol (1+k);
  728.     _pllsty (j);
  729.  
  730.     _plline (3, xl, yl);
  731.     plptex (num2str (llevel[i]), xl[3]+.1, yl[3], , , 0);
  732.   }
  733.  
  734.   # Flush  and go back to text mode.
  735.   _plflush ();
  736.   _pltext ();
  737.   
  738.   #
  739.   # Increment the plot no. so that next time
  740.   # we use the correct settings.
  741.   #
  742.   
  743.   WIN.[P].subplot = WIN.[P].subplot + 1;
  744.   
  745.   return P;
  746. };
  747.  
  748. #
  749. # error bar plot
  750. #
  751.  
  752. plerry = function (x, y, y_low, y_high)
  753. {
  754.   local (hscale, xmin, xmax, ymin, ymax, i, p)
  755.   
  756.   check_plot_object ();
  757.   
  758.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  759.   WIN.[P].desc.[p] = inf(); // use own legend
  760.   
  761.   if (x.nr != y.nr || x.nr != y_low.nr || x.nr != y_high.nr) 
  762.   {
  763.     error(" Size inconsistent in plerry.");
  764.   }
  765.   _plgra ();
  766.   _plcol (1);
  767.   _pllsty (1);
  768.   _plfont (WIN.[P].font[p]);
  769.   _plwid (WIN.[P].width[p]);
  770.   xy_scales ( real([x,y,y_low,y_high]), p, xmin, xmax, ymin, ymax );
  771.   
  772.   _pladv (0);
  773.   _plvsta ();
  774.   _plwind (xmin, xmax, ymin, ymax);
  775.   _plbox (WIN.[P].gridx[p], 0, 0, WIN.[P].gridy[p], 0, 0);
  776.  
  777.   if (plot_matrix ( [x,y], p, 0, xmin, xmax, ymin, ymax, ymax-ymin ) < 0) 
  778.   { 
  779.     return -1; 
  780.   }
  781.  
  782.   _plcol (3);
  783.   _plerry(x.nr, x, y_low, y_high);
  784.   _plcol (1);
  785.   _pllab (WIN.[P].xlabel[p], WIN.[P].ylabel[p], WIN.[P].title[p]);  
  786.   _plflush ();
  787.   _pltext ();
  788.   
  789.   #
  790.   # Increment the plot no. so that next time
  791.   # we use the correct settings.
  792.   #
  793.   
  794.   WIN.[P].subplot = WIN.[P].subplot + 1;
  795.   return P;      
  796. };
  797.  
  798. ##############################################################################
  799. #
  800. # Plot a Histogram(s), from the columns of a matrix.
  801. #
  802. ##############################################################################
  803.  
  804. plhist = function ( M , nbin )
  805. {
  806.   local (hscale, i, k, np, p, ymin, ymax, xmax, v, xl, yl, desc)
  807.  
  808.   check_plot_object ();
  809.   
  810.   if (!exist (nbin)) { nbin = 10; }
  811.   
  812.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  813.   np = M.nr;
  814.   
  815.   # Compute max/min values of data
  816.   
  817.   ymin = min (min (real (M)));
  818.   ymax = max (max (real (M)));
  819.   
  820.   #
  821.   # Check computed scale limits against user's
  822.   #
  823.   
  824.   if (!isinf (WIN.[P].ymin[p])) { ymin = WIN.[P].ymin[p]; }
  825.   if (!isinf (WIN.[P].ymax[p])) { ymax = WIN.[P].ymax[p]; }
  826.  
  827.   _plgra ();
  828.   _plcol (1);
  829.   _plfont (WIN.[P].font[p]);
  830.   _plwid (WIN.[P].width[p]);
  831.   
  832.   for (i in 1:M.nc) 
  833.   { 
  834.     hscale[i] = hist_scales (M[;i], nbin);
  835.   }  
  836.  
  837.   _pladv (0);
  838.   _plvsta ();
  839.   _plwind (ymin, ymax, 0, max (hscale));
  840.   _plbox (WIN.[P].gridx[p], 0, 0, WIN.[P].gridy[p], 0, 0);
  841.  
  842.   v = max (hscale);
  843.   xmax = ymax;
  844.   for (i in 1:M.nc)
  845.   {
  846.     k = mod (i, 14) + 1;
  847.     _plcol (WIN.[P].color[p;k]);
  848.     _plhist (np, real(M[;i]), ymin, ymax, nbin, 1);
  849.     
  850.     if (!any (any (isinf (WIN.[P].desc.[p]))))
  851.     {
  852.       # Use the default if necessary
  853.       if (WIN.[P].desc.[p][1] == "default") 
  854.       {
  855.     desc = "c1";
  856.     else if (WIN.[P].desc.[p].n >= i) {
  857.       desc = WIN.[P].desc.[p][i];
  858.     else
  859.       # Not sure what to do, user has messed up.
  860.       desc = "";
  861.       } }
  862.       
  863.       v = v - max(hscale)/11;
  864.       xl = (ymax-ymin)*[10.5/12, 11/12, 11.5/12]' + ymin;
  865.       yl = [v, v, v]';
  866.  
  867.       _plline (3, xl, yl);
  868.       plptex(desc, xl[1]-(ymax-ymin)/25, yl[3], , , 1);
  869.     }
  870.   }
  871.  
  872.   _plcol (1);
  873.   _pllab (WIN.[P].xlabel[p], WIN.[P].ylabel[p], WIN.[P].title[p]);  
  874.   _plflush ();
  875.   _pltext ();
  876.  
  877.   #
  878.   # Increment the plot no. so that next time
  879.   # we use the correct settings.
  880.   #
  881.   
  882.   WIN.[P].subplot = WIN.[P].subplot + 1;
  883.  
  884.   return 1;
  885. };
  886.  
  887. ##############################################################################
  888. #
  889. # Various support functions for the WIN list
  890. #
  891. ##############################################################################
  892.  
  893. #
  894. # Replot
  895. #
  896.  
  897. replot = function ( )
  898. {
  899.   check_plot_object ();
  900.   _replot ();
  901. };
  902.  
  903. #
  904. # Set the X-axis label
  905. #
  906.  
  907. xlabel = function ( xstr )
  908. {
  909.   local (i);
  910.   
  911.   check_plot_object ();
  912.   if (!exist (xstr)) { xstr = ""; }
  913.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  914.   WIN.[P].xlabel[i] = xstr;
  915. };
  916.  
  917. #
  918. # Set the Y-axis label
  919. #
  920.  
  921. ylabel = function ( xstr )
  922. {
  923.   local (i);
  924.   
  925.   check_plot_object ();
  926.   if (!exist (xstr)) { xstr = ""; }
  927.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  928.   WIN.[P].ylabel[i] = xstr;
  929. };
  930.  
  931. #
  932. # Set the Z-axis label
  933. #
  934.  
  935. zlabel = function ( xstr )
  936. {
  937.   local (i);
  938.   
  939.   check_plot_object ();
  940.   if (!exist (xstr)) { xstr = ""; }
  941.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  942.   WIN.[P].zlabel[i] = xstr;
  943. };
  944.  
  945. #
  946. # Set the plot-title
  947. #
  948.  
  949. ptitle = function ( xstr )
  950. {
  951.   local (p);
  952.   
  953.   check_plot_object ();
  954.   if (!exist (xstr)) { xstr = ""; }
  955.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  956.   WIN.[P].title[p] = xstr;
  957. };
  958.  
  959. plimits = function ( xmin, xmax, ymin, ymax, zmin, zmax )
  960. {
  961.   local (i);
  962.  
  963.   check_plot_object ();
  964.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  965.   if (exist (xmin)) {
  966.     WIN.[P].xmin[i] = xmin;
  967.   else
  968.     WIN.[P].xmin[i] = inf ();
  969.   }
  970.   if (exist (xmax)) { 
  971.     WIN.[P].xmax[i] = xmax;
  972.   else
  973.     WIN.[P].xmax[i] = inf ();
  974.   }
  975.   if (exist (ymin)) {
  976.     WIN.[P].ymin[i] = ymin;
  977.   else
  978.     WIN.[P].ymin[i] = inf ();
  979.   }
  980.   if (exist (ymax)) {
  981.     WIN.[P].ymax[i] = ymax;
  982.   else
  983.     WIN.[P].ymax[i] = inf ();
  984.   }
  985.   if (exist (zmin)) {
  986.     WIN.[P].zmin[i] = zmin;
  987.   else
  988.     WIN.[P].zmin[i] = inf ();
  989.   }
  990.   if (exist (zmax)) {
  991.     WIN.[P].zmax[i] = zmax;
  992.   else
  993.     WIN.[P].zmax[i] = inf ();
  994.   }
  995. };
  996.  
  997. plgrid = function ( sty_x, sty_y )
  998. {
  999.   local (i)
  1000.  
  1001.   check_plot_object ();
  1002.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  1003.  
  1004.   if (exist (sty_x)) 
  1005.   { 
  1006.     if (class (sty_x) == "string")
  1007.     {
  1008.       WIN.[P].gridx[i] = sty_x;
  1009.     else
  1010.       error ("plgrid: requires string argument GRID_STY_X");
  1011.     }
  1012.   else
  1013.     WIN.[P].gridx[i] = grid_x_default;
  1014.   }
  1015.   if (exist (sty_y)) 
  1016.   { 
  1017.     if (class (sty_y) == "string")
  1018.     {
  1019.       WIN.[P].gridy[i] = sty_y;
  1020.     else
  1021.       error ("plgrid: requires string argument GRID_STY_Y");
  1022.     }
  1023.   else
  1024.     WIN.[P].gridy[i] = grid_y_default;
  1025.   }
  1026. };
  1027.  
  1028.  plgrid3 = function ( sty_x, sty_y, sty_z )
  1029.  {
  1030.    local (i)
  1031.  
  1032.    check_plot_object ();
  1033.    i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  1034.    if (exist (sty_x)) 
  1035.    { 
  1036.      if (class (sty_x) == "string")
  1037.      {
  1038.        WIN.[P].grid3x[i] = sty_x;
  1039.      else
  1040.        error ("plgrid3: requires string argument GRID_STY_X");
  1041.      }
  1042.    else
  1043.      WIN.[P].grid3x[i] = grid_3x_default;
  1044.    }
  1045.    if (exist (sty_y)) 
  1046.    { 
  1047.      if (class (sty_y) == "string")
  1048.      {
  1049.        WIN.[P].grid3y[i] = sty_y;
  1050.      else
  1051.        error ("plgrid3: requires string argument GRID_STY_Y");
  1052.      }
  1053.    else
  1054.      WIN.[P].grid3y[i] = grid_3y_default;
  1055.    }
  1056.    if (exist (sty_z)) 
  1057.    { 
  1058.      if (class (sty_z) == "string")
  1059.      {
  1060.        WIN.[P].grid3z[i] = sty_z;
  1061.      else
  1062.        error ("plgrid3: requires string argument GRID_STY_Z");
  1063.      }
  1064.    else
  1065.      WIN.[P].grid3z[i] = grid_3z_default;
  1066.    }
  1067.  };
  1068.  
  1069.  #
  1070.  # A friendlier interface to changing grid/axis
  1071.  # styles.
  1072.  #
  1073.  
  1074.  plaxis = function ( X_STR, Y_STR )
  1075.  {
  1076.    local (i)
  1077.  
  1078.    check_plot_object ();
  1079.    i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  1080.  
  1081.    if (exist (X_STR))
  1082.    {
  1083.      if (X_STR == "log") { WIN.[P].gridx[i] = "bcngstl"; }
  1084.    else
  1085.      WIN.[P].gridx[i] = grid_x_default;
  1086.    }
  1087.  
  1088.    if (exist (Y_STR))
  1089.    {
  1090.      if (Y_STR == "log") { WIN.[P].gridy[i] = "bcngstlv"; }
  1091.    else
  1092.      WIN.[P].gridy[i] = grid_y_default;
  1093.    }
  1094.    return P;
  1095.  };
  1096.  
  1097.  #
  1098.  # Change plot aspect ratio
  1099.  #
  1100.  
  1101.  plaspect = function ( style )
  1102.  {
  1103.    local (i)
  1104.  
  1105.    check_plot_object ();
  1106.    i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  1107.    if (exist (style)) { WIN.[P].aspect[i] = style; }
  1108.  };
  1109.  
  1110.  #
  1111.  # Change plot line style
  1112.  #
  1113.  
  1114.  plstyle = function ( style )
  1115.  {
  1116.    local (i)
  1117.  
  1118.    check_plot_object ();
  1119.    i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  1120.    if (exist (style)) 
  1121.    {
  1122.      if (class (style) == "string") 
  1123.      {
  1124.        WIN.[P].style.[i] = style;
  1125.      }
  1126.      return 1;
  1127.    }
  1128.    WIN.[P].style[i] = "line";  
  1129.  };
  1130.  
  1131.  #
  1132.  # Get the right value of line-style
  1133.  #
  1134.  
  1135.  get_style = function ( STY, K )
  1136.  {
  1137.    local (sty);
  1138.    sty = mod(K, STY.n);
  1139.    if(sty == 0) { sty = STY.n; }
  1140.    return STY[sty];
  1141.  };
  1142.  
  1143.  #
  1144.  # Change fonts
  1145.  #
  1146.  
  1147.  plfont = function ( font )
  1148.  {
  1149.    local (i)
  1150.  
  1151.    check_plot_object ();
  1152.    i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  1153.  
  1154.    if (!exist (font)) { font = 1; }
  1155.  
  1156.    if (WIN.[P].fontld == 0)
  1157.    {
  1158.      _plfontld (1);
  1159.      WIN.[P].fontld = 1;
  1160.    }
  1161.  
  1162.    WIN.[P].font[i] = font;
  1163.    return P;
  1164.  };
  1165.  
  1166.  #
  1167.  # Change pen width
  1168.  #
  1169.  
  1170.  plwid = function ( width )
  1171.  {
  1172.    local (i)
  1173.  
  1174.    check_plot_object ();
  1175.    i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  1176.  
  1177.    if (!exist (width)) { width = 1; }
  1178.    WIN.[P].width[i] = width;
  1179.    return P;
  1180.  };
  1181.  
  1182.  #
  1183.  # Place some text on the plot
  1184.  #
  1185.  
  1186.  plptex = function ( text, x , y , dx , dy , just )
  1187.  {
  1188.    if (!check_plot_object ()) {
  1189.      printf ("Must use plot() before plptex()\n");
  1190.      return 0;
  1191.    }
  1192.  
  1193.    if (!exist (x)) { x = 0; }
  1194.    if (!exist (y)) { y = 0; }
  1195.    if (!exist (dx)) { dx = abs(x)+1; }
  1196.    if (!exist (dy)) { dy = 0; }
  1197.    if (!exist (just)) { just = 0; }
  1198.  
  1199.    _plptex (x, y, dx, dy, just, text);
  1200.  };
  1201.  
  1202.  #
  1203.  # Set up the viewing altitude for 3-D plots
  1204.  #
  1205.  
  1206.  plalt = function ( ALT )
  1207.  {
  1208.    local (i)
  1209.  
  1210.    check_plot_object ();
  1211.    i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  1212.    if (exist (ALT)) 
  1213.    { 
  1214.      WIN.[P].alt[i] = ALT; 
  1215.    else
  1216.      WIN.[P].alt[i] = 60;
  1217.    }
  1218.    return P;
  1219.  };
  1220.  
  1221.  #
  1222.  # Set the viewing azimuth for 3-D plots
  1223.  #
  1224.  
  1225.  plaz = function ( AZ )
  1226.  {
  1227.    local (i)
  1228.  
  1229.    check_plot_object ();
  1230.    i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  1231.    if (exist (AZ)) 
  1232.    { 
  1233.      WIN.[P].az[i] = AZ; 
  1234.    else
  1235.      WIN.[P].az[i] = 45;
  1236.    }
  1237.    return P;
  1238.  };
  1239.  
  1240.  ##############################################################################
  1241.  #
  1242.  # Various internal support functions. Eventually these will be static.
  1243.  #
  1244.  ##############################################################################
  1245.  
  1246.  #
  1247.  # Find the X and Y scales for a single matrix.
  1248.  #
  1249.  
  1250.  xy_scales = function ( M, p, xmin, xmax, ymin, ymax )
  1251.  {
  1252.  
  1253.    #
  1254.    # 1st check for un-plottable data
  1255.    #
  1256.  
  1257.    if (any (any (isinf (M))))
  1258.      { error ("plot: cannot plot infs"); }
  1259.    if (any (any (isnan (M))))
  1260.      { error ("plot: cannot plot NaNs"); }
  1261.  
  1262.    if (M.nc == 1)
  1263.    {
  1264.      xmin = 1;
  1265.      xmax = M.nr;
  1266.      ymin = min (M);
  1267.      ymax = max (M);
  1268.    else
  1269.      xmin = min (M[;1]);
  1270.      xmax = max (M[;1]);
  1271.      ymin = min (min (M[;2:M.nc]));
  1272.      ymax = max (max (M[;2:M.nc]));
  1273.    }
  1274.  
  1275.    #
  1276.    # Check computed scale limits against user's
  1277.    #
  1278.  
  1279.    if (!isinf (WIN.[P].xmin[p])) { xmin = WIN.[P].xmin[p]; }
  1280.    if (!isinf (WIN.[P].xmax[p])) { xmax = WIN.[P].xmax[p]; }
  1281.    if (!isinf (WIN.[P].ymin[p])) { ymin = WIN.[P].ymin[p]; }
  1282.    if (!isinf (WIN.[P].ymax[p])) { ymax = WIN.[P].ymax[p]; }
  1283.  
  1284.    #
  1285.    # Check for potential errors
  1286.    #
  1287.  
  1288.    if (xmin == xmax) 
  1289.    { 
  1290.      xmin = xmin - 1;
  1291.      xmax = xmax + 1;
  1292.    }
  1293.  
  1294.    if (ymin == ymax)
  1295.    {
  1296.      ymin = ymin - 1;
  1297.      ymax = ymax + 1;
  1298.    }
  1299.  
  1300.    #
  1301.    # Finally, adjust if log-scales
  1302.    #
  1303.  
  1304.    if (find_char (WIN.[P].gridx[p], "l"))
  1305.    {
  1306.      if (xmin <= 0 || xmax <= 0) { error ("plot: cannot plot log <= 0"); }
  1307.      xmin = log10 (xmin);
  1308.      xmax = log10 (xmax);
  1309.    }
  1310.    if (find_char (WIN.[P].gridy[p], "l"))
  1311.    {
  1312.      if (ymin <= 0 || ymax <= 0) { error ("plot: cannot plot log <= 0"); }
  1313.      ymin = log10 (ymin);
  1314.      ymax = log10 (ymax);
  1315.    }
  1316.  
  1317.    return 1;
  1318.  };
  1319.  
  1320.  #
  1321.  # Find the X, Y and Z scales for a single matrix.
  1322.  #
  1323.  
  1324.  xyz_scales = function ( L, p, xmin, xmax, ymin, ymax, zmin, zmax )
  1325.  {
  1326.    # X - scale
  1327.    if (any (any (isinf (L.x))))
  1328.      { error ("plot3: cannot plot infs"); }
  1329.    if (any (any (isnan (L.x))))
  1330.      { error ("plot3: cannot plot NaNs"); }
  1331.    
  1332.    xmin = min (real (L.x));
  1333.    xmax = max (real (L.x));
  1334.  
  1335.    # Y - scale
  1336.    if (any (any (isinf (L.y))))
  1337.      { error ("plot3: cannot plot infs"); }
  1338.    if (any (any (isnan (L.y))))
  1339.      { error ("plot3: cannot plot NaNs"); }
  1340.    
  1341.    ymin = min (real (L.y));
  1342.    ymax = max (real (L.y));
  1343.  
  1344.    # Z - scale
  1345.    if (any (any (isinf (L.z))))
  1346.      { error ("plot3: cannot plot infs"); }
  1347.    if (any (any (isnan (L.z))))
  1348.      { error ("plot3: cannot plot NaNs"); }
  1349.  
  1350.    zmin = min (min (real (L.z)));
  1351.    zmax = max (max (real (L.z)));
  1352.  
  1353.    #
  1354.    # Check computed scale limits against user's
  1355.    #
  1356.  
  1357.    if (!isinf (WIN.[P].xmin[p])) { xmin = WIN.[P].xmin[p]; }
  1358.    if (!isinf (WIN.[P].xmax[p])) { xmax = WIN.[P].xmax[p]; }
  1359.    if (!isinf (WIN.[P].ymin[p])) { ymin = WIN.[P].ymin[p]; }
  1360.    if (!isinf (WIN.[P].ymax[p])) { ymax = WIN.[P].ymax[p]; }
  1361.    if (!isinf (WIN.[P].zmin[p])) { zmin = WIN.[P].zmin[p]; }
  1362.    if (!isinf (WIN.[P].zmax[p])) { zmax = WIN.[P].zmax[p]; }
  1363.  
  1364.    return 1;
  1365.  };
  1366.  
  1367.  #
  1368.  # Find the X and Y scales for a list of matrices
  1369.  #
  1370.  
  1371.  list_scales = function ( data, p, Xmin, Xmax, Ymin, Ymax )
  1372.  {
  1373.    local (M, i, xmin, xmax, ymin, ymax, once)
  1374.    once = 1;
  1375.  
  1376.    for (i in members (data))
  1377.    {
  1378.      M = real (data.[i]);
  1379.      if (class (M) != "num") { continue; }
  1380.  
  1381.      #
  1382.      # 1st check for un-plottable data
  1383.      #
  1384.  
  1385.      if (any (any (isinf (M))))
  1386.        { error ("plot: cannot plot infs"); }
  1387.      if (any (any (isnan (M))))
  1388.        { error ("plot: cannot plot NaNs"); }
  1389.  
  1390.      if (M.nc == 1)
  1391.      {
  1392.        xmin = 1;
  1393.        xmax = M.nr;
  1394.        ymin = min (M);
  1395.        ymax = max (M);
  1396.      else
  1397.        xmin = min (M[;1]);
  1398.        xmax = max (M[;1]);
  1399.        ymin = min (min (M[;2:M.nc]));
  1400.        ymax = max (max (M[;2:M.nc]));
  1401.      }
  1402.      if (once) { 
  1403.        Xmin = xmin; Xmax = xmax; Ymin = ymin; Ymax = ymax; 
  1404.        once = 0; 
  1405.      }
  1406.      if (xmin < Xmin) { Xmin = xmin; }
  1407.      if (xmax > Xmax) { Xmax = xmax; }
  1408.      if (ymin < Ymin) { Ymin = ymin; }
  1409.      if (ymax > Ymax) { Ymax = ymax; }
  1410.    }
  1411.  
  1412.    #
  1413.    # Check computed scale limits against user's
  1414.    #
  1415.  
  1416.    if (!isinf (WIN.[P].xmin[p])) { Xmin = WIN.[P].xmin[p]; }
  1417.    if (!isinf (WIN.[P].xmax[p])) { Xmax = WIN.[P].xmax[p]; }
  1418.    if (!isinf (WIN.[P].ymin[p])) { Ymin = WIN.[P].ymin[p]; }
  1419.    if (!isinf (WIN.[P].ymax[p])) { Ymax = WIN.[P].ymax[p]; }
  1420.  
  1421.    #
  1422.    # Finally, adjust if log-scales
  1423.    #
  1424.  
  1425.    if (find_char (WIN.[P].gridx[p], "l"))
  1426.    {
  1427.      if (Xmin <= 0 || Xmax <= 0) { error ("plot: cannot plot log x <= 0"); }
  1428.      Xmin = log10 (Xmin);
  1429.      Xmax = log10 (Xmax);
  1430.    }
  1431.    if (find_char (WIN.[P].gridy[p], "l"))
  1432.    {
  1433.      if (Ymin <= 0 || Ymax <= 0) { error ("plot: cannot plot log y <= 0"); }
  1434.      Ymin = log10 (Ymin);
  1435.      Ymax = log10 (Ymax);
  1436.    }
  1437.  
  1438.    return 1;
  1439.  };
  1440.    
  1441.  #
  1442.  # Find the maximum number of elements in a bin for a single 
  1443.  # column matrix.
  1444.  #
  1445.  
  1446.  hist_scales = function ( data, nbin )
  1447.  {
  1448.    local (i, binval, dmin, dmax, dbin)
  1449.  
  1450.    dmin = min (real (data));
  1451.    dmax = max (real (data));
  1452.    dbin = linspace (dmin, dmax, nbin+1);
  1453.    binval = zeros (nbin, 1);
  1454.  
  1455.    for (i in 1:nbin)
  1456.    {
  1457.      binval[i] = length (find (data >= dbin[i] && data < dbin[i+1]));
  1458.    }
  1459.  
  1460.    return max (binval);
  1461.  };
  1462.  
  1463.  #
  1464.  # Plot the columns of a matrix (core function)
  1465.  #
  1466.  # Notes: This is the core function for plotting a matrix. If the
  1467.  # matrix is a single column, then the matrix elements are plotted
  1468.  # versus the row numbers. If it is a multi-column matrix, then
  1469.  # columns 2:N are plotted versus column 1.
  1470.  #
  1471.  # p, K, k and l are indices for plot features.
  1472.  # p: the current plot index (the plot #)
  1473.  # K: usually 0. This index is used to start of the line style and
  1474.  # color index (k = color index, l = line-style index). This is mostly
  1475.  # used by plot_list, which may call plot_matrix repeatedly.
  1476.  # k: the line color index. This value determines the line color used
  1477.  # for each column of data. If K = 0, then k goes like 2:14, then
  1478.  # flops back to 1:14.
  1479.  # l: the line style inex. This value determines the line style used
  1480.  # for each column of data - not the line-type (points, or lines). If
  1481.  # K = 0, then l goes like 2:8, then flops back to 1:8.
  1482.  #
  1483.  
  1484.  plot_matrix = function ( M, p, K, xmin, xmax, ymin, ymax, v )
  1485.  {
  1486.    local (ans, i, k, l, np, x, y, xl, yl, desc)
  1487.  
  1488.    np = M.nr;
  1489.   
  1490.    if (M.nc == 1)
  1491.    {
  1492.      x = 1:M.nr;
  1493.      y = real (M);
  1494.      k = mod (1+K, 14) + 1;
  1495.      l = mod (1+K, 8) + 1;
  1496.  
  1497.      if (find_char (WIN.[P].gridx[p], "l"))
  1498.        { x = log10 (x); }
  1499.      if (find_char (WIN.[P].gridy[p], "l"))
  1500.        { y = log10 (y); }
  1501.  
  1502.      _plcol (WIN.[P].color[p;k]);
  1503.      _pllsty (WIN.[P].lstyle[p;l]);
  1504.  
  1505.      if (get_style (WIN.[P].style.[p], k-1) == "line") 
  1506.      {
  1507.        _plline (M.nr, x, y);
  1508.      else if (get_style (WIN.[P].style.[p], k-1) == "point") {
  1509.        _plpoin (M.nr, x, y, WIN.[P].pstyle[p]+k);
  1510.      else if (get_style (WIN.[P].style.[p], k-1) == "line-point") {
  1511.        _plline (M.nr, x, y);
  1512.        _plpoin (M.nr, x, y, WIN.[P].pstyle[p]+k);
  1513.      else {
  1514.        _plline (M.nr, x, y);
  1515.      }}}}
  1516.  
  1517.      #
  1518.      # Now do the legend 
  1519.      #
  1520.  
  1521.      if (!any (any (isinf (WIN.[P].desc.[p]))))
  1522.      {
  1523.        # Use the default if necessary
  1524.        if (WIN.[P].desc.[p][1] == "default") 
  1525.        {
  1526.          desc = "c1";
  1527.        else if (WIN.[P].desc.[p].n >= k-1) {
  1528.          desc = WIN.[P].desc.[p][k-1];
  1529.        else
  1530.          # Not sure what to do, user has messed up.
  1531.          desc = "";
  1532.        }}
  1533.              
  1534.        v = v - (ymax-ymin)/11;
  1535.        xl = (xmax-xmin)*[10.5/12, 11/12, 11.5/12]' + xmin;
  1536.        yl = [v, v, v]' + ymin;
  1537.  
  1538.        if (get_style (WIN.[P].style.[p], k-1) == "line") 
  1539.        {
  1540.          _plline (3, xl, yl);
  1541.        else if (get_style (WIN.[P].style.[p], k-1) == "point") {
  1542.          _plpoin (3, xl, yl, WIN.[P].pstyle[p]+k);
  1543.        else if (get_style (WIN.[P].style.[p]) == "line-point") {
  1544.          _plline (3, xl, yl);
  1545.          _plpoin (3, xl, yl, WIN.[P].pstyle[p]+k);
  1546.        }}}
  1547.  
  1548.        plptex(desc, xl[1]-(xmax-xmin)/25, yl[3], , , 1);
  1549.  
  1550.      }
  1551.  
  1552.    else
  1553.  
  1554.      #
  1555.      # Check for large column dimension
  1556.      #
  1557.  
  1558.      if (M.nc > 3*M.nr)
  1559.      {
  1560.        printf (" Plot %i columns and %i rows, are you sure (y/n) ? "...
  1561.                , M.nc, M.nr);
  1562.        ans = getline ("stdin");
  1563.        if (ans.[1] != "y") { return -1; }
  1564.      }
  1565.  
  1566.      for (i in 2:M.nc)
  1567.      {
  1568.        x = real (M[;1]);
  1569.        y = real (M[;i]);
  1570.        if (find_char (WIN.[P].gridx[p], "l"))
  1571.          { x = log10 (x); }
  1572.        if (find_char (WIN.[P].gridy[p], "l"))
  1573.          { y = log10 (y); }
  1574.  
  1575.        k = mod (i-1 + K, 14) + 1;
  1576.        l = mod (i-1 + K, 8) + 1;
  1577.  
  1578.        _plcol (WIN.[P].color[p;k]);
  1579.        _pllsty (WIN.[P].lstyle[p;l]);
  1580.  
  1581.        if (get_style (WIN.[P].style.[p], k-1) == "line") 
  1582.        {
  1583.          _plline (np, x, y);
  1584.        else if (get_style (WIN.[P].style.[p], k-1) == "point") {
  1585.          _plpoin (np, x, y, WIN.[P].pstyle[p]+k);
  1586.        else if (get_style (WIN.[P].style.[p], k-1) == "line-point") {
  1587.          _plline (np, x, y);
  1588.          _plpoin (np, x, y, WIN.[P].pstyle[p]+k);
  1589.        else {
  1590.          _plline (np, x, y);
  1591.        }}}}
  1592.  
  1593.        #
  1594.        # Now do the legend 
  1595.        #
  1596.  
  1597.        if (!any (any (isinf (WIN.[P].desc.[p]))))
  1598.        {
  1599.          # Use the default if necessary
  1600.          if (WIN.[P].desc.[p][1] == "default") 
  1601.          {
  1602.            desc = "c" + num2str (i);
  1603.          else if (WIN.[P].desc.[p].n >= k-1) {
  1604.            desc = WIN.[P].desc.[p][k-1];
  1605.          else
  1606.            # Not sure what to do, user has messed up.
  1607.            desc = "";
  1608.          }}
  1609.              
  1610.          v = v - (ymax-ymin)/11;
  1611.          xl = (xmax-xmin)*[10.5/12, 11/12, 11.5/12]' + xmin;
  1612.          yl = [v, v, v]' + ymin;
  1613.  
  1614.          if (get_style (WIN.[P].style.[p], k-1) == "line") 
  1615.          {
  1616.            _plline (3, xl, yl);
  1617.          else if (get_style (WIN.[P].style.[p], k-1) == "point") {
  1618.            _plpoin (3, xl, yl, WIN.[P].pstyle[p]+k);
  1619.          else if (get_style (WIN.[P].style.[p]) == "line-point") {
  1620.            _plline (3, xl, yl);
  1621.            _plpoin (3, xl, yl, WIN.[P].pstyle[p]+k);
  1622.          }}}
  1623.  
  1624.          plptex(desc, xl[1]-(xmax-xmin)/25, yl[3], , , 1);
  1625.  
  1626.        }
  1627.      }
  1628.    }
  1629.  
  1630.   return k-1;
  1631.  };
  1632.  
  1633.  #
  1634.  # Plot all of the matrices in a list on the same plot
  1635.  #
  1636.  
  1637.  plot_list = function ( L , p, xmin, xmax, ymin, ymax )
  1638.  {
  1639.    local (M, i, k, v)
  1640.    k = 0;
  1641.    v = ymax - ymin;
  1642.    for (i in members (L))
  1643.    {
  1644.      M = L.[i];
  1645.      if (class (M) != "num") { continue; }
  1646.      if ((k = plot_matrix (M, p, k, xmin, xmax, ymin, ymax, v)) < 0) { return k; }
  1647.    }
  1648.    return 1;
  1649.  };
  1650.  
  1651.  #
  1652.  # Check the elements of LIST.
  1653.  # LIST must contain elements `x', `y',
  1654.  # and `z'
  1655.  #
  1656.  
  1657.  check_3d_list = function ( LIST )
  1658.  {
  1659.    #
  1660.    # Check existence and types
  1661.    #
  1662.  
  1663.    if (class (LIST) != "list") {
  1664.      error ("plot3: argument must be a list");
  1665.    }
  1666.    if (!exist (LIST.x)) {
  1667.      error ("plot3: arg must contain `x' member");
  1668.    else if (class (LIST.x) != "num") {
  1669.      error ("plot3: x must be numeric");
  1670.    }}
  1671.    if (!exist (LIST.y)) {
  1672.      error ("plot3: arg must contain `y' member");
  1673.    else if (class (LIST.y) != "num") {
  1674.      error ("plot3: y must be numeric");
  1675.    }}
  1676.    if (!exist (LIST.z)) {
  1677.      error ("plot3: arg must contain `z' member");
  1678.    else if (class (LIST.z) != "num") {
  1679.      error ("plot3: z must be numeric");
  1680.    }}
  1681.  
  1682.    #
  1683.    # Check sizes
  1684.    #
  1685.  
  1686.    if (LIST.x.n != LIST.z.nr) {
  1687.      error ("plot3: x.n != z.nr");
  1688.    }
  1689.  
  1690.    if (LIST.y.n != LIST.z.nc) {
  1691.      error ("plot3: y.n != z.nc");
  1692.    }
  1693.  
  1694.  };
  1695.  
  1696.   find_char = function ( str , char )
  1697.   {
  1698.     local (i, tmp)
  1699.  
  1700.     tmp = strsplt (str);
  1701.     for (i in 1:tmp.n)
  1702.     {
  1703.       if (tmp[i] == char) 
  1704.       {
  1705.         return i;
  1706.       }
  1707.     }
  1708.     return 0;
  1709.   };
  1710.  
  1711.  plhistx = function ( M , nbin )
  1712.  {
  1713.    local (i, j, k, l, dbin, np, p, binval, xt, yt, v, desc)
  1714.  
  1715.    check_plot_object ();
  1716.  
  1717.    if (!exist (nbin)) { nbin = 10; }
  1718.  
  1719.    p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  1720.    np = M.nr;
  1721.  
  1722.    # Compute max/min values of data
  1723.  
  1724.    ymin = min (min(real (M)));
  1725.    ymax = max (max(real (M)));
  1726.  
  1727.    #
  1728.    # Check computed scale limits against user's
  1729.    #
  1730.  
  1731.    if (!isinf (WIN.[P].ymin[p])) { ymin = WIN.[P].ymin[p]; }
  1732.    if (!isinf (WIN.[P].ymax[p])) { ymax = WIN.[P].ymax[p]; }
  1733.  
  1734.    _plgra ();
  1735.    _plcol (15);
  1736.    _plfont (WIN.[P].font[p]);
  1737.    _plwid (WIN.[P].width[p]);
  1738.  
  1739.    dbin = (linspace (ymin, ymax, nbin+1))';
  1740.    for (j in 1:M.nc) {
  1741.      // counting
  1742.      for (i in 1:nbin) {
  1743.        binval[i;j] = length (find (M[;j] >= dbin[i] && M[;j] < dbin[i+1]));
  1744.      }
  1745.    }
  1746.  
  1747.    _pladv (0);
  1748.    _plvsta ();
  1749.    xmin = 0;
  1750.    xmax =  max(max(binval));
  1751.    _plwind (ymin, ymax, xmin, xmax);
  1752.    _plbox (WIN.[P].gridx[p], 0, 0, WIN.[P].gridy[p], 0, 0);
  1753.    dbin = (linspace (ymin, ymax, nbin))';
  1754.  
  1755.    v = xmax;
  1756.    for (i in 1:M.nc)
  1757.    {
  1758.      k = mod (i, 14) + 1;
  1759.      l = mod (i,  8) + 1;
  1760.      _plcol (WIN.[P].color[p;k]);
  1761.      _pllsty (WIN.[P].lstyle[p;l]);
  1762.  
  1763.      if (get_style (WIN.[P].style.[p], k-1) == "line") 
  1764.      {
  1765.        _plline (nbin, dbin, binval[;i]);
  1766.      else if (get_style (WIN.[P].style.[p], k-1) == "point") {
  1767.        _plpoin (nbin, dbin, binval[;i], WIN.[P].pstyle[p]+k);
  1768.      else if (get_style (WIN.[P].style.[p], k-1) == "line-point") {
  1769.        _plline (nbin, dbin, binval[;i]);
  1770.        _plpoin (nbin, dbin, binval[;i], WIN.[P].pstyle[p]+k);     
  1771.      }}}
  1772.  
  1773.      // write legend around upper-right corner.
  1774.      // it is better to have user to choose location for legend.
  1775.  
  1776.      if (!any (any (isinf (WIN.[P].desc.[p]))))
  1777.      {
  1778.        # Use the default if necessary
  1779.        if (WIN.[P].desc.[p][1] == "default") 
  1780.        {
  1781.          desc = "c"+num2str(i);
  1782.        else if (WIN.[P].desc.[p].n >= i) {
  1783.          desc = WIN.[P].desc.[p][i];
  1784.        else
  1785.          # Not sure what to do, user has messed up.
  1786.          desc = "";
  1787.        }}
  1788.  
  1789.        v = v - (xmax)/11;
  1790.        xt = (ymax-ymin)*[10.5/12, 11/12, 11.5/12]' + ymin;
  1791.        yt = [v, v, v]';
  1792.  
  1793.        if (get_style (WIN.[P].style.[p], k-1) == "line") 
  1794.        {
  1795.          _plline (3, xt, yt);
  1796.        else if (get_style (WIN.[P].style.[p], k-1) == "point") {
  1797.          _plpoin (3, xt, yt, WIN.[P].pstyle[p]+k);
  1798.        else if (get_style (WIN.[P].style.[p]) == "line-point") {
  1799.          _plline (3, xt, yt);
  1800.          _plpoin (3, xt, yt, WIN.[P].pstyle[p]+k);
  1801.        }}}
  1802.  
  1803.        plptex(desc, xt[1]-(ymax-ymin)/25, yt[3], , , 1);
  1804.      }
  1805.    }
  1806.  
  1807.    _plcol (15);
  1808.    _pllab (WIN.[P].xlabel[p], WIN.[P].ylabel[p], WIN.[P].title[p]);  
  1809.    _plflush ();
  1810.    _pltext ();
  1811.  
  1812.    #
  1813.    # Increment the plot no. so that next time
  1814.    # we use the correct settings.
  1815.    #
  1816.  
  1817.    WIN.[P].subplot = WIN.[P].subplot + 1;
  1818.  
  1819.    return 1;
  1820.  };
  1821.  
  1822.  #
  1823.  # Create a legend in the current plot window
  1824.  #
  1825.  # if pobj.desc.[p] = inf()        no legend
  1826.  # if pobj.desc.[p] = "default"        default ("c1", "c2", ...)
  1827.  # if pobj.desc.[p] = "string"        use "string" as description
  1828.  #
  1829.  
  1830.  #
  1831.  # Set the current plot legend string
  1832.  #
  1833.  plegend = function ( LEGEND )
  1834.  {
  1835.    local (p)
  1836.  
  1837.    check_plot_object ();
  1838.  
  1839.    p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  1840.  
  1841.    if (!exist (LEGEND)) 
  1842.    {
  1843.      WIN.[P].desc.[p] = inf();
  1844.      return P;
  1845.    }
  1846.  
  1847.    if (class (LEGEND) == "string")
  1848.    {
  1849.      WIN.[P].desc.[p] = LEGEND;
  1850.    }
  1851.  
  1852.    return P;
  1853.  };
  1854.