c----------------------------------------------------------------------
      subroutine master_draw
c----------------------------------------------------------------------
c  
c          This subroutine is the master drawing driver.  Depending on
c     the values of the logical variables used to specify the type
c     of drawing, it will divert the flow to the specific drawing
c     subroutines.
c
c          It does, in addition, provide some needed services like
c     sorting if needed, etc.
c
      parameter (natom = 10000)
      parameter (nbond = 50000)
c
      real max_z,min_z
      logical wireframe,ball_and_stick,need_sort,label_num,label_cha,
     &        need_lim,perspective,round_bonds,autoscl,showhb,havehb
c
      common/draw_log/wireframe,ball_and_stick,need_sort,label_num,
     &                label_cha,need_lim,perspective,round_bonds,
     &                autoscl,showhb,havehb
      common/mol_num/numat,num(natom),co(3,natom),atrad(natom)
      common/frame/scale,xori,yori,max_z,min_z
c
      call draw_grid
      call draw_labels
      call draw_all_arcs
      call draw_lines
      if (ball_and_stick) then
         if (need_sort) then
            call rank_seq(numat)
            need_sort=.false.
         endif
c        call calc_trap
         call draw_ball_stick
      elseif (wireframe) then
         call draw_wire
      endif
      return
      end
c----------------------------------------------------------------------
      subroutine draw_wire
c----------------------------------------------------------------------
c
c          This subroutine draws wireframe-type models.
c
      parameter (natom = 10000) 
      parameter (nbond = 50000) 
c
      integer bndtyp, reversed
      logical wireframe,ball_and_stick,need_sort,label_num,label_cha
      logical full_bond,need_lim,perspective,round_bonds,autoscl
      logical showhb,havehb
      logical show_atom
c
      common/bonds/numbonds,ibond(nbond),jbond(nbond),nbonds(nbond),
     &             lastbond(natom),numpoint
      common/bond_type/bndtyp(nbond)
      common/draw_log/wireframe,ball_and_stick,need_sort,label_num,
     &                label_cha,need_lim,perspective,round_bonds,
     &                autoscl,showhb,havehb
      common/mol_num/numat,num(natom),co(3,natom),atrad(natom)
      common/show/show_atom(natom)
      common/per_color/iatcolor(105), reversed
c
      inibond=1
      do i=1,numat
         if (show_atom(i) .and. (label_num .or. label_cha))
     &                       call label_atom(i)
         x_i = co(1,i)
         y_i = co(2,i)
         z_i = co(3,i)
         do k=inibond,lastbond(i)
            l=nbonds(k)
            full_bond = bndtyp(l).eq.0
            if (i.eq.ibond(l)) then
               j=jbond(l)
            else
               j=ibond(l)
            endif
            if (j.gt.i .and. (show_atom(i).and.show_atom(j)) ) then
               if(bndtyp(l).ne.2 .or. (bndtyp(l).eq.2 .and. showhb))then
                  x_j = co(1,j)
                  y_j = co(2,j)
                  z_j = co(3,j)
                  if (perspective) then
                     call plot_md(persp_x(x_i,z_i),persp_y(
     &                            y_i,z_i),3)
                  else
                     call plot_md(x_i,y_i,3)
                  endif
                  xhalf = 0.5*(x_i+x_j)
                  yhalf = 0.5*(y_i+y_j)
                  zhalf = 0.5*(z_i+z_j)
                  mol1 = num(i)
                  mol2 = num(j)
c
c                 Figure out the colors for half-bonds.
c
                  call newpen(iatcolor(mol1))

                  if (.not.full_bond) call setdash(1) 
                  if (perspective) then
                     call plot_md(persp_x(xhalf,zhalf),persp_y(
     &                            yhalf,zhalf),2)
                  else
                     call plot_md(xhalf,yhalf,2)
                  endif
c
c second half of bond
c
                  call newpen(iatcolor(mol2))

                  if (.not.full_bond) call setdash(1) 
                  if (perspective) then
                     call plot_md(persp_x(x_j,z_j),persp_y(
     &                            y_j,z_j),2)
                  else
                     call plot_md(x_j,y_j,2)
                  endif
                  call newpen(1)
                  if (.not.full_bond) call setdash(0)
               endif
            endif
         enddo
         inibond = lastbond(i)+1
      enddo
      return
      end
c-----------------------------------------------------------------------
      subroutine draw_ball_stick
c-----------------------------------------------------------------------
c
c          This subroutine draws ball and stick models.
c
      parameter (natom = 10000) 
      parameter (nbond = 50000) 
c 
      integer atom
      integer bndtyp
      logical wireframe,ball_and_stick,need_sort,label_num,
     &        label_cha,need_lim,perspective,round_bonds,
     &        autoscl,showhb,havehb
      logical show_atom
c
      common/draw_log/wireframe,ball_and_stick,need_sort,label_num,
     &                label_cha,need_lim,perspective,round_bonds,
     &                autoscl,showhb,havehb
      common/mol_num/numat,num(natom),co(3,natom),atrad(natom)
      common/sortseq/iseq(natom),irank(natom)
      common/bonds/numbonds,ibond(nbond),jbond(nbond),nbonds(nbond),
     &             lastbond(natom),numpoint
      common/bond_type/bndtyp(nbond)
      common/show/show_atom(natom)
c
      do atom=1,numat
         i=iseq(atom)
         if (show_atom(i)) then
c
c           plot the atom numbered i
c
            call draw_at(i)
            if(label_num .or. label_cha) call label_atom(i)
c
c           i = number of atom
c
            inibond=1
            if (i.gt.1) inibond=lastbond(i-1)+1
            do k=inibond,lastbond(i)
               l=nbonds(k)
               if (i.eq.ibond(l)) then
                  j=jbond(l)
               else
                  j=ibond(l)
               endif
               if(bndtyp(l).ne.2 .or. (bndtyp(l).eq.2 .and. showhb))then
                  if (show_atom(j)) then
                     if (irank(j).gt.irank(i)) call draw_bond(l)
                   endif
                endif
             enddo
         endif
      enddo
      return
      end
c----------------------------------------------------------------------
      subroutine draw_at(atom)
c----------------------------------------------------------------------
c
c          This subroutine draws a colored circle for the ATOM atom.
c
      parameter (natom = 10000) 
      parameter (nbond = 50000) 
c 
      integer atom
c
      logical wireframe,ball_and_stick,need_sort,label_num,
     &                label_cha,need_lim,perspective,round_bonds,
     &                autoscl,showhb,havehb

c
      common/draw_log/wireframe,ball_and_stick,need_sort,label_num,
     &                label_cha,need_lim,perspective,round_bonds,
     &                autoscl,showhb,havehb
      common/mol_num/numat,num(natom),co(3,natom),atrad(natom)
c
      if (perspective) then
         z = co(3,atom)
         call poly_gon(persp_x(co(1,atom),z),persp_y(co(2,atom),z),
     &                   persp(atrad(atom),z),num(atom))
      else
         call poly_gon(co(1,atom),co(2,atom),atrad(atom),num(atom))
      endif
      return
      end
c-----------------------------------------------------------------------
      subroutine draw_bond(bond)
c-----------------------------------------------------------------------
c
      parameter (natom = 10000) 
      parameter (nbond = 50000) 
c 
      integer bond
      logical wireframe,ball_and_stick,need_sort,label_num,round_bonds,
     &                label_cha,need_lim,perspective,full_bond,
     &                autoscl,showhb,havehb
c
      dimension xb(20),yb(20)
c
      common/bonds/numbonds,ibond(nbond),jbond(nbond),nbonds(nbond),
     &             lastbond(natom),numpoint
      common/bond_type/bndtyp(nbond)
      common/draw_log/wireframe,ball_and_stick,need_sort,label_num,
     &                label_cha,need_lim,perspective,round_bonds,
     &                autoscl,showhb,havehb
c
c
c     plot the bond numbered "BOND"
c
      full_bond = (bndtyp(bond).eq.0)
      call calc_bond(bond,xb,yb,numpbond)
      call fill_bond(xb,yb,full_bond,numpbond)
c      
      return
      end
c----------------------------------------------------------------------
      subroutine label_atom(iat)
c----------------------------------------------------------------------
c
c          This subroutine writes the atomic label for atom IAT.  It
c     can be either the sequential number in the input or the string
c     label given in the input file, depending on the logical
c     variables LABEL_NUM and LABEL_CHA.
c
      parameter (natom = 10000) 
      parameter (nbond = 50000) 
c 
      real max_z,min_z
c
      integer reversed
      logical wireframe,ball_and_stick,need_sort,label_num,label_cha
      logical dark_atom,need_lim,perspective,round_bonds,autoscl
      logical showhb,havehb
      logical color_ws,mono_ws,PostScript
      character*4 str,atlabel
c
      common/draw_log/wireframe,ball_and_stick,need_sort,label_num,
     &                label_cha,need_lim,perspective,round_bonds,
     &                autoscl,showhb,havehb
      common/frame/scale,xori,yori,max_z,min_z
      common/mol_char/atlabel(natom)
      common/mol_num/numat,num(natom),co(3,natom),atrad(natom)
      common /device_type/ color_ws,mono_ws,PostScript
      common/per_tbl_num/radii(105)
      common/draw_spec/rad_fact,bond_fact,top_z,bottom_z
      common/per_color/iatcolor(105), reversed
c
c     define a label size for dummy atoms
c
      def_radius = radii(6)*rad_fact
c
c     In the case of dark atoms change line color to white, unless
c     reversed was specified.
c
      if (reversed .eq. 0) then
      if (color_ws .and. ball_and_stick) then
         if (num(iat).eq.6 .or. num(iat).eq.7 .or. num(iat).eq.8) then
            dark_atom = .true.
            call char_color(8)
         else
            dark_atom = .false.
         endif
      elseif (PostScript .and. ball_and_stick) then
         if (num(iat).eq.6 .or. num(iat).eq.8) then
            dark_atom = .true.
            call char_color(8)
         else
            dark_atom = .false.
         endif
      endif 
      endif
c
c     calculate the screen position and dimension of the string :
c
      x=co(1,iat)
      y=co(2,iat)
      tmprad = atrad(iat)
      if (num(iat).eq.99) tmprad = def_radius
      ahi=0.666666*tmprad / scale
c
c
c     then set the string :
c
      if (label_num) write(str,'(i4)') iat
      if (label_cha) str=atlabel(iat)
      call limits(str,is,ie)
      nchar=ie-is+1
c
      if (perspective) then
           x=persp_x(x,co(3,iat))
           y=persp_y(y,co(3,iat))
           ahi= persp(ahi*scale,co(3,iat))/scale
      endif
c
      x=(x-xori)/scale-0.4e0*float(nchar)*ahi
      y=(y-yori)/scale-0.5*ahi
c
c     plot the string :
c
      call plot(x,y,3)
      call symbol(x,y,ahi,str(is:ie),0.0,nchar)
c
c     return line color to black, if needed :
c
      if (dark_atom) call char_color(1)
      return 
      end
c-----------------------------------------------------------------------
      subroutine poly_gon(x0,y0,r,mol)
c-----------------------------------------------------------------------
c
c     Draws a MOL-dependent color-filled circle of radius R at X0-Y0
c
      real max_z,min_z
      integer reversed
c
      common/frame/scale,xori,yori,max_z,min_z
      common/per_color/iatcolor(105), reversed
c
      if (reversed .eq. 0) then
         call polycolor(iatcolor(mol))
      else if (reversed .eq. 1) then
         call polycolor(8)
      else if (reversed .eq. 2) then
         call polycolor(8)
      endif
c
      ax0 = (x0 - xori) / scale
      ay0 = (y0 - yori) / scale
      ar  = r / scale
      call circle(ax0,ay0,ar)
      return
      end
c-----------------------------------------------------------------------
      subroutine fill_bond(xb,yb,full_bond,numpbond)
c-----------------------------------------------------------------------
c
c     Draws a white-filled trapeze for a bond as specified by the
c     coordinates of its numpbond verteces
c
      real max_z,min_z
      logical full_bond
      integer reversed
c
      dimension xb(20),yb(20),axp(20),ayp(20)
c
      common/frame/scale,xori,yori,max_z,min_z
      common/per_color/iatcolor(105), reversed
c
      if (reversed .eq. 0) then
         call polycolor(8)
      else if (reversed .eq. 1) then
         call polycolor(7)
      else if (reversed .eq. 2) then
         call polycolor(8)
      endif
      if (.not.full_bond) then
         call setdash(1)
      endif
      do i =1,numpbond
         axp(i) = (xb(i) - xori) / scale
         ayp(i) = (yb(i) - yori) / scale
      enddo
      call polygon(axp,ayp,numpbond)
      if (.not. full_bond) then
         call setdash(0)
      endif
      return
      end
c----------------------------------------------------------------------
      subroutine auto_scale
c----------------------------------------------------------------------
c
c          When this subroutine is called it sets all the parameters needed
c     for the correct plotting of a structure in the screen and plotter.
c    
c          It initially finds the maximum X and Y distances, from
c     them and the screen dimensions it evaluates two separate scale
c     factors XSCALE and YSCALE.  It takes the largest of the two, and
c     uses it to figure the X and Y origin coordinates needed to center
c     the figure.
c
c                                              Apr 16, 1987
c
      parameter (natom = 10000) 
      parameter (nbond = 50000) 
c 
      real max_z,min_z
c
      common/mol_num/numat,num(natom),co(3,natom),atrad(natom)
      common/frame/scale,xori,yori,max_z,min_z
      common/useratt/user_scale
c
      scr_x = 10.0e0
      scr_y = 8.15e0
      xscale=0.0
      yscale=0.0
      scale=0.0
      xmin=co(1,1)-atrad(1)
      xmax=co(1,1)+atrad(1)
      ymin=co(2,1)-atrad(1)
      ymax=co(2,1)+atrad(1)
      zmin=co(3,1)-atrad(1)
      zmax=co(3,1)+atrad(1)
c
      do i=2,numat
         xmr=co(1,i)-atrad(i)
         xpr=co(1,i)+atrad(i)
         ymr=co(2,i)-atrad(i)
         ypr=co(2,i)+atrad(i)
         zmr=co(3,i)-atrad(i)
         zpr=co(3,i)+atrad(i)
         xmin=min(xmin,xmr)
         xmax=max(xmax,xpr)
         ymin=min(ymin,ymr)
         ymax=max(ymax,ypr)
         zmin=min(zmin,zmr)
         zmax=max(zmax,zpr)
      enddo
c
      xscale=1.05*(xmax-xmin)/scr_x
      yscale=1.05*(ymax-ymin)/scr_y
      scale=max(xscale,yscale)
c
      scale = scale / user_scale
      xori=xmin+(xmax-xmin-(scr_x*scale))*0.5e0
      yori=ymin+(ymax-ymin-(scr_y*scale))*0.5e0
      return
      end
c----------------------------------------------------------------------
      subroutine plot_md(x,y,mode)
c----------------------------------------------------------------------
c
c          This subroutine will act as an intermediate between this
c     program and the graphics subroutine.  Its main function is
c     to take the absolute coordinates given by the program and scale
c     them to the range required for screen and/or paper plotting by
c     the plotting library.
c
      real max_z,min_z
c
      common/frame/scale,xori,yori,max_z,min_z
c
      call plot((x-xori)/scale,(y-yori)/scale,mode)
      return
      end