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
      subroutine restore_plot
      implicit none
c
      include 'sysdep.cmn'
      include 'io.cmn'
c
      character*1  ch
      character*(PM_MAX_PATHLEN) savefile
      logical err
c
      call clear_screen
      call get_str_with_default(savefilename,
     1        'Enter Postmini command file name: ',0)
      savefile = savefilename
      call do_restore_plot(savefile,err)
      if( err .AND. .NOT. lbatch ) then
         call get_char(ch,'Hit return to continue... ',0)
      endif
c
      return
      end
c
      subroutine do_restore_plot(savefile,err)
      implicit none
c
      include 'sysdep.cmn'
      include 'const.h'
      include 'io.cmn'
      include 'plot.cmn'
      include 'units.cmn'
      include 'defaults.cmn'
      include 'compare.cmn'
      include 'annotations.cmn'
      include 'save.cmn'
      include 'gui.cmn'
c
      character*(*) savefile
      logical err
c
c     GUI error level: 0=success, 1=warning, 2=error
c     Locals
      integer nrecl, ierr, ioerr, istat, pm_put_error_msg,
     1        pm_put_info_msg
c
      err = .FALSE.
      gui_error_level = 0
c
c     Read savefile into save memory, including annotations
      call copy_init_axis_to_save
      call copy_init_curve_to_save
      call copy_init_params_to_save
      call copy_init_notes_to_save
      call read_savefile(savefile,.FALSE.,err)
      if( err ) then
         if( lgui ) then
            istat = pm_put_error_msg(
     1              'Could not open PMI file, or it had syntax errors')
            gui_error_level = 2
         endif
         return
      endif
c
c     Store annotations in POSTMINI.NOTE file for later recall,
c     if there are any...
      if( num_boxes .ne. 0 .or. num_lines .ne. 0 .or. num_text .ne. 0
     1    .or. num_markers .ne. 0 .or. num_ellipses .ne. 0 ) then
         istat = pm_put_info_msg('Annotations from command file ' //
     1           'will be placed in file POSTMINI.NOTE')
c
         nrecl=0
         call pm_open_file(notefil, 'POSTMINI.NOTE',
     1                     'UNKNOWN','FORMATTED','READWRITE','LIST',
     2                     'NONSHARED',nrecl, ioerr, ierr)
         if( ierr .ne. 0 ) goto 510
         rewind(notefil)
c
c        write annotations to save file
         call write_save_annotations(notefil)
c
         close(unit=notefil)
         goto 520
510      istat = pm_put_error_msg(
     1        'Warning: could not open file POSTMINI.NOTE')
         istat = pm_put_error_msg(
     1        'Annotations will not be saved in that file')
         gui_error_level = 1
520      continue
      endif
c
c     Do some rudimentary checking of plot parameters for validity
      call check_save_plot_params(err)
      if( err ) then
         gui_error_level = 2
         istat = pm_put_error_msg(
     1           'Error: PMI file had invalid parameters')
         return
      endif
c
c     Switch on type of plot
      if( save_plot_type .eq. PLOT_NULL_COMMAND ) then
         istat = pm_put_error_msg(
     1        'Error: NULL plot type specified, cannot restore')
         err = .TRUE.
         gui_error_level = 2
      else if( save_plot_type .eq. PLOT_COMPARE_COMMAND ) then
         call copy_save_to_compare(err)
         if( err ) then
            istat = pm_put_error_msg('Error while restoring plot')
            gui_error_level = 2
            return
         endif
         if( lgui ) then
            call pm_copy_save_axis_to_compare
         else
            call restore_compare_plot(err,save_plot_subtype)
         endif
      else if( save_plot_type .eq. PLOT_1D_COMMAND ) then
         if( lgui ) then
            istat = pm_put_error_msg(
     1           'PLOT 1D not supported, cannot restore')
            gui_error_level = 2
         else
            call load_datafile_from_save(err)
            if( err ) return
            call restore_1d_plot(err)
         endif
      else if( save_plot_type .eq. PLOT_2D_COMMAND ) then
         if( lgui ) then
            istat = pm_put_error_msg(
     1           'PLOT 2D not supported yet, cannot restore')
            gui_error_level = 2
         else
            call load_datafile_from_save(err)
            if( err ) return
            call restore_2d_plot(err)
         endif
      else if( save_plot_type .eq. PLOT_OVERLAY_COMMAND ) then
         if( lgui ) then
            istat = pm_put_error_msg(
     1           'PLOT OVERLAY not supported yet, cannot restore')
            gui_error_level = 2
         else
            call copy_save_to_multi_2d(err)
            if( err ) then
               istat = pm_put_error_msg('Error while restoring plot')
               gui_error_level = 2
               return
            endif
            call restore_multi_2d_plot(err)
         endif
      else if( save_plot_type .eq. PLOT_3D_COMMAND ) then
         if( lgui ) then
            istat = pm_put_error_msg(
     1           'PLOT 3D not yet supported, cannot restore')
            gui_error_level = 2
         else
            call load_datafile_from_save(err)
            if( err ) return
            call restore_3d_plot(err)
         endif
      else
         istat = pm_put_error_msg('Unknown plot type, cannot restore')
         err = .TRUE.
         gui_error_level = 2
      endif
c
      return
      end
c
      subroutine copy_save_to_compare(err)
      implicit none
c
      include 'sysdep.cmn'
      include 'const.h'
      include 'io.cmn'
      include 'plot.cmn'
      include 'mesh.cmn'
      include 'data.cmn'
      include 'units.cmn'
      include 'defaults.cmn'
      include 'compare.cmn'
      include 'user.cmn'
      include 'analytical.cmn'
      include 'annotations.cmn'
      include 'save.cmn'
c
      logical err
c
      integer i, j, k, l, entry, ierr, index, menu_id, np, ip,
     1        itype, find_new_entry, delete_entry
      logical read_error, load_error, bad_entry, update_datafile,
     1        lfound
c     Triangular grid cut line extraction can take up lots of storage
c     temporarily
      integer TMPSIZE
      parameter (TMPSIZE = 4*COMP_ARRAYSIZE)
      real    tmp(TMPSIZE,3)
      character*40 datatype
c
c     Mark any existing compare data as deleted, 
c     unless it has already been loaded via the callable interface.
c
      j = 0
      do 1000 i=1,MULTI_SIZE
         k = slot_index(i)
         if( k .ne. -1 ) then
            if( multi_1d(k)%datatype .ne. DT_CALLABLE_XY ) then
               in_use(i) = .FALSE.
               slot_index(i) = -1
            else
               j = j + 1
            endif
         endif
1000  continue
c
c     Store curve information
      n_cuts = 0
      do 2000 i=1,num_curves
c
c        First, look for XY data already loaded via callable interface
         if( curve_dat(i)%datafile_type .eq.
     1       SAVE_DATATYPE_CALLABLE_XY ) then
c           Look for matching entry ID string in data already loaded
            lfound = .FALSE.
            do 2001 k=1,MULTI_SIZE
               if( in_use(k) ) then
                  l = slot_index(k)
                  if( multi_1d(l)%datatype .eq. DT_CALLABLE_XY .and.
     1                multi_1d(l)%label .eq. curve_dat(i)%label ) then
                     lfound = .TRUE.
                  endif
               endif
2001        continue
            if( .NOT. lfound ) then
               write(luttyo,2002) curve_dat(i)%label
2002           format(' ','Error: callable XY data not found'/
     1                ' ','       label = ',a)
               err = .TRUE.
            endif
            goto 2000
         endif
c
c        Else, read in data files
         entry = find_new_entry()
         if( entry .eq. -1 ) then
            write(luttyo,2010)
2010        format(' ','Error: all data slots used up')
            err = .TRUE.
         else
            bad_entry = .FALSE.
            datatype = curve_dat(i)%datafile_type
c           Switch between ASCII/ANALYTICAL and simulator datatypes
c           (for 1D only)
            if( datatype .eq. SAVE_DATATYPE_ASCII .and.
     1          curve_dat(i)%datafile_dimension .eq. 1 ) then
               j = j + 1
c              Set up user common area
               user_1d_asciifile = curve_dat(i)%datafile
               user_1d_col_x = curve_dat(i)%colx
               user_1d_n_ycols = curve_dat(i)%n_ycols
               do l=1,user_1d_n_ycols
                  user_1d_col_y(l) = curve_dat(i)%coly(l)
               end do
               user_1d_col_z = curve_dat(i)%colz
               user_1d_skip  = curve_dat(i)%skip
               user_1d_nlines = curve_dat(i)%lines_to_read
               user_1d_ax = curve_dat(i)%ax
               user_1d_ay = curve_dat(i)%ay
               user_1d_az = curve_dat(i)%az
               user_1d_mx = curve_dat(i)%mx
               user_1d_my = curve_dat(i)%my
               user_1d_mz = curve_dat(i)%mz
               user_1d_x_expression = curve_dat(i)%x_expression
               user_1d_y_expression = curve_dat(i)%y_expression
               user_1d_xabs = curve_dat(i)%xabs
               user_1d_yabs = curve_dat(i)%yabs
               user_1d_zabs = curve_dat(i)%zabs
               user_1d_sort = curve_dat(i)%sort
               user_1d_sample_freq = curve_dat(i)%sample_freq
               user_1d_x_char_data = curve_dat(i)%x_char_data
               user_1d_rect_grid = .NOT. curve_dat(i)%data_is_trigrid
               user_1d_merge_ycols = curve_dat(i)%merge_ycols
               user_1d_sweep = curve_dat(i)%sweep
               user_1d_contour = curve_dat(i)%contour
               user_1d_errorbars = curve_dat(i)%errorbars
               user_1d_dim  = curve_dat(i)%datafile_dimension
               user_1d_dataset_number  = curve_dat(i)%datafile_record
               if( user_1d_rect_grid ) then
                  call do_read_asciifile
               else
                  call do_read_asciifile_2d_tri
               endif
c
               if( have_user_1d_data ) then
                    call copy_data_to_slot(entry,
     1                 arr2d(1,3), ni,
     2                 arr2d(1,2), 0.0, '?',
     3                 data_name(1), data_nickname(1),
     4                 data_units(1), COMP_USER_TYPE,
     4                 curve_dat(i)%int_limits,
     5                 curve_dat(i)%label,
     6                 curve_dat(i)%color, curve_dat(i)%linetype,
     7                 curve_dat(i)%linewidth,
     8                 curve_dat(i)%markertype,
     9                 curve_dat(i)%markercolor,
     A                 curve_dat(i)%markersize,
     B                 curve_dat(i)%xaxis_tag,
     C                 curve_dat(i)%yaxis_tag,
     D                 curve_dat(i)%suppress,
     E                 curve_dat(i)%mx, curve_dat(i)%ax,
     F                 curve_dat(i)%my, curve_dat(i)%ay,
     G                 curve_dat(i)%xabs, curve_dat(i)%yabs,
     H                 .FALSE.,
     I                 curve_dat(i)%x_expression,
     J                 curve_dat(i)%y_expression)
               else
                  bad_entry = .TRUE.
               endif
            else if( datatype .eq. SAVE_DATATYPE_ANALYTICAL .and.
     1          curve_dat(i)%datafile_dimension .eq. 1 ) then
               j = j + 1
c              Set up analytical common area
               binfiletype = DT_ANALYTICAL
               analytical_expression = 
     1            curve_dat(i)%analytical_expression
               analytical_xstart = curve_dat(i)%analytical_xstart
               analytical_xend   = curve_dat(i)%analytical_xend
               analytical_xstep  = curve_dat(i)%analytical_xstep
               analytical_pstart = curve_dat(i)%analytical_pstart
               analytical_pend   = curve_dat(i)%analytical_pend
               analytical_pstep  = curve_dat(i)%analytical_pstep
               analytical_type   = curve_dat(i)%analytical_type
               analytical_sweep  = curve_dat(i)%sweep
               analytical_dim    = curve_dat(i)%datafile_dimension
               call do_eval_analytical
c
               if( have_analytical_data ) then
                    call copy_data_to_slot(entry,
     1                 arr2d(1,3), ni,
     2                 arr2d(1,2), 0.0, '?',
     3                 data_name(1), data_nickname(1),
     4                 data_units(1), COMP_CUT_TYPE,
     4                 curve_dat(i)%int_limits,
     5                 curve_dat(i)%label,
     6                 curve_dat(i)%color, curve_dat(i)%linetype,
     7                 curve_dat(i)%linewidth,
     8                 curve_dat(i)%markertype,
     9                 curve_dat(i)%markercolor,
     A                 curve_dat(i)%markersize,
     B                 curve_dat(i)%xaxis_tag,
     C                 curve_dat(i)%yaxis_tag,
     D                 curve_dat(i)%suppress,
     E                 curve_dat(i)%mx, curve_dat(i)%ax,
     F                 curve_dat(i)%my, curve_dat(i)%ay,
     G                 curve_dat(i)%xabs, curve_dat(i)%yabs,
     H                 .FALSE.,
     I                 curve_dat(i)%x_expression,
     J                 curve_dat(i)%y_expression)
               else
                  bad_entry = .TRUE.
               endif
            elseif( datatype .eq. SAVE_DATATYPE_MINIMOS  .or.
     1              datatype .eq. SAVE_DATATYPE_PROMIS   .or.
     2              datatype .eq. SAVE_DATATYPE_VLSICAP  .or.
     3              datatype .eq. SAVE_DATATYPE_BAMBI    .or.
     4              datatype .eq. SAVE_DATATYPE_USEOUT   .or.
     5              datatype .eq. SAVE_DATATYPE_2DOP     .or.
     6              datatype .eq. SAVE_DATATYPE_SUPREM3S .or.
     7              datatype .eq. SAVE_DATATYPE_SUPREM3P .or.
     8              datatype .eq. SAVE_DATATYPE_SUPREM4  .or.
     9              datatype .eq. SAVE_DATATYPE_PISCES   .or.
     A              datatype .eq. SAVE_DATATYPE_TIF      .or.
     B              datatype .eq. SAVE_DATATYPE_FAIM     .or.
     C              datatype .eq. SAVE_DATATYPE_MCP2D    .or.
     C              datatype .eq. SAVE_DATATYPE_MEDICI   .or.
     D              datatype .eq. SAVE_DATATYPE_CALLABLE .or.
     F              datatype .eq. SAVE_DATATYPE_ASCII .or.
     G              datatype .eq. SAVE_DATATYPE_GRAPES .or.
     H              datatype .eq. SAVE_DATATYPE_SPICE3 ) then
               j = j + 1
c              Set up data area and read file
c              Optimization: if file is loaded, and name is same, don't re-read
               if( update_datafile(i) ) then
                  prev_binfilename = binfilename
                  binfilename     = curve_dat(i)%datafile
                  plane_value     = curve_dat(i)%cut_plane_coord
                  cut_plane_axis  = curve_dat(i)%cut_plane_axis
                  datafile_record = curve_dat(i)%datafile_record
                  binfile_convert_type = curve_dat(i)%bin_datatype
                  if( curve_dat(i)%datafile_dimension .eq. 3 ) then
                     three_d = .TRUE.
                  else
                     three_d = .FALSE.
                  endif
                  if( datatype .eq. SAVE_DATATYPE_MINIMOS     ) then
                     binfiletype = DT_MINIMOS
                     call read_mini(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_PROMIS ) then
                     binfiletype = DT_PROMIS
                     call read_promis(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_SUPREM4 ) then
                     binfiletype = DT_SUPREM4
                     call read_sup4(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_PISCES ) then
                     binfiletype = DT_PISCES
                     medici_solution_file = curve_dat(i)%misc_ch(1)
                     call read_pisces(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_MEDICI ) then
                     binfiletype = DT_MEDICI
                     medici_solution_file = curve_dat(i)%misc_ch(1)
                     call read_medici(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_VLSICAP) then
                     binfiletype = DT_VLSICAP
                     vlsi_edge = curve_dat(i)%misc(1)
                     vlsi_yint = curve_dat(i)%misc(2)
                     call read_vlsicap(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_BAMBI  ) then
                     binfiletype = DT_BAMBI
                     call read_bambi(read_error)
                  else if( datatype .eq. SAVE_DATATYPE_USEOUT ) then
                     binfiletype = DT_USEOUT
                     call read_useout(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_2DOP   ) then
                     binfiletype = DT_2DOP
                     tdop_gate_length = curve_dat(i)%misc(1)
                     tdop_interface_file = curve_dat(i)%misc_ch(1)
                     call read_2dop(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_SUPREM3S) then
                     binfiletype = DT_SUPREM3S
                     call read_suprem3(2,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_SUPREM3P) then
                     binfiletype = DT_SUPREM3P
                     call read_suprem3(1,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_CALLABLE ) then
                     binfiletype = DT_CALLABLE
                     call load_callable_data(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_TIF ) then
                     tif_ixcol = curve_dat(i)%colx
                     tif_iycol = curve_dat(i)%coly(1)
                     binfiletype = DT_TIF
                     call read_tif(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_FAIM ) then
                     binfiletype = DT_FAIM
                     call read_faim(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_MCP2D ) then
                     binfiletype = DT_MCP2D
                     call read_mcp2d(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_GRAPES ) then
                     binfiletype = DT_GRAPES
                     call read_grapesfile(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_SPICE3 ) then
                     binfiletype = DT_SPICE3
                     call do_read_spice3file(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_ASCII ) then
                     binfiletype = DT_ASCII
c                    Set up user common area
                     user_1d_asciifile = curve_dat(i)%datafile
                     user_1d_col_x = curve_dat(i)%colx
                     user_1d_n_ycols = curve_dat(i)%n_ycols
                     user_1d_col_y(1) = curve_dat(i)%coly(1)
                     user_1d_col_z = curve_dat(i)%colz
                     user_1d_skip  = curve_dat(i)%skip
                     user_1d_nlines = curve_dat(i)%lines_to_read
                     user_1d_ax = curve_dat(i)%ax
                     user_1d_ay = curve_dat(i)%ay
                     user_1d_az = curve_dat(i)%az
                     user_1d_mx = curve_dat(i)%mx
                     user_1d_my = curve_dat(i)%my
                     user_1d_mz = curve_dat(i)%mz
                     user_1d_x_expression = curve_dat(i)%x_expression
                     user_1d_y_expression = curve_dat(i)%y_expression
                     user_1d_xabs = curve_dat(i)%xabs
                     user_1d_yabs = curve_dat(i)%yabs
                     user_1d_zabs = curve_dat(i)%zabs
                     user_1d_sort = curve_dat(i)%sort
                     user_1d_x_char_data = curve_dat(i)%x_char_data
                     user_1d_sample_freq = curve_dat(i)%sample_freq
                     user_1d_rect_grid =
     1                   .NOT. curve_dat(i)%data_is_trigrid
                     user_1d_sweep = curve_dat(i)%sweep
                     user_1d_merge_ycols = curve_dat(i)%merge_ycols
                     user_1d_contour = curve_dat(i)%contour
                     user_1d_errorbars = curve_dat(i)%errorbars
                     user_1d_dim  = curve_dat(i)%datafile_dimension
                     user_1d_dataset_number=curve_dat(i)%datafile_record
                     if( user_1d_rect_grid ) then
                        call do_read_asciifile_2d
                     else
                        call do_read_asciifile_2d_tri
                     endif
                     read_error = .NOT. have_user_1d_data
                  else if( datatype .eq. SAVE_DATATYPE_ANALYTICAL ) then
c                   Set up analytical common area
                    binfiletype = DT_ANALYTICAL
                    analytical_expression = 
     1                curve_dat(i)%analytical_expression
                     analytical_xstart = curve_dat(i)%analytical_xstart
                     analytical_xend   = curve_dat(i)%analytical_xend
                     analytical_xstep  = curve_dat(i)%analytical_xstep
                     analytical_ystart = curve_dat(i)%analytical_ystart
                     analytical_yend   = curve_dat(i)%analytical_yend
                     analytical_ystep  = curve_dat(i)%analytical_ystep
                     analytical_type   = curve_dat(i)%analytical_type
                     analytical_sweep  = curve_dat(i)%sweep
                     analytical_dim    = curve_dat(i)%datafile_dimension
                     call do_eval_analytical_2d
                     read_error = .NOT. have_analytical_data
                  endif
               endif
               if(.NOT. read_error) then
                  if( binfiletype .eq. DT_GRAPES ) then
                     grapes_x_quantity = curve_dat(i)%misc_ch(1)
                     grapes_y_quantity = curve_dat(i)%datafile_quantity
                     grapes_use_dictionary=curve_dat(i)%use_grapes_dict
                     call load_grapesfile_data(read_error)
                  endif
                  if( binfiletype .eq. DT_SPICE3 ) then
                     spice3_x_quantity = curve_dat(i)%misc_ch(1)
                     spice3_y_quantity = curve_dat(i)%datafile_quantity
                     call load_spice3file_data(read_error)
                  endif
                  call lookup_menuid_from_nickname(
     1               curve_dat(i)%datafile_quantity, menu_id)
                  if( menu_id .eq. -1 ) then
                     write(luttyo,2020) curve_dat(i)%datafile_quantity
2020                 format(' ','Error: unknown quantity: ',a)
                     bad_entry = .TRUE.
                  else
                     index = menu_to_arr(menu_id)
                     if( is_data_1d(menu_id) ) then
                        call copy_data_to_slot(entry,
     1                     arr2d(1,index), ni,
     2                     xyz(1,idisti), 0.0, 'Y',
     3                     data_name(menu_id), data_nickname(menu_id),
     4                     data_units(menu_id), COMP_CUT_TYPE,
     5                     curve_dat(i)%int_limits,
     6                     curve_dat(i)%label,
     7                     curve_dat(i)%color, curve_dat(i)%linetype,
     8                     curve_dat(i)%linewidth,
     9                     curve_dat(i)%markertype,
     A                     curve_dat(i)%markercolor,
     B                     curve_dat(i)%markersize,
     C                     curve_dat(i)%xaxis_tag,
     D                     curve_dat(i)%yaxis_tag,
     E                     curve_dat(i)%suppress,
     F                     curve_dat(i)%mx, curve_dat(i)%ax,
     G                     curve_dat(i)%my, curve_dat(i)%ay,
     H                     curve_dat(i)%xabs, curve_dat(i)%yabs,
     I                     is_data_log(menu_id),
     J                     curve_dat(i)%x_expression,
     K                     curve_dat(i)%y_expression)
                     else
                        if( curve_dat(i)%int ) then
                        call extract_int_line(
     1                     arrdat(arrind(index)), ni,
     2                     does_data_jump_at_oxide(menu_id),
     3                     is_data_log(menu_id),
     4                     is_data_in_oxide(menu_id),
     5                     curve_dat(i)%cut_line_axis,
     6                     curve_dat(i)%int_limits,
     7                     tmp, TMPSIZE, 3, np, ip,
     8                     load_error)
                        if( .NOT. load_error ) then
                           call copy_data_to_slot(entry,
     1                     tmp(1,1), np, xyz(1,ip),
     2                     curve_dat(i)%cut_line_coord,
     3                     curve_dat(i)%cut_line_axis,
     4                     data_name(menu_id), data_nickname(menu_id),
     5                     data_units(menu_id), COMP_INT_TYPE,
     6                     curve_dat(i)%int_limits,
     7                     curve_dat(i)%label,
     8                     curve_dat(i)%color, curve_dat(i)%linetype,
     9                     curve_dat(i)%linewidth,
     A                     curve_dat(i)%markertype,
     B                     curve_dat(i)%markercolor,
     C                     curve_dat(i)%markersize,
     D                     curve_dat(i)%xaxis_tag,
     E                     curve_dat(i)%yaxis_tag,
     F                     curve_dat(i)%suppress,
     G                     curve_dat(i)%mx, curve_dat(i)%ax,
     H                     curve_dat(i)%my, curve_dat(i)%ay,
     I                     curve_dat(i)%xabs, curve_dat(i)%yabs,
     J                     is_data_log(menu_id),
     K                     curve_dat(i)%x_expression,
     L                     curve_dat(i)%y_expression)
                        else
                           bad_entry = .TRUE.
                        endif
                        else
                        call extract_cut_line(
     1                     arrdat(arrind(index)), ni, nj,
     2                     does_data_jump_at_oxide(menu_id),
     2                     is_data_log(menu_id),
     2                     is_data_in_oxide(menu_id),
     3                     curve_dat(i)%cut_line_coord,
     4                     curve_dat(i)%cut_line_axis,
     5                     curve_dat(i)%oxide, data_is_trigrid,
     6                     tmp, TMPSIZE, 3, np, ip,
     7                     load_error)
                        if( .NOT. load_error ) then
                           if( curve_dat(i)%oxide ) then
                              itype = COMP_INTERFACE_TYPE
                           else
                              itype = COMP_CUT_TYPE
                           endif
                           call copy_data_to_slot(entry,
     1                     tmp(1,1), np, tmp(1,2),
     2                     curve_dat(i)%cut_line_coord,
     3                     curve_dat(i)%cut_line_axis,
     4                     data_name(menu_id), data_nickname(menu_id),
     5                     data_units(menu_id), itype,
     6                     curve_dat(i)%int_limits,
     7                     curve_dat(i)%label,
     8                     curve_dat(i)%color, curve_dat(i)%linetype,
     9                     curve_dat(i)%linewidth,
     A                     curve_dat(i)%markertype,
     B                     curve_dat(i)%markercolor,
     C                     curve_dat(i)%markersize,
     D                     curve_dat(i)%xaxis_tag,
     E                     curve_dat(i)%yaxis_tag,
     F                     curve_dat(i)%suppress,
     G                     curve_dat(i)%mx, curve_dat(i)%ax,
     H                     curve_dat(i)%my, curve_dat(i)%ay,
     I                     curve_dat(i)%xabs, curve_dat(i)%yabs,
     J                     is_data_log(menu_id),
     K                     curve_dat(i)%x_expression,
     L                     curve_dat(i)%y_expression)
                        else
                           bad_entry = .TRUE.
                        endif
                        endif
                     endif
                  endif
               else
                  bad_entry = .TRUE.
               endif
            else
               write(luttyo,2050) datatype
2050           format(' ','Error: unknown data type: ',a)
               bad_entry = .TRUE.
            endif
            if( bad_entry ) then
               ierr = delete_entry(entry)
               j = j - 1
               err = .TRUE.
            endif
         endif
2000  continue
      n_cuts = j
c
      return
      end
c      
      subroutine restore_compare_plot(err,mode)
      implicit none
c
      logical err
      integer mode
c
      call pm_copy_save_axis_to_compare
c
      call compare_plot(.TRUE.,.FALSE.,mode)
c
      err = .FALSE.
c
      return
      end
c      
      subroutine pm_copy_save_axis_to_compare
      implicit none
c
      include 'sysdep.cmn'
      include 'const.h'
      include 'plot.cmn'
      include 'units.cmn'
      include 'defaults.cmn'
      include 'data.cmn'
      include 'compare.cmn'
      include 'annotations.cmn'
      include 'save.cmn'
      include 'io.cmn'
c
c     Locals
      logical lcurv, lnolabels
      integer i,j
c
c     init_common will only allow one intialization, so it's safe to
c     call anywhere.
      call init_common
c
c     Set up sticky data
      comp_lfirst = .FALSE.
c     Temporary hack: save l_curv,comp_nolabels and restore after call
c     to pm_init_comp_sticky
      lcurv = l_curv
      lnolabels = comp_nolabels
      call pm_init_comp_sticky
      l_curv = lcurv
      comp_nolabels = lnolabels
c
      j = 0
      do 100 i=1,4,3
         j = j + 1
         comp_horiz_label(j)   = axis_dat(i)%label
         comp_hunits(j)        = axis_dat(i)%units
         comp_h_min(j)         = axis_dat(i)%min
         comp_h_max(j)         = axis_dat(i)%max
         comp_h_tic_major(j)   = axis_dat(i)%major_tic
         comp_h_minor_tic_freq(j) = axis_dat(i)%minor_tic_freq
         comp_h_exponent(j)    = axis_dat(i)%exp
         comp_h_abs(j)         = axis_dat(i)%abs
         comp_h_log(j)         = axis_dat(i)%log
         comp_h_used(j)        = axis_dat(i)%used
         comp_h_autoscale_min(j)   = axis_dat(i)%autoscale_min
         comp_h_autoscale_max(j)   = axis_dat(i)%autoscale_max
c
c        Curv related parameters
         comp_h_omit(j)        = axis_dat(i)%omit
         comp_h_nonumbering(j) = axis_dat(i)%nonumbering
         comp_h_color(j)       = axis_dat(i)%color
         comp_h_minortic_size(j)    = axis_dat(i)%minor_tic_size
         comp_h_majortic_factor(j)  = axis_dat(i)%major_tic_factor
         comp_h_thickness(j)   = axis_dat(i)%thickness
         comp_h_format(j)      = axis_dat(i)%format
         comp_h_major_grid(j)  = axis_dat(i)%major_grid
         comp_h_minor_grid(j)  = axis_dat(i)%minor_grid
         comp_h_major_grid_linetype(j) = axis_dat(i)%major_grid_linetype
         comp_h_minor_grid_linetype(j) = axis_dat(i)%minor_grid_linetype
c        Note: log tic size/scale sticky variables are kept at default values
100   continue
c
      j = 0
      do 200 i=2,5,3
         j = j + 1
         comp_vert_label(j)    = axis_dat(i)%label
         comp_vunits(j)        = axis_dat(i)%units
         comp_v_min(j)         = axis_dat(i)%min
         comp_v_max(j)         = axis_dat(i)%max
         comp_v_tic_major(j)   = axis_dat(i)%major_tic
         comp_v_minor_tic_freq(j) = axis_dat(i)%minor_tic_freq
         comp_v_exponent(j)    = axis_dat(i)%exp
         comp_v_abs(j)         = axis_dat(i)%abs
         comp_v_log(j)         = axis_dat(i)%log
         comp_v_used(j)        = axis_dat(i)%used
         comp_v_autoscale_min(j)   = axis_dat(i)%autoscale_min
         comp_v_autoscale_max(j)   = axis_dat(i)%autoscale_max
c
c        Curv related parameters
c
c        Curv related parameters
         comp_v_omit(j)        = axis_dat(i)%omit
         comp_v_nonumbering(j) = axis_dat(i)%nonumbering
         comp_v_color(j)       = axis_dat(i)%color
         comp_v_minortic_size(j)    = axis_dat(i)%minor_tic_size
         comp_v_majortic_factor(j)  = axis_dat(i)%major_tic_factor
         comp_v_thickness(j)   = axis_dat(i)%thickness
         comp_v_format(j)      = axis_dat(i)%format
         comp_v_major_grid(j)  = axis_dat(i)%major_grid
         comp_v_minor_grid(j)  = axis_dat(i)%minor_grid
         comp_v_major_grid_linetype(j) = axis_dat(i)%major_grid_linetype
         comp_v_minor_grid_linetype(j) = axis_dat(i)%minor_grid_linetype
c        Note: log tic size/scale sticky variables are kept at default values
200   continue
c
      lsavefile_read = .TRUE.
c
      return
      end
c
      subroutine extract_cut_line(data_array, isize, jsize,
     1                     discontinuity,
     2                     data_is_log, value_in_oxide,
     3                     cut_line_coord, cut_line_axis,
     4                     along_oxide, data_is_trigrid,
     5                     work_vec, irows, icols, np, ip, err)
      implicit none
c
      include 'const.h'
      include 'sysdep.cmn'
      include 'mesh.cmn'
      include 'units.cmn'
      include 'io.cmn'
c
      integer isize, jsize, irows, icols, np, ip
      real  data_array(isize,jsize), work_vec(irows, icols)
      real  cut_line_coord
      logical discontinuity, along_oxide, data_is_log, value_in_oxide,
     1        data_is_trigrid, err
      character*1 cut_line_axis
c
c     Extract data values from a 2D array of data along a vertical or
c     horizontal cut line. Return values in work_vec(*,1) and coordinates in
c     work_vec(*,2)
c
      integer V_CUT, H_CUT
      parameter (V_CUT = 1, H_CUT = 2)
c
c     locals
c
      integer type, index
      integer i, j
      character*1 discont_type
c
      integer getindex  ! function
c
      err = .FALSE.
c
      if( icols .lt. 3 ) then
         write(luttyo,30)
30       format(' ','Error: (extract_cut_line) ',
     1              'insufficient temporary storage')
         err = .TRUE.
      endif
c
      if( along_oxide ) then
         if( data_is_trigrid ) then
            write(luttyo,40)
40          format(' ','Error: Cut along oxide not supported yet on ',
     1             'triangular mesh')
            err = .TRUE.
         else
            call find_values_along_oxide(data_array, isize, jsize,
     1      cut_line_coord, discontinuity, data_is_log, value_in_oxide,
     2      work_vec(1,1), irows)
            np = ni
            do 41 i=1,np
               work_vec(1,2) = xyz(i,idisti)
41          continue
         endif
         return
      endif
c
      if( data_is_trigrid ) then
         call extract_cut_line_tri(data_array, tri_npts,
     1                     discontinuity,
     2                     data_is_log, value_in_oxide,
     3                     cut_line_coord, cut_line_axis,
     4                     along_oxide,
     5                     work_vec, irows, icols, np, err)
         ip = 0
         return
      endif
c
c     Determine if a horizontal or vertical cut
      if ( cut_line_axis .eq. charv ) then
         type = H_CUT
      else if ( cut_line_axis .eq. charh ) then
         type = V_CUT
      else
         write(luttyo,50) cut_line_axis
50       format(' ','Error: Invalid cut line axis: ',a)
         err = .TRUE.
         return
      endif
c
c     Note well: In the code this was modeled after (store_line),
c     we ask the user for help if the quantity has a discontinuity
c     at the oxide and semiconductor interface as to which side s/he wants.
c     We don't have this info from the save file (yet), so always select
c     the semiconductor side.
      discont_type = 'S'
c
c     get equivalent index in array
      if( type .eq. V_CUT ) then     ! vertical cut
         index = getindex(cut_line_coord, xyz(1,idisti), ni)
      else
         index = getindex(cut_line_coord, xyz(1,idistj), nj) ! horizontal cut
      endif
      if( index .eq. 0 ) then
         write(luttyo,100) cut_line_coord*pm_scale_cm2plot
100      format(' ','Error: Cut line value ',1pe10.3,
     1              ' is outside structure')
         err = .TRUE.
         return
      endif
c
      if( type .eq. V_CUT ) then
         do 200 j=1,nj
            work_vec(j,2) = data_array(index,j)
            work_vec(j,3) = data_array(index+1,j)
200         continue
         call interp_line(work_vec(1,2), work_vec(1,3), work_vec(1,1),
     1      nj,
     1      xyz(index,idisti), xyz(index+1,idisti), cut_line_coord)
         np = nj
         ip = idistj
         do 201 i=1,np
            work_vec(i,2) = xyz(i,ip)
201      continue
      else
         if( discont_type .eq. 'O' ) then
c           special handling of interface...
c           User wants oxide cut_line_coord, take one step up into oxide
            index = index - 1
            cut_line_coord = xyz(index,idistj)
         endif
         do 300 i=1,ni
            work_vec(i,2) = data_array(i,index)
            work_vec(i,3) = data_array(i,index+1)
300         continue
         call interp_line(work_vec(1,2), work_vec(1,3), work_vec(1,1),
     1      ni,
     1      xyz(index,idistj), xyz(index+1,idistj), cut_line_coord)
         np = ni
         ip = idisti
         do 301 i=1,np
            work_vec(i,2) = xyz(i,ip)
301      continue
      endif
c
      return
      end
c
      subroutine extract_int_line(data_array, isize,
     1                     discontinuity,
     2                     data_is_log, value_in_oxide,
     3                     int_axis,
     4                     int_limits,
     5                     work_vec, irows, icols, np, ip, err)
      implicit none
c
      include 'const.h'
      include 'sysdep.cmn'
      include 'mesh.cmn'
      include 'units.cmn'
      include 'io.cmn'
c
      integer isize, irows, icols, np, ip
      real  data_array(isize,isize), work_vec(irows, icols)
      real  int_limits(2)
      logical discontinuity, along_oxide, data_is_log, value_in_oxide,
     1        err
      character*1 int_axis
c
      integer V_INT, H_INT
      parameter (V_INT = 1, H_INT = 2)
c
c     locals
c
      real value, lower, upper
      integer type, ilower, iupper
      integer i
      logical novalue_in_oxide
      integer getindex  ! function
c
      err = .FALSE.
c
c     Determine if a horizontal or vertical integration
      if ( int_axis .eq. charv ) then
         type = V_INT
      else if ( int_axis .eq. charh ) then
         type = H_INT
      else
         write(luttyo,50) int_axis
50       format(' ','Error: Invalid integration direction: ',a)
         err = .TRUE.
         return
      endif
c
c     Check integration limits for validity
      lower = int_limits(1)
      upper = int_limits(2)
      if( lower .ge. upper ) then
         write(luttyo,201)
201      format(' ',
     1      'Error: Integration minimum must be less than maximum')
         err = .TRUE.
         return
      endif
c
      novalue_in_oxide = .NOT. value_in_oxide
      if( type .eq. V_INT .and. plane_type .ne. XZ_PLANE .and.
     1    novalue_in_oxide .and. 
     1   (lower .lt. 0.0 .or. upper .lt. 0.0) ) then
         write(luttyo,202)
202      format(' ','Error: Integration in oxide (y < 0) not defined ',
     1              'for this variable')
         err = .TRUE.
         return
      endif
c
c     Find closest mesh point to the given limits
      if( type .eq. H_INT ) then
         ilower = getindex(lower, xyz(1,idisti), ni)
      else
         ilower = getindex(lower, xyz(1,idistj), nj)
      endif
      if( ilower .eq. 0 ) then
         write(luttyo,1000) pm_scale_cm2plot*lower
1000     format(' ','Error: integration limit ',1pe10.3,
     1              ' outside bounds')
         err = .TRUE.
         return
      endif
c
      if( type .eq. H_INT ) then
         iupper = getindex(upper, xyz(1,idisti), ni)
      else
         iupper = getindex(upper, xyz(1,idistj), nj)
      endif
      if( iupper .eq. 0 ) then
         write(luttyo,1000) pm_scale_cm2plot*upper
         err = .TRUE.
         return
      endif
c
c     Adjust index to get closest to specified value.
c     getindex() function always gets mesh index of left (lower)
c     side of interval.
      if( type .eq. H_INT ) then
c        Round lower limit
         if( (xyz(ilower+1,idisti) - lower) .lt.
     1       (lower - xyz(ilower,idisti)) ) then
            ilower = ilower + 1
         endif
c        Round upper limit
         if( iupper .ne. ni ) then
            if( (xyz(iupper+1,idisti) - upper) .lt.
     1          (upper - xyz(iupper,idisti)) ) then
               iupper = iupper + 1
            endif
         endif
      else
c        Round lower limit
         if( (xyz(ilower+1,idistj) - lower) .lt.
     1       (lower - xyz(ilower,idistj)) ) then
            ilower = ilower + 1
         endif
c
c        Round upper limit
         if( iupper .ne. nj ) then
            if( (xyz(iupper+1,idistj) - upper) .lt.
     1          (upper - xyz(iupper,idistj)) ) then
               iupper = iupper + 1
            endif
         endif
      endif
c
      if( type .eq. V_INT ) then
c        Integrate sections
         do 300 i=1,ni
            call do_1d_integrate(data_array, ni, nj,
     1         lower, upper, ilower, iupper, i,
     2         novalue_in_oxide, 'V', value)
            work_vec(i,1) = value
300      continue
         np = ni
         ip = idisti
      else
         do 350 i=1,nj
            call do_1d_integrate(data_array, ni, nj,
     1         lower, upper, ilower, iupper, i,
     2         novalue_in_oxide, 'H', value)
            work_vec(i,1) = value
350      continue
         np = nj
         ip = idistj
      endif
c
      return
      end
c
      logical function update_datafile(i)
      implicit none
c
      include 'const.h'
      include 'sysdep.cmn'
      include 'io.cmn'
      include 'data.cmn'
      include 'mesh.cmn'
      include 'annotations.cmn'
      include 'compare.cmn'
      include 'save.cmn'
c
      integer i
      logical ltmp
      character*40 ch_tmp
c
c     Check a number of conditions to see if we need to read
c     the next datafile.
c
      ltmp = .FALSE.
c     If no data file loaded
      if( .NOT. datafile_loaded ) ltmp = .TRUE.
c
c     File names don't match
      if( binfilename .ne. curve_dat(i)%datafile ) ltmp = .TRUE.
c
c     File record numbers don't match
      if( datafile_record .ne. curve_dat(i)%datafile_record )
     1   ltmp = .TRUE.
c
c     File types don't match
      call datatype_to_name(binfiletype,ch_tmp)
      if( ch_tmp .ne. curve_dat(i)%datafile_type ) ltmp = .TRUE.
c
c     3D plane types and coord don't match (MINIMOS and 2DOP)
      if( datafile_loaded .and. three_d .and.
     1   (binfiletype .eq. DT_MINIMOS .or. binfiletype .eq. DT_2DOP) )
     2   then
         if( curve_dat(i)%cut_plane_axis .eq. 'Z' .and.
     1       plane_type .ne. XY_PLANE ) ltmp = .TRUE.
         if( curve_dat(i)%cut_plane_axis .eq. 'X' .and.
     1       plane_type .ne. YZ_PLANE ) ltmp = .TRUE.
         if( curve_dat(i)%cut_plane_axis .eq. 'Y' .and.
     1       plane_type .ne. XZ_PLANE ) ltmp = .TRUE.
c        Compare cut plane coordinates. Since this comparison may be
c        sensitive to round off, use a relative error check to within
c        1 angstrom (coordinates are in microns).
         if( ABS(curve_dat(i)%cut_plane_coord - plane_actual_value)
     1       .ge. 0.0001 )
     1       ltmp = .TRUE.
      endif
c
c     If 2DOP file, check gate length and interface file
      if( curve_dat(i)%datafile_type .eq. SAVE_DATATYPE_2DOP .and.
     1    (tdop_gate_length .ne. curve_dat(i)%misc(1) .or.
     2     tdop_interface_file .ne. curve_dat(i)%misc_ch(1)) )
     3    ltmp = .TRUE.
c
c     If MEDICI file, check solution file
      if( curve_dat(i)%datafile_type .eq. SAVE_DATATYPE_MEDICI .and.
     1    medici_solution_file .ne. curve_dat(i)%misc_ch(1) )
     2    ltmp = .TRUE.
c
c     If PISCES file, check solution file
      if( curve_dat(i)%datafile_type .eq. SAVE_DATATYPE_PISCES .and.
     1    medici_solution_file .ne. curve_dat(i)%misc_ch(1) )
     2    ltmp = .TRUE.
c
c     If VLSICAP, check misc parameters
      if( curve_dat(i)%datafile_type .eq. SAVE_DATATYPE_VLSICAP .and.
     1    (vlsi_edge .ne. curve_dat(i)%misc(1) .or.
     1     vlsi_yint .ne. curve_dat(i)%misc(2)) ) ltmp = .TRUE.
c
c     If TIF file, check for 1d data file
      if( tif_1d_data ) ltmp = .TRUE.
c
c     Always 'update' analytical datatype
      if( curve_dat(i)%datafile_type .eq. SAVE_DATATYPE_ANALYTICAL )
     1    ltmp = .TRUE.
c
      update_datafile = ltmp
      return
      end
c
      subroutine load_datafile_from_save(err)
      implicit none
c
      include 'sysdep.cmn'
      include 'const.h'
      include 'io.cmn'
      include 'plot.cmn'
      include 'mesh.cmn'
      include 'data.cmn'
      include 'units.cmn'
      include 'defaults.cmn'
      include 'compare.cmn'
      include 'user.cmn'
      include 'analytical.cmn'
      include 'annotations.cmn'
      include 'save.cmn'
c
      logical err
c
      integer menu_id, i
      logical read_error, update_datafile
      character*40 datatype
c
c     Set up data area and read file
      datatype = curve_dat(1)%datafile_type
      if( datatype .eq. SAVE_DATATYPE_MINIMOS  .or.
     1    datatype .eq. SAVE_DATATYPE_PROMIS   .or.
     2    datatype .eq. SAVE_DATATYPE_VLSICAP  .or.
     3    datatype .eq. SAVE_DATATYPE_BAMBI    .or.
     4    datatype .eq. SAVE_DATATYPE_USEOUT   .or.
     5    datatype .eq. SAVE_DATATYPE_2DOP     .or.
     6    datatype .eq. SAVE_DATATYPE_SUPREM3S .or.
     7    datatype .eq. SAVE_DATATYPE_SUPREM3P .or.
     8    datatype .eq. SAVE_DATATYPE_SUPREM4  .or.
     8    datatype .eq. SAVE_DATATYPE_PISCES   .or.
     9    datatype .eq. SAVE_DATATYPE_CALLABLE .or.
     A    datatype .eq. SAVE_DATATYPE_TIF      .or.
     B    datatype .eq. SAVE_DATATYPE_FAIM     .or.
     C    datatype .eq. SAVE_DATATYPE_MCP2D    .or.
     D    datatype .eq. SAVE_DATATYPE_MEDICI   .or.
     F    datatype .eq. SAVE_DATATYPE_ASCII    .or.
     G    datatype .eq. SAVE_DATATYPE_GRAPES .or.
     H    datatype .eq. SAVE_DATATYPE_SPICE3 .or.
     I    datatype .eq. SAVE_DATATYPE_ANALYTICAL ) then
c        Optimization: if file is loaded, and name is same, don't re-read
         if( update_datafile(1) ) then
            prev_binfilename = binfilename
            binfilename     = curve_dat(1)%datafile
            plane_value     = curve_dat(1)%cut_plane_coord
            cut_plane_axis  = curve_dat(1)%cut_plane_axis
            datafile_record = curve_dat(1)%datafile_record
            binfile_convert_type = curve_dat(1)%bin_datatype
            if( curve_dat(1)%datafile_dimension .eq. 3 ) then
               three_d = .TRUE.
            else
               three_d = .FALSE.
            endif
            prev_binfiletype = binfiletype
            if( datatype .eq. SAVE_DATATYPE_MINIMOS     ) then
               binfiletype = DT_MINIMOS
               call read_mini(.TRUE.,read_error)
            else if( datatype .eq. SAVE_DATATYPE_PROMIS ) then
               binfiletype = DT_PROMIS
               call read_promis(.TRUE.,read_error)
            else if( datatype .eq. SAVE_DATATYPE_SUPREM4 ) then
               binfiletype = DT_SUPREM4
               call read_sup4(.TRUE.,read_error)
            else if( datatype .eq. SAVE_DATATYPE_PISCES ) then
               binfiletype = DT_PISCES
               medici_solution_file = curve_dat(1)%misc_ch(1)
               call read_pisces(.TRUE.,read_error)
            else if( datatype .eq. SAVE_DATATYPE_MEDICI   ) then
               binfiletype = DT_MEDICI
               medici_solution_file = curve_dat(1)%misc_ch(1)
               call read_medici(.TRUE.,read_error)
            else if( datatype .eq. SAVE_DATATYPE_VLSICAP) then
               binfiletype = DT_VLSICAP
               vlsi_edge = curve_dat(1)%misc(1)
               vlsi_yint = curve_dat(1)%misc(2)
               call read_vlsicap(.TRUE.,read_error)
            else if( datatype .eq. SAVE_DATATYPE_BAMBI  ) then
               binfiletype = DT_BAMBI
               call read_bambi(read_error)
            else if( datatype .eq. SAVE_DATATYPE_USEOUT ) then
               binfiletype = DT_USEOUT
               call read_useout(.TRUE.,read_error)
            else if( datatype .eq. SAVE_DATATYPE_2DOP   ) then
               binfiletype = DT_2DOP
               tdop_gate_length = curve_dat(1)%misc(1)
               tdop_interface_file = curve_dat(1)%misc_ch(1)
               call read_2dop(.TRUE.,read_error)
            else if( datatype .eq. SAVE_DATATYPE_SUPREM3S) then
               binfiletype = DT_SUPREM3S
               call read_suprem3(2,read_error)
            else if( datatype .eq. SAVE_DATATYPE_SUPREM3P) then
               binfiletype = DT_SUPREM3P
               call read_suprem3(1,read_error)
            else if( datatype .eq. SAVE_DATATYPE_TIF ) then
               tif_ixcol = curve_dat(1)%colx
               tif_iycol = curve_dat(1)%coly(1)
               binfiletype = DT_TIF
               call read_tif(.TRUE.,read_error)
            else if( datatype .eq. SAVE_DATATYPE_GRAPES ) then
               binfiletype = DT_GRAPES
               call read_grapesfile(.TRUE.,read_error)
            else if( datatype .eq. SAVE_DATATYPE_SPICE3 ) then
               binfiletype = DT_SPICE3
               call do_read_spice3file(.TRUE.,read_error)
            else if( datatype .eq. SAVE_DATATYPE_FAIM ) then
               binfiletype = DT_FAIM
               call read_faim(.TRUE.,read_error)
            else if( datatype .eq. SAVE_DATATYPE_MCP2D ) then
               binfiletype = DT_MCP2D
               call read_mcp2d(.TRUE.,read_error)
            else if( datatype .eq. SAVE_DATATYPE_CALLABLE ) then
               binfiletype = DT_CALLABLE
               call load_callable_data(.TRUE.,read_error)
            else if( datatype .eq. SAVE_DATATYPE_ASCII ) then
               binfiletype = DT_ASCII
c              Set up user area
               user_1d_asciifile = curve_dat(1)%datafile
               user_1d_col_x = curve_dat(1)%colx
               user_1d_n_ycols = curve_dat(1)%n_ycols
               do i=1,user_1d_n_ycols
                  user_1d_col_y(i) = curve_dat(1)%coly(i)
               end do
               user_1d_col_z = curve_dat(1)%colz
               user_1d_skip  = curve_dat(1)%skip
               user_1d_nlines = curve_dat(1)%lines_to_read
               user_1d_ax = curve_dat(1)%ax
               user_1d_ay = curve_dat(1)%ay
               user_1d_az = curve_dat(1)%az
               user_1d_mx = curve_dat(1)%mx
               user_1d_my = curve_dat(1)%my
               user_1d_mz = curve_dat(1)%mz
               user_1d_x_expression = curve_dat(1)%x_expression
               user_1d_y_expression = curve_dat(1)%y_expression
               user_1d_xabs = curve_dat(1)%xabs
               user_1d_yabs = curve_dat(1)%yabs
               user_1d_zabs = curve_dat(1)%zabs
               user_1d_sort = curve_dat(1)%sort
               user_1d_x_char_data = curve_dat(1)%x_char_data
               user_1d_sample_freq = curve_dat(1)%sample_freq
               user_1d_rect_grid = .NOT. curve_dat(1)%data_is_trigrid
               user_1d_sweep = curve_dat(1)%sweep
               user_1d_merge_ycols = curve_dat(1)%merge_ycols
               user_1d_contour = curve_dat(1)%contour
               user_1d_errorbars = curve_dat(1)%errorbars
               user_1d_dim  = curve_dat(1)%datafile_dimension
               user_1d_dataset_number  = curve_dat(1)%datafile_record
               if( user_1d_dim .eq. 1 ) then
                  call do_read_asciifile
               elseif ( user_1d_rect_grid ) then
                  call do_read_asciifile_2d
               else
                  call do_read_asciifile_2d_tri
               endif
               read_error = .NOT. have_user_1d_data
            else if( datatype .eq. SAVE_DATATYPE_ANALYTICAL ) then
c              Set up analytical common area
               binfiletype = DT_ANALYTICAL
               analytical_expression = 
     1             curve_dat(1)%analytical_expression
               analytical_xstart = curve_dat(1)%analytical_xstart
               analytical_xend   = curve_dat(1)%analytical_xend
               analytical_xstep  = curve_dat(1)%analytical_xstep
               analytical_pstart = curve_dat(1)%analytical_pstart
               analytical_pend   = curve_dat(1)%analytical_pend
               analytical_pstep  = curve_dat(1)%analytical_pstep
               analytical_ystart = curve_dat(1)%analytical_ystart
               analytical_yend   = curve_dat(1)%analytical_yend
               analytical_ystep  = curve_dat(1)%analytical_ystep
               analytical_type   = curve_dat(1)%analytical_type
               analytical_sweep  = curve_dat(1)%sweep
               analytical_dim    = curve_dat(1)%datafile_dimension
               if( analytical_dim .eq. 1 ) then
                  call do_eval_analytical
               else
                  call do_eval_analytical_2d
               endif
               read_error = .NOT. have_analytical_data
            endif
         endif
         if(.NOT. read_error) then
            if( binfiletype .eq. DT_GRAPES ) then
               grapes_x_quantity = curve_dat(1)%misc_ch(1)
               grapes_y_quantity = curve_dat(1)%datafile_quantity
               grapes_use_dictionary = curve_dat(1)%use_grapes_dict
               call load_grapesfile_data(read_error)
            endif
            if( binfiletype .eq. DT_SPICE3 ) then
               spice3_x_quantity = curve_dat(1)%misc_ch(1)
               spice3_y_quantity = curve_dat(1)%datafile_quantity
               call load_spice3file_data(read_error)
            endif
            call lookup_menuid_from_nickname(
     1         curve_dat(1)%datafile_quantity, menu_id)
            if( menu_id .eq. -1 ) then
               write(luttyo,2020) curve_dat(1)%datafile_quantity
2020           format(' ','Error: unknown quantity: ',a)
               err = .TRUE.
            else
            endif
         else
            err = .TRUE.
         endif
      else
         write(luttyo,2040) datatype
2040     format(' ','Error: unknown data type: ',a)
         err = .TRUE.
      endif
c
      return
      end
c      
      subroutine restore_2d_plot(err)
      implicit none
c
      include 'sysdep.cmn'
      include 'const.h'
      include 'plot.cmn'
      include 'units.cmn'
      include 'defaults.cmn'
      include 'mesh.cmn'
      include 'data.cmn'
      include 'compare.cmn'
      include 'annotations.cmn'
      include 'save.cmn'
      include 'io.cmn'
      include 'title.cmn'
c
      logical err
c
      integer i, menu_id
c
c     load_datafile_from_save has already checked for valid quantity
      call lookup_menuid_from_nickname(
     1         curve_dat(1)%datafile_quantity, menu_id)
      if( is_data_1d(menu_id) ) then
         write(luttyo,50)
50       format(' ','Error: Cannot plot contours of a 1d quantity!')
         err = .TRUE.
         return
      endif
c
c     Check save quantities for validity in this context
      if( net_dop_array_index .eq. 0 .and. contour_junctions ) then
         write(luttyo,51)
51       format(' ','Error: Cannot plot junctions, no doping!')
         err = .TRUE.
         return
      endif
c
      if( contour_color .and. contour_clabel ) then
         write(luttyo,52)
52       format(' ',
     1          'Error: Cannot plot labels on color fill contour plot!')
         err = .TRUE.
         return
      endif
c
c     Set up sticky information
      twod_lfirst = .FALSE.
      prev_plane_type = plane_type
      twod_h_min  = pm_scale_plot2cm*axis_dat(1)%min
      twod_h_max  = pm_scale_plot2cm*axis_dat(1)%max
      twod_h_tic_major = pm_scale_plot2cm*axis_dat(1)%major_tic
      twod_h_minor_tic_freq = axis_dat(1)%minor_tic_freq
      twod_h_label = axis_dat(1)%label
      twod_h_units = axis_dat(1)%units
      twod_h_autoscale = axis_dat(1)%autoscale_min .or.
     1                   axis_dat(1)%autoscale_max
c
      twod_v_min  = pm_scale_plot2cm*axis_dat(2)%min
      twod_v_max  = pm_scale_plot2cm*axis_dat(2)%max
      twod_v_tic_major = pm_scale_plot2cm*axis_dat(2)%major_tic
      twod_v_minor_tic_freq = axis_dat(2)%minor_tic_freq
      twod_v_label = axis_dat(2)%label
      twod_v_units = axis_dat(2)%units
      twod_v_autoscale = axis_dat(2)%autoscale_min .or.
     1                   axis_dat(2)%autoscale_max
c
      twod_color     = curve_dat(1)%contour_color
      twod_junctions = curve_dat(1)%contour_junctions
      twod_grid      = curve_dat(1)%contour_grid
      twod_clabel    = curve_dat(1)%contour_clabel
      twod_negate    = contour_negate
      twod_invert    = .NOT. axis_dat(2)%invert
      twod_solid     = curve_dat(1)%contour_solid
      twod_log_interp = curve_dat(1)%contour_log
      twod_expression = curve_dat(1)%z_expression
      twod_colormat  = contour_colormat
      twod_num_contours = curve_dat(1)%num_contours
      do 100 i=1,twod_num_contours
         twod_contours(i) = curve_dat(1)%contours(i)
100   continue
c
c     Set up plot title
      plot_title = title_dat(1)
c
      lsavefile_read = .TRUE.
c
      call ship_contours(arrdat(arrind(menu_to_arr(menu_id))),
     1        ni, nj, data_name(menu_id),
     2        data_units(menu_id), data_nickname(menu_id),
     3        is_data_log(menu_id), is_data_in_oxide(menu_id),
     4        .TRUE.)
c
c      return
      end
c
      subroutine restore_3d_plot(err)
      implicit none
c
      include 'const.h'
      include 'sysdep.cmn'
      include 'plot.cmn'
      include 'units.cmn'
      include 'defaults.cmn'
      include 'mesh.cmn'
      include 'data.cmn'
      include 'compare.cmn'
      include 'annotations.cmn'
      include 'save.cmn'
      include 'io.cmn'
      include 'title.cmn'
c
      logical err
c
      integer menu_id
c
c     load_datafile_from_save has already checked for valid quantity
      call lookup_menuid_from_nickname(
     1         curve_dat(1)%datafile_quantity, menu_id)
      if( is_data_1d(menu_id) ) then
         write(luttyo,50)
50       format(' ','Error: Cannot plot surface of a 1d quantity!')
         err = .TRUE.
         return
      endif
c
c     Set sticky information
      threed_lfirst = .FALSE.
      prev_plane_type = plane_type
      threed_h_min       = pm_scale_plot2cm*axis_dat(1)%min
      threed_h_max       = pm_scale_plot2cm*axis_dat(1)%max
      threed_h_tic_major = pm_scale_plot2cm*axis_dat(1)%major_tic
      threed_h_minor_tic_freq  = axis_dat(1)%minor_tic_freq
      threed_h_label     = axis_dat(1)%label
      threed_h_units     = axis_dat(1)%units
      threed_h_autoscale = axis_dat(1)%autoscale_min .or.
     1                     axis_dat(1)%autoscale_max
c
      threed_v_min        = pm_scale_plot2cm*axis_dat(2)%min
      threed_v_max        = pm_scale_plot2cm*axis_dat(2)%max
      threed_v_tic_major  = pm_scale_plot2cm*axis_dat(2)%major_tic
      threed_v_minor_tic_freq  = axis_dat(2)%minor_tic_freq
      threed_v_label      = axis_dat(2)%label
      threed_v_units      = axis_dat(2)%units
      threed_v_autoscale  = axis_dat(2)%autoscale_min .or.
     1                      axis_dat(2)%autoscale_max
c
      threed_z_min        = axis_dat(3)%min
      threed_z_max        = axis_dat(3)%max
      threed_z_tic_major  = axis_dat(3)%major_tic
      threed_z_minor_tic_freq  = axis_dat(3)%minor_tic_freq
      threed_z_exponent   = axis_dat(3)%exp
      threed_z_abs        = axis_dat(3)%abs
      threed_z_log        = axis_dat(3)%log
      threed_z_autoscale  = axis_dat(3)%autoscale_min .or.
     1                      axis_dat(3)%autoscale_max
      threed_z_label      = axis_dat(3)%label
      threed_z_units      = axis_dat(3)%units
c
      threed_eye_position(1) = surface_eye_position(1)
      threed_eye_position(2) = surface_eye_position(2)
      threed_eye_position(3) = surface_eye_position(3)
c
      threed_swapxy = surface_swapxy
      threed_expression = surface_expression
c
c     Set up plot title
      plot_title = title_dat(1)
c
c     End sticky information
c
      lsavefile_read = .TRUE.
c
      call ship_surface(arrdat(arrind(menu_to_arr(menu_id))), ni, nj,
     1                  data_name(menu_id), data_units(menu_id),
     2                  data_nickname(menu_id),
     3                  is_data_log(menu_id),.TRUE.)
c
      return
      end
c      
      subroutine restore_1d_plot(err)
      implicit none
c
      include 'const.h'
      include 'sysdep.cmn'
      include 'plot.cmn'
      include 'units.cmn'
      include 'defaults.cmn'
      include 'mesh.cmn'
      include 'data.cmn'
      include 'compare.cmn'
      include 'annotations.cmn'
      include 'save.cmn'
      include 'io.cmn'
      include 'title.cmn'
c
      logical err
c
      integer i, menu_id, index, ip, np
      logical load_error, luser
c
      real tmp(COMP_ARRAYSIZE), tmp2(COMP_ARRAYSIZE), clc
      character*1 cla
      character*20 dfq
c
c     load_datafile_from_save has already check for valid quantity
      call lookup_menuid_from_nickname(
     1         curve_dat(1)%datafile_quantity, menu_id)
      index = menu_to_arr(menu_id)
c
      if( is_data_1d(menu_id) ) then
c        Load into temporary arrays because plot_line re-uses work_vec
         do 100 i=1,ni
            tmp(i)  = arr2d(i,index)
            tmp2(i) = xyz(i,idisti)
100      continue
         np = ni
         ip = idisti
         load_error = .FALSE.
      else
         call extract_cut_line(arrdat(arrind(index)), ni, nj,
     1        does_data_jump_at_oxide(menu_id),
     2        is_data_log(menu_id),
     3        is_data_in_oxide(menu_id),
     4        curve_dat(1)%cut_line_coord,
     5        curve_dat(1)%cut_line_axis,
     6        curve_dat(1)%oxide, data_is_trigrid,
     7        work_vec, COMP_ARRAYSIZE, MULTI_SIZE, np, ip,
     8        load_error)
c        Load into temporary arrays because plot_line re-uses work_vec
         do 200 i=1,np
            tmp(i) = work_vec(i,1)
            tmp2(i) = work_vec(i,2)
200      continue
      endif
c
      if( load_error ) then
         err = .TRUE.
         return
      endif
c
c     Set up 1d sticky information
      oned_lfirst = .FALSE.
      prev_plane_type = plane_type
      oned_h_label = axis_dat(1)%label
      oned_h_units = axis_dat(1)%units
      oned_h_min  = axis_dat(1)%min
      oned_h_max  = axis_dat(1)%max
      oned_h_tic_major = axis_dat(1)%major_tic
      oned_h_minor_tic_freq = axis_dat(1)%minor_tic_freq
      oned_h_exponent = axis_dat(1)%exp
      oned_h_log = axis_dat(1)%log
      oned_h_abs = axis_dat(1)%abs
      oned_h_autoscale = axis_dat(1)%autoscale_min .or.
     1                   axis_dat(1)%autoscale_max
c
      oned_v_min  = axis_dat(2)%min
      oned_v_max  = axis_dat(2)%max
      oned_v_tic_major = axis_dat(2)%major_tic
      oned_v_minor_tic_freq = axis_dat(2)%minor_tic_freq
      oned_v_exponent = axis_dat(2)%exp
      oned_v_log = axis_dat(2)%log
      oned_v_abs = axis_dat(2)%abs
      oned_v_autoscale = axis_dat(2)%autoscale_min .or.
     1                   axis_dat(2)%autoscale_max
      if( curve_dat(1)%markercolor .ne. 0 .and.
     1    curve_dat(1)%markertype  .ne. 0 ) then
         oned_do_markers = .TRUE.
      else
         oned_do_markers = .FALSE.
      endif
c
c     Set up plot title
      plot_title = title_dat(1)
c
      lsavefile_read = .TRUE.
c
c     Store these variables in temporaries, because the 1D plot routine
c     will initialize curve_dat(*)%*. Axis_dat quatities are stored
c     locally in plot_line.
      clc = curve_dat(1)%cut_line_coord
      cla = curve_dat(1)%cut_line_axis
      dfq = curve_dat(1)%datafile_quantity
      luser = dfq .eq. 'ASCII'
c
      call plot_line(tmp, tmp2, np,
     1     clc, cla, dfq,
     2     axis_dat(2)%label, axis_dat(2)%units, is_data_log(menu_id),
     3     .TRUE., luser)
c
      return
      end
c
      subroutine check_save_plot_params(err)
      implicit none
c
c     Check plot parameters for validity
c
c
      include 'const.h'
      include 'sysdep.cmn'
      include 'plot.cmn'
      include 'units.cmn'
      include 'defaults.cmn'
      include 'mesh.cmn'
      include 'data.cmn'
      include 'compare.cmn'
      include 'annotations.cmn'
      include 'save.cmn'
      include 'io.cmn'
c
      logical err
c
c     locals
      real rdum1, rdum2
      integer i, iscale
      character*1 axes(5)
      data axes/'X','Y','Z','X','Y'/
c
c     Check for min less than max
      do 100 i=1,MAX_AXES
         if( axis_dat(i)%used .and.
     1       axis_dat(i)%max .le. axis_dat(i)%min .and.
     2       .NOT. (axis_dat(i)%autoscale_min .or.
     3              axis_dat(i)%autoscale_max) ) then
            write(luttyo,101) axes(i)
101         format(' ','Error: ',a1,' axis max is less than min')
            err = .TRUE.
         endif
100   continue
c
c     Check for missing axes
      if( .not. axis_dat(1)%used .and. .not. axis_dat(4)%used ) then
         write(luttyo,201) 'X'
201      format(' ','Error: ',a1,' axis was not specified')
         err = .TRUE.
      endif
      if( .not. axis_dat(2)%used .and. .not. axis_dat(5)%used ) then
         write(luttyo,201) 'Y'
         err = .TRUE.
      endif
      if( save_plot_type .eq. PLOT_3D_COMMAND ) then
         if( .not. axis_dat(3)%used ) then
            write(luttyo,201) 'Z'
            err = .TRUE.
          endif
      endif
c
c     Check for min/max less than zero on log scale
      do 300 i=1,MAX_AXES
         if( axis_dat(i)%used .and. axis_dat(i)%log .and.
     1       .NOT. (axis_dat(i)%autoscale_min .or.
     2              axis_dat(i)%autoscale_max) ) then
            if (axis_dat(i)%min .le. 0.0 ) then
               write(luttyo,301) axes(i)
301            format(' ','Error: Log scale specified for ',a1,' axis'/
     1             ' ','       and min is less than or equal to zero')
               err = .TRUE.
            endif
            if (axis_dat(i)%max .le. 0.0 ) then
               write(luttyo,302) axes(i)
302            format(' ','Error: Log scale specified for ',a1,' axis'/
     1             ' ','       and max is less than or equal to zero')
               err = .TRUE.
            endif
         endif
300   continue
c
c     Check for invalid tic mark parameters, try to set them to
c     reasonable values
      do 400 i=1,MAX_AXES
         if( axis_dat(i)%used .and. .NOT. axis_dat(i)%log .and. .NOT.
     1      (axis_dat(i)%autoscale_min .or. axis_dat(i)%autoscale_max)
     2     ) then
c           
            if( axis_dat(i)%minor_tic_freq .le. 0 .or.
     1          axis_dat(i)%major_tic .le. 0.0    ) then
               write(luttyo,401) axes(i)
401            format(' ','Warning: ',a1,
     1                ' axis tic mark parameters are invalid, ',
     2                'using default values')
               call get_default_limits(axis_dat(i)%min,axis_dat(i)%max,
     1              rdum1, rdum2, axis_dat(i)%major_tic,
     2              axis_dat(i)%minor_tic_freq, iscale, .FALSE.)
               if( iscale .ne. 0 ) then
                  axis_dat(i)%major_tic = axis_dat(i)%major_tic*
     1                                    10.0**iscale
               endif
            endif
         endif
400   continue
c
c     In comparison plots, check for data file attached to an axis
c     that has not been defined yet...
      if( save_plot_type .eq. PLOT_COMPARE_COMMAND ) then
         do i=1,num_curves
c           X bottom
            if( curve_dat(i)%xaxis_tag .eq. 1 .and.
     1          .not. axis_dat(1)%used ) then
               write(luttyo,501) i
501            format(' ','Error: Lower X axis for curve number ',i5,
     1                ' was never defined')
               err = .TRUE.
            endif
c           X top
            if( curve_dat(i)%xaxis_tag .eq. 2 .and.
     1          .not. axis_dat(4)%used ) then
               write(luttyo,502) i
502            format(' ','Error: Upper X axis for curve number ',i5,
     1                ' was never defined')
               err = .TRUE.
            endif
c           Y left
            if( curve_dat(i)%yaxis_tag .eq. 1 .and.
     1          .not. axis_dat(2)%used ) then
               write(luttyo,503) i
503            format(' ','Error: Left Y axis for curve number ',i5,
     1                ' was never defined')
               err = .TRUE.
            endif
c           Y right
            if( curve_dat(i)%yaxis_tag .eq. 2 .and.
     1          .not. axis_dat(5)%used ) then
               write(luttyo,504) i
504            format(' ','Error: Right Y axis for curve number ',i5,
     1                ' was never defined')
               err = .TRUE.
            endif
         end do
      endif
c
c     Some sanity checks for ASCII file data
c     If user specifies loading multiple columns, then make sure that
c     datafile dimension is one...
      do i=1,num_curves
         if( curve_dat(i)%datafile_type .eq. SAVE_DATATYPE_ASCII ) then
            if( curve_dat(i)%datafile_dimension .gt. 2 ) then
c              Second test is our quirky way of telling Postmini
c              to load all remaining Y columns...
               if( curve_dat(i)%n_ycols .gt. 1 .or.
     1            (curve_dat(i)%n_ycols .eq. 1. and.
     2             curve_dat(i)%coly(1) .eq. -1) ) then
                  write(luttyo,505) i
505               format(' ','Error: Curve number ',i5,
     1               ': cannot load multiple Y columns for X-Y-Z data')
                  err = .TRUE.
               endif
            endif
         endif
      end do
c
c     Some general sanity checks...
      do i=1,num_curves
         if( curve_dat(i)%datafile_record .le. 0 ) then
            write(luttyo,600) i
600         format(' ','Error: Curve number ',i5,
     1      ': datafile record must be greater than zero')
            err = .TRUE.
         endif
      end do
c
      return
      end
c
      subroutine extract_cut_line_tri(data_array, isize,
     1                     discontinuity,
     2                     data_is_log, value_in_oxide,
     3                     cut_line_coord, cut_line_axis,
     4                     along_oxide,
     5                     work_vec, irows, icols, np, err)
      implicit none
c
      include 'const.h'
      include 'sysdep.cmn'
      include 'mesh.cmn'
      include 'units.cmn'
      include 'io.cmn'
c
      integer isize, irows, icols, np
      real  data_array(isize), work_vec(irows, icols)
      real  cut_line_coord
      logical discontinuity, along_oxide, data_is_log, value_in_oxide,
     1        err
      character*1 cut_line_axis
c
c     locals
      integer i, ie, type, inext(3), itmp
      integer kreg(4*COMP_ARRAYSIZE), ireg
      real x(3), y(3), z(3), s, dz, pt, t1, t2, tmp
      data inext/2,3,1/
c
      integer V_CUT, H_CUT
      parameter (V_CUT = 1, H_CUT = 2)
c
      err = .FALSE.

c     Check for possible local data overrun
      if( irows .gt. 4*COMP_ARRAYSIZE ) then
         write(luttyo,30)
30       format(' ','Error: (extract_cut_line_tri) ',
     1              'insufficient temporary storage')
         np = 0
         err = .TRUE.
         return
      endif
c
c     Extract data along a vertical or horizontal line from a 2D
c     dataset defined on a triangular grid and return data
c     in work_vec(*,1) and coordinates in work_vec(*,2)
c
c     Determine if a horizontal or vertical cut
      if ( cut_line_axis .eq. charv ) then
         type = H_CUT
      else if ( cut_line_axis .eq. charh ) then
         type = V_CUT
      else
         write(luttyo,50) cut_line_axis
50       format(' ','Error: Invalid cut line axis: ',a)
         err = .TRUE.
         return
      endif
c
c     Is cut line outside bounding box?
      if( type .eq. H_CUT ) then
         if( cut_line_coord .lt. xyzlimits(1,2) .or.
     1       cut_line_coord .gt. xyzlimits(2,2) ) then
            write(luttyo,40) cut_line_coord*pm_scale_cm2plot
40          format(' ','Error: Cut line value ',1pe10.3,
     1                 ' is outside structure')
            err = .TRUE.
            return
         endif
      else
         if( cut_line_coord .lt. xyzlimits(1,1) .or.
     1       cut_line_coord .gt. xyzlimits(2,1) ) then
            write(luttyo,40) cut_line_coord*pm_scale_cm2plot
            err = .TRUE.
            return
         endif
      endif
c
      np = 0
c
      do 100 ie=1,tri_ne
         do 101 i=1,3
            x(i) = cord(1,nop(i,ie))
            y(i) = cord(2,nop(i,ie))
101      continue
c        Check for element outside cut line
         if( type .eq. H_CUT ) then
            if( cut_line_coord .lt. MIN(y(1),y(2),y(3)) .or.
     1          cut_line_coord .gt. MAX(y(1),y(2),y(3)) ) goto 100
         else
            if( cut_line_coord .lt. MIN(x(1),x(2),x(3)) .or.
     1          cut_line_coord .gt. MAX(x(1),x(2),x(3)) ) goto 100
         endif
c        Cut line intersects element. Check each side for intersection
         do 102 i=1,3
            z(i) = data_array(nop(i,ie))
102      continue
         ireg = tri2reg(ie)
         do 103 i=1,3
            if( type .eq. H_CUT ) then
c              Horizontal cut
               dz = y(inext(i)) - y(i)
               if( dz .eq. 0.0 ) then
                  if( cut_line_coord .eq. y(i) ) then
c                    Add both endpoints
CC                     np = np + 1
CC                     if( np .gt. irows ) goto 1000
CC                     work_vec(np,1) = z(i)
CC                     work_vec(np,2) = x(i)
CC                     kreg(np) = ireg
CC                     np = np + 1
CC                     if( np .gt. irows ) goto 1000
CC                     work_vec(np,1) = z(inext(i))
CC                     work_vec(np,2) = x(inext(i))
CC                     kreg(np) = ireg
                  endif
               else
                  s = (x(inext(i)) - x(i))/dz
                  pt = s*(cut_line_coord - y(i)) + x(i)
                  t1 = MIN(x(i),x(inext(i)))
                  t2 = MAX(x(i),x(inext(i)))
                  if( pt .ge. t1 .and. pt .le. t2 ) then
                     np = np + 1
                     if( np .gt. irows ) goto 1000
                     work_vec(np,1) = (cut_line_coord - y(i))/dz*
     1                  (z(inext(i)) - z(i)) + z(i)
                     work_vec(np,2) = pt
                     kreg(np) = ireg
                  endif
               endif
            else
c              Vertical cut
               dz = x(inext(i)) - x(i)
               if( dz .eq. 0.0 ) then
                  if( cut_line_coord .eq. x(i) ) then
c                    Add both endpoints
CC                     np = np + 1
CC                     if( np .gt. irows ) goto 1000
CC                     work_vec(np,1) = z(i)
CC                     work_vec(np,2) = y(i)
CC                     kreg(np) = ireg
CC                     np = np + 1
CC                     if( np .gt. irows ) goto 1000
CC                     work_vec(np,1) = z(inext(i))
CC                     work_vec(np,2) = y(inext(i))
CC                     kreg(np) = ireg
                  endif
               else
                  s = (y(inext(i)) - y(i))/dz
                  pt = s*(cut_line_coord - x(i)) + y(i)
                  t1 = MIN(y(i),y(inext(i)))
                  t2 = MAX(y(i),y(inext(i)))
                  if( pt .ge. t1 .and. pt .le. t2 ) then
                     np = np + 1
                     if( np .gt. irows ) goto 1000
                     work_vec(np,1) = (cut_line_coord - x(i))/dz*
     1                  (z(inext(i)) - z(i)) + z(i)
                     work_vec(np,2) = pt
                     kreg(np) = ireg
                  endif
               endif
            endif
103      continue
100   continue
c
c     Sort in increasing order and remove any duplicate points
      call sort_uniq_1d(work_vec(1,2), work_vec(1,1), kreg, np)
c
c     Adjust co-incident X points so that similar regions are contiguous
      do 200 i=2,np-1
         if( work_vec(i,2) .eq. work_vec(i+1,2) .and.
     1       kreg(i+1) .eq. kreg(i-1) .and. kreg(i) .ne. kreg(i-1)) then
            tmp = work_vec(i+1,1)
            work_vec(i+1,1) = work_vec(i,1)
            work_vec(i,1)   = tmp
            tmp = work_vec(i+1,2)
            work_vec(i+1,2) = work_vec(i,2)
            work_vec(i,2)   = tmp
            itmp = kreg(i+1)
            kreg(i+1) = kreg(i)
            kreg(i)   = itmp
         endif
200   continue
c
      if( np .eq. 0 ) then
         write(luttyo,60)
60       format(' ','Error: No points extracted along cut line')
         err = .TRUE.
      endif
c
      return
c
1000  continue
c     Tried to overrun point storage
      write(luttyo,70)
70    format(' ','Error: Too many points on cut line')
      err = .TRUE.
      return
c
      end
c
      subroutine sort_uniq_1d(x,y,info,npts)
      implicit none
c
      integer npts, info(npts)
      real x(npts), y(npts)
c
      real tmp
      integer i, j, k, itmp
c
c     Simple (read slow) in place sort
      do 100 j=npts-1,1,-1
         do 100 i=1,j
            if( x(i+1) .le. x(i) ) then
               tmp    = x(i+1)
               x(i+1) = x(i)
               x(i)   = tmp
               tmp    = y(i+1)
               y(i+1) = y(i)
               y(i)   = tmp
               itmp    = info(i+1)
               info(i+1) = info(i)
               info(i)   = itmp
            endif
100   continue
c
c     Cull duplicates
      j = npts
      i = 1
      do while( i .lt. j )
         do while( x(i) .eq. x(i+1) .and. y(i) .eq. y(i+1)
     1             .and. i .lt. j)
            do 201 k=i+1,j-1
               x(k) = x(k+1)
               y(k) = y(k+1)
               info(k) = info(k+1)
201         continue
            j = j - 1
         end do
         i = i + 1
      end do
      npts = j
c
      return
      end
c
      subroutine copy_save_to_multi_2d(err)
      implicit none
c
      include 'sysdep.cmn'
      include 'const.h'
      include 'io.cmn'
      include 'plot.cmn'
      include 'mesh.cmn'
      include 'data.cmn'
      include 'units.cmn'
      include 'defaults.cmn'
      include 'compare_2d.cmn'
      include 'user.cmn'
      include 'analytical.cmn'
      include 'annotations.cmn'
      include 'save.cmn'
c
      logical err
c
      integer i, j, k, l, entry, ierr, index, menu_id, np, ip,
     1        itype, num_callable, find_new_entry_2d, delete_entry_2d
      logical read_error, load_error, bad_entry, update_datafile,
     1        lfound, init_contours
      character*40 datatype
c
c     Delete any existing compare_2d data, 
c     unless it has already been loaded via the callable interface.
c
      num_callable = 0
      do i=1,MULTI_SIZE
         k = slot_index_2d(i)
         if( k .ne. -1 ) then
            if( multi_2d(k)%datatype .ne. DT_CALLABLE ) then
               num_callable = num_callable + 1
            endif
         endif
      end do
c
      do while( n_data_2d .gt. num_callable )
         ierr = delete_entry_2d(1)
      end do
!!      j = n_data_2d
c
c     Store 2D dataset information
!!      n_data_2d = 0
      do 2000 i=1,num_curves
c
c        First, look for XY data already loaded via callable interface
         if( curve_dat(i)%datafile_type .eq.
     1       SAVE_DATATYPE_CALLABLE ) then
c           Look for matching entry ID string in data already loaded
            lfound = .FALSE.
            do 2001 k=1,MULTI_SIZE
               if( in_use_2d(k) ) then
                  l = slot_index_2d(k)
                  if( multi_2d(l)%datatype .eq. DT_CALLABLE .and.
     1                multi_2d(l)%label .eq. curve_dat(i)%label ) then
                     lfound = .TRUE.
                  endif
               endif
2001        continue
            if( .NOT. lfound ) then
               write(luttyo,2002) curve_dat(i)%label
2002           format(' ','Error: callable data not found'/
     1                ' ','       label = ',a)
               err = .TRUE.
            endif
            goto 2000
         endif
c
c        Else, read in data files
         entry = find_new_entry_2d()
         if( entry .eq. -1 ) then
            write(luttyo,2010)
2010        format(' ','Error: all data slots used up')
            err = .TRUE.
         else
            bad_entry = .FALSE.
            datatype = curve_dat(i)%datafile_type
            if( datatype .eq. SAVE_DATATYPE_GRAPES .or.
     1          datatype .eq. SAVE_DATATYPE_SPICE3 .or.
     2          datatype .eq. SAVE_DATATYPE_SUPREM3S .or.
     3          datatype .eq. SAVE_DATATYPE_SUPREM3P .or.
     4          datatype .eq. SAVE_DATATYPE_CALLABLE_XY ) then
               write(luttyo,2011)
2011           format(' ','Error: Cannot plot contours of 1D data')
               bad_entry = .TRUE.
            else if( datatype .eq. SAVE_DATATYPE_MINIMOS  .or.
     1          datatype .eq. SAVE_DATATYPE_PROMIS   .or.
     2          datatype .eq. SAVE_DATATYPE_VLSICAP  .or.
     3          datatype .eq. SAVE_DATATYPE_BAMBI    .or.
     4          datatype .eq. SAVE_DATATYPE_USEOUT   .or.
     5          datatype .eq. SAVE_DATATYPE_2DOP     .or.
     8          datatype .eq. SAVE_DATATYPE_SUPREM4  .or.
     9          datatype .eq. SAVE_DATATYPE_PISCES   .or.
     A          datatype .eq. SAVE_DATATYPE_TIF      .or.
     B          datatype .eq. SAVE_DATATYPE_FAIM     .or.
     C          datatype .eq. SAVE_DATATYPE_MCP2D    .or.
     C          datatype .eq. SAVE_DATATYPE_MEDICI   .or.
     D          datatype .eq. SAVE_DATATYPE_CALLABLE .or.
     F          datatype .eq. SAVE_DATATYPE_ASCII .or.
     G          datatype .eq. SAVE_DATATYPE_ANALYTICAL ) then
               j = j + 1
c              Set up data area and read file
c              Optimization: if file is loaded, and name is same, don't re-read
               if( update_datafile(i) ) then
                  prev_binfilename = binfilename
                  binfilename     = curve_dat(i)%datafile
                  plane_value     = curve_dat(i)%cut_plane_coord
                  cut_plane_axis  = curve_dat(i)%cut_plane_axis
                  datafile_record = curve_dat(i)%datafile_record
                  binfile_convert_type = curve_dat(i)%bin_datatype
                  if( curve_dat(i)%datafile_dimension .eq. 3 ) then
                     three_d = .TRUE.
                  else
                     three_d = .FALSE.
                  endif
                  if( datatype .eq. SAVE_DATATYPE_MINIMOS     ) then
                     binfiletype = DT_MINIMOS
                     call read_mini(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_PROMIS ) then
                     binfiletype = DT_PROMIS
                     call read_promis(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_SUPREM4 ) then
                     binfiletype = DT_SUPREM4
                     call read_sup4(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_PISCES ) then
                     binfiletype = DT_PISCES
                     medici_solution_file = curve_dat(i)%misc_ch(1)
                     call read_pisces(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_MEDICI ) then
                     binfiletype = DT_MEDICI
                     medici_solution_file = curve_dat(i)%misc_ch(1)
                     call read_medici(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_VLSICAP) then
                     binfiletype = DT_VLSICAP
                     vlsi_edge = curve_dat(i)%misc(1)
                     vlsi_yint = curve_dat(i)%misc(2)
                     call read_vlsicap(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_BAMBI  ) then
                     binfiletype = DT_BAMBI
                     call read_bambi(read_error)
                  else if( datatype .eq. SAVE_DATATYPE_USEOUT ) then
                     binfiletype = DT_USEOUT
                     call read_useout(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_2DOP   ) then
                     binfiletype = DT_2DOP
                     tdop_gate_length = curve_dat(i)%misc(1)
                     tdop_interface_file = curve_dat(i)%misc_ch(1)
                     call read_2dop(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_SUPREM3S) then
                     binfiletype = DT_SUPREM3S
                     call read_suprem3(2,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_SUPREM3P) then
                     binfiletype = DT_SUPREM3P
                     call read_suprem3(1,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_CALLABLE ) then
                     binfiletype = DT_CALLABLE
                     call load_callable_data(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_TIF ) then
                     tif_ixcol = curve_dat(i)%colx
                     tif_iycol = curve_dat(i)%coly(1)
                     binfiletype = DT_TIF
                     call read_tif(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_FAIM ) then
                     binfiletype = DT_FAIM
                     call read_faim(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_MCP2D ) then
                     binfiletype = DT_MCP2D
                     call read_mcp2d(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_GRAPES ) then
                     binfiletype = DT_GRAPES
                     call read_grapesfile(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_SPICE3 ) then
                     binfiletype = DT_SPICE3
                     call do_read_spice3file(.TRUE.,read_error)
                  else if( datatype .eq. SAVE_DATATYPE_ASCII ) then
                     binfiletype = DT_ASCII
c                    Set up user common area
                     user_1d_asciifile = curve_dat(i)%datafile
                     user_1d_col_x = curve_dat(i)%colx
                     user_1d_n_ycols = curve_dat(i)%n_ycols
                     user_1d_col_y(1) = curve_dat(i)%coly(1)
                     user_1d_col_z = curve_dat(i)%colz
                     user_1d_skip  = curve_dat(i)%skip
                     user_1d_nlines = curve_dat(i)%lines_to_read
                     user_1d_ax = curve_dat(i)%ax
                     user_1d_ay = curve_dat(i)%ay
                     user_1d_az = curve_dat(i)%az
                     user_1d_mx = curve_dat(i)%mx
                     user_1d_my = curve_dat(i)%my
                     user_1d_mz = curve_dat(i)%mz
                     user_1d_x_expression = curve_dat(i)%x_expression
                     user_1d_y_expression = curve_dat(i)%y_expression
                     user_1d_rect_grid =
     1                   .NOT. curve_dat(i)%data_is_trigrid
                     user_1d_dim  = curve_dat(i)%datafile_dimension
                     user_1d_dataset_number=curve_dat(i)%datafile_record
c                    These items are not used for 2D data,
c                    but we set them so they have some value...
                     user_1d_errorbars = curve_dat(i)%errorbars
                     user_1d_contour = curve_dat(i)%contour
                     user_1d_sweep = curve_dat(i)%sweep
                     user_1d_merge_ycols = curve_dat(i)%merge_ycols
                     user_1d_xabs = curve_dat(i)%xabs
                     user_1d_yabs = curve_dat(i)%yabs
                     user_1d_zabs = curve_dat(i)%zabs
                     user_1d_sort = curve_dat(i)%sort
                     user_1d_x_char_data = curve_dat(i)%x_char_data
                     user_1d_sample_freq = curve_dat(i)%sample_freq
c
                     if( user_1d_rect_grid ) then
                        call do_read_asciifile_2d
                     else
                        call do_read_asciifile_2d_tri
                     endif
                     read_error = .NOT. have_user_1d_data
                  else if( datatype .eq. SAVE_DATATYPE_ANALYTICAL ) then
                     binfiletype = DT_ANALYTICAL
                     analytical_expression = 
     1                  curve_dat(i)%analytical_expression
                     analytical_xstart = curve_dat(i)%analytical_xstart
                     analytical_xend   = curve_dat(i)%analytical_xend
                     analytical_xstep  = curve_dat(i)%analytical_xstep
                     analytical_ystart = curve_dat(i)%analytical_ystart
                     analytical_yend   = curve_dat(i)%analytical_yend
                     analytical_ystep  = curve_dat(i)%analytical_ystep
                     analytical_type   = curve_dat(i)%analytical_type
                     analytical_sweep  = curve_dat(i)%sweep
                     analytical_dim    = curve_dat(i)%datafile_dimension
                     call do_eval_analytical_2d
                     read_error = .NOT. have_analytical_data
                  endif
               endif
               if(.NOT. read_error) then
                  call lookup_menuid_from_nickname(
     1               curve_dat(i)%datafile_quantity, menu_id)
                  if( menu_id .eq. -1 ) then
                     write(luttyo,2020) curve_dat(i)%datafile_quantity
2020                 format(' ','Error: unknown quantity: ',a)
                     bad_entry = .TRUE.
                  else
                     index = menu_to_arr(menu_id)
                     if( is_data_1d(menu_id) ) then
                        write(luttyo,2012)
2012                    format(' ',
     1                     'Error: Cannot plot contours of 1D data')
                        bad_entry = .TRUE.
                     else
                        init_contours = .FALSE.
                        if( data_is_trigrid ) then
                           np = tri_npts
                        else
                           np = ni*nj
                        endif
                        call copy_2d_data_to_slot(entry,
     1                     arrdat(arrind(index)), np,
     1                     data_name(menu_id), data_nickname(menu_id),
     1                     data_units(menu_id),
     1                     curve_dat(i)%contours,
     1                     curve_dat(i)%num_contours, init_contours,
     1                     curve_dat(i)%label,
     1                     curve_dat(i)%contour_color,
     1                     curve_dat(i)%contour_junctions,
     1                     curve_dat(i)%contour_grid,
     1                     curve_dat(i)%contour_clabel,
     1                     curve_dat(i)%contour_solid,
     1                     curve_dat(i)%contour_colormat,
     1                     is_data_log(menu_id),
     1                     is_data_in_oxide(menu_id),
     1                     curve_dat(i)%ax, curve_dat(i)%mx,
     1                     curve_dat(i)%ay, curve_dat(i)%my,
     1                     curve_dat(i)%az, curve_dat(i)%mz,
     1                     curve_dat(i)%x_expression,
     1                     curve_dat(i)%y_expression,
     1                     curve_dat(i)%z_expression)
                     endif
                  endif
               else
                  bad_entry = .TRUE.
               endif
            else
               write(luttyo,2050) datatype
2050           format(' ','Error: unknown data type: ',a)
               bad_entry = .TRUE.
            endif
            if( bad_entry ) then
               ierr = delete_entry_2d(entry)
               j = j - 1
               err = .TRUE.
            endif
         endif
2000  continue
!!      n_data_2d = j
c
      return
      end
c      
      subroutine restore_multi_2d_plot(err)
      implicit none
c
      logical err
      integer mode
c
      call pm_copy_save_axis_to_multi_2d
c
      call compare_plot_2d(.TRUE.,.FALSE.)
c
      err = .FALSE.
c
      return
      end
c      
      subroutine pm_copy_save_axis_to_multi_2d
      implicit none
c
      include 'sysdep.cmn'
      include 'const.h'
      include 'plot.cmn'
      include 'units.cmn'
      include 'defaults.cmn'
      include 'data.cmn'
      include 'compare_2d.cmn'
      include 'annotations.cmn'
      include 'save.cmn'
      include 'io.cmn'
      include 'title.cmn'
c
c     Locals
      integer i
c
c     init_common will only allow one intialization, so it's safe to
c     call anywhere.
      call init_common
c
c     Set up sticky information
      twod_lfirst = .FALSE.
      twod_h_min  = pm_scale_plot2cm*axis_dat(1)%min
      twod_h_max  = pm_scale_plot2cm*axis_dat(1)%max
      twod_h_tic_major = pm_scale_plot2cm*axis_dat(1)%major_tic
      twod_h_minor_tic_freq = axis_dat(1)%minor_tic_freq
      twod_h_label = axis_dat(1)%label
      twod_h_units = axis_dat(1)%units
      twod_h_autoscale = axis_dat(1)%autoscale_min .or.
     1                   axis_dat(1)%autoscale_max
c
      twod_v_min  = pm_scale_plot2cm*axis_dat(2)%min
      twod_v_max  = pm_scale_plot2cm*axis_dat(2)%max
      twod_v_tic_major = pm_scale_plot2cm*axis_dat(2)%major_tic
      twod_v_minor_tic_freq = axis_dat(2)%minor_tic_freq
      twod_v_label = axis_dat(2)%label
      twod_v_units = axis_dat(2)%units
      twod_v_autoscale = axis_dat(2)%autoscale_min .or.
     1                   axis_dat(2)%autoscale_max
c
      twod_invert    = .NOT. axis_dat(2)%invert
c
c     Set the active contours to dataset 1
      twod_num_contours = curve_dat(1)%num_contours
      do 100 i=1,twod_num_contours
         twod_contours(i) = curve_dat(1)%contours(i)
100   continue
c
c     Overlay plot uses title_dat array, so no need to set it here
c
      lsavefile_read = .TRUE.
c
      return
      end
