Index: LMDZ4/branches/LMDZ4_AR5/libf/dyn3dpar/filtreg_p.F
===================================================================
--- LMDZ4/branches/LMDZ4_AR5/libf/dyn3dpar/filtreg_p.F	(revision 1599)
+++ LMDZ4/branches/LMDZ4_AR5/libf/dyn3dpar/filtreg_p.F	(revision 1601)
@@ -208,24 +208,39 @@
                IF( ifiltre.EQ.-2 )   THEN
                   DO j = jdfil,jffil
+#ifdef BLAS
                      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)
+#else
+                     champ_fft(:,j-jdfil+1,:)
+     &                    =matmul(matrinvn(:,:,j),champ_loc(:iim,j,:))
+#endif
                   ENDDO
                   
                ELSE IF ( griscal )     THEN
                   DO j = jdfil,jffil
+#ifdef BLAS
                      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)
+#else
+                     champ_fft(:,j-jdfil+1,:)
+     &                    =matmul(matriceun(:,:,j),champ_loc(:iim,j,:))
+#endif
                   ENDDO
                   
                ELSE 
                   DO j = jdfil,jffil
+#ifdef BLAS
                      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)
+#else
+                     champ_fft(:,j-jdfil+1,:)
+     &                    =matmul(matricevn(:,:,j),champ_loc(:iim,j,:))
+#endif
                   ENDDO
                   
@@ -236,8 +251,14 @@
                IF( ifiltre.EQ.-2 )   THEN
                   DO j = jdfil,jffil
+#ifdef BLAS
                      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)
+#else
+                     champ_fft(:,j-jdfil+1,:)
+     &                    =matmul(matrinvs(:,:,j-jfiltsu+1),
+     &                            champ_loc(:iim,j,:))
+#endif
                   ENDDO
                   
@@ -245,8 +266,14 @@
                   
                   DO j = jdfil,jffil
+#ifdef BLAS
                      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)
+#else
+                     champ_fft(:,j-jdfil+1,:)
+     &                    =matmul(matriceus(:,:,j-jfiltsu+1),
+     &                            champ_loc(:iim,j,:))
+#endif
                   ENDDO
                   
@@ -254,8 +281,14 @@
                   
                   DO j = jdfil,jffil
+#ifdef BLAS
                      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)
+#else
+                     champ_fft(:,j-jdfil+1,:)
+     &                    =matmul(matricevs(:,:,j-jfiltsv+1),
+     &                            champ_loc(:iim,j,:))
+#endif
                   ENDDO
                   
Index: LMDZ4/branches/LMDZ4_AR5/libf/filtrez/coefils.h
===================================================================
--- LMDZ4/branches/LMDZ4_AR5/libf/filtrez/coefils.h	(revision 1599)
+++ LMDZ4/branches/LMDZ4_AR5/libf/filtrez/coefils.h	(revision 1601)
@@ -1,4 +1,4 @@
 !
-! $Header$
+! $Id $
 !
       COMMON/coefils/jfiltnu,jfiltsu,jfiltnv,jfiltsv,sddu(iim),sddv(iim)&
@@ -7,5 +7,10 @@
      & ,coefilu2(iim,jjm),coefilv2(iim,jjm)
 !c
-      INTEGER jfiltnu,jfiltsu,jfiltnv,jfiltsv,modfrstu,modfrstv
+      INTEGER jfiltnu ! index of the last lat line filtered in NH (U grid)
+      INTEGER jfiltsu ! index of the first lat line filtered in SH (U grid)
+      INTEGER jfiltnv ! index of the last lat line filtered in NH (V grid)
+      INTEGER jfiltsv ! index of the first lat line filtered in SH (V grid)
+      INTEGER modfrstu ! number of retained (ie: unfiltered) modes on U grid
+      INTEGER modfrstv ! number of retained (ie: unfiltered) modes on V grid
       REAL    sddu,sddv,unsddu,unsddv,coefilu,coefilv,eignfnu,eignfnv
       REAL    coefilu2,coefilv2
Index: LMDZ4/branches/LMDZ4_AR5/libf/filtrez/filtreg_mod.F90
===================================================================
--- LMDZ4/branches/LMDZ4_AR5/libf/filtrez/filtreg_mod.F90	(revision 1599)
+++ LMDZ4/branches/LMDZ4_AR5/libf/filtrez/filtreg_mod.F90	(revision 1601)
@@ -1,2 +1,5 @@
+!
+! $Id $
+!
 MODULE filtreg_mod
 
@@ -42,5 +45,4 @@
     INTEGER ixmineq
 #endif
-    EXTERNAL  inifgn
     !
     ! ------------------------------------------------------------
@@ -71,7 +73,7 @@
     CALL inifgn(eignvl)
     !
-    PRINT *,' EIGNVL '
+    PRINT *,'inifilr: EIGNVL '
     PRINT 250,eignvl
-250 FORMAT( 1x,5e13.6)
+250 FORMAT( 1x,5e14.6)
     !
     ! compute eigenvalues and eigenfunctions
@@ -113,4 +115,10 @@
 #endif
     !
+    ! For a regular grid, we want the filter to start at latitudes
+    ! corresponding to lengths dx of the same size as dy (in terms
+    ! of angles: dx=2*dy) => at colat0=0.5 (i.e. colatitude=30 degrees
+    !  <=> latitude=60 degrees).
+    ! Same idea for the zoomed grid: start filtering polewards as soon
+    ! as length dx becomes of the same size as dy 
     !
     colat0  =  MIN( 0.5, dymin/dxmin )
@@ -158,16 +166,25 @@
     imx  = iim
     !
-    PRINT *,' TRUNCATION AT ',imx
-    !
+    PRINT *,'inifilr: TRUNCATION AT ',imx
+    !
+! Ehouarn: set up some defaults
+    jfiltnu=2 ! avoid north pole
+    jfiltsu=jjm ! avoid south pole (which is at jjm+1)
+    jfiltnv=1 ! NB: no poles on the V grid
+    jfiltsv=jjm
+
     DO j = 2, jjm/2+1
        cof = COS( rlatu(j) )/ colat0
        IF ( cof .LT. 1. ) THEN
-          IF( rlamda(imx) * COS(rlatu(j) ).LT.1. ) jfiltnu= j
+          IF( rlamda(imx) * COS(rlatu(j) ).LT.1. ) THEN
+            jfiltnu= j
+          ENDIF
        ENDIF
 
        cof = COS( rlatu(jjp1-j+1) )/ colat0
        IF ( cof .LT. 1. ) THEN
-          IF( rlamda(imx) * COS(rlatu(jjp1-j+1) ).LT.1. ) &
+          IF( rlamda(imx) * COS(rlatu(jjp1-j+1) ).LT.1. ) THEN
                jfiltsu= jjp1-j+1
+          ENDIF
        ENDIF
     ENDDO
@@ -176,16 +193,18 @@
        cof = COS( rlatv(j) )/ colat0
        IF ( cof .LT. 1. ) THEN
-          IF( rlamda(imx) * COS(rlatv(j) ).LT.1. ) jfiltnv= j
+          IF( rlamda(imx) * COS(rlatv(j) ).LT.1. ) THEN
+            jfiltnv= j
+          ENDIF
        ENDIF
 
        cof = COS( rlatv(jjm-j+1) )/ colat0
        IF ( cof .LT. 1. ) THEN
-          IF( rlamda(imx) * COS(rlatv(jjm-j+1) ).LT.1. ) &
+          IF( rlamda(imx) * COS(rlatv(jjm-j+1) ).LT.1. ) THEN
                jfiltsv= jjm-j+1
+          ENDIF
        ENDIF
     ENDDO
     !                                 
 
-    IF ( jfiltnu.LE.0 ) jfiltnu=1
     IF( jfiltnu.GT. jjm/2 +1 )  THEN
        PRINT *,' jfiltnu en dehors des valeurs acceptables ' ,jfiltnu
@@ -193,5 +212,4 @@
     ENDIF
 
-    IF( jfiltsu.LE.0) jfiltsu=1
     IF( jfiltsu.GT.  jjm  +1 )  THEN
        PRINT *,' jfiltsu en dehors des valeurs acceptables ' ,jfiltsu
@@ -199,5 +217,4 @@
     ENDIF
 
-    IF( jfiltnv.LE.0) jfiltnv=1
     IF( jfiltnv.GT. jjm/2    )  THEN
        PRINT *,' jfiltnv en dehors des valeurs acceptables ' ,jfiltnv
@@ -205,5 +222,4 @@
     ENDIF
 
-    IF( jfiltsv.LE.0) jfiltsv=1
     IF( jfiltsv.GT.     jjm  )  THEN
        PRINT *,' jfiltsv en dehors des valeurs acceptables ' ,jfiltsv
@@ -211,14 +227,14 @@
     ENDIF
 
-    PRINT *,' jfiltnv jfiltsv jfiltnu jfiltsu ' , &
+    PRINT *,'inifilr: jfiltnv jfiltsv jfiltnu jfiltsu ' , &
          jfiltnv,jfiltsv,jfiltnu,jfiltsu
 
     IF(first_call_inifilr) THEN
        ALLOCATE(matriceun(iim,iim,jfiltnu))
-       ALLOCATE(matriceus(iim,iim,jfiltsu))
+       ALLOCATE(matriceus(iim,iim,jjm-jfiltsu+1))
        ALLOCATE(matricevn(iim,iim,jfiltnv))
-       ALLOCATE(matricevs(iim,iim,jfiltsv))
+       ALLOCATE(matricevs(iim,iim,jjm-jfiltsv+1))
        ALLOCATE( matrinvn(iim,iim,jfiltnu))
-       ALLOCATE( matrinvs(iim,iim,jfiltsu))
+       ALLOCATE( matrinvs(iim,iim,jjm-jfiltsu+1))
        first_call_inifilr = .FALSE.
     ENDIF
@@ -230,4 +246,5 @@
     !
     DO j = 1,jjm
+    !default initialization: all modes are retained (i.e. no filtering)
        modfrstu( j ) = iim
        modfrstv( j ) = iim
@@ -306,5 +323,6 @@
 
     IF(jfiltnv.GE.jjm/2 .OR. jfiltnu.GE.jjm/2)THEN
-
+! Ehouarn: and what are these for??? Trying to handle a limit case
+!          where filters extend to and meet at the equator?
        IF(jfiltnv.EQ.jfiltsv)jfiltsv=1+jfiltnv
        IF(jfiltnu.EQ.jfiltsu)jfiltsu=1+jfiltnu
@@ -334,5 +352,5 @@
              eignft(i,k) = eignfnv(k,i) * coff
           ENDDO
-       ENDDO
+       ENDDO ! of DO i=1,iim
 #ifdef CRAY
        CALL MXM( eignfnv,iim,eignft,iim,matriceun(1,1,j),iim )
@@ -350,9 +368,9 @@
              ENDDO
           ENDDO
-       ENDDO
-#endif
-#endif
-
-    ENDDO
+       ENDDO ! of DO k = 1, iim
+#endif
+#endif
+
+    ENDDO ! of DO j = 2, jfiltnu
 
     DO j = jfiltsu, jjm
@@ -364,5 +382,5 @@
              eignft(i,k) = eignfnv(k,i) * coff
           ENDDO
-       ENDDO
+       ENDDO ! of DO i=1,iim
 #ifdef CRAY
        CALL MXM(eignfnv,iim,eignft,iim,matriceus(1,1,j-jfiltsu+1),iim)
@@ -381,9 +399,9 @@
              ENDDO
           ENDDO
-       ENDDO
-#endif
-#endif
-
-    ENDDO
+       ENDDO ! of DO k = 1, iim
+#endif
+#endif
+
+    ENDDO ! of DO j = jfiltsu, jjm
 
     !   ...................................................................
@@ -421,5 +439,5 @@
 #endif
 
-    ENDDO
+    ENDDO ! of DO j = 1, jfiltnv
 
     DO j = jfiltsv, jjm
@@ -452,5 +470,5 @@
 #endif
 
-    ENDDO
+    ENDDO ! of DO j = jfiltsv, jjm
 
     !   ...................................................................
@@ -488,5 +506,5 @@
 #endif
 
-    ENDDO
+    ENDDO ! of DO j = 2, jfiltnu
 
     DO j = jfiltsu, jjm
@@ -518,5 +536,5 @@
 #endif
 
-    ENDDO
+    ENDDO ! of DO j = jfiltsu, jjm
 
     IF (use_filtre_fft) THEN
Index: LMDZ4/branches/LMDZ4_AR5/libf/filtrez/mod_fft_fftw.F90
===================================================================
--- LMDZ4/branches/LMDZ4_AR5/libf/filtrez/mod_fft_fftw.F90	(revision 1599)
+++ LMDZ4/branches/LMDZ4_AR5/libf/filtrez/mod_fft_fftw.F90	(revision 1601)
@@ -1,39 +1,84 @@
+!
+! $Id$
+!
+
 MODULE mod_fft_fftw
 
 #ifdef FFT_FFTW
 
-  REAL,SAVE,ALLOCATABLE    :: Table_forward(:)
-  REAL,SAVE,ALLOCATABLE    :: Table_backward(:)
-  REAL,SAVE                :: scale_factor
-  INTEGER,SAVE             :: vsize
-  INTEGER,PARAMETER        :: inc=1
+  REAL, SAVE                   :: scale_factor
+  INTEGER, SAVE                :: vsize
+  INTEGER, PARAMETER           :: inc=1
   
-  INTEGER,SAVE             :: plan_forward
-  INTEGER,SAVE             :: plan_backward
+  INTEGER*8, ALLOCATABLE, DIMENSION(:), SAVE :: plan_forward
+  INTEGER*8, ALLOCATABLE, DIMENSION(:), SAVE :: plan_backward
   
 CONTAINS
   
-  SUBROUTINE Init_fft(iim)
+  SUBROUTINE Init_fft(iim,nvectmax)
   IMPLICIT NONE
-#include <rfftw.h>
+#include <fftw3.f>
     INTEGER :: iim
-    REAL    :: rtmp=1.
-    COMPLEX*16 :: ctmp
-    INTEGER :: itmp=1
-    INTEGER :: isign=0
-    INTEGER :: ierr
+    INTEGER :: nvectmax
+
+    INTEGER :: itmp
+
+    INTEGER               :: rank
+    INTEGER               :: howmany
+    INTEGER               :: istride, idist
+    INTEGER               :: ostride, odist
+    INTEGER, DIMENSION(1) :: n_array, inembed, onembed
+
+    REAL,    DIMENSION(iim+1,nvectmax) :: dbidon
+    COMPLEX, DIMENSION(iim/2+1,nvectmax) :: cbidon
+
+    vsize = iim
+    scale_factor = 1./SQRT(1.*vsize)
+
+    dbidon = 0
+    cbidon = 0
+
+    ALLOCATE(plan_forward(nvectmax))
+    ALLOCATE(plan_backward(nvectmax))
     
-    vsize=iim
-    scale_factor=1./SQRT(1.*vsize)
-    ALLOCATE(Table_forward(2*vsize+64))
-    ALLOCATE(Table_backward(2*vsize+64))
+    WRITE(*,*)"!---------------------!"
+    WRITE(*,*)"!                     !"
+    WRITE(*,*)"! INITIALISATION FFTW !"
+    WRITE(*,*)"!                     !"
+    WRITE(*,*)"!---------------------!"
     
-!    CALL DZFFTM(isign,vsize,itmp,scale_factor,rtmp,vsize+inc,ctmp,vsize/2+1,table_forward,rtmp,ierr)
-    
-!    CALL ZDFFTM(isign,vsize,itmp,scale_factor,ctmp,vsize/2+1,rtmp,vsize+inc,table_backward,rtmp,ierr)
+! On initialise tous les plans 
+    DO itmp = 1, nvectmax
+       rank       = 1
+       n_array(1) = iim
+       howmany    = itmp
+       inembed(1) = iim + 1 ; onembed(1) = iim/2 + 1
+       istride    = 1       ; ostride    = 1
+       idist      = iim + 1 ; odist      = iim/2 + 1
 
-    CALL rfftw_f77_create_plan(plan_forward,iim,FFTW_REAL_TO_COMPLEX,FFTW_ESTIMATE)
-    CALL rfftw_f77_create_plan(plan_backward,iim,FFTW_COMPLEX_TO_REAL,FFTW_ESTIMATE)
-    
+       CALL dfftw_plan_many_dft_r2c(plan_forward(itmp), rank, n_array, howmany, &
+            & dbidon, inembed, istride, idist, &
+            & cbidon, onembed, ostride, odist, &
+            & FFTW_ESTIMATE)
+
+       rank       = 1
+       n_array(1) = iim
+       howmany    = itmp
+       inembed(1) = iim/2 + 1 ; onembed(1) = iim + 1
+       istride    = 1         ; ostride    = 1
+       idist      = iim/2 + 1 ; odist      = iim + 1
+       CALL dfftw_plan_many_dft_c2r(plan_backward(itmp), rank, n_array, howmany, &
+            & cbidon, inembed, istride, idist, &
+            & dbidon, onembed, ostride, odist, &
+            & FFTW_ESTIMATE)
+
+    ENDDO
+
+    WRITE(*,*)"!-------------------------!"
+    WRITE(*,*)"!                         !"
+    WRITE(*,*)"! FIN INITIALISATION FFTW !"
+    WRITE(*,*)"!                         !"
+    WRITE(*,*)"!-------------------------!"
+
   END SUBROUTINE Init_fft
   
@@ -41,29 +86,25 @@
   SUBROUTINE fft_forward(vect,TF_vect,nb_vect)
     IMPLICIT NONE
-#include <rfftw.h>
-    INTEGER,INTENT(IN)  :: nb_vect
-    REAL,INTENT(IN)     :: vect(vsize+inc,nb_vect)
-    COMPLEX*16,INTENT(OUT) :: TF_vect(vsize/2+1,nb_vect)
-    REAL                :: work(4*vsize*nb_vect)
-    INTEGER             :: ierr
-    INTEGER, PARAMETER :: isign=-1
-  
-!    CALL DZFFTM(isign,vsize,nb_vect,scale_factor,vect,vsize+inc,TF_vect,vsize/2+1,table_forward,work,ierr)
-     CALL rfftwnd_f77_real_to_complex(plan_forward,nb_vect,vect, 1, vsize+inc , TF_vect, 1, vsize/2+1);  
-    
+#include <fftw3.f>
+    INTEGER,INTENT(IN)     :: nb_vect
+    REAL,INTENT(IN)        :: vect(vsize+inc,nb_vect)
+    COMPLEX,INTENT(OUT) :: TF_vect(vsize/2+1,nb_vect)
+
+    CALL dfftw_execute_dft_r2c(plan_forward(nb_vect),vect,TF_vect)
+
+    TF_vect = scale_factor * TF_vect
+
   END SUBROUTINE fft_forward
   
   SUBROUTINE fft_backward(TF_vect,vect,nb_vect)
     IMPLICIT NONE
-#include <rfftw.h>
-    INTEGER,INTENT(IN)  :: nb_vect
-    REAL,INTENT(OUT)    :: vect(vsize+inc,nb_vect)
-    COMPLEX*16,INTENT(IN ) :: TF_vect(vsize/2+1,nb_vect)
-    REAL                :: work(4*vsize*nb_vect)
-    INTEGER             :: ierr
-    INTEGER, PARAMETER :: isign=1
-  
-!    CALL ZDFFTM(isign,vsize,nb_vect,scale_factor,TF_vect,vsize/2+1,vect,vsize+inc,table_backward,work,ierr)
-    CALL rfftwnd_f77_complex_to_real(plan_forward,nb_vect,TF_vect, 1, vsize/2+1 , vect, 1, vsize+inc);  
+#include <fftw3.f>
+    INTEGER,INTENT(IN)     :: nb_vect
+    REAL,INTENT(OUT)       :: vect(vsize+inc,nb_vect)
+    COMPLEX,INTENT(IN ) :: TF_vect(vsize/2+1,nb_vect)
+
+    CALL dfftw_execute_dft_c2r(plan_backward(nb_vect),TF_vect,vect)
+
+    vect = scale_factor * vect
 
   END SUBROUTINE fft_backward
@@ -72,3 +113,2 @@
   
 END MODULE mod_fft_fftw
-
Index: LMDZ4/branches/LMDZ4_AR5/libf/filtrez/mod_fft_mathkeisan.F90
===================================================================
--- LMDZ4/branches/LMDZ4_AR5/libf/filtrez/mod_fft_mathkeisan.F90	(revision 1599)
+++ LMDZ4/branches/LMDZ4_AR5/libf/filtrez/mod_fft_mathkeisan.F90	(revision 1601)
@@ -15,5 +15,5 @@
     INTEGER :: nb_vect_max
     REAL    :: rtmp=1.
-    COMPLEX*16 :: ctmp
+    COMPLEX :: ctmp
     INTEGER :: itmp=1
     INTEGER :: isign=0
@@ -37,5 +37,5 @@
     INTEGER,INTENT(IN)  :: nb_vect
     REAL,INTENT(IN)     :: vect(vsize+inc,nb_vect)
-    COMPLEX*16,INTENT(OUT) :: TF_vect(vsize/2+1,nb_vect)
+    COMPLEX,INTENT(OUT) :: TF_vect(vsize/2+1,nb_vect)
     REAL                :: work(4*vsize*nb_vect)
     INTEGER             :: ierr
@@ -51,5 +51,5 @@
     INTEGER,INTENT(IN)  :: nb_vect
     REAL,INTENT(OUT)    :: vect(vsize+inc,nb_vect)
-    COMPLEX*16,INTENT(IN ) :: TF_vect(vsize/2+1,nb_vect)
+    COMPLEX,INTENT(IN ) :: TF_vect(vsize/2+1,nb_vect)
     REAL                :: work(4*vsize*nb_vect)
     INTEGER             :: ierr
Index: LMDZ4/branches/LMDZ4_AR5/libf/filtrez/mod_fft_mkl.F90
===================================================================
--- LMDZ4/branches/LMDZ4_AR5/libf/filtrez/mod_fft_mkl.F90	(revision 1599)
+++ LMDZ4/branches/LMDZ4_AR5/libf/filtrez/mod_fft_mkl.F90	(revision 1601)
@@ -24,5 +24,5 @@
     INTEGER :: nb_vect_max
     REAL    :: rtmp=1.
-    COMPLEX*16 :: ctmp
+    COMPLEX :: ctmp
     INTEGER :: itmp=1
     INTEGER :: isign=0
@@ -60,5 +60,5 @@
     INTEGER,INTENT(IN)  :: nb_vect
     REAL,INTENT(IN)     :: vect((vsize+inc)*nb_vect)
-    COMPLEX*16,INTENT(OUT) :: TF_vect((vsize/2+1)*nb_vect)
+    COMPLEX,INTENT(OUT) :: TF_vect((vsize/2+1)*nb_vect)
     REAL                :: work(4*vsize*nb_vect)
     INTEGER             :: ierr
@@ -102,5 +102,5 @@
     INTEGER,INTENT(IN)  :: nb_vect
     REAL,INTENT(OUT)    :: vect((vsize+inc)*nb_vect)
-    COMPLEX*16,INTENT(IN ) :: TF_vect((vsize/2+1)*nb_vect)
+    COMPLEX,INTENT(IN ) :: TF_vect((vsize/2+1)*nb_vect)
     REAL                :: work(4*vsize*nb_vect)
     INTEGER             :: ierr
Index: LMDZ4/branches/LMDZ4_AR5/libf/filtrez/mod_fft_wrapper.F90
===================================================================
--- LMDZ4/branches/LMDZ4_AR5/libf/filtrez/mod_fft_wrapper.F90	(revision 1599)
+++ LMDZ4/branches/LMDZ4_AR5/libf/filtrez/mod_fft_wrapper.F90	(revision 1601)
@@ -19,5 +19,5 @@
     INTEGER,INTENT(IN)  :: nb_vect
     REAL,INTENT(IN)     :: vect(vsize+inc,nb_vect)
-    COMPLEX*16,INTENT(INOUT) :: TF_vect(vsize/2+1,nb_vect)
+    COMPLEX,INTENT(INOUT) :: TF_vect(vsize/2+1,nb_vect)
     
     STOP "wrapper fft : une FFT doit etre specifiee a l'aide d'une clee CPP, sinon utiliser le filtre classique"
@@ -29,5 +29,5 @@
     INTEGER,INTENT(IN)  :: nb_vect
     REAL,INTENT(INOUT)    :: vect(vsize+inc,nb_vect)
-    COMPLEX*16,INTENT(IN ) :: TF_vect(vsize/2+1,nb_vect)
+    COMPLEX,INTENT(IN ) :: TF_vect(vsize/2+1,nb_vect)
   
     STOP "wrapper fft : une FFT doit etre specifiee a l'aide d'une clee CPP, sinon utiliser le filtre classique"
Index: LMDZ4/branches/LMDZ4_AR5/libf/filtrez/mod_filtre_fft.F90
===================================================================
--- LMDZ4/branches/LMDZ4_AR5/libf/filtrez/mod_filtre_fft.F90	(revision 1599)
+++ LMDZ4/branches/LMDZ4_AR5/libf/filtrez/mod_filtre_fft.F90	(revision 1601)
@@ -1,2 +1,6 @@
+!
+! $Id$
+!
+
 MODULE mod_filtre_fft
 
@@ -23,5 +27,6 @@
     INTEGER            :: index_vp(iim)
     INTEGER            :: i,j
-    
+    INTEGER            :: l,ll_nb
+
     index_vp(1)=1
     DO i=1,iim/2
@@ -98,7 +103,13 @@
     ENDDO
     
-    
+#ifdef FFT_FFTW
+
+    WRITE (*,*)"COTH jfiltnu,jfiltsu,jfiltnv,jjm-jfiltsv"
+    WRITE (*,*)jfiltnu,jfiltsu,jfiltnv,jjm-jfiltsv
+    WRITE (*,*)MAX(jfiltnu-2,jjm-jfiltsu,jfiltnv-2,jjm-jfiltsv)+1
+    CALL Init_FFT(iim,(llm+1)*(MAX(jfiltnu-2,jjm-jfiltsu,jfiltnv-2,jjm-jfiltsv)+1))
+#else    
     CALL Init_FFT(iim,(jjm+1)*(llm+1))
-        
+#endif        
     
   END SUBROUTINE Init_filtre_fft
@@ -118,11 +129,8 @@
 
     REAL               :: vect(iim+inc,jj_end-jj_begin+1,nbniv)
-!    REAL               :: vect_test(iim+inc,jj_end-jj_begin+1,nbniv)
-    COMPLEX*16         :: TF_vect(iim/2+1,jj_end-jj_begin+1,nbniv)
-!    COMPLEX*16         :: TF_vect_test(iim/2+1,jj_end-jj_begin+1,nbniv)
+    COMPLEX         :: TF_vect(iim/2+1,jj_end-jj_begin+1,nbniv)
     INTEGER            :: nb_vect
     INTEGER :: i,j,l
     INTEGER :: ll_nb
-!    REAL               :: vect_tmp(iim+inc,jj_end-jj_begin+1,nbniv)
     
     ll_nb=0
@@ -140,15 +148,5 @@
     nb_vect=(jj_end-jj_begin+1)*ll_nb
 
-!    vect_tmp=vect
-
     CALL FFT_forward(vect,TF_vect,nb_vect)
-
-!    CALL FFT_forward(vect,TF_vect_test,nb_vect)
-!      PRINT *,"XXXXXXXXXXXXX Filtre_u_FFT xxxxxxxxxxxx"
-!      DO j=1,jj_end-jj_begin+1
-!      DO i=1,iim/2+1
-!         PRINT *,"====",i,j,"----->",TF_vect_test(i,j,1)
-!       ENDDO
-!      ENDDO
 
     DO l=1,ll_nb
@@ -159,15 +157,8 @@
       ENDDO
     ENDDO
-       
+  
     CALL FFT_backward(TF_vect,vect,nb_vect)
-!    CALL FFT_backward(TF_vect_test,vect_test,nb_vect)
-          
-!      PRINT *,"XXXXXXXXXXXXX Filtre_u_FFT xxxxxxxxxxxx"
-!      DO j=1,jj_end-jj_begin+1
-!         DO i=1,iim
-!           PRINT *,"====",i,j,"----->",vect_test(i,j,1)
-!         ENDDO
-!      ENDDO
-
+      
+      
     ll_nb=0
 !$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
@@ -199,5 +190,5 @@
 
     REAL               :: vect(iim+inc,jj_end-jj_begin+1,nbniv)
-    COMPLEX*16         :: TF_vect(iim/2+1,jj_end-jj_begin+1,nbniv)
+    COMPLEX            :: TF_vect(iim/2+1,jj_end-jj_begin+1,nbniv)
     INTEGER            :: nb_vect
     INTEGER :: i,j,l
@@ -260,6 +251,6 @@
     REAL,INTENT(INOUT) :: vect_inout(iim+1,nlat,nbniv)
 
-    REAL               :: vect(iim+inc,jj_end-jj_begin+1,nbniv)
-    COMPLEX*16         :: TF_vect(iim/2+1,jj_end-jj_begin+1,nbniv)
+     REAL               :: vect(iim+inc,jj_end-jj_begin+1,nbniv)
+    COMPLEX            :: TF_vect(iim/2+1,jj_end-jj_begin+1,nbniv)
     INTEGER            :: nb_vect
     INTEGER :: i,j,l
@@ -305,37 +296,4 @@
 
   END SUBROUTINE Filtre_inv_fft  
-  
-  
-!  SUBROUTINE get_ll_index(nbniv,ll_index,ll_nb)
-!  IMPLICIT NONE
-!    INTEGER,INTENT(IN)  :: nbniv
-!    INTEGER,INTENT(OUT) :: ll_index(nbniv)
-!    INTEGER,INTENT(OUT) :: ll_nb
-!
-!    INTEGER :: l,ll_begin, ll_end
-!   INTEGER :: omp_rank,omp_size
-!   INTEGER :: OMP_GET_NUM_THREADS
-!   INTEGER :: omp_chunk
-!   EXTERNAL OMP_GET_NUM_THREADS
-!   INTEGER :: OMP_GET_THREAD_NUM
-!   EXTERNAL OMP_GET_THREAD_NUM
-!
-!   
-!   omp_size=OMP_GET_NUM_THREADS()
-!   omp_rank=OMP_GET_THREAD_NUM()    
-!   omp_chunk=nbniv/omp_size+min(1,MOD(nbniv,omp_size))
-!   
-!   ll_begin=omp_rank*OMP_CHUNK+1
-!   ll_nb=0
-!   DO WHILE (ll_begin<=nbniv)
-!     ll_end=min(ll_begin+OMP_CHUNK-1,nbniv)
-!     DO l=ll_begin,ll_end
-!       ll_nb=ll_nb+1
-!       ll_index(ll_nb)=l
-!     ENDDO
-!     ll_begin=ll_begin+omp_size*OMP_CHUNK
-!   ENDDO
-!  
-!  END SUBROUTINE get_ll_index
    
 END MODULE mod_filtre_fft
