spgradx.f

Go to the documentation of this file.

C> @file
C> @brief 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>
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>
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>
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> @note This subprogram is thread-safe.
C>
C> @param M Fourier wavenumber truncation
C> @param INCW first dimension of the complex amplitude array
C> (INCW >= M+1)
C> @param KMAX number of Fourier fields
C> @param MP identifiers
C> (0 or 10 for scalar, 1 for vector)
C> @param CLAT cosine of latitude
C> @param[out] W Fourier amplitudes corrected when MP=10 and CLAT=0
C> @param[out] WX complex amplitudes of x-gradients
C>
C> @author IREDELL @date 96-02-20
      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
 
        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
 
      END SUBROUTINE