subroutine trwave(zz,rr,wez,wer,wbphi,i,wtp)
c	nn is used in this subroutine, so common/misc/
c	cannot be used here.
c---get field of traveling wave at r,z, and time wtp
c----------------------------------------------------------------------
c
      include 'param_sz.h'
      include 'cfldscom.h'
      include 'constcom.h'
      include 'coordcom.h'
      include 'debug.h'
      include 'pcordcom.h'
      include 'syscom.h'
      include 'tstepcom.h'
      include 'wavescom.h'
      include 'ucom.h'
c       
c--------------------------------------------------------------------------
c*
      wer=0.
      wez=0.
      wbphi=0.
c------------------------- if nprint gt 0 ---------------------------------
      if (rne.eq.0.)then
        do 100 i=1,nel
        if(zz.ge.zloc(i-1).and.zz.le.zloc(i).and.ntype(i).eq.9)goto 200
c---if zz is in range of a cell card output the field of that cell
        if(zz.ge.zloc(i-1).and.zz.le.zloc(i).and.ntype(i).eq.7)then
           ne=i
           rne=ne
           ph=wtp*radian*el(4,ne)/freq
           sp=sin(ph)
           cp=cos(ph)
           nc=el(6,ne)
           zc=zz+hcll(nc)-zloc(ne)
           if (el(8,ne).gt.0.) zc=zc-hcll(nc)
           call cfield(rr,zc,nc,ez,er,bt)
           s=sp*el(10,ne)+cp*el(9,ne)
           wez=el(5,ne)*100.*ez*s
           c=cp*el(10,ne)-sp*el(9,ne)
           wer=el(5,ne)*100.*er*s
           wbphi=el(5,ne)*100.*bt*c
           rne=0.
        endif ! cell
 100    continue
        i=0
        return
c---the z position is within the defined region of the i-th traveling wave---
 200    continue
        ne=i
        nc=el(6,ne)
        i=ne+1-nelw(nc)
c--------------------- end of if nprint gt 0 --------------------------------
      else
        ne=rne
        nc=el(6,ne)
        i=ne+1-nelw(nc)
      endif
c---------now add up all the spatial components from nmin to nmax-----------
      amp0=el(5,ne)
      do 300 n=nmin,nmax
        nn=n-nmin+1
        cay=twopi/wavel*el(8,ne)/freq
        cay3=wk3(nn,i,nc)
        cay1=wk1(nn,i,nc)
cliu    cay1rsq=(cay1*rr)**2
        cay1rsq=-cay1*rr**2         
        if(ne.ne.nelw(nc))then
          phase1=wtp*radian*el(8,ne)/freq
c	------------------------------------------
          if(ne.ne.nw(nc)) then
		phase2=-cay3*(zz-zloc(ne-1)-el(1,ne)/2.)
          else
		phase2=-cay3*(zz-zloc(ne-1)-el(1,ne))
          endif
c	----------------------------------------------
          phase3=wpoff(nn,i,nc)
          phase4=el(4,ne)
	  phase=phase1+phase2+phase3+phase4
c	  phase=wtp*radian*el(8,ne)/freq-cay3*(zz-zloc(ne-1)
C     &		-el(1,ne)/2.)+wpoff(nn,i,nc)+el(4,ne)
	else
          phase1=wtp*radian*el(8,ne)/freq
	  phase2=-cay3*(zz-zloc(ne-1))
          phase3=wpoff(nn,i,nc)
          phase4=el(4,ne)
          phase=phase1+phase2+phase3+phase4
c         phase=wtp*radian*el(8,ne)/freq-cay3*(ZZ-zloc(ne-1))
C     &		+wpoff(nn,i,nc)+el(4,ne)
	endif
	 amplitd=amp0*wa(nn,i,nc)
	 wez= wez+amplitd*xi0(cay1rsq)*cos(phase)
         dern= -amplitd*xi1ox(cay1rsq)*rr
cliu     wer=wer + dern*cos(phase)*cay3
         wer=wer + dern*sin(phase)*cay3    
         wbphi= wbphi+ dern*sin(phase)*cay
 300  continue
      return
      end
c++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++*