Head

GCC Code Coverage Report

Directory:	./		Exec	Total	Coverage
File:	filtrez/filtreg.F	Lines:	60	74	81.1 %
Date:	2023-06-30 12:56:34	Branches:	156	200	78.0 %


!
! $Header$
!
      SUBROUTINE filtreg ( champ, nlat, nbniv, ifiltre,iaire,
     &     griscal ,iter)

      USE filtreg_mod

      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 "dimensions.h"
      INCLUDE "paramet.h"
      INCLUDE "coefils.h"

      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

c            IF( iaire.EQ.1 )  THEN
c               CALL SCOPY(  iim,    sddv, 1,  sdd1, 1 )
c               CALL SCOPY(  iim,  unsddv, 1,  sdd2, 1 )
c            ELSE
c               CALL SCOPY(  iim,  unsddv, 1,  sdd1, 1 )
c               CALL SCOPY(  iim,    sddv, 1,  sdd2, 1 )
c            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

c            IF( iaire.EQ.1 )  THEN
c               CALL SCOPY(  iim,    sddu, 1,  sdd1, 1 )
c               CALL SCOPY(  iim,  unsddu, 1,  sdd2, 1 )
c            ELSE
c               CALL SCOPY(  iim,  unsddu, 1,  sdd1, 1 )
c               CALL SCOPY(  iim,    sddu, 1,  sdd2, 1 )
c            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


Generated by: GCOVR (Version 4.2)

Line	Branch	Exec	Source
1			!
2			! $Header$
3			!
4		77664	SUBROUTINE filtreg ( champ, nlat, nbniv, ifiltre,iaire,
5			& griscal ,iter)
6
7			USE filtreg_mod
8
9			IMPLICIT NONE
10			c=======================================================================
11			c
12			c Auteur: P. Le Van 07/10/97
13			c ------
14			c
15			c Objet: filtre matriciel longitudinal ,avec les matrices precalculees
16			c pour l'operateur Filtre .
17			c ------
18			c
19			c Arguments:
20			c ----------
21			c
22			c nblat nombre de latitudes a filtrer
23			c nbniv nombre de niveaux verticaux a filtrer
24			c champ(iip1,nblat,nbniv) en entree : champ a filtrer
25			c en sortie : champ filtre
26			c ifiltre +1 Transformee directe
27			c -1 Transformee inverse
28			c +2 Filtre directe
29			c -2 Filtre inverse
30			c
31			c iaire 1 si champ intensif
32			c 2 si champ extensif (pondere par les aires)
33			c
34			c iter 1 filtre simple
35			c
36			c=======================================================================
37			c
38			c
39			c Variable Intensive
40			c ifiltre = 1 filtre directe
41			c ifiltre =-1 filtre inverse
42			c
43			c Variable Extensive
44			c ifiltre = 2 filtre directe
45			c ifiltre =-2 filtre inverse
46			c
47			c
48			INCLUDE "dimensions.h"
49			INCLUDE "paramet.h"
50			INCLUDE "coefils.h"
51
52			INTEGER nlat,nbniv,ifiltre,iter
53			INTEGER i,j,l,k
54			INTEGER iim2,immjm
55			INTEGER jdfil1,jdfil2,jffil1,jffil2,jdfil,jffil
56
57			REAL champ( iip1,nlat,nbniv)
58
59		25888	REAL eignq(iim,nlat,nbniv), sdd1(iim),sdd2(iim)
60			LOGICAL griscal
61			INTEGER hemisph, iaire
62
63			LOGICAL,SAVE :: first=.TRUE.
64
65			REAL, SAVE :: sdd12(iim,4)
66
67			INTEGER, PARAMETER :: type_sddu=1
68			INTEGER, PARAMETER :: type_sddv=2
69			INTEGER, PARAMETER :: type_unsddu=3
70			INTEGER, PARAMETER :: type_unsddv=4
71
72			INTEGER :: sdd1_type, sdd2_type
73
74			if ( iim == 1 ) return ! no filtre in 2D y-z
75
76	✓✓	12944	IF (first) THEN
77	✓✓	33	sdd12(1:iim,type_sddu) = sddu(1:iim)
78	✓✓	33	sdd12(1:iim,type_sddv) = sddv(1:iim)
79	✓✓	33	sdd12(1:iim,type_unsddu) = unsddu(1:iim)
80	✓✓	33	sdd12(1:iim,type_unsddv) = unsddv(1:iim)
81
82		1	first=.FALSE.
83			ENDIF
84
85	✗✓	12944	IF(ifiltre.EQ.1.or.ifiltre.EQ.-1)
86			& STOP'Pas de transformee simple dans cette version'
87
88	✗✓	12944	IF( iter.EQ. 2 ) THEN
89			PRINT *,' Pas d iteration du filtre dans cette version !'
90			& , ' Utiliser old_filtreg et repasser !'
91			STOP
92			ENDIF
93
94	✓✓✗✓	12944	IF( ifiltre.EQ. -2 .AND..NOT.griscal ) THEN
95			PRINT *,' Cette routine ne calcule le filtre inverse que '
96			& , ' sur la grille des scalaires !'
97			STOP
98			ENDIF
99
100	✗✓	12944	IF( ifiltre.NE.2 .AND.ifiltre.NE. - 2 ) THEN
101			PRINT *,' Probleme dans filtreg car ifiltre NE 2 et NE -2'
102			& , ' corriger et repasser !'
103			STOP
104			ENDIF
105
106			iim2 = iim * iim
107			immjm = iim * jjm
108
109	✓✓	12944	IF( griscal ) THEN
110	✗✓	9861	IF( nlat. NE. jjp1 ) THEN
111			PRINT 1111
112			STOP
113			ELSE
114
115	✓✓	9861	IF( iaire.EQ.1 ) THEN
116			sdd1_type = type_sddv
117			sdd2_type = type_unsddv
118			ELSE
119			sdd1_type = type_unsddv
120			sdd2_type = type_sddv
121			ENDIF
122
123			c IF( iaire.EQ.1 ) THEN
124			c CALL SCOPY( iim, sddv, 1, sdd1, 1 )
125			c CALL SCOPY( iim, unsddv, 1, sdd2, 1 )
126			c ELSE
127			c CALL SCOPY( iim, unsddv, 1, sdd1, 1 )
128			c CALL SCOPY( iim, sddv, 1, sdd2, 1 )
129			c END IF
130
131			jdfil1 = 2
132		9861	jffil1 = jfiltnu
133		9861	jdfil2 = jfiltsu
134			jffil2 = jjm
135			END IF
136			ELSE
137	✗✓	3083	IF( nlat.NE.jjm ) THEN
138			PRINT 2222
139			STOP
140			ELSE
141
142	✓✓	3083	IF( iaire.EQ.1 ) THEN
143			sdd1_type = type_sddu
144			sdd2_type = type_unsddu
145			ELSE
146			sdd1_type = type_unsddu
147			sdd2_type = type_sddu
148			ENDIF
149
150			c IF( iaire.EQ.1 ) THEN
151			c CALL SCOPY( iim, sddu, 1, sdd1, 1 )
152			c CALL SCOPY( iim, unsddu, 1, sdd2, 1 )
153			c ELSE
154			c CALL SCOPY( iim, unsddu, 1, sdd1, 1 )
155			c CALL SCOPY( iim, sddu, 1, sdd2, 1 )
156			c END IF
157
158			jdfil1 = 1
159		3083	jffil1 = jfiltnv
160		3083	jdfil2 = jfiltsv
161			jffil2 = jjm
162			END IF
163			END IF
164
165	✓✓	38832	DO hemisph = 1, 2
166
167	✓✓	25888	IF ( hemisph.EQ.1 ) THEN
168			jdfil = jdfil1
169			jffil = jffil1
170			ELSE
171			jdfil = jdfil2
172			jffil = jffil2
173			END IF
174
175	✓✓	997140	DO l = 1, nbniv
176	✓✓	5853400	DO j = jdfil,jffil
177	✓✓	161227832	DO i = 1, iim
178		160256580	champ(i,j,l) = champ(i,j,l) * sdd12(i,sdd1_type) ! sdd1(i)
179			END DO
180			END DO
181			END DO
182
183	✓✓	25888	IF( hemisph. EQ. 1 ) THEN
184
185	✓✓	12944	IF( ifiltre. EQ. -2 ) THEN
186
187	✓✓	2022	DO j = jdfil,jffil
188			#ifdef BLAS
189			CALL SGEMM("N", "N", iim, nbniv, iim, 1.0,
190			& matrinvn(1,1,j),
191			& iim, champ(1,j,1), iip1*nlat, 0.0,
192			& eignq(1,j-jdfil+1,1), iim*nlat)
193			#else
194			eignq(:,j-jdfil+1,:)
195	✗✓✗✗ ✗✗✗✗ ✗✗✗✗ ✗✗✗✗ ✓✓✓✓	2170617	$ = matmul(matrinvn(:,:,j), champ(:iim,j,:))
196			#endif
197			END DO
198
199	✓✓	12607	ELSE IF ( griscal ) THEN
200
201	✓✓	57144	DO j = jdfil,jffil
202			#ifdef BLAS
203			CALL SGEMM("N", "N", iim, nbniv, iim, 1.0,
204			& matriceun(1,1,j),
205			& iim, champ(1,j,1), iip1*nlat, 0.0,
206			& eignq(1,j-jdfil+1,1), iim*nlat)
207			#else
208			eignq(:,j-jdfil+1,:)
209	✓✓✓✓ ✓✓✓✓ ✓✓✓✓ ✓✗✓✗ ✓✓✓✓	62598504	$ = matmul(matriceun(:,:,j), champ(:iim,j,:))
210			#endif
211			END DO
212
213			ELSE
214
215	✓✓	18498	DO j = jdfil,jffil
216			#ifdef BLAS
217			CALL SGEMM("N", "N", iim, nbniv, iim, 1.0,
218			& matricevn(1,1,j),
219			& iim, champ(1,j,1), iip1*nlat, 0.0,
220			& eignq(1,j-jdfil+1,1), iim*nlat)
221			#else
222			eignq(:,j-jdfil+1,:)
223	✓✓✓✓ ✓✓✓✓ ✓✓✓✓ ✓✗✓✗ ✓✓✓✓	20692873	$ = matmul(matricevn(:,:,j), champ(:iim,j,:))
224			#endif
225			END DO
226
227			ENDIF
228
229			ELSE
230
231	✓✓	12944	IF( ifiltre. EQ. -2 ) THEN
232
233	✓✓	2022	DO j = jdfil,jffil
234			#ifdef BLAS
235			CALL SGEMM("N", "N", iim, nbniv, iim, 1.0,
236			& matrinvs(1,1,j-jfiltsu+1),
237			& iim, champ(1,j,1), iip1*nlat, 0.0,
238			& eignq(1,j-jdfil+1,1), iim*nlat)
239			#else
240			eignq(:,j-jdfil+1,:)
241		1685	$ = matmul(matrinvs(:,:,j-jfiltsu+1),
242	✗✓✗✗ ✗✗✗✗ ✗✗✗✗ ✗✗✗✗ ✓✓✓✓	2170617	$ champ(:iim,j,:))
243			#endif
244			END DO
245
246
247	✓✓	12607	ELSE IF ( griscal ) THEN
248
249	✓✓	57144	DO j = jdfil,jffil
250			#ifdef BLAS
251			CALL SGEMM("N", "N", iim, nbniv, iim, 1.0,
252			& matriceus(1,1,j-jfiltsu+1),
253			& iim, champ(1,j,1), iip1*nlat, 0.0,
254			& eignq(1,j-jdfil+1,1), iim*nlat)
255			#else
256			eignq(:,j-jdfil+1,:)
257		47620	$ = matmul(matriceus(:,:,j-jfiltsu+1),
258	✓✓✓✓ ✓✓✓✓ ✓✓✓✓ ✓✗✓✗ ✓✓✓✓	62598504	$ champ(:iim,j,:))
259			#endif
260			END DO
261
262			ELSE
263
264	✓✓	18498	DO j = jdfil,jffil
265			#ifdef BLAS
266			CALL SGEMM("N", "N", iim, nbniv, iim, 1.0,
267			& matricevs(1,1,j-jfiltsv+1),
268			& iim, champ(1,j,1), iip1*nlat, 0.0,
269			& eignq(1,j-jdfil+1,1), iim*nlat)
270			#else
271			eignq(:,j-jdfil+1,:)
272		15415	$ = matmul(matricevs(:,:,j-jfiltsv+1),
273	✓✓✓✓ ✓✓✓✓ ✓✓✓✓ ✓✗✓✗ ✓✓✓✓	20693883	$ champ(:iim,j,:))
274			#endif
275			END DO
276
277			ENDIF
278
279			ENDIF
280
281	✓✓	25888	IF( ifiltre.EQ. 2 ) THEN
282
283	✓✓	970180	DO l = 1, nbniv
284	✓✓	5695010	DO j = jdfil,jffil
285	✓✓	156864356	DO i = 1, iim
286			champ( i,j,l ) =
287			& (champ(i,j,l) + eignq(i,j-jdfil+1,l))
288		155919390	& * sdd12(i,sdd2_type) ! sdd2(i)
289			END DO
290			END DO
291			END DO
292
293			ELSE
294
295	✓✓	26960	DO l = 1, nbniv
296	✓✓	158390	DO j = jdfil,jffil
297	✓✓	4363476	DO i = 1, iim
298			champ( i,j,l ) =
299			& (champ(i,j,l) - eignq(i,j-jdfil+1,l))
300		4337190	& * sdd12(i,sdd2_type) ! sdd2(i)
301			END DO
302			END DO
303			END DO
304
305			ENDIF
306
307	✓✓	1010084	DO l = 1, nbniv
308	✓✓	5853400	DO j = jdfil,jffil
309		5827512	champ( iip1,j,l ) = champ( 1,j,l )
310			END DO
311			END DO
312
313
314			ENDDO
315
316			1111 FORMAT(//20x,'ERREUR dans le dimensionnement du tableau CHAMP a
317			& filtrer, sur la grille des scalaires'/)
318			2222 FORMAT(//20x,'ERREUR dans le dimensionnement du tableau CHAMP a fi
319			& ltrer, sur la grille de V ou de Z'/)
320		12944	RETURN
321			END