      subroutine azplot(string,tutm,itower,az,ntower,coor,area)
      implicit double precision (a-h,o-z)
      double precision tutm(2,*),az(ntower),coor(2),area
      integer itower(ntower),ntower
      character*(*) string
c   local variables
c   maxchr - maximum number of characters to use
c            (determined by screen width)
c   maxlin - maximum number of lines to use
c            (determined by screen height)
      integer maxchr,maxlin
      parameter (maxchr=79,maxlin=24)
      character plot(maxchr,maxlin)
      double precision minx,miny,maxx,maxy,xsin,xcos,delta
c
c   define statement functions to convert utm coordinates
c   to plot coordinates
c
      ix(x)=int(((x-minx)/(maxx-minx))*dble(maxchr)+1.5)
      iy(y)=int(((maxy-y)/(maxy-miny))*dble(maxlin)+1.5)
c
c   clear plot area
c
       do 5 i=1,maxlin
         do 5 j=1,maxchr
            plot(j,i)=' '
 5          continue
c
c   get smallest and largest coordinates of plot
      minx=tutm(1,itower(1))
      miny=tutm(2,itower(1))
      maxx=tutm(1,itower(1))
      maxy=tutm(2,itower(1))
      do 10 i=2,ntower
         minx=min(minx,tutm(1,itower(i)))
         miny=min(miny,tutm(2,itower(i)))
         maxx=max(maxx,tutm(1,itower(i)))
         maxy=max(maxy,tutm(2,itower(i)))
 10      continue
      if (coor(1).gt.0.d0 .and. coor(2).gt.0.d0) then
         minx=min(minx,coor(1))
         miny=min(miny,coor(2))
         maxx=max(maxx,coor(1))
         maxy=max(maxy,coor(2))
      endif
c   set boundaries to +/- 10% of mins and maxs.
      minx=minx-(maxx-minx)*0.1
      miny=miny-(maxy-miny)*0.1
      maxx=maxx+(maxx-minx)*0.1
      maxy=maxy+(maxy-miny)*0.1
c
c   plot azimuths as string of *'s
c
      delta=min((maxy-miny)/dble(maxlin),(maxx-minx)/dble(maxchr))
      do 30 i=1,ntower
         xsin=sin(az(itower(i)))
         xcos=cos(az(itower(i)))
         ip=0
 40      ip=ip+1
         iiy=iy(tutm(2,itower(i))+dble(ip)*delta*xsin)
         if (iiy.lt.1 .or. iiy.gt.maxlin) go to 30
         iix=ix(tutm(1,itower(i))+dble(ip)*delta*xcos)
         if (iix.lt.1 .or. iix.gt.maxchr) go to 30
         plot(iix,iiy)=char(249)
         go to 40
 30      continue
c
c   plot towers as ^ symbols
c
      do 20 i=1,ntower
         plot(ix(tutm(1,itower(i))),iy(tutm(2,itower(i))))=char(30)
 20      continue
c
c   plot estimated location as X
c
      if (coor(1).gt.0.d0 .and. coor(2).gt.0.d0) then
         plot(ix(coor(1)),iy(coor(2)))=char(219)
      endif
c
c   put border on plot
c
      do 56 j=2,maxchr-1
         plot(j,1)=char(205)
 56      plot(j,maxlin)=char(205)
      do 57 i=2,maxlin-1
         plot(1,i)=char(186)
 57      plot(maxchr,i)=char(186)
      plot(1,1)=char(201)
      plot(1,maxlin)=char(200)
      plot(maxchr,1)=char(187)
      plot(maxchr,maxlin)=char(188)
c
c   print plot on screen
c
      do 60 i=1,maxlin
         write(*,'(1x,80a1)') (plot(j,i),j=1,maxchr)
 60      continue
      write(*,'(1x,a,2i8.7,f6.2,a)') string(1:46),int(coor(1)+0.5),
     1   int(coor(2)+0.5),area,' Continue?'
      read(*,*)
      return
      end
