Index: LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/filtreg_p.F
===================================================================
--- LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/filtreg_p.F	(revision 1074)
+++ LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/filtreg_p.F	(revision 1086)
@@ -2,10 +2,13 @@
 
       SUBROUTINE filtreg_p ( champ, ibeg, iend, nlat, nbniv, 
-     .                       ifiltre, iaire, griscal ,iter)
-      USE Parallel, only : OMP_CHUNK 
+     &     ifiltre, iaire, griscal ,iter)
+      USE Parallel, only : OMP_CHUNK
       USE mod_filtre_fft
       USE timer_filtre
+      
+      USE filtre
+      
       IMPLICIT NONE
-
+      
 c=======================================================================
 c
@@ -50,5 +53,4 @@
 #include "dimensions.h"
 #include "paramet.h"
-#include "parafilt.h"
 #include "coefils.h"
 c
@@ -57,14 +59,6 @@
       INTEGER iim2,immjm
       INTEGER jdfil1,jdfil2,jffil1,jffil2,jdfil,jffil
-
+      
       REAL  champ( iip1,nlat,nbniv)
-      REAL matriceun,matriceus,matricevn,matricevs,matrinvn,matrinvs
-      COMMON/matrfil/matriceun(iim,iim,nfilun),matriceus(iim,iim,nfilus)
-     ,             , matricevn(iim,iim,nfilvn),matricevs(iim,iim,nfilvs)
-     ,             ,  matrinvn(iim,iim,nfilun),matrinvs (iim,iim,nfilus)
-cym      REAL  eignq(iim), sdd1(iim),sdd2(iim)
-
-      REAL  eignq(iim)
-      REAL :: sdd1(iim),sdd2(iim)
       
       LOGICAL    griscal
@@ -74,23 +68,27 @@
       REAL :: champ_in(iip1,nlat,nbniv)
       
-      REAL,SAVE,TARGET :: sddu_loc(iim)
-      REAL,SAVE,TARGET :: sddv_loc(iim)
-      REAL,SAVE,TARGET :: unsddu_loc(iim)
-      REAL,SAVE,TARGET :: unsddv_loc(iim)
-c$OMP THREADPRIVATE(sddu_loc,sddv_loc,unsddu_loc,unsddv_loc)
       LOGICAL,SAVE     :: first=.TRUE.
 c$OMP THREADPRIVATE(first) 
 
+      REAL, DIMENSION(iip1,nlat,nbniv) :: champ_loc
+      INTEGER :: ll_nb, nbniv_loc
+      REAL, SAVE :: sdd12(iim,4)
+c$OMP THREADPRIVATE(sdd12) 
+
+      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 (first) THEN
-        sddu_loc(1:iim)=sddu(1:iim)
-        sddv_loc(1:iim)=sddv(1:iim)
-        unsddu_loc(1:iim)=unsddu(1:iim)
-        unsddv_loc(1:iim)=unsddv(1:iim)
-        CALL Init_timer
-        first=.FALSE.
-c	PRINT *,"----> sddu_loc=",sddu_loc
-c	PRINT *,"----> sddv_loc=",sddv_loc
-c	PRINT *,"----> unsddu_loc=",unsddu_loc
-c	PRINT *,"----> unsddv_loc=",unsddv_loc
+         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)
+
+         CALL Init_timer
+         first=.FALSE.
       ENDIF
 
@@ -99,23 +97,25 @@
 c$OMP END MASTER
 
+c-------------------------------------------------------c
+
       IF(ifiltre.EQ.1.or.ifiltre.EQ.-1) 
-     *    STOP'Pas de transformee simple dans cette version'
-
+     &     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
+         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
+         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
+         PRINT *,' Probleme dans filtreg car ifiltre NE 2 et NE -2'
+     &        , ' corriger et repasser !'
+         STOP
       ENDIF
 c
@@ -127,77 +127,53 @@
       IF( griscal )   THEN
          IF( nlat. NE. jjp1 )  THEN
-             PRINT  1111
-             STOP
+            PRINT  1111
+            STOP
          ELSE
-c
-             IF( iaire.EQ.1 )  THEN
-cym                CALL SCOPY(  iim,    sddv, 1,  sdd1, 1 ) 
-cym                CALL SCOPY(  iim,  unsddv, 1,  sdd2, 1 )
-cym               sdd1=>sddv_loc
-cym               sdd2=>unsddv_loc
-               sdd1(1:iim)=sddv_loc(1:iim)
-               sdd2(1:iim)=unsddv_loc(1:iim)
-             ELSE
-cym                CALL SCOPY(  iim,  unsddv, 1,  sdd1, 1 )
-cym                CALL SCOPY(  iim,    sddv, 1,  sdd2, 1 )
-               sdd1(1:iim)=unsddv_loc(1:iim)
-               sdd2(1:iim)=sddv_loc(1:iim)
-             END IF
-c
-             jdfil1 = 2
-             jffil1 = jfiltnu
-             jdfil2 = jfiltsu
-             jffil2 = jjm
-          END IF
+c     
+            IF( iaire.EQ.1 )  THEN
+               sdd1_type = type_sddv
+               sdd2_type = type_unsddv
+            ELSE
+               sdd1_type = type_unsddv
+               sdd2_type = type_sddv
+            ENDIF
+c
+            jdfil1 = 2
+            jffil1 = jfiltnu
+            jdfil2 = jfiltsu
+            jffil2 = jjm
+         ENDIF
       ELSE
-          IF( nlat.NE.jjm )  THEN
-             PRINT  2222
-             STOP
-          ELSE
-c
-             IF( iaire.EQ.1 )  THEN
-cym                CALL SCOPY(  iim,    sddu, 1,  sdd1, 1 ) 
-cym                CALL SCOPY(  iim,  unsddu, 1,  sdd2, 1 )
-cym                sdd1=>sddu_loc
-cym                sdd2=>unsddu_loc
-                sdd1(1:iim)=sddu_loc(1:iim)
-                sdd2(1:iim)=unsddu_loc(1:iim)
-
-             ELSE
-cym                CALL SCOPY(  iim,  unsddu, 1,  sdd1, 1 )
-cym                CALL SCOPY(  iim,    sddu, 1,  sdd2, 1 )
-cym               sdd1=>unsddu_loc
-cym               sdd2=>sddu_loc
-	       sdd1(1:iim)=unsddu_loc(1:iim)
-               sdd2(1:iim)=sddu_loc(1:iim)
-             END IF
-c
-             jdfil1 = 1
-             jffil1 = jfiltnv
-             jdfil2 = jfiltsv
-             jffil2 = jjm
-          END IF
-      END IF
-
-c      PRINT *,"APPEL a filtreg --> sdd1=",sdd1
-c      PRINT *,"APPEL a filtreg --> sdd2=",sdd2
-c      PRINT *,"----> sddu_loc=",sddu_loc
-c	PRINT *,"----> sddv_loc=",sddv_loc
-c	PRINT *,"----> unsddu_loc=",unsddu_loc
-c	PRINT *,"----> unsddv_loc=",unsddv_loc
- 
-c
-c
-      DO 100  hemisph = 1, 2
-c
-      IF ( hemisph.EQ.1 )  THEN
-c ym
-          jdfil = max(jdfil1,ibeg)
-          jffil = min(jffil1,iend)
-      ELSE
-c ym
-          jdfil = max(jdfil2,ibeg)
-          jffil = min(jffil2,iend)
-      END IF
+         IF( nlat.NE.jjm )  THEN
+            PRINT  2222
+            STOP
+         ELSE
+c
+            IF( iaire.EQ.1 )  THEN
+               sdd1_type = type_sddu
+               sdd2_type = type_unsddu
+            ELSE
+               sdd1_type = type_unsddu
+               sdd2_type = type_sddu
+            ENDIF
+c     
+            jdfil1 = 1
+            jffil1 = jfiltnv
+            jdfil2 = jfiltsv
+            jffil2 = jjm
+         ENDIF
+      ENDIF
+c      
+      DO hemisph = 1, 2
+c     
+         IF ( hemisph.EQ.1 )  THEN
+cym
+            jdfil = max(jdfil1,ibeg)
+            jffil = min(jffil1,iend)
+         ELSE
+cym
+            jdfil = max(jdfil2,ibeg)
+            jffil = min(jffil2,iend)
+         ENDIF
 
 
@@ -206,144 +182,199 @@
 cccccccccccccccccccccccccccccccccccccccccccc
 
-      IF (.NOT. use_filtre_fft) THEN
-      
-c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)    
-      DO 50  l = 1, nbniv
-        DO 30  j = jdfil,jffil
- 
- 
-          DO  5  i = 1, iim
-            champ(i,j,l) = champ(i,j,l) * sdd1(i)
-   5      CONTINUE
-c
-
-          IF( hemisph. EQ. 1 )      THEN
-
-            IF( ifiltre. EQ. -2 )   THEN
-
-
-              CALL SGEMV("N", iim,iim, 1.0, matrinvn(1,1,j),iim,
-     .                     champ(1,j,l), 1, 0.0, eignq, 1)
-
-
-            ELSE IF ( griscal )     THEN
-
-              CALL SGEMV("N", iim,iim, 1.0, matriceun(1,1,j),iim,
-     .                    champ(1,j,l), 1, 0.0, eignq, 1)
-
-            ELSE 
-
-              CALL SGEMV("N", iim,iim, 1.0, matricevn(1,1,j),iim,
-     .                   champ(1,j,l), 1, 0.0, eignq, 1)
-            ENDIF
-
-          ELSE
-
-            IF( ifiltre. EQ. -2 )   THEN
-      
-              CALL SGEMV("N",iim,iim,1.0, matrinvs(1,1,j-jfiltsu+1),iim,
-     .                   champ(1,j,l), 1, 0.0, eignq, 1)
-      
-            ELSE IF ( griscal )     THEN
-      
-              CALL SGEMV("N",iim,iim,1.0,matriceus(1,1,j-jfiltsu+1),iim,
-     .                   champ(1,j,l), 1, 0.0, eignq, 1)
-            ELSE 
-          
-	      CALL SGEMV("N",iim,iim,1.0,matricevs(1,1,j-jfiltsv+1),iim,
-     .                    champ(1,j,l), 1, 0.0, eignq, 1)
-            ENDIF
-
-          ENDIF
-
-
-c
-          IF( ifiltre.EQ. 2 )  THEN
-         
-	    DO 15 i = 1, iim
-              champ( i,j,l ) = ( champ(i,j,l) + eignq(i) ) * sdd2(i)
-  15        CONTINUE
-         
-	  ELSE
-        
-	    DO 16 i=1,iim
-               champ( i,j,l ) = ( champ(i,j,l) - eignq(i) ) * sdd2(i)
-16          CONTINUE
-          
-	  ENDIF
-c
-          champ( iip1,j,l ) = champ( 1,j,l )
-c
-  30    CONTINUE
-c
-  50  CONTINUE
+         IF (.NOT. use_filtre_fft) THEN
+      
+c     !---------------------------------!
+c     ! Agregation des niveau verticaux !
+c     ! uniquement necessaire pour une  !
+c     ! execution OpenMP                !
+c     !---------------------------------!
+            ll_nb = 0
+c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
+            DO l = 1, nbniv
+               ll_nb = ll_nb+1
+               DO j = jdfil,jffil
+                  DO i = 1, iim
+                     champ_loc(i,j,ll_nb) = 
+     &                    champ(i,j,l) * sdd12(i,sdd1_type)
+                  ENDDO
+               ENDDO
+            ENDDO
 c$OMP END DO NOWAIT
 
+            nbniv_loc = ll_nb
+
+            IF( hemisph.EQ.1 )      THEN
+               
+               IF( ifiltre.EQ.-2 )   THEN
+                  DO j = jdfil,jffil
+                     CALL DGEMM("N", "N", iim, nbniv_loc, iim, 1.0, 
+     &                    matrinvn(1,1,j), iim, 
+     &                    champ_loc(1,j,1), iip1*nlat, 0.0,
+     &                    champ_fft(1,j-jdfil+1,1), iip1*nlat)
+                  ENDDO
+                  
+               ELSE IF ( griscal )     THEN
+                  DO j = jdfil,jffil
+                     CALL DGEMM("N", "N", iim, nbniv_loc, iim, 1.0, 
+     &                    matriceun(1,1,j), iim, 
+     &                    champ_loc(1,j,1), iip1*nlat, 0.0,
+     &                    champ_fft(1,j-jdfil+1,1), iip1*nlat)
+                  ENDDO
+                  
+               ELSE 
+                  DO j = jdfil,jffil
+                     CALL DGEMM("N", "N", iim, nbniv_loc, iim, 1.0, 
+     &                    matricevn(1,1,j), iim, 
+     &                    champ_loc(1,j,1), iip1*nlat, 0.0, 
+     &                    champ_fft(1,j-jdfil+1,1), iip1*nlat)
+                  ENDDO
+                  
+               ENDIF
+               
+            ELSE
+               
+               IF( ifiltre.EQ.-2 )   THEN
+                  DO j = jdfil,jffil
+                     CALL DGEMM("N", "N", iim, nbniv_loc, iim, 1.0, 
+     &                    matrinvs(1,1,j-jfiltsu+1), iim, 
+     &                    champ_loc(1,j,1), iip1*nlat, 0.0, 
+     &                    champ_fft(1,j-jdfil+1,1), iip1*nlat)
+                  ENDDO
+                  
+               ELSE IF ( griscal )     THEN
+                  
+                  DO j = jdfil,jffil
+                     CALL DGEMM("N", "N", iim, nbniv_loc, iim, 1.0, 
+     &                    matriceus(1,1,j-jfiltsu+1), iim, 
+     &                    champ_loc(1,j,1), iip1*nlat, 0.0, 
+     &                    champ_fft(1,j-jdfil+1,1), iip1*nlat)
+                  ENDDO
+                  
+               ELSE 
+                  
+                  DO j = jdfil,jffil
+                     CALL DGEMM("N", "N", iim, nbniv_loc, iim, 1.0, 
+     &                    matricevs(1,1,j-jfiltsv+1), iim, 
+     &                    champ_loc(1,j,1), iip1*nlat, 0.0, 
+     &                    champ_fft(1,j-jdfil+1,1), iip1*nlat)
+                  ENDDO
+                  
+               ENDIF
+               
+            ENDIF
+!     c     
+            IF( ifiltre.EQ.2 )  THEN
+               
+c     !-------------------------------------!
+c     ! Dés-agregation des niveau verticaux !
+c     ! uniquement necessaire pour une      !
+c     ! execution OpenMP                    !
+c     !-------------------------------------!
+               ll_nb = 0
+c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
+               DO l = 1, nbniv
+                  ll_nb = ll_nb + 1
+                  DO j = jdfil,jffil
+                     DO i = 1, iim
+                        champ( i,j,l ) = (champ_loc(i,j,ll_nb) 
+     &                       + champ_fft(i,j-jdfil+1,ll_nb))
+     &                       * sdd12(i,sdd2_type)
+                     ENDDO
+                  ENDDO
+               ENDDO
+c$OMP END DO NOWAIT
+               
+            ELSE
+               
+               ll_nb = 0
+c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
+               DO l = 1, nbniv_loc
+                  ll_nb = ll_nb + 1
+                  DO j = jdfil,jffil
+                     DO i = 1, iim
+                        champ( i,j,l ) = (champ_loc(i,j,ll_nb) 
+     &                       - champ_fft(i,j-jdfil+1,ll_nb))
+     &                       * sdd12(i,sdd2_type)
+                     ENDDO
+                  ENDDO
+               ENDDO
+c$OMP END DO NOWAIT
+               
+            ENDIF
+            
+c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
+            DO l = 1, nbniv
+               DO j = jdfil,jffil
+                  champ( iip1,j,l ) = champ( 1,j,l )
+               ENDDO
+            ENDDO
+c$OMP END DO NOWAIT
+            
 ccccccccccccccccccccccccccccccccccccccccccccc
 c Utilisation du filtre FFT
 ccccccccccccccccccccccccccccccccccccccccccccc
         
-       ELSE
+         ELSE
        
 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
-	  DO l=1,nbniv
-	    DO j=jdfil,jffil
-	      DO  i = 1, iim
-                champ( i,j,l)= champ(i,j,l)*sdd1(i)
-		champ_fft( i,j,l) = champ(i,j,l)
-              ENDDO
+            DO l=1,nbniv
+               DO j=jdfil,jffil
+                  DO  i = 1, iim
+                     champ( i,j,l)= champ(i,j,l)*sdd12(i,sdd1_type)
+                     champ_fft( i,j,l) = champ(i,j,l)
+                  ENDDO
+               ENDDO
             ENDDO
-	  ENDDO
 c$OMP END DO NOWAIT
 
-      IF (jdfil<=jffil) THEN
-        IF( ifiltre. EQ. -2 )   THEN
-          CALL Filtre_inv_fft(champ_fft,nlat,jdfil,jffil,nbniv) 
-        ELSE IF ( griscal )     THEN
-	  CALL Filtre_u_fft(champ_fft,nlat,jdfil,jffil,nbniv)
-	ELSE
-	  CALL Filtre_v_fft(champ_fft,nlat,jdfil,jffil,nbniv)
-	ENDIF
-      ENDIF
-
-
-        IF( ifiltre.EQ. 2 )  THEN
+            IF (jdfil<=jffil) THEN
+               IF( ifiltre. EQ. -2 )   THEN
+                  CALL Filtre_inv_fft(champ_fft,nlat,jdfil,jffil,nbniv) 
+               ELSE IF ( griscal )     THEN
+                  CALL Filtre_u_fft(champ_fft,nlat,jdfil,jffil,nbniv)
+               ELSE
+                  CALL Filtre_v_fft(champ_fft,nlat,jdfil,jffil,nbniv)
+               ENDIF
+            ENDIF
+
+
+            IF( ifiltre.EQ. 2 )  THEN
 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)          
-	  DO l=1,nbniv
-	    DO j=jdfil,jffil
-	      DO  i = 1, iim
-                champ( i,j,l)=(champ(i,j,l)+champ_fft(i,j,l))
-     .                             *sdd2(i)
-              ENDDO
-            ENDDO
-	  ENDDO
+               DO l=1,nbniv
+                  DO j=jdfil,jffil
+                     DO  i = 1, iim
+                        champ( i,j,l)=(champ(i,j,l)+champ_fft(i,j,l))
+     &                       *sdd12(i,sdd2_type)
+                     ENDDO
+                  ENDDO
+               ENDDO
 c$OMP END DO NOWAIT	  
-	ELSE
+            ELSE
         
 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)  
-	  DO l=1,nbniv
-	    DO j=jdfil,jffil
-	      DO  i = 1, iim
-                champ(i,j,l)=(champ(i,j,l)-champ_fft(i,j,l))
-     .                            *sdd2(i)
-              ENDDO
+               DO l=1,nbniv
+                  DO j=jdfil,jffil
+                     DO  i = 1, iim
+                        champ(i,j,l)=(champ(i,j,l)-champ_fft(i,j,l))
+     &                       *sdd12(i,sdd2_type)
+                     ENDDO
+                  ENDDO
+               ENDDO
+c$OMP END DO NOWAIT          
+            ENDIF
+c
+c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 
+            DO l=1,nbniv
+               DO j=jdfil,jffil
+!            champ_FFT( iip1,j,l ) = champ_FFT( 1,j,l )
+                  champ( iip1,j,l ) = champ( 1,j,l )
+               ENDDO
             ENDDO
-	  ENDDO
-c$OMP END DO NOWAIT          
-        ENDIF
-c
-c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 
-	DO l=1,nbniv
-	  DO j=jdfil,jffil
-!            champ_FFT( iip1,j,l ) = champ_FFT( 1,j,l )
-	    champ( iip1,j,l ) = champ( 1,j,l )
-          ENDDO
-	ENDDO
 c$OMP END DO NOWAIT          	
-      ENDIF 
+         ENDIF 
 c Fin de la zone de filtrage
 
 	
- 100  CONTINUE
+      ENDDO
 
 !      DO j=1,nlat
@@ -359,8 +390,8 @@
       
 c
-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'/)
+ 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'/)
 c$OMP MASTER      
       CALL stop_timer
