c
c     Copyright (c) 1999,2000, Compaq Computer Corporation
c
c     See the file "license.terms" for information on usage and redistribution
c     of this file, and for a DISCLAIMER OF ALL WARRANTIES.
c
c   Subroutine to read promis solution files
c
c
      subroutine read_promis(batchmode,read_error)
      implicit none
c
      include 'const.h'
      include 'sysdep.cmn'
      include 'mesh.cmn'
      include 'data.cmn'
      include 'io.cmn'
      include 'plot.cmn'
c
c     read PROMIS binary solution file
c
      logical batchmode, read_error
c
c     locals
c
      integer i
      logical lpromis_tri
c
      call check_promis_version(lpromis_tri, read_error)
c
      if( lpromis_tri ) then
         call read_promis_tri_binary(batchmode, read_error)
      else
         call read_promis_binary(batchmode, read_error)
      endif
c
      if ( read_error ) goto 10000
c
c     compute some stuff related to mesh...
c
c     set up boundary limits information
      call coord_limits
c
      if( lpromis_tri ) then
         call load_promis_tri_specific
      else
         call load_promis_specific
      endif
c
      if( lfirst_read ) then
         lfirst_read = .FALSE.
         prev_plane_type = plane_type
         prev_binfilename = binfilename
         prev_binfiletype = binfiletype
      endif
      if( prev_binfiletype .ne. binfiletype ) then
         oned_lfirst = .TRUE.
         twod_lfirst = .TRUE.
         threed_lfirst = .TRUE.
      endif
      datafile_loaded = .TRUE.
c
      if( .NOT. lpromis_tri ) then
c        Set up some mesh spacings for rectagular mesh
         nxm = nx - 1
         nym = ny - 1
         nytotm = nytot - 1
c
c        inverse mesh spacing...
c
         do 210 i=1,nxm
210         edxyz(i,1) = 1.0/(xyz(i+1,1) - xyz(i,1))
         do 211 i=1,nytotm
211        edxyz(i,2) = 1.0/(xyz(i+1,2) - xyz(i,2))
         do 212 i=1,nz-1
212         edxyz(i,3) = 1.0/(xyz(i+1,3) - xyz(i,3))
c
c       average mesh spacing...
c
         dxyzm(1,1)  = (xyz(2,1)  - xyz(1,1))/2.0
         dxyzm(nx,1) = (xyz(nx,1) - xyz(nx-1,1))/2.0
         do 220 i=2,nxm
220         dxyzm(i,1) = (xyz(i+1,1) - xyz(i-1,1))/2.0
c
         dxyzm(1,2)     = (xyz(2,2)  - xyz(1,2))/2.0
         dxyzm(nytot,2) = (xyz(nytot,2) - xyz(nytot-1,2))/2.0
         do 221 i=2,nytotm
221         dxyzm(i,2) = (xyz(i+1,2) - xyz(i-1,2))/2.0
c
         dxyzm(1,3)     = (xyz(2,3)  - xyz(1,3))/2.0
         dxyzm(nz,3) = (xyz(nz,3) - xyz(nz-1,3))/2.0
         do 222 i=2,nz-1
222         dxyzm(i,3) = (xyz(i+1,3) - xyz(i-1,3))/2.0
      endif
c
10000 return
      end
c
      subroutine read_promis_tri_binary(batchmode, read_error)
      implicit none
c
      include 'const.h'
      include 'sysdep.cmn'
      include 'mesh.cmn'
      include 'data.cmn'
      include 'io.cmn'
      include 'plot.cmn'
      include 'bincodes.cmn'
c
c     read PROMIS binary postprocessing file
c
      logical batchmode, read_error
c
c     locals
c
      integer i, j, k, ilen, ioerr, string_len,
     1        local_record, nequ, read_type, int_vec(1), ierr, nrecl
      character*80 search_string
      character*70 save_text
      character*(PM_MAX_PATHLEN) expanded_name
      character*25 vernum
      character*1  answer
c
      write(luttyo,10)
10    format(' ',
     1       'Reading PROMIS triangular grid postprocessing file...')
c
c     Open PROMIS postprocessing file
      nrecl=0
      call pm_open_file(binfil, binfilename,
     1               'OLD','UNFORMATTED','READ','NONE','NONSHARED',
     2               nrecl, ioerr, ierr)

      if( ierr .ne. 0 ) goto 10000
      local_record = 0
c
c     get stem of binary file name
      call get_filestem(binfilename, binfilestem)
c
      do 51 j=1,NUM_DATA_ARRAYS
         in_binary_file(j) = .FALSE.
51    continue
c
c     Split up data storage into chunks
      call allocate_tri_arrays(ARR2DSIZE,NUM_DATA_ARRAYS,ierr)
c
c     Clear out arrays. Not really necessary.
      call zero_real_array(arrdat,ARR2DSIZE,NUM_DATA_ARRAYS)
c
      do i=1,MAX_REGIONS
         imat(i) = 0
      end do
      do i=1,MAX_TRI
         do j=1,3
           nop(j,i) = 0
         end do
      end do
c
      tri_npts = 0
      tri_ne   = 0
      tri_nreg = 0
c
500   continue
      if( .not. batchmode ) then
         write(luttyo,120)
120      format(' ','You may extract a PROMIS data record by:'/
     1      ' ','   R -- Record number'/
     1      ' ','   D -- Display record label and ask user whether ',
     1          'to extract it'/
     1      ' ','   F -- Find record containing a specified string ',
     1          ' in the label')
         call get_char(answer,'Enter read option (R/D/F) [D] ',0)
         call upcase(answer)
         if( answer .eq. 'R' ) then
            read_type = 1
126         continue
            call get_ints(int_vec,1,'Enter record number: ',0)
            if( int_vec(1) .le. 0 ) then
               write(luttyo,127)
127            format(' ','Invalid record number')
               goto 126
            endif
            datafile_record = int_vec(1)
         else if (answer .eq. 'F' ) then
            read_type = 2
            search_string = 'FINAL'
            call get_str_with_default(search_string,
     1           'Enter search string ',0)
            call upcase(search_string)
         else
            read_type = 3
         endif
      endif
c
c     Read header
      read(binfil,end=10010,err=10010,iostat=ioerr) vernum
      if( vernum .ne. 'PROMIS-PM-VERSION-1.000  ') goto 10020
c
525   read(binfil,end=10010,err=10010,iostat=ioerr) save_text
      local_record = local_record + 1
c
      if( batchmode ) then
         if( local_record .eq. datafile_record ) then
            answer = 'E'
            write(luttyo,125) datafile_record, save_text
125         format(' ','Extracting PROMIS record ',i5,' with label:'/
     1             ' ',2x,a)
         else
            answer = 'S'
         endif
      else
         if( read_type .eq. 3 ) then
            write(luttyo,130) local_record, save_text
130         format(' ','PROMIS record (',i5,') found with label:'/3x,a)
            call get_char(answer,'Extract or Skip? (E/S) [S] ',0)
            call upcase(answer)
         elseif( read_type .eq. 2 ) then
            expanded_name = save_text
            call upcase(expanded_name)
            ilen = string_len(search_string)
            if( INDEX(expanded_name,search_string(1:ilen)) .ge. 1 ) then
               write(luttyo,130) local_record, save_text
               call get_char(answer,'Extract or Skip? (E/S) [S] ',0)
               call upcase(answer)
            else
               answer = 'S'
            endif
         else
            if( local_record .eq. datafile_record ) then
               answer = 'E'
               write(luttyo,125) datafile_record, save_text
            else
               answer = 'S'
            endif
         endif
      endif
c
      if( answer .eq. 'E' ) then
         goto 550
      else
         call skip_promis_tri_data(binfil,
     1        nelem, cord, nop, tri2reg, imat, tri_nextel, xyz,
     2        arrdat, ARR2DSIZE, COMP_ARRAYSIZE, MAX_TRI, MAX_REGIONS,
     3        read_error)
         if( read_error ) then
            close(unit=binfil)
            return
         endif
      endif
      goto 525
c
550   continue
      read(binfil,end=10010,err=10010,iostat=ioerr)
     1    tri_npts,tri_ne,tri_nreg,nequ
      read(binfil,end=10010,err=10010,iostat=ioerr) (nelem(i),i=1,nequ)
      if( tri_npts .gt. ARR2DSIZE ) goto 10040
      if( tri_ne   .gt. MAX_TRI ) goto 10080
      if( tri_nreg .gt. MAX_REGIONS ) goto 10050
c     
c     Coordinates
      read(binfil,end=10010,err=10010,iostat=ioerr)
     1             ((cord(j,i),j=1,2),i=1,tri_npts)
      do i=1,tri_npts
         do j=1,2
            cord(j,i) = cord(j,i)*MU2CM
         end do
      end do
c     Triangles
      read(binfil,end=10010,err=10010,iostat=ioerr)
     1             ((nop(i,j),i=1,3),j=1,tri_ne)
c     Triangle to region
      read(binfil,end=10010,err=10010,iostat=ioerr)
     1             (tri2reg(i),i=1,tri_ne)
c     Region to material array
      read(binfil,end=10010,err=10010,iostat=ioerr)
     1             (imat(i),i=1,tri_nreg)
c     Triangle adjacency (neighbor) information
      read(binfil,end=10010,err=10010,iostat=ioerr)
     1             ((tri_nextel(j,i),j=1,3),i=1,tri_ne)
c
      tri_ndata = 0
c
      do j=1,nequ
         tri_ndata = tri_ndata + 1
         read(binfil,end=10010,err=10010,iostat=ioerr)
     1                (arrdat(i+arrind(tri_ndata)-1),i=1,tri_npts)
         in_binary_file(tri_ndata) = .TRUE.
      end do
      promis_nequ = nequ
c
c     Read 1D interface information
      read(binfil,end=10010,err=10010,iostat=ioerr) nx
      if( nx .ne. 0 ) then
         read(binfil,end=10010,err=10010,iostat=ioerr)
     1                   (xyz(j,1),j=1,nx)
         do i=1,2
            read(binfil,end=10010,err=10010,iostat=ioerr)
     1                   (arr2d(j,i),j=1,nx)
            in_binary_file(tri_ndata+i) = .TRUE.
         end do
      endif
c
c     Compute net doping
      tri_ndata = tri_ndata + 1
      do k=1,nequ
         if( nelem(k) .eq. 5 ) then
c           Element is Boron, negate sign
            do i=1,tri_npts
               arrdat(i+arrind(tri_ndata)-1) =
     1         arrdat(i+arrind(tri_ndata)-1) - arrdat(i+arrind(k)-1)
            end do
         else if( nelem(k) .eq. 15 .or. nelem(k) .eq. 33 .or.
     1      nelem(k) .eq. 51 ) then
c           Phos, Arsenic or Antimony
            do i=1,tri_npts
               arrdat(i+arrind(tri_ndata)-1) =
     1         arrdat(i+arrind(tri_ndata)-1) + arrdat(i+arrind(k)-1)
            end do
         endif
      end do
c
c     Set placeholder values for old grid info. X info required for
c     1D data.
105   continue
c     Set nx to 1 to make sure that if it is used as an adjustable
c     array bound, nothing fails.
      if( nx .eq. 0 ) nx = 1
      ny = 1
      nt = ny
      nyoxm = 0
      nytot= ny
      nxsp = 1
      nxsoxm = 1
      plane_index = 1
      prev_plane_type = plane_type
      plane_type = XY_PLANE
      xyz(1,3) = -xw/2.0
      nxy=nx*ny
      nxytot=nx*nytot
c
c     Load generic mesh pointers, depending on plane type
c
      call set_3d_index_limits
      cut_plane_axis = charb
c
c     This variable must be set by all read routines...
      data_is_trigrid = .TRUE.
c
      write(luttyo,110)
110   format(' ','...done')
      datafile_record = local_record
      read_error = .FALSE.
      close(unit=binfil)
      goto 10200
c
c     errors of commission/omission
c
c     I/O error on open of file
10000 continue
      call expand_filename(binfilename,expanded_name)
      ilen = string_len(expanded_name)
      write(luttyo,10001) ioerr, expanded_name(1:ilen)
10001 format(/' ','FORTRAN I/O error #',i5,' on open of file:'/
     1       ' ',a)
      write(luttyo,10002)
10002 format(' ','Check for these errors:'/
     1       ' ','   - Incorrect file, directory or disk name'/
     1       ' ','   - You must have read permission on the file'/)
      read_error = .TRUE.  ! go back and get new name...
      goto 10200
c
c     I/O error on read of file
10010 write(luttyo,10011) ioerr
10011 format(/' ','FORTRAN I/O error #',i5,' on read of file')
      read_error = .TRUE.
      close(unit=binfil)
      goto 10200
c
10020 continue
      write(luttyo,10021)
10021 format(' ','Error: Invalid header in PROMIS save file')
      read_error = .TRUE.
      close(unit=binfil)
      goto 10200
c
10040 continue
      write(luttyo,10041)
10041 format(' ','Error: Number of points exceeds internal storage')
      read_error = .TRUE.
      close(unit=binfil)
      goto 10200
c
10050 continue
      write(luttyo,10051)
10051 format(' ','Error: Number of regions exceeds internal storage')
      read_error = .TRUE.
      close(unit=binfil)
      goto 10200
c
10080 continue
      write(luttyo,10081)
10081 format(' ','Error: Number of triangles exceeds internal storage')
      read_error = .TRUE.
      close(unit=binfil)
      goto 10200
c
10200 continue
      return
      end
c
      subroutine load_promis_tri_specific
      implicit none
c
c     Load descriptive arrays with PROMIS specific data
c
      include 'const.h'
      include 'sysdep.cmn'
      include 'mesh.cmn'
      include 'data.cmn'
c
c     Data varies linear or logarithmic
      logical LOG, LIN
      parameter (LOG = .TRUE., LIN = .FALSE.)
c     Data has value in oxide, or not
      logical OX, NO_OX
      parameter (OX = .TRUE., NO_OX = .FALSE.)
c     Data jumps at oxide interfaces
      logical JUMP, NO_JUMP
      parameter (JUMP = .TRUE., NO_JUMP = .FALSE.)
c     Data is 1D or 2D
      logical ONED, TWOD
      parameter (ONED = .TRUE., TWOD = .FALSE.)
c
      integer i
      character*25 eltype
      character*20 nickname
      character*10 ext
c
c     Number of menu items
      if( nx .gt. 1 ) then
         max_menu_num = tri_ndata + 2
      else
         max_menu_num = tri_ndata
      endif
c
c     (2D) array of net doping
      net_dop_array_index = tri_ndata
c     File has only 2D data
      three_d = .FALSE.
c
      do i=1,promis_nequ
         if( nelem(i) .eq. 5 ) then
            eltype = 'Boron'
            ext = '.boron'
            nickname = 'B'
         else if( nelem(i) .eq. 15 ) then
            eltype = 'Phosphorous'
            ext = '.phos'
            nickname = 'P'
         else if( nelem(i) .eq. 33 ) then
            eltype = 'Arsenic'
            ext = '.as'
            nickname = 'AS'
         else if( nelem(i) .eq. 51 ) then
            eltype = 'Antimony'
            ext = '.sb'
            nickname = 'SB'
         else if( nelem(i) .eq. 94 ) then
            eltype = 'Clustered As'
            ext = '.as_cl'
            nickname = 'AS_INACTIVE'
         else if( nelem(i) .eq. 95 ) then
            eltype = 'Interstitials'
            ext = '.inter'
            nickname = 'INTERSTITIALS'
         else if( nelem(i) .eq. 96 ) then
            eltype = 'Vacancies'
            ext = '.vac'
            nickname = 'VACANCIES'
         else if( nelem(i) .eq. 97 ) then
            eltype = 'Electron concentration'
            ext = '.e_conc'
            nickname = 'ELECTRONS'
         else if( nelem(i) .eq. 98 ) then
            eltype = 'Hole concentration'
            ext = '.h_conc'
            nickname = 'HOLES'
         else if( nelem(i) .eq. 99 ) then
            eltype = 'Potential'
            ext = '.psi'
            nickname = 'POTENTIAL'
         else if( nelem(i) .eq. 101 ) then
            ext = '.oxi'
         else
            write(eltype,100) nelem(i)
100         format('Atomic number ',i5)
            if( nelem(i) .lt. 10 ) then
               write(nickname,101) nelem(i)
101            format('ATOMIC_NO_',i1)
            elseif( nelem(i) .lt. 100) then
               write(nickname,102) nelem(i)
102            format('ATOMIC_NO_',i2)
            else
               write(nickname,103) nelem(i)
103            format('ATOMIC_NO_',i3)
            endif
            ext = '.unk'
         endif
         if( nelem(i) .eq. 99 ) then
            call load_desc_arr(i, eltype, nickname,
     1                        'V','V cm','cm^{2}',ext,
     2                         LIN,OX,NO_JUMP,TWOD,i)
c        Special case for element "100" (arbitrary)
         else if( nelem(i) .eq. 100 ) then
            call load_desc_arr(i, eltype, nickname,
     1                        ' ',' ',' ',ext,
     2                         LIN,NO_OX,NO_JUMP,TWOD,i)
         else
            call load_desc_arr(i, eltype, nickname,
     1                        'cm^{-3}','cm^{-2}','cm^{-1}',ext,
     2                         LOG,NO_OX,NO_JUMP,TWOD,i)
         endif
      end do
c
      if( nx .gt. 1 ) then
      eltype = 'Si/SiO2 boundary'
      call load_desc_arr(promis_nequ+1, eltype, 'SI_SIO2_INTERFACE',
     1                  'microns','micron cm','micron cm^{2}',ext,
     2                   LIN,OX,NO_JUMP,ONED,1)
      eltype = 'SiO2 Thickness'
      call load_desc_arr(promis_nequ+2, eltype, 'SIO2_THICKNESS',
     1                   'microns','micron cm','micron cm^{2}',ext,
     2                   LIN,OX,NO_JUMP,ONED,2)
      endif
c
      call load_desc_arr(max_menu_num, 'Net doping', 'NET_DOPING',
     1                   'cm^{-3}','cm^{-2}','cm^{-1}','.net',
     2                   LOG,NO_OX,NO_JUMP,TWOD,tri_ndata)
c
      return
      end
c
      subroutine skip_promis_tri_data(luprom,
     1           nelem, cord, nop, tri2reg, imat, tri_nextel, xyz,
     2           arrdat, arr2dsize, comp_arraysize,
     3           max_tri, max_regions,
     4           read_error)
      implicit none
c
      integer luprom
      integer arr2dsize, comp_arraysize, max_tri, max_regions
      real    cord(2,arr2dsize), xyz(comp_arraysize,3),
     1        arrdat(arr2dsize)
      integer nelem(10), nop(3,max_tri), tri2reg(max_tri),
     1        imat(max_regions), tri_nextel(3,max_tri)
      logical read_error
c
      include 'sysdep.cmn'
      include 'io.cmn'
c
      integer i, j, tri_npts, tri_ne, tri_nreg, nequ, nx, ioerr
c
      read_error = .FALSE.
c     Number of points, triangles, etc.
      read(luprom,end=10010,err=10010,iostat=ioerr)
     1    tri_npts,tri_ne,tri_nreg,nequ
      if( tri_npts .gt. ARR2DSIZE .or.
     1    tri_ne   .gt. MAX_TRI   .or.
     2    tri_nreg .gt. MAX_REGIONS ) then
          write(luttyo,100)
100       format(' ','Error: file too large for POSTMINI arrays')
          read_error = .TRUE.
          return
      endif
c
c     Element types
      read(luprom,end=10010,err=10010,iostat=ioerr)
     1    (nelem(i),i=1,nequ)
c     Coordinates
      read(luprom,end=10010,err=10010,iostat=ioerr)
     1             ((cord(j,i),j=1,2),i=1,tri_npts)
c     Triangles
      read(luprom,end=10010,err=10010,iostat=ioerr)
     1             ((nop(i,j),i=1,3),j=1,tri_ne)
c     Triangle to region array
      read(luprom,end=10010,err=10010,iostat=ioerr)
     1             (tri2reg(i),i=1,tri_ne)
c     Region to material array
      read(luprom,end=10010,err=10010,iostat=ioerr)
     1             (imat(i),i=1,tri_nreg)
c     Triangle adjacency (neighbor) information
      read(luprom,end=10010,err=10010,iostat=ioerr)
     1             ((tri_nextel(j,i),j=1,3),i=1,tri_ne)
c
c     Read 2D data
      do i=1,nequ
         read(luprom,end=10010,err=10010,iostat=ioerr)
     1                (arrdat(j),j=1,tri_npts)
      end do
c
c     Read 1D interface information
      read(luprom,end=10010,err=10010,iostat=ioerr) nx
      do i=1,3
         read(luprom,end=10010,err=10010,iostat=ioerr)
     1                (xyz(j,1),j=1,nx)
      end do
c
      return
c
10010 continue
      write(luttyo,10011) ioerr
10011 format(' ','Error: error during read of PROMIS binary file'/
     1       ' ','FORTRAN I/O error number was ',i5)
      read_error = .TRUE.
      end
c
      subroutine check_promis_version(lpromis_tri, read_error)
      implicit none
c
      include 'sysdep.cmn'
      include 'io.cmn'
c
      logical lpromis_tri, read_error
c
      integer nrecl, ioerr, ierr
      character*25 version_string
c
      read_error = .FALSE.
      nrecl=0
      call pm_open_file(binfil, binfilename,
     1               'OLD','UNFORMATTED','READ','NONE','NONSHARED',
     2               nrecl, ioerr, ierr)
      if( ierr .ne. 0 ) goto 10000
      read(binfil,err=10000,end=10000) version_string
      if( version_string(1:17) .eq. 'PROMIS-PM-VERSION' ) then
         lpromis_tri = .TRUE.
      else
         lpromis_tri = .FALSE.
      endif
      close(unit=binfil)
c
      return
c
10000 continue
      read_error = .TRUE.
      lpromis_tri = .FALSE.
      end


c
      subroutine read_promis_binary(batchmode, read_error)
      implicit none
c
      include 'const.h'
      include 'sysdep.cmn'
      include 'mesh.cmn'
      include 'data.cmn'
      include 'io.cmn'
      include 'plot.cmn'
      include 'bincodes.cmn'
c
c     read double precision PROMIS binary file
c     Note that PROMIS defines X and Y  *REVERSED* from MINIMOS (sigh)
c
      logical batchmode, read_error
c
c     locals
c
      double precision vbndlo(10)
      double precision vbndup(10)
      integer intmod(10), nptx, npty, nequ
      integer i, j, k, kk, ilen, ioerr, num_2d_items, local_record,
     1        string_len, int_vec(1), read_type, ioff1, ioff2, ierr,
     2        nrecl
      character*1 answer
      character*80 search_string
      character*70 save_text, com_text
      character*(PM_MAX_PATHLEN) expanded_name
      logical linear(10), have_oxide
c
      write(luttyo,10)
10    format(' ','Reading PROMIS save file...',$)
c
c     Append .SAVE if no extension
      call add_extension(binfilename,'.save')
c
      nrecl=0
      call pm_open_file(binfil, binfilename,
     1               'OLD','UNFORMATTED','READ','NONE','NONSHARED',
     2               nrecl, ioerr, ierr)
      if( ierr .ne. 0 ) goto 10000
      local_record = 0
c
c     get stem of binary file name
      call get_filestem(binfilename, binfilestem)
c
c     clear out arrays
      do 45 j=1,3
         do 45 i=1,COMP_ARRAYSIZE
            xyz(i,j) = 0.0
            edxyz(i,j) = 0.0
            dxyzm(i,j) = 0.0
45    continue
      do 51 j=1,MENUSIZE
         in_binary_file(j) = .FALSE.
51       continue
c
500   continue
      if( .not. batchmode ) then
         write(luttyo,120)
120      format(' ','You may extract a PROMIS data record by:'/
     1      ' ','   R -- Record number'/
     1      ' ','   D -- Display record label and ask user whether ',
     1          'to extract it'/
     1      ' ','   F -- Find record containing a specified string ',
     1          ' in the label')
         call get_char(answer,'Enter read option (R/D/F) [D] ',0)
         call upcase(answer)
         if( answer .eq. 'R' ) then
            read_type = 1
126         continue
            call get_ints(int_vec,1,'Enter record number: ',0)
            if( int_vec(1) .le. 0 ) then
               write(luttyo,127)
127            format(' ','Invalid record number')
               goto 126
            endif
            datafile_record = int_vec(1)
         else if (answer .eq. 'F' ) then
            read_type = 2
            search_string = 'FINAL'
            call get_str_with_default(search_string,
     1           'Enter search string ',0)
            call upcase(search_string)
         else
            read_type = 3
         endif
      endif
c
c     1. first read title from binary file
525   read(binfil,err=10010,iostat=ioerr) com_text
      read(binfil,err=10010,iostat=ioerr) save_text
      local_record = local_record + 1
c
      if( batchmode ) then
         if( local_record .eq. datafile_record ) then
            answer = 'E'
            write(luttyo,125) datafile_record, save_text
125         format(' ','Extracting PROMIS record ',i5,' with label:'/
     1             ' ',2x,a)
         else
            answer = 'S'
         endif
      else
         if( read_type .eq. 3 ) then
            write(luttyo,130) local_record, save_text
130         format(' ','PROMIS record (',i5,') found with label:'/3x,a)
            call get_char(answer,'Extract or Skip? (E/S) [S] ',0)
            call upcase(answer)
         elseif( read_type .eq. 2 ) then
            expanded_name = save_text
            call upcase(expanded_name)
            ilen = string_len(search_string)
            if( INDEX(expanded_name,search_string(1:ilen)) .ge. 1 ) then
               write(luttyo,130) local_record, save_text
               call get_char(answer,'Extract or Skip? (E/S) [S] ',0)
               call upcase(answer)
            else
               answer = 'S'
            endif
         else
            if( local_record .eq. datafile_record ) then
               answer = 'E'
               write(luttyo,125) datafile_record, save_text
            else
               answer = 'S'
            endif
         endif
      endif
c
      if( answer .eq. 'E' ) then
         goto 550
      else
         call skip_promis_data(nequ,nptx,npty,
     1      nelem,linear,intmod,vbndlo,vbndup,10,work_dbl,work_dbl,
     1      work_dbl,COMP_ARRAYSIZE,ARR2DSIZE,binfil,read_error)
         if( read_error ) then
            close(unit=binfil)
            return
         endif
      endif
      goto 525
c
550   title = save_text
      read(binfil,end=10010,err=10010) nequ, nptx, npty
      nx = npty
      nytot = nptx
      nz = 1
      promis_nequ = nequ
c
c     Is array storage large enough?
      if( MAX(nx,nytot,nz) .gt. COMP_ARRAYSIZE .or.
     1    (nx*nytot) .gt. ARR2DSIZE .or.
     1    (nequ+1) .gt. NUM_DATA_ARRAYS) then
         write(luttyo,10050)
10050    format(' ','Error: binary file too large for POSTMINI arrays')
         read_error = .TRUE.  ! go back
         close(unit=binfil)
         return
      endif
c
      read(binfil,end=10010,err=10010) (nelem(i), i=1,nequ)
      read(binfil,end=10010,err=10010) (linear(i), i=1,nequ)
      read(binfil,end=10010,err=10010) (intmod(i), i=1,nequ)
      read(binfil,end=10010,err=10010) (vbndlo(i), i=1,nequ)
      read(binfil,end=10010,err=10010) (vbndup(i), i=1,nequ)
      read(binfil,end=10010,err=10010) (work_dbl(i), i=1,nptx)
      do i=1,nptx
         xyz(i,2) = MU2CM*SNGL(work_dbl(i))
      end do
      read(binfil,end=10010,err=10010) (work_dbl(i), i=1,npty)
      do i=1,npty
         xyz(i,1) = MU2CM*SNGL(work_dbl(i))
      end do
c
      plane_index = 1
      prev_plane_type = plane_type
      if( save_text(47:50) .eq. 'TLCS' ) then
         plane_type = XT_PLANE
      else if( save_text(47:50) .eq. 'TVCS' ) then
         plane_type = YT_PLANE
      else
         plane_type = XY_PLANE
      endif
      ny = nytot
      nt    = ny
c     No "oxide"
      nyoxm = 0
c     These values for nxsp and nxsoxm are innocuous (used by
c     POSTMINI to draw in source/drain contacts)
      nxsp = 2
      nxsoxm = 1
      nxy=nx*ny
      nxytot=nx*nytot
c
      call allocate_rect_arrays(nx,ny,NUM_DATA_ARRAYS,ierr)
      if( ierr .ne. 0 ) then
c        Failure during array allocation
         write(luttyo,10050)
         read_error = .TRUE.
         close(unit=binfil)
         return
      endif
c
      have_oxide = .FALSE.
      do i=1,nequ
         if( nelem(i) .eq. 101 ) have_oxide = .TRUE.
      end do
c
c     Increment number of menu items by net doping and (1D) oxide elements
      if( have_oxide ) then
         max_menu_num = nequ + 2
         num_2d_items = nequ
      else
         max_menu_num = nequ + 1
         num_2d_items = nequ + 1
      endif
c
      do i=1,num_2d_items
         call zero_real_array(arrdat(arrind(i)),nx*ny,1)
      end do
      if( have_oxide ) then
         do i=1,2
            call zero_real_array(arr2d(1,i),COMP_ARRAYSIZE,1)
         end do
      endif
c
c     Read in data elements
      
      kk = 1
      do k=1,nequ
c        Special case for oxide. Load it into two 2D array slots.
c        One for silicon surface, the other for oxide thickness
c        (half of array is lower oxide shape, the other half is upper
c        oxide shape)
         if( nelem(k) .eq. 101 ) then
            do i=1,nptx
               read(binfil,end=10010,err=10010,iostat=ioerr)
     1            (work_dbl(j), j=1,npty)
               if( i .eq. 1 ) then
                  do j=1,npty
                     arr2d(j,1) = SNGL(work_dbl(j))
                  end do
               else if( i .eq. nptx ) then
                  do j=1,npty
                     arr2d(j,2) = SNGL(work_dbl(j))
                  end do
c                 Compute oxide thickness
                  do j=1,npty
                     arr2d(j,2) = arr2d(j,1) - arr2d(j,2)
                  end do
               endif
            end do
            promis_2d_data_map(k) = 0
         else
            do i=1,nptx
               read(binfil,end=10010,err=10010,iostat=ioerr)
     1            (work_dbl(j), j=1,npty)
               ioff1 = arrind(kk) - 1 + (i-1)*npty
               do j=1,npty
                  arrdat(ioff1+j) = SNGL(work_dbl(j))
               end do
            end do
            promis_2d_data_map(k) = kk
            kk = kk + 1
         end if
      end do
c
c     Compute net doping
      do k=1,nequ
         kk = promis_2d_data_map(k)
         if( nelem(k) .eq. 5 ) then
c           Element is Boron, negate sign
            do i=1,nptx
               ioff1 = arrind(num_2d_items) - 1 + (i-1)*npty
               ioff2 = arrind(kk) - 1 + (i-1)*npty
               do j=1,npty
                  arrdat(ioff1+j) = arrdat(ioff1+j) - arrdat(ioff2+j)
               end do
            end do
         else if( nelem(k) .eq. 15 .or. nelem(k) .eq. 33 .or.
     1      nelem(k) .eq. 51 ) then
c           Phos, Arsenic or Antimony
            do i=1,nptx
               ioff1 = arrind(num_2d_items) - 1 + (i-1)*npty
               ioff2 = arrind(kk) - 1 + (i-1)*npty
               do j=1,npty
                  arrdat(ioff1+j) = arrdat(ioff1+j) + arrdat(ioff2+j)
               end do
            end do
         endif
      end do
c
      do i=1,max_menu_num
         in_binary_file(i) = .TRUE.
      end do
c
c     Load generic mesh pointers, depending on plane type
c
      call set_3d_index_limits
      cut_plane_axis = charb
c     This variable must be set by all read routines...
      data_is_trigrid = .FALSE.
c
      datafile_record = local_record
      write(luttyo,110)
110   format(' ','...done')
      read_error = .FALSE.
      close(unit=binfil)
      goto 10100
c
c     errors of commission/omission
c
c     I/O error on open of binary file
10000 continue
      call expand_filename(binfilename,expanded_name)
      ilen = string_len(expanded_name)
      write(luttyo,10001) ioerr, expanded_name(1:ilen)
10001 format(/' ','FORTRAN I/O error #',i5,' on open of binary file:'/
     1       ' ',a)
      write(luttyo,10002)
10002 format(' ','Check for these errors:'/
     1       ' ','   - Incorrect file, directory or disk name'/
     1       ' ','   - You must have read permission on the file'/)
      read_error = .TRUE.  ! go back and get new name...
      goto 10100
c
c
c     I/O error on read of binary file
10010 write(luttyo,10011) ioerr
10011 format(/' ','FORTRAN I/O error #',i5,' on read of binary file')
      if( ioerr .eq. 24 ) then       ! End-of-file during read
         write(luttyo,10012)
10012    format(' ','Unexpected end-of-file during read...'/
     1          ' ','File may not be a PROMIS binary file'/)
      else if( ioerr .eq. 35 ) then  ! Not binary file format
         write(luttyo,10013)
10013    format(' ','File is probably not a FORTRAN binary file...'/)
      else if( ioerr .eq. 67 ) then  ! Not enough data in record
         write(luttyo,10014)
10014    format(' ','Insufficient data in record'/)
      else if( ioerr .eq. -1 ) then  ! End of file at start of read
         write(luttyo,10015)
10015    format(' ','End-of-file reached')
         if( batchmode ) then
            write(luttyo,10016) datafile_record
10016       format(' ','Error: could not find record number: ',i5)
         else
            call get_char(answer,
     1         'Read file over from beginning? (Y/N) [Y] ',0)
            call upcase(answer)
            if( answer .ne. 'N' ) then
               rewind(binfil)
               local_record = 0
               goto 500
            endif
         endif
      endif
      read_error = .TRUE.  ! go back and get new name...
      close(unit=binfil)
      goto 10100
c
10100 continue
      return
      end
c
      subroutine load_promis_specific
      implicit none
c
c     Load descriptive arrays with PROMIS specific data
c
      include 'const.h'
      include 'sysdep.cmn'
      include 'data.cmn'
c
c     Data varies linear or logarithmic
      logical LOG, LIN
      parameter (LOG = .TRUE., LIN = .FALSE.)
c     Data has value in oxide, or not
      logical OX, NO_OX
      parameter (OX = .TRUE., NO_OX = .FALSE.)
c     Data jumps at oxide interfaces
      logical JUMP, NO_JUMP
      parameter (JUMP = .TRUE., NO_JUMP = .FALSE.)
c     Data is 1D or 2D
      logical ONED, TWOD
      parameter (ONED = .TRUE., TWOD = .FALSE.)
c
c     Locals
      integer i, j, k
      character*25 eltype
      character*20 nickname
      character*10 ext
c
c     File has 2d data
      three_d = .FALSE.
c
      j = 1
      k = 1
      do i=1,promis_nequ
         if( nelem(i) .eq. 5 ) then
            eltype = 'Boron'
            ext = '.boron'
            nickname = 'B'
         else if( nelem(i) .eq. 15 ) then
            eltype = 'Phosphorous'
            ext = '.phos'
            nickname = 'P'
         else if( nelem(i) .eq. 33 ) then
            eltype = 'Arsenic'
            ext = '.as'
            nickname = 'AS'
         else if( nelem(i) .eq. 51 ) then
            eltype = 'Antimony'
            ext = '.sb'
            nickname = 'SB'
         else if( nelem(i) .eq. 94 ) then
            eltype = 'Clustered As'
            ext = '.as_cl'
            nickname = 'AS_INACTIVE'
         else if( nelem(i) .eq. 95 ) then
            eltype = 'Interstitials'
            ext = '.inter'
            nickname = 'INTERSTITIALS'
         else if( nelem(i) .eq. 96 ) then
            eltype = 'Vacancies'
            ext = '.vac'
            nickname = 'VACANCIES'
         else if( nelem(i) .eq. 97 ) then
            eltype = 'Electron concentration'
            ext = '.e_conc'
            nickname = 'ELECTRONS'
         else if( nelem(i) .eq. 98 ) then
            eltype = 'Hole concentration'
            ext = '.h_conc'
            nickname = 'HOLES'
         else if( nelem(i) .eq. 99 ) then
            eltype = 'Potential'
            ext = '.psi'
            nickname = 'POTENTIAL'
         else if( nelem(i) .eq. 101 ) then
            ext = '.oxi'
         else
            write(eltype,100) nelem(i)
100         format('Atomic number ',i5)
            if( nelem(i) .lt. 10 ) then
               write(nickname,101) nelem(i)
101            format('ATOMIC_NO_',i1)
            elseif( nelem(i) .lt. 100) then
               write(nickname,102) nelem(i)
102            format('ATOMIC_NO_',i2)
            else
               write(nickname,103) nelem(i)
103            format('ATOMIC_NO_',i3)
            endif
            ext = '.unk'
         endif
         if( nelem(i) .eq. 99 ) then
            call load_desc_arr(j, eltype, nickname,
     1                        'V','V cm','cm^{2}',ext,
     2                         LIN,OX,NO_JUMP,TWOD,k)
            k = k + 1
            j = j + 1
         else if(nelem(i) .eq. 101 ) then
            eltype = 'Si/SiO2 boundary'
            call load_desc_arr(j, eltype, 'SI_SIO2_INTERFACE',
     1                        'microns','micron cm','micron cm^{2}',ext,
     2                         LIN,OX,NO_JUMP,ONED,1)
            eltype = 'SiO2 Thickness'
            call load_desc_arr(j+1, eltype, 'SIO2_THICKNESS',
     1                        'microns','micron cm','micron cm^{2}',ext,
     2                         LIN,OX,NO_JUMP,ONED,2)
            j = j + 2
c        Special case for element "100" (arbitrary)
         else if( nelem(i) .eq. 100 ) then
            call load_desc_arr(j, eltype, nickname,
     1                        ' ',' ',' ',ext,
     2                         LIN,NO_OX,NO_JUMP,TWOD,k)
            k = k + 1
            j = j + 1
         else
            call load_desc_arr(j, eltype, nickname,
     1                        'cm^{-3}','cm^{-2}','cm^{-1}',ext,
     2                         LOG,NO_OX,NO_JUMP,TWOD,k)
            k = k + 1
            j = j + 1
         endif
      end do
c
      call load_desc_arr(max_menu_num, 'Net doping', 'NET_DOPING',
     1                   'cm^{-3}','cm^{-2}','cm^{-1}','.net',
     2                   LOG,NO_OX,NO_JUMP,TWOD,k)
c
c     Array index of net doping
      net_dop_array_index = k
c
      return
      end
c
      subroutine skip_promis_data(nequ,nptx,npty,
     1      nelem,linear,intmod,vbndlo,vbndup,isize,x,y,
     1      array,ixysize,iarrsize,luprom,read_error)
      implicit none
c
      include 'sysdep.cmn'
      include 'io.cmn'
c
      integer nequ,nptx,npty,isize,ixysize,iarrsize,luprom
      logical linear(isize), read_error
      integer nelem(isize), intmod(isize)
      double precision vbndlo(isize), vbndup(isize)
      double precision x(ixysize), y(ixysize), array(iarrsize)
c
      integer ioerr, ioff, i, j, k
c
      read(luprom,end=4000,err=3000,iostat=ioerr) nequ, nptx, npty
      if( MAX(nptx,npty) .gt. ixysize .or.
     1    (nptx*npty) .gt. iarrsize ) then
         write(luttyo,100)
100      format(' ','Error: file too large for POSTMINI arrays')
         read_error = .TRUE.
         return
      endif
      read(luprom,end=4000,err=3000,iostat=ioerr) (nelem(i), i=1,nequ)
      read(luprom,end=4000,err=3000,iostat=ioerr) (linear(i), i=1,nequ)
      read(luprom,end=4000,err=3000,iostat=ioerr) (intmod(i), i=1,nequ)
      read(luprom,end=4000,err=3000,iostat=ioerr) (vbndlo(i), i=1,nequ)
      read(luprom,end=4000,err=3000,iostat=ioerr) (vbndup(i), i=1,nequ)
      read(luprom,end=4000,err=3000,iostat=ioerr) (x(i), i=1,nptx)
      read(luprom,end=4000,err=3000,iostat=ioerr) (y(i), i=1,npty)
c
      do k=1,nequ
         do i=1,nptx
            ioff = (i-1)*npty
            read(luprom,end=4000,err=3000,iostat=ioerr)
     1         (array(ioff+j), j=1,npty)
         end do
      end do
      read_error = .FALSE.
      return
c
3000  write(luttyo,3001) ioerr
3001  format(' ','Error: error during read of PROMIS binary file'/
     1       ' ','FORTRAN I/O error number was ',i5)
      read_error = .TRUE.
      return
c
4000  write(luttyo,4001) ioerr
4001  format(' ','Error: unexpected end of file during read ',
     1           'of PROMIS binary file'/
     1       ' ','FORTRAN I/O error number was ',i5)
      read_error = .TRUE.
      return
      end
