/[lmdze]/trunk/filtrez/filtreg_scal.f
ViewVC logotype

Contents of /trunk/filtrez/filtreg_scal.f

Parent Directory Parent Directory | Revision Log Revision Log


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

  ViewVC Help
Powered by ViewVC 1.1.21