c--meanfft.f
c--computing the average amplitude of 5 quiet days 
c--modified from mfft.f by Su Xu 97/7/14
c---------------------------------------------------------------
      Program MEANFFT 
      parameter(ndim=1000000)
      implicit real*8(a-h,o-y),integer*4(i-n),logical(z)
      dimension amp(ndim),phase(ndim)
      character*50 title
      character*40 fname1,fname2
      character*20 fmt,dunit,funit
      common // x,y
c--
c--get filenames
      write(*,'(/a)') ' MEANFFT'
      write(*,'(a)')  ' Finds mean spectrum of 5 quiet days'
      write(*,'(a)')  '------------------------------------'
      write(*,'(a$)') ' Filename for input files --> '
      read(*,'(a)') fname2
      open(2,file=fname2)
c
c--zero the arrays
      do 10 i=1,1000000
        amp(i) = 0.d0
        phase(i) = 0.d0
  10  continue
c
c--now read each fft filename 
      nf = 0
  20  read(2,'(a)',end=50) fname1
c--make sure filename has no leading blanks
      len = lstr(fname1,-40)
c
c--read in its data
      open(1,file = fname1)
      write(*,'(2a)') '   reading ',fname1
      read(1,'(a)') title 
      read(1,'(3i10,g12.5,2i10)') inorm,no,nyq,delf,iw,lcb
      read(1,'(3a20)') fmt,dunit,funit
c
c--accumulate amplitude and phase from previous files
      do 30 i=1,nyq
	read(1,fmt) xa,xp
        amp(i) = amp(i)+xa
        phase(i) = phase(i)+xp
  30  continue
      nf = nf+1
      close(1)
      goto 20
c
c--do the averaging
  50  xnf = float(nf)
      do 60 i=1,nyq
        amp(i) = amp(i)/xnf
        phase(i) = phase(i)/xnf
  60  continue
c
c--write result        
      write(*,'(a$)') ' Filename for output file  --> '
      read(*,'(a)') fname2
      open(2,file=fname2)
      write(2,'(a)') title 
      write(2,'(3i10,g12.5,2i10)') inorm,no,nyq,delf,iw,lcb
      write(2,'(3a20)') fmt,dunit,funit
      write(2,fmt) (amp(i),phase(i),i=1,nyq)
c
      stop
      end
c----------------------------------------------------------------------------
      function lstr(name,n)
c
c function to return the length of a string. On input the string
c 'name' has a nominal length n; if n is negative, blank characters
c at the beginning of 'name' are first truncated.
c
      character*1 name(1)
c
c--zero length case first
      if(n.eq.0) then
        lstr = 0
        return
      endif
      k = 0
      ld = n
c
c--ld is negative, we truncate leading blanks
      if(ld.lt.0) then
        ld = -ld
        j = 0
  10    j = j+1
        if(name(j).eq.' ') goto 10
        nb = j-1
c
c--there are nb blanks to be erased in front of 'name'
        if(nb.gt.0) then
          do 20 i = nb+1,ld
            name(i-nb) = name(i)
            name(i) = ' '
  20      continue
c
c--reduce length of string by number of blanks
          ld = ld-nb
        endif
      endif
c
c--ld is non-zero, find length of string
c--first set true length equal to nominal length
      k = ld+1
c--and decrease k until a non-blank character appears or beginning
c--of string is reached
  40  k = k-1
      if(k.gt.0.and.name(k).eq.' ') goto 40
      lstr = k
c
      return
      end
  
