pkg/zonal_filt/zonal_filt_nofill.F

42c525bfb4 Alis*0001 #include "ZONAL_FILT_OPTIONS.h"
                0002
                0003       SUBROUTINE ZONAL_FILT_NOFILL(
                0004      U           field,
                0005      I           jMin, jMax, kMin, kMax, bi, bj, gridLoc, myThid )
                0006 C     /==========================================================\
                0007 C     | S/R ZONAL_FILT_NOFILL                                    |
                0008 C     | o Apply FFT filter to a latitude circle.                 |
                0009 C     |   *No fill*                                              |
                0010 C     \==========================================================/
                0011       IMPLICIT NONE
                0012
                0013 C     == Global data ==
                0014 #include "SIZE.h"
                0015 #include "ZONAL_FILT.h"
                0016 #include "FFTPACK.h"
                0017
                0018 C     == Routine arguments ==
                0019 C     jMin - Range of points to filter
                0020 C     jMax
                0021 C     kMin
                0022 C     kMax
                0023 C     bi
                0024 C     bj
                0025 C     myThid  - Thread number of this instance of FILTER_LATCIRC_FFT_APPLY
                0026 C     field   - Field to filter
                0027 C     gridLoc - Orientation (U or V) of field.
                0028       INTEGER myThid
                0029       INTEGER gridLoc
                0030 C     Real*8 field(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
                0031       Real*8 field(1-OLx:sNx+OLx,1-OLy:sNy+OLy,1,nSx,nSy)
                0032       INTEGER jMin, jMax, kMin, kMax, bi, bj
                0033
                0034 #ifdef ALLOW_ZONAL_FILT
                0035
                0036 C     == Local data ==
                0037       Real*8 phi(Nx)
                0038       INTEGER I, J, K
                0039
                0040       DO k=kMin, kMax
                0041        DO j=jMin, jMax
                0042
                0043 C     o Copy zonal line of field into local workspace
                0044         DO i=1,sNx
                0045          phi(I) = field(i,j,k,bi,bj)
                0046         ENDDO
                0047
                0048 C     o Forward transform (using specific FFT package)
                0049 C       CALL R8FFTF( Nx, phi, FFTPACKWS(1,bj) )
                0050         CALL R8FFTF1( Nx, phi,
                0051      &    FFTPACKWS1(1,bj), FFTPACKWS2(1,bj),FFTPACKWS3(1,bj) )
                0052
                0053 C     o Apply amplitude filter and normalize
                0054         IF (gridLoc .EQ. 1) THEN
                0055          DO i=1, Nx
                0056           phi(i)=phi(i)*ampFactor(i,j,bi,bj)/float(Nx)
                0057          ENDDO
                0058         ELSEIF (gridLoc .EQ. 2) THEN
                0059          DO i=1, Nx
                0060           phi(i)=phi(i)*ampFactorV(i,j,bi,bj)/float(Nx)
                0061          ENDDO
                0062         ELSE
2cfc9d59a2 Patr*0063          WRITE(*,*) 'Error: gridLoc = ',gridLoc
42c525bfb4 Alis*0064          STOP 'Error: gridLoc has illegal value'
                0065         ENDIF
                0066
                0067 C     o Backward transform (using specific FFT package)
                0068 C       CALL R8FFTB( Nx, phi, FFTPACKWS(1,bj) )
                0069         CALL R8FFTB1( Nx, phi,
                0070      &    FFTPACKWS1(1,bj), FFTPACKWS2(1,bj),FFTPACKWS3(1,bj) )
                0071
                0072 C       o Do periodic wrap around by hand
                0073         DO i=1-OLx,0
                0074          field(i,j,k,bi,bj) = phi(sNx+i)
                0075         ENDDO
                0076         DO i=1,sNx
                0077          field(i,j,k,bi,bj) = phi(I)
                0078         ENDDO
                0079         DO i=sNx+1,sNx+OLx
                0080          field(i,j,k,bi,bj) = phi(i-sNx)
                0081         ENDDO
                0082
                0083        ENDDO
                0084       ENDDO
                0085
                0086 #endif /* ALLOW_ZONAL_FILT */
                0087
                0088       RETURN
                0089       END