C     P R O G R A M   V E L P L
C     *************************
C
C     PROGRAM VELPL IS DESIGNED FOR THE PLOTTING OF PLANE SECTIONS
C     OF SLOWNESS, PHASE VELOCITY AND GROUP VELOCITY  SURFACES 
C     FOR GENERAL ANISOTROPIC MEDIA
C
C     ************************************************************
C
      character*80 text
C
      mode=0
      call serv(mode,lin,lou,llu,lu1,lu2)
      if(mode.eq.0)call plots(ldum1,ldum2,7)
      if(mode.eq.0)lin=5
      if(mode.eq.0)lou=6
C
      IRUN=0
      read(lin,100)text
      write(lou,100)text
      READ(lin,107)SC,ISHIFT,IPRINT,LU
      if(mode.ne.0)lu=llu
      WRITE(lou,107)SC,ISHIFT,IPRINT,LU
      IF(ISHIFT.EQ.0)ISHIFT=8
      SHIFT=FLOAT(ISHIFT)
C
      xlen=10.
      ylen=10.
      call plot(xlen,ylen,-3)
c
    1 continue
      do 3 i=1,3
      iplot=0
      read(lu,101)itype,npar,inum
      if(itype.eq.0)return
      if(npar.eq.1)write(lou,103)
      if(npar.eq.2)write(lou,104)
      if(npar.eq.3)write(lou,105)
      do 2 k=1,inum
      read(lu,102)ddelta,xv,yv,zv
      v=xv*xv+yv*yv+zv*zv
      if(npar.le.2)then
        vv=1./v
        v=sqrt(vv)
      end if
      if(npar.eq.1)then
        if(k.eq.1.and.i.eq.1)vconst=vv
        xv=xv*vconst
        yv=yv*vconst
        zv=zv*vconst
      end if
      if(npar.eq.2)then
        xv=xv*vv
        yv=yv*vv
        zv=zv*vv
      end if
      if(npar.eq.3)then
        v=sqrt(v)
      end if
      if(iprint.ne.0)write(lou,106)ddelta,v,xv,yv,zv
      call symbol(xv,zv,.05,char(1),0.,-1)
      iplot=1
      vaux=v
    2 continue
    3 continue
      call plot(shift,-yv-ylen,-3)
      go to 1
c
  100 FORMAT(A)
  101 FORMAT(16I5)
  102 FORMAT(4E15.5)
  103 FORMAT(5x,'ANGLE',5x,'PHVEL',8x,'PX',8x,'PY',8x,'PZ')
  104 FORMAT(5x,'ANGLE',5x,'PHVEL',7x,'VPX',7x,'VPY',7x,'VPZ')
  105 FORMAT(5x,'ANGLE',5x,'GRVEL',7x,'VGX',7x,'VGY',7x,'VGZ')
  106 format(5f10.5)
  107 format(f10.5,4i5)
c
   99 call plot(0.,0.,999)
      stop
      end