trunk/filtrez/filtreg_scal.f

!
! $Header: /home/cvsroot/LMDZ4/libf/filtrez/filtreg.F,v 1.1.1.1 2004/05/19 12:53:09 lmdzadmin Exp $
!
      SUBROUTINE filtreg ( champ, nlat, nbniv, ifiltre,iaire,
     .   griscal ,iter)

      use dimens_m
      use paramet_m
             use parafilt
      IMPLICIT NONE
c=======================================================================
c
c   Auteur: P. Le Van        07/10/97
c   ------
c
c   Objet: filtre matriciel longitudinal ,avec les matrices precalculees
c                     pour l'operateur  Filtre    .
c   ------
c
c   Arguments:
c   ----------
c
c      nblat                 nombre de latitudes a filtrer
c      nbniv                 nombre de niveaux verticaux a filtrer
c      champ(iip1,nblat,nbniv)  en entree : champ a filtrer
c                            en sortie : champ filtre
c      ifiltre               +1  Transformee directe
c                            -1  Transformee inverse
c                            +2  Filtre directe
c                            -2  Filtre inverse
c
c      iaire                 1   si champ intensif
c                            2   si champ extensif (pondere par les aires)
c
c      iter                  1   filtre simple
c
c=======================================================================
c
c
c                      Variable Intensive
c                ifiltre = 1     filtre directe
c                ifiltre =-1     filtre inverse
c
c                      Variable Extensive
c                ifiltre = 2     filtre directe
c                ifiltre =-2     filtre inverse
c
c
      include "coefils.h"
c
      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 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)
      REAL  eignq(iim), sdd1(iim),sdd2(iim)
      LOGICAL    griscal
      INTEGER    hemisph, iaire
c

      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
c

      iim2   = iim * iim
      immjm  = iim * jjm
c
c
      IF( griscal )   THEN
         IF( nlat. NE. jjp1 )  THEN
             PRINT  1111
             STOP
         ELSE
c
             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
c
             jdfil1 = 2
             jffil1 = jfiltnu
             jdfil2 = jfiltsu
             jffil2 = jjm
          END IF
      ELSE
          IF( nlat.NE.jjm )  THEN
             PRINT  2222
             STOP
          ELSE
c
             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
c
             jdfil1 = 1
             jffil1 = jfiltnv
             jdfil2 = jfiltsv
             jffil2 = jjm
          END IF
      END IF
c
c
      DO 100  hemisph = 1, 2
c
      IF ( hemisph.EQ.1 )  THEN
          jdfil = jdfil1
          jffil = jffil1
      ELSE
          jdfil = jdfil2
          jffil = jffil2
      END IF

 
      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
      DO k = 1, iim
         eignq(k) = 0.0
      ENDDO
      DO k = 1, iim
      DO i = 1, iim
         eignq(k) = eignq(k) + matrinvn(k,i,j)*champ(i,j,l)
      ENDDO
      ENDDO
        ELSE IF ( griscal )     THEN
      DO k = 1, iim
         eignq(k) = 0.0
      ENDDO
      DO i = 1, iim
      DO k = 1, iim
         eignq(k) = eignq(k) + matriceun(k,i,j)*champ(i,j,l)
      ENDDO
      ENDDO
        ELSE 
      DO k = 1, iim
         eignq(k) = 0.0
      ENDDO
      DO i = 1, iim
      DO k = 1, iim
         eignq(k) = eignq(k) + matricevn(k,i,j)*champ(i,j,l)
      ENDDO
      ENDDO
        ENDIF

      ELSE

        IF( ifiltre. EQ. -2 )   THEN
      DO k = 1, iim
         eignq(k) = 0.0
      ENDDO
      DO i = 1, iim
      DO k = 1, iim
         eignq(k) = eignq(k) + matrinvs(k,i,j-jfiltsu+1)*champ(i,j,l)
      ENDDO
      ENDDO
        ELSE IF ( griscal )     THEN
      DO k = 1, iim
         eignq(k) = 0.0
      ENDDO
      DO i = 1, iim
      DO k = 1, iim
         eignq(k) = eignq(k) + matriceus(k,i,j-jfiltsu+1)*champ(i,j,l)
      ENDDO
      ENDDO
        ELSE 
      DO k = 1, iim
         eignq(k) = 0.0
      ENDDO
      DO i = 1, iim
      DO k = 1, iim
         eignq(k) = eignq(k) + matricevs(k,i,j-jfiltsv+1)*champ(i,j,l)
      ENDDO
      ENDDO
        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
c    
 100  CONTINUE
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'/)
      RETURN
      END
1	guez	3	!
2			! $Header: /home/cvsroot/LMDZ4/libf/filtrez/filtreg.F,v 1.1.1.1 2004/05/19 12:53:09 lmdzadmin Exp $
3			!
4			SUBROUTINE filtreg ( champ, nlat, nbniv, ifiltre,iaire,
5			. griscal ,iter)
6
7			use dimens_m
8			use paramet_m
9	guez	5	use parafilt
10	guez	3	IMPLICIT NONE
11			c=======================================================================
12			c
13			c Auteur: P. Le Van 07/10/97
14			c ------
15			c
16			c Objet: filtre matriciel longitudinal ,avec les matrices precalculees
17			c pour l'operateur Filtre .
18			c ------
19			c
20			c Arguments:
21			c ----------
22			c
23			c nblat nombre de latitudes a filtrer
24			c nbniv nombre de niveaux verticaux a filtrer
25			c champ(iip1,nblat,nbniv) en entree : champ a filtrer
26			c en sortie : champ filtre
27			c ifiltre +1 Transformee directe
28			c -1 Transformee inverse
29			c +2 Filtre directe
30			c -2 Filtre inverse
31			c
32			c iaire 1 si champ intensif
33			c 2 si champ extensif (pondere par les aires)
34			c
35			c iter 1 filtre simple
36			c
37			c=======================================================================
38			c
39			c
40			c Variable Intensive
41			c ifiltre = 1 filtre directe
42			c ifiltre =-1 filtre inverse
43			c
44			c Variable Extensive
45			c ifiltre = 2 filtre directe
46			c ifiltre =-2 filtre inverse
47			c
48			c
49			include "coefils.h"
50			c
51			INTEGER nlat,nbniv,ifiltre,iter
52			INTEGER i,j,l,k
53			INTEGER iim2,immjm
54			INTEGER jdfil1,jdfil2,jffil1,jffil2,jdfil,jffil
55
56			REAL champ( iip1,nlat,nbniv)
57			REAL matriceun,matriceus,matricevn,matricevs,matrinvn,matrinvs
58			COMMON/matrfil/matriceun(iim,iim,nfilun),matriceus(iim,iim,nfilus)
59			, , matricevn(iim,iim,nfilvn),matricevs(iim,iim,nfilvs)
60			, , matrinvn(iim,iim,nfilun),matrinvs (iim,iim,nfilus)
61			REAL eignq(iim), sdd1(iim),sdd2(iim)
62			LOGICAL griscal
63			INTEGER hemisph, iaire
64			c
65
66			IF(ifiltre.EQ.1.or.ifiltre.EQ.-1)
67			* STOP'Pas de transformee simple dans cette version'
68
69			IF( iter.EQ. 2 ) THEN
70			PRINT *,' Pas d iteration du filtre dans cette version !'
71			* , ' Utiliser old_filtreg et repasser !'
72			STOP
73			ENDIF
74
75			IF( ifiltre.EQ. -2 .AND..NOT.griscal ) THEN
76			PRINT *,' Cette routine ne calcule le filtre inverse que ',
77			* ' sur la grille des scalaires !'
78			STOP
79			ENDIF
80
81			IF( ifiltre.NE.2 .AND.ifiltre.NE. - 2 ) THEN
82			PRINT *,' Probleme dans filtreg car ifiltre NE 2 et NE -2'
83			*,' corriger et repasser !'
84			STOP
85			ENDIF
86			c
87
88			iim2 = iim * iim
89			immjm = iim * jjm
90			c
91			c
92			IF( griscal ) THEN
93			IF( nlat. NE. jjp1 ) THEN
94			PRINT 1111
95			STOP
96			ELSE
97			c
98			IF( iaire.EQ.1 ) THEN
99			CALL SCOPY( iim, sddv, 1, sdd1, 1 )
100			CALL SCOPY( iim, unsddv, 1, sdd2, 1 )
101			ELSE
102			CALL SCOPY( iim, unsddv, 1, sdd1, 1 )
103			CALL SCOPY( iim, sddv, 1, sdd2, 1 )
104			END IF
105			c
106			jdfil1 = 2
107			jffil1 = jfiltnu
108			jdfil2 = jfiltsu
109			jffil2 = jjm
110			END IF
111			ELSE
112			IF( nlat.NE.jjm ) THEN
113			PRINT 2222
114			STOP
115			ELSE
116			c
117			IF( iaire.EQ.1 ) THEN
118			CALL SCOPY( iim, sddu, 1, sdd1, 1 )
119			CALL SCOPY( iim, unsddu, 1, sdd2, 1 )
120			ELSE
121			CALL SCOPY( iim, unsddu, 1, sdd1, 1 )
122			CALL SCOPY( iim, sddu, 1, sdd2, 1 )
123			END IF
124			c
125			jdfil1 = 1
126			jffil1 = jfiltnv
127			jdfil2 = jfiltsv
128			jffil2 = jjm
129			END IF
130			END IF
131			c
132			c
133			DO 100 hemisph = 1, 2
134			c
135			IF ( hemisph.EQ.1 ) THEN
136			jdfil = jdfil1
137			jffil = jffil1
138			ELSE
139			jdfil = jdfil2
140			jffil = jffil2
141			END IF
142
143
144			DO 50 l = 1, nbniv
145			DO 30 j = jdfil,jffil
146
147
148			DO 5 i = 1, iim
149			champ(i,j,l) = champ(i,j,l) * sdd1(i)
150			5 CONTINUE
151			c
152
153			IF( hemisph. EQ. 1 ) THEN
154
155			IF( ifiltre. EQ. -2 ) THEN
156			DO k = 1, iim
157			eignq(k) = 0.0
158			ENDDO
159			DO k = 1, iim
160			DO i = 1, iim
161			eignq(k) = eignq(k) + matrinvn(k,i,j)*champ(i,j,l)
162			ENDDO
163			ENDDO
164			ELSE IF ( griscal ) THEN
165			DO k = 1, iim
166			eignq(k) = 0.0
167			ENDDO
168			DO i = 1, iim
169			DO k = 1, iim
170			eignq(k) = eignq(k) + matriceun(k,i,j)*champ(i,j,l)
171			ENDDO
172			ENDDO
173			ELSE
174			DO k = 1, iim
175			eignq(k) = 0.0
176			ENDDO
177			DO i = 1, iim
178			DO k = 1, iim
179			eignq(k) = eignq(k) + matricevn(k,i,j)*champ(i,j,l)
180			ENDDO
181			ENDDO
182			ENDIF
183
184			ELSE
185
186			IF( ifiltre. EQ. -2 ) THEN
187			DO k = 1, iim
188			eignq(k) = 0.0
189			ENDDO
190			DO i = 1, iim
191			DO k = 1, iim
192			eignq(k) = eignq(k) + matrinvs(k,i,j-jfiltsu+1)*champ(i,j,l)
193			ENDDO
194			ENDDO
195			ELSE IF ( griscal ) THEN
196			DO k = 1, iim
197			eignq(k) = 0.0
198			ENDDO
199			DO i = 1, iim
200			DO k = 1, iim
201			eignq(k) = eignq(k) + matriceus(k,i,j-jfiltsu+1)*champ(i,j,l)
202			ENDDO
203			ENDDO
204			ELSE
205			DO k = 1, iim
206			eignq(k) = 0.0
207			ENDDO
208			DO i = 1, iim
209			DO k = 1, iim
210			eignq(k) = eignq(k) + matricevs(k,i,j-jfiltsv+1)*champ(i,j,l)
211			ENDDO
212			ENDDO
213			ENDIF
214
215			ENDIF
216			c
217			IF( ifiltre.EQ. 2 ) THEN
218			DO 15 i = 1, iim
219			champ( i,j,l ) = ( champ(i,j,l) + eignq(i) ) * sdd2(i)
220			15 CONTINUE
221			ELSE
222			DO 16 i=1,iim
223			champ( i,j,l ) = ( champ(i,j,l) - eignq(i) ) * sdd2(i)
224			16 CONTINUE
225			ENDIF
226			c
227			champ( iip1,j,l ) = champ( 1,j,l )
228			c
229			30 CONTINUE
230			c
231			50 CONTINUE
232			c
233			100 CONTINUE
234			c
235			1111 FORMAT(//20x,'ERREUR dans le dimensionnement du tableau CHAMP a
236			*filtrer, sur la grille des scalaires'/)
237			2222 FORMAT(//20x,'ERREUR dans le dimensionnement du tableau CHAMP a fi
238			*ltrer, sur la grille de V ou de Z'/)
239			RETURN
240			END