spgradx.f

Go to the documentation of this file.

C> @file
C>
C> Compute x-gradient in fourier space
C> @author IREDELL @date 96-02-20
 
C> THIS SUBPROGRAM COMPUTES THE X-GRADIENT OF FIELDS
C> IN COMPLEX FOURIER SPACE.
C> THE X-GRADIENT OF A VECTOR FIELD W IS
C> WX=CONJG(W)*L/RERTH
C> WHERE L IS THE WAVENUMBER AND RERTH IS THE EARTH RADIUS,
C> SO THAT THE RESULT IS THE X-GRADIENT OF THE PSEUDO-VECTOR.
C> THE X-GRADIENT OF A SCALAR FIELD W IS
C> WX=CONJG(W)*L/(RERTH*CLAT)
C> WHERE CLAT IS THE COSINE OF LATITUDE.
C> AT THE POLE THIS IS UNDEFINED, SO THE WAY TO GET
C> THE X-GRADIENT AT THE POLE IS BY PASSING BOTH
C> THE WEIGHTED WAVENUMBER 0 AND THE UNWEIGHTED WAVENUMBER 1 
C> AMPLITUDES AT THE POLE AND SETTING MP=10.
C> IN THIS CASE, THE WAVENUMBER 1 AMPLITUDES ARE USED
C> TO COMPUTE THE X-GRADIENT AND THEN ZEROED OUT.
C>
C> @param M        - INTEGER FOURIER WAVENUMBER TRUNCATION
C> @param INCW     - INTEGER FIRST DIMENSION OF THE COMPLEX AMPLITUDE ARRAY
C>                (INCW >= M+1)
C> @param KMAX     - INTEGER NUMBER OF FOURIER FIELDS
C> @param MP       - INTEGER (KM) IDENTIFIERS
C>                (0 OR 10 FOR SCALAR, 1 FOR VECTOR)
C> @param CLAT     - REAL COSINE OF LATITUDE
C> @param[out] W        - COMPLEX(INCW,KMAX) FOURIER AMPLITUDES
C>                CORRECTED WHEN MP=10 AND CLAT=0
C> @param[out] WX       - COMPLEX(INCW,KMAX) COMPLEX AMPLITUDES OF X-GRADIENTS
C>
C> @note THIS SUBPROGRAM IS THREAD-SAFE.
      SUBROUTINE spgradx(M,INCW,KMAX,MP,CLAT,W,WX)
 
        IMPLICIT NONE
        INTEGER,INTENT(IN):: M,INCW,KMAX,MP(KMAX)
        REAL,INTENT(IN):: CLAT
        REAL,INTENT(INOUT):: W(2*INCW,KMAX)
        REAL,INTENT(OUT):: WX(2*INCW,KMAX)
        INTEGER K,L
        REAL,PARAMETER:: RERTH=6.3712e6
C - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
        DO k=1,kmax
          IF(mp(k).EQ.1) THEN
            DO l=0,m
              wx(2*l+1,k)=-w(2*l+2,k)*(l/rerth)
              wx(2*l+2,k)=+w(2*l+1,k)*(l/rerth)
            ENDDO
          ELSEIF(clat.EQ.0.) THEN
            DO l=0,m
              wx(2*l+1,k)=0
              wx(2*l+2,k)=0
            ENDDO
            IF(mp(k).EQ.10.AND.m.GE.2) THEN
              wx(3,k)=-w(4,k)/rerth
              wx(4,k)=+w(3,k)/rerth
              w(3,k)=0
              w(4,k)=0
            ENDIF
          ELSE
            DO l=0,m
              wx(2*l+1,k)=-w(2*l+2,k)*(l/(rerth*clat))
              wx(2*l+2,k)=+w(2*l+1,k)*(l/(rerth*clat))
            ENDDO
          ENDIF
        ENDDO
C - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
      END SUBROUTINE