trunk/filtrez/filtreg_scal.f

module filtreg_m

  IMPLICIT NONE

contains

  SUBROUTINE filtreg(champ, direct, intensive)

    ! From filtrez/filtreg.F, version 1.1.1.1, 2004/05/19 12:53:09
    ! Author: P. Le Van
    ! Objet : filtre matriciel longitudinal, avec les matrices pr\'ecalcul\'ees
    ! pour l'op\'erateur filtre. 

    USE coefils, ONLY: sddu, sddv, unsddu, unsddv
    USE dimens_m, ONLY: iim, jjm
    use inifilr_m, only: jfiltnu, jfiltnv, jfiltsu, jfiltsv, matriceun, &
         matriceus, matricevn, matricevs, matrinvn, matrinvs
    use nr_util, only: assert

    REAL, intent(inout):: champ(:, :, :) ! (iim + 1, nlat, nbniv)
    ! en entr\'ee : champ \`a filtrer, en sortie : champ filtr\'e

    logical, intent(in):: direct ! filtre direct ou inverse 

    logical, intent(in):: intensive
    ! champ intensif ou extensif (pond\'er\'e par les aires)

    ! Local:
    LOGICAL griscal
    INTEGER nlat ! nombre de latitudes \`a filtrer
    integer nbniv ! nombre de niveaux verticaux \`a filtrer 
    INTEGER jdfil1, jdfil2, jffil1, jffil2, jdfil, jffil
    INTEGER j, l
    REAL eignq(iim), sdd1(iim), sdd2(iim)
    INTEGER hemisph

    !-----------------------------------------------------------

    call assert(size(champ, 1) == iim + 1, "filtreg iim + 1")
    nlat = size(champ, 2)
    nbniv = size(champ, 3)
    call assert(nlat == jjm .or. nlat == jjm + 1, "filtreg nlat")
    griscal = nlat == jjm + 1

    IF (.not. direct .AND. nlat == jjm) THEN
       PRINT *, 'filtreg: inverse filter on scalar grid only'
       STOP 1
    END IF

    IF (griscal) THEN
       IF (intensive) THEN
          sdd1 = sddv
          sdd2 = unsddv
       ELSE
          sdd1 = unsddv
          sdd2 = sddv
       END IF

       jdfil1 = 2
       jffil1 = jfiltnu
       jdfil2 = jfiltsu
       jffil2 = jjm
    ELSE
       IF (intensive) THEN
          sdd1 = sddu
          sdd2 = unsddu
       ELSE
          sdd1 = unsddu
          sdd2 = sddu
       END IF

       jdfil1 = 1
       jffil1 = jfiltnv
       jdfil2 = jfiltsv
       jffil2 = jjm
    END IF

    DO hemisph = 1, 2
       IF (hemisph==1) THEN
          jdfil = jdfil1
          jffil = jffil1
       ELSE
          jdfil = jdfil2
          jffil = jffil2
       END IF

       DO l = 1, nbniv
          DO j = jdfil, jffil
             champ(:iim, j, l) = champ(:iim, j, l) * sdd1

             IF (hemisph==1) THEN
                IF (.not. direct) THEN
                   eignq = matmul(matrinvn(:, :, j), champ(:iim, j, l))
                ELSE IF (griscal) THEN
                   eignq = matmul(matriceun(:, :, j), champ(:iim, j, l))
                ELSE
                   eignq = matmul(matricevn(:, :, j), champ(:iim, j, l))
                END IF
             ELSE
                IF (.not. direct) THEN
                   eignq = matmul(matrinvs(:, :, j - jfiltsu + 1), &
                        champ(:iim, j, l))
                ELSE IF (griscal) THEN
                   eignq = matmul(matriceus(:, :, j - jfiltsu + 1), &
                        champ(:iim, j, l))
                ELSE
                   eignq = matmul(matricevs(:, :, j - jfiltsv + 1), &
                        champ(:iim, j, l))
                END IF
             END IF

             IF (direct) THEN
                champ(:iim, j, l) = (champ(:iim, j, l) + eignq) * sdd2
             ELSE
                champ(:iim, j, l) = (champ(:iim, j, l) - eignq) * sdd2
             END IF

             champ(iim + 1, j, l) = champ(1, j, l)
          END DO
       END DO
    end DO

  END SUBROUTINE filtreg

end module filtreg_m
1	guez	27	module filtreg_m
2	guez	3
3	guez	27	IMPLICIT NONE
4	guez	3
5	guez	27	contains
6	guez	3
7	guez	107	SUBROUTINE filtreg(champ, direct, intensive)
8	guez	3
9	guez	27	! From filtrez/filtreg.F, version 1.1.1.1, 2004/05/19 12:53:09
10			! Author: P. Le Van
11	guez	133	! Objet : filtre matriciel longitudinal, avec les matrices pr\'ecalcul\'ees
12			! pour l'op\'erateur filtre.
13	guez	3
14	guez	71	USE coefils, ONLY: sddu, sddv, unsddu, unsddv
15	guez	64	USE dimens_m, ONLY: iim, jjm
16	guez	54	use inifilr_m, only: jfiltnu, jfiltnv, jfiltsu, jfiltsv, matriceun, &
17			matriceus, matricevn, matricevs, matrinvn, matrinvs
18	guez	107	use nr_util, only: assert
19	guez	3
20	guez	107	REAL, intent(inout):: champ(:, :, :) ! (iim + 1, nlat, nbniv)
21	guez	133	! en entr\'ee : champ \`a filtrer, en sortie : champ filtr\'e
22	guez	3
23	guez	107	logical, intent(in):: direct ! filtre direct ou inverse
24	guez	3
25	guez	107	logical, intent(in):: intensive
26	guez	133	! champ intensif ou extensif (pond\'er\'e par les aires)
27	guez	3
28	guez	107	! Local:
29			LOGICAL griscal
30	guez	133	INTEGER nlat ! nombre de latitudes \`a filtrer
31			integer nbniv ! nombre de niveaux verticaux \`a filtrer
32	guez	27	INTEGER jdfil1, jdfil2, jffil1, jffil2, jdfil, jffil
33	guez	135	INTEGER j, l
34	guez	27	REAL eignq(iim), sdd1(iim), sdd2(iim)
35			INTEGER hemisph
36	guez	25
37	guez	27	!-----------------------------------------------------------
38	guez	25
39	guez	107	call assert(size(champ, 1) == iim + 1, "filtreg iim + 1")
40			nlat = size(champ, 2)
41			nbniv = size(champ, 3)
42			call assert(nlat == jjm .or. nlat == jjm + 1, "filtreg nlat")
43			griscal = nlat == jjm + 1
44	guez	25
45	guez	107	IF (.not. direct .AND. nlat == jjm) THEN
46			PRINT *, 'filtreg: inverse filter on scalar grid only'
47	guez	54	STOP 1
48	guez	27	END IF
49	guez	25
50	guez	27	IF (griscal) THEN
51	guez	107	IF (intensive) THEN
52			sdd1 = sddv
53			sdd2 = unsddv
54	guez	27	ELSE
55	guez	107	sdd1 = unsddv
56			sdd2 = sddv
57			END IF
58	guez	25
59	guez	107	jdfil1 = 2
60			jffil1 = jfiltnu
61			jdfil2 = jfiltsu
62			jffil2 = jjm
63	guez	27	ELSE
64	guez	107	IF (intensive) THEN
65			sdd1 = sddu
66			sdd2 = unsddu
67	guez	27	ELSE
68	guez	107	sdd1 = unsddu
69			sdd2 = sddu
70			END IF
71	guez	25
72	guez	107	jdfil1 = 1
73			jffil1 = jfiltnv
74			jdfil2 = jfiltsv
75			jffil2 = jjm
76	guez	27	END IF
77	guez	25
78	guez	135	DO hemisph = 1, 2
79	guez	27	IF (hemisph==1) THEN
80			jdfil = jdfil1
81			jffil = jffil1
82			ELSE
83			jdfil = jdfil2
84			jffil = jffil2
85			END IF
86	guez	25
87	guez	135	DO l = 1, nbniv
88			DO j = jdfil, jffil
89			champ(:iim, j, l) = champ(:iim, j, l) * sdd1
90	guez	25
91	guez	27	IF (hemisph==1) THEN
92	guez	107	IF (.not. direct) THEN
93	guez	135	eignq = matmul(matrinvn(:, :, j), champ(:iim, j, l))
94	guez	27	ELSE IF (griscal) THEN
95	guez	135	eignq = matmul(matriceun(:, :, j), champ(:iim, j, l))
96	guez	27	ELSE
97	guez	135	eignq = matmul(matricevn(:, :, j), champ(:iim, j, l))
98	guez	27	END IF
99			ELSE
100	guez	107	IF (.not. direct) THEN
101	guez	135	eignq = matmul(matrinvs(:, :, j - jfiltsu + 1), &
102			champ(:iim, j, l))
103	guez	27	ELSE IF (griscal) THEN
104	guez	135	eignq = matmul(matriceus(:, :, j - jfiltsu + 1), &
105			champ(:iim, j, l))
106	guez	27	ELSE
107	guez	135	eignq = matmul(matricevs(:, :, j - jfiltsv + 1), &
108			champ(:iim, j, l))
109	guez	27	END IF
110			END IF
111	guez	25
112	guez	107	IF (direct) THEN
113	guez	135	champ(:iim, j, l) = (champ(:iim, j, l) + eignq) * sdd2
114	guez	27	ELSE
115	guez	135	champ(:iim, j, l) = (champ(:iim, j, l) - eignq) * sdd2
116	guez	27	END IF
117	guez	25
118	guez	27	champ(iim + 1, j, l) = champ(1, j, l)
119	guez	135	END DO
120			END DO
121			end DO
122	guez	27
123			END SUBROUTINE filtreg
124
125			end module filtreg_m