c***********************************************************************

      subroutine ufidl_data

c***********************************************************************

      INCLUDE 'ssbox.inc'

      real FAMPROO(6000),FPHASEROO(6000)
      real FAMPVX(6000),FPHASEVX(6000)
      real FAMPVY(6000),FPHASEVY(6000)
      real FAMPTT(6000),FPHASETT(6000)
      real vxt(nstar),vyt(nstar)
      real dispvx,dispvy,ttmean

      real vxm(12000),vym(12000),vxxm(12000),vyym(12000)
      real vzzm(12000),nm(12000),ttm(12000)

      if(i60.eq.0) return


      dispvx=0.
      dispvy=0.
      do 2001 I=1,lkm
         vxt(i)=rp(i,1)
         vyt(i)=rp(i,2)+shear*r(i,1)
         dispvx=dispvx+vxt(i)**2
         dispvy=dispvy+vyt(i)**2
 2001 continue
      dispvx=sqrt(dispvx/lkm)
      dispvy=sqrt(dispvy/lkm)


      do 2012 ii=1,i60
         famproo(ii)=0.
         fampvx(ii)=0.
         fampvy(ii)=0.
         famptt(ii)=0.
         fphaseroo(ii)=0.
         fphasevx(ii)=0.
         fphasevy(ii)=0.
         fphasett(ii)=0.
 2012 continue

c     divide system into i60*2 zones

      i120=2*i60

      do 2013 ii=1,i120
         vxm(ii)=0.
         vym(ii)=0.
         vxxm(ii)=0.
         vyym(ii)=0.
         vzzm(ii)=0.
         nm(ii)=0.
         ttm(ii)=0.
 2013 continue

      do 2014 i=1,lkm
         izone=int((r(i,1)+xl)/2./xl*i120)+1
         if(izone.le.0) then izone=izone+i120
         if(izone.gt.i120) then izone=izone-i120
         if(izone.ge.1.and.izone.le.i120) then
            vxm(izone)=vxm(izone)+vxt(i)
            vym(izone)=vym(izone)+vyt(i)
            vxxm(izone)=vxxm(izone)+vxt(i)**2
            vyym(izone)=vyym(izone)+vyt(i)**2
            vzzm(izone)=vzzm(izone)+rp(i,3)**2
            nm(izone)=nm(izone)+1.
         endif
 2014 continue

      ttmean=0.
      nttmean=0
      do 2015 ii=1,i120
         if(nm(ii).ne.0) then
            vxm(ii)=vxm(ii)/nm(ii)
            vym(ii)=vym(ii)/nm(ii)
            vxxm(ii)=vxxm(ii)/nm(ii)
            vyym(ii)=vyym(ii)/nm(ii)
            vzzm(ii)=vzzm(ii)/nm(ii)
            ttm(ii)=(vxxm(ii)-vxm(ii)**2+vyym(ii)-vym(ii)**2+vzzm(ii))/3.
            ttmean=ttmean+ttm(ii)
            nttmean=nttmean+1
         endif
 2015 continue
      IF(nttmean.ne.0) ttmean=ttmean/nttmean

      do 2023 II=1,i60

         fampcr=0.
         fampsr=0.
         fampcx=0.
         fampsx=0.
         fampcy=0.
         fampsy=0.
         fampct=0.
         fampst=0.

         do 2024 i=1,i120
            fiit=(i-0.5)/i120*2.*pii-pii
            sinf=sin(fiit*II)
            cosf=cos(fiit*II)
            fampcr=fampcr+cosf*nm(i)
            fampsr=fampsr+sinf*nm(i)
            fampcx=fampcx+cosf*vxm(i)
            fampsx=fampsx+sinf*vxm(i)
            fampcy=fampcy+cosf*vym(i)
            fampsy=fampsy+sinf*vym(i)
            fampct=fampct+cosf*ttm(i)
            fampst=fampst+sinf*ttm(i)
 2024    continue

         famproo(ii)=2.*sqrt((fampcr/i120)**2+(fampsr/i120)**2)/
     +        (lkm/i120)
         fphaseroo(ii)=atan2(fampsr,fampcr)/ii
         fampvx(ii)=2.*sqrt((fampcx/i120)**2+(fampsx/i120)**2)
         fphasevx(ii)=atan2(fampsx,fampcx)/ii
         fampvy(ii)=2.*sqrt((fampcy/i120)**2+(fampsy/i120)**2)
         fphasevy(ii)=atan2(fampsy,fampcy)/ii
         famptt(ii)=2.*sqrt((fampct/i120)**2+(fampst/i120)**2)
         fphasett(ii)=atan2(fampst,fampct)/ii
 2023 continue

c     normalization
c     famproo     relative perturbation in surface density
c     fampvx,vy   absolute perturbation in m/sec
c     famptt      absolute perturbation in (m/sec)^2

c     dispvx,dispvy  m/sec
c     ttmean (m/sec)^2


      write(71) orb,dispvx,dispvy,ttmean
      write(71) i60
      write(71) (famproo(ii),ii=1,i60)
      write(71) (fampvx(ii),ii=1,i60)
      write(71) (fampvy(ii),ii=1,i60)
      write(71) (famptt(ii),ii=1,i60)
      write(71) (fphaseroo(ii),ii=1,i60)
      write(71) (fphasevx(ii),ii=1,i60)
      write(71) (fphasevy(ii),ii=1,i60)
      write(71) (fphasett(ii),ii=1,i60)

      RETURN
      END