c
c     Copyright (c) 1999,2000, Compaq Computer Corporation
c
c     See the file "license.terms" for information on usage and redistribution
c     of this file, and for a DISCLAIMER OF ALL WARRANTIES.
c
c     historgram plot utilities
c
      subroutine plot_histogram
      implicit none
c
c     Asks user what array to plot, then plot it!
c
c     Programmer: John Faricelli
c     Date: 25-Jun-85
c
      include 'const.h'
      include 'sysdep.cmn'
      include 'mesh.cmn'
      include 'data.cmn'
      include 'io.cmn'
      include 'bincodes.cmn'
c
c     locals
c
      logical luser
      integer answer, string_len
      character*10 buffer
      character*80 errfmt
c
      errfmt = '('' '')'
c
c     Short circuit this menu if there is only one data item to select from
c
      if( max_menu_num .eq. 1 ) then
         if( is_data_1d(1) ) then
            luser = data_nickname(1) .eq. 'ASCII'
            call plot_histo(arr2d(1,menu_to_arr(1)), xyz(1,idisti),
     1         ni, 0.0, charv,
     2         data_nickname(1), data_name(1),
     3         data_units(1), is_data_log(1), .FALSE., luser)
         else
            write(luttyo,100)
100         format(' ','Histograms are for 1D data only')
         endif
         return
      endif
c
200   call clear_screen
      write(luttyo,errfmt)      ! print any error messages
      errfmt = '('' '')'   ! clear format
      call write_data_menu('1D histogram plot')
      write(luttyo,105)
105   format(/' ','Request? ',$)
      read(luttyi,106,err=300,end=300) buffer
106   format(a)
      call upcase(buffer)
      if( string_len(buffer) .eq. 0 ) then
         answer = EXIT_MENU
      else if( buffer(1:1) .eq. 'Q' .or.  buffer(1:1) .eq. 'E' ) then
         answer = EXIT_MENU
      else
         read(buffer,*,err=300) answer
      endif
c
c    Error check:
c       1) Exit menu
c       2) Answer in range
c       3) Is requested data in binary file?
c
      if( answer .eq. EXIT_MENU ) then
         goto 1000
      else if ( answer .lt. 1 .or. answer .gt. max_menu_num ) then
         errfmt = '('' '',''Invalid input: Please try again...'')'
         call beep_terminal
      else if( .NOT. in_binary_file(answer) ) then
         errfmt = '('' '',''Item is not present in data file'')'
         call beep_terminal
      else if( is_data_1d(answer) ) then
         luser = data_nickname(answer) .eq. 'ASCII'
         call plot_histo(arr2d(1,menu_to_arr(answer)), xyz(1,idisti),
     1        ni, 0.0, charv,
     2        data_nickname(answer), data_name(answer),
     3        data_units(answer), is_data_log(answer), .FALSE., luser)
      else
         errfmt = '('' '',''Histograms are for 1D data only'')'
         call beep_terminal
      endif
      goto 200
c
c     illegal input
300   errfmt = '('' '',''Invalid input: Please try again...'')'
      call beep_terminal
      goto 200
c
1000  continue
      return
      end
c
      subroutine plot_histo(data_vector, dist_vector, ipts, 
     1                     cut_line_coord, cut_line_axis,
     1                     nickname, vlabel, vunits, vlog, lsave,
     1                     luser)
      implicit none
c
      include 'const.h'
      include 'sysdep.cmn'
      include 'mesh.cmn'
      include 'data.cmn'
      include 'io.cmn'
      include 'plot.cmn'
      include 'units.cmn'
      include 'defaults.cmn'
c
      logical lsave, vlog, luser
      integer ipts
      real data_vector(ipts), dist_vector(ipts), cut_line_coord
      character*(*) nickname, vlabel, vunits
      character*1 cut_line_axis
c
c     locals
c
      logical do_plot_title, do_oxide
      logical do_markers, fake_ln03
      character*40 vert_label, vert_units, horiz_label, horiz_units
      character*30 cut_line_string
      integer i, ilen, ulen, string_len
      integer v_exponent, v_minor_tic_freq
      real real_vec(10)
      real v_tic_major, delta
      real v_min, v_max
      real width, hwidth
c
CC BEGIN MENU
      real hmin_full, hmax_full, h_min, h_max
      real prev_hmin_plot, prev_hmax_plot,
     1     prev_h_tic_major
      integer prev_h_minor_tic_freq, prev_h_exponent, iuser_hscale
      real vmin_full, vmax_full, vmin_plot, vmax_plot
      real prev_vmin_plot, prev_vmax_plot,
     1     prev_v_tic_major
      integer item_number, int_vec(2),
     1     prev_v_minor_tic_freq
      integer iuser_vscale, prev_iuser_vscale,
     1        prev_v_exponent
      character*80 string
      character*6 str_hlog, str_vlog
      character*3 str_habs, str_vabs, markers
      character*1 prev_cut_line_axis
      real h
      integer ierr
      logical eof, do_vlog, do_vabs, do_hlog, do_habs,
     1        do_vuscale, prev_vabs, abort, is_vtxxx, lopen,
     2        delayed_message, prev_v_saved, prev_h_saved, lannotate,
     3        do_huscale, prev_habs
CC END MENU
CC ZOOM
      real h1, h2, h3, h4, x1, x2, y1, y2
      logical do_zoom, do_unzoom
CC END ZOOM
c
      real rdum1, rdum2, rdum3
      integer idum, idum2
      character*1 cdum
c
      data prev_cut_line_axis/' '/
      save prev_cut_line_axis
c
      lopen = .FALSE.
      lannotate = .FALSE.
      delayed_message = .FALSE.
      do_plot_title = .TRUE.
      cdum = ' '
      do_huscale = .FALSE.
c
c     Form cut line string for second line of header...
      if( luser ) then
         cut_line_string = ' '
      else
         write(cut_line_string,160) cut_line_axis,
     1                              cut_line_coord*pm_scale_cm2plot
160      format('at ',a1,' = ',1pg12.5)
      endif
c
      do_oxide = .FALSE.
c
      do 100 i=1,ipts
         work_vec(i,1) = data_vector(i)
         work_vec(i,2) = dist_vector(i)*pm_scale_cm2plot
100   continue
c
      if( lsave ) then
         horiz_units = oned_h_units
         horiz_label = oned_h_label
      elseif ( luser ) then
         horiz_units = ' '
         horiz_label = ' '
      else
         ulen = string_len(pm_length_unit)
         horiz_units = pm_length_unit(1:ulen)
         if( cut_line_axis .eq. charv ) then
           if( plane_type .eq. XY_PLANE ) then
c             Special case for SUPREM3, which flips X/Y axes
               if( charv .eq. 'X' ) then
                  horiz_label = 'Depth'
               else
                  horiz_label = 'Length'
               endif
            else if( plane_type .eq. YZ_PLANE ) then
               horiz_label = 'Width'
            else if( plane_type .eq. XZ_PLANE ) then
               horiz_label = 'Length'
            else if( plane_type .eq. XT_PLANE ) then
               horiz_label = 'Time'
               horiz_units = 's'
            else if( plane_type .eq. YT_PLANE ) then
               horiz_label = 'Time'
               horiz_units = 's'
            endif
         else
            if( plane_type .eq. XY_PLANE ) then
               horiz_label = 'Depth'
            else if( plane_type .eq. YZ_PLANE ) then
               horiz_label = 'Depth'
            else if( plane_type .eq. XZ_PLANE ) then
               horiz_label = 'Width'
            else if( plane_type .eq. XT_PLANE ) then
               horiz_label = 'Length'
            else if( plane_type .eq. YT_PLANE ) then
               horiz_label = 'Depth'
            endif
         endif
      endif
c
      if( lsave ) then
         do_habs = oned_h_abs
         do_hlog = oned_h_log
c        Get min/max of X axis data only
         iuser_hscale = 0.0
         call recompute_minmax(work_vec(1,2), ipts,
     1                         hmin_full, hmax_full,
     2                         rdum1, rdum2,
     3                         rdum3, idum,
     4                         iuser_hscale, idum2,
     5                         do_habs, do_hlog, do_huscale,
     6                         work_vec(1,3))
         if( idum2 .ne. oned_h_exponent ) then
c
c           We'd have code here if we supported user scaling of x axis
c
         endif
      else if( oned_lfirst .or. l_autoscale .or. luser .or.
     1        (prev_cut_line_axis .ne. cut_line_axis) ) then
c        Actually, we should reset bounds every time we read in a new binary
c        file.
         oned_lfirst = .FALSE.
         do_habs  = .FALSE.
         do_hlog  = .FALSE.
         call recompute_minmax(work_vec(1,2), ipts,
     1                         hmin_full, hmax_full,
     2                         oned_h_min, oned_h_max,
     3                         oned_h_tic_major, oned_h_minor_tic_freq,
     4                         iuser_hscale, oned_h_exponent,
     5                         do_habs, do_hlog, do_huscale,
     6                         work_vec(1,3))
c
c        Assumes uniform spacing of histogram
         width = work_vec(2,2) - work_vec(1,2)
         hwidth = width*10.0**-oned_h_exponent
         oned_h_min = oned_h_min - hwidth
         oned_h_max = oned_h_max + hwidth
c        Get updated tic mark info
         call get_default_limits(oned_h_min, oned_h_max,
     1              rdum1, rdum2,
     2              oned_h_tic_major, oned_h_minor_tic_freq,
     3              idum, .FALSE.)
c        Unscale absissa
CC why? jvf
CC         if( oned_h_exponent .ne. 0 ) then
CC            oned_h_tic_major = oned_h_tic_major*10.0**oned_h_exponent
CC            oned_h_exponent = 0
CC         endif
      else
c        Reset hmin and hmax, in case bounds changed
         iuser_hscale = 0.0
         call recompute_minmax(work_vec(1,2), ipts,
     1                         hmin_full, hmax_full,
     2                         rdum1, rdum2,
     3                         rdum3, idum,
     4                         iuser_hscale, idum2,
     5                         .FALSE., .FALSE., .FALSE.,
     6                         work_vec(1,3))
      endif
c
c     Compute defaults for "prev" quantities
      if( do_hlog ) then
         prev_h_saved = .FALSE.
      else
         prev_h_saved = .TRUE.
         prev_hmin_plot = oned_h_min
         prev_hmax_plot = oned_h_max
         prev_h_tic_major = oned_h_tic_major
         prev_h_minor_tic_freq = oned_h_minor_tic_freq
         prev_h_exponent = oned_h_exponent
      endif
c
      if( lsave ) then
         do_vlog = oned_v_log
         do_vabs = oned_v_abs
         do_markers = oned_do_markers
         vmin_plot = oned_v_min
         vmax_plot = oned_v_max
         v_tic_major = oned_v_tic_major
         v_minor_tic_freq = oned_v_minor_tic_freq
         iuser_vscale = 0
         do_vuscale = .FALSE.
         call recompute_minmax(work_vec(1,1), ipts,
     1                         vmin_full, vmax_full,
     2                         rdum1, rdum2,
     3                         rdum3, idum,
     4                         iuser_vscale, idum2,
     5                         do_vabs, do_vlog, do_vuscale,
     6                         work_vec(1,3))
         if( idum2 .ne. oned_v_exponent ) then
            do_vuscale = .TRUE.
            iuser_vscale = oned_v_exponent
            h = 10.0**oned_v_exponent/10.0**idum2
         endif
         v_exponent = oned_v_exponent
      else
         if( vlog ) then
            do_vlog  = .TRUE.
            do_vabs  = .TRUE.
         else
            do_vlog  = .FALSE.
            do_vabs  = .FALSE.
         endif
         prev_vabs  = .FALSE.
         do_markers = .FALSE.
         v_exponent = 0
         iuser_vscale = 0
         do_vuscale = .FALSE.
c
         call recompute_minmax(work_vec(1,1), ipts,
     1                         vmin_full, vmax_full,
     2                         vmin_plot, vmax_plot,
     3                         v_tic_major, v_minor_tic_freq,
     4                         iuser_vscale, v_exponent,
     5                         do_vabs, do_vlog, do_vuscale,
     6                         work_vec(1,3))
c
c        For histograms, it's clearer if the Y limit is expanded to zero
         if( .NOT. vlog ) then
            if( vmin_plot .gt. 0.0 .or. vmax_plot .lt. 0.0 ) then
               if( vmin_plot .gt. 0.0 ) vmin_plot = 0.0
               if( vmax_plot .lt. 0.0 ) vmax_plot = 0.0
c              Get updated tic mark info
               call get_default_limits(vmin_plot, vmax_plot,
     1              rdum1, rdum2,
     2              v_tic_major, v_minor_tic_freq,
     3              idum, .FALSE.)
            endif
         endif
      endif
c
      if( do_vlog ) then
         prev_v_saved = .FALSE.
      else
         prev_v_saved = .TRUE.
         prev_iuser_vscale = iuser_vscale
         prev_vmin_plot = vmin_plot
         prev_vmax_plot = vmax_plot
         prev_v_tic_major   = v_tic_major
         prev_v_minor_tic_freq  = v_minor_tic_freq
         prev_v_exponent = v_exponent
      endif
      prev_cut_line_axis = cut_line_axis
c
      vert_label = vlabel
      vert_units = vunits
c
      if(do_hlog) then
         str_hlog = 'log'
      else
         str_hlog = 'linear'
      endif
      if(do_habs) then
         str_habs = 'yes'
      else
         str_habs = 'no'
      endif
c
      if(do_vlog) then
         str_vlog = 'log'
      else
         str_vlog = 'linear'
      endif
      if(do_vabs) then
         str_vabs = 'yes'
      else
         str_vabs = 'no'
      endif
c
      if(do_markers) then
         markers = 'yes'
      else
         markers = 'no'
      endif
c
c     If batch processing, perform hardcopy and return
      if( lbatch .or. lrestore ) then
         previous_plot_device = plot_device
         if( lbatch ) then
            plot_device = hardcopy_device
            if( lplotfile ) then
               plot_file_name = given_plotfilename
            endif
         endif
c
c        Stuff axis information into save structures
         call copy_axis_to_save(PLOT_1D_COMMAND,
     1        oned_h_min, oned_h_max,
     2        oned_h_exponent, oned_h_tic_major, oned_h_minor_tic_freq,
     3        horiz_units,
     4        horiz_label, do_hlog, do_habs,
     5        vmin_plot, vmax_plot,
     6        v_exponent, v_tic_major, v_minor_tic_freq, vert_units,
     7        vert_label, do_vlog, do_vabs,
     8        rdum1, rdum1,
     9        idum, rdum1, idum, cdum,
     A        cdum, .FALSE., .FALSE.)
c
c        Stuff data information into save structures
         call copy_1d_to_save(cut_line_axis,
     1        cut_line_coord,
     2        nickname, do_markers)
c        Copy data into temp store, take absolute value, log if necessary...
         do 120 i=1,ipts
            if( do_vabs .or. do_vlog ) then
               work_vec(i,3) = ABS(data_vector(i))
               if( do_vlog ) then
c                 check for log(0.0), fudge
                  if( work_vec(i,3) .eq. 0.0 ) work_vec(i,3) =
     1                                         PM_MIN_FLOAT
                  work_vec(i,3) = ALOG10(work_vec(i,3))
               endif
            else
               work_vec(i,3) = data_vector(i)
            endif
120      continue
         if( do_vlog ) then
            v_min = ALOG10(vmin_plot)
            v_max = ALOG10(vmax_plot)
         else
            v_min = vmin_plot
            v_max = vmax_plot
         endif
c
         if( do_hlog ) then
            h_min = ALOG10(oned_h_min)
            h_max = ALOG10(oned_h_max)
         else
            h_min = oned_h_min
            h_max = oned_h_max
         endif
c
         call do_histo_plot(work_vec(1,3), work_vec(1,2), ipts,
     1        h_min, h_max,
     2        oned_h_exponent, oned_h_tic_major, oned_h_minor_tic_freq,
     3        v_min, v_max,
     4        v_exponent, v_tic_major, v_minor_tic_freq,
     5        horiz_label, horiz_units, vert_label, vert_units,
     6        do_vlog, do_vabs, do_hlog, do_habs,
     7        cut_line_string, do_oxide, do_plot_title, do_markers,
     8        lsave, lopen, lannotate, ierr)
         plot_device = previous_plot_device
         return
      endif
c
c     Main return point for menu
c
4000  continue
c     Write out menu form:
      call clear_screen
      if( .NOT. luser ) then
         call display_string(1,1 ,'Cut line selected:     =',24,0)
      endif
      call display_string(3,1 ,'1)  Title:',10, 0)
c
      call display_string(5 ,1 ,'HORIZONTAL AXIS PARAMETERS:',27,0)
      call display_string(6 ,1 ,'Data range:',11,0)
      call display_string(6 ,24,' to ',4,0)
      call display_string(8 ,1 ,'2)  Plot range:',15,0)
      call display_string(8 ,24,' to ',4,0)
      call display_string(9 ,1 ,'3)  Major tics at:',18,0)
      call display_string(10,1 ,'4)  Minor tic freq:',19,0)
      call display_string(11,1 ,'5)  Scale factor:       10**',28,0)
      call display_string(6 ,43,'6)  Linear/log: ',16, 0)
      call display_string(7 ,43,'7)  Absolute value: ',20, 0)
      call display_string(8 ,43,'8)  Recompute default axis parameters',
     1                    37,0)
      call display_string(9 ,43,'9)  Label: ',11, 0)
      call display_string(10,43,'10) Units: ',11, 0)
c
      call display_string(13,1 ,'VERTICAL AXIS PARAMETERS:',25,0)
      call display_string(14,1 ,'Data range:',11,0)
      call display_string(14,24,' to',3,0)
      call display_string(16,1 ,'11) Plot range:',15,0)
      call display_string(16,24,' to',3,0)
      call display_string(17,1 ,'12) Major tics at:',18,0)
      call display_string(18,1 ,'13) Minor tic freq:',19,0)
      call display_string(19,1 ,'14) Scale factor:       10**',28,0)
      call display_string(14,43,'15) Linear/log: ',16, 0)
      call display_string(15,43,'16) Absolute value: ',20, 0)
      call display_string(16,43,'17) Recompute default axis parameters',
     1                    37,0)
      call display_string(17,43,'18) Label: ',11, 0)
      call display_string(18,43,'19) Units: ',11, 0)
      call display_string(19,43,'20) Symbols: ',13, 0)
c
      call display_string(23,1 ,'Messages:',9,0)
c
c     Cut line
      if( .NOT. luser ) then
         call display_string(1 , 22, cut_line_axis, 1, 0)
         call display_real(1 , 26, cut_line_coord*pm_scale_cm2plot, 0)
      endif
c     Title
      call display_string(3 , 11, plot_title, 69, 1)
c
c     X Range
      call display_real  (6 , 15, hmin_full, 1)
      call display_real  (6 , 28, hmax_full, 1)
      call display_real  (8 , 15, oned_h_min, 1)
      call display_real  (8 , 28, oned_h_max, 1)
c     X major tics
      call display_real  (9 , 23, oned_h_tic_major, 1)
c     X minor tic freq
      call display_int   (10, 27, oned_h_minor_tic_freq, 5)
c     X scale factor
      call display_int  (11, 29, oned_h_exponent, 3)
c     Linear/log
      call display_string(6 , 59, str_hlog, 6, 1)
c     Absolute value
      call display_string(7 , 63, str_habs, 3, 1)
c     Horizontal label
      call display_string(9 , 54, horiz_label, 27, 1)
c     Horizontal Units
      call display_string(10, 54, horiz_units, 20, 1)
c
c     Y Range
      call display_real  (14, 15, vmin_full, 1)
      call display_real  (14, 28, vmax_full, 1)
      call display_real  (16, 15, vmin_plot, 1)
      call display_real  (16, 28, vmax_plot, 1)
      if( do_vlog ) then
         call display_string(17, 23, '       NA', 9, 1)
         call display_string(18, 27, '   NA', 5, 1)
      else
c        Y major tics
         call display_real  (17, 23, v_tic_major, 1)
c        Y minor tic freq
         call display_int   (18, 27, v_minor_tic_freq, 5)
      endif
c     Y scale factor
      call display_int  (19, 29, v_exponent, 3)
c     Linear/log
      call display_string(14, 59, str_vlog, 6, 1)
c     Absolute value
      call display_string(15, 63, str_vabs, 3, 1)
c     Vertical label
      call display_string(17, 54, vert_label, 27, 1)
c     Vertical Units
      call display_string(18, 54, vert_units, 20, 1)
      call display_string(19, 56, markers, 3, 1)
c
c     Display any messaged from before window refresh
      if( delayed_message ) then
         call display_message(string)
         delayed_message = .FALSE.
      endif
c
2000  call display_prompt(
     1 'Enter number of item, or '//
     2 'Annotate Plot Hardcopy Keepwin Repaint Zoom Unzoom Exit')
      call read_from_input_area(string, eof)
      if( eof ) goto 2000
      call upcase(string)
      if( string(1:1) .eq. 'P' ) then
         call display_prompt('Making plot')
c
c        Stuff axis information into save structures
         call copy_axis_to_save(PLOT_1D_COMMAND,
     1        oned_h_min, oned_h_max,
     2        oned_h_exponent, oned_h_tic_major, oned_h_minor_tic_freq,
     3        horiz_units,
     4        horiz_label, do_hlog, do_habs,
     5        vmin_plot, vmax_plot,
     6        v_exponent, v_tic_major, v_minor_tic_freq, vert_units,
     7        vert_label, do_vlog, do_vabs,
     8        rdum1, rdum1,
     9        idum, rdum1, idum, cdum,
     A        cdum, .FALSE., .FALSE.)
c
c        Stuff data information into save structures
         call copy_1d_to_save(cut_line_axis,
     1        cut_line_coord,
     2        nickname, do_markers)
c        Copy data into temp store, take absolute value, log if necessary...
         do 165 i=1,ipts
            if( do_vabs .or. do_vlog ) then
               work_vec(i,3) = ABS(data_vector(i))
               if( do_vlog ) then
c                 check for log(0.0), fudge
                  if( work_vec(i,3) .eq. 0.0 ) work_vec(i,3) =
     1                                         PM_MIN_FLOAT
                  work_vec(i,3) = ALOG10(work_vec(i,3))
               endif
            else
               work_vec(i,3) = data_vector(i)
            endif
165      continue
         if( do_vlog ) then
            v_min = ALOG10(vmin_plot)
            v_max = ALOG10(vmax_plot)
         else
            v_min = vmin_plot
            v_max = vmax_plot
         endif
c
         if( do_hlog ) then
            h_min = ALOG10(oned_h_min)
            h_max = ALOG10(oned_h_max)
         else
            h_min = oned_h_min
            h_max = oned_h_max
         endif
c
         call do_histo_plot(work_vec(1,3), work_vec(1,2), ipts,
     1        h_min, h_max,
     2        oned_h_exponent, oned_h_tic_major, oned_h_minor_tic_freq,
     3        v_min, v_max,
     4        v_exponent, v_tic_major, v_minor_tic_freq,
     5        horiz_label, horiz_units, vert_label, vert_units,
     6        do_vlog, do_vabs, do_hlog, do_habs,
     7        cut_line_string, do_oxide, do_plot_title, do_markers,
     8        lsave, lopen, lannotate, ierr)
         if( ierr .eq. PLT$NORMAL ) then
            lopen = .TRUE.
         else
            call plt$get_error_message(ierr,string)
            delayed_message = .TRUE.
         endif
c
         goto 4000
c
      else if( string(1:1) .eq. 'Z' .or. string(1:1) .eq. 'U') then
         if( lopen ) then
            if( is_vtxxx(current_window_id) ) then
               call display_message(
     1           'Sorry: cannot zoom on VTxxx terminals')
               goto 2000
            else
               if( string(1:1) .eq. 'Z' ) then
                  call display_prompt(
     1     'Enter lower left, upper right. ' // PM_CANCEL_BUTTON //
     1     ' button cancels input')
                  do_unzoom = .FALSE.
                  if( do_hlog ) then
                     h1 = ALOG10(ABS(MAX(PM_MIN_FLOAT,oned_h_min)))
                     h2 = ALOG10(ABS(MAX(PM_MIN_FLOAT,oned_h_max)))
                  else
                     h1 = oned_h_min
                     h2 = oned_h_max
                  endif
                  if( do_vlog ) then
                     h3 = ALOG10(ABS(MAX(PM_MIN_FLOAT,vmin_plot)))
                     h4 = ALOG10(ABS(MAX(PM_MIN_FLOAT,vmax_plot)))
                  else
                     h3 = vmin_plot
                     h4 = vmax_plot
                  endif
c                 Magic constants are values to GKS$SET_VIEWPORT in
c                 subroutine do_histo_plot.
                  call window_zoom(current_window_id,
     1               h1, h2, h3, h4,
     2               0.175, 0.825, 0.10, 0.80, .FALSE.,
     3               x1, x2, y1, y2, do_zoom)
                  if( do_zoom ) then
                     if( do_hlog ) then
                        x1 = 10.0**(MIN(x1,38.0))
                        x2 = 10.0**(MIN(x2,38.0))
                     endif
                     if( do_vlog ) then
                        y1 = 10.0**(MIN(y1,38.0))
                        y2 = 10.0**(MIN(y2,38.0))
                     endif
                  endif
               else
                  x1 = hmin_full
                  x2 = hmax_full
                  y1 = vmin_full
                  y2 = vmax_full
                  do_unzoom = .TRUE.
               endif
               if( do_zoom .or. do_unzoom ) then
                  oned_h_min  = x1
                  oned_h_max  = x2
                  vmin_plot  = y1
                  vmax_plot  = y2
c
c                 Get tic mark parameters
                  if( .NOT. do_hlog ) then
                     call get_default_limits(oned_h_min, oned_h_max,
     1                    rdum1, rdum2,
     2                    oned_h_tic_major, oned_h_minor_tic_freq,
     3                    idum, .FALSE.)
c                    Handle the case of unzoom
                     if( idum .ne. 0 ) then
                        oned_h_min = oned_h_min * 10.0**(-idum)
                        oned_h_max = oned_h_max * 10.0**(-idum)
CC                        oned_h_tic_major=oned_h_tic_major*10.0**(idum)
                        oned_h_exponent = idum
                     endif
                  endif
                  if( do_unzoom ) then
c                    Recompute all vertical axis parameters
                     call recompute_minmax(work_vec(1,1), ipts,
     1                  rdum1, rdum2,
     2                  vmin_plot, vmax_plot,
     3                  v_tic_major, v_minor_tic_freq,
     4                  iuser_vscale, v_exponent,
     5                  do_vabs, do_vlog, do_vuscale, work_vec(1,3))
                  else if( .NOT. do_vlog ) then
c                    Just get updated tic mark info
                     call get_default_limits(vmin_plot, vmax_plot,
     1                  rdum1, rdum2,
     2                  v_tic_major, v_minor_tic_freq,
     3                  idum, .FALSE.)
                  endif
c
c                 Stuff axis information into save structures
                  call copy_axis_to_save(PLOT_1D_COMMAND,
     1               oned_h_min, oned_h_max,
     2               oned_h_exponent, oned_h_tic_major,
     3               oned_h_minor_tic_freq, horiz_units,
     4               horiz_label, do_hlog, do_habs,
     5               vmin_plot, vmax_plot,
     6               v_exponent, v_tic_major,
     7               v_minor_tic_freq, vert_units,
     8               vert_label, do_vlog, do_vabs,
     9               rdum1, rdum2,
     A               idum, rdum3, idum2, cdum,
     B               cdum, .FALSE., .FALSE.)
c
c                 Stuff data information into save structures
                  call copy_1d_to_save(cut_line_axis,
     1                 cut_line_coord, nickname, do_markers)
c
c                 Copy data into temp store, take absolute value,
c                 log if necessary...
                  do 140 i=1,ipts
                     if( do_vabs .or. do_vlog ) then
                        work_vec(i,3) = ABS(data_vector(i))
                        if( do_vlog ) then
c                          check for log(0.0), fudge
                           if( work_vec(i,3) .eq. 0.0 )
     1                         work_vec(i,3) = PM_MIN_FLOAT
                           work_vec(i,3) = ALOG10(work_vec(i,3))
                        endif
                     else
                        work_vec(i,3) = data_vector(i)
                     endif
140               continue
                  if( do_vlog ) then
                     v_min = ALOG10(vmin_plot)
                     v_max = ALOG10(vmax_plot)
                  else
                     v_min = vmin_plot
                     v_max = vmax_plot
                  endif
c
                 if( do_hlog ) then
                    h_min = ALOG10(oned_h_min)
                    h_max = ALOG10(oned_h_max)
                 else
                    h_min = oned_h_min
                    h_max = oned_h_max
                 endif
c
                  call do_histo_plot(work_vec(1,3), work_vec(1,2),
     1               ipts, h_min, h_max,
     2               oned_h_exponent, oned_h_tic_major,
     3               oned_h_minor_tic_freq,
     4               v_min, v_max,
     5               v_exponent, v_tic_major, v_minor_tic_freq,
     6               horiz_label, horiz_units, vert_label, vert_units,
     7               do_vlog, do_vabs, do_hlog, do_habs,
     8               cut_line_string,
     9               do_oxide, do_plot_title, do_markers,
     A               lsave, lopen, lannotate, ierr)
                  endif
               endif
            goto 4000
         else
            call display_message(
     1         'Plot must be displayed before you can zoom it')
         endif
         goto 2000
c
      else if( string(1:1) .eq. 'H' ) then
         fake_ln03 = .TRUE.
         call inquire_hardcopy_menu(fake_ln03, abort)
         if( abort ) goto 4000
         previous_plot_device = plot_device
         plot_device = hardcopy_device
         call display_prompt('Making hardcopy')
c
c        Stuff axis information into save structures
         call copy_axis_to_save(PLOT_1D_COMMAND,
     1        oned_h_min, oned_h_max,
     2        oned_h_exponent, oned_h_tic_major, oned_h_minor_tic_freq,
     3        horiz_units,
     4        horiz_label, do_hlog, do_habs,
     5        vmin_plot, vmax_plot,
     6        v_exponent, v_tic_major, v_minor_tic_freq, vert_units,
     7        vert_label, do_vlog, do_vabs,
     8        rdum1, rdum1,
     9        idum, rdum1, idum, cdum,
     A        cdum, .FALSE., .FALSE.)
c
c        Stuff data information into save structures
         call copy_1d_to_save(cut_line_axis,
     1        cut_line_coord,
     2        nickname, do_markers)
c        Copy data into temp store, take absolute value, log if necessary...
         do 175 i=1,ipts
            if( do_vabs .or. do_vlog ) then
               work_vec(i,3) = ABS(data_vector(i))
               if( do_vlog ) then
c                 check for log(0.0), fudge
                  if( work_vec(i,3) .eq. 0.0 ) work_vec(i,3) =
     1                                         PM_MIN_FLOAT
                  work_vec(i,3) = ALOG10(work_vec(i,3))
               endif
            else
               work_vec(i,3) = data_vector(i)
            endif
175      continue
         if( do_vlog ) then
            v_min = ALOG10(vmin_plot)
            v_max = ALOG10(vmax_plot)
         else
            v_min = vmin_plot
            v_max = vmax_plot
         endif
c
         if( do_hlog ) then
            h_min = ALOG10(oned_h_min)
            h_max = ALOG10(oned_h_max)
         else
            h_min = oned_h_min
            h_max = oned_h_max
         endif
c
         call do_histo_plot(work_vec(1,3), work_vec(1,2), ipts,
     1        h_min, h_max,
     2        oned_h_exponent, oned_h_tic_major, oned_h_minor_tic_freq,
     3        v_min, v_max,
     4        v_exponent, v_tic_major, v_minor_tic_freq,
     5        horiz_label, horiz_units, vert_label, vert_units,
     6        do_vlog, do_vabs, do_hlog, do_habs,
     7        cut_line_string, do_oxide, do_plot_title, do_markers,
     8        lsave, lopen, lannotate, ierr)
         plot_device = previous_plot_device
         if( ierr .eq. PLT$NORMAL ) then
            ilen = string_len(plot_file_name)
            write(string,200,err=201) plot_file_name(1:ilen)
200         format('Plot file ',a' created')
201         continue
         else
            call plt$get_error_message(ierr,string)
         endif
         delayed_message = .TRUE.
         goto 4000
c
      else if( string(1:1) .eq. 'A' ) then
         if( lopen ) then
            if( is_vtxxx(current_window_id) ) then
               call display_message(
     1           'Sorry: cannot annotate on VTxxx terminals')
            else
               call anno_annotate_graph(current_window_id)
               lannotate = .TRUE.
               goto 4000
            endif
         else
            call display_message(
     1         'Plot must be displayed before you can annotate it')
         endif
         goto 2000
c     Special case to catch E or Q (which are "legal" floating point
c     numbers!)
      else if( string(1:1) .eq. 'E' .or. string(1:1) .eq. 'Q') then
         goto 5000
      else if( string(1:1) .eq. 'K') then
         if( lopen ) then
            if( is_vtxxx(current_window_id) ) then
               call display_message(
     1           'Sorry: no window functions on VTxxx terminals')
               goto 2000
            else
c              Set "save this window" flag
c              (only for workstations and when using retained window style)
               call retain_current_window
            endif
         else
            call display_message(
     1      'Plot must be displayed in a window before you can keep it')
         endif
         goto 2000
      else if( string(1:1) .eq. 'R') then ! Repaint screen
         goto 4000
      else
         read(string,*,err=1000) item_number
         if( item_number .eq. 1 ) then
            call get_string_with_prompt(plot_title, 
     1              'Enter title (may use mathmode)')
            call display_string( 3, 12, plot_title, 69, 1)
c
         elseif ( item_number .eq. 2 ) then
c           X axis min/max
            real_vec(1) = oned_h_min
            real_vec(2) = oned_h_max
            call get_reals_with_prompt(real_vec, 2,
     1              'Enter min and max values for horizontal axis')
            if( real_vec(1) .ge. real_vec(2) ) then
               call display_message(
     1                 'Min value must be less than max')
            else if( do_hlog .and.
     1         (real_vec(1) .le. 0.0 .or. real_vec(2) .le. 0.0) ) then
               call display_message(
     1                 'Log limits must be greater than zero')
            else
               oned_h_min = real_vec(1)
               call display_real(8 , 15, oned_h_min, 1)
               oned_h_max = real_vec(2)
               call display_real(8 , 28, oned_h_max, 1)
c              Generate compatible tic mark parameters
               if( .NOT. do_hlog ) then
                  call get_default_limits(oned_h_min, oned_h_max,
     1               rdum1, rdum2,
     2               oned_h_tic_major, oned_h_minor_tic_freq,
     3               idum, .FALSE.)
c                 Unscale
                  oned_h_tic_major = oned_h_tic_major*10.0**idum
                  call display_real(9 , 23, oned_h_tic_major, 1)
                  call display_int (10, 27, oned_h_minor_tic_freq,5)
               endif
            endif
c
         elseif ( item_number .eq. 3 ) then
c           X axis major tic increment
            real_vec(1) = oned_h_tic_major
            if( do_hlog ) then
               call display_message(
     1            'Not used when plotting on log scale')
            else
               call get_reals_with_prompt(real_vec, 1,
     1                 'Enter major tic increment for horizontal axis')
               delta = oned_h_max - oned_h_min
               if( real_vec(1) .gt. delta .or.
     1            real_vec(1) .le. 0.0 ) then
                  call display_message(
     1                   'Increment < 0 or larger than plot range')
               else
                  oned_h_tic_major = real_vec(1)
                  call display_real  (9 , 23, oned_h_tic_major, 1)
               endif
            endif
c
         elseif ( item_number .eq. 4 ) then
c           Minor tic frequency
            int_vec(1) = oned_h_minor_tic_freq
            if( do_hlog ) then
               call display_message(
     1            'Not used when plotting on log scale')
            else
               call get_ints_with_prompt(int_vec, 1,
     1                 'Enter minor tic frequency for horizontal axis')
               if( int_vec(1) .lt. 1 ) then
                  call display_message('Illegal value')
               else
                  oned_h_minor_tic_freq = int_vec(1)
                  call display_int(10, 27, oned_h_minor_tic_freq,5)
               endif
            endif
c
         elseif ( item_number .eq. 5 ) then
c           X axis scale factor
            call display_message('Sorry: not implemented yet')
c
         elseif ( item_number .eq. 6 ) then
c           X axis log
            do_hlog = .NOT. do_hlog
            if( do_hlog ) then
               prev_habs = do_habs
               do_habs = .TRUE.
               str_hlog = 'log'
               str_habs = 'yes'
c              Save previous values
               prev_h_exponent = oned_h_exponent
               prev_hmin_plot  = oned_h_min
               prev_hmax_plot  = oned_h_max
               prev_h_tic_major  = oned_h_tic_major
               prev_h_minor_tic_freq = oned_h_minor_tic_freq
               prev_h_saved = .TRUE.
               call recompute_minmax(work_vec(1,2), ipts,
     1            hmin_full, hmax_full,
     2            oned_h_min, oned_h_max,
     3            oned_h_tic_major, oned_h_minor_tic_freq,
     4            iuser_hscale, oned_h_exponent,
     5            do_habs, do_hlog, do_huscale, work_vec(1,3))
               oned_h_exponent = 0
c              Mark X axis tics as not applicable
               call display_string(9 , 23, '       NA', 9, 1)
               call display_string(10, 27, '   NA', 5, 1)
            else
               do_habs = prev_habs
               str_hlog = 'linear'
               if( do_habs ) then
                  str_habs = 'yes'
               else
                  str_habs = 'no'
               endif
c              Recompute hmin_full/hmax_full
               call recompute_minmax(work_vec(1,2), ipts,
     1            hmin_full, hmax_full,
     2            oned_h_min, oned_h_max,
     3            oned_h_tic_major, oned_h_minor_tic_freq,
     4            iuser_hscale, oned_h_exponent,
     5            do_habs, do_hlog, do_huscale, work_vec(1,3))
c              Restore previous axis values
               if( prev_h_saved ) then
                  oned_h_exponent = prev_h_exponent
                  oned_h_min = prev_hmin_plot
                  oned_h_max = prev_hmax_plot
                  oned_h_tic_major = prev_h_tic_major
                  oned_h_minor_tic_freq = prev_h_minor_tic_freq
               endif
c              Re-display tic/scale data
               call display_real  (9 , 23, oned_h_tic_major, 1)
               call display_int   (10, 27, oned_h_minor_tic_freq,5)
            endif
            call display_real  (6 , 15, hmin_full, 1)
            call display_real  (6 , 28, hmax_full, 1)
            call display_real  (8 , 15, oned_h_min, 1)
            call display_real  (8 , 28, oned_h_max, 1)
            call display_int   (11, 29, oned_h_exponent,3)
            call display_string(6 , 59, str_hlog, 6, 1)
            call display_string(7 , 63, str_habs, 3, 1)
c
         elseif ( item_number .eq. 7 ) then
c           X axis absolute value (not implented)
            call display_message('Sorry: not implemented yet')
c
         elseif ( item_number .eq. 8 ) then
c           X axis recompute defaults (not implented)
            call display_message('Sorry: not implemented yet')
c
         elseif ( item_number .eq. 9 ) then
            call get_string_with_prompt(horiz_label,
     1              'Enter horizontal label (may use mathmode)')
            call display_string( 9, 54, horiz_label, 27, 1)
c
         elseif ( item_number .eq. 10 ) then
            call get_string_with_prompt(horiz_units,
     1              'Enter units (may use mathmode)')
            call display_string(10, 54, horiz_units, 20, 1)
c
         elseif ( item_number .eq. 11 ) then
c           Y axis min/max
            real_vec(1) = vmin_plot
            real_vec(2) = vmax_plot
            call get_reals_with_prompt(real_vec, 2,
     1              'Enter min and max values for vertical axis')
            if( real_vec(1) .ge. real_vec(2) ) then
               call display_message(
     1                 'Min value must be less than max')
            else if( do_vlog .and.
     1         (real_vec(1) .le. 0.0 .or. real_vec(2) .le. 0.0) ) then
               call display_message(
     1                 'Log limits must be greater than zero')
            else
               vmin_plot = real_vec(1)
               call display_real(16, 15, vmin_plot, 1)
               vmax_plot = real_vec(2)
               call display_real(16, 28, vmax_plot, 1)
c              Generate compatible tic mark parameters
               if( .NOT. do_vlog ) then
                  call get_default_limits(vmin_plot, vmax_plot,
     1               rdum1, rdum2,
     2               v_tic_major, v_minor_tic_freq,
     3               idum, .FALSE.)
c                 Unscale absissa
                  v_tic_major = v_tic_major*10.0**idum
                  call display_real(17, 23, v_tic_major, 1)
                  call display_int (18, 27, v_minor_tic_freq,5)
               endif
            endif
c
         elseif ( item_number .eq. 12 ) then
c           Major Y tic increment
            real_vec(1) = v_tic_major
            if( do_vlog ) then
               call display_message(
     1            'Not used when plotting on log scale')
            else
               call get_reals_with_prompt(real_vec, 1,
     1              'Enter major tic increment for vertical axis')
               delta = vmax_plot - vmin_plot
               if( real_vec(1) .gt. delta .or.
     1             real_vec(1) .le. 0.0 ) then
                   call display_message('Illegal value')
               else
                  v_tic_major = real_vec(1)
                  call display_real  (17, 23, v_tic_major, 1)
               endif
            endif
c
         elseif ( item_number .eq. 13 ) then
c           Minor tic frequency
            int_vec(1) = v_minor_tic_freq
            if( do_vlog ) then
               call display_message(
     1            'Not used when plotting on log scale')
            else
               call get_ints_with_prompt(int_vec, 1,
     1              'Enter minor tic frequency for vertical axis')
               if( int_vec(1) .lt. 1 ) then
                 call display_message('Illegal value')
               else
                  v_minor_tic_freq = int_vec(1)
                  call display_int(18, 27, v_minor_tic_freq,5)
               endif
            endif
c
         elseif ( item_number .eq. 14 ) then
c           Y axis scale factor
            if( do_vlog ) then
               call display_message(
     1            'Cannot use scale factor on log plot')
            else
               int_vec(1) = v_exponent
               call get_ints_with_prompt(int_vec, 1,
     1      'Enter integer power of ten that will scale (divide) data')
               if( int_vec(1) .gt. 38 .or.
     1             int_vec(1) .lt. -38 ) then
                  call display_message('Illegal value')
               else
                  iuser_vscale = int_vec(1)
                  do_vuscale = .TRUE.
                  h = 10.0**iuser_vscale/10.0**v_exponent
                  vmin_plot = vmin_plot/h
                  vmax_plot = vmax_plot/h
                  v_tic_major = v_tic_major/h
                  v_exponent = iuser_vscale
c                 Re-display tic/scale data
                  call display_real(16, 15, vmin_plot, 1)
                  call display_real(16, 28, vmax_plot, 1)
                  call display_int (19, 29, v_exponent, 3)
                  call display_real(17, 23, v_tic_major, 1)
                  call display_int (18, 27, v_minor_tic_freq, 5)
               endif
            endif
c
         elseif ( item_number .eq. 15 ) then
c           Logarithm
            do_vlog = .NOT. do_vlog
            if( do_vlog ) then
               prev_vabs = do_vabs
               do_vabs = .TRUE.
               str_vlog = 'log'
               str_vabs = 'yes'
c              Save previous values
               prev_v_exponent = v_exponent
               prev_iuser_vscale = iuser_vscale
               prev_vmin_plot = vmin_plot
               prev_vmax_plot = vmax_plot
               prev_v_tic_major  = v_tic_major
               prev_v_minor_tic_freq = v_minor_tic_freq
               prev_v_saved = .TRUE.
               iuser_vscale = 0
c
               call recompute_minmax(work_vec(1,1), ipts,
     1                   vmin_full, vmax_full,
     2                   vmin_plot, vmax_plot,
     3                   v_tic_major, v_minor_tic_freq,
     4                   iuser_vscale, v_exponent,
     5                   do_vabs, do_vlog, do_vuscale, work_vec(1,3))
               if( iuser_vscale .ne. 0 ) v_exponent = iuser_vscale
c              Mark Y axis tics as not applicable
               call display_string(17, 23, '       NA', 9, 1)
               call display_string(18, 27, '   NA', 5, 1)
            else
               do_vabs = prev_vabs
               str_vlog = 'linear'
               if( do_vabs ) then
                  str_vabs = 'yes'
               else
                  str_vabs = 'no'
               endif
c              Recompute vmin_full/vmax_full
               call recompute_minmax(work_vec(1,1), ipts,
     1                   vmin_full, vmax_full,
     2                   vmin_plot, vmax_plot,
     3                   v_tic_major, v_minor_tic_freq,
     4                   iuser_vscale, v_exponent,
     5                   do_vabs, do_vlog, do_vuscale, work_vec(1,3))
c              Restore previous axis values
               if( prev_v_saved ) then
                  iuser_vscale = prev_iuser_vscale
                  v_exponent = prev_v_exponent
                  vmin_plot = prev_vmin_plot
                  vmax_plot = prev_vmax_plot
                  v_tic_major = prev_v_tic_major
                  v_minor_tic_freq = prev_v_minor_tic_freq
               endif
               if( iuser_vscale .ne. 0 ) v_exponent = iuser_vscale
c              Re-display tic/scale data
               call display_real  (17, 23, v_tic_major, 1)
               call display_int   (18, 27, v_minor_tic_freq, 5)
            endif
            call display_real  (14, 15, vmin_full, 1)
            call display_real  (14, 28, vmax_full, 1)
            call display_real  (16, 15, vmin_plot, 1)
            call display_real  (16, 28, vmax_plot, 1)
            call display_int   (19, 29, v_exponent, 3)
            call display_string(14, 59, str_vlog, 6, 1)
            call display_string(15, 63, str_vabs, 3, 1)
c
c           Absolute value
         elseif ( item_number .eq. 16 ) then
            if( do_vlog .AND. do_vabs ) then
               call display_message(
     1            'Absolute value required for log scale')
            else
               prev_vabs = do_vabs
               do_vabs = .NOT. do_vabs
               if( do_vabs ) then
                  str_vabs = 'yes'
               else
                  str_vabs = 'no'
               endif
            endif
            call recompute_minmax(work_vec(1,1), ipts,
     1                vmin_full, vmax_full,
     2                vmin_plot, vmax_plot,
     3                v_tic_major, v_minor_tic_freq,
     4                iuser_vscale, v_exponent,
     5                do_vabs, do_vlog, do_vuscale, work_vec(1,3))
            if( do_vuscale ) v_exponent = iuser_vscale
c           Re-display tic/scale data
            call display_real  (17, 23, v_tic_major, 1)
            call display_int   (18, 27, v_minor_tic_freq, 5)
            call display_int   (19, 29, v_exponent, 3)
            call display_string(15, 63, str_vabs, 3, 1)
            call display_real  (14, 15, vmin_full, 1)
            call display_real  (14, 28, vmax_full, 1)
            call display_real  (16, 15, vmin_plot, 1)
            call display_real  (16, 28, vmax_plot, 1)
c
         elseif ( item_number .eq. 17 ) then
c           Autoscale horizontal axis
            do_vuscale = .FALSE.
            iuser_vscale = 0
            call recompute_minmax(work_vec(1,1), ipts,
     1                rdum1, rdum2,
     2                vmin_plot, vmax_plot,
     3                v_tic_major, v_minor_tic_freq,
     4                iuser_vscale, v_exponent,
     5                do_vabs, do_vlog, do_vuscale, work_vec(1,3))
c           Re-display tic/scale data
            if( .NOT. do_vlog ) then
               call display_real  (17, 23, v_tic_major, 1)
               call display_int   (18, 27, v_minor_tic_freq, 5)
            endif
            call display_int   (19, 29, v_exponent, 3)
            call display_string(15, 63, str_vabs, 3, 1)
            call display_real  (16, 15, vmin_plot, 1)
            call display_real  (16, 28, vmax_plot, 1)
c
         elseif ( item_number .eq. 18 ) then
            call get_string_with_prompt(vert_label,
     1              'Enter vertical label (may use mathmode)')
            call display_string(17, 54, vert_label, 27, 1)
c
         elseif ( item_number .eq. 19 ) then
            call get_string_with_prompt(vert_units,
     1              'Enter units (may use mathmode)')
            call display_string(18, 54, vert_units, 20, 1)
c
         elseif ( item_number .eq. 20 ) then
c           Markers
            do_markers = .NOT. do_markers
            if( do_markers ) then
               markers = 'yes'
            else
               markers = 'no'
            endif
            call display_string(19, 56, markers, 3, 1)
         else
            call display_message('Item number out of range')
         endif
         goto 2000
c
1000     call display_message('Illegal input')
         goto 2000
      endif
c
5000  continue
c
      return
      end
c
      subroutine do_histo_plot(data_vec, dist, num_points,
     1     h_win_min, h_win_max,
     2     h_exponent, h_tic_major, h_minor_tic_freq,
     3     v_win_min, v_win_max,
     4     v_exponent, v_tic_major, v_minor_tic_freq,
     5     horiz_label, hunits, vert_label, vunits,
     6     v_log, v_abs, h_log, h_abs,
     7     cut_line_string, do_oxide, do_header, do_markers,
     8     add_notes_from_save, ws_window_open, add_notes_from_file,
     9     ierror)
c
      implicit none
c
c     Subroutine which actually calls graphics routines.
c     NOTE: values in calls to tune_charsize were painfully picked to give
c     readable results with VTxxx terminals. If you change the
c     size of windows where text is written, expect the characters to be
c     screwed up, clipped outside the plot window, have terrible
c     aspect ratios, etc. etc. etc.
c
      logical do_oxide, do_header, do_markers, v_log, v_abs,
     1        add_notes_from_save, ws_window_open, add_notes_from_file,
     2        h_log, h_abs
      character*(*) cut_line_string
      character*(*) horiz_label, vert_label, hunits, vunits
      integer num_points
      real data_vec(1)
      real dist(1)
      real h_win_min, h_win_max, v_win_min, v_win_max
      real h_tic_major, v_tic_major
      integer h_exponent, v_exponent
      integer h_minor_tic_freq, v_minor_tic_freq
      integer ierror
c
      include 'const.h'
      include 'sysdep.cmn'
      include 'plot.cmn'
      include 'io.cmn'
      include 'user.cmn'
c
c     locals
c
      logical exp_flag, lerr, is_vtxxx
      integer i, j, ilen, ilen2, ilen3, window_id, return_stat
      integer limit, string_len
      integer get_new_window_id
      character*1 response
      character*8 time_string
      character*12 num_label, num_label1, num_label2
      character*85 work_string
      real fuzz, delta, h_cur, v_cur, v_minor, ticvalue
      real h_scale, v_scale, v_max, h_max, h_minor, scaled_zero
c
      real temp_vec(COMP_ARRAYSIZE), temp_coord(COMP_ARRAYSIZE)
      real hscl, vscl, ZEROO, ONEE
      integer k, l, npts
      real rleft, rrt
c     Two parameters close to zero and one, but within the range (0,1)
      parameter( ZEROO = 0.000001 )
      parameter( ONEE  = 0.999999 )
c
      fuzz = 0.999
c
c     Scale plot to 0.0 to 1.0 (avoid bug in GKS V5.0)
c     Compute horizontal and vertical scale factors
      hscl = h_win_max - h_win_min
      vscl = v_win_max - v_win_min
      if( hscl .le. 0.0 .or. vscl .le. 0.0 ) then
         ierror = PLT$INVALID_MINMAX
         if( lbatch .or. lrestore ) then
            call plt$get_error_message(return_stat,work_string)
            ilen = string_len(work_string)
            call plt$error_message('do_histo_plot', work_string(1:ilen))
         endif
         return
      endif
c
c     Meaning of {h,v}_scale is now inverse of pre 29-Apr-93 code
      h_scale = 10.0**(-h_exponent)
      v_scale = 10.0**(-v_exponent)
      hscl = 1.0/hscl
      vscl = 1.0/vscl
c
c     open plot device
c
      window_id = get_new_window_id(PLT$1D_WINDOW_ID)
      if( window_id .eq. -1 ) then
c        This should only happen when we are using the retained
c        style, and we run out of plot window ids.
         ierror = PLT$TOO_MANY_WINDOWS
         if( lbatch .or. lrestore ) then
            call plt$get_error_message(ierror,work_string)
            ilen = string_len(work_string)
            call plt$error_message('do_histo_plot', work_string(1:ilen))
         endif
         return
      endif
c
      call open_plot_device(window_id, plot_device, 'Postmini',
     1                      return_stat)
      if( return_stat .ne. PLT$NORMAL ) then
         if( lbatch .or. lrestore ) then
            call plt$get_error_message(return_stat,work_string)
            ilen = string_len(work_string)
            call plt$error_message('do_histo_plot', work_string(1:ilen))
         endif
         ierror = return_stat
         return
      endif
c
c     Activate plot device, clear it, start segment on VAXstations
      call activate_plot_device(window_id, .TRUE., .TRUE.)
      call pop_vs_window(window_id)
      call init_device_depend
      call GKS$SET_WINDOW(1, 0.0,1.0, 0.0,1.0)
c     create virtual plot window at 65% of the screen (NDC space)
      call GKS$SET_VIEWPORT(1, 0.175*ratio_x, 0.825*ratio_x,
     1                         0.1*ratio_y, 0.80*ratio_y)
      call GKS$SET_CLIPPING(GKS$K_NOCLIP)
      call GKS$SELECT_XFORM(1)
c
c    Draw box around plot.
c
      call set_line_index(DEFAULT_LINE)
      call move_abs(ZEROO, ZEROO)
      call draw_abs(ONEE,  ZEROO)
      call draw_abs(ONEE,  ONEE)
      call draw_abs(ZEROO, ONEE)
      call draw_abs(ZEROO, ZEROO)
c
      call GKS$SET_CLIPPING(1, GKS$K_CLIP)
      call GKS$SELECT_XFORM(1)
c
c     plot data
      call GKS$SELECT_XFORM(1)
      call set_line_index(POLY_LINE)
c
c     Since underlying window system may not clip correctly when
c     zooming in on a small portion of data, try to limit polyline
c     to minimum. This happens most often in practive when zooming on absissa.
c
c     To avoid bug in GKS V5.0, scale data to 0.0 to 1.0
      do 50 i=1,num_points
         if( h_log ) then
            temp_coord(i) = (ALOG10(MAX(PM_MIN_FLOAT,ABS(dist(i))))
     1                      *h_scale - h_win_min)*hscl
         else
            temp_coord(i) = (dist(i)*h_scale - h_win_min)*hscl
         endif
         temp_vec(i) = (data_vec(i)*v_scale - v_win_min)*vscl
50    continue
      rleft = 0.0
      k = num_points+1
      do 100 i=1,num_points
         if( temp_coord(i) .ge. rleft ) then
            k = MAX(1,i-1)
            goto 101
         endif
100   continue
101   continue
      rrt = 1.0
      l = 0
      do 102 i=num_points,1,-1
         if( temp_coord(i) .le. rrt ) then
            l = MIN(num_points,i+1)
            goto 103
         endif
102   continue
103   continue
      npts = l-k+1
c
      if( npts .gt. 0 ) then
         scaled_zero = -v_win_min*vscl
         call draw_histogram(npts, temp_coord(k), temp_vec(k),
     1                       scaled_zero)
      endif
c
      if( do_markers ) then
         if( npts .gt. 0 ) then
            call GKS$SET_PMARK_TYPE(GKS$K_MARKERTYPE_CIRCLE)
            call GKS$POLYMARKER(npts, temp_coord(k), temp_vec(k))
         endif
      endif
c
c     do_oxide no longer supported. jvf
c
      call set_line_index(DEFAULT_LINE)
c
c     plot x-axis tic marks...
      if( .not. h_log ) then
         call GKS$SET_WINDOW(1, 0.0,1.0, 0.0,1.0)
         call GKS$SELECT_XFORM(1)
         call set_line_index(DEFAULT_LINE)
c
c        bottom of axis marks...
c        Find first major tic mark below h_win_min
         if( h_win_min .lt. 0.0 ) then
            h_cur = FLOAT(INT(h_win_min/h_tic_major) - 1)*h_tic_major
         else
            h_cur = FLOAT(INT(h_win_min/h_tic_major))*h_tic_major
         endif
c        Normalize to 0.0 to 1.0 scale
         h_cur = (h_cur - h_win_min)/(h_win_max - h_win_min)
         delta = h_tic_major/
     1           (FLOAT(h_minor_tic_freq)*(h_win_max - h_win_min))
         i = 0
c        We use 1.0001 instead of 1.0 as a safeguard agains fp rounding
         do while(h_cur .le. 1.0001)
c          Plot tic marks only if within axis bounds
c           We use -0.0001 instead of 0.0 as a safeguard agains fp rounding
            if(h_cur .ge. -0.0001) then
               call move_abs(h_cur,0.0)
c              Determine where a big tic mark should go
               if( MOD(i,h_minor_tic_freq) .eq. 0 ) then
                  call draw_abs(h_cur,0.05)
               else
                  call draw_abs(h_cur,0.025)
              endif
            endif
            h_cur = h_cur + delta
            i = i + 1
         enddo
c        top of axis marks...
         if( h_win_min .lt. 0.0 ) then
            h_cur = FLOAT(INT(h_win_min/h_tic_major) - 1)*h_tic_major
         else
            h_cur = FLOAT(INT(h_win_min/h_tic_major))*h_tic_major
         endif
c        Normalize to 0.0 to 1.0 scale
         h_cur = (h_cur - h_win_min)/(h_win_max - h_win_min)
         i = 0
c        We use 1.0001 instead of 1.0 as a safeguard agains fp rounding
         do while(h_cur .le. 1.0001)
c           Plot tic marks only if within axis bounds
c           We use -0.0001 instead of 0.0 as a safeguard agains fp rounding
            if(h_cur .ge. -0.0001) then
               call move_abs(h_cur,1.0)
c              Determine where a big tic mark should go
               if( MOD(i,h_minor_tic_freq) .eq. 0 ) then
                  call draw_abs(h_cur,0.95)
               else
                  call draw_abs(h_cur,0.975)
               endif
            endif
            h_cur = h_cur + delta
            i = i + 1
         enddo
      else
c        log plot; tics every decade
         call GKS$SET_WINDOW(1, h_win_min,h_win_max,0.0,1.0)
         call GKS$SELECT_XFORM(1)
c        bottom of plot...
         h_cur = FLOAT(INT(h_win_min)) - 1.0
         h_max = FLOAT(INT(h_win_max)) + 1.0
         limit = INT((h_max - h_cur)/fuzz)
         do i=1,limit
            if( h_cur .ge. h_win_min .or. h_cur .le. h_win_max) then
               call move_abs(h_cur, 0.0)
               call draw_abs(h_cur, 0.025)
            endif
c           minor tic marks...
            do j=2,9
               h_minor = h_cur + ALOG10(FLOAT(j))
               if( h_minor .ge. h_win_min .or. h_minor .le. h_win_max)
     1            then
                  call move_abs(h_minor,0.0)
                  call draw_abs(h_minor,0.0125)
               endif
            end do
            h_cur = h_cur + 1.0
         end do
         if( h_cur .le. h_win_max) then
            call move_abs(h_cur, 0.0)
            call draw_abs(h_cur, 0.025)
         endif
c
c        top of plot...
         h_cur = FLOAT(INT(h_win_min)) - 1.0
         h_max = FLOAT(INT(h_win_max)) + 1.0
         limit = INT((h_max - h_cur)/fuzz)
         do i=1,limit
            if( h_cur .ge. h_win_min .or. h_cur .le. h_win_max) then
               call move_abs(h_cur, 1.0)
               call draw_abs(h_cur, 0.975)
            endif
c           minor tic marks...
            do j=2,9
               h_minor = h_cur + ALOG10(FLOAT(j))
               if( h_minor .ge. h_win_min .or. h_minor .le. h_win_max)
     1            then
                  call move_abs(h_minor, 1.0)
                  call draw_abs(h_minor, 0.9875)
               endif
            end do
            h_cur = h_cur + 1.0
         end do
         if( h_cur .ge. h_win_min .or. h_cur .le. h_win_max) then
            call move_abs(h_cur,1.0)
            call draw_abs(h_cur,0.975)
         endif
      endif
c
c     plot y-axis tic marks...
      if ( .not. v_log ) then
c
         call GKS$SET_WINDOW(1, 0.0,1.0,0.0,1.0)
         call GKS$SELECT_XFORM(1)
c        left side of plot...
c        Find first tic mark
         if( v_win_min .lt. 0.0 ) then
            v_cur = FLOAT(INT(v_win_min/v_tic_major) - 1)*v_tic_major
         else
            v_cur = FLOAT(INT(v_win_min/v_tic_major))*v_tic_major
         endif
c        Normalize to 0.0 to 1.0 scale
         v_cur = (v_cur - v_win_min)/(v_win_max - v_win_min)
         delta = v_tic_major/
     1           (FLOAT(v_minor_tic_freq)*(v_win_max - v_win_min))
         i = 0
c        We use 1.0001 instead of 1.0 as a safeguard agains fp rounding
         do while(v_cur .le. 1.0001)
c           Plot tic marks only if within axis bounds
c           We use -0.0001 instead of 0.0 as a safeguard agains fp rounding
            if( v_cur .ge. -0.0001 ) then
               call move_abs(0.0,v_cur)
c              Determine where a big tic mark should go
               if( MOD(i,v_minor_tic_freq) .eq. 0 ) then
                  call draw_abs(0.05,v_cur)
               else
                  call draw_abs(0.025,v_cur)
               endif
            endif
            v_cur = v_cur + delta
            i = i + 1
         enddo
c        right side of plot...
         if( v_win_min .lt. 0.0 ) then
            v_cur = FLOAT(INT(v_win_min/v_tic_major) - 1)*v_tic_major
         else
            v_cur = FLOAT(INT(v_win_min/v_tic_major))*v_tic_major
         endif
c        Normalize to 0.0 to 1.0 scale
         v_cur = (v_cur - v_win_min)/(v_win_max - v_win_min)
         i = 0
c        We use 1.0001 instead of 1.0 as a safeguard agains fp rounding
         do while(v_cur .le. 1.0001)
c           Plot tic marks only if within axis bounds
c           We use -0.0001 instead of 0.0 as a safeguard agains fp rounding
            if( v_cur .ge. -0.0001 ) then
               call move_abs(1.0,v_cur)
c              Determine where a big tic mark should go
               if( MOD(i,v_minor_tic_freq) .eq. 0 ) then
                  call draw_abs(0.95,v_cur)
               else
                  call draw_abs(0.975,v_cur)
               endif
            endif
            v_cur = v_cur + delta
            i = i + 1
         enddo
      else
c
c        log plot; tics every decade
         call GKS$SET_WINDOW(1, 0.0,1.0,v_win_min,v_win_max)
         call GKS$SELECT_XFORM(1)
c        left side of plot...
         v_cur = FLOAT(INT(v_win_min)) - 1.0
         v_max = FLOAT(INT(v_win_max)) + 1.0
         limit = INT((v_max - v_cur)/fuzz)
         do 620 i=1,limit
            if( v_cur .ge. v_win_min .or. v_cur .le. v_win_max) then
               call move_abs(0.0, v_cur)
               call draw_abs(0.025, v_cur)
            endif
c           minor tic marks...
            do 621 j=2,9
               v_minor = v_cur + ALOG10(FLOAT(j))
               if( v_minor .ge. v_win_min .or. v_minor .le. v_win_max)
     1            then
                  call move_abs(0.0,v_minor)
                  call draw_abs(0.0125,v_minor)
               endif
621         continue
            v_cur = v_cur + 1.0
620      continue
         if( v_cur .le. v_win_max) then
            call move_abs(0.0, v_cur)
            call draw_abs(0.025, v_cur)
         endif
c
c        right side of plot...
         v_cur = FLOAT(INT(v_win_min)) - 1.0
         v_max = FLOAT(INT(v_win_max)) + 1.0
         limit = INT((v_max - v_cur)/fuzz)
         do 625 i=1,limit
            if( v_cur .ge. v_win_min .or. v_cur .le. v_win_max) then
               call move_abs(1.0,v_cur)
               call draw_abs(0.975,v_cur)
            endif
c           minor tic marks...
            do 626 j=2,9
               v_minor = v_cur + ALOG10(FLOAT(j))
               if( v_minor .ge. v_win_min .or. v_minor .le. v_win_max)
     1            then
                  call move_abs(1.0,v_minor)
                  call draw_abs(0.9875,v_minor)
               endif
626         continue
            v_cur = v_cur + 1.0
625      continue
         if( v_cur .ge. v_win_min .or. v_cur .le. v_win_max) then
            call move_abs(1.0,v_cur)
            call draw_abs(0.975,v_cur)
         endif
      endif
c
c     horizontal label...
c
c     Magic constants are amount window must be expanded on each side
c     to keep in proportion with previous plotting window,
c     but with full viewport in x. They scale with the ratio 13/20 =
c     (0.825 - 0.175)/1.0 for the current plot window setting.
c        The left x window:  13/20 = (0.0 - 0.175)/x  x = -7/26
c        The right x window: 13/20 = (1.0 - 0.825)/x  x =  7/26
      call GKS$SET_WINDOW(1, -0.269, 1.269, 0.0, 0.2)
      call GKS$SET_VIEWPORT(1, 0.0, ratio_x,
     1                         0.0*ratio_y,  0.1*ratio_y)
      call GKS$SET_CLIPPING(1, GKS$K_NOCLIP)
      call GKS$SELECT_XFORM(1)
      call tune_charsize(0.035, 0.06, -0.5, 0.025, 0.0458, 0.0)
      if ( h_log ) then
         call GKS$SET_TEXT_ALIGN(GKS$K_TEXT_HALIGN_CENTER,
     1                           GKS$K_TEXT_VALIGN_HALF)
c        Log axis. Label each major tic mark, adjust left and
c        right label to align with plot frame.
         h_cur = FLOAT(INT(h_win_min)) - 1
         h_minor = h_cur
c        Normalize to 0.0 to 1.0 scale
         h_cur = (h_cur - h_win_min)/(h_win_max - h_win_min)
         delta = 1.0/(h_win_max - h_win_min)
         do while( h_cur .le. 1.0001 )
            if( h_cur .ge. -0.0001) then
c              Format into mathmode string 10^{tic_value}
               call format_mathmode_exponent(h_minor,num_label1)
               call mathmode(ws_id, h_cur, 0.11, num_label1,
     1                       .FALSE., lerr)
            endif
            h_minor = h_minor + 1.0
            h_cur = h_cur + delta
         end do
      else
c        Write out labels at major tic intervals
         call GKS$SET_TEXT_ALIGN(GKS$K_TEXT_HALIGN_CENTER,
     1                           GKS$K_TEXT_VALIGN_HALF)
c        Find first major tic mark below h_win_min
         if( h_win_min .lt. 0.0 ) then
            h_cur = FLOAT(INT(h_win_min/h_tic_major) - 1)*h_tic_major
         else
            h_cur = FLOAT(INT(h_win_min/h_tic_major))*h_tic_major
         endif
c        Normalize to 0.0 to 1.0 scale
         h_cur = (h_cur - h_win_min)/(h_win_max - h_win_min)
         delta = h_tic_major/
     1           (FLOAT(h_minor_tic_freq)*(h_win_max - h_win_min))
         i = 0
c        We use 1.0001 instead of 1.0 as a safeguard agains fp rounding
         do while(h_cur .le. 1.0001)
c           Plot tic marks only if within axis bounds
c           We use -0.0001 instead of 0.0 as a safeguard agains fp rounding
            if(h_cur .ge. -0.0001) then
c              Determine when to plot a tic mark label
               if( MOD(i,h_minor_tic_freq) .eq. 0 ) then
                  ticvalue = (h_cur*(h_win_max - h_win_min) + h_win_min)
                  call format_axis_label(ticvalue, h_tic_major,
     1                                   num_label, exp_flag)
                  ilen = string_len(num_label)
                  call GKS$TEXT(h_cur, 0.12, num_label(1:ilen))
               endif
            endif
            h_cur = h_cur + delta
            i = i + 1
         enddo
      endif
c
c     Write horizontal label
      call tune_charsize(0.035, 0.06 , -0.5, 0.03, 0.05, 0.0)
      call GKS$SET_TEXT_ALIGN(GKS$K_TEXT_HALIGN_CENTER,
     1                        GKS$K_TEXT_VALIGN_BOTTOM)
      ilen = string_len(horiz_label)
      ilen3 = string_len(hunits)
      if( h_exponent .eq. 0 ) then
         if( ilen3 .eq. 0 ) then
            work_string = horiz_label(1:ilen)
         else
            work_string = horiz_label(1:ilen) // ' (' //
     1                    hunits(1:ilen3) // ')'
            ilen = ilen + ilen3 + 3
         endif
      else
         call format_mathmode_exponent(FLOAT(h_exponent),num_label1)
         ilen2 = string_len(num_label1)
         if( ilen3 .ne. 0 ) then
            work_string =  horiz_label(1:ilen) //
     1                     ' (' //  num_label1(1:ilen2) // ' ' //
     2                     hunits(1:ilen3) // ')'
            ilen = ilen + ilen2 + ilen3 + 4
         else
            work_string =  horiz_label(1:ilen) //
     1                     ' (' //  num_label1(1:ilen2) // ')'
            ilen = ilen + ilen2 + 3
         endif
      endif
      call mathmode(ws_id, 0.5, 0.0, work_string(1:ilen),.FALSE., lerr)
c
c     vertical tic mark numbering...
c
c     Magic constants are amount window must be expanded on each side
c     to keep in proportion with previous plotting window,
c     but with full viewport in y. They scale with the ratio 7/10 =
c     (0.8 - 0.1)/1.0 for the current plot window setting.
c        The lower y window:  7/10 = (0.0 - 0.1)/x  x = -1/7
c        The upper y window:  7/10 = (1.0 - 0.8)/x  x =  2/7
      call GKS$SET_WINDOW  (1, 0.0, 0.15, -0.1429, 1.2857)
      call GKS$SET_VIEWPORT(1, 0.0*ratio_x, 0.175*ratio_x,
     1                         0.0*ratio_y, ratio_y)
      call GKS$SET_CLIPPING(1, GKS$K_NOCLIP)
      call GKS$SELECT_XFORM(1)
      if ( v_log ) then
c
c     Log axis. Label each major tic mark, adjust bottom and
c     top label to align with plot frame.
c
c         call tune_charsize(0.015, 0.04, -0.2, 0.015, 0.0225, 0.0)
         call tune_charsize(0.015, 0.057, -0.2, 0.015, 0.032, 0.0)
         call GKS$SET_TEXT_ALIGN(GKS$K_TEXT_HALIGN_RIGHT,
     1                           GKS$K_TEXT_VALIGN_HALF)
         v_cur = FLOAT(INT(v_win_min)) - 1
         v_minor = v_cur
c        Normalize to 0.0 to 1.0 scale
         v_cur = (v_cur - v_win_min)/(v_win_max - v_win_min)
         delta = 1.0/(v_win_max - v_win_min)
         do while( v_cur .le. 1.0001 )
            if( v_cur .ge. -0.0001) then
c              Format into mathmode string 10^{tic_value}
               call format_mathmode_exponent(v_minor,num_label1)
               call mathmode(ws_id, 0.14, v_cur, num_label1,
     1                       .FALSE., lerr)
            endif
            v_minor = v_minor + 1.0
            v_cur = v_cur + delta
         end do
      else
c
c     Linear axis
c
         if( exp_flag ) then
c           Use smaller font so we don't run into vertical label
c            call tune_charsize(0.015, 0.04, -0.2, 0.012, 0.02, 0.0)
            call tune_charsize(0.015, 0.057, -0.2, 0.012, 0.029, 0.0)
         else
c            call tune_charsize(0.015, 0.04, -0.2, 0.015, 0.0225, 0.0)
            call tune_charsize(0.015, 0.057, -0.2, 0.015, 0.032, 0.0)
         endif
         call GKS$SET_TEXT_ALIGN(GKS$K_TEXT_HALIGN_RIGHT,
     1                           GKS$K_TEXT_VALIGN_HALF)
c        Write labels on major tic marks...
c        Find first tic mark
         if( v_win_min .lt. 0.0 ) then
            v_cur = FLOAT(INT(v_win_min/v_tic_major) - 1)*v_tic_major
         else
            v_cur = FLOAT(INT(v_win_min/v_tic_major))*v_tic_major
         endif
c        Normalize to 0.0 to 1.0 scale
         v_cur = (v_cur - v_win_min)/(v_win_max - v_win_min)
         delta = v_tic_major/
     1           (FLOAT(v_minor_tic_freq)*(v_win_max - v_win_min))
         i = 0
c        We use 1.0001 instead of 1.0 as a safeguard agains fp rounding
         do while(v_cur .le. 1.0001)
c           Plot tic marks only if within axis bounds
c           We use -0.0001 instead of 0.0 as a safeguard agains fp rounding
            if( v_cur .ge. -0.0001 ) then
c              Determine where a big tic mark should go
               if( MOD(i,v_minor_tic_freq) .eq. 0 ) then
                  ticvalue = (v_cur*(v_win_max - v_win_min) + v_win_min)
                  call format_axis_label(ticvalue, v_tic_major,
     1                 num_label2, exp_flag)
                  ilen = string_len(num_label2)
                  call GKS$TEXT(0.13, v_cur, num_label2(1:ilen))
               endif
            endif
            v_cur = v_cur + delta
            i = i + 1
         enddo
      endif
c
c     vertical label
c
c      call tune_charsize(0.024, 0.02, 0.0, 0.03, 0.02, 0.0)
      call tune_charsize(0.034, 0.02, 0.0, 0.043, 0.02, 0.0)
      call GKS$SET_TEXT_UPVEC(-1.0,0.0)
      call GKS$SET_TEXT_ALIGN(GKS$K_TEXT_HALIGN_CENTER,
     1                        GKS$K_TEXT_VALIGN_HALF)
      ilen = string_len(vert_label)
      ilen3 = string_len(vunits)
      if( v_exponent .eq. 0 ) then
         if( ilen3 .ne. 0 ) then
            work_string = vert_label(1:ilen) // ' (' // vunits(1:ilen3)
     1                    // ')'
            ilen = ilen + ilen3 + 3
         else
            work_string = vert_label(1:ilen)
         endif
      else
         call format_mathmode_exponent(FLOAT(v_exponent),num_label1)
         ilen2 = string_len(num_label1)
         if( ilen3 .ne. 0 ) then
            work_string =  vert_label(1:ilen) //
     1                     ' (' //  num_label1(1:ilen2) // ' ' //
     2                     vunits(1:ilen3) // ')'
            ilen = ilen + ilen2 + ilen3 + 4
         else
            work_string =  vert_label(1:ilen) //
     1                     ' (' //  num_label1(1:ilen2) // ')'
            ilen = ilen + ilen2 + 3
         endif
      endif
      call mathmode(ws_id, 0.05, 0.45, work_string(1:ilen),
     1                 .FALSE., lerr)
      call GKS$SET_TEXT_UPVEC(0.0,1.0)
      call GKS$SET_TEXT_ALIGN(GKS$K_TEXT_HALIGN_NORMAL,
     1                        GKS$K_TEXT_VALIGN_NORMAL)
c
c     titling...
c
      call GKS$SET_WINDOW(1, 0.0,1.0,0.0,0.2)
      call GKS$SET_VIEWPORT(1, 0.0*ratio_x, 1.0*ratio_x,
     1                         0.80*ratio_y, 1.0*ratio_y)
      call GKS$SET_CLIPPING(1, GKS$K_NOCLIP)
      call GKS$SELECT_XFORM(1)
c     Gag: yet another set of magic constants :-( jvf
      call tune_charsize(0.025, 0.04, -0.5, 0.02, 0.025, 0.0)
      call GKS$SET_TEXT_ALIGN(GKS$K_TEXT_HALIGN_NORMAL,
     1                        GKS$K_TEXT_VALIGN_NORMAL)
      if( do_header ) then
         call mathmode(ws_id, 0.0, 0.15, plot_title,
     1                 .FALSE., lerr)
      endif
c
c     Print auxiliary information
c
      call TIME(time_string)
      ilen = string_len(vert_label)
      if( v_abs ) then
         work_string = vert_label(1:ilen) // ' (abs) ' //
     1         cut_line_string // ';  ' //
     1         date_string  // ' ' // time_string
      else
         work_string = vert_label(1:ilen) // ' ' //
     1         cut_line_string // ';  ' //
     1         date_string  // ' ' // time_string
      endif
      ilen = string_len(work_string)
      call GKS$TEXT(0.0, 0.10, work_string(1:ilen))
      ilen = string_len(binfilename)
      work_string =  'Data source: ' // binfilename(1:ilen)
      ilen = string_len(work_string)
      call GKS$TEXT(0.0, 0.05, work_string(1:ilen))
c
      call deactivate_plot_device
c
      if( plot_device .eq. hardcopy_device ) then
c        Add in annotations to hardcopy plot
         if( ws_window_open ) then
            call anno_replay_annotations(current_window_id, window_id,
     1                 .TRUE.)
         else
            if( add_notes_from_save ) then
               call anno_restore_notes_from_save(window_id, lerr)
            endif
         endif
         call close_plot_device(window_id)
      else
         if( add_notes_from_save ) then
            call anno_restore_notes_from_save(window_id, lerr)
         else if ( add_notes_from_file ) then
c           Restore any previous annotations
            call read_savefile('POSTMINI.NOTE',.TRUE.,lerr)
            if( .NOT. lerr )
     1         call anno_restore_notes_from_save(window_id, lerr)
         endif
         if( is_vtxxx(window_id) ) then
            call move_cursor_bottom_of_screen
            call beep_terminal
            call get_char(response,'Hit return to continue ',0)
c           Restore default color map
            call reset_vt_colormap
         endif
      endif
c
      ierror = PLT$NORMAL
c
      end
c
      subroutine draw_histogram(npts, coords, data, scaled_zero)
      implicit none
c
      integer npts
      real coords(npts), data(npts), scaled_zero
c
      include 'const.h'
      include 'sysdep.cmn'
      include 'plot.cmn'
c
      integer i
      real x(4), y(4), width, hwidth
c
      width = coords(2) - coords(1) !  Assumes uniform spacing of histogram
      hwidth = 0.25*width
c
      call GKS$SET_FILL_INT_STYLE(GKS$K_INTSTYLE_SOLID)
      call GKS$SET_FILL_COLOR_INDEX(1)
      do i=1,npts
         x(1) = coords(i) - hwidth
         x(2) = coords(i) + hwidth
         x(3) = x(2)
         x(4) = x(1)
         y(1) = scaled_zero
         y(2) = scaled_zero
         y(3) = data(i)
         y(4) = data(i)
         call GKS$FILL_AREA(4,x,y)
      end do
c
      return
      end

