1 | subroutine interp_horiz (varo,varn,imo,jmo,imn,jmn,lm, |
---|
2 | & rlonuo,rlatvo,rlonun,rlatvn) |
---|
3 | |
---|
4 | c=========================================================== |
---|
5 | c Interpolation Horizontales des variables d'une grille LMDZ |
---|
6 | c (des points SCALAIRES au point SCALAIRES) |
---|
7 | c dans une autre grille LMDZ en conservant la quantite |
---|
8 | c totale pour les variables intensives (/m2) : ex : Pression au sol |
---|
9 | c |
---|
10 | c Francois Forget (01/1995) |
---|
11 | c=========================================================== |
---|
12 | |
---|
13 | IMPLICIT NONE |
---|
14 | |
---|
15 | c Declarations: |
---|
16 | c ============== |
---|
17 | c |
---|
18 | c ARGUMENTS |
---|
19 | c """"""""" |
---|
20 | |
---|
21 | INTEGER,INTENT(IN) :: imo, jmo ! dimensions ancienne grille (input) |
---|
22 | INTEGER,INTENT(IN) :: imn,jmn ! dimensions nouvelle grille (input) |
---|
23 | |
---|
24 | REAL,INTENT(IN) :: rlonuo(imo+1) ! Latitude et |
---|
25 | REAL,INTENT(IN) :: rlatvo(jmo) ! longitude des |
---|
26 | REAL,INTENT(IN) :: rlonun(imn+1) ! bord des |
---|
27 | REAL,INTENT(IN) :: rlatvn(jmn) ! cases "scalaires" (input) |
---|
28 | |
---|
29 | INTEGER,INTENT(IN) :: lm ! dimension verticale (input) |
---|
30 | REAL,INTENT(IN) :: varo (imo+1, jmo+1,lm) ! var dans l'ancienne grille (input) |
---|
31 | REAL,INTENT(OUT) :: varn (imn+1,jmn+1,lm) ! var dans la nouvelle grille (output) |
---|
32 | |
---|
33 | c Autres variables |
---|
34 | c """""""""""""""" |
---|
35 | INTEGER imnmx2,jmnmx2 |
---|
36 | cc parameter (imnmx2=190,jmnmx2=100) |
---|
37 | c parameter (imnmx2=360,jmnmx2=190) |
---|
38 | parameter (imnmx2=720,jmnmx2=380) |
---|
39 | REAL airetest(imnmx2+1,jmnmx2+1) |
---|
40 | INTEGER ii,jj,l |
---|
41 | |
---|
42 | REAL,SAVE :: airen ((imnmx2+1)*(jmnmx2+1)) ! aire dans la nouvelle grille |
---|
43 | REAL airentotn ! aire totale pole nord dans la nouvelle grille |
---|
44 | REAL airentots ! aire totale pole sud dans la nouvelle grille |
---|
45 | c Info sur les ktotal intersection entre les cases new/old grille |
---|
46 | |
---|
47 | c kmax: le nombre max d'intersections entre les 2 grilles horizontales |
---|
48 | c On fixe kmax a la taille de la grille des donnees martiennes (360x179) |
---|
49 | c + des pouiemes (cas ou une maille est a cheval sur 2 ou 4 mailles) |
---|
50 | c Il y a un test dans iniinterp_h pour s'assurer que ktotal < kmax |
---|
51 | INTEGER kmax, k |
---|
52 | integer,save :: ktotal |
---|
53 | c parameter (kmax = 360*179 + 200000) |
---|
54 | cc parameter (kmax = 360*179 + 40000) |
---|
55 | parameter (kmax = 360*179 + 360000) |
---|
56 | |
---|
57 | INTEGER,SAVE :: iik(kmax), jjk(kmax),jk(kmax),ik(kmax) |
---|
58 | REAL,SAVE :: intersec(kmax) |
---|
59 | REAL r |
---|
60 | REAL totn, tots |
---|
61 | integer,save :: prev_sumdim=0 |
---|
62 | |
---|
63 | logical,save :: firsttest=.true. , aire_ok=.true. |
---|
64 | |
---|
65 | integer,save :: imoS,jmoS,imnS,jmnS |
---|
66 | |
---|
67 | c Test dimensions imnmx2 jmnmx2 |
---|
68 | c"""""""""""""""""""""""""""""" |
---|
69 | c test dimensionnement tableau airetest |
---|
70 | if (imn.GT.imnmx2.OR.jmn.GT.jmnmx2) then |
---|
71 | write(*,*) 'STOP pb dimensionnement tableau airetest' |
---|
72 | write(*,*) 'il faut imn < imnmx2 et jmn < jmnmx2' |
---|
73 | write(*,*) 'imn imnmx2', imn,imnmx2 |
---|
74 | write(*,*) 'jmn jmnmx2', jmn,jmnmx2 |
---|
75 | call exit(1) |
---|
76 | endif |
---|
77 | |
---|
78 | c initialisation |
---|
79 | c -------------- |
---|
80 | c Si c'est le premier appel, on prepare l'interpolation |
---|
81 | c en calculant pour chaque case autour d'un point scalaire de la |
---|
82 | c nouvelle grille, la surface de intersection avec chaque |
---|
83 | c case de l'ancienne grille. |
---|
84 | |
---|
85 | c This must also be done if we change the dimension |
---|
86 | if (imo+jmo+imn+jmn.ne.prev_sumdim) then |
---|
87 | firsttest=.true. |
---|
88 | prev_sumdim=imo+jmo+imn+jmn |
---|
89 | end if |
---|
90 | |
---|
91 | if (firsttest) then |
---|
92 | call iniinterp_h(imo,jmo,imn,jmn ,kmax, |
---|
93 | & rlonuo,rlatvo,rlonun,rlatvn, |
---|
94 | & ktotal,iik,jjk,jk,ik,intersec,airen) |
---|
95 | imoS=imo |
---|
96 | jmoS=jmo |
---|
97 | imnS=imn |
---|
98 | jmnS=jmn |
---|
99 | else |
---|
100 | if(imo.NE.imoS.OR.jmo.NE.jmoS.OR.imn.NE.imnS.OR.jmn.NE.jmnS) then |
---|
101 | call iniinterp_h(imo,jmo,imn,jmn ,kmax, |
---|
102 | & rlonuo,rlatvo,rlonun,rlatvn, |
---|
103 | & ktotal,iik,jjk,jk,ik,intersec,airen) |
---|
104 | imoS=imo |
---|
105 | jmoS=jmo |
---|
106 | imnS=imn |
---|
107 | jmnS=jmn |
---|
108 | end if |
---|
109 | end if |
---|
110 | |
---|
111 | ! initialize varn() to zero |
---|
112 | varn(1:imn+1,1:jmn+1,1:lm)=0. |
---|
113 | |
---|
114 | c Interpolation horizontale |
---|
115 | c ------------------------- |
---|
116 | c boucle sur toute les ktotal intersections entre les cases |
---|
117 | c de l'ancienne et la nouvelle grille |
---|
118 | c |
---|
119 | ! Ehouarn 2012: for some strange reason, with ifort v12.x, |
---|
120 | ! when the order of the loop below is changed |
---|
121 | ! values of varn(:,:,l=2...) are then sometimes remain zero! |
---|
122 | do l=1,lm |
---|
123 | do k=1,ktotal |
---|
124 | varn(iik(k),jjk(k),l) = varn(iik(k),jjk(k),l) |
---|
125 | & + varo(ik(k), jk(k),l)*intersec(k)/airen(iik(k) |
---|
126 | & +(jjk(k)-1)*(imn+1)) |
---|
127 | end do |
---|
128 | end do |
---|
129 | |
---|
130 | c Une seule valeur au pole pour les variables ! : |
---|
131 | c ----------------------------------------------- |
---|
132 | DO l=1, lm |
---|
133 | totn =0. |
---|
134 | tots =0. |
---|
135 | |
---|
136 | |
---|
137 | c moyenne du champ au poles (ponderee par les aires) |
---|
138 | c""""""""""""""""""""""""""""""" |
---|
139 | airentotn=0. |
---|
140 | airentots=0. |
---|
141 | |
---|
142 | do ii =1, imn+1 |
---|
143 | totn = totn + varn(ii,1,l)*airen(ii) |
---|
144 | tots = tots + varn (ii,jmn+1,l)*airen(jmn*(imn+1)+ii) |
---|
145 | airentotn=airentotn + airen(ii) |
---|
146 | airentots=airentots + airen(jmn*(imn+1)+ii) |
---|
147 | end do |
---|
148 | |
---|
149 | do ii =1, imn+1 |
---|
150 | varn(ii,1,l) = totn/airentotn |
---|
151 | varn(ii,jmn+1,l) = tots/airentots |
---|
152 | end do |
---|
153 | |
---|
154 | ENDDO ! of DO l=1, lm |
---|
155 | |
---|
156 | |
---|
157 | c--------------------------------------------------------------- |
---|
158 | c TEST TEST TEST TEST TEST TEST TEST TEST TEST TEST |
---|
159 | if (firsttest) then |
---|
160 | firsttest = .false. |
---|
161 | write (*,*) 'INTERP. HORIZ. : TEST SUR LES AIRES:' |
---|
162 | |
---|
163 | do jj =1 , jmn+1 |
---|
164 | do ii=1, imn+1 |
---|
165 | airetest(ii,jj) =0. |
---|
166 | end do |
---|
167 | end do |
---|
168 | do k=1,ktotal |
---|
169 | airetest(iik(k),jjk(k))= airetest(iik(k),jjk(k)) +intersec(k) |
---|
170 | end do |
---|
171 | do jj =1 , jmn+1 |
---|
172 | do ii=1, imn+1 |
---|
173 | r = airen(ii+(jj-1)*(imn+1))/airetest(ii,jj) |
---|
174 | if ((r.gt.1.001).or.(r.lt.0.999)) then |
---|
175 | write (*,*) '********** PROBLEME D'' AIRES !!!', |
---|
176 | & ' DANS L''INTERPOLATION HORIZONTALE' |
---|
177 | write(*,*)'ii,jj,airen,airetest', |
---|
178 | & ii,jj,airen(ii+(jj-1)*(imn+1)),airetest(ii,jj) |
---|
179 | aire_ok = .false. |
---|
180 | end if |
---|
181 | end do |
---|
182 | end do |
---|
183 | if (aire_ok) write(*,*) 'INTERP. HORIZ. : AIRES OK' |
---|
184 | endif ! of if (firsttest) |
---|
185 | |
---|
186 | c FIN TEST FIN TEST FIN TEST FIN TEST FIN TEST FIN TEST FIN TEST |
---|
187 | c -------------------------------------------------------------------- |
---|
188 | |
---|
189 | |
---|
190 | end |
---|