! ! $Header$ ! SUBROUTINE filtreg ( champ, nlat, nbniv, ifiltre,iaire, & griscal ,iter) USE filtreg_mod USE dimensions_mod, ONLY: iim, jjm, llm, ndm USE paramet_mod_h USE coefils_mod_h IMPLICIT NONE !======================================================================= ! ! Auteur: P. Le Van 07/10/97 ! ------ ! ! Objet: filtre matriciel longitudinal ,avec les matrices precalculees ! pour l'operateur Filtre . ! ------ ! ! Arguments: ! ---------- ! ! nblat nombre de latitudes a filtrer ! nbniv nombre de niveaux verticaux a filtrer ! champ(iip1,nblat,nbniv) en entree : champ a filtrer ! en sortie : champ filtre ! ifiltre +1 Transformee directe ! -1 Transformee inverse ! +2 Filtre directe ! -2 Filtre inverse ! ! iaire 1 si champ intensif ! 2 si champ extensif (pondere par les aires) ! ! iter 1 filtre simple ! !======================================================================= ! ! ! Variable Intensive ! ifiltre = 1 filtre directe ! ifiltre =-1 filtre inverse ! ! Variable Extensive ! ifiltre = 2 filtre directe ! ifiltre =-2 filtre inverse ! INTEGER :: nlat,nbniv,ifiltre,iter INTEGER :: i,j,l,k INTEGER :: iim2,immjm INTEGER :: jdfil1,jdfil2,jffil1,jffil2,jdfil,jffil REAL :: champ( iip1,nlat,nbniv) REAL :: eignq(iim,nlat,nbniv), sdd1(iim),sdd2(iim) LOGICAL :: griscal INTEGER :: hemisph, iaire LOGICAL,SAVE :: first=.TRUE. REAL, SAVE :: sdd12(iim,4) INTEGER, PARAMETER :: type_sddu=1 INTEGER, PARAMETER :: type_sddv=2 INTEGER, PARAMETER :: type_unsddu=3 INTEGER, PARAMETER :: type_unsddv=4 INTEGER :: sdd1_type, sdd2_type if ( iim == 1 ) return ! no filtre in 2D y-z IF (first) THEN sdd12(1:iim,type_sddu) = sddu(1:iim) sdd12(1:iim,type_sddv) = sddv(1:iim) sdd12(1:iim,type_unsddu) = unsddu(1:iim) sdd12(1:iim,type_unsddv) = unsddv(1:iim) first=.FALSE. ENDIF IF(ifiltre.EQ.1.or.ifiltre.EQ.-1) & STOP 'Pas de transformee simple dans cette version' IF( iter.EQ. 2 ) THEN PRINT *,' Pas d iteration du filtre dans cette version !'& & , ' Utiliser old_filtreg et repasser !' STOP ENDIF IF( ifiltre.EQ. -2 .AND..NOT.griscal ) THEN PRINT *,' Cette routine ne calcule le filtre inverse que ' & , ' sur la grille des scalaires !' STOP ENDIF IF( ifiltre.NE.2 .AND.ifiltre.NE. - 2 ) THEN PRINT *,' Probleme dans filtreg car ifiltre NE 2 et NE -2' & , ' corriger et repasser !' STOP ENDIF iim2 = iim * iim immjm = iim * jjm IF( griscal ) THEN IF( nlat.NE. jjp1 ) THEN PRINT 1111 STOP ELSE IF( iaire.EQ.1 ) THEN sdd1_type = type_sddv sdd2_type = type_unsddv ELSE sdd1_type = type_unsddv sdd2_type = type_sddv ENDIF ! IF( iaire.EQ.1 ) THEN ! CALL SCOPY( iim, sddv, 1, sdd1, 1 ) ! CALL SCOPY( iim, unsddv, 1, sdd2, 1 ) ! ELSE ! CALL SCOPY( iim, unsddv, 1, sdd1, 1 ) ! CALL SCOPY( iim, sddv, 1, sdd2, 1 ) ! END IF jdfil1 = 2 jffil1 = jfiltnu jdfil2 = jfiltsu jffil2 = jjm END IF ELSE IF( nlat.NE.jjm ) THEN PRINT 2222 STOP ELSE IF( iaire.EQ.1 ) THEN sdd1_type = type_sddu sdd2_type = type_unsddu ELSE sdd1_type = type_unsddu sdd2_type = type_sddu ENDIF ! IF( iaire.EQ.1 ) THEN ! CALL SCOPY( iim, sddu, 1, sdd1, 1 ) ! CALL SCOPY( iim, unsddu, 1, sdd2, 1 ) ! ELSE ! CALL SCOPY( iim, unsddu, 1, sdd1, 1 ) ! CALL SCOPY( iim, sddu, 1, sdd2, 1 ) ! END IF jdfil1 = 1 jffil1 = jfiltnv jdfil2 = jfiltsv jffil2 = jjm END IF END IF DO hemisph = 1, 2 IF ( hemisph.EQ.1 ) THEN jdfil = jdfil1 jffil = jffil1 ELSE jdfil = jdfil2 jffil = jffil2 END IF DO l = 1, nbniv DO j = jdfil,jffil DO i = 1, iim champ(i,j,l) = champ(i,j,l) * sdd12(i,sdd1_type) ! sdd1(i) END DO END DO END DO IF( hemisph.EQ. 1 ) THEN IF( ifiltre.EQ. -2 ) THEN DO j = jdfil,jffil #ifdef BLAS CALL SGEMM("N", "N", iim, nbniv, iim, 1.0, & matrinvn(1,1,j), & iim, champ(1,j,1), iip1*nlat, 0.0, & eignq(1,j-jdfil+1,1), iim*nlat) #else eignq(:,j-jdfil+1,:) & = matmul(matrinvn(:,:,j), champ(:iim,j,:)) #endif END DO ELSE IF ( griscal ) THEN DO j = jdfil,jffil #ifdef BLAS CALL SGEMM("N", "N", iim, nbniv, iim, 1.0, & matriceun(1,1,j), & iim, champ(1,j,1), iip1*nlat, 0.0, & eignq(1,j-jdfil+1,1), iim*nlat) #else eignq(:,j-jdfil+1,:) & = matmul(matriceun(:,:,j), champ(:iim,j,:)) #endif END DO ELSE DO j = jdfil,jffil #ifdef BLAS CALL SGEMM("N", "N", iim, nbniv, iim, 1.0, & matricevn(1,1,j), & iim, champ(1,j,1), iip1*nlat, 0.0, & eignq(1,j-jdfil+1,1), iim*nlat) #else eignq(:,j-jdfil+1,:) & = matmul(matricevn(:,:,j), champ(:iim,j,:)) #endif END DO ENDIF ELSE IF( ifiltre.EQ. -2 ) THEN DO j = jdfil,jffil #ifdef BLAS CALL SGEMM("N", "N", iim, nbniv, iim, 1.0, & matrinvs(1,1,j-jfiltsu+1), & iim, champ(1,j,1), iip1*nlat, 0.0, & eignq(1,j-jdfil+1,1), iim*nlat) #else eignq(:,j-jdfil+1,:) & = matmul(matrinvs(:,:,j-jfiltsu+1), & champ(:iim,j,:)) #endif END DO ELSE IF ( griscal ) THEN DO j = jdfil,jffil #ifdef BLAS CALL SGEMM("N", "N", iim, nbniv, iim, 1.0, & matriceus(1,1,j-jfiltsu+1), & iim, champ(1,j,1), iip1*nlat, 0.0, & eignq(1,j-jdfil+1,1), iim*nlat) #else eignq(:,j-jdfil+1,:) & = matmul(matriceus(:,:,j-jfiltsu+1), & champ(:iim,j,:)) #endif END DO ELSE DO j = jdfil,jffil #ifdef BLAS CALL SGEMM("N", "N", iim, nbniv, iim, 1.0, & matricevs(1,1,j-jfiltsv+1), & iim, champ(1,j,1), iip1*nlat, 0.0, & eignq(1,j-jdfil+1,1), iim*nlat) #else eignq(:,j-jdfil+1,:) & = matmul(matricevs(:,:,j-jfiltsv+1), & champ(:iim,j,:)) #endif END DO ENDIF ENDIF IF( ifiltre.EQ. 2 ) THEN DO l = 1, nbniv DO j = jdfil,jffil DO i = 1, iim champ( i,j,l ) = & (champ(i,j,l) + eignq(i,j-jdfil+1,l)) & * sdd12(i,sdd2_type) ! sdd2(i) END DO END DO END DO ELSE DO l = 1, nbniv DO j = jdfil,jffil DO i = 1, iim champ( i,j,l ) = & (champ(i,j,l) - eignq(i,j-jdfil+1,l)) & * sdd12(i,sdd2_type) ! sdd2(i) END DO END DO END DO ENDIF DO l = 1, nbniv DO j = jdfil,jffil champ( iip1,j,l ) = champ( 1,j,l ) END DO END DO ENDDO 1111 FORMAT(//20x,'ERREUR dans le dimensionnement du tableau CHAMP a& & filtrer, sur la grille des scalaires'/) 2222 FORMAT(//20x,'ERREUR dans le dimensionnement du tableau CHAMP a fi& & ltrer, sur la grille de V ou de Z'/) RETURN END SUBROUTINE filtreg