1 | MODULE wrf_lmdz_mod |
---|
2 | ! |
---|
3 | ! Module for the interface to include LMDZ (http://lmdz.lmd.jussieu.fr/) physical packages in WRF |
---|
4 | ! |
---|
5 | ! L. Fita, Laboratoire Meterologie Dynamique, IPSL, UPMC, CNRS, Jussieu, Paris, France |
---|
6 | ! November 2013 |
---|
7 | ! |
---|
8 | |
---|
9 | ! WRF modules |
---|
10 | ! USE module_domain_type |
---|
11 | |
---|
12 | !!!!!!! Subroutines |
---|
13 | ! wrf_varKsoil: Subroutine to provide WRF variable values from the 4 LMDZ soil types |
---|
14 | ! lmdz_varKsoil: Subroutine to provide variable values to the 4 LMDZ soil types according to different methods |
---|
15 | ! LMDZmoist_WRFmoist: Subroutine to transform from LMDZ moisture array to the WRF moisture array |
---|
16 | ! WRFmoist_LMDZmoist: Subroutine to transform from WRF moisture array to the LMDZ moisture array |
---|
17 | ! eta_to_pressure: Subroutine to transform from eta levels to pressure levels |
---|
18 | ! temp_to_theta: Subroutine to transform from temperatures to WRF potential temperatures |
---|
19 | ! theta_to_temp: Subroutine to transform from WRF potential temperatures to temperature |
---|
20 | ! temp_to_theta1: Subroutine to transform from temperature to potential temperature |
---|
21 | ! theta_to_temp1: Subroutine to transform from potential temperature to temperature |
---|
22 | ! def_get_scalar_value: Subroutine to get a scalar value from a .def file |
---|
23 | ! NOTE: For that scalar values with scientific notation (EN, or e-N) it will not work... |
---|
24 | ! next version! |
---|
25 | ! value_kind: Subroutine to determine which kind of value is given |
---|
26 | ! string_split: Subtoutine to split a string according to a character and returns the 'Nsec' of the split |
---|
27 | ! load_lmdz: Subroutine to load LMDZ variables with values from WRF |
---|
28 | ! load_lmdz_limit: Subroutine to load LMDZ limit variables with values from WRF |
---|
29 | ! get_lmdz: Subroutine to get LMDZ values and fill WRF variables |
---|
30 | ! get_lowbdy: Subroutine to load LMDZ limit variables with values from lowbdy WRF |
---|
31 | ! vect_mat: Subroutine to transform a LMDZ physic 1D vector to a 3D matrix |
---|
32 | ! vect_mat_zstagg: Subroutine to transform a LMDZ physic 1D vector to a 3D matrix |
---|
33 | ! mat_vect: Subroutine to transform a 3D matrix to a LMDZ physic 1D vector |
---|
34 | ! mat_vect_zstagg: Subroutine to transform a 3D matrix z-staggered to a LMDZ physic 1D vector |
---|
35 | ! WRFlanduse_LMDZKsoil: Subroutine to retrieve LMDZ Ksoil classification using WRF landuse |
---|
36 | ! interpolate_layers: Subroutine to interpolate values between different layers |
---|
37 | ! (constant values along the depth of a layer) |
---|
38 | ! interpolate_lin: Program to interpolate values y=a+b*x with: (from clWRF modifications) |
---|
39 | ! table_characteristics: Subroutine to read values from a WRF .TBL |
---|
40 | ! read_table: Subroutine to read values from a WRF .TBL |
---|
41 | ! compute_soil_mass: Subroutine to compute the total soil mass (kg) from a point with |
---|
42 | ! different percentages of soil types |
---|
43 | ! compute_TOTsoil_mass_water: Subroutine to compute total soil mass and water from a |
---|
44 | ! given grid point using a given data-set from SOILPARM.TBL |
---|
45 | ! compute_TOTsoil_mass_water_values: Subroutine to compute total soil mass and water |
---|
46 | ! from a given grid point using the values from a given data-set of SOILPARM.TBL |
---|
47 | ! get_lmdz_out: Subroutine to get all the LMDZ output variables into the WRF |
---|
48 | ! Registry structure |
---|
49 | ! put_lmdz_in_WRFout: Subroutine to get LMDZ output and get into WRF standard output variables |
---|
50 | |
---|
51 | ! LMDZ modules |
---|
52 | ! USE indice_sol_mod |
---|
53 | |
---|
54 | ! IMPLICIT NONE |
---|
55 | |
---|
56 | ! INCLUDE "dimsoil.h" |
---|
57 | |
---|
58 | CONTAINS |
---|
59 | |
---|
60 | SUBROUTINE wrf_varKsoil(is,ie,js,je,dx,dy,wbdy,kmethod,ter,lic,oce,sic,vter,vlic, & |
---|
61 | voce,vsic,var) |
---|
62 | ! Subroutine to provide WRF variable values from the 4 LMDZ soil types |
---|
63 | |
---|
64 | IMPLICIT NONE |
---|
65 | |
---|
66 | INTEGER, INTENT(IN) :: is,ie,js,je,dx,dy,wbdy |
---|
67 | CHARACTER(LEN=50), INTENT(IN) :: kmethod |
---|
68 | REAL, DIMENSION(is:ie,js:je), INTENT(IN) :: ter,lic,oce,sic |
---|
69 | REAL, DIMENSION(is:ie,js:je), INTENT(IN) :: vter,vlic,voce,vsic |
---|
70 | REAL, DIMENSION(is:ie,js:je), INTENT(OUT) :: var |
---|
71 | |
---|
72 | ! Local |
---|
73 | INTEGER :: i,j,ddx,ddy |
---|
74 | CHARACTER(LEN=50) :: errmsg, fname |
---|
75 | REAL*8, DIMENSION(is:ie,js:je) :: prop |
---|
76 | |
---|
77 | !!!!!!! Variables |
---|
78 | ! d[x/y]: dimension of the matrices |
---|
79 | ! kmethod: method to use |
---|
80 | ! 'prod': product, variable is distributed by percentage of the soil type vsoil=var*soil |
---|
81 | ! ter: percentage of continent |
---|
82 | ! lic: percentage of frozen continent |
---|
83 | ! oce: percentage of ocean |
---|
84 | ! sic: percentage of frozen ocean |
---|
85 | ! vter: variable on continent |
---|
86 | ! vlic: variable on frozen continent |
---|
87 | ! voce: variable on ocean |
---|
88 | ! vsic: variable on frozen ocean |
---|
89 | ! prop: proportion between soil fractions |
---|
90 | ! var: variable to be reinterpreted from LMDZ soil types |
---|
91 | |
---|
92 | fname="'wrf_varksoil'" |
---|
93 | errmsg='ERROR -- error -- ERROR -- error' |
---|
94 | |
---|
95 | prop=-1.d0 |
---|
96 | var=0. |
---|
97 | ddx=dx-2*wbdy |
---|
98 | ddy=dy-2*wbdy |
---|
99 | |
---|
100 | methodkind: IF (TRIM(kmethod) == 'prod') THEN |
---|
101 | DO i=1,ddx |
---|
102 | DO j=1,ddy |
---|
103 | var(i,j)=vter(i,j)*ter(i,j)+vlic(i,j)*lic(i,j)+voce(i,j)*oce(i,j)+vsic(i,j)*sic(i,j) |
---|
104 | END DO |
---|
105 | END DO |
---|
106 | |
---|
107 | ELSE |
---|
108 | PRINT *,TRIM(errmsg) |
---|
109 | PRINT *,' ' // TRIM(fname) // ": method '" // TRIM(kmethod) // "' does not exist !!!!!" |
---|
110 | STOP |
---|
111 | |
---|
112 | END IF methodkind |
---|
113 | |
---|
114 | END SUBROUTINE wrf_varKsoil |
---|
115 | |
---|
116 | SUBROUTINE lmdz_varKsoil(is,ie,js,je,dx,dy,wbdy,var,kmethod,ter,lic,oce,sic,vter, & |
---|
117 | vlic,voce,vsic) |
---|
118 | ! Subroutine to provide variable values to the 4 LMDZ soil types according to |
---|
119 | ! different methods |
---|
120 | |
---|
121 | IMPLICIT NONE |
---|
122 | |
---|
123 | INTEGER, INTENT(IN) :: is,ie,js,je,dx,dy,wbdy |
---|
124 | REAL, DIMENSION(is:ie,js:je), INTENT(IN) :: var |
---|
125 | CHARACTER(LEN=50), INTENT(IN) :: kmethod |
---|
126 | REAL, DIMENSION(is:ie,js:je), INTENT(IN) :: ter,lic,oce,sic |
---|
127 | REAL, DIMENSION(is:ie,js:je), INTENT(OUT) :: vter,vlic,voce,vsic |
---|
128 | |
---|
129 | ! Local |
---|
130 | INTEGER :: i,j,ddx,ddy |
---|
131 | CHARACTER(LEN=50) :: errmsg, fname |
---|
132 | REAL*8, DIMENSION(is:ie,js:je) :: prop |
---|
133 | |
---|
134 | !!!!!!! Variables |
---|
135 | ! d[x/y]: dimension of the matrices |
---|
136 | ! var: variable to be reinterpreted in LMDZ soil types |
---|
137 | ! kmethod: method to use |
---|
138 | ! 'NOlic': NO lic, variable is not located in ice continent soil type: lic |
---|
139 | ! 'NOlnd': NO land, variable is distributed only to sea soil types: oce, sic |
---|
140 | ! 'NOneg': NO negative, variable is distributed only to free ice soil types: ter, oce |
---|
141 | ! 'NOoce': NO oce, variable is not located in ocean soil type: oce |
---|
142 | ! 'NOpos': NO positive, variable is distributed only to iced soil types: lic, sic |
---|
143 | ! 'NOsea': NO sea, variable is distributed only to land soil types: ter, lic |
---|
144 | ! 'NOsic': NO sic, variable is not located in ice ocean soil type: sic |
---|
145 | ! 'NOter': NO ter, variable is not located in land soil type: ter |
---|
146 | ! 'prod': product, variable is distributed by percentage of the soil type vsoil=var*soil |
---|
147 | ! 'lic:' lic, variable only is located in ice land soil type: lic |
---|
148 | ! 'oce:' oce, variable only is located in ocean soil type: oce |
---|
149 | ! 'sic:' sic, variable only is located in ice ocean soil type: sic |
---|
150 | ! 'ter:' ter, variable only is located in land soil type: ter |
---|
151 | ! ter: percentage of continent |
---|
152 | ! lic: percentage of frozen continent |
---|
153 | ! oce: percentage of ocean |
---|
154 | ! sic: percentage of frozen ocean |
---|
155 | ! vter: variable on continent |
---|
156 | ! vlic: variable on frozen continent |
---|
157 | ! voce: variable on ocean |
---|
158 | ! vsic: variable on frozen ocean |
---|
159 | ! prop: proportion between soil fractions |
---|
160 | |
---|
161 | fname="'lmdz_varksoil'" |
---|
162 | errmsg='ERROR -- error -- ERROR -- error' |
---|
163 | |
---|
164 | ddx=dx-2*wbdy |
---|
165 | ddy=dy-2*wbdy |
---|
166 | |
---|
167 | prop=-1.d0 |
---|
168 | vter=0. |
---|
169 | vlic=0. |
---|
170 | voce=0. |
---|
171 | vsic=0. |
---|
172 | |
---|
173 | methodkind: IF (TRIM(kmethod) == 'NOlic') THEN |
---|
174 | WHERE (ter+oce+sic /= 0.) prop=1.d0/(ter+oce+sic) |
---|
175 | DO i=1,ddx |
---|
176 | DO j=1,ddy |
---|
177 | IF (prop(i,j) >= 0.) THEN |
---|
178 | vter(i,j)=var(i,j)*ter(i,j)*prop(i,j) |
---|
179 | voce(i,j)=var(i,j)*oce(i,j)*prop(i,j) |
---|
180 | vsic(i,j)=var(i,j)*sic(i,j)*prop(i,j) |
---|
181 | END IF |
---|
182 | END DO |
---|
183 | END DO |
---|
184 | |
---|
185 | ELSEIF (TRIM(kmethod) == 'NOlnd') THEN |
---|
186 | WHERE (oce+sic /= 0.) prop=oce*1.d0/(oce+sic) |
---|
187 | DO i=1,ddx |
---|
188 | DO j=1,ddy |
---|
189 | IF (prop(i,j) >= 0.) THEN |
---|
190 | voce(i,j)=var(i,j)*prop(i,j) |
---|
191 | vsic(i,j)=var(i,j)*(1.-prop(i,j)) |
---|
192 | END IF |
---|
193 | END DO |
---|
194 | END DO |
---|
195 | |
---|
196 | ELSEIF (TRIM(kmethod) == 'NOneg') THEN |
---|
197 | WHERE (ter+oce /= 0.) prop=ter*1.d0/(ter+oce) |
---|
198 | DO i=1,ddx |
---|
199 | DO j=1,ddy |
---|
200 | IF (prop(i,j) >= 0.) THEN |
---|
201 | vter(i,j)=var(i,j)*prop(i,j) |
---|
202 | voce(i,j)=var(i,j)*(1.-prop(i,j)) |
---|
203 | END IF |
---|
204 | END DO |
---|
205 | END DO |
---|
206 | |
---|
207 | ELSEIF (TRIM(kmethod) == 'NOoce') THEN |
---|
208 | WHERE (ter+lic+sic /= 0.) prop=1.d0/(ter+lic+sic) |
---|
209 | DO i=1,ddx |
---|
210 | DO j=1,ddy |
---|
211 | IF (prop(i,j) >= 0.) THEN |
---|
212 | vter(i,j)=var(i,j)*ter(i,j)*prop(i,j) |
---|
213 | vlic(i,j)=var(i,j)*lic(i,j)*prop(i,j) |
---|
214 | vsic(i,j)=var(i,j)*sic(i,j)*prop(i,j) |
---|
215 | END IF |
---|
216 | END DO |
---|
217 | END DO |
---|
218 | |
---|
219 | ELSEIF (TRIM(kmethod) == 'NOpos') THEN |
---|
220 | WHERE (lic+sic /=0.) prop=lic*1.d0/(lic+sic) |
---|
221 | DO i=1,ddx |
---|
222 | DO j=1,ddy |
---|
223 | IF (prop(i,j) >= 0.) THEN |
---|
224 | vlic(i,j)=var(i,j)*prop(i,j) |
---|
225 | vsic(i,j)=var(i,j)*(1.-prop(i,j)) |
---|
226 | END IF |
---|
227 | END DO |
---|
228 | END DO |
---|
229 | |
---|
230 | ELSEIF (TRIM(kmethod) == 'NOsea') THEN |
---|
231 | WHERE (ter+lic /= 0.) prop=ter*1.d0/(ter+lic) |
---|
232 | DO i=1,ddx |
---|
233 | DO j=1,ddy |
---|
234 | IF (prop(i,j) >= 0.) THEN |
---|
235 | vter(i,j)=var(i,j)*prop(i,j) |
---|
236 | vlic(i,j)=var(i,j)*(1.-prop(i,j)) |
---|
237 | END IF |
---|
238 | END DO |
---|
239 | END DO |
---|
240 | |
---|
241 | ELSEIF (TRIM(kmethod) == 'NOsic') THEN |
---|
242 | WHERE (ter+lic+oce /= 0.) prop=1.d0/(ter+lic+oce) |
---|
243 | DO i=1,ddx |
---|
244 | DO j=1,ddy |
---|
245 | IF (prop(i,j) >= 0.) THEN |
---|
246 | vter(i,j)=var(i,j)*ter(i,j)*prop(i,j) |
---|
247 | vlic(i,j)=var(i,j)*lic(i,j)*prop(i,j) |
---|
248 | voce(i,j)=var(i,j)*oce(i,j)*prop(i,j) |
---|
249 | END IF |
---|
250 | END DO |
---|
251 | END DO |
---|
252 | |
---|
253 | ELSEIF (TRIM(kmethod) == 'NOter') THEN |
---|
254 | WHERE (lic+sic+oce /= 0.) prop=1.d0/(lic+sic+oce) |
---|
255 | DO i=1,ddx |
---|
256 | DO j=1,ddy |
---|
257 | IF (prop(i,j) >= 0.) THEN |
---|
258 | vlic(i,j)=var(i,j)*lic(i,j)*prop(i,j) |
---|
259 | vsic(i,j)=var(i,j)*sic(i,j)*prop(i,j) |
---|
260 | voce(i,j)=var(i,j)*oce(i,j)*prop(i,j) |
---|
261 | END IF |
---|
262 | END DO |
---|
263 | END DO |
---|
264 | |
---|
265 | ELSEIF (TRIM(kmethod) == 'lic') THEN |
---|
266 | DO i=1,ddx |
---|
267 | DO j=1,ddy |
---|
268 | IF (lic(i,j) /= 0.) THEN |
---|
269 | vlic(i,j)=var(i,j) |
---|
270 | END IF |
---|
271 | END DO |
---|
272 | END DO |
---|
273 | |
---|
274 | ELSEIF (TRIM(kmethod) == 'oce') THEN |
---|
275 | DO i=1,ddx |
---|
276 | DO j=1,ddy |
---|
277 | IF (oce(i,j) /= 0.) THEN |
---|
278 | voce(i,j)=var(i,j) |
---|
279 | END IF |
---|
280 | END DO |
---|
281 | END DO |
---|
282 | |
---|
283 | ELSEIF (TRIM(kmethod) == 'prod') THEN |
---|
284 | DO i=1,ddx |
---|
285 | DO j=1,ddy |
---|
286 | vter(i,j)=var(i,j)*ter(i,j)*1.d0 |
---|
287 | vlic(i,j)=var(i,j)*lic(i,j)*1.d0 |
---|
288 | voce(i,j)=var(i,j)*oce(i,j)*1.d0 |
---|
289 | vsic(i,j)=var(i,j)*sic(i,j)*1.d0 |
---|
290 | END DO |
---|
291 | END DO |
---|
292 | |
---|
293 | ELSEIF (TRIM(kmethod) == 'sic') THEN |
---|
294 | DO i=1,ddx |
---|
295 | DO j=1,ddy |
---|
296 | IF (sic(i,j) /= 0.) THEN |
---|
297 | vsic(i,j)=var(i,j) |
---|
298 | END IF |
---|
299 | END DO |
---|
300 | END DO |
---|
301 | |
---|
302 | ELSEIF (TRIM(kmethod) == 'ter') THEN |
---|
303 | DO i=1,ddx |
---|
304 | DO j=1,ddy |
---|
305 | IF (ter(i,j) /= 0.) THEN |
---|
306 | vter(i,j)=var(i,j) |
---|
307 | END IF |
---|
308 | END DO |
---|
309 | END DO |
---|
310 | |
---|
311 | ELSE |
---|
312 | PRINT *,TRIM(errmsg) |
---|
313 | PRINT *,' ' // TRIM(fname) // ": method '" // TRIM(kmethod) // "' does not exist !!!!!" |
---|
314 | STOP |
---|
315 | END IF methodkind |
---|
316 | |
---|
317 | END SUBROUTINE lmdz_varKsoil |
---|
318 | |
---|
319 | |
---|
320 | SUBROUTINE LMDZmoist_WRFmoist(lmdzMixRat, dx, dy, dz, Nwmoist, wbdy, wqvid, wqcid, & |
---|
321 | & wqrid, Nlmdzq, wMOIST) |
---|
322 | ! Subroutine to transform from LMDZ moisture array to the WRF moisture array |
---|
323 | |
---|
324 | IMPLICIT NONE |
---|
325 | |
---|
326 | INTEGER, INTENT(IN) :: dx, dy, dz, wbdy |
---|
327 | INTEGER, INTENT(IN) :: Nwmoist, wqvid, wqcid, & |
---|
328 | wqrid |
---|
329 | INTEGER, INTENT(IN) :: Nlmdzq |
---|
330 | REAL, DIMENSION((dx-2*wbdy)*(dy-2*wbdy),dz,Nlmdzq), INTENT(IN) :: & |
---|
331 | lmdzMixRat |
---|
332 | REAL, DIMENSION(1-wbdy:dx-wbdy+1,dz+1,1-wbdy:dy-wbdy+1,Nwmoist), INTENT(OUT) :: wMOIST |
---|
333 | |
---|
334 | ! Local |
---|
335 | INTEGER :: i,j,k, iq |
---|
336 | INTEGER :: ddx, ddy |
---|
337 | INTEGER, ALLOCATABLE, DIMENSION(:) :: widlmdz |
---|
338 | REAL, DIMENSION(1-wbdy:dx-wbdy+1,1-wbdy:dy-wbdy+1) :: qcr_ratio, wcrMOIST |
---|
339 | |
---|
340 | !!!!!!! Variables |
---|
341 | ! widlmdz: vector with the id of the WRF tracers to be also included in the LMDZ moist array |
---|
342 | ! d[x,y,z]: Spatial dimension of the WRF moisture array |
---|
343 | ! Nwmoist: number of species in wMOIST |
---|
344 | ! wqvid: index for the water vapor mixing ratio |
---|
345 | ! wqcid: index for the cloud water mixing ratio |
---|
346 | ! wqrid: index for the rain water mixing ratio |
---|
347 | ! Nlmdzq: number of species in LMDZ |
---|
348 | ! lmdxMixRat: LMDZ mixing ratio array |
---|
349 | ! qcr_ratio: ratio between cloud/rain liquid water in WRF moisture array |
---|
350 | ! wcrMOIST: WRF array combination of cloud and rain liquid water |
---|
351 | ! wMOIST: WRF moisture array |
---|
352 | |
---|
353 | !!!!!! L. Fita, July 203 |
---|
354 | ! At this stage, we only care about water vapor, cloud liquid water and rain liquid water. |
---|
355 | |
---|
356 | ddx=dx-2*wbdy |
---|
357 | ddy=dy-2*wbdy |
---|
358 | |
---|
359 | wMOIST=0. |
---|
360 | wcrMOIST=0. |
---|
361 | |
---|
362 | DO k=1,dz |
---|
363 | DO j=1,ddy |
---|
364 | DO i=1,ddx |
---|
365 | wMOIST(i,k,j,wqvid) = lmdzMixRat(ddx*(j-1)+i,k,1) |
---|
366 | wMOIST(i,k,j,wqcid) = lmdzMixRat(ddx*(j-1)+i,k,2) |
---|
367 | !! wcrMOIST(i,j) = lmdzMixRat(ddx*(j-1)+1,k,2) |
---|
368 | END DO |
---|
369 | END DO |
---|
370 | !! L. Fita, July 2013 |
---|
371 | !! We can not know how to split liquid water between cloud and rain. Thus, we use the previous |
---|
372 | !! time-step proportion |
---|
373 | !! qcr_ratio=wMOIST(:,k,:,wqcid)/(wMOIST(:,k,:,wqcid) + wMOIST(:,k,:,wqrid)) |
---|
374 | !! wMOIST(:,k,:,wqcid) = wcrMOIST*qcr_ratio |
---|
375 | !! wMOIST(:,k,:,wqrid) = wcrMOIST*(1.-qcr_ratio) |
---|
376 | ENDDO |
---|
377 | |
---|
378 | END SUBROUTINE LMDZmoist_WRFmoist |
---|
379 | |
---|
380 | |
---|
381 | SUBROUTINE WRFmoist_LMDZmoist(wMOIST, dx, dy, dz, Nwmoist, wbdy, wqvid, wqcid, & |
---|
382 | & wqrid, Nlmdzq, lmdzMixRat) |
---|
383 | ! Subroutine to transform from WRF moisture array to the LMDZ moisture array |
---|
384 | |
---|
385 | IMPLICIT NONE |
---|
386 | |
---|
387 | INTEGER, INTENT(IN) :: dx, dy, dz, wbdy |
---|
388 | INTEGER, INTENT(IN) :: Nwmoist, wqvid, wqcid, wqrid |
---|
389 | INTEGER, INTENT(IN) :: Nlmdzq |
---|
390 | REAL, DIMENSION(1-wbdy:dx-wbdy+1,dz+1,1-wbdy:dy-wbdy+1,Nwmoist), INTENT(IN) :: wMOIST |
---|
391 | REAL, DIMENSION((dx-2*wbdy)*(dy-2*wbdy),dz,Nlmdzq), INTENT(OUT) :: lmdzMixRat |
---|
392 | |
---|
393 | ! Local |
---|
394 | INTEGER :: i, j, k, iq, dx2, dy2 |
---|
395 | INTEGER :: ddx, ddy |
---|
396 | ! INTEGER, ALLOCATABLE, DIMENSION(:) :: widlmdz |
---|
397 | ! REAL, DIMENSION(dx,dy) :: wcrMOIST |
---|
398 | REAL, DIMENSION(1-wbdy:dx-wbdy+1,1-wbdy:dy-wbdy+1) :: qcr_ratio, wcrMOIST |
---|
399 | |
---|
400 | |
---|
401 | !!!!!!! Variables |
---|
402 | ! wMOIST: WRF moisture array |
---|
403 | ! d[x,y,z]: Spatial dimension of the WRF moisture array |
---|
404 | ! Nwmoist: number of species in wMOIST |
---|
405 | ! wqvid: index for the water vapor mixing ratio |
---|
406 | ! wqcid: index for the cloud water mixing ratio |
---|
407 | ! wqrid: index for the rain water mixing ratio |
---|
408 | ! wcrMOIST: WRF array combination of cloud and rain liquid water |
---|
409 | ! Nlmdzq: number of species in LMDZ |
---|
410 | ! lmdxMixRat: LMDZ mixing ratio array |
---|
411 | ! widlmdz: vector with the id of the WRF tracers to be also included in the LMDZ moist array |
---|
412 | |
---|
413 | ddx=dx-2*wbdy |
---|
414 | ddy=dy-2*wbdy |
---|
415 | |
---|
416 | lmdzMixRat=0. |
---|
417 | |
---|
418 | DO k=1,dz |
---|
419 | ! Only taking water vapor and cloud water mixing ratios |
---|
420 | DO j=1,ddy |
---|
421 | DO i=1,ddx |
---|
422 | lmdzMixRat(ddx*(j-1)+i,k,1) = wMOIST(i,k,j,wqvid) |
---|
423 | lmdzMixRat(ddx*(j-1)+i,k,2) = wMOIST(i,k,j,wqcid) |
---|
424 | END DO |
---|
425 | END DO |
---|
426 | ENDDO |
---|
427 | |
---|
428 | !! IF (Nlmdzq > 2) THEN |
---|
429 | ! Whenever we would like to use more WRF water spicies some sort of mapping should be needed. |
---|
430 | ! Maybe something like this: |
---|
431 | ! IF (ALLOCATED(widlmdz)) DEALLOCATE(widlmdz) |
---|
432 | ! ALLOCATE(widlmdz(Nlmdzq-2)) |
---|
433 | |
---|
434 | !! DO iq=3,Nlmdzq |
---|
435 | !! DO k=1,dz |
---|
436 | !! lmdzMixRat(:,k,iq) = (iq-3.+(dz-k-1.)*1./(dz+1.))*1. |
---|
437 | !! END DO |
---|
438 | !! END DO |
---|
439 | !! END IF |
---|
440 | |
---|
441 | ! DEALLOCATE(widlmdz) |
---|
442 | |
---|
443 | END SUBROUTINE WRFmoist_LMDZmoist |
---|
444 | |
---|
445 | |
---|
446 | SUBROUTINE eta_to_pressure(etaP, psfc, ptop, dz, presP) |
---|
447 | ! Subroutine to transform from eta levels to pressure levels |
---|
448 | |
---|
449 | IMPLICIT NONE |
---|
450 | |
---|
451 | INTEGER, INTENT(IN) :: dz |
---|
452 | REAL, INTENT(IN) :: psfc, ptop |
---|
453 | REAL, DIMENSION(dz), INTENT(IN) :: etaP |
---|
454 | REAL, DIMENSION(dz), INTENT(OUT) :: presP |
---|
455 | |
---|
456 | !!!!!!Variables |
---|
457 | ! etaP: vector with the eta values |
---|
458 | ! psfc: pressure at the surface (in Pa) |
---|
459 | ! ptop: pressure at the top (in Pa) |
---|
460 | ! dz: number of vertical levels |
---|
461 | ! presP: equivalent vector of pressures (in Pa) |
---|
462 | |
---|
463 | presP = etaP*(psfc - ptop) + ptop |
---|
464 | |
---|
465 | END SUBROUTINE eta_to_pressure |
---|
466 | |
---|
467 | |
---|
468 | SUBROUTINE temp_to_theta(tempValues, pressValues, thetaValues, & |
---|
469 | & ids,ide, jds,jde, kds,kde, & |
---|
470 | & ims,ime, jms,jme, kms,kme, & |
---|
471 | & ips,ipe, jps,jpe, kps,kpe, & ! patch dims |
---|
472 | & i_start,i_end,j_start,j_end,kts,kte,num_tiles & |
---|
473 | ) |
---|
474 | ! Subroutine to transform from temperatures to WRF potential temperatures |
---|
475 | |
---|
476 | USE module_model_constants |
---|
477 | |
---|
478 | IMPLICIT NONE |
---|
479 | |
---|
480 | !-- thetaValues potential temperature values from WRF (300. to be added) |
---|
481 | !-- pressValues pressure |
---|
482 | !-- tempValues temperature values in K |
---|
483 | ! |
---|
484 | !-- ids start index for i in domain |
---|
485 | !-- ide end index for i in domain |
---|
486 | !-- jds start index for j in domain |
---|
487 | !-- jde end index for j in domain |
---|
488 | !-- kds start index for k in domain |
---|
489 | !-- kde end index for k in domain |
---|
490 | !-- ims start index for i in memory |
---|
491 | !-- ime end index for i in memory |
---|
492 | !-- jms start index for j in memory |
---|
493 | !-- jme end index for j in memory |
---|
494 | !-- ips start index for i in patch |
---|
495 | !-- ipe end index for i in patch |
---|
496 | !-- jps start index for j in patch |
---|
497 | !-- jpe end index for j in patch |
---|
498 | !-- kms start index for k in memory |
---|
499 | !-- kme end index for k in memory |
---|
500 | !-- i_start start indices for i in tile |
---|
501 | !-- i_end end indices for i in tile |
---|
502 | !-- j_start start indices for j in tile |
---|
503 | !-- j_end end indices for j in tile |
---|
504 | !-- kts start index for k in tile |
---|
505 | !-- kte end index for k in tile |
---|
506 | !-- num_tiles number of tiles |
---|
507 | ! |
---|
508 | !====================================================================== |
---|
509 | |
---|
510 | INTEGER, INTENT(IN) :: & |
---|
511 | ids,ide, jds,jde, kds,kde, & |
---|
512 | ims,ime, jms,jme, kms,kme, & |
---|
513 | ips,ipe, jps,jpe, kps,kpe, & |
---|
514 | kts,kte, & |
---|
515 | num_tiles |
---|
516 | |
---|
517 | INTEGER, DIMENSION(num_tiles), INTENT(IN) :: & |
---|
518 | & i_start,i_end,j_start,j_end |
---|
519 | |
---|
520 | |
---|
521 | REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(IN) :: tempValues, & |
---|
522 | & pressValues |
---|
523 | REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(OUT) :: thetaValues |
---|
524 | |
---|
525 | ! Local |
---|
526 | REAL, PARAMETER :: pres0 = 100000. |
---|
527 | |
---|
528 | thetaValues = tempValues * ( pres0 / pressValues ) ** (r_d / cp) |
---|
529 | |
---|
530 | thetaValues = thetaValues - 300. |
---|
531 | |
---|
532 | END SUBROUTINE temp_to_theta |
---|
533 | |
---|
534 | SUBROUTINE theta_to_temp(thetaValues, pressValues, tempValues, & |
---|
535 | & ids,ide, jds,jde, kds,kde, & |
---|
536 | & ims,ime, jms,jme, kms,kme, & |
---|
537 | & ips,ipe, jps,jpe, kps,kpe, & ! patch dims |
---|
538 | & i_start,i_end,j_start,j_end,kts,kte,num_tiles & |
---|
539 | ) |
---|
540 | ! Subroutine to transform from WRF potential temperatures to temperature |
---|
541 | |
---|
542 | USE module_model_constants |
---|
543 | |
---|
544 | IMPLICIT NONE |
---|
545 | |
---|
546 | !-- thetaValues potential temperature values from WRF (300. to be added) |
---|
547 | !-- pressValues pressure |
---|
548 | !-- tempValues temperature values in K |
---|
549 | ! |
---|
550 | !-- ids start index for i in domain |
---|
551 | !-- ide end index for i in domain |
---|
552 | !-- jds start index for j in domain |
---|
553 | !-- jde end index for j in domain |
---|
554 | !-- kds start index for k in domain |
---|
555 | !-- kde end index for k in domain |
---|
556 | !-- ims start index for i in memory |
---|
557 | !-- ime end index for i in memory |
---|
558 | !-- jms start index for j in memory |
---|
559 | !-- jme end index for j in memory |
---|
560 | !-- ips start index for i in patch |
---|
561 | !-- ipe end index for i in patch |
---|
562 | !-- jps start index for j in patch |
---|
563 | !-- jpe end index for j in patch |
---|
564 | !-- kms start index for k in memory |
---|
565 | !-- kme end index for k in memory |
---|
566 | !-- i_start start indices for i in tile |
---|
567 | !-- i_end end indices for i in tile |
---|
568 | !-- j_start start indices for j in tile |
---|
569 | !-- j_end end indices for j in tile |
---|
570 | !-- kts start index for k in tile |
---|
571 | !-- kte end index for k in tile |
---|
572 | !-- num_tiles number of tiles |
---|
573 | ! |
---|
574 | !====================================================================== |
---|
575 | |
---|
576 | INTEGER, INTENT(IN) :: & |
---|
577 | ids,ide, jds,jde, kds,kde, & |
---|
578 | ims,ime, jms,jme, kms,kme, & |
---|
579 | ips,ipe, jps,jpe, kps,kpe, & |
---|
580 | kts,kte, & |
---|
581 | num_tiles |
---|
582 | |
---|
583 | INTEGER, DIMENSION(num_tiles), INTENT(IN) :: & |
---|
584 | & i_start,i_end,j_start,j_end |
---|
585 | |
---|
586 | |
---|
587 | REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(IN) :: thetaValues, & |
---|
588 | & pressvalues |
---|
589 | REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(OUT) :: tempValues |
---|
590 | |
---|
591 | ! Local |
---|
592 | REAL, PARAMETER :: pres0 = 100000. |
---|
593 | INTEGER :: i,j,k |
---|
594 | |
---|
595 | !! PRINT *,' thetaValues: ', LBOUND(thetaValues),' x ', UBOUND(thetaValues) |
---|
596 | |
---|
597 | !! DO k=kms,kme-2 |
---|
598 | !! DO j=jms,jme |
---|
599 | !! DO i=ims,ime |
---|
600 | !! IF (thetaValues(i,k,j) == 0. .AND. i > 0 .AND. i < 32 & |
---|
601 | !! .AND. j > 0 .AND. j < 32) THEN |
---|
602 | !! PRINT *,' Zth: T ',i,k,j,thetaValues(i,k,j) |
---|
603 | !! END IF |
---|
604 | !! END DO |
---|
605 | !! END DO |
---|
606 | !! END DO |
---|
607 | |
---|
608 | tempValues = ( thetaValues + 300.)*( pressValues / pres0 ) ** (r_d / cp) |
---|
609 | |
---|
610 | !! DO k=kms,kme-2 |
---|
611 | END SUBROUTINE theta_to_temp |
---|
612 | |
---|
613 | SUBROUTINE temp_to_theta1(tempValue, pressValue, thetaValue) |
---|
614 | ! Subroutine to transform from temperature to potential temperature |
---|
615 | |
---|
616 | USE module_model_constants |
---|
617 | |
---|
618 | IMPLICIT NONE |
---|
619 | |
---|
620 | !!!!!!! Variables |
---|
621 | ! thetaValue potential temperature value |
---|
622 | ! pressValue pressure |
---|
623 | ! tempValue temperature value in K |
---|
624 | ! |
---|
625 | REAL, INTENT(IN) :: tempValue, pressValue |
---|
626 | REAL, INTENT(OUT) :: thetaValue |
---|
627 | |
---|
628 | ! Local |
---|
629 | REAL, PARAMETER :: pres0 = 100000. |
---|
630 | |
---|
631 | thetaValue = tempValue * ( pres0 / pressValue ) ** (r_d / cp) |
---|
632 | |
---|
633 | END SUBROUTINE temp_to_theta1 |
---|
634 | |
---|
635 | SUBROUTINE theta_to_temp1(thetaValue, pressValue, tempValue) |
---|
636 | ! Subroutine to transform from potential temperature to temperature |
---|
637 | |
---|
638 | USE module_model_constants |
---|
639 | |
---|
640 | IMPLICIT NONE |
---|
641 | |
---|
642 | !!!!!!! Variables |
---|
643 | ! thetaValue potential temperature |
---|
644 | ! pressValue pressure |
---|
645 | ! tempValue temperature value in K |
---|
646 | |
---|
647 | REAL, INTENT(IN) :: thetaValue, pressvalue |
---|
648 | REAL, INTENT(OUT) :: tempValue |
---|
649 | |
---|
650 | ! Local |
---|
651 | REAL, PARAMETER :: pres0 = 100000. |
---|
652 | |
---|
653 | tempValue = thetaValue*( pressValue / pres0 ) ** (r_d / cp) |
---|
654 | |
---|
655 | END SUBROUTINE theta_to_temp1 |
---|
656 | |
---|
657 | SUBROUTINE def_get_scalar_value(filen, valn, ds, di, dr, db) |
---|
658 | ! Subroutine to get a scalar value from a .def file |
---|
659 | ! NOTE: For that scalar values with scientific notation (EN, or e-N) it will not work... |
---|
660 | ! next version! |
---|
661 | |
---|
662 | IMPLICIT NONE |
---|
663 | |
---|
664 | CHARACTER(LEN=100), INTENT(IN) :: filen, valn |
---|
665 | CHARACTER(LEN=100), INTENT(OUT) :: ds |
---|
666 | INTEGER, INTENT(OUT) :: di |
---|
667 | REAL, INTENT(OUT) :: dr |
---|
668 | LOGICAL, INTENT(OUT) :: db |
---|
669 | |
---|
670 | ! Local |
---|
671 | INTEGER :: funit, ios, Lvaln |
---|
672 | INTEGER :: iline, lineval |
---|
673 | CHARACTER(LEN=1) :: kindvalue |
---|
674 | CHARACTER(LEN=100) :: varname,varname0 |
---|
675 | CHARACTER(LEN=50) :: subname, errmsg, value, & |
---|
676 | value0 |
---|
677 | LOGICAL :: fexists, is_used |
---|
678 | |
---|
679 | !!!!!!! Variables |
---|
680 | ! filen: Name of the def file to get the value |
---|
681 | ! valn: name of the value to get |
---|
682 | ! ds: string value (up to 100 characters) |
---|
683 | ! di: integer value |
---|
684 | ! dr: real value |
---|
685 | ! db: logical value |
---|
686 | ! value: value get from the file |
---|
687 | ! kindvalue: kind of the value ('s': string, 'i': integer, 'r': real, 'b':boolean ) |
---|
688 | |
---|
689 | subname = 'def_get_scalar_value' |
---|
690 | errmsg = 'ERROR -- error -- ERROR -- error' |
---|
691 | |
---|
692 | INQUIRE(FILE=TRIM(filen), EXIST=fexists) |
---|
693 | |
---|
694 | IF (fexists) THEN |
---|
695 | Lvaln = LEN_TRIM(valn) |
---|
696 | |
---|
697 | ds = 'None' |
---|
698 | di = 0 |
---|
699 | dr = 0. |
---|
700 | db = .FALSE. |
---|
701 | |
---|
702 | DO funit=10,100 |
---|
703 | INQUIRE(UNIT=funit, OPENED=is_used) |
---|
704 | IF (.NOT. is_used) EXIT |
---|
705 | END DO |
---|
706 | OPEN(funit, FILE=TRIM(filen), STATUS='old', FORM='formatted', IOSTAT=ios) |
---|
707 | IF ( ios /= 0 ) THEN |
---|
708 | PRINT *,errmsg |
---|
709 | PRINT *, ' ' // TRIM(subname) // ": ERROR opening '" // TRIM(filen) // "'" |
---|
710 | PRINT *,'ios: ',ios |
---|
711 | STOP |
---|
712 | END IF |
---|
713 | ! iunit = get_unused_unit() |
---|
714 | ! IF ( iunit <= 0 ) THEN |
---|
715 | ! IF ( wrf_dm_on_monitor() ) THEN |
---|
716 | ! CALL wrf_error_fatal('Error in module_ra_cam_support: could not find a free Fortran unit.') |
---|
717 | ! END IF |
---|
718 | ! END IF |
---|
719 | lineval = 0 |
---|
720 | DO iline=1, 10000 |
---|
721 | READ(funit,*,END=125)varname0 |
---|
722 | varname = TRIM(varname0) |
---|
723 | IF (varname(1:Lvaln) == TRIM(valn)) THEN |
---|
724 | lineval = iline |
---|
725 | EXIT |
---|
726 | END IF |
---|
727 | END DO |
---|
728 | |
---|
729 | 125 CONTINUE |
---|
730 | IF (lineval == 0) THEN |
---|
731 | PRINT *,errmsg |
---|
732 | PRINT *,' ' // TRIM(subname) // ": In file '" // TRIM(filen) // & |
---|
733 | & "' does not exist a variable called '" // TRIM(valn) // '" !!!' |
---|
734 | STOP |
---|
735 | ELSE |
---|
736 | CLOSE(funit) |
---|
737 | OPEN(funit, FILE=TRIM(filen), STATUS='old', FORM='formatted', IOSTAT=ios) |
---|
738 | DO iline=1,lineval-1 |
---|
739 | READ(funit,*)varname |
---|
740 | END DO |
---|
741 | READ(funit,*)value0 |
---|
742 | CLOSE(UNIT=funit) |
---|
743 | |
---|
744 | CALL string_split(value0, '=', 2, value) |
---|
745 | CALL string_split(value, ':', 2, value) |
---|
746 | CALL value_kind(value, kindvalue) |
---|
747 | |
---|
748 | ! PRINT *,'value: ',TRIM(value), ' kind: ',kindvalue |
---|
749 | |
---|
750 | IF (kindvalue == 'i') THEN |
---|
751 | READ(value, '(I50)')di |
---|
752 | ELSE IF (kindvalue == 'r') THEN |
---|
753 | READ(value, '(F50.25)')dr |
---|
754 | ELSE IF (kindvalue == 'e') THEN |
---|
755 | READ (value, '(E50.25)')dr |
---|
756 | ELSE IF (kindvalue == 's') THEN |
---|
757 | ds = TRIM(value) |
---|
758 | ELSE |
---|
759 | IF (TRIM(value) == 'y' .OR. TRIM(value) == 'yes' .OR. TRIM(value) == 'Y' & |
---|
760 | .OR. TRIM(value) == 'Yes' .OR. TRIM(value) == 'YES') db = .TRUE. |
---|
761 | |
---|
762 | IF (TRIM(value) == 'n' .OR. TRIM(value) == 'no' .OR. TRIM(value) == 'N' & |
---|
763 | .OR. TRIM(value) == 'No' .OR. TRIM(value) == 'NO') db = .FALSE. |
---|
764 | END IF |
---|
765 | END IF |
---|
766 | ELSE |
---|
767 | PRINT *,errmsg |
---|
768 | PRINT *,' ' // TRIM(subname) // ": File '" // TRIM(filen) // "' does not exist !!" |
---|
769 | STOP |
---|
770 | END IF |
---|
771 | |
---|
772 | END SUBROUTINE def_get_scalar_value |
---|
773 | |
---|
774 | SUBROUTINE value_kind(value, kindval) |
---|
775 | ! Subroutine to determine which kind of value is given |
---|
776 | |
---|
777 | IMPLICIT NONE |
---|
778 | |
---|
779 | CHARACTER(LEN=50), INTENT(IN) :: value |
---|
780 | CHARACTER(LEN=1), INTENT(OUT) :: kindval |
---|
781 | |
---|
782 | ! Local |
---|
783 | INTEGER :: ic |
---|
784 | INTEGER :: Lval, Nc, NnumVal, NstrVal |
---|
785 | LOGICAL :: is_point, is_exp |
---|
786 | |
---|
787 | !!!!!!! Variables |
---|
788 | ! value: value to determine the cahracteristics |
---|
789 | ! kindval: kind of value. 's': string, 'i': integer, 'r': real, 'e': real scientific notation, |
---|
790 | ! 'b':boolean |
---|
791 | ! Lval: length of characters of the value |
---|
792 | ! NnumVal: number of numeric characers |
---|
793 | ! NstrVal: number of string characters |
---|
794 | ! is_point: boolean variable indicating the presence of a '.' |
---|
795 | ! is_exp: boolean variable indicating the presence of a 'E' or 'e' (exponential) |
---|
796 | |
---|
797 | |
---|
798 | IF (TRIM(value) == 'y' .OR. TRIM(value) == 'n' .OR. TRIM(value) == 'yes' .OR. & |
---|
799 | TRIM(value) == 'no' .OR. TRIM(value) == 'Y' .OR. TRIM(value) == 'N' .OR. & |
---|
800 | TRIM(value) == 'Yes' .OR. TRIM(value) == 'No' .OR. TRIM(value) == 'YES' .OR. & |
---|
801 | TRIM(value) == 'NO') THEN |
---|
802 | ! PRINT *,' boolean!' |
---|
803 | kindval = 'b' |
---|
804 | ELSE |
---|
805 | Lval = LEN_TRIM(value) |
---|
806 | NnumVal = 0 |
---|
807 | NstrVal = 0 |
---|
808 | is_point = .FALSE. |
---|
809 | |
---|
810 | DO ic=1,Lval |
---|
811 | Nc = ICHAR(value(ic:ic)) |
---|
812 | ! PRINT *, 'char: ',value(ic:ic), ' Nc:',Nc |
---|
813 | IF (Nc >= 48 .AND. Nc <= 57) THEN |
---|
814 | Nnumval = Nnumval + 1 |
---|
815 | ELSE |
---|
816 | NstrVal = NstrVal + 1 |
---|
817 | IF (value(ic:ic) == '.') is_point = .TRUE. |
---|
818 | IF ((value(ic:ic) == 'e') .OR. (value(ic:ic) == 'E')) is_exp = .TRUE. |
---|
819 | END IF |
---|
820 | END DO |
---|
821 | |
---|
822 | ! PRINT *,'NnumVal:', NnumVal, 'NstrVal: ', NstrVal, 'point: ', is_point |
---|
823 | IF (Nnumval == Lval) THEN |
---|
824 | ! PRINT *,' integer!' |
---|
825 | kindval = 'i' |
---|
826 | ELSE IF ( (Nnumval == Lval -1) .AND. (is_point .EQV. .TRUE.) .AND. & |
---|
827 | (NstrVal == 1)) THEN |
---|
828 | ! PRINT *,' real!' |
---|
829 | kindval = 'r' |
---|
830 | ELSE IF ( ((Nnumval == Lval -2) .OR.(Nnumval == Lval -3)) .AND. & |
---|
831 | (is_point .EQV. .TRUE.) .AND. (is_exp .EQV. .TRUE.) .AND. & |
---|
832 | ((NstrVal == 2) .OR. (NstrVal == 3)) ) THEN |
---|
833 | ! PRINT *,' real exponential!' |
---|
834 | kindval = 'e' |
---|
835 | ELSE |
---|
836 | ! PRINT *,' string!' |
---|
837 | kindval = 's' |
---|
838 | END IF |
---|
839 | END IF |
---|
840 | |
---|
841 | END SUBROUTINE value_kind |
---|
842 | |
---|
843 | SUBROUTINE string_split(string, charsplit, Nsec, sec) |
---|
844 | ! Subtoutine to split a string according to a character and returns the 'Nsec' of the split |
---|
845 | |
---|
846 | IMPLICIT NONE |
---|
847 | |
---|
848 | CHARACTER(LEN=50), INTENT(IN) :: string |
---|
849 | CHARACTER(LEN=1), INTENT(IN) :: charsplit |
---|
850 | INTEGER, INTENT(IN) :: Nsec |
---|
851 | CHARACTER(LEN=50), INTENT(OUT) :: sec |
---|
852 | |
---|
853 | ! Local |
---|
854 | INTEGER :: ic |
---|
855 | INTEGER :: poschar, Lstring |
---|
856 | INTEGER :: isec, inisec |
---|
857 | CHARACTER(LEN=50) :: str, secchar |
---|
858 | |
---|
859 | !!!!!!! Variables |
---|
860 | ! string: string to split |
---|
861 | ! charsplit: character to use to split |
---|
862 | ! Nsec: number of section to return |
---|
863 | ! sec: returned section of the string |
---|
864 | |
---|
865 | ! PRINT *,'string: ',string |
---|
866 | |
---|
867 | str = TRIM(string) |
---|
868 | poschar=0 |
---|
869 | |
---|
870 | Lstring = LEN_TRIM(str) |
---|
871 | isec = 1 |
---|
872 | poschar = INDEX(str, charsplit) |
---|
873 | inisec = 1 |
---|
874 | |
---|
875 | IF (poschar == 0) THEN |
---|
876 | sec = str |
---|
877 | |
---|
878 | ELSE |
---|
879 | DO ic=poschar,50 |
---|
880 | secchar = str(inisec:poschar-1) |
---|
881 | IF (isec == Nsec) sec = secchar |
---|
882 | inisec = poschar + 1 |
---|
883 | poschar = INDEX(str(inisec:Lstring), charsplit) |
---|
884 | IF (poschar == 0) poschar = Lstring + 2 |
---|
885 | isec = isec + 1 |
---|
886 | END DO |
---|
887 | |
---|
888 | END IF |
---|
889 | |
---|
890 | END SUBROUTINE string_split |
---|
891 | |
---|
892 | SUBROUTINE get_lmdz(ddx,ddy,ddz,bdy,dlmdz,Nks,NzO,plon,plat,pzmasq,ppctsrf,pftsol, & |
---|
893 | pftsoil,pqsol,pqsurf,pfalb1,pfalb2,pevap,psnow,pradsol,psolsw,psollw,pfder, & |
---|
894 | prain_fall,psnow_fall,pagesno,pzmea,pzstd,pzgam,pzthe,pzpic,pzval,prugoro, & |
---|
895 | prugos,pzmax0,pf0,pwake_s,pwake_cstar,pwake_pe,pwake_fip,psst, & |
---|
896 | palbedo,pt_ancien,pq_ancien,pu_ancien,pv_ancien,pclwcon,prnebcon,pratqs, & |
---|
897 | pema_work1,pema_work2,pwake_deltat,pwake_deltaq,pfm_therm,pentr_therm, & |
---|
898 | pdetr_therm,pnat,pwrf_grid) |
---|
899 | ! Subroutine to get LMDZ variables with values from WRF |
---|
900 | |
---|
901 | ! WRF modules |
---|
902 | USE module_domain_type |
---|
903 | |
---|
904 | ! LMDZ modules |
---|
905 | USE indice_sol_mod |
---|
906 | |
---|
907 | IMPLICIT NONE |
---|
908 | |
---|
909 | TYPE(domain) , TARGET :: pwrf_grid |
---|
910 | INTEGER, INTENT(IN) :: ddx, ddy, ddz, bdy, & |
---|
911 | dlmdz, Nks, NzO |
---|
912 | REAL, DIMENSION(dlmdz), INTENT(IN) :: plon, plat, pzmasq, & |
---|
913 | pqsol, pradsol, psolsw, psollw, pfder, prain_fall, psnow_fall, pzmea, pzstd, & |
---|
914 | pzgam, pzthe, pzpic, pzval, prugoro, pzmax0, pf0, pwake_s, & |
---|
915 | pwake_cstar, pwake_pe, pwake_fip, psst, palbedo, pnat |
---|
916 | ! What to do with these variables? ,pnat,pbils,prug |
---|
917 | REAL, DIMENSION(dlmdz,Nks), INTENT(IN) :: ppctsrf,pftsol,pqsurf, & |
---|
918 | pfalb1,pfalb2,pevap,psnow,prugos,pagesno |
---|
919 | REAL, DIMENSION(dlmdz,ddz), INTENT(IN) :: pt_ancien,pq_ancien, & |
---|
920 | pu_ancien,pv_ancien,pclwcon,prnebcon,pratqs,pema_work1,pema_work2,pwake_deltat,& |
---|
921 | pwake_deltaq,pentr_therm,pdetr_therm |
---|
922 | REAL, DIMENSION(dlmdz,NzO,Nks), INTENT(IN) :: pftsoil |
---|
923 | REAL, DIMENSION(dlmdz,ddz+1), INTENT(IN) :: pfm_therm |
---|
924 | |
---|
925 | ! Local |
---|
926 | INTEGER :: ix,iy,iz,ixy,dx,dy |
---|
927 | CHARACTER(LEN=50) :: LMDZvarmethod |
---|
928 | |
---|
929 | ! Filling all 2D variables |
---|
930 | !! |
---|
931 | PRINT *,' Lluis before get_lmdz: 5 11 ltksoil: ',pwrf_grid%ltksoil(5,:,11), & |
---|
932 | ' ftsol: ', pftsol(315,:) |
---|
933 | |
---|
934 | DO iy=1,ddy |
---|
935 | DO ix=1,ddx |
---|
936 | ixy=ddx*(iy-1)+ix |
---|
937 | ! It does not change on time: pzmask |
---|
938 | ! It might change due to ice dynamics? |
---|
939 | pwrf_grid%lter(ix,iy) = ppctsrf(ixy,is_ter) |
---|
940 | pwrf_grid%llic(ix,iy) = ppctsrf(ixy,is_lic) |
---|
941 | pwrf_grid%loce(ix,iy) = ppctsrf(ixy,is_oce) |
---|
942 | pwrf_grid%lsic(ix,iy) = ppctsrf(ixy,is_sic) |
---|
943 | DO iz=1,Nks |
---|
944 | pwrf_grid%ltksoil(ix,iz,iy) = pftsol(ixy,iz) |
---|
945 | END DO |
---|
946 | DO iz=1,NzO |
---|
947 | pwrf_grid%lotter(ix,iz,iy) = pftsoil(ixy,iz,1) |
---|
948 | pwrf_grid%lotlic(ix,iz,iy) = pftsoil(ixy,iz,2) |
---|
949 | pwrf_grid%lotoce(ix,iz,iy) = pftsoil(ixy,iz,3) |
---|
950 | pwrf_grid%lotsic(ix,iz,iy) = pftsoil(ixy,iz,4) |
---|
951 | END DO |
---|
952 | DO iz=1,Nks |
---|
953 | pwrf_grid%lqksoil(ix,iz,iy) = pqsurf(ixy,iz) |
---|
954 | END DO |
---|
955 | pwrf_grid%lwsol(ix,iy) = pqsol(ixy) |
---|
956 | DO iz=1,Nks |
---|
957 | pwrf_grid%lalbksoil(ix,iz,iy) = pfalb1(ixy,iz) |
---|
958 | pwrf_grid%llwalbksoil(ix,iz,iy) = pfalb2(ixy,iz) |
---|
959 | pwrf_grid%levapksoil(ix,iz,iy) = pevap(ixy,iz) |
---|
960 | pwrf_grid%lsnowksoil(ix,iz,iy) = psnow(ixy,iz) |
---|
961 | END DO |
---|
962 | pwrf_grid%lrads(ix,iy) = pradsol(ixy) |
---|
963 | pwrf_grid%lsolsw(ix,iy) = psolsw(ixy) |
---|
964 | pwrf_grid%lsollw(ix,iy) = psollw(ixy) |
---|
965 | pwrf_grid%lfder(ix,iy) = pfder(ixy) |
---|
966 | pwrf_grid%lrain(ix,iy) = prain_fall(ixy) |
---|
967 | pwrf_grid%lsnow(ix,iy) = psnow_fall(ixy) |
---|
968 | DO iz=1,Nks |
---|
969 | pwrf_grid%lrugksoil(ix,iz,iy) = prugos(ixy,iz) |
---|
970 | pwrf_grid%lagesnoksoil(ix,iz,iy) = pagesno(ixy,iz) |
---|
971 | END DO |
---|
972 | ! These variables do not change on time! |
---|
973 | !! |
---|
974 | ! pzmea, pzstd, pzsig, pzgam, pzthe, pzpic, pzval, prugsrel |
---|
975 | pwrf_grid%lrugsea(ix,iy) = prugos(ixy,is_oce) |
---|
976 | ! pwrf_grid%lrunofflic(ix,iy)= prun_off_lic(ixy) |
---|
977 | pwrf_grid%lzmax0(ix,iy) = pzmax0(ixy) |
---|
978 | pwrf_grid%lf0(ix,iy) = pf0(ixy) |
---|
979 | pwrf_grid%lwake_s(ix,iy) = pwake_s(ixy) |
---|
980 | pwrf_grid%lwake_cstar(ix,iy) = pwake_cstar(ixy) |
---|
981 | pwrf_grid%lwake_pe(ix,iy) = pwake_pe(ixy) |
---|
982 | pwrf_grid%lwake_fip(ix,iy) = pwake_fip(ixy) |
---|
983 | pwrf_grid%lnat(ix,iy) = pnat(ixy) |
---|
984 | pwrf_grid%sst(ix,iy) = psst(ixy) |
---|
985 | ! pwrf_grid%lbils(ix,iy) = pbils(ixy) |
---|
986 | pwrf_grid%albedo(ix,iy) = palbedo(ixy) |
---|
987 | ! pwrf_grid%lrug(ix,iy) = prug(ixy) |
---|
988 | END DO |
---|
989 | END DO |
---|
990 | |
---|
991 | dx=ddx+2*bdy |
---|
992 | dy=ddy+2*bdy |
---|
993 | |
---|
994 | CALL vect_mat(pclwcon, dx, dy, ddz, bdy, pwrf_grid%lclwcon) |
---|
995 | CALL vect_mat(prnebcon, dx, dy, ddz, bdy, pwrf_grid%lrnebcon) |
---|
996 | CALL vect_mat(pratqs, dx, dy, ddz, bdy, pwrf_grid%lratqs) |
---|
997 | CALL vect_mat(pema_work1, dx, dy, ddz, bdy, pwrf_grid%lema_work1) |
---|
998 | CALL vect_mat(pema_work2, dx, dy, ddz, bdy, pwrf_grid%lema_work2) |
---|
999 | CALL vect_mat(pwake_deltat, dx, dy, ddz, bdy, pwrf_grid%lwake_deltat) |
---|
1000 | CALL vect_mat(pwake_deltaq, dx, dy, ddz, bdy, pwrf_grid%lwake_deltaq) |
---|
1001 | CALL vect_mat(pfm_therm, dx, dy, ddz+1, bdy, pwrf_grid%lfm_therm) |
---|
1002 | CALL vect_mat(pentr_therm, dx, dy, ddz, bdy, pwrf_grid%lentr_therm) |
---|
1003 | CALL vect_mat(pdetr_therm, dx, dy, ddz, bdy, pwrf_grid%ldetr_therm) |
---|
1004 | |
---|
1005 | PRINT *,' Lluis after get_lmdz: 5 11 ltksoil: ',pwrf_grid%ltksoil(5,:,11), & |
---|
1006 | ' ftsol: ', pftsol(315,:) |
---|
1007 | |
---|
1008 | RETURN |
---|
1009 | |
---|
1010 | END SUBROUTINE get_lmdz |
---|
1011 | |
---|
1012 | SUBROUTINE get_lowbdy(ddx, ddy, ddz, bdy, dlmdz, Nks, Nz0, pkglo, pwrf_grid, & |
---|
1013 | ppctsrf, pftsol, prugos, palbedo, psst, ptsurf_lim, pz0_lim, palb_lim, prug_glo, & |
---|
1014 | ppct_glo) |
---|
1015 | ! Subroutine to load LMDZ limit variables with values from lowbdy WRF |
---|
1016 | |
---|
1017 | ! WRF modules |
---|
1018 | ! USE module_model_constants |
---|
1019 | USE module_domain_type |
---|
1020 | |
---|
1021 | ! LMDZ modules |
---|
1022 | USE indice_sol_mod |
---|
1023 | |
---|
1024 | IMPLICIT NONE |
---|
1025 | |
---|
1026 | TYPE(domain) , TARGET :: pwrf_grid |
---|
1027 | INTEGER, INTENT(IN) :: ddx, ddy, ddz, bdy, & |
---|
1028 | dlmdz, Nks, Nz0, pkglo |
---|
1029 | REAL, DIMENSION(dlmdz), INTENT(INOUT) :: psst, palbedo, & |
---|
1030 | ptsurf_lim, pz0_lim, palb_lim |
---|
1031 | ! What to do with these variables? ,pnat,pbils |
---|
1032 | REAL, DIMENSION(dlmdz,Nks), INTENT(INOUT) :: ppctsrf, pftsol, prugos |
---|
1033 | REAL, DIMENSION(pkglo), INTENT(OUT) :: prug_glo |
---|
1034 | REAL, DIMENSION(pkglo,Nks), INTENT(OUT) :: ppct_glo |
---|
1035 | |
---|
1036 | ! Local |
---|
1037 | INTEGER :: ix,iy,iz,ixy,dx,dy,ixx,iyy |
---|
1038 | CHARACTER(LEN=50) :: LMDZvarmethod,fname,errmsg |
---|
1039 | |
---|
1040 | fname='get_lowbdy' |
---|
1041 | errmsg='ERROR -- error -- ERROR -- error' |
---|
1042 | |
---|
1043 | ! Filling all 2D variables |
---|
1044 | !! |
---|
1045 | DO iy=1,ddy |
---|
1046 | DO ix=1,ddx |
---|
1047 | ixy=ddx*(iy-1)+ix |
---|
1048 | ppctsrf(ixy,is_ter) = pwrf_grid%lter(ix,iy) |
---|
1049 | ppctsrf(ixy,is_lic) = pwrf_grid%llic(ix,iy) |
---|
1050 | IF ( pwrf_grid%xice(ix,iy) /= pwrf_grid%lsic(ix,iy) ) THEN |
---|
1051 | ! No ocean/frozen ocean fraction before |
---|
1052 | IF (pwrf_grid%loce(ix,iy) == 0.) pwrf_grid%lrugksoil(ix,3,iy) = 0.0001 |
---|
1053 | IF (pwrf_grid%lsic(ix,iy) == 0.) pwrf_grid%lrugksoil(ix,4,iy) = 0.001 |
---|
1054 | pwrf_grid%loce(ix,iy) = 1. - pwrf_grid%lsic(ix,iy) |
---|
1055 | pwrf_grid%lsic(ix,iy) = pwrf_grid%lsic(ix,iy) |
---|
1056 | ppctsrf(ixy,is_oce) = pwrf_grid%loce(ix,iy) |
---|
1057 | ppctsrf(ixy,is_sic) = pwrf_grid%lsic(ix,iy) |
---|
1058 | ! Recomputing total roughness |
---|
1059 | !! |
---|
1060 | pwrf_grid%lrug(ix,iy)=pwrf_grid%lter(ix,iy)*pwrf_grid%lrugksoil(ix,1,iy) +& |
---|
1061 | pwrf_grid%llic(ix,iy)*pwrf_grid%lrugksoil(ix,2,iy) + & |
---|
1062 | pwrf_grid%loce(ix,iy)*pwrf_grid%lrugksoil(ix,3,iy) + & |
---|
1063 | pwrf_grid%lsic(ix,iy)*pwrf_grid%lrugksoil(ix,4,iy) |
---|
1064 | END IF |
---|
1065 | ! Full SST grid points must have a SST value from the wrfinput (WRF 1/0 landmask) |
---|
1066 | ! but that LMDZ points with fractions of land/mask could not. Using the closest one |
---|
1067 | ! if it is an isolated grid point, using TSK instead |
---|
1068 | IF (pwrf_grid%sst(ix,iy) < 200.15) THEN |
---|
1069 | pwrf_grid%ltksoil(ix,3,iy) = -9. |
---|
1070 | IF ( (ix > 1).AND.(ix < ddx) .AND. (iy > 1).AND.(iy < ddy) ) THEN |
---|
1071 | DO ixx=-1,1 |
---|
1072 | DO iyy=-1,1 |
---|
1073 | ! PRINT *,ixx,iyy,wrf_grid%sst(ix+ixx,iy+iyy) |
---|
1074 | IF (pwrf_grid%sst(ix+ixx,iy+iyy) > 200.15) THEN |
---|
1075 | pwrf_grid%ltksoil(ix,3,iy) = pwrf_grid%sst(ix+ixx,iy+iyy) |
---|
1076 | EXIT |
---|
1077 | END IF |
---|
1078 | END DO |
---|
1079 | IF ( pwrf_grid%ltksoil(ix,3,iy) > 200.15) EXIT |
---|
1080 | END DO |
---|
1081 | IF ( pwrf_grid%ltksoil(ix,3,iy) == -9.) pwrf_grid%ltksoil(ix,3,iy) = & |
---|
1082 | pwrf_grid%tsk(ix,iy) |
---|
1083 | ELSE |
---|
1084 | pwrf_grid%ltksoil(ix,3,iy) = pwrf_grid%tsk(ix,iy) |
---|
1085 | END IF |
---|
1086 | ELSE |
---|
1087 | pwrf_grid%ltksoil(ix,3,iy) = pwrf_grid%sst(ix,iy) |
---|
1088 | END IF |
---|
1089 | IF ( (pwrf_grid%ltksoil(ix,3,iy) < 200.15) .OR. & |
---|
1090 | (pwrf_grid%ltksoil(ix,3,iy) > 370.) ) THEN |
---|
1091 | PRINT *,TRIM(errmsg) |
---|
1092 | ! WRITE(wrf_err_message,*) ' ' // TRIM(fname) // & |
---|
1093 | ! ' SST inconsistency at point: ', ix, ', ', iy,' with loce: ', & |
---|
1094 | ! wrf_grid%loce(ix,iy),' has a tsoil(oce)= ', & |
---|
1095 | ! wrf_grid%ltksoil(ix,3,iy),' sst: ',wrf_grid%sst(ix,iy),' tsk: ', & |
---|
1096 | ! wrf_grid%tsk(ix,iy) |
---|
1097 | ! CALL wrf_error_fatal(TRIM(wrf_err_message)) |
---|
1098 | PRINT *, ' ' // TRIM(fname) // & |
---|
1099 | ' SST inconsistency at point: ', ix, ', ', iy,' with loce: ', & |
---|
1100 | pwrf_grid%loce(ix,iy),' has a tsoil(oce)= ', & |
---|
1101 | pwrf_grid%ltksoil(ix,3,iy),' sst: ', pwrf_grid%sst(ix,iy),' tsk: ', & |
---|
1102 | pwrf_grid%tsk(ix,iy) |
---|
1103 | END IF |
---|
1104 | psst(ixy) = pwrf_grid%ltksoil(ix,3,iy) |
---|
1105 | palbedo(ixy) = pwrf_grid%albbck(ix,iy) |
---|
1106 | DO iz=1,Nks |
---|
1107 | prugos(ixy,iz) = pwrf_grid%lrugksoil(ix,iz,iy) |
---|
1108 | END DO |
---|
1109 | |
---|
1110 | ! In ocean_forced_mod is used as tsurf_lim = sst(knindex(1:knon)). Let's complicate it |
---|
1111 | ! a little bit... |
---|
1112 | ptsurf_lim(ixy)= psst(ixy) |
---|
1113 | ! In surf_land_bucket_mod is used as z0_new=rugos(knindex(1:knon),1) [z0_new == z0_limit] |
---|
1114 | pz0_lim(ixy) = prugos(ixy,1) |
---|
1115 | palb_lim(ixy) = palbedo(ixy) |
---|
1116 | |
---|
1117 | prug_glo(ixy) = prugos(ixy,is_ter) |
---|
1118 | ppct_glo(ixy,:) = ppctsrf(ixy,:) |
---|
1119 | |
---|
1120 | END DO |
---|
1121 | END DO |
---|
1122 | |
---|
1123 | RETURN |
---|
1124 | |
---|
1125 | END SUBROUTINE get_lowbdy |
---|
1126 | |
---|
1127 | SUBROUTINE load_lmdz(pwrf_grid,ddx,ddy,ddz,bdy,dlmdz,Nks,NzO,prlon,prlat,pzmasq, & |
---|
1128 | ppctsrf,pftsol,pftsoil,pqsurf,pqsol,pfalb1,pfalb2,pevap,psnow,pradsol,psolsw, & |
---|
1129 | psollw,pfder,prain_fall,psnow_fall,pagesno,pzmea,pzstd,pzgam,pzthe,pzpic,pzval, & |
---|
1130 | prugoro,prugos,prun_off_lic,pzmax0,pf0,pwake_s,pwake_cstar,pwake_pe,pwake_fip, & |
---|
1131 | psst,palbedo,pt_ancien,pq_ancien,pu_ancien,pv_ancien,pclwcon,prnebcon,pratqs, & |
---|
1132 | pema_work1,pema_work2,pwake_deltat,pwake_deltaq,pfm_therm,pentr_therm, & |
---|
1133 | pdetr_therm,pnat) |
---|
1134 | ! Subroutine to load LMDZ variables with values from WRF |
---|
1135 | |
---|
1136 | ! WRF modules |
---|
1137 | ! USE module_model_constants |
---|
1138 | USE module_domain_type |
---|
1139 | |
---|
1140 | ! LMDZ modules |
---|
1141 | USE indice_sol_mod |
---|
1142 | |
---|
1143 | IMPLICIT NONE |
---|
1144 | |
---|
1145 | TYPE(domain) , TARGET :: pwrf_grid |
---|
1146 | INTEGER, INTENT(IN) :: ddx, ddy, ddz, bdy, & |
---|
1147 | dlmdz, Nks, NzO |
---|
1148 | REAL, DIMENSION(dlmdz), INTENT(OUT) :: prlon, prlat, pzmasq, & |
---|
1149 | pqsol, pradsol, psolsw, psollw, pfder, prain_fall, psnow_fall, pzmea, pzstd, & |
---|
1150 | pzgam, pzthe, pzpic, pzval, prugoro, prun_off_lic, pzmax0, pf0, pwake_s, & |
---|
1151 | pwake_cstar, pwake_pe, pwake_fip, psst, palbedo,pnat |
---|
1152 | REAL, DIMENSION(dlmdz,Nks), INTENT(OUT) :: ppctsrf,pftsol,pqsurf, & |
---|
1153 | pfalb1,pfalb2,pevap,psnow,prugos,pagesno |
---|
1154 | REAL, DIMENSION(dlmdz,ddz), INTENT(OUT) :: pt_ancien,pq_ancien, & |
---|
1155 | pu_ancien,pv_ancien,pclwcon,prnebcon,pratqs,pema_work1,pema_work2,pwake_deltat,& |
---|
1156 | pwake_deltaq,pentr_therm,pdetr_therm |
---|
1157 | REAL, DIMENSION(dlmdz,NzO,Nks), INTENT(OUT) :: pftsoil |
---|
1158 | REAL, DIMENSION(dlmdz,ddz+1), INTENT(OUT) :: pfm_therm |
---|
1159 | |
---|
1160 | |
---|
1161 | ! Local |
---|
1162 | INTEGER :: ix,iy,iz,ixy,dx,dy |
---|
1163 | CHARACTER(LEN=50) :: LMDZvarmethod,fname |
---|
1164 | |
---|
1165 | fname='load_lmdz' |
---|
1166 | |
---|
1167 | ! Filling all 2D variables |
---|
1168 | !! |
---|
1169 | PRINT *,' ' // TRIM(fname) // '; LUBOUND(grid%lmsk): ',LBOUND(pwrf_grid%lmsk), UBOUND(pwrf_grid%lmsk) |
---|
1170 | PRINT *,' ' // TRIM(fname) // '; ddx ddy: ',ddy,ddx |
---|
1171 | |
---|
1172 | DO iy=1,ddy |
---|
1173 | DO ix=1,ddx |
---|
1174 | ixy=ddx*(iy-1)+ix |
---|
1175 | prlon(ixy) = pwrf_grid%xlong(ix,iy) |
---|
1176 | prlat(ixy) = pwrf_grid%xlat(ix,iy) |
---|
1177 | pzmasq(ixy) = pwrf_grid%lmsk(ix,iy) |
---|
1178 | ppctsrf(ixy,is_ter) = pwrf_grid%lter(ix,iy) |
---|
1179 | ppctsrf(ixy,is_lic) = pwrf_grid%llic(ix,iy) |
---|
1180 | ppctsrf(ixy,is_oce) = pwrf_grid%loce(ix,iy) |
---|
1181 | ppctsrf(ixy,is_sic) = pwrf_grid%lsic(ix,iy) |
---|
1182 | DO iz=1,Nks |
---|
1183 | pftsol(ixy,iz) = pwrf_grid%ltksoil(ix,iz,iy) |
---|
1184 | END DO |
---|
1185 | DO iz=1,NzO |
---|
1186 | pftsoil(ixy,iz,1) = pwrf_grid%lotter(ix,iz,iy) |
---|
1187 | pftsoil(ixy,iz,2) = pwrf_grid%lotlic(ix,iz,iy) |
---|
1188 | pftsoil(ixy,iz,3) = pwrf_grid%lotoce(ix,iz,iy) |
---|
1189 | pftsoil(ixy,iz,4) = pwrf_grid%lotsic(ix,iz,iy) |
---|
1190 | END DO |
---|
1191 | DO iz=1,Nks |
---|
1192 | pqsurf(ixy,iz) = pwrf_grid%lqksoil(ix,iz,iy) |
---|
1193 | END DO |
---|
1194 | pqsol(ixy) = pwrf_grid%lwsol(ix,iy) |
---|
1195 | DO iz=1,Nks |
---|
1196 | pfalb1(ixy,iz) = pwrf_grid%lalbksoil(ix,iz,iy) |
---|
1197 | pfalb2(ixy,iz) = pwrf_grid%llwalbksoil(ix,iz,iy) |
---|
1198 | pevap(ixy,iz) = pwrf_grid%levapksoil(ix,iz,iy) |
---|
1199 | psnow(ixy,iz) = pwrf_grid%lsnowksoil(ix,iz,iy) |
---|
1200 | END DO |
---|
1201 | ! No way to get back the independent radiative values (unless using LMDZ 'hist*' |
---|
1202 | ! files). Use it in the LMDZ's added way also from grid (no need to initialize it, |
---|
1203 | ! no radiation on t=0) |
---|
1204 | pradsol(ixy) = pwrf_grid%lrads(ix,iy) |
---|
1205 | psolsw(ixy) = pwrf_grid%lsolsw(ix,iy) |
---|
1206 | psollw(ixy) = pwrf_grid%lsollw(ix,iy) |
---|
1207 | pfder(ixy) = pwrf_grid%lfder(ix,iy) |
---|
1208 | ! LMDZ does not provide a way (in the restart) to distingusih between convective, |
---|
1209 | ! non-convective. Thus Use it as a whole |
---|
1210 | prain_fall(ixy) = pwrf_grid%lrain(ix,iy) |
---|
1211 | psnow_fall(ixy) = pwrf_grid%lsnow(ix,iy) |
---|
1212 | DO iz=1,Nks |
---|
1213 | prugos(ixy,iz) = pwrf_grid%lrugksoil(ix,iz,iy) |
---|
1214 | pagesno(ixy,iz) = pwrf_grid%lagesnoksoil(ix,iz,iy) |
---|
1215 | END DO |
---|
1216 | pzmea(ixy) = pwrf_grid%lzmea(ix,iy) |
---|
1217 | pzstd(ixy) = pwrf_grid%lzstd(ix,iy) |
---|
1218 | pzgam(ixy) = pwrf_grid%lzgam(ix,iy) |
---|
1219 | pzthe(ixy) = pwrf_grid%lzthe(ix,iy) |
---|
1220 | pzpic(ixy) = pwrf_grid%lzpic(ix,iy) |
---|
1221 | pzval(ixy) = pwrf_grid%lzval(ix,iy) |
---|
1222 | prugoro(ixy) = pwrf_grid%lzrugsrel(ix,iy) |
---|
1223 | ! prugos(ixy,is_oce) = pwrf_grid%lrugsea(ix,iy) |
---|
1224 | prun_off_lic(ixy) = pwrf_grid%lrunofflic(ix,iy) |
---|
1225 | pzmax0(ixy) = pwrf_grid%lzmax0(ix,iy) |
---|
1226 | pf0(ixy) = pwrf_grid%lf0(ix,iy) |
---|
1227 | pwake_s(ixy) = pwrf_grid%lwake_s(ix,iy) |
---|
1228 | pwake_cstar(ixy) = pwrf_grid%lwake_cstar(ix,iy) |
---|
1229 | pwake_pe(ixy) = pwrf_grid%lwake_pe(ix,iy) |
---|
1230 | pwake_fip(ixy) = pwrf_grid%lwake_fip(ix,iy) |
---|
1231 | !! limit.nc |
---|
1232 | ! nat == pctsrf, not needed ? |
---|
1233 | pnat(ixy) = pwrf_grid%lnat(ix,iy) |
---|
1234 | psst(ixy) = pwrf_grid%sst(ix,iy) |
---|
1235 | ! not used, not defined.... |
---|
1236 | ! pbils(ixy) = pwrf_grid%lbils(ix,iy) |
---|
1237 | palbedo(ixy) = pwrf_grid%albedo(ix,iy) |
---|
1238 | ! already done in previous steps |
---|
1239 | ! prugos(ixy) = pwrf_grid%lrug(ix,iy) |
---|
1240 | |
---|
1241 | END DO |
---|
1242 | END DO |
---|
1243 | |
---|
1244 | dx=ddx+2*bdy |
---|
1245 | dy=ddy+2*bdy |
---|
1246 | |
---|
1247 | CALL mat_vect(pwrf_grid%t_2, dx, dy, ddz, bdy, pt_ancien) |
---|
1248 | ! New variable added in the Registry qv_2 |
---|
1249 | CALL mat_vect(pwrf_grid%qv_2, dx, dy, ddz, bdy, pq_ancien) |
---|
1250 | CALL mat_vect(pwrf_grid%u_2, dx, dy, ddz, bdy, pu_ancien) |
---|
1251 | CALL mat_vect(pwrf_grid%v_2, dx, dy, ddz, bdy, pv_ancien) |
---|
1252 | CALL mat_vect(pwrf_grid%lclwcon, dx, dy, ddz, bdy, pclwcon) |
---|
1253 | CALL mat_vect(pwrf_grid%lrnebcon, dx, dy, ddz, bdy, prnebcon) |
---|
1254 | CALL mat_vect(pwrf_grid%lratqs, dx, dy, ddz, bdy, pratqs) |
---|
1255 | CALL mat_vect(pwrf_grid%lema_work1, dx, dy, ddz, bdy, pema_work1) |
---|
1256 | CALL mat_vect(pwrf_grid%lema_work2, dx, dy, ddz, bdy, pema_work2) |
---|
1257 | CALL mat_vect(pwrf_grid%lwake_deltat, dx, dy, ddz, bdy, pwake_deltat) |
---|
1258 | CALL mat_vect(pwrf_grid%lwake_deltaq, dx, dy, ddz, bdy, pwake_deltaq) |
---|
1259 | CALL mat_vect(pwrf_grid%lfm_therm, dx, dy, ddz+1, bdy, pfm_therm) |
---|
1260 | CALL mat_vect(pwrf_grid%lentr_therm, dx, dy, ddz, bdy, pentr_therm) |
---|
1261 | CALL mat_vect(pwrf_grid%ldetr_therm, dx, dy, ddz, bdy, pdetr_therm) |
---|
1262 | |
---|
1263 | RETURN |
---|
1264 | |
---|
1265 | END SUBROUTINE load_lmdz |
---|
1266 | |
---|
1267 | SUBROUTINE vect_mat(vect, dx, dy, dz, wbdy, mat) |
---|
1268 | ! Subroutine to transform a LMDZ physic 1D vector to a 3D matrix |
---|
1269 | |
---|
1270 | IMPLICIT NONE |
---|
1271 | |
---|
1272 | INTEGER, INTENT(IN) :: dx, dy, dz, wbdy |
---|
1273 | REAL, DIMENSION((dx-2*wbdy)*(dy-2*wbdy),dz), INTENT(IN) :: vect |
---|
1274 | REAL, DIMENSION(1-wbdy:dx-wbdy+1,dz+1,1-wbdy:dy-wbdy+1), INTENT(OUT) :: mat |
---|
1275 | |
---|
1276 | ! Local |
---|
1277 | INTEGER :: i,j,k |
---|
1278 | INTEGER :: ddx, ddy |
---|
1279 | INTEGER :: lp, il, jl |
---|
1280 | |
---|
1281 | ! Variables |
---|
1282 | ! d[x/y/z]: dimensions |
---|
1283 | ! mat: matrix to transform |
---|
1284 | ! vect: vector to produce |
---|
1285 | mat=0. |
---|
1286 | |
---|
1287 | ddx=dx-2*wbdy |
---|
1288 | ddy=dy-2*wbdy |
---|
1289 | ! lp=321 |
---|
1290 | ! jl=INT(lp/ddx)+1 |
---|
1291 | ! il=lp-ddx*(jl-1) |
---|
1292 | |
---|
1293 | DO k=1,dz |
---|
1294 | DO j=1,ddy |
---|
1295 | DO i=1,ddx |
---|
1296 | mat(i,k,j)=vect(ddx*(j-1)+i,k) |
---|
1297 | END DO |
---|
1298 | END DO |
---|
1299 | ! PRINT *,k,' Lluis vect: ',vect(lp,k),mat(il,k,jl),' lp: ',lp,' il jl :',il,jl |
---|
1300 | END DO |
---|
1301 | |
---|
1302 | END SUBROUTINE vect_mat |
---|
1303 | |
---|
1304 | |
---|
1305 | SUBROUTINE vect_mat_zstagg(vect, dx, dy, dz, wbdy, mat) |
---|
1306 | ! Subroutine to transform a LMDZ physic 1D vector to a 3D matrix |
---|
1307 | |
---|
1308 | IMPLICIT NONE |
---|
1309 | |
---|
1310 | INTEGER, INTENT(IN) :: dx, dy, dz, wbdy |
---|
1311 | REAL, DIMENSION((dx-2*wbdy)*(dy-2*wbdy),dz), INTENT(IN) :: vect |
---|
1312 | REAL, DIMENSION(1-wbdy:dx-wbdy+1,dz,1-wbdy:dy-wbdy+1), INTENT(OUT) :: mat |
---|
1313 | |
---|
1314 | ! Local |
---|
1315 | INTEGER :: i,j,k |
---|
1316 | INTEGER :: ddx, ddy |
---|
1317 | |
---|
1318 | ! Variables |
---|
1319 | ! d[x/y/z]: dimensions |
---|
1320 | ! mat: matrix to transform |
---|
1321 | ! vect: vector to produce |
---|
1322 | |
---|
1323 | ddx=dx-2*wbdy |
---|
1324 | ddy=dy-2*wbdy |
---|
1325 | |
---|
1326 | mat=0. |
---|
1327 | |
---|
1328 | DO k=1,dz |
---|
1329 | DO j=1,ddy |
---|
1330 | DO i=1,ddx |
---|
1331 | mat(i,k,j)=vect(ddx*(j-1)+i,k) |
---|
1332 | END DO |
---|
1333 | END DO |
---|
1334 | END DO |
---|
1335 | |
---|
1336 | END SUBROUTINE vect_mat_zstagg |
---|
1337 | |
---|
1338 | |
---|
1339 | SUBROUTINE mat_vect(mat, dx, dy, dz, wbdy, vect) |
---|
1340 | ! Subroutine to transform a 3D matrix to a LMDZ physic 1D vector |
---|
1341 | |
---|
1342 | IMPLICIT NONE |
---|
1343 | |
---|
1344 | INTEGER, INTENT(IN) :: dx, dy, dz, wbdy |
---|
1345 | REAL, DIMENSION(1-wbdy:dx-wbdy+1,dz+1,1-wbdy:dy-wbdy+1), INTENT(IN) :: mat |
---|
1346 | REAL, DIMENSION((dx-2*wbdy)*(dy-2*wbdy),dz), INTENT(OUT) :: & |
---|
1347 | vect |
---|
1348 | |
---|
1349 | ! Local |
---|
1350 | INTEGER :: i,j,k |
---|
1351 | INTEGER :: ddx,ddy |
---|
1352 | |
---|
1353 | ! Variables |
---|
1354 | ! d[x/y/z]: dimensions |
---|
1355 | ! mat: matrix to transform |
---|
1356 | ! vect: vector to produce |
---|
1357 | |
---|
1358 | ! Removing HALO |
---|
1359 | ddx=dx-2*wbdy |
---|
1360 | ddy=dy-2*wbdy |
---|
1361 | |
---|
1362 | vect=0. |
---|
1363 | |
---|
1364 | DO k=1,dz |
---|
1365 | DO j=1,ddy |
---|
1366 | DO i=1,ddx |
---|
1367 | vect(ddx*(j-1)+i,k)= mat(i,k,j) |
---|
1368 | END DO |
---|
1369 | END DO |
---|
1370 | END DO |
---|
1371 | |
---|
1372 | END SUBROUTINE mat_vect |
---|
1373 | |
---|
1374 | SUBROUTINE mat_vect_zstagg(mat, dx, dy, dz, wbdy, vect) |
---|
1375 | ! Subroutine to transform a 3D matrix z-staggered to a LMDZ physic 1D vector |
---|
1376 | |
---|
1377 | IMPLICIT NONE |
---|
1378 | |
---|
1379 | INTEGER, INTENT(IN) :: dx, dy, dz, wbdy |
---|
1380 | REAL, DIMENSION(1-wbdy:dx-wbdy+1,dz,1-wbdy:dy-wbdy+1), INTENT(IN) :: mat |
---|
1381 | REAL, DIMENSION((dx-2*wbdy)*(dy-2*wbdy),dz), INTENT(OUT) :: & |
---|
1382 | vect |
---|
1383 | |
---|
1384 | ! Local |
---|
1385 | INTEGER :: i,j,k |
---|
1386 | INTEGER :: ddx,ddy |
---|
1387 | |
---|
1388 | ! Variables |
---|
1389 | ! d[x/y/z]: dimensions |
---|
1390 | ! mat: matrix to transform |
---|
1391 | ! vect: vector to produce |
---|
1392 | |
---|
1393 | ! Removing HALO |
---|
1394 | ddx=dx-2*wbdy |
---|
1395 | ddy=dy-2*wbdy |
---|
1396 | |
---|
1397 | DO k=1,dz |
---|
1398 | DO j=1,ddy |
---|
1399 | DO i=1,ddx |
---|
1400 | vect(ddx*(j-1)+i,k)= mat(i,k,j) |
---|
1401 | END DO |
---|
1402 | END DO |
---|
1403 | END DO |
---|
1404 | |
---|
1405 | END SUBROUTINE mat_vect_zstagg |
---|
1406 | |
---|
1407 | SUBROUTINE WRFlanduse_LMDZKsoil(lndcn,is,ie,js,je,dx,dy,ncat,lndu,mask,kt,kl,ko,ks) |
---|
1408 | ! Subroutine to retrieve LMDZ Ksoil classification using WRF landuse |
---|
1409 | |
---|
1410 | IMPLICIT NONE |
---|
1411 | |
---|
1412 | CHARACTER(LEN=50), INTENT(IN) :: lndcn |
---|
1413 | INTEGER, INTENT(IN) :: is,ie,js,je,dx,dy,ncat |
---|
1414 | REAL, DIMENSION(is:ie,ncat,js:je), INTENT(IN) :: lndu |
---|
1415 | REAL, DIMENSION(is:ie,js:je), INTENT(IN) :: mask |
---|
1416 | REAL, DIMENSION(is:ie,js:je), INTENT(OUT) :: kt,kl,ko,ks |
---|
1417 | |
---|
1418 | ! Local |
---|
1419 | INTEGER :: i,j,k,ij,Nnones,Ncattot |
---|
1420 | CHARACTER(LEN=50) :: fname, errmsg, wrnmsg |
---|
1421 | REAL, DIMENSION(is:ie,js:je) :: totlndu |
---|
1422 | INTEGER, ALLOCATABLE, DIMENSION(:) :: ter_wk, lic_wk, oce_wk, & |
---|
1423 | sic_wk |
---|
1424 | LOGICAL, ALLOCATABLE, DIMENSION(:,:) :: wk |
---|
1425 | INTEGER :: Nter_wk, Nlic_wk, & |
---|
1426 | Noce_wk, Nsic_wk, klks_equal |
---|
1427 | REAL*8, DIMENSION(is:ie,js:je) :: Ksoilval |
---|
1428 | |
---|
1429 | |
---|
1430 | !!!!!!!! Variables |
---|
1431 | ! lndcn: land-use data-set (same as in VEGPARM.TBL) |
---|
1432 | ! dx,dy: horizontal dimension of the domain |
---|
1433 | ! ncat: number of categories of the data-set |
---|
1434 | ! lndu: WRF landuse matrix |
---|
1435 | ! mask: land-sea mask from WRF |
---|
1436 | ! k[t/l/o/s]: LMDZ ter, lic, oce, sic matrices |
---|
1437 | ! Ncattot: theoretical number of laduses for a given data-base |
---|
1438 | ! [ter/lic/oce/sic]_wk: vectors with the indices of the landuse for LMDZ Ksoil |
---|
1439 | ! N[ter/lic/oce/sic]_wk: number of indices of the landuse for LMDZ Ksoil |
---|
1440 | ! wk: boolean matrix for each Ksoil indicating 'true' if the WRF landuse is for Ksoil |
---|
1441 | |
---|
1442 | fname="'WRFlanduse_LMDZKsoil'" |
---|
1443 | errmsg='ERROR -- error -- ERROR -- error' |
---|
1444 | wrnmsg='WARNING -- warning -- WARNING -- warning' |
---|
1445 | |
---|
1446 | ! Checking if at all the grid points SUM(landuse) == 1. |
---|
1447 | !! |
---|
1448 | PRINT *,' ' // TRIM(fname) // " using '" // TRIM(lndcn) // "' WRF land types" |
---|
1449 | |
---|
1450 | totlndu = SUM(lndu, 2) |
---|
1451 | |
---|
1452 | Nnones=0 |
---|
1453 | DO i=1,dx |
---|
1454 | DO j=1,dy |
---|
1455 | IF (totlndu(i,j) /= 1. ) Nnones = Nnones+1 |
---|
1456 | END DO |
---|
1457 | END DO |
---|
1458 | |
---|
1459 | IF (Nnones /= 0) THEN |
---|
1460 | PRINT *,TRIM(errmsg) |
---|
1461 | PRINT *,' ' // TRIM(fname) // ': ', Nnones,' grid points do not have a ' // & |
---|
1462 | 'TOTlanduse == 1. from WRF !!!!!' |
---|
1463 | DO i=1,dx |
---|
1464 | DO j=1,dy |
---|
1465 | IF (totlndu(i,j) /= 1.) PRINT *,' ',i,', ',j,' total land use: ', & |
---|
1466 | totlndu(i,j),' indiv. values :', lndu(i,:,j) |
---|
1467 | END DO |
---|
1468 | END DO |
---|
1469 | END IF |
---|
1470 | |
---|
1471 | SELECT CASE (lndcn) |
---|
1472 | CASE ('OLD') |
---|
1473 | Ncattot = 13 |
---|
1474 | Nter_wk = 11 |
---|
1475 | Nlic_wk = 1 |
---|
1476 | Noce_wk = 1 |
---|
1477 | Nsic_wk = 1 |
---|
1478 | |
---|
1479 | IF (ALLOCATED(ter_wk)) DEALLOCATE(ter_wk) |
---|
1480 | ALLOCATE(ter_wk(Nter_wk)) |
---|
1481 | IF (ALLOCATED(lic_wk)) DEALLOCATE(lic_wk) |
---|
1482 | ALLOCATE(lic_wk(Nlic_wk)) |
---|
1483 | IF (ALLOCATED(oce_wk)) DEALLOCATE(oce_wk) |
---|
1484 | ALLOCATE(oce_wk(Noce_wk)) |
---|
1485 | IF (ALLOCATED(sic_wk)) DEALLOCATE(sic_wk) |
---|
1486 | ALLOCATE(sic_wk(Nsic_wk)) |
---|
1487 | |
---|
1488 | ! It could be done much more 'elegantly' ... reading values from LANDUSE.TBL (v3.3) |
---|
1489 | !! |
---|
1490 | ! 1,'Urban land' 2,'Agriculture' |
---|
1491 | ! 3,'Range-grassland' 4,'Deciduous forest' |
---|
1492 | ! 5,'Coniferous forest' 6,'Mixed forest and wet land' |
---|
1493 | ! 7,'Water' 8,'Marsh or wet land' |
---|
1494 | ! 9,'Desert' 10,'Tundra' |
---|
1495 | ! 11,'Permanent ice' 12,'Tropical or subtropical forest' |
---|
1496 | ! 13,'Savannah' |
---|
1497 | |
---|
1498 | ter_wk = (/ 1,2,3,4,5,6,8,9,10,12,13 /) |
---|
1499 | lic_wk = (/ 11 /) |
---|
1500 | oce_wk = (/ 7 /) |
---|
1501 | sic_wk = (/ 11 /) |
---|
1502 | |
---|
1503 | CASE ('USGS') |
---|
1504 | Ncattot = 33 |
---|
1505 | Nter_wk = 22 |
---|
1506 | Nlic_wk = 1 |
---|
1507 | Noce_wk = 1 |
---|
1508 | Nsic_wk = 1 |
---|
1509 | |
---|
1510 | IF (ALLOCATED(ter_wk)) DEALLOCATE(ter_wk) |
---|
1511 | ALLOCATE(ter_wk(Nter_wk)) |
---|
1512 | IF (ALLOCATED(lic_wk)) DEALLOCATE(lic_wk) |
---|
1513 | ALLOCATE(lic_wk(Nlic_wk)) |
---|
1514 | IF (ALLOCATED(oce_wk)) DEALLOCATE(oce_wk) |
---|
1515 | ALLOCATE(oce_wk(Noce_wk)) |
---|
1516 | IF (ALLOCATED(sic_wk)) DEALLOCATE(sic_wk) |
---|
1517 | ALLOCATE(sic_wk(Nsic_wk)) |
---|
1518 | |
---|
1519 | ! It could be done much more 'elegantly' ... reading values from LANDUSE.TBL (v3.3) |
---|
1520 | !! |
---|
1521 | ! 1,'Urban and Built-Up Land' 2,'Dryland Cropland and Pasture' |
---|
1522 | ! 3,'Irrigated Cropland and Pasture' 4,'Mixed Dryland/Irrigated Cropland and Pasture' |
---|
1523 | ! 5,'Cropland/Grassland Mosaic' 6,'Cropland/Woodland Mosaic' |
---|
1524 | ! 7,'Grassland' 8,'Shrubland' |
---|
1525 | ! 9,'Mixed Shrubland/Grassland' 10,'Savanna' |
---|
1526 | !11,'Deciduous Broadleaf Forest' 12,'Deciduous Needleleaf Forest' |
---|
1527 | !13,'Evergreen Broadleaf Forest' 14,'Evergreen Needleleaf Forest' |
---|
1528 | !15,'Mixed Forest' 16,'Water Bodies' |
---|
1529 | !17,'Herbaceous Wetland' 18,'Wooded Wetland' |
---|
1530 | !19,'Barren or Sparsely Vegetated' 20,'Herbaceous Tundra' |
---|
1531 | !21,'Wooded Tundra' 22,'Mixed Tundra' |
---|
1532 | !23,'Bare Ground Tundra' 24,'Snow or Ice' |
---|
1533 | !25,'Playa' 26,'Lava' |
---|
1534 | !27,'White Sand' 28,'Unassigned' |
---|
1535 | !29,'Unassigned' 30,'Unassigned' |
---|
1536 | !31,'Low Intensity Residential ' 32,'High Intensity Residential' |
---|
1537 | !33,'Industrial or Commercial' |
---|
1538 | |
---|
1539 | ! Correct values |
---|
1540 | ! ter_wk = (/ 1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,17,18,19,20,21,22,23,25,26, & |
---|
1541 | ! 27,31,32,33 /) |
---|
1542 | ! How ever, USGS by default has 24 so... |
---|
1543 | ter_wk = (/ 1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,17,18,19,20,21,22,23 /) |
---|
1544 | lic_wk = (/ 24 /) |
---|
1545 | oce_wk = (/ 16 /) |
---|
1546 | sic_wk = (/ 24 /) |
---|
1547 | |
---|
1548 | CASE ('MODIFIED_IGBP_MODIS_NOAH') |
---|
1549 | Ncattot = 33 |
---|
1550 | Nter_wk = 28 |
---|
1551 | Nlic_wk = 1 |
---|
1552 | Noce_wk = 1 |
---|
1553 | Nsic_wk = 1 |
---|
1554 | |
---|
1555 | IF (ALLOCATED(ter_wk)) DEALLOCATE(ter_wk) |
---|
1556 | ALLOCATE(ter_wk(Nter_wk)) |
---|
1557 | IF (ALLOCATED(lic_wk)) DEALLOCATE(lic_wk) |
---|
1558 | ALLOCATE(lic_wk(Nlic_wk)) |
---|
1559 | IF (ALLOCATED(oce_wk)) DEALLOCATE(oce_wk) |
---|
1560 | ALLOCATE(oce_wk(Noce_wk)) |
---|
1561 | IF (ALLOCATED(sic_wk)) DEALLOCATE(sic_wk) |
---|
1562 | ALLOCATE(sic_wk(Nsic_wk)) |
---|
1563 | |
---|
1564 | ! It could be done much more 'elegantly' ... reading values from LANDUSE.TBL (v3.3) |
---|
1565 | !! |
---|
1566 | ! 1,'Evergreen Needleleaf Forest' 2,'Evergreen Broadleaf Forest' |
---|
1567 | ! 3,'Deciduous Needleleaf Forest' 4,'Deciduous Broadleaf Forest' |
---|
1568 | ! 5,'Mixed Forests' 6,'Closed Shrublands' |
---|
1569 | ! 7,'Open Shrublands' 8,'Woody Savannas' |
---|
1570 | ! 9,'Savannas' 10,'Grasslands' |
---|
1571 | ! 11,'Permanent wetlands' 12,'Croplands' |
---|
1572 | ! 13,'Urban and Built-Up' 14,'cropland/natural vegetation mosaic' |
---|
1573 | ! 15,'Snow and Ice' 16,'Barren or Sparsely Vegetated' |
---|
1574 | ! 17,'Water' 18,'Wooded Tundra' |
---|
1575 | ! 19,'Mixed Tundra' 20,'Barren Tundra' |
---|
1576 | ! 21,'Unassigned' 22,'Unassigned' |
---|
1577 | ! 23,'Unassigned' 24,'Unassigned' |
---|
1578 | ! 25,'Unassigned' 26,'Unassigned' |
---|
1579 | ! 27,'Unassigned' 28,'Unassigned' |
---|
1580 | ! 29,'Unassigned' 30,'Unassigned' |
---|
1581 | ! 31,'Low Intensity Residential ' 32,'High Intensity Residential' |
---|
1582 | ! 33,'Industrial or Commercial' |
---|
1583 | |
---|
1584 | ter_wk = (/ 1,2,3,4,5,6,7,8,9,10,11,12,13,14,16,18,19,20,31,32,33 /) |
---|
1585 | lic_wk = (/ 15 /) |
---|
1586 | oce_wk = (/ 17 /) |
---|
1587 | sic_wk = (/ 15 /) |
---|
1588 | |
---|
1589 | CASE ('SiB') |
---|
1590 | Ncattot = 16 |
---|
1591 | Nter_wk = 14 |
---|
1592 | Nlic_wk = 1 |
---|
1593 | Noce_wk = 1 |
---|
1594 | Nsic_wk = 1 |
---|
1595 | |
---|
1596 | IF (ALLOCATED(ter_wk)) DEALLOCATE(ter_wk) |
---|
1597 | ALLOCATE(ter_wk(Nter_wk)) |
---|
1598 | IF (ALLOCATED(lic_wk)) DEALLOCATE(lic_wk) |
---|
1599 | ALLOCATE(lic_wk(Nlic_wk)) |
---|
1600 | IF (ALLOCATED(oce_wk)) DEALLOCATE(oce_wk) |
---|
1601 | ALLOCATE(oce_wk(Noce_wk)) |
---|
1602 | IF (ALLOCATED(sic_wk)) DEALLOCATE(sic_wk) |
---|
1603 | ALLOCATE(sic_wk(Nsic_wk)) |
---|
1604 | |
---|
1605 | ! It could be done much more 'elegantly' ... reading values from LANDUSE.TBL (v3.3) |
---|
1606 | !! |
---|
1607 | ! 1,'Evergreen Broadleaf Trees' 2,'Broadleaf Deciduous Trees' |
---|
1608 | ! 3,'Deciduous and Evergreen Trees' 4,'Evergreen Needleleaf Trees' |
---|
1609 | ! 5,'Deciduous Needleleaf Trees' 6,'Ground Cover with Trees and Shrubs' |
---|
1610 | ! 7,'Groundcover Only' 8,'Broadleaf Shrubs with Perennial Ground Cover' |
---|
1611 | ! 9,'Broadleaf Shrubs with Bare Soil' 10,'Groundcover with Dwarf Trees and Shrubs' |
---|
1612 | ! 11,'Bare Soil' 12,'Agriculture or C3 Grassland' |
---|
1613 | ! 13,'Persistent Wetland' 14,'Dry Coastal Complexes' |
---|
1614 | ! 15,'Water' 16,'Ice Cap and Glacier' |
---|
1615 | |
---|
1616 | ter_wk = (/ 1,2,3,4,5,6,7,8,9,10,11,12,13,14 /) |
---|
1617 | lic_wk = (/ 16 /) |
---|
1618 | oce_wk = (/ 15 /) |
---|
1619 | sic_wk = (/ 16 /) |
---|
1620 | |
---|
1621 | CASE ('LW12') |
---|
1622 | Ncattot = 3 |
---|
1623 | Nter_wk = 1 |
---|
1624 | Nlic_wk = 1 |
---|
1625 | Noce_wk = 1 |
---|
1626 | Nsic_wk = 1 |
---|
1627 | |
---|
1628 | IF (ALLOCATED(ter_wk)) DEALLOCATE(ter_wk) |
---|
1629 | ALLOCATE(ter_wk(Nter_wk)) |
---|
1630 | IF (ALLOCATED(lic_wk)) DEALLOCATE(lic_wk) |
---|
1631 | ALLOCATE(lic_wk(Nlic_wk)) |
---|
1632 | IF (ALLOCATED(oce_wk)) DEALLOCATE(oce_wk) |
---|
1633 | ALLOCATE(oce_wk(Noce_wk)) |
---|
1634 | IF (ALLOCATED(sic_wk)) DEALLOCATE(sic_wk) |
---|
1635 | ALLOCATE(sic_wk(Nsic_wk)) |
---|
1636 | |
---|
1637 | ! It could be done much more 'elegantly' ... reading values from LANDUSE.TBL (v3.3) |
---|
1638 | !! |
---|
1639 | ! 1,'Land' 2,'Water' |
---|
1640 | ! 3,'Snow or Ice' |
---|
1641 | |
---|
1642 | ter_wk = (/ 1 /) |
---|
1643 | lic_wk = (/ 3 /) |
---|
1644 | oce_wk = (/ 2 /) |
---|
1645 | sic_wk = (/ 3 /) |
---|
1646 | |
---|
1647 | CASE DEFAULT |
---|
1648 | PRINT *,TRIM(errmsg) |
---|
1649 | PRINT *,' ' // TRIM(fname) // ": landuse data-base '" // TRIM(lndcn) // & |
---|
1650 | "' not ready !!!" |
---|
1651 | STOP |
---|
1652 | |
---|
1653 | END SELECT |
---|
1654 | |
---|
1655 | ! Checking the correct number of categories |
---|
1656 | IF (Ncattot /= ncat) THEN |
---|
1657 | PRINT *,TRIM(wrnmsg) |
---|
1658 | PRINT *,' ' // TRIM(fname) // ": '" // TRIM(lndcn) // "' landuse data-base " & |
---|
1659 | // ' has ',Ncattot,' land-uses types ',ncat,' passed !!!!' |
---|
1660 | ! STOP |
---|
1661 | END IF |
---|
1662 | |
---|
1663 | ! Mimic WRF landuse with 4 LMDZ Ksoils |
---|
1664 | !! |
---|
1665 | IF (ALLOCATED(wk)) DEALLOCATE(wk) |
---|
1666 | ! ALLOCATE(wk(4,ncattot)) |
---|
1667 | ALLOCATE(wk(4,ncat)) |
---|
1668 | |
---|
1669 | wk = .FALSE. |
---|
1670 | |
---|
1671 | DO i=1,Nter_wk |
---|
1672 | wk(1,ter_wk(i)) = .TRUE. |
---|
1673 | END DO |
---|
1674 | DO i=1,Nlic_wk |
---|
1675 | wk(2,lic_wk(i)) = .TRUE. |
---|
1676 | END DO |
---|
1677 | DO i=1,Noce_wk |
---|
1678 | wk(3,oce_wk(i)) = .TRUE. |
---|
1679 | END DO |
---|
1680 | DO i=1,Nsic_wk |
---|
1681 | wk(4,sic_wk(i)) = .TRUE. |
---|
1682 | END DO |
---|
1683 | |
---|
1684 | klks_equal = 0 |
---|
1685 | DO i=1,dx |
---|
1686 | DO j=1,dy |
---|
1687 | ! Using double precission numbers to avoid truncation errors... |
---|
1688 | Ksoilval(1,1) = SUM(lndu(i,:,j),MASK=wk(1,:)) |
---|
1689 | kt(i,j) = Ksoilval(1,1) |
---|
1690 | Ksoilval(1,1) = SUM(lndu(i,:,j),MASK=wk(2,:)) |
---|
1691 | kl(i,j) = Ksoilval(1,1) |
---|
1692 | Ksoilval(1,1) = SUM(lndu(i,:,j),MASK=wk(3,:)) |
---|
1693 | ko(i,j) = Ksoilval(1,1) |
---|
1694 | Ksoilval(1,1) = SUM(lndu(i,:,j),MASK=wk(4,:)) |
---|
1695 | ks(i,j) = Ksoilval(1,1) |
---|
1696 | |
---|
1697 | IF (kl(i,j) == ks(i,j)) klks_equal = klks_equal+1 |
---|
1698 | END DO |
---|
1699 | END DO |
---|
1700 | |
---|
1701 | ! Assigning lic/sic to all that ice/snow according to the percentage of ter/oce |
---|
1702 | ! Is there any other better way? Satellites, do not know what is below the ice ?! |
---|
1703 | |
---|
1704 | IF (klks_equal == dx*dy) THEN |
---|
1705 | ! Assuming that kl=ks=ice |
---|
1706 | DO i=1,dx |
---|
1707 | DO j=1,dy |
---|
1708 | IF ((kt(i,j)+ko(i,j)) /= 0.) THEN |
---|
1709 | Ksoilval(1,1) = kl(i,j)*kt(i,j)/(kt(i,j)+ko(i,j)) |
---|
1710 | kl(i,j) = Ksoilval(1,1) |
---|
1711 | Ksoilval(1,1) = ks(i,j)*ko(i,j)/(kt(i,j)+ko(i,j)) |
---|
1712 | ks(i,j) = Ksoilval(1,1) |
---|
1713 | ELSE |
---|
1714 | Ksoilval(1,1) = kl(i,j)/(kl(i,j)+ks(i,j)) |
---|
1715 | Ksoilval(1,2) = ks(i,j)/(kl(i,j)+ks(i,j)) |
---|
1716 | kl(i,j) = Ksoilval(1,1) |
---|
1717 | ks(i,j) = Ksoilval(1,2) |
---|
1718 | END IF |
---|
1719 | ! Introducing WRF land-sea mask |
---|
1720 | IF (mask(i,j) == 1.) THEN |
---|
1721 | kl(i,j)=kl(i,j)+ks(i,j) |
---|
1722 | ks(i,j)=0. |
---|
1723 | ELSE |
---|
1724 | ks(i,j)=kl(i,j)+ks(i,j) |
---|
1725 | kl(i,j)=0. |
---|
1726 | END IF |
---|
1727 | END DO |
---|
1728 | END DO |
---|
1729 | ELSE |
---|
1730 | PRINT *,' ' // TRIM(fname) // ': land-use data provided different lic & sic !' |
---|
1731 | END IF |
---|
1732 | |
---|
1733 | ! Checking if at all the grid points kt+kl+ko+ks == 1. |
---|
1734 | !! |
---|
1735 | totlndu = kt+kl+ko+ks |
---|
1736 | |
---|
1737 | Nnones=0 |
---|
1738 | DO i=1,dx |
---|
1739 | DO j=1,dy |
---|
1740 | IF (totlndu(i,j) /= 1. ) Nnones = Nnones+1 |
---|
1741 | IF (kt(i,j) + ko(i,j) == 0.)THEN |
---|
1742 | PRINT *,TRIM(wrnmsg) |
---|
1743 | PRINT *,' ' // TRIM(fname) // ': point without land or water category !!' |
---|
1744 | PRINT *,'i j WRF_ter ter ______________' |
---|
1745 | DO k=1,Nter_wk |
---|
1746 | PRINT *,i,j,lndu(i,ter_wk(k),j), kt(i,j) |
---|
1747 | END DO |
---|
1748 | PRINT *,'i j WRF_oce oce ______________' |
---|
1749 | DO k=1,Noce_wk |
---|
1750 | PRINT *,i,j,lndu(i,oce_wk(k),j), ko(i,j) |
---|
1751 | END DO |
---|
1752 | PRINT *,' wrf%landusef: ', lndu(i,:,j) |
---|
1753 | ! STOP |
---|
1754 | END IF |
---|
1755 | END DO |
---|
1756 | END DO |
---|
1757 | |
---|
1758 | ij=1 |
---|
1759 | IF (Nnones /= 0) THEN |
---|
1760 | PRINT *,TRIM(errmsg) |
---|
1761 | PRINT *,' ' // TRIM(fname) // ': ', Nnones,' grid points do not have a ' // & |
---|
1762 | 'TOTlanduse == 1.!!!!!' |
---|
1763 | PRINT *,'pt i j ter+lic+oce+sic : ter lic oce sic ______________' |
---|
1764 | DO i=1,dx |
---|
1765 | DO j=1,dy |
---|
1766 | IF (totlndu(i,j) /= 1.) THEN |
---|
1767 | PRINT *,ij,'; ',i,j,totlndu(i,j),' : ',kt(i,j),kl(i,j),ko(i,j),ks(i,j) |
---|
1768 | ij=ij+1 |
---|
1769 | END IF |
---|
1770 | END DO |
---|
1771 | END DO |
---|
1772 | END IF |
---|
1773 | |
---|
1774 | DEALLOCATE(ter_wk) |
---|
1775 | DEALLOCATE(lic_wk) |
---|
1776 | DEALLOCATE(oce_wk) |
---|
1777 | DEALLOCATE(sic_wk) |
---|
1778 | DEALLOCATE(wk) |
---|
1779 | |
---|
1780 | RETURN |
---|
1781 | |
---|
1782 | END SUBROUTINE WRFlanduse_LMDZKsoil |
---|
1783 | |
---|
1784 | SUBROUTINE interpolate_layers(NvalsA,posA,valsA,NvalsB,posB,valsB) |
---|
1785 | ! Subroutine to interpolate values between different layers (constant values along the depth of a layer) |
---|
1786 | |
---|
1787 | IMPLICIT NONE |
---|
1788 | |
---|
1789 | INTEGER, INTENT(IN) :: NvalsA,NvalsB |
---|
1790 | REAL, DIMENSION(NvalsA), INTENT(IN) :: posA,valsA |
---|
1791 | REAL, DIMENSION(NvalsB), INTENT(IN) :: posB |
---|
1792 | REAL, DIMENSION(NvalsB), INTENT(OUT) :: valsB |
---|
1793 | |
---|
1794 | ! Local |
---|
1795 | INTEGER :: i,k |
---|
1796 | REAL :: xA,yA,xB,yB |
---|
1797 | |
---|
1798 | !!!!!!! Variables |
---|
1799 | ! A: original values used to interpolate |
---|
1800 | ! B: values to obtain |
---|
1801 | ! Nvals[A/B]: number of values |
---|
1802 | ! pos[A/B]: depths of each layer |
---|
1803 | ! vals[A/B]: values of each layer |
---|
1804 | |
---|
1805 | !!!!!!! Functions/Subroutines |
---|
1806 | ! interpolate_lin: subroutine to linearly interpolate |
---|
1807 | |
---|
1808 | i=1 |
---|
1809 | DO k=1,NvalsB |
---|
1810 | IF (posB(k) > posA(i)) i = i+1 |
---|
1811 | IF (i-1 < 1) THEN |
---|
1812 | xA = posA(i) |
---|
1813 | yA = valsA(i) |
---|
1814 | xB = posA(i+1) |
---|
1815 | yB = valsA(i+1) |
---|
1816 | ELSEIF ( i > NvalsA) THEN |
---|
1817 | xA = posA(i-2) |
---|
1818 | yA = valsA(i-2) |
---|
1819 | xB = posA(i-1) |
---|
1820 | yB = valsA(i-1) |
---|
1821 | ELSE |
---|
1822 | xA = posA(i-1) |
---|
1823 | yA = valsA(i-1) |
---|
1824 | xB = posA(i) |
---|
1825 | yB = valsA(i) |
---|
1826 | END IF |
---|
1827 | IF (posB(k) <= xB .AND. i-1 > 1) THEN |
---|
1828 | CALL interpolate_lin(xA,yA,xB,yB,posB(k),valsB(k)) |
---|
1829 | ELSEIF (posB(k) <= xA .AND. i-1 < 1) THEN |
---|
1830 | valsB(k) = yA + (yB - yA)*(posB(k) - xA)/(xB - xA) |
---|
1831 | ELSE |
---|
1832 | valsB(k) = yB + (yB - yA)*(posB(k) - xB)/(xB - xA) |
---|
1833 | END IF |
---|
1834 | ! PRINT *,'i ',i,' Interpolating for: ',posB(k),' with xA: ', xA,' yA: ',yA,' & ', & |
---|
1835 | ! ' xB: ', xB,' yB: ',yB,' --> ',valsB(k) |
---|
1836 | END DO |
---|
1837 | |
---|
1838 | END SUBROUTINE interpolate_layers |
---|
1839 | |
---|
1840 | SUBROUTINE interpolate_lin(x1,y1,x2,y2,x0,y) |
---|
1841 | ! Program to interpolate values y=a+b*x with: (from clWRF modifications) |
---|
1842 | ! a=y1 |
---|
1843 | ! b=(y2-y1)/(x2-x1) |
---|
1844 | ! x=(x2-x0) |
---|
1845 | |
---|
1846 | IMPLICIT NONE |
---|
1847 | |
---|
1848 | REAL, INTENT (IN) :: x1,x2,x0 |
---|
1849 | ! REAL(r8), INTENT (IN) :: y1,y2 |
---|
1850 | ! REAL(r8), INTENT (OUT) :: y |
---|
1851 | ! REAL(r8) :: a,b,x |
---|
1852 | REAL, INTENT (IN) :: y1,y2 |
---|
1853 | REAL, INTENT (OUT) :: y |
---|
1854 | REAL :: a,b,x |
---|
1855 | |
---|
1856 | a=y1 |
---|
1857 | b=(y2-y1)/(x2-x1) |
---|
1858 | |
---|
1859 | IF (x0 .GE. x1) THEN |
---|
1860 | x=x0-x1 |
---|
1861 | ELSE |
---|
1862 | x=x1-x0 |
---|
1863 | b=-b |
---|
1864 | ENDIF |
---|
1865 | |
---|
1866 | y=a+b*x |
---|
1867 | |
---|
1868 | ! PRINT *,'y=a+b*x' |
---|
1869 | ! IF (b .LT. 0.) THEN |
---|
1870 | ! PRINT *,'y=',y1,' ',b,'*x' |
---|
1871 | ! ELSE |
---|
1872 | ! PRINT *,'y=',y1,'+',b,'*x' |
---|
1873 | ! ENDIF |
---|
1874 | ! PRINT *,'y(',x,')=',y |
---|
1875 | |
---|
1876 | END SUBROUTINE interpolate_lin |
---|
1877 | |
---|
1878 | SUBROUTINE table_characteristics(datasetn,tablen,datasetline,Nvar,Nr,Nc) |
---|
1879 | ! Subroutine to read values from a WRF .TBL |
---|
1880 | |
---|
1881 | IMPLICIT NONE |
---|
1882 | |
---|
1883 | CHARACTER(LEN=50), INTENT(IN) :: datasetn |
---|
1884 | CHARACTER(LEN=100), INTENT(IN) :: tablen |
---|
1885 | INTEGER, INTENT(OUT) :: datasetline, Nvar, Nr, Nc |
---|
1886 | |
---|
1887 | ! Local |
---|
1888 | INTEGER :: il, iline, ic |
---|
1889 | CHARACTER(LEN=50) :: fname, errmsg |
---|
1890 | CHARACTER(LEN=200) :: lineval |
---|
1891 | LOGICAL :: existsfile, is_used |
---|
1892 | INTEGER :: funit, Llineval, icoma |
---|
1893 | |
---|
1894 | !!!!!! Variables |
---|
1895 | ! datasetn: name of the data set to use |
---|
1896 | ! tablen: name of the table to read |
---|
1897 | ! Nvar: number of variations (g.e. seasons) |
---|
1898 | ! Nr: table row number |
---|
1899 | ! Nc: table column number |
---|
1900 | ! datasetline: line in the table where the dataset starts |
---|
1901 | |
---|
1902 | fname='table_characteristics' |
---|
1903 | errmsg='ERROR -- error -- ERROR -- error' |
---|
1904 | |
---|
1905 | INQUIRE(FILE=TRIM(tablen), EXIST=existsfile) |
---|
1906 | |
---|
1907 | IF (existsfile) THEN |
---|
1908 | DO funit=10,100 |
---|
1909 | INQUIRE(UNIT=funit, OPENED=is_used) |
---|
1910 | IF (.NOT. is_used) EXIT |
---|
1911 | END DO |
---|
1912 | |
---|
1913 | ELSE |
---|
1914 | PRINT *,TRIM(errmsg) |
---|
1915 | PRINT *,TRIM(fname),': file "'// TRIM(tablen) // '" does not exist !!!' |
---|
1916 | STOP |
---|
1917 | END IF |
---|
1918 | |
---|
1919 | OPEN(UNIT=funit, FILE=TRIM(tablen), STATUS='OLD') |
---|
1920 | |
---|
1921 | ! Looking for the place where data-set starts |
---|
1922 | !! |
---|
1923 | datasetline = -1 |
---|
1924 | iline=1 |
---|
1925 | DO il=1,100000 |
---|
1926 | READ(funit,*,END=100)lineval |
---|
1927 | IF (TRIM(lineval) == TRIM(datasetn)) datasetline = il |
---|
1928 | END DO |
---|
1929 | |
---|
1930 | 100 CONTINUE |
---|
1931 | IF (datasetline == -1) THEN |
---|
1932 | PRINT *,TRIM(errmsg) |
---|
1933 | PRINT *,TRIM(fname),': in file "'// TRIM(tablen) // '" dataset: "'// & |
---|
1934 | TRIM(datasetn) // '" does not exist !!!' |
---|
1935 | STOP |
---|
1936 | END IF |
---|
1937 | |
---|
1938 | ! Looking for number of values |
---|
1939 | !! |
---|
1940 | REWIND(UNIT=funit) |
---|
1941 | DO ic=1,datasetline |
---|
1942 | READ(funit,*)lineval |
---|
1943 | END DO |
---|
1944 | READ(funit,*)Nr,Nvar |
---|
1945 | |
---|
1946 | IF (Nvar > 1) READ(funit,*)lineval |
---|
1947 | |
---|
1948 | READ(funit,'(200A)')lineval |
---|
1949 | Llineval=LEN(TRIM(lineval)) |
---|
1950 | |
---|
1951 | ! PRINT *,' lineval: ******' // TRIM(lineval) // '******', Llineval |
---|
1952 | |
---|
1953 | Nc = 0 |
---|
1954 | ic = 1 |
---|
1955 | DO WHILE (ic <= Llineval) |
---|
1956 | ! PRINT *,lineval(ic:Llineval) |
---|
1957 | icoma = SCAN(lineval(ic:Llineval),',') |
---|
1958 | IF (icoma > 1) THEN |
---|
1959 | Nc = Nc + 1 |
---|
1960 | ! PRINT *,' Nc: ',Nc,' ic ',ic,' icoma: ',icoma |
---|
1961 | ic = ic + icoma |
---|
1962 | ELSE |
---|
1963 | ic = ic + Llineval |
---|
1964 | END IF |
---|
1965 | END DO |
---|
1966 | PRINT *,' dataset :"' // TRIM(datasetn) // '" at line ',datasetline, & |
---|
1967 | ' in table :"'// TRIM(tablen) // '" has N variations: ',Nvar, & |
---|
1968 | ' Number of types: ',Nr, ' Number of values: ',Nc |
---|
1969 | |
---|
1970 | CLOSE(funit) |
---|
1971 | |
---|
1972 | RETURN |
---|
1973 | |
---|
1974 | END SUBROUTINE table_characteristics |
---|
1975 | |
---|
1976 | SUBROUTINE read_table(datasetn,tablen,datasetline,Nvar,Nr,Nc,vals,valns) |
---|
1977 | ! Subroutine to read values from a WRF .TBL |
---|
1978 | |
---|
1979 | IMPLICIT NONE |
---|
1980 | |
---|
1981 | CHARACTER(LEN=50), INTENT(IN) :: datasetn |
---|
1982 | CHARACTER(LEN=100), INTENT(IN) :: tablen |
---|
1983 | INTEGER, INTENT(IN) :: datasetline, Nvar, Nr, Nc |
---|
1984 | REAL, DIMENSION(Nvar,Nr,Nc), INTENT(OUT) :: vals |
---|
1985 | CHARACTER(LEN=50), DIMENSION(Nr), INTENT(OUT) :: valns |
---|
1986 | |
---|
1987 | |
---|
1988 | ! Local |
---|
1989 | INTEGER :: iv, ir, ic |
---|
1990 | CHARACTER(LEN=50) :: fname, errmsg, line |
---|
1991 | LOGICAL :: existsfile, is_used |
---|
1992 | INTEGER :: funit |
---|
1993 | |
---|
1994 | !!!!!!! Variables |
---|
1995 | ! datasetn: name of the data set to use |
---|
1996 | ! tablen: name of the table to read |
---|
1997 | ! Nvar: number of variations (g.e. seasons) |
---|
1998 | ! Nr: table row number |
---|
1999 | ! Nc: table column number |
---|
2000 | ! datasetline: line in the table where the dataset starts |
---|
2001 | ! vals: values of the data-set |
---|
2002 | ! valns: assigned name of each value in the dataset |
---|
2003 | |
---|
2004 | fname='read_table' |
---|
2005 | errmsg='ERROR -- error -- ERROR -- error' |
---|
2006 | |
---|
2007 | INQUIRE(FILE=TRIM(tablen), EXIST=existsfile) |
---|
2008 | |
---|
2009 | IF (existsfile) THEN |
---|
2010 | DO funit=10,100 |
---|
2011 | INQUIRE(UNIT=funit, OPENED=is_used) |
---|
2012 | IF (.NOT. is_used) EXIT |
---|
2013 | END DO |
---|
2014 | |
---|
2015 | ELSE |
---|
2016 | PRINT *,TRIM(errmsg) |
---|
2017 | PRINT *,TRIM(fname),': file "'// TRIM(tablen) // '" does not exist !!!' |
---|
2018 | STOP |
---|
2019 | END IF |
---|
2020 | |
---|
2021 | OPEN(UNIT=funit, FILE=TRIM(tablen), STATUS='OLD') |
---|
2022 | |
---|
2023 | DO ic=1,datasetline+1 |
---|
2024 | READ(funit,*)line |
---|
2025 | END DO |
---|
2026 | |
---|
2027 | DO iv=1, Nvar |
---|
2028 | IF (Nvar > 1) READ(funit,*)line |
---|
2029 | DO ir=1, Nr |
---|
2030 | READ(funit,*)(vals(iv,ir,ic), ic=1,Nc),valns(ir) |
---|
2031 | END DO |
---|
2032 | END DO |
---|
2033 | |
---|
2034 | ! PRINT *,' values dataset :"' // TRIM(datasetn) // '" in table :"' //TRIM(tablen) |
---|
2035 | ! DO iv=1,Nvar |
---|
2036 | ! PRINT *,' variation: ',iv |
---|
2037 | ! DO ir=1,Nr |
---|
2038 | ! PRINT *,' ',(vals(iv,ir,ic), ic=1,Nc),valns(ir) |
---|
2039 | ! END DO |
---|
2040 | ! END DO |
---|
2041 | |
---|
2042 | CLOSE(funit) |
---|
2043 | |
---|
2044 | RETURN |
---|
2045 | |
---|
2046 | END SUBROUTINE read_table |
---|
2047 | |
---|
2048 | |
---|
2049 | SUBROUTINE compute_soil_mass(Nvar,Nr,Nc,vals,valns,soilt,Nnost,NOsoilt,smois,ddx, & |
---|
2050 | ddy,ddz,soilm,soilq) |
---|
2051 | ! Subroutine to compute the total soil mass (kg) from a point with different |
---|
2052 | ! percentages of soil types |
---|
2053 | |
---|
2054 | IMPLICIT NONE |
---|
2055 | |
---|
2056 | INTEGER, INTENT(IN) :: Nvar, Nr, Nc, Nnost |
---|
2057 | REAL, DIMENSION(Nvar,Nr,Nc), INTENT(IN) :: vals |
---|
2058 | CHARACTER(LEN=50), DIMENSION(Nr), INTENT(IN) :: valns |
---|
2059 | REAL, DIMENSION(Nr), INTENT(IN) :: soilt |
---|
2060 | CHARACTER(LEN=50), DIMENSION(Nnost), INTENT(IN) :: NOsoilt |
---|
2061 | REAL, INTENT(IN) :: smois, ddx, ddy, ddz |
---|
2062 | REAL, INTENT(OUT) :: soilm,soilq |
---|
2063 | |
---|
2064 | ! Local |
---|
2065 | INTEGER :: iv, ir, ic |
---|
2066 | CHARACTER(LEN=50) :: fname, errmsg |
---|
2067 | LOGICAL :: is_land |
---|
2068 | |
---|
2069 | !!!!!!! Variables |
---|
2070 | ! Nvar: number of variations (g.e. seasons) |
---|
2071 | ! Nr: table row number |
---|
2072 | ! Nc: table column number |
---|
2073 | ! vals: values of the data-set (col2: density [g cm-3]) |
---|
2074 | ! valns: assigned name of each value in the dataset |
---|
2075 | ! soilt: vector with the percentages of soil types |
---|
2076 | ! Nnost: number of soil types which are not land |
---|
2077 | ! NOsoilt: labels of the non-land soil types |
---|
2078 | ! smois: soil moisture [m3 m-3] |
---|
2079 | ! ddx,ddy,ddz: volume of the grid cell [m] |
---|
2080 | ! soilm: total soil mass [kg] |
---|
2081 | ! soilq: total soil water [kg kg-1] |
---|
2082 | |
---|
2083 | fname='read_table' |
---|
2084 | errmsg='ERROR -- error -- ERROR -- error' |
---|
2085 | |
---|
2086 | soilm = 0. |
---|
2087 | soilq = 0. |
---|
2088 | |
---|
2089 | ! PRINT *,' soilmoisture: ',smois |
---|
2090 | DO ir=1,Nr |
---|
2091 | is_land = .TRUE. |
---|
2092 | IF (soilt(ir) /= 0.) THEN |
---|
2093 | DO iv=1, Nnost |
---|
2094 | IF (TRIM(valns(ir)) == TRIM(NOsoilt(iv))) is_land = .FALSE. |
---|
2095 | END DO |
---|
2096 | IF (is_land) THEN |
---|
2097 | ! PRINT *,' soil moisture range: ', vals(1,ir,3), vals(1,ir,5) |
---|
2098 | soilm = soilm + soilt(ir)*vals(1,ir,2)*ddx*ddy*ddz*1000. |
---|
2099 | soilq = soilq + soilt(ir)*vals(1,ir,2)*smois*ddx*ddy*ddz*1000. |
---|
2100 | ! PRINT *,' Percentage of soil "' // TRIM(valns(ir)) //'": ',soilt(ir)*100.,& |
---|
2101 | ! ' % density: ',vals(1,ir,2),' partial mass: ', & |
---|
2102 | ! soilt(ir)*vals(1,ir,2)*ddx*ddy*ddz*1000. , ' partial water: ', & |
---|
2103 | ! soilt(ir)*vals(1,ir,2)*smois*ddx*ddy*ddz*1000. |
---|
2104 | END IF |
---|
2105 | END IF |
---|
2106 | END DO |
---|
2107 | ! Remember 1 kg H20 = 1 l = 1mm / m2 |
---|
2108 | !! PRINT *,' TOTAL. soil mass: ',soilm,' [kg] soil water: ',soilq,' [kg]', & |
---|
2109 | !! ' equiv. water depth in soil ',soilq/(ddx*ddy*1000.),' [m]' |
---|
2110 | |
---|
2111 | END SUBROUTINE compute_soil_mass |
---|
2112 | |
---|
2113 | SUBROUTINE compute_TOTsoil_mass_water(datname, tabname, Nsoiltypes, soiltypest, & |
---|
2114 | soiltypesb, Nsoillayers, smois, soillayers, sizeX, sizeY, soilTOTmass, & |
---|
2115 | soilTOThumidity) |
---|
2116 | ! Subroutine to compute total soil mass and water from a given grid point using a |
---|
2117 | ! given data-set from SOILPARM.TBL |
---|
2118 | |
---|
2119 | IMPLICIT NONE |
---|
2120 | |
---|
2121 | CHARACTER(LEN=50), INTENT(IN) :: datname |
---|
2122 | CHARACTER(LEN=100), INTENT(IN) :: tabname |
---|
2123 | INTEGER, INTENT(IN) :: Nsoiltypes, Nsoillayers |
---|
2124 | REAL, DIMENSION(Nsoiltypes), INTENT(IN) :: soiltypest, soiltypesb |
---|
2125 | REAL, DIMENSION(Nsoillayers), INTENT(IN) :: smois, soillayers |
---|
2126 | REAL, INTENT(IN) :: sizeX, sizeY |
---|
2127 | REAL, INTENT(OUT) :: soilTOTmass, & |
---|
2128 | soilTOThumidity |
---|
2129 | |
---|
2130 | ! Local |
---|
2131 | INTEGER :: ix,iy,iz |
---|
2132 | REAL, ALLOCATABLE, DIMENSION(:,:,:) :: values |
---|
2133 | CHARACTER(LEN=50), ALLOCATABLE, DIMENSION(:) :: namevalues,NOsoiltypes |
---|
2134 | INTEGER :: linedataset,Nvariations, & |
---|
2135 | Nrows, Ncols |
---|
2136 | INTEGER :: NNOsoiltypes |
---|
2137 | REAL, ALLOCATABLE, DIMENSION(:) :: soiltypes |
---|
2138 | REAL :: dz,soilmass,soilhumidity |
---|
2139 | CHARACTER(LEN=50) :: fname, errmsg |
---|
2140 | |
---|
2141 | !!!!!!! Variables |
---|
2142 | ! datname: name of the data set to use |
---|
2143 | ! tabname: name of the table to read |
---|
2144 | ! values: mmatrix with the name of the values: variation (Nvariations, g.e. seasons), |
---|
2145 | ! table row (Nrows), table column (Ncols) |
---|
2146 | ! linedataset: line inside data-set where values start |
---|
2147 | ! namevalues: string assigned to each table row |
---|
2148 | ! Nsoiltypes: number of soil types |
---|
2149 | ! soiltypes[t/b]: percentage of the soil types at a given grid point on top/bottom |
---|
2150 | ! soiltypes: percentage of the soil types at a given grid point (lineray interpolated) |
---|
2151 | ! NOsoiltypes: vector with the names of soil types which are not land |
---|
2152 | |
---|
2153 | fname='"compute_TOTsoil_mass_water"' |
---|
2154 | errmsg='ERROR -- error -- ERROR -- error' |
---|
2155 | |
---|
2156 | CALL table_characteristics(datname, tabname, linedataset, Nvariations, & |
---|
2157 | Nrows, Ncols) |
---|
2158 | |
---|
2159 | IF (Nrows /= Nsoiltypes) THEN |
---|
2160 | PRINT *,errmsg |
---|
2161 | PRINT ' ' // TRIM(fname) // ': Number of soil-types: ',Nsoiltypes, & |
---|
2162 | ' does not concide with the number of values: ', Nrows, & |
---|
2163 | ' from the data-set "' // TRIM(datname) // '" !!!' |
---|
2164 | STOP |
---|
2165 | END IF |
---|
2166 | |
---|
2167 | IF (ALLOCATED(values)) DEALLOCATE(values) |
---|
2168 | ALLOCATE(values(Nvariations,Nrows,Ncols)) |
---|
2169 | |
---|
2170 | IF (ALLOCATED(namevalues)) DEALLOCATE(namevalues) |
---|
2171 | ALLOCATE(namevalues(Nrows)) |
---|
2172 | |
---|
2173 | CALL read_table(datname, tabname, linedataset, Nvariations, Nrows, Ncols, & |
---|
2174 | values, namevalues) |
---|
2175 | |
---|
2176 | IF (ALLOCATED(soiltypes)) DEALLOCATE(soiltypes) |
---|
2177 | ALLOCATE(soiltypes(Nrows)) |
---|
2178 | |
---|
2179 | ! On WRF3.3 SOILPARM.TBL |
---|
2180 | SELECT CASE (TRIM(datname)) |
---|
2181 | |
---|
2182 | CASE ('STAS') |
---|
2183 | NNOsoiltypes = 1 |
---|
2184 | IF (ALLOCATED(NOsoiltypes)) DEALLOCATE(NOsoiltypes) |
---|
2185 | ALLOCATE(NOsoiltypes(NNOsoiltypes)) |
---|
2186 | |
---|
2187 | NOsoiltypes=(/ 'WATER' /) |
---|
2188 | CASE ('STAS-RUC') |
---|
2189 | NNOsoiltypes = 1 |
---|
2190 | IF (ALLOCATED(NOsoiltypes)) DEALLOCATE(NOsoiltypes) |
---|
2191 | ALLOCATE(NOsoiltypes( NNOsoiltypes)) |
---|
2192 | |
---|
2193 | NOsoiltypes=(/ 'WATER' /) |
---|
2194 | END SELECT |
---|
2195 | |
---|
2196 | soilTOTmass = 0. |
---|
2197 | soilTOThumidity = 0. |
---|
2198 | |
---|
2199 | DO iz=1,Nsoillayers |
---|
2200 | IF (iz == 1 ) THEN |
---|
2201 | dz = soillayers(iz) |
---|
2202 | ELSE |
---|
2203 | dz = soillayers(iz) - soillayers(iz - 1) |
---|
2204 | END IF |
---|
2205 | IF (( iz > 1) .AND. (iz < Nsoillayers)) THEN |
---|
2206 | DO ix=1,Nrows |
---|
2207 | CALL interpolate_lin(soillayers(1), soiltypest(ix), & |
---|
2208 | soillayers(Nsoillayers), soiltypesb(ix), soillayers(iz), soiltypes(ix)) |
---|
2209 | END DO |
---|
2210 | ELSEIF (iz == 1) THEN |
---|
2211 | soiltypes = soiltypest |
---|
2212 | ELSE |
---|
2213 | soiltypes = soiltypesb |
---|
2214 | END IF |
---|
2215 | ! PRINT *,'lev: ',iz,' dz: ',dz,' soil types: ' ,soiltypes |
---|
2216 | |
---|
2217 | CALL compute_soil_mass(Nvariations,Nrows,Ncols,values,namevalues,soiltypes, & |
---|
2218 | NNOsoiltypes,NOsoiltypes,smois(iz),1.,1.,dz,soilmass,soilhumidity) |
---|
2219 | soilTOTmass = soilTOTmass + soilmass |
---|
2220 | soilTOThumidity = soilTOThumidity + soilhumidity |
---|
2221 | END DO |
---|
2222 | |
---|
2223 | DEALLOCATE(values) |
---|
2224 | DEALLOCATE(namevalues) |
---|
2225 | DEALLOCATE(soiltypes) |
---|
2226 | DEALLOCATE(NOsoiltypes) |
---|
2227 | |
---|
2228 | END SUBROUTINE compute_TOTsoil_mass_water |
---|
2229 | |
---|
2230 | SUBROUTINE compute_TOTsoil_mass_water_values(datname, Nvariations, Nrows, Ncols, & |
---|
2231 | values, namevalues, ip, jp, Nsoiltypes, soiltypest, soiltypesb, Nsoillayers, & |
---|
2232 | smois, soillayers, sizeX, sizeY, soilTOTmass, soilTOThumidity) |
---|
2233 | ! Subroutine to compute total soil mass and water from a given grid point using the |
---|
2234 | ! values from a given data-set of SOILPARM.TBL |
---|
2235 | |
---|
2236 | IMPLICIT NONE |
---|
2237 | |
---|
2238 | CHARACTER(LEN=50), INTENT(IN) :: datname |
---|
2239 | INTEGER, INTENT(IN) :: Nvariations, Nrows, Ncols |
---|
2240 | REAL, DIMENSION(Nvariations, Nrows, Ncols), INTENT(IN) :: values |
---|
2241 | CHARACTER(LEN=50), DIMENSION(Nrows), INTENT(IN) :: namevalues |
---|
2242 | INTEGER, INTENT(IN) :: ip, jp |
---|
2243 | INTEGER, INTENT(IN) :: Nsoiltypes, Nsoillayers |
---|
2244 | REAL, DIMENSION(Nsoiltypes), INTENT(IN) :: soiltypest, soiltypesb |
---|
2245 | REAL, DIMENSION(Nsoillayers), INTENT(IN) :: smois, soillayers |
---|
2246 | REAL, INTENT(IN) :: sizeX, sizeY |
---|
2247 | REAL, INTENT(OUT) :: soilTOTmass, & |
---|
2248 | soilTOThumidity |
---|
2249 | |
---|
2250 | ! Local |
---|
2251 | INTEGER :: ix,iy,iz |
---|
2252 | CHARACTER(LEN=50), ALLOCATABLE, DIMENSION(:) :: NOsoiltypes |
---|
2253 | INTEGER :: linedataset |
---|
2254 | INTEGER :: NNOsoiltypes |
---|
2255 | REAL, ALLOCATABLE, DIMENSION(:) :: soiltypes |
---|
2256 | REAL :: dz,soilmass,soilhumidity |
---|
2257 | CHARACTER(LEN=50) :: fname, errmsg |
---|
2258 | |
---|
2259 | !!!!!!! Variables |
---|
2260 | ! datname: name of the data set to use |
---|
2261 | ! values: mmatrix with the name of the values: variation (Nvariations, g.e. seasons), |
---|
2262 | ! table row (Nrows), table column (Ncols) |
---|
2263 | ! linedataset: line inside data-set where values start |
---|
2264 | ! namevalues: string assigned to each table row |
---|
2265 | ! ip, jp: coordinates of the grid point |
---|
2266 | ! Nsoiltypes: number of soil types |
---|
2267 | ! soiltypes[t/b]: percentage of the soil types at a given grid point on top/bottom |
---|
2268 | ! soiltypes: percentage of the soil types at a given grid point (lineray interpolated) |
---|
2269 | ! NOsoiltypes: vector with the names of soil types which are not land |
---|
2270 | |
---|
2271 | fname='"compute_TOTsoil_mass_water_values"' |
---|
2272 | errmsg='ERROR -- error -- ERROR -- error' |
---|
2273 | |
---|
2274 | IF (Nrows /= Nsoiltypes) THEN |
---|
2275 | PRINT *,errmsg |
---|
2276 | PRINT ' ' // TRIM(fname) // ': Number of soil-types: ',Nsoiltypes, & |
---|
2277 | ' does not concide with the number of values: ', Nrows, & |
---|
2278 | ' from the data-set "' // TRIM(datname) // '" !!!' |
---|
2279 | STOP |
---|
2280 | END IF |
---|
2281 | |
---|
2282 | IF (ALLOCATED(soiltypes)) DEALLOCATE(soiltypes) |
---|
2283 | ALLOCATE(soiltypes(Nrows)) |
---|
2284 | |
---|
2285 | ! On WRF3.3 SOILPARM.TBL |
---|
2286 | SELECT CASE (TRIM(datname)) |
---|
2287 | |
---|
2288 | CASE ('STAS') |
---|
2289 | NNOsoiltypes = 1 |
---|
2290 | IF (ALLOCATED(NOsoiltypes)) DEALLOCATE(NOsoiltypes) |
---|
2291 | ALLOCATE(NOsoiltypes(NNOsoiltypes)) |
---|
2292 | |
---|
2293 | NOsoiltypes=(/ 'WATER' /) |
---|
2294 | CASE ('STAS-RUC') |
---|
2295 | NNOsoiltypes = 1 |
---|
2296 | IF (ALLOCATED(NOsoiltypes)) DEALLOCATE(NOsoiltypes) |
---|
2297 | ALLOCATE(NOsoiltypes( NNOsoiltypes)) |
---|
2298 | |
---|
2299 | NOsoiltypes=(/ 'WATER' /) |
---|
2300 | END SELECT |
---|
2301 | |
---|
2302 | soilTOTmass = 0. |
---|
2303 | soilTOThumidity = 0. |
---|
2304 | |
---|
2305 | DO iz=1,Nsoillayers |
---|
2306 | IF (iz == 1 ) THEN |
---|
2307 | dz = soillayers(iz) |
---|
2308 | ELSE |
---|
2309 | dz = soillayers(iz) - soillayers(iz - 1) |
---|
2310 | END IF |
---|
2311 | IF (( iz > 1) .AND. (iz < Nsoillayers)) THEN |
---|
2312 | DO ix=1,Nrows |
---|
2313 | CALL interpolate_lin(soillayers(1), soiltypest(ix), & |
---|
2314 | soillayers(Nsoillayers), soiltypesb(ix), soillayers(iz), soiltypes(ix)) |
---|
2315 | END DO |
---|
2316 | ELSEIF (iz == 1) THEN |
---|
2317 | soiltypes = soiltypest |
---|
2318 | ELSE |
---|
2319 | soiltypes = soiltypesb |
---|
2320 | END IF |
---|
2321 | ! PRINT *,'lev: ',iz,' dz: ',dz,' soil types: ' ,soiltypes |
---|
2322 | |
---|
2323 | CALL compute_soil_mass(Nvariations,Nrows,Ncols,values,namevalues,soiltypes, & |
---|
2324 | NNOsoiltypes,NOsoiltypes,smois(iz),sizeX,sizeY,dz,soilmass,soilhumidity) |
---|
2325 | soilTOTmass = soilTOTmass + soilmass |
---|
2326 | soilTOThumidity = soilTOThumidity + soilhumidity |
---|
2327 | END DO |
---|
2328 | |
---|
2329 | IF ( (soilTOTmass < 0.) .OR. (soilTOThumidity < 0.) ) THEN |
---|
2330 | PRINT *,errmsg |
---|
2331 | PRINT *,' ' // TRIM(fname) // ' at point ',ip, ', ', jp, & |
---|
2332 | ' negative total soil mass: ',soilTOTmass, ' or soil water content', & |
---|
2333 | soilTOThumidity |
---|
2334 | PRINT *,' Soil_name density soil_types_________' |
---|
2335 | DO ix=1, Nrows |
---|
2336 | PRINT *,ix, TRIM(namevalues(ix)), values(1,ix,2), soiltypes(ix) |
---|
2337 | END DO |
---|
2338 | PRINT *,' depth_layer partial_mass soil+moisure partial_moisture__________' |
---|
2339 | DO iz=1,Nsoillayers |
---|
2340 | IF (iz == 1 ) THEN |
---|
2341 | dz = soillayers(iz) |
---|
2342 | ELSE |
---|
2343 | dz = soillayers(iz) - soillayers(iz - 1) |
---|
2344 | END IF |
---|
2345 | IF (( iz > 1) .AND. (iz < Nsoillayers)) THEN |
---|
2346 | DO ix=1,Nrows |
---|
2347 | CALL interpolate_lin(soillayers(1), soiltypest(ix), & |
---|
2348 | soillayers(Nsoillayers), soiltypesb(ix), soillayers(iz),soiltypes(ix)) |
---|
2349 | END DO |
---|
2350 | ELSEIF (iz == 1) THEN |
---|
2351 | soiltypes = soiltypest |
---|
2352 | ELSE |
---|
2353 | soiltypes = soiltypesb |
---|
2354 | END IF |
---|
2355 | PRINT *,' ',iz, dz, values(1,iz,2)*sizeX*sizeY*dz*1000., smois(iz), & |
---|
2356 | smois(iz)*values(1,iz,2)*sizeX*sizeY*dz*1000. |
---|
2357 | END DO |
---|
2358 | STOP |
---|
2359 | END IF |
---|
2360 | |
---|
2361 | |
---|
2362 | DEALLOCATE(soiltypes) |
---|
2363 | DEALLOCATE(NOsoiltypes) |
---|
2364 | |
---|
2365 | END SUBROUTINE compute_TOTsoil_mass_water_values |
---|
2366 | |
---|
2367 | ! cat Registry_out.LMDZ | awk '{print " pwrf_grid%"tolower($3)" "$4" = p"tolower($3)" "$4}' >> get_lmdz_out.f90 |
---|
2368 | ! cat get_lmdz_out.f90n' ',' '{print $3}' | tr '\n |
---|
2369 | SUBROUTINE get_lmdz_out(ddx, ddy, ddz, bdy, dlmdz, Nks, Nz0, nqtot, ivap, iliq, & |
---|
2370 | iflag_pbl, iflag_con, iflag_thermals, iflag_wake, read_climoz, itap, & |
---|
2371 | pdtphys, paprs, pplay, pphi, pphis, presnivs, qx, & |
---|
2372 | d_u, d_v, d_t, d_qx, d_ps, & |
---|
2373 | pwrf_grid) |
---|
2374 | ! Subroutine to get all the LMDZ output variables into the WRF Registry structure |
---|
2375 | |
---|
2376 | ! WRF modules |
---|
2377 | USE module_domain_type |
---|
2378 | |
---|
2379 | ! LMDZ modules |
---|
2380 | USE indice_sol_mod |
---|
2381 | USE phys_state_var_mod |
---|
2382 | USE phys_output_var_mod |
---|
2383 | ! USE pbl_surface_mod |
---|
2384 | USE phys_local_var_mod |
---|
2385 | USE comgeomphy |
---|
2386 | |
---|
2387 | ! WRF_LMDZ modules |
---|
2388 | USE output_lmdz_NOmodule |
---|
2389 | |
---|
2390 | ! netcdf |
---|
2391 | USE netcdf, ONLY: nf90_fill_real |
---|
2392 | |
---|
2393 | IMPLICIT NONE |
---|
2394 | |
---|
2395 | TYPE(domain) , TARGET :: pwrf_grid |
---|
2396 | INTEGER, INTENT(IN) :: ddx, ddy, ddz, bdy, & |
---|
2397 | dlmdz, Nks, Nz0, nqtot, ivap, iliq, iflag_pbl, iflag_con, iflag_thermals, & |
---|
2398 | iflag_wake, read_climoz, itap |
---|
2399 | REAL, INTENT(IN) :: pdtphys |
---|
2400 | REAL, DIMENSION(dlmdz), INTENT(IN) :: pphis, d_ps |
---|
2401 | REAL, DIMENSION(ddz), INTENT(IN) :: presnivs |
---|
2402 | REAL, DIMENSION(dlmdz,ddz+1), INTENT(IN) :: paprs |
---|
2403 | REAL, DIMENSION(dlmdz,ddz), INTENT(IN) :: pplay, pphi, d_u, d_v, d_t |
---|
2404 | REAL, DIMENSION(dlmdz,ddz,nqtot), INTENT(IN) :: qx, d_qx |
---|
2405 | |
---|
2406 | ! EXTERNAL suphel |
---|
2407 | |
---|
2408 | ! Local |
---|
2409 | INTEGER :: dx, dy, ix, iy, ixy, iz |
---|
2410 | INTEGER :: i,j,k,l,n |
---|
2411 | INTEGER :: nsrf |
---|
2412 | REAL, DIMENSION(dlmdz,ddz) :: var3d |
---|
2413 | REAL, DIMENSION(dlmdz,ddz+1) :: var3dstagg |
---|
2414 | REAL :: RDAY, RG, RD, RMO3 |
---|
2415 | CHARACTER(LEN=50) :: fname, errmsg |
---|
2416 | CHARACTER(LEN=4) :: bb2 |
---|
2417 | REAL, DIMENSION(dlmdz) :: zx_tmp_fi2d ! variable temporaire grille physique |
---|
2418 | REAL, DIMENSION(dlmdz,ddz) :: zx_tmp_fi3d ! variable temporaire pour champs 3D |
---|
2419 | REAL, DIMENSION(dlmdz,ddz+1) :: zx_tmp_fi3d1 !variable temporaire pour champs 3D (kelvp1) |
---|
2420 | REAL(KIND=8), DIMENSION(dlmdz,ddz) :: zx_tmp2_fi3d ! variable temporaire pour champs 3D |
---|
2421 | REAL, PARAMETER :: dobson_u = 2.1415e-05 |
---|
2422 | |
---|
2423 | ! klon=dldmz |
---|
2424 | ! klev=ddz |
---|
2425 | |
---|
2426 | |
---|
2427 | !! To think about / discuss with Frederic |
---|
2428 | |
---|
2429 | INCLUDE "declare_STDlev.h" |
---|
2430 | ! INCLUDE "calcul_STDlev.h" |
---|
2431 | |
---|
2432 | !!!!!!! Variables |
---|
2433 | ! pwrf_grid: WRF grid structure |
---|
2434 | |
---|
2435 | ! L. Fita, LMD January 2014. It should work using suphel, but not? |
---|
2436 | ! CALL suphel |
---|
2437 | |
---|
2438 | fname='"get_lmdz_out"' |
---|
2439 | errmsg='ERROR -- error -- ERROR -- error' |
---|
2440 | |
---|
2441 | RDAY = 86400. |
---|
2442 | RG = 9.81 |
---|
2443 | RD = 287. |
---|
2444 | RMO3=47.9942 |
---|
2445 | |
---|
2446 | PRINT *,' Lluis in ' //TRIM(fname)//'!' |
---|
2447 | |
---|
2448 | ! Remove variable assignation already done in get_lmdz? These variables are not used |
---|
2449 | ! in restart/boundaries |
---|
2450 | |
---|
2451 | DO ix=1,ddx |
---|
2452 | DO iy=1,ddy |
---|
2453 | ixy=ddx*(iy-1)+ix |
---|
2454 | |
---|
2455 | pwrf_grid%lages_ter(ix,iy) = agesno(ixy,1) |
---|
2456 | pwrf_grid%lages_lic(ix,iy) = agesno(ixy,2) |
---|
2457 | pwrf_grid%lages_oce(ix,iy) = agesno(ixy,3) |
---|
2458 | pwrf_grid%lages_sic(ix,iy) = agesno(ixy,4) |
---|
2459 | !! pwrf_grid%lahyb(misc??) = pahyb misc |
---|
2460 | pwrf_grid%laire(ix,iy) = airephy(ixy) |
---|
2461 | pwrf_grid%laireter(ix,iy) = paire_ter(ixy) |
---|
2462 | pwrf_grid%lalb1(ix,iy) = albsol1(ixy) |
---|
2463 | pwrf_grid%lalb2(ix,iy) = albsol2(ixy) |
---|
2464 | pwrf_grid%lalbe_ter(ix,iy) = falb1(ixy,1) |
---|
2465 | pwrf_grid%lalbe_lic(ix,iy) = falb1(ixy,2) |
---|
2466 | pwrf_grid%lalbe_oce(ix,iy) = falb1(ixy,3) |
---|
2467 | pwrf_grid%lalbe_sic(ix,iy) = falb1(ixy,4) |
---|
2468 | pwrf_grid%lale(ix,iy) = ale(ixy) |
---|
2469 | pwrf_grid%lale_bl(ix,iy) = ale_bl(ixy) |
---|
2470 | pwrf_grid%lale_wk(ix,iy) = ale_wake(ixy) |
---|
2471 | pwrf_grid%lalp(ix,iy) = alp(ixy) |
---|
2472 | pwrf_grid%lalp_bl(ix,iy) = alp_bl(ixy) |
---|
2473 | pwrf_grid%lalp_wk(ix,iy) = alp_wake(ixy) |
---|
2474 | !! pwrf_grid%lalt(misc??) = palt misc |
---|
2475 | !! pwrf_grid%lbhyb(misc??) = pbhyb misc |
---|
2476 | pwrf_grid%lbils(ix,iy) = bils(ixy) |
---|
2477 | pwrf_grid%lbils_diss(ix,iy) = bils_diss(ixy) |
---|
2478 | pwrf_grid%lbils_ec(ix,iy) = bils_ec(ixy) |
---|
2479 | pwrf_grid%lbils_enthalp(ix,iy) = bils_enthalp(ixy) |
---|
2480 | pwrf_grid%lbils_kinetic(ix,iy) = bils_kinetic(ixy) |
---|
2481 | pwrf_grid%lbils_latent(ix,iy) = bils_latent(ixy) |
---|
2482 | pwrf_grid%lbils_tke(ix,iy) = bils_tke(ixy) |
---|
2483 | pwrf_grid%lcape(ix,iy) = cape(ixy) |
---|
2484 | pwrf_grid%lcape_max(ix,iy) = cape(ixy) |
---|
2485 | pwrf_grid%lcdrh(ix,iy) = cdragh(ixy) |
---|
2486 | pwrf_grid%lcdrm(ix,iy) = cdragm(ixy) |
---|
2487 | pwrf_grid%lcldh(ix,iy) = cldh(ixy) |
---|
2488 | pwrf_grid%lcldl(ix,iy) = cldl(ixy) |
---|
2489 | pwrf_grid%lcldm(ix,iy) = cldm(ixy) |
---|
2490 | pwrf_grid%lcldq(ix,iy) = cldq(ixy) |
---|
2491 | pwrf_grid%lcldt(ix,iy) = cldt(ixy) |
---|
2492 | pwrf_grid%lcontfracatm(ix,iy) = pctsrf(ixy,is_ter)+pctsrf(ixy,is_lic) |
---|
2493 | pwrf_grid%lcontfracor(ix,iy) = pctsrf(ixy,is_ter) |
---|
2494 | !NOTknown! pwrf_grid%lcumpb(ix,iy) = cumpb(ixy) |
---|
2495 | !NOTknown! pwrf_grid%lcumrn(ix,iy) = cumrn(ixy) |
---|
2496 | pwrf_grid%ldthmin(ix,iy) = dthmin(ixy) |
---|
2497 | pwrf_grid%ldtsvdft(ix,iy) = d_ts(ixy,is_ter) |
---|
2498 | pwrf_grid%ldtsvdfo(ix,iy) = d_ts(ixy,is_oce) |
---|
2499 | pwrf_grid%ldtsvdfg(ix,iy) = d_ts(ixy,is_lic) |
---|
2500 | pwrf_grid%ldtsvdfi(ix,iy) = d_ts(ixy,is_sic) |
---|
2501 | pwrf_grid%levap(ix,iy) = fevap(ixy,1) |
---|
2502 | pwrf_grid%levap_ter(ix,iy) = fevap(ixy,1) |
---|
2503 | pwrf_grid%levap_lic(ix,iy) = fevap(ixy,2) |
---|
2504 | pwrf_grid%levap_oce(ix,iy) = fevap(ixy,3) |
---|
2505 | pwrf_grid%levap_sic(ix,iy) = fevap(ixy,4) |
---|
2506 | pwrf_grid%levappot_ter(ix,iy) = evap_pot(ixy,1) |
---|
2507 | pwrf_grid%levappot_lic(ix,iy) = evap_pot(ixy,2) |
---|
2508 | pwrf_grid%levappot_oce(ix,iy) = evap_pot(ixy,3) |
---|
2509 | pwrf_grid%levappot_sic(ix,iy) = evap_pot(ixy,4) |
---|
2510 | pwrf_grid%lfbase(ix,iy) = ema_cbmf(ixy) |
---|
2511 | pwrf_grid%lfder(ix,iy) = fder(ixy) |
---|
2512 | pwrf_grid%lffonte(ix,iy) = zxffonte(ixy) |
---|
2513 | pwrf_grid%lflat(ix,iy) = zxfluxlat(ixy) |
---|
2514 | pwrf_grid%lflw_ter(ix,iy) = fsollw(ixy,1) |
---|
2515 | pwrf_grid%lflw_lic(ix,iy) = fsollw(ixy,2) |
---|
2516 | pwrf_grid%lflw_oce(ix,iy) = fsollw(ixy,3) |
---|
2517 | pwrf_grid%lflw_sic(ix,iy) = fsollw(ixy,4) |
---|
2518 | pwrf_grid%lfqcalving(ix,iy) = zxfqcalving(ixy) |
---|
2519 | pwrf_grid%lfqfonte(ix,iy) = zxfqfonte(ixy) |
---|
2520 | pwrf_grid%lfract_ter(ix,iy) = pctsrf(ixy,1) |
---|
2521 | pwrf_grid%lfract_lic(ix,iy) = pctsrf(ixy,2) |
---|
2522 | pwrf_grid%lfract_oce(ix,iy) = pctsrf(ixy,3) |
---|
2523 | pwrf_grid%lfract_sic(ix,iy) = pctsrf(ixy,4) |
---|
2524 | pwrf_grid%lfsnow(ix,iy) = zfra_o(ixy) |
---|
2525 | pwrf_grid%lfsw_ter(ix,iy) = fsolw(ixy,1) |
---|
2526 | pwrf_grid%lfsw_lic(ix,iy) = fsolw(ixy,2) |
---|
2527 | pwrf_grid%lfsw_oce(ix,iy) = fsolw(ixy,3) |
---|
2528 | pwrf_grid%lfsw_sic(ix,iy) = fsolw(ixy,4) |
---|
2529 | !NOtnow! pwrf_grid%lftime_con(ix,iy) = float(itau_con)/float(itap) |
---|
2530 | pwrf_grid%lftime_th(ix,iy) = 0. |
---|
2531 | pwrf_grid%liwp(ix,iy) = fiwp(ixy) |
---|
2532 | !! pwrf_grid%llat j = pllat j |
---|
2533 | pwrf_grid%llat_ter(ix,iy) = fluxlat(ixy,1) |
---|
2534 | pwrf_grid%llat_lic(ix,iy) = fluxlat(ixy,2) |
---|
2535 | pwrf_grid%llat_oce(ix,iy) = fluxlat(ixy,3) |
---|
2536 | pwrf_grid%llat_sic(ix,iy) = fluxlat(ixy,4) |
---|
2537 | pwrf_grid%llmaxth(ix,iy) = lmax_th(ixy) |
---|
2538 | !! pwrf_grid%llon i = pllon i |
---|
2539 | pwrf_grid%llwdn200(ix,iy) = lwdn200(ixy) |
---|
2540 | pwrf_grid%llwdn200clr(ix,iy) = lwdn200clr(ixy) |
---|
2541 | pwrf_grid%llwdnsfc(ix,iy) = sollwdown(ixy) |
---|
2542 | pwrf_grid%llwdnsfcclr(ix,iy) = sollwdownclr(ixy) |
---|
2543 | pwrf_grid%llwdownor(ix,iy) = sollwdown(ixy) |
---|
2544 | pwrf_grid%llwp(ix,iy) = flwp(ixy) |
---|
2545 | pwrf_grid%llwup200(ix,iy) = lwup200(ixy) |
---|
2546 | pwrf_grid%llwup200clr(ix,iy) = lwup200clr(ixy) |
---|
2547 | pwrf_grid%llwupsfc(ix,iy) = sollwdown(ixy)-sollw(ixy) |
---|
2548 | pwrf_grid%llwupsfcclr(ix,iy) = -1.*lwdn0(ixy,1)-sollw0(ixy) |
---|
2549 | pwrf_grid%lmsnow(ix,iy) = snow_o(ixy) |
---|
2550 | pwrf_grid%lndayrain(ix,iy) = nday_rain(ixy) |
---|
2551 | pwrf_grid%lnettop(ix,iy) = topsw(ixy)-toplw(ixy) |
---|
2552 | pwrf_grid%lpbase(ix,iy) = ema_pcb(ixy) |
---|
2553 | pwrf_grid%lphis(ix,iy) = pphis(ixy) |
---|
2554 | pwrf_grid%lplcl(ix,iy) = plcl(ixy) |
---|
2555 | pwrf_grid%lplfc(ix,iy) = plfc(ixy) |
---|
2556 | pwrf_grid%lpluc(ix,iy) = rain_con(ixy) + snow_con(ixy) |
---|
2557 | pwrf_grid%lplul(ix,iy) = rain_lsc(ixy) + snow_lsc(ixy) |
---|
2558 | pwrf_grid%lplulst(ix,iy) = plul_st(ixy) |
---|
2559 | pwrf_grid%lplulth(ix,iy) = plul_th(ixy) |
---|
2560 | pwrf_grid%lpourc_ter(ix,iy) = pctsrf(ixy,1)*100. |
---|
2561 | pwrf_grid%lpourc_lic(ix,iy) = pctsrf(ixy,2)*100. |
---|
2562 | pwrf_grid%lpourc_oce(ix,iy) = pctsrf(ixy,3)*100. |
---|
2563 | pwrf_grid%lpourc_sic(ix,iy) = pctsrf(ixy,4)*100. |
---|
2564 | pwrf_grid%lprecip(ix,iy) = rain_fall(ixy) + snow_fall(ixy) |
---|
2565 | !! pwrf_grid%lpresnivs k = plpresnivs k |
---|
2566 | pwrf_grid%lprw(ix,iy) = prw(ixy) |
---|
2567 | pwrf_grid%lpsol(ix,iy) = paprs(ixy,1) |
---|
2568 | pwrf_grid%lptop(ix,iy) = ema_pct(ixy) |
---|
2569 | pwrf_grid%lq2m(ix,iy) = zq2m(ixy) |
---|
2570 | pwrf_grid%lqsat2m(ix,iy) = qsat2m(ixy) |
---|
2571 | pwrf_grid%lqsol(ix,iy) = qsol(ixy) |
---|
2572 | pwrf_grid%lqsurf(ix,iy) = zxqsurf(ixy) |
---|
2573 | pwrf_grid%lradsol(ix,iy) = radsol(ixy) |
---|
2574 | pwrf_grid%lrh2m(ix,iy) = MIN(100.,rh2m(ixy)*100.) |
---|
2575 | pwrf_grid%lrh2m_max(ix,iy) = MIN(100.,rh2m(ixy)*100.) |
---|
2576 | pwrf_grid%lrh2m_min(ix,iy) = MIN(100.,rh2m(ixy)*100.) |
---|
2577 | pwrf_grid%lrugs(ix,iy) = zxrugs(ixy) |
---|
2578 | pwrf_grid%lrugs_ter(ix,iy) = frugs(ixy,1) |
---|
2579 | pwrf_grid%lrugs_lic(ix,iy) = frugs(ixy,2) |
---|
2580 | pwrf_grid%lrugs_oce(ix,iy) = frugs(ixy,3) |
---|
2581 | pwrf_grid%lrugs_sic(ix,iy) = frugs(ixy,4) |
---|
2582 | pwrf_grid%lsens(ix,iy) = -1.*sens(ixy) |
---|
2583 | pwrf_grid%lsicf(ix,iy) = pctsrf(ixy,is_sic) |
---|
2584 | pwrf_grid%ls_lcl(ix,iy) = s_lcl(ixy) |
---|
2585 | pwrf_grid%lslp(ix,iy) = slp(ixy) |
---|
2586 | pwrf_grid%lsnow(ix,iy) = snow_fall(ixy) |
---|
2587 | pwrf_grid%lsnowl(ix,iy) = snow_lsc(ixy) |
---|
2588 | pwrf_grid%lsoll(ix,iy) = sollw(ixy) |
---|
2589 | pwrf_grid%lsoll0(ix,iy) = sollw0(ixy) |
---|
2590 | pwrf_grid%lsollwdown(ix,iy) = sollwdown(ixy) |
---|
2591 | pwrf_grid%lsols(ix,iy) = solsw(ixy) |
---|
2592 | pwrf_grid%lsols0(ix,iy) = solsw0(ixy) |
---|
2593 | pwrf_grid%ls_pblh(ix,iy) = s_pblh(ixy) |
---|
2594 | pwrf_grid%ls_pblt(ix,iy) = s_pblt(ixy) |
---|
2595 | pwrf_grid%ls_therm(ix,iy) = s_therm(ixy) |
---|
2596 | pwrf_grid%lswdn200(ix,iy) = swdn200(ixy) |
---|
2597 | pwrf_grid%lswdn200clr(ix,iy) = swdn200clr(ixy) |
---|
2598 | pwrf_grid%lswdnsfc(ix,iy) = swdn(ixy,1) |
---|
2599 | pwrf_grid%lswdnsfcclr(ix,iy) = swdn0(ixy,1) |
---|
2600 | pwrf_grid%lswdntoa(ix,iy) = swdn(ixy,klev+1) |
---|
2601 | pwrf_grid%lswdntoaclr(ix,iy) = swdn0(ixy,klev+1) |
---|
2602 | pwrf_grid%lswdownor(ix,iy) = solsw(ixy)/(1.-albsol1(ixy)) |
---|
2603 | pwrf_grid%lswnetor(ix,iy) = fsolsw(ixy,is_ter) |
---|
2604 | pwrf_grid%lswup200(ix,iy) = swup200(ixy) |
---|
2605 | pwrf_grid%lswup200clr(ix,iy) = swup200clr(ixy) |
---|
2606 | pwrf_grid%lswupsfc(ix,iy) = swup(ixy,1) |
---|
2607 | pwrf_grid%lswupsfcclr(ix,iy) = swup0(ixy,1) |
---|
2608 | pwrf_grid%lswuptoa(ix,iy) = swup(ixy,klev+1) |
---|
2609 | pwrf_grid%lswuptoaclr(ix,iy) = swup0(ixy,klev+1) |
---|
2610 | pwrf_grid%lt2m(ix,iy) = zt2m(ixy) |
---|
2611 | pwrf_grid%lt2m_max(ix,iy) = zt2m(ixy) |
---|
2612 | pwrf_grid%lt2m_min(ix,iy) = zt2m(ixy) |
---|
2613 | pwrf_grid%lt2m_ter(ix,iy) = t2m(ixy,1) |
---|
2614 | pwrf_grid%lt2m_lic(ix,iy) = t2m(ixy,2) |
---|
2615 | pwrf_grid%lt2m_oce(ix,iy) = t2m(ixy,3) |
---|
2616 | pwrf_grid%lt2m_sic(ix,iy) = t2m(ixy,4) |
---|
2617 | pwrf_grid%ltaux(ix,iy) = 0. |
---|
2618 | pwrf_grid%ltauy(ix,iy) = 0. |
---|
2619 | DO iz=1,4 |
---|
2620 | pwrf_grid%ltaux(ix,iy) = pwrf_grid%ltauy(ix,iy) + pctsrf(ixy,iz)* & |
---|
2621 | fluxu(ixy,1,iz) |
---|
2622 | pwrf_grid%ltauy(ix,iy) = pwrf_grid%ltauy(ix,iy) + pctsrf(ixy,iz)* & |
---|
2623 | fluxv(ixy,1,iz) |
---|
2624 | END DO |
---|
2625 | pwrf_grid%ltaux_ter(ix,iy) = fluxu(ixy,1,1) |
---|
2626 | pwrf_grid%ltaux_lic(ix,iy) = fluxu(ixy,1,2) |
---|
2627 | pwrf_grid%ltaux_oce(ix,iy) = fluxu(ixy,1,3) |
---|
2628 | pwrf_grid%ltaux_sic(ix,iy) = fluxu(ixy,1,4) |
---|
2629 | pwrf_grid%ltauy_ter(ix,iy) = fluxv(ixy,1,1) |
---|
2630 | pwrf_grid%ltauy_lic(ix,iy) = fluxv(ixy,1,2) |
---|
2631 | pwrf_grid%ltauy_oce(ix,iy) = fluxv(ixy,1,3) |
---|
2632 | pwrf_grid%ltauy_sic(ix,iy) = fluxv(ixy,1,4) |
---|
2633 | !! pwrf_grid%ltime - = pltime - |
---|
2634 | !! pwrf_grid%ltime_bounds {lbnds} = pltime_bounds {lbnds} |
---|
2635 | !! pwrf_grid%lti86400 - = plti86400 - |
---|
2636 | IF ( (pctsrf(ixy,is_oce).GT.epsfra).OR.(pctsrf(ixy,is_sic).GT.epsfra) ) THEN |
---|
2637 | pwrf_grid%lt_oce_sic(ix,iy) = (ftsol(ixy, is_oce) * pctsrf(ixy,is_oce)+ & |
---|
2638 | ftsol(ixy, is_sic)*pctsrf(ixy,is_sic))/(pctsrf(ixy,is_oce)+ & |
---|
2639 | pctsrf(ixy,is_sic)) |
---|
2640 | ELSE |
---|
2641 | pwrf_grid%lt_oce_sic(ix,iy) = 273.15 |
---|
2642 | END IF |
---|
2643 | !! pwrf_grid%lt_op_00001800(misc??) = pt_op_00001800 misc |
---|
2644 | !! pwrf_grid%lt_op_00001800_bnds(misc??) = pt_op_00001800_bnds misc |
---|
2645 | pwrf_grid%ltopl(ix,iy) = toplw(ixy) |
---|
2646 | pwrf_grid%ltopl0(ix,iy) = toplw0(ixy) |
---|
2647 | pwrf_grid%ltops(ix,iy) = topsw(ixy) |
---|
2648 | pwrf_grid%ltops0(ix,iy) = topsw0(ixy) |
---|
2649 | pwrf_grid%ltpot(ix,iy) = tpot(ixy) |
---|
2650 | pwrf_grid%ltpote(ix,iy) = tpote(ixy) |
---|
2651 | pwrf_grid%ltsol(ix,iy) = zxtsol(ixy) |
---|
2652 | pwrf_grid%ltsol_ter(ix,iy) = ftsol(ixy,1) |
---|
2653 | pwrf_grid%ltsol_lic(ix,iy) = ftsol(ixy,2) |
---|
2654 | pwrf_grid%ltsol_oce(ix,iy) = ftsol(ixy,3) |
---|
2655 | pwrf_grid%ltsol_sic(ix,iy) = ftsol(ixy,4) |
---|
2656 | pwrf_grid%lu10m(ix,iy) = zu10m(ixy) |
---|
2657 | pwrf_grid%lu10m_ter(ix,iy) = u10m(ixy,1) |
---|
2658 | pwrf_grid%lu10m_lic(ix,iy) = u10m(ixy,2) |
---|
2659 | pwrf_grid%lu10m_oce(ix,iy) = u10m(ixy,3) |
---|
2660 | pwrf_grid%lu10m_sic(ix,iy) = u10m(ixy,4) |
---|
2661 | pwrf_grid%lue(ix,iy) = ue(ixy) |
---|
2662 | pwrf_grid%luq(ix,iy) = uq(ixy) |
---|
2663 | pwrf_grid%lustar(ix,iy) = zustar(ixy) |
---|
2664 | pwrf_grid%lustar_ter(ix,iy) = ustar(ixy,1) |
---|
2665 | pwrf_grid%lustar_lic(ix,iy) = ustar(ixy,2) |
---|
2666 | pwrf_grid%lustar_oce(ix,iy) = ustar(ixy,3) |
---|
2667 | pwrf_grid%lustar_sic(ix,iy) = ustar(ixy,4) |
---|
2668 | pwrf_grid%lv10m(ix,iy) = zv10m(ixy) |
---|
2669 | pwrf_grid%lv10m_ter(ix,iy) = v10m(ixy,1) |
---|
2670 | pwrf_grid%lv10m_lic(ix,iy) = v10m(ixy,2) |
---|
2671 | pwrf_grid%lv10m_oce(ix,iy) = v10m(ixy,3) |
---|
2672 | pwrf_grid%lv10m_sic(ix,iy) = v10m(ixy,4) |
---|
2673 | pwrf_grid%lve(ix,iy) = ve(ixy) |
---|
2674 | pwrf_grid%lvq(ix,iy) = vq(ixy) |
---|
2675 | pwrf_grid%lwake_h(ix,iy) = wake_h(ixy) |
---|
2676 | pwrf_grid%lwake_s(ix,iy) = wake_s(ixy) |
---|
2677 | pwrf_grid%lwape(ix,iy) = wake_pe(ixy) |
---|
2678 | pwrf_grid%lwbeff(ix,iy) = wbeff(ixy) |
---|
2679 | pwrf_grid%lwbilo_ter(ix,iy) = wfbilo(ixy,1) |
---|
2680 | pwrf_grid%lwbilo_lic(ix,iy) = wfbilo(ixy,2) |
---|
2681 | pwrf_grid%lwbilo_oce(ix,iy) = wfbilo(ixy,3) |
---|
2682 | pwrf_grid%lwbilo_sic(ix,iy) = wfbilo(ixy,4) |
---|
2683 | pwrf_grid%lwbils_ter(ix,iy) = wfbils(ixy,1) |
---|
2684 | pwrf_grid%lwbils_lic(ix,iy) = wfbils(ixy,2) |
---|
2685 | pwrf_grid%lwbils_oce(ix,iy) = wfbils(ixy,3) |
---|
2686 | pwrf_grid%lwbils_sic(ix,iy) = wfbils(ixy,4) |
---|
2687 | pwrf_grid%lweakinv(ix,iy) = weak_inversion(ixy) |
---|
2688 | pwrf_grid%lwind10m(ix,iy) = SQRT(zu10m(ixy)*zu10m(ixy)+zv10m(ixy)*zv10m(ixy)) |
---|
2689 | pwrf_grid%lwind10max(ix,iy) = SQRT(zu10m(ixy)*zu10m(ixy)+zv10m(ixy)*zv10m(ixy)) |
---|
2690 | pwrf_grid%lzmax_th(ix,iy) = zmax_th(ixy) |
---|
2691 | ! Output at the standard levels |
---|
2692 | !! |
---|
2693 | ! DO iz=1, nlevSTD |
---|
2694 | ! bb2=clevSTD(iz) |
---|
2695 | |
---|
2696 | ! IF (TRIM(bb2) == "850") THEN |
---|
2697 | ! pwrf_grid%lq850(ix,iy) = q850(ixy) |
---|
2698 | ! pwrf_grid%lt850(ix,iy) = t850(ixy) |
---|
2699 | ! pwrf_grid%lu850(ix,iy) = u850(ixy) |
---|
2700 | ! pwrf_grid%lv850(ix,iy) = v850(ixy) |
---|
2701 | ! pwrf_grid%lw850(ix,iy) = w850(ixy) |
---|
2702 | ! pwrf_grid%lz850(ix,iy) = z850(ixy) |
---|
2703 | ! ELSE IF (TRIM(bb2) == "700") THEN |
---|
2704 | ! pwrf_grid%lq700(ix,iy) = q700(ixy) |
---|
2705 | ! pwrf_grid%lt700(ix,iy) = t700(ixy) |
---|
2706 | ! pwrf_grid%lu700(ix,iy) = u700(ixy) |
---|
2707 | ! pwrf_grid%lv700(ix,iy) = v700(ixy) |
---|
2708 | ! pwrf_grid%lw700(ix,iy) = w700(ixy) |
---|
2709 | ! pwrf_grid%lz700(ix,iy) = z700(ixy) |
---|
2710 | ! ELSE IF (TRIM(bb2) == "500") THEN |
---|
2711 | ! pwrf_grid%lq500(ix,iy) = q500(ixy) |
---|
2712 | ! pwrf_grid%lt500(ix,iy) = t500(ixy) |
---|
2713 | ! pwrf_grid%lu500(ix,iy) = u500(ixy) |
---|
2714 | ! pwrf_grid%lv500(ix,iy) = v500(ixy) |
---|
2715 | ! pwrf_grid%lw500(ix,iy) = w500(ixy) |
---|
2716 | ! pwrf_grid%lz500(ix,iy) = z500(ixy) |
---|
2717 | ! ELSE IF (TRIM(bb2) == "200")THEN |
---|
2718 | ! pwrf_grid%lq200(ix,iy) = q200(ixy) |
---|
2719 | ! pwrf_grid%lt200(ix,iy) = t200(ixy) |
---|
2720 | ! pwrf_grid%lu200(ix,iy) = u200(ixy) |
---|
2721 | ! pwrf_grid%lv200(ix,iy) = v200(ixy) |
---|
2722 | ! pwrf_grid%lw200(ix,iy) = w200(ixy) |
---|
2723 | ! pwrf_grid%lz200(ix,iy) = z200(ixy) |
---|
2724 | ! ELSE IF (TRIM(bb2) == "100") THEN |
---|
2725 | ! pwrf_grid%lq100(ix,iy) = q100(ixy) |
---|
2726 | ! pwrf_grid%lt100(ix,iy) = t100(ixy) |
---|
2727 | ! pwrf_grid%lu100(ix,iy) = u100(ixy) |
---|
2728 | ! pwrf_grid%lv100(ix,iy) = v100(ixy) |
---|
2729 | ! pwrf_grid%lw100(ix,iy) = w100(ixy) |
---|
2730 | ! pwrf_grid%lz100(ix,iy) = z100(ixy) |
---|
2731 | ! ELSE IF (TRIM(bb2) == "50") THEN |
---|
2732 | ! pwrf_grid%lq50(ix,iy) = q50(ixy) |
---|
2733 | ! pwrf_grid%lt50(ix,iy) = t50(ixy) |
---|
2734 | ! pwrf_grid%lu50(ix,iy) = u50(ixy) |
---|
2735 | ! pwrf_grid%lv50(ix,iy) = v50(ixy) |
---|
2736 | ! pwrf_grid%lw50(ix,iy) = w50(ixy) |
---|
2737 | ! pwrf_grid%lz50(ix,iy) = z50(ixy) |
---|
2738 | ! ELSE IF (TRIM(bb2) == "10") THEN |
---|
2739 | ! pwrf_grid%lq10(ix,iy) = q10(ixy) |
---|
2740 | ! pwrf_grid%lt10(ix,iy) = t10(ixy) |
---|
2741 | ! pwrf_grid%lu10(ix,iy) = u10(ixy) |
---|
2742 | ! pwrf_grid%lv10(ix,iy) = v10(ixy) |
---|
2743 | ! pwrf_grid%lw10(ix,iy) = w10(ixy) |
---|
2744 | ! pwrf_grid%lz10(ix,iy) = z10(ixy) |
---|
2745 | ! END IF |
---|
2746 | ! END DO |
---|
2747 | |
---|
2748 | END DO |
---|
2749 | END DO |
---|
2750 | |
---|
2751 | dx=ddx+2*bdy |
---|
2752 | dy=ddy+2*bdy |
---|
2753 | |
---|
2754 | CALL vect_mat_zstagg(fraca, dx, dy, ddz+1, bdy, pwrf_grid%la_th) |
---|
2755 | CALL vect_mat(beta_prec, dx, dy, ddz, bdy, pwrf_grid%lbeta_prec) |
---|
2756 | var3d=upwd + dnwd+ dnwd0 |
---|
2757 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldmc) |
---|
2758 | CALL vect_mat(dnwd, dx, dy, ddz, bdy, pwrf_grid%ldnwd) |
---|
2759 | CALL vect_mat(dnwd0, dx, dy, ddz, bdy, pwrf_grid%ldnwd0) |
---|
2760 | ! Originally is var3d = d_q_ajsb(1:ddx*ddy,1:ddz)/pdtphys |
---|
2761 | var3d = d_q_ajsb/pdtphys |
---|
2762 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldqajs) |
---|
2763 | var3d = d_q_con/pdtphys |
---|
2764 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldqcon) |
---|
2765 | CALL vect_mat(d_q_dyn, dx, dy, ddz, bdy, pwrf_grid%ldqdyn) |
---|
2766 | var3d = d_q_eva/pdtphys |
---|
2767 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldqeva) |
---|
2768 | var3d = d_q_lsc/pdtphys |
---|
2769 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldqlsc) |
---|
2770 | var3d = d_q_lscst/pdtphys |
---|
2771 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldqlscst) |
---|
2772 | var3d = d_q_lscth/pdtphys |
---|
2773 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldqlscth) |
---|
2774 | ! L. Fita, LMD January 2014. d_qx is the tendency of moisture due to physiq! |
---|
2775 | ! MOISTtend in the WRF_LMDZ coupling |
---|
2776 | CALL vect_mat(d_qx(:,:,ivap), dx, dy, ddz, bdy, pwrf_grid%ldqphy) |
---|
2777 | var3d = d_q_ajs/pdtphys - d_q_ajsb/pdtphys |
---|
2778 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldqthe) |
---|
2779 | var3d = d_q_vdf/pdtphys |
---|
2780 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldqvdf) |
---|
2781 | var3d = d_q_wake/pdtphys |
---|
2782 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldqwak) |
---|
2783 | var3d = d_t_ajsb/pdtphys |
---|
2784 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldtajs) |
---|
2785 | var3d = d_t_con/pdtphys |
---|
2786 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldtcon) |
---|
2787 | var3d = d_t_diss/pdtphys |
---|
2788 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldtdis) |
---|
2789 | CALL vect_mat(d_t_dyn, dx, dy, ddz, bdy, pwrf_grid%ldtdyn) |
---|
2790 | var3d = d_t_ec/pdtphys |
---|
2791 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldtec) |
---|
2792 | var3d = d_t_eva/pdtphys |
---|
2793 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldteva) |
---|
2794 | CALL vect_mat(detr_therm, dx, dy, ddz, bdy, pwrf_grid%ld_th) |
---|
2795 | CALL vect_mat(cldemi, dx, dy, ddz, bdy, pwrf_grid%lcldemi) |
---|
2796 | CALL vect_mat(cldtau, dx, dy, ddz, bdy, pwrf_grid%lcldtau) |
---|
2797 | CALL vect_mat(clwcon, dx, dy, ddz, bdy, pwrf_grid%lclwcon) |
---|
2798 | var3d = d_t_lif/pdtphys |
---|
2799 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldtlif) |
---|
2800 | var3d = d_t_lsc/pdtphys |
---|
2801 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldtlsc) |
---|
2802 | var3d = (d_t_lsc+d_t_eva)/pdtphys |
---|
2803 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldtlschr) |
---|
2804 | var3d = d_t_lscst/pdtphys |
---|
2805 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldtlscst) |
---|
2806 | var3d = d_t_lscth/pdtphys |
---|
2807 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldtlscth) |
---|
2808 | var3d=-1.*cool0/RDAY |
---|
2809 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldtlw0) |
---|
2810 | var3d = -1.*cool/RDAY |
---|
2811 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldtlwr) |
---|
2812 | var3d=d_t_oro/pdtphys |
---|
2813 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldtoro) |
---|
2814 | CALL vect_mat(d_t, dx, dy, ddz, bdy, pwrf_grid%ldtphy) |
---|
2815 | var3d=heat0/RDAY |
---|
2816 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldtsw0) |
---|
2817 | var3d=heat/RDAY |
---|
2818 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldtswr) |
---|
2819 | var3d = d_t_ajs/pdtphys - d_t_ajsb/pdtphys |
---|
2820 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldtthe) |
---|
2821 | var3d=d_t_vdf/pdtphys |
---|
2822 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldtvdf) |
---|
2823 | var3d=d_t_wake/pdtphys |
---|
2824 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldtwak) |
---|
2825 | var3d=d_u_con/pdtphys |
---|
2826 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%lducon) |
---|
2827 | CALL vect_mat(d_u_dyn, dx, dy, ddz, bdy, pwrf_grid%ldudyn) |
---|
2828 | var3d=d_u_lif/pdtphys |
---|
2829 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldulif) |
---|
2830 | var3d=d_u_oro/pdtphys |
---|
2831 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%lduoro) |
---|
2832 | var3d=d_u_vdf/pdtphys |
---|
2833 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%lduvdf) |
---|
2834 | var3d=d_v_con/pdtphys |
---|
2835 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldvcon) |
---|
2836 | CALL vect_mat(d_v_dyn, dx, dy, ddz, bdy, pwrf_grid%ldvdyn) |
---|
2837 | var3d=d_v_lif/pdtphys |
---|
2838 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldvlif) |
---|
2839 | var3d=d_v_oro/pdtphys |
---|
2840 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldvoro) |
---|
2841 | var3d=d_v_vdf/pdtphys |
---|
2842 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ldvvdf) |
---|
2843 | !NOTfound! CALL vect_mat(ec550aer, dx, dy, ddz, bdy, pwrf_grid%lec550aer) |
---|
2844 | CALL vect_mat(entr_therm, dx, dy, ddz, bdy, pwrf_grid%le_th) |
---|
2845 | CALL vect_mat(coefm(:,:,is_ave), dx, dy, ddz+1, bdy, pwrf_grid%levu) |
---|
2846 | CALL vect_mat(fl, dx, dy, ddz, bdy, pwrf_grid%lfl) |
---|
2847 | CALL vect_mat(zphi, dx, dy, ddz, bdy, pwrf_grid%lgeop) |
---|
2848 | var3d = q_seri + ql_seri |
---|
2849 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%lh2o) |
---|
2850 | CALL vect_mat(fiwc, dx, dy, ddz, bdy, pwrf_grid%liwcon) |
---|
2851 | CALL vect_mat(coefh(:,:,is_ave), dx, dy, ddz+1, bdy, pwrf_grid%lkz) |
---|
2852 | CALL vect_mat(coefh(:,:,is_ave), dx, dy, ddz+1, bdy, pwrf_grid%lkz_max) |
---|
2853 | ! CALL vect_mat(lambda_th, dx, dy, ddz, bdy, pwrf_grid%llambda_th) |
---|
2854 | CALL vect_mat(flwc, dx, dy, ddz, bdy, pwrf_grid%llwcon) |
---|
2855 | !NOTknown! CALL vect_mat(ma, dx, dy, ddz, bdy, pwrf_grid%lma) |
---|
2856 | CALL vect_mat(zmasse, dx, dy, ddz, bdy, pwrf_grid%lmass) |
---|
2857 | !NOTknown! CALL vect_mat(mc, dx, dy, ddz, bdy, pwrf_grid%lmc) |
---|
2858 | !NOTknown! CALL vect_mat(mcd, dx, dy, ddz, bdy, pwrf_grid%lmcd) |
---|
2859 | CALL vect_mat(ql_seri, dx, dy, ddz, bdy, pwrf_grid%loliq) |
---|
2860 | CALL vect_mat(q_seri, dx, dy, ddz, bdy, pwrf_grid%lovap) |
---|
2861 | ! L. Fita, LMD January 2014. d_qx is the tendency of moisture due to physiq! |
---|
2862 | ! MOISTtend in the WRF_LMDZ coupling |
---|
2863 | CALL vect_mat(qx(:,:,ivap), dx, dy, ddz, bdy, pwrf_grid%lovapinit) |
---|
2864 | !Notnow! var3d = wo(:,:,1)*dobson_u*1e3 / zmasse / rmo3 * rmd) |
---|
2865 | !Notnow! CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%lozone) |
---|
2866 | CALL vect_mat_zstagg(paprs, dx, dy, ddz+1, bdy, pwrf_grid%lpaprs) |
---|
2867 | !NOTknown! CALL vect_mat(pb, dx, dy, ddz, bdy, pwrf_grid%lpb) |
---|
2868 | CALL vect_mat_zstagg(pmflxs, dx, dy, ddz+1, bdy, pwrf_grid%lpr_con_i) |
---|
2869 | CALL vect_mat_zstagg(pmflxr, dx, dy, ddz+1, bdy, pwrf_grid%lpr_con_l) |
---|
2870 | CALL vect_mat(pplay, dx, dy, ddz, bdy, pwrf_grid%lpres) |
---|
2871 | CALL vect_mat_zstagg(psfl, dx, dy, ddz+1, bdy, pwrf_grid%lpr_lsc_i) |
---|
2872 | CALL vect_mat_zstagg(prfl, dx, dy, ddz+1, bdy, pwrf_grid%lpr_lsc_l) |
---|
2873 | var3d = 0. |
---|
2874 | WHERE (ptconv) var3d = 1. |
---|
2875 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%lptconv) |
---|
2876 | CALL vect_mat(zqasc, dx, dy, ddz, bdy, pwrf_grid%lq_th) |
---|
2877 | CALL vect_mat(ratqs, dx, dy, ddz, bdy, pwrf_grid%lratqs) |
---|
2878 | CALL vect_mat(re, dx, dy, ddz, bdy, pwrf_grid%lre) |
---|
2879 | CALL vect_mat(ref_ice, dx, dy, ddz, bdy, pwrf_grid%lref_ice) |
---|
2880 | CALL vect_mat(ref_liq, dx, dy, ddz, bdy, pwrf_grid%lref_liq) |
---|
2881 | CALL vect_mat(zx_rh, dx, dy, ddz, bdy, pwrf_grid%lrhum) |
---|
2882 | CALL vect_mat(lwdn, dx, dy, ddz, bdy, pwrf_grid%lrld) |
---|
2883 | CALL vect_mat(lwdn0, dx, dy, ddz, bdy, pwrf_grid%lrldcs) |
---|
2884 | CALL vect_mat(lwup, dx, dy, ddz, bdy, pwrf_grid%lrlu) |
---|
2885 | CALL vect_mat(lwup0, dx, dy, ddz, bdy, pwrf_grid%lrlucs) |
---|
2886 | !NOTknown! CALL vect_mat(rn, dx, dy, ddz, bdy, pwrf_grid%lrn) |
---|
2887 | CALL vect_mat(cldfra, dx, dy, ddz, bdy, pwrf_grid%lrneb) |
---|
2888 | CALL vect_mat(rnebcon, dx, dy, ddz, bdy, pwrf_grid%lrnebcon) |
---|
2889 | CALL vect_mat(rneb, dx, dy, ddz, bdy, pwrf_grid%lrnebls) |
---|
2890 | CALL vect_mat(swdn, dx, dy, ddz, bdy, pwrf_grid%lrsd) |
---|
2891 | CALL vect_mat(swdn0, dx, dy, ddz, bdy, pwrf_grid%lrsdcs) |
---|
2892 | CALL vect_mat(swup, dx, dy, ddz, bdy, pwrf_grid%lrsu) |
---|
2893 | CALL vect_mat(swup0, dx, dy, ddz, bdy, pwrf_grid%lrsucs) |
---|
2894 | CALL vect_mat(fluxt(:,:,1), dx, dy, ddz, bdy, pwrf_grid%lsens_ter) |
---|
2895 | CALL vect_mat(fluxt(:,:,2), dx, dy, ddz, bdy, pwrf_grid%lsens_lic) |
---|
2896 | CALL vect_mat(fluxt(:,:,3), dx, dy, ddz, bdy, pwrf_grid%lsens_oce) |
---|
2897 | CALL vect_mat(fluxt(:,:,4), dx, dy, ddz, bdy, pwrf_grid%lsens_sic) |
---|
2898 | CALL vect_mat(t_seri, dx, dy, ddz, bdy, pwrf_grid%ltemp) |
---|
2899 | CALL vect_mat(theta, dx, dy, ddz, bdy, pwrf_grid%ltheta) |
---|
2900 | var3d=0. |
---|
2901 | DO nsrf=1,nbsrf |
---|
2902 | DO iz=1,ddz+1 |
---|
2903 | var3dstagg(:,iz)=var3dstagg(:,iz)+pctsrf(:,nsrf)*pbl_tke(:,iz,nsrf) |
---|
2904 | ENDDO |
---|
2905 | ENDDO |
---|
2906 | CALL vect_mat_zstagg(var3dstagg, dx, dy, ddz+1, bdy, pwrf_grid%ltke) |
---|
2907 | CALL vect_mat_zstagg(var3dstagg, dx, dy, ddz+1, bdy, pwrf_grid%ltke_max) |
---|
2908 | CALL vect_mat(pbl_tke(:,:,1), dx, dy, ddz+1, bdy, pwrf_grid%ltke_ter) |
---|
2909 | CALL vect_mat(pbl_tke(:,:,2), dx, dy, ddz+1, bdy, pwrf_grid%ltke_lic) |
---|
2910 | CALL vect_mat(pbl_tke(:,:,3), dx, dy, ddz+1, bdy, pwrf_grid%ltke_oce) |
---|
2911 | CALL vect_mat(pbl_tke(:,:,4), dx, dy, ddz+1, bdy, pwrf_grid%ltke_sic) |
---|
2912 | !NOTinRegisty! CALL vect_mat(pbl_tke(:,:,1), dx, dy, ddz+1, bdy, pwrf_grid%ltke_ter_max) |
---|
2913 | !NOTinRegisty! CALL vect_mat(pbl_tke(:,:,2), dx, dy, ddz+1, bdy, pwrf_grid%ltke_lic_max) |
---|
2914 | !NOTinRegisty! CALL vect_mat(pbl_tke(:,:,3), dx, dy, ddz+1, bdy, pwrf_grid%ltke_oce_max) |
---|
2915 | !NOTinRegisty! CALL vect_mat(pbl_tke(:,:,4), dx, dy, ddz+1, bdy, pwrf_grid%ltke_sic_max) |
---|
2916 | var3d = d_qx(:,:,ivap)+d_q_dyn |
---|
2917 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ltnhus) |
---|
2918 | IF (iflag_thermals == 1) THEN |
---|
2919 | var3d = d_q_con/pdtphys |
---|
2920 | ELSE IF ((iflag_thermals > 1) .AND. (iflag_wake == 1)) THEN |
---|
2921 | var3d = d_q_con/pdtphys + d_q_ajs/pdtphys + d_q_wake/pdtphys |
---|
2922 | END IF |
---|
2923 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ltnhusc) |
---|
2924 | var3d = d_q_lsc/pdtphys + d_q_eva/pdtphys |
---|
2925 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ltnhusscpbl) |
---|
2926 | var3d = d_t + d_t_dyn |
---|
2927 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ltnt) |
---|
2928 | IF (iflag_thermals == 1) THEN |
---|
2929 | var3d = d_t_con/pdtphys + d_t_ajsb/pdtphys |
---|
2930 | ELSE IF ((iflag_thermals > 1) .AND. (iflag_wake == 1)) THEN |
---|
2931 | var3d=d_t_con/pdtphys + d_t_ajs/pdtphys + d_t_wake/pdtphys |
---|
2932 | END IF |
---|
2933 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ltntc) |
---|
2934 | var3d = heat/RDAY - cool/RDAY |
---|
2935 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ltntr) |
---|
2936 | var3d = (d_t_lsc + d_t_eva + d_t_vdf)/pdtphys |
---|
2937 | CALL vect_mat(var3d, dx, dy, ddz, bdy, pwrf_grid%ltntscpbl) |
---|
2938 | CALL vect_mat(upwd, dx, dy, ddz, bdy, pwrf_grid%lupwd) |
---|
2939 | CALL vect_mat(u_seri, dx, dy, ddz, bdy, pwrf_grid%lvitu) |
---|
2940 | CALL vect_mat(v_seri, dx, dy, ddz, bdy, pwrf_grid%lvitv) |
---|
2941 | CALL vect_mat(omega, dx, dy, ddz, bdy, pwrf_grid%lvitw) |
---|
2942 | CALL vect_mat_zstagg(Vprecip, dx, dy, ddz+1, bdy, pwrf_grid%lvprecip) |
---|
2943 | CALL vect_mat(wake_deltaq, dx, dy, ddz, bdy, pwrf_grid%lwake_deltaq) |
---|
2944 | CALL vect_mat(wake_deltat, dx, dy, ddz, bdy, pwrf_grid%lwake_deltat) |
---|
2945 | CALL vect_mat(wake_omg, dx, dy, ddz, bdy, pwrf_grid%lwake_omg) |
---|
2946 | CALL vect_mat(wdtrainA, dx, dy, ddz, bdy, pwrf_grid%lwdtraina) |
---|
2947 | CALL vect_mat(wdtrainM, dx, dy, ddz, bdy, pwrf_grid%lwdtrainm) |
---|
2948 | ! L. Fita, LMD, January 2014. pphis is the given value of the pressure at the surface |
---|
2949 | var3dstagg(:,1) = pphis(:)/RG |
---|
2950 | DO iz=1,ddz |
---|
2951 | var3dstagg(:,iz+1)= var3dstagg(:,iz)-(t_seri(:,iz)*RD*(paprs(:,iz+1)- & |
---|
2952 | paprs(:,iz)) )/ (pplay(:,iz)*RG) |
---|
2953 | ENDDO |
---|
2954 | CALL vect_mat(var3dstagg, dx, dy, ddz, bdy, pwrf_grid%lzfull) |
---|
2955 | var3dstagg(:,1) = pphis(:)/RG-( (t_seri(:,1)+zxtsol)/2.*RD*(pplay(:,1)- & |
---|
2956 | paprs(:,1)) )/ ((paprs(:,1)+pplay(:,1))/2.*RG) |
---|
2957 | DO iz=1, ddz-1 |
---|
2958 | var3dstagg(:,iz+1)= var3dstagg(:,iz)-( (t_seri(:,iz)+t_seri(:,iz+1))/ & |
---|
2959 | 2.*RD*(pplay(:,iz+1)-pplay(:,iz)))/( paprs(:,iz)*RG) |
---|
2960 | ENDDO |
---|
2961 | CALL vect_mat_zstagg(var3dstagg, dx, dy, ddz+1, bdy, pwrf_grid%lzhalf) |
---|
2962 | |
---|
2963 | RETURN |
---|
2964 | |
---|
2965 | END SUBROUTINE get_lmdz_out |
---|
2966 | |
---|
2967 | SUBROUTINE put_lmdz_in_WRFout(is, ie, js, je, ks, ke, ddx, ddy, ddz, bdy, dlmdz,xp,& |
---|
2968 | yp, lmdp, Nks, Nz0, nqtot, ivap, iliq, iflag_pbl, iflag_con, iflag_thermals, & |
---|
2969 | iflag_wake, itap, pdtphys, paprs, pplay, pphi, pphis, presnivs, qx, d_u, d_v, & |
---|
2970 | d_t, d_qx, d_ps, orchidlev, tsoilorchid, nslayers, wslay, wzs, wdzs, wnest_pos, & |
---|
2971 | wq2, wt2, wth2, wpsfc, wu10, wv10, wresm, wzetatop, wtslb, wsmois, & |
---|
2972 | wsh2o, wsmcrel, wsfroff, wudroff, wgrdflx, wacgrdflx, wsnow, wsnowh, & |
---|
2973 | wcanwat, wsst, wsstsk, wlai, wh_diabatic, wsina, wtsk, wtiso, & |
---|
2974 | wmax_msftx, wmax_msfty, wrainc, wraincv, wrainsh, wrainshv, wrainnc, wrainncv, & |
---|
2975 | wsnownc, wsnowncv, wgraupelnc, wgraupelncv, whailnc, whailncv, wstepave_count, & |
---|
2976 | wcldfra, wswdown, wglw, wswnorm, wolr, wemiss, wtmn, wust, wpblh, whfx, & |
---|
2977 | wqfx, wlh, wachfx, waclhf, wsnowc, wsr, wpotevp, wsnopcx, wsoiltb) |
---|
2978 | ! wsave_topo_from_real, wavg_fuel_frac, wuah, wvah, wseed1, wseed2) |
---|
2979 | ! Subroutine to get LMDZ output and get into WRF standard output variables |
---|
2980 | ! (as they appear in the wrfout) |
---|
2981 | |
---|
2982 | ! LMDZ modules |
---|
2983 | USE indice_sol_mod |
---|
2984 | USE phys_state_var_mod |
---|
2985 | USE phys_output_var_mod |
---|
2986 | ! USE pbl_surface_mod |
---|
2987 | USE phys_local_var_mod |
---|
2988 | USE comgeomphy |
---|
2989 | |
---|
2990 | ! WRF_LMDZ modules |
---|
2991 | USE output_lmdz_NOmodule |
---|
2992 | USE lmdz_wrf_variables_mod |
---|
2993 | |
---|
2994 | IMPLICIT NONE |
---|
2995 | |
---|
2996 | INTEGER, INTENT(IN) :: is, ie, js, je, ks, ke, & |
---|
2997 | ddx, ddy, ddz, bdy, dlmdz, xp, yp, lmdp, Nks, Nz0, nqtot, ivap, iliq, & |
---|
2998 | iflag_pbl, iflag_con, iflag_thermals, iflag_wake, itap, nslayers |
---|
2999 | REAL, INTENT(IN) :: pdtphys |
---|
3000 | REAL, DIMENSION(dlmdz), INTENT(IN) :: pphis, d_ps |
---|
3001 | REAL, DIMENSION(ddz), INTENT(IN) :: presnivs |
---|
3002 | REAL, DIMENSION(dlmdz,ddz+1), INTENT(IN) :: paprs |
---|
3003 | REAL, DIMENSION(dlmdz,ddz), INTENT(IN) :: pplay, pphi, d_u, d_v, d_t |
---|
3004 | REAL, DIMENSION(dlmdz,ddz,nqtot), INTENT(IN) :: qx, d_qx |
---|
3005 | REAL, DIMENSION(Nz0), INTENT(IN) :: orchidlev |
---|
3006 | REAL, DIMENSION(dlmdz,Nz0,Nks), INTENT(IN) :: tsoilorchid |
---|
3007 | INTEGER, INTENT(INOUT) :: wstepave_count |
---|
3008 | REAL, INTENT(INOUT) :: wresm, wzetatop, wtiso, & |
---|
3009 | wmax_msftx, wmax_msfty |
---|
3010 | REAL, DIMENSION(1:nslayers), INTENT(INOUT) :: wzs, wdzs |
---|
3011 | REAL, DIMENSION(1:nslayers), INTENT(IN) :: wslay |
---|
3012 | REAL, DIMENSION(is:ie,js:je), INTENT(INOUT) :: wnest_pos, wq2, wt2, & |
---|
3013 | wth2, wpsfc, wu10, wv10, wsfroff, wudroff, wgrdflx, wacgrdflx, wsnow, wsnowh, & |
---|
3014 | wcanwat, wsst, wsstsk, wlai, wsina, wtsk, wrainc, wraincv, wrainsh, wrainshv, & |
---|
3015 | wrainnc, wrainncv, wsnownc, wsnowncv, wgraupelnc, wgraupelncv, whailnc, & |
---|
3016 | whailncv, wswdown, wglw, wswnorm, wolr, wemiss, wtmn, wust, wpblh, whfx, & |
---|
3017 | wqfx, wlh, wachfx, waclhf, wsnowc, wsr, wpotevp, wsnopcx, wsoiltb |
---|
3018 | REAL, DIMENSION(is:ie,1:nslayers,js:je), INTENT(INOUT) :: wtslb, wsmois, wsh2o, & |
---|
3019 | wsmcrel |
---|
3020 | REAL, DIMENSION(is:ie,ks:ke,js:je), INTENT(INOUT) :: wh_diabatic, wcldfra |
---|
3021 | ! EXTERNAL suphel |
---|
3022 | |
---|
3023 | ! Local |
---|
3024 | INTEGER :: dx, dy, ix, iy, ixy, iz |
---|
3025 | INTEGER :: nsrf |
---|
3026 | REAL :: varixy |
---|
3027 | REAL, DIMENSION(Nz0) :: yvals1 |
---|
3028 | REAL :: RDAY, RG, RD |
---|
3029 | CHARACTER(LEN=50) :: fname |
---|
3030 | |
---|
3031 | !!!!!!! Variables |
---|
3032 | ! pwrf_grid: WRF grid structure |
---|
3033 | |
---|
3034 | fname = "'put_lmdz_in_WRFout'" |
---|
3035 | |
---|
3036 | ! L. Fita, LMD January 2014. It should work using suphel, but not? |
---|
3037 | ! CALL suphel |
---|
3038 | |
---|
3039 | RDAY = 86400. |
---|
3040 | RG = 9.81 |
---|
3041 | RD = 287. |
---|
3042 | |
---|
3043 | wdzs(1)=wslay(nslayers) |
---|
3044 | DO iz=1,nslayers-1 |
---|
3045 | wdzs(iz+1) = wslay(nslayers-iz+1) - wslay(nslayers-iz) |
---|
3046 | END DO |
---|
3047 | ! wmax_msftx = |
---|
3048 | ! wmax_msfty = |
---|
3049 | ! wresm = |
---|
3050 | ! wstepave_count = |
---|
3051 | ! wtiso = |
---|
3052 | ! wzetatop = |
---|
3053 | wzs = wslay |
---|
3054 | |
---|
3055 | DO ix=1,ddx |
---|
3056 | DO iy=1,ddy |
---|
3057 | ixy=ddx*(iy-1)+ix |
---|
3058 | wacgrdflx(ix,iy) = wacgrdflx(ix,iy) + 0. |
---|
3059 | wachfx(ix,iy) = wachfx(ix,iy) + sens(ixy) |
---|
3060 | waclhf(ix,iy) = waclhf(ix,iy) + zxfluxlat(ixy) |
---|
3061 | DO iz=1,ddz |
---|
3062 | wcldfra(ix,iz,iy) = cldfra(ixy,iz) |
---|
3063 | END DO |
---|
3064 | wcanwat(ix,iy) = 0. |
---|
3065 | wglw(ix,iy) = sollwdown(ixy) |
---|
3066 | ! These are precipitation fluxes !! |
---|
3067 | ! wgraupelnc(ix,iy) = wgraupelnc(ix,iy) + prmflxs(ixy) + prfs(ixy,iz) |
---|
3068 | ! wgraupelncv(ix,iy) = prmflxs(ixy) + prfs(ixy,iz) |
---|
3069 | wgrdflx(ix,iy) = 0. |
---|
3070 | whailnc(ix,iy) = whailnc(ix,iy) + 0. |
---|
3071 | whailncv(ix,iy) = 0. |
---|
3072 | wh_diabatic(ix,ks:ke,iy) = 0. |
---|
3073 | whfx(ix,iy) = sens(ixy) |
---|
3074 | wlai(ix,iy) = 0. |
---|
3075 | wlh(ix,iy) = zxfluxlat(ixy) |
---|
3076 | wnest_pos(ix,iy) = 0. |
---|
3077 | !! wnoahres(ix,iy) = 0. |
---|
3078 | wolr(ix,iy) = toplw(ixy) |
---|
3079 | wpblh(ix,iy) = s_pblh(ixy) |
---|
3080 | ! Should be something from evappot_srf ? |
---|
3081 | wpotevp(ix,iy) = wpotevp(ix,iy) + 0. |
---|
3082 | wpsfc(ix,iy) = paprs(ixy,1) |
---|
3083 | wq2(ix,iy) = zq2m(ixy) |
---|
3084 | wqfx(ix,iy) = 0. |
---|
3085 | wrainc(ix,iy) = wrainc(ix,iy) + rain_fall(ixy)*pdtphys |
---|
3086 | wraincv(ix,iy) = rain_fall(ixy)*pdtphys |
---|
3087 | wrainnc(ix,iy) = wrainnc(ix,iy) + rain_lsc(ixy)*pdtphys |
---|
3088 | wrainncv(ix,iy) = rain_lsc(ixy)*pdtphys |
---|
3089 | wrainsh(ix,iy) = wrainsh(ix,iy) + 0. |
---|
3090 | wrainshv(ix,iy) = 0. |
---|
3091 | ! L. Fita, LMD. February 2014 |
---|
3092 | ! LMDZ run_off_lic is for continental ice. We could try to get run_off_ter, but.... |
---|
3093 | ! wsfroff(ix,iy) = 0. |
---|
3094 | DO iz=1,Nks |
---|
3095 | yvals1(iz) = tsoilorchid(ixy,iz,1)*pctsrf(ixy,is_ter) + & |
---|
3096 | tsoilorchid(ixy,iz,2)*pctsrf(ixy,is_lic) + & |
---|
3097 | tsoilorchid(ixy,iz,3)*pctsrf(ixy,is_oce) + & |
---|
3098 | tsoilorchid(ixy,iz,4)*pctsrf(ixy,is_sic) |
---|
3099 | END DO |
---|
3100 | CALL interpolate_layers(Nz0,orchidlev,yvals1,nslayers,wslay,wtslb(ix,:,iy)) |
---|
3101 | ! No way ??? |
---|
3102 | ! wsh2o(is:ie,iz,js:je) = wsh2o(is:ie,iz,js:je) + |
---|
3103 | ! wsmcrel(is:ie,iz,js:je) = 0. |
---|
3104 | |
---|
3105 | ! wsinalpha(ix,iy) = |
---|
3106 | wsina(ix,iy) = 0. |
---|
3107 | wsnopcx(ix,iy) = wsnopcx(ix,iy) + 0. |
---|
3108 | !! wsnowh(ix,iy) = snowght(ixy) |
---|
3109 | wsnownc(ix,iy) = wsnownc(ix,iy) + snow_fall(ixy)*pdtphys |
---|
3110 | wsnowncv(ix,iy) = snow_fall(ixy)*pdtphys |
---|
3111 | wsoiltb(ix,iy) = wtslb(ix,nslayers,iy) |
---|
3112 | IF (rain_fall(ixy) /= 0.) THEN |
---|
3113 | wsr(ix,iy) = ( snow_fall(ixy) ) / ( rain_fall(ixy)) |
---|
3114 | ELSE |
---|
3115 | wsr(ix,iy) = 0. |
---|
3116 | END IF |
---|
3117 | wsstsk(ix,iy) = ftsol(ixy,is_oce) |
---|
3118 | wsst(ix,iy) = sst(ixy) |
---|
3119 | wswdown(ix,iy) = swdn(ixy,1)*pctsrf(ixy,is_ter) + & |
---|
3120 | swdn(ixy,2)*pctsrf(ixy,is_lic) + swdn(ixy,3)*pctsrf(ixy,is_oce) + & |
---|
3121 | swdn(ixy,4)*pctsrf(ixy,is_sic) |
---|
3122 | wswnorm(ix,iy) = 0. |
---|
3123 | wt2(ix,iy) = zt2m(ixy) |
---|
3124 | ! This is not exact!! |
---|
3125 | CALL temp_to_theta1(zt2m(ixy), paprs(ixy,1), wth2(ix,iy)) |
---|
3126 | wtsk(ix,iy) = ftsol(ixy, is_ter)*pctsrf(ixy,is_ter) + ftsol(ixy, is_lic)* & |
---|
3127 | pctsrf(ixy,is_lic) + ftsol(ixy, is_oce)*pctsrf(ixy,is_oce) + & |
---|
3128 | ftsol(ixy, is_sic)*pctsrf(ixy,is_sic) |
---|
3129 | wu10(ix,iy) = zu10m(ixy) |
---|
3130 | wudroff(ix,iy) = wudroff(ix,iy) + 0. |
---|
3131 | wv10(ix,iy) = zv10m(ixy) |
---|
3132 | END DO |
---|
3133 | END DO |
---|
3134 | |
---|
3135 | END SUBROUTINE put_lmdz_in_WRFout |
---|
3136 | |
---|
3137 | END MODULE wrf_lmdz_mod |
---|
3138 | |
---|