c-----*----------------------------------------------------6---------7--
c     Calculating lagged and leading correlation coefficients 
c       rt(-nt:nt) between two anomaly series x(i) and y(i).
c     input: n,x(n),y(n),and nt
c       n: number of data
c       x(n) and y(n): raw anomaly series
c       nt: input parameter, maximum lag or lead time
c     output: rt(-nt:nt)
c       rt(j): lagged and leading correlation coefficients between x and y
c              j<0, lagged correlations, x trails y;
c              j=0, simultaneous correltion bewteen x and y;
c              j>0, leaing correlations, x lead y. 
c     By Dr. LI Jianping, April 3, 2001.
      subroutine llcorrelation(n,x,y,nt,rt)
      dimension x(n),y(n),rt(-nt:nt)
      dimension x1(n),y1(n)               ! Work Array
c     Lagged correlations: x(i) trails y(i). 
      do 20 ll=1,nt
        do 10 i=1,n-ll
          x1(i)=x(i+ll)
          y1(i)=y(i)
  10    continue
        call correlation(n-ll,x1,y1,r)
        rt(-ll)=r
  20  continue
      call correlation(n,x,y,r)
      rt(0)=r
c     Leading correlations: x(i) leads y(i).
      do 40 ll=1,nt
        do 30 i=1,n-ll
          x1(i)=x(i)
          y1(i)=y(i+ll)
  30    continue
        call correlation(n-ll,x1,y1,r)
        rt(ll)=r
  40  continue
      return
      end
c-----*----------------------------------------------------6---------7--
c     For the correlation coefficient r between two series 
c       x(i) and y(i), where i=1,...,n.
c     input: n,x(n),y(n)
c       n: number of time series
c       x(n): raw series 
c       y(n): raw series 
c     output: r
c       r: correlation coefficient between x and y
c     By Dr. LI Jianping, January 5, 2000.
      subroutine correlation(n,x,y,r)
      dimension x(n),y(n)
      call meanvar(n,x,ax,sx,vx)
      call meanvar(n,y,ay,sy,vy)
      sxy=0.
      do 10 i=1,n
        sxy=sxy+(x(i)-ax)*(y(i)-ay)
  10  continue
      sxy=sxy/float(n)
      r=sxy/(sx*sy)
      return
      end
c-----*----------------------------------------------------6---------7--
c     Computing the mean ax, standard deviation sx
c       and variance vx of a series x(i) (i=1,...,n).
c     input: n and x(n)
c       n: number of raw series
c       x(n): raw series
c     output: ax, sx and vx
c       ax: the mean value of x(n)
c       sx: the standard deviation of x(n)
c       vx: the variance of x(n)
c     By Dr. LI Jianping, May 6, 1999.
      subroutine meanvar(n,x,ax,sx,vx)
      dimension x(n)
      ax=0.
      vx=0.
      sx=0.
      do 10 i=1,n
        ax=ax+x(i)
  10  continue
      ax=ax/float(n)
      do 20 i=1,n
        vx=vx+(x(i)-ax)**2
  20  continue
      vx=vx/float(n)
      sx=sqrt(vx)
      return
      end
