module_soil_pre.F
References to this file elsewhere.
1 #if ( ! NMM_CORE == 1 )
2
3 MODULE module_soil_pre
4
5 USE module_date_time
6 USE module_state_description
7
8 CONTAINS
9
10 SUBROUTINE adjust_for_seaice_pre ( xice , landmask , tsk , ivgtyp , vegcat , lu_index , &
11 xland , landusef , isltyp , soilcat , soilctop , &
12 soilcbot , tmn , &
13 seaice_threshold , &
14 num_veg_cat , num_soil_top_cat , num_soil_bot_cat , &
15 iswater , isice , &
16 scheme , &
17 ids , ide , jds , jde , kds , kde , &
18 ims , ime , jms , jme , kms , kme , &
19 its , ite , jts , jte , kts , kte )
20
21 IMPLICIT NONE
22
23 INTEGER , INTENT(IN) :: ids , ide , jds , jde , kds , kde , &
24 ims , ime , jms , jme , kms , kme , &
25 its , ite , jts , jte , kts , kte , &
26 iswater , isice
27 INTEGER , INTENT(IN) :: num_veg_cat , num_soil_top_cat , num_soil_bot_cat , scheme
28
29 REAL , DIMENSION(ims:ime,1:num_veg_cat,jms:jme) , INTENT(INOUT):: landusef
30 REAL , DIMENSION(ims:ime,1:num_soil_top_cat,jms:jme) , INTENT(INOUT):: soilctop
31 REAL , DIMENSION(ims:ime,1:num_soil_bot_cat,jms:jme) , INTENT(INOUT):: soilcbot
32 INTEGER , DIMENSION(ims:ime,jms:jme), INTENT(OUT) :: isltyp , ivgtyp
33 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(INOUT) :: landmask , xice , tsk , lu_index , &
34 vegcat, xland , soilcat , tmn
35 REAL , INTENT(IN) :: seaice_threshold
36
37 INTEGER :: i , j , num_seaice_changes , loop
38 CHARACTER (LEN=132) :: message
39
40 num_seaice_changes = 0
41 fix_seaice : SELECT CASE ( scheme )
42
43 CASE ( SLABSCHEME )
44 DO j = jts , MIN(jde-1,jte)
45 DO i = its , MIN(ide-1,ite)
46 IF ( xice(i,j) .GT. 200.0 ) THEN
47 xice(i,j) = 0.
48 num_seaice_changes = num_seaice_changes + 1
49 END IF
50 END DO
51 END DO
52 IF ( num_seaice_changes .GT. 0 ) THEN
53 WRITE ( message , FMT='(A,I6)' ) &
54 'Total pre number of sea ice locations removed (due to FLAG values) = ', &
55 num_seaice_changes
56 CALL wrf_debug ( 0 , message )
57 END IF
58 num_seaice_changes = 0
59 DO j = jts , MIN(jde-1,jte)
60 DO i = its , MIN(ide-1,ite)
61 IF ( ( xice(i,j) .GE. 0.5 ) .OR. &
62 ( ( landmask(i,j) .LT. 0.5 ) .AND. ( tsk(i,j) .LT. seaice_threshold ) ) ) THEN
63 xice(i,j) = 1.
64 num_seaice_changes = num_seaice_changes + 1
65 if(landmask(i,j) .LT. 0.5 )tmn(i,j) = 271.4
66 vegcat(i,j)=isice
67 ivgtyp(i,j)=isice
68 lu_index(i,j)=isice
69 landmask(i,j)=1.
70 xland(i,j)=1.
71 DO loop=1,num_veg_cat
72 landusef(i,loop,j)=0.
73 END DO
74 landusef(i,ivgtyp(i,j),j)=1.
75
76 isltyp(i,j) = 16
77 soilcat(i,j)=isltyp(i,j)
78 DO loop=1,num_soil_top_cat
79 soilctop(i,loop,j)=0
80 END DO
81 DO loop=1,num_soil_bot_cat
82 soilcbot(i,loop,j)=0
83 END DO
84 soilctop(i,isltyp(i,j),j)=1.
85 soilcbot(i,isltyp(i,j),j)=1.
86 END IF
87 END DO
88 END DO
89 IF ( num_seaice_changes .GT. 0 ) THEN
90 WRITE ( message , FMT='(A,I6)' ) &
91 'Total pre number of sea ice location changes (water to land) = ', num_seaice_changes
92 CALL wrf_debug ( 0 , message )
93 END IF
94
95 CASE ( LSMSCHEME , RUCLSMSCHEME )
96 num_seaice_changes = 0
97 DO j = jts , MIN(jde-1,jte)
98 DO i = its , MIN(ide-1,ite)
99 IF ( landmask(i,j) .GT. 0.5 ) THEN
100 if (xice(i,j).gt.0) num_seaice_changes = num_seaice_changes + 1
101 xice(i,j) = 0.
102 END IF
103 END DO
104 END DO
105 IF ( num_seaice_changes .GT. 0 ) THEN
106 WRITE ( message , FMT='(A,I6)' ) &
107 'Total pre number of land location changes (seaice set to zero) = ', num_seaice_changes
108 CALL wrf_debug ( 0 , message )
109 END IF
110
111 END SELECT fix_seaice
112
113 END SUBROUTINE adjust_for_seaice_pre
114
115 SUBROUTINE adjust_for_seaice_post ( xice , landmask , tsk_old , tsk , ivgtyp , vegcat , lu_index , &
116 xland , landusef , isltyp , soilcat , soilctop , &
117 soilcbot , tmn , vegfra , &
118 tslb , smois , sh2o , &
119 seaice_threshold , &
120 num_veg_cat , num_soil_top_cat , num_soil_bot_cat , &
121 num_soil_layers , &
122 iswater , isice , &
123 scheme , &
124 ids , ide , jds , jde , kds , kde , &
125 ims , ime , jms , jme , kms , kme , &
126 its , ite , jts , jte , kts , kte )
127
128 IMPLICIT NONE
129
130 INTEGER , INTENT(IN) :: ids , ide , jds , jde , kds , kde , &
131 ims , ime , jms , jme , kms , kme , &
132 its , ite , jts , jte , kts , kte , &
133 iswater , isice
134 INTEGER , INTENT(IN) :: num_veg_cat , num_soil_top_cat , num_soil_bot_cat , scheme
135 INTEGER , INTENT(IN) :: num_soil_layers
136
137 REAL , DIMENSION(ims:ime,1:num_veg_cat,jms:jme) , INTENT(INOUT):: landusef
138 REAL , DIMENSION(ims:ime,1:num_soil_top_cat,jms:jme) , INTENT(INOUT):: soilctop
139 REAL , DIMENSION(ims:ime,1:num_soil_bot_cat,jms:jme) , INTENT(INOUT):: soilcbot
140 REAL , DIMENSION(ims:ime,1:num_soil_layers,jms:jme) , INTENT(INOUT):: tslb , smois , sh2o
141 INTEGER , DIMENSION(ims:ime,jms:jme), INTENT(OUT) :: isltyp , ivgtyp
142 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(INOUT) :: landmask , xice , tsk , lu_index , &
143 vegcat, xland , soilcat , tmn , &
144 tsk_old , vegfra
145 REAL , INTENT(IN) :: seaice_threshold
146 REAL :: total_depth , mid_point_depth
147
148 INTEGER :: i , j , num_seaice_changes , loop
149 CHARACTER (LEN=132) :: message
150
151 num_seaice_changes = 0
152 fix_seaice : SELECT CASE ( scheme )
153
154 CASE ( SLABSCHEME )
155
156 CASE ( LSMSCHEME , RUCLSMSCHEME )
157 DO j = jts , MIN(jde-1,jte)
158 DO i = its , MIN(ide-1,ite)
159 IF ( xice(i,j) .GT. 200.0 ) THEN
160 xice(i,j) = 0.
161 num_seaice_changes = num_seaice_changes + 1
162 END IF
163 END DO
164 END DO
165 IF ( num_seaice_changes .GT. 0 ) THEN
166 WRITE ( message , FMT='(A,I6)' ) &
167 'Total post number of sea ice locations removed (due to FLAG values) = ', &
168 num_seaice_changes
169 CALL wrf_debug ( 0 , message )
170 END IF
171 num_seaice_changes = 0
172 DO j = jts , MIN(jde-1,jte)
173 DO i = its , MIN(ide-1,ite)
174 IF ( ( ( tsk(i,j) .LT. 170 ) .OR. ( tsk(i,j) .GT. 400 ) ) .AND. &
175 ( ( tsk_old(i,j) .GT. 170 ) .AND. ( tsk_old(i,j) .LT. 400 ) ) )THEN
176 tsk(i,j) = tsk_old(i,j)
177 END IF
178 IF ( ( ( tsk(i,j) .LT. 170 ) .OR. ( tsk(i,j) .GT. 400 ) ) .AND. &
179 ( ( tsk_old(i,j) .LT. 170 ) .OR. ( tsk_old(i,j) .GT. 400 ) ) )THEN
180 print *,'TSK woes in seaice post, i,j=',i,j,' tsk = ',tsk(i,j), tsk_old(i,j)
181 CALL wrf_error_fatal ( 'TSK is unrealistic, problems for seaice post')
182 ELSE IF ( ( xice(i,j) .GE. 0.5 ) .OR. &
183 ( ( landmask(i,j) .LT. 0.5 ) .AND. ( tsk(i,j) .LT. seaice_threshold ) ) ) THEN
184 xice(i,j) = 1.
185 num_seaice_changes = num_seaice_changes + 1
186 if(landmask(i,j) .LT. 0.5 )tmn(i,j) = 271.4
187 vegcat(i,j)=isice
188 ivgtyp(i,j)=isice
189 lu_index(i,j)=isice
190 landmask(i,j)=1.
191 xland(i,j)=1.
192 vegfra(i,j)=0.
193 DO loop=1,num_veg_cat
194 landusef(i,loop,j)=0.
195 END DO
196 landusef(i,ivgtyp(i,j),j)=1.
197
198 tsk_old(i,j) = tsk(i,j)
199
200 isltyp(i,j) = 16
201 soilcat(i,j)=isltyp(i,j)
202 DO loop=1,num_soil_top_cat
203 soilctop(i,loop,j)=0
204 END DO
205 DO loop=1,num_soil_bot_cat
206 soilcbot(i,loop,j)=0
207 END DO
208 soilctop(i,isltyp(i,j),j)=1.
209 soilcbot(i,isltyp(i,j),j)=1.
210
211 total_depth = 3. ! ice is 3 m deep, num_soil_layers equispaced layers
212 DO loop = 1,num_soil_layers
213 mid_point_depth=(total_depth/num_soil_layers)/2. + &
214 (loop-1)*(total_depth/num_soil_layers)
215 tslb(i,loop,j) = ( (total_depth-mid_point_depth)*tsk(i,j) + &
216 mid_point_depth*tmn(i,j) ) / total_depth
217 END DO
218
219 DO loop=1,num_soil_layers
220 smois(i,loop,j) = 1.0
221 sh2o(i,loop,j) = 0.0
222 END DO
223 ELSE IF ( xice(i,j) .LT. 0.5 ) THEN
224 xice(i,j) = 0.
225 END IF
226 END DO
227 END DO
228 IF ( num_seaice_changes .GT. 0 ) THEN
229 WRITE ( message , FMT='(A,I6)' ) &
230 'Total post number of sea ice location changes (water to land) = ', num_seaice_changes
231 CALL wrf_debug ( 0 , message )
232 END IF
233
234 END SELECT fix_seaice
235
236 END SUBROUTINE adjust_for_seaice_post
237
238 SUBROUTINE process_percent_cat_new ( landmask , &
239 landuse_frac , soil_top_cat , soil_bot_cat , &
240 isltyp , ivgtyp , &
241 num_veg_cat , num_soil_top_cat , num_soil_bot_cat , &
242 ids , ide , jds , jde , kds , kde , &
243 ims , ime , jms , jme , kms , kme , &
244 its , ite , jts , jte , kts , kte , &
245 iswater )
246
247 IMPLICIT NONE
248
249 INTEGER , INTENT(IN) :: ids , ide , jds , jde , kds , kde , &
250 ims , ime , jms , jme , kms , kme , &
251 its , ite , jts , jte , kts , kte , &
252 iswater
253 INTEGER , INTENT(IN) :: num_veg_cat , num_soil_top_cat , num_soil_bot_cat
254 REAL , DIMENSION(ims:ime,1:num_veg_cat,jms:jme) , INTENT(INOUT):: landuse_frac
255 REAL , DIMENSION(ims:ime,1:num_soil_top_cat,jms:jme) , INTENT(IN):: soil_top_cat
256 REAL , DIMENSION(ims:ime,1:num_soil_bot_cat,jms:jme) , INTENT(IN):: soil_bot_cat
257 INTEGER , DIMENSION(ims:ime,jms:jme), INTENT(OUT) :: isltyp , ivgtyp
258 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(INOUT) :: landmask
259
260 INTEGER :: i , j , l , ll, dominant_index
261 REAL :: dominant_value
262
263 #ifdef WRF_CHEM
264 ! REAL :: lwthresh = .99
265 REAL :: lwthresh = .50
266 #else
267 REAL :: lwthresh = .50
268 #endif
269
270 INTEGER , PARAMETER :: iswater_soil = 14
271 INTEGER :: iforce
272 CHARACTER (LEN=132) :: message
273 integer :: change_water , change_land
274 change_water = 0
275 change_land = 0
276
277 ! Sanity check on the 50/50 points
278
279 DO j = jts , MIN(jde-1,jte)
280 DO i = its , MIN(ide-1,ite)
281 dominant_value = landuse_frac(i,iswater,j)
282 IF ( dominant_value .EQ. lwthresh ) THEN
283 DO l = 1 , num_veg_cat
284 IF ( l .EQ. iswater ) CYCLE
285 IF ( ( landuse_frac(i,l,j) .EQ. lwthresh ) .AND. ( landmask(i,j) .LT. 0.5 ) ) THEN
286 PRINT *,i,j,' water and category ',l,' both at 50%, landmask is ',landmask(i,j)
287 landuse_frac(i,l,j) = lwthresh - .01
288 landuse_frac(i,iswater,j) = lwthresh + 0.01
289 ELSE IF ( ( landuse_frac(i,l,j) .EQ. lwthresh ) .AND. ( landmask(i,j) .GT. 0.5 ) ) THEN
290 PRINT *,i,j,' water and category ',l,' both at 50%, landmask is ',landmask(i,j)
291 landuse_frac(i,l,j) = lwthresh + .01
292 landuse_frac(i,iswater,j) = lwthresh - 0.01
293 END IF
294 END DO
295 END IF
296 END DO
297 END DO
298
299 ! Compute the dominant VEGETATION INDEX.
300
301 DO j = jts , MIN(jde-1,jte)
302 DO i = its , MIN(ide-1,ite)
303 dominant_value = landuse_frac(i,1,j)
304 dominant_index = 1
305 DO l = 2 , num_veg_cat
306 IF ( l .EQ. iswater ) THEN
307 ! wait a bit
308 ELSE IF ( ( l .NE. iswater ) .AND. ( landuse_frac(i,l,j) .GT. dominant_value ) ) THEN
309 dominant_value = landuse_frac(i,l,j)
310 dominant_index = l
311 END IF
312 END DO
313 IF ( landuse_frac(i,iswater,j) .GT. lwthresh ) THEN
314 dominant_value = landuse_frac(i,iswater,j)
315 dominant_index = iswater
316 ELSE IF ( ( landuse_frac(i,iswater,j) .EQ. lwthresh) .AND. &
317 ( landmask(i,j) .LT. 0.5) .AND. &
318 ( dominant_value .EQ. lwthresh) ) THEN
319 dominant_value = landuse_frac(i,iswater,j)
320 dominant_index = iswater
321 ELSE IF ( ( landuse_frac(i,iswater,j) .EQ. lwthresh) .AND. &
322 ( landmask(i,j) .GT. 0.5) .AND. &
323 ( dominant_value .EQ. lwthresh) ) THEN
324 !no op
325 ELSE IF ( ( landuse_frac(i,iswater,j) .EQ. lwthresh ) .AND. &
326 ( dominant_value .LT. lwthresh ) ) THEN
327 dominant_value = landuse_frac(i,iswater,j)
328 dominant_index = iswater
329 END IF
330 IF ( dominant_index .EQ. iswater ) THEN
331 if(landmask(i,j).gt.lwthresh) then
332 !print *,'changing to water at point ',i,j
333 !print '(24(i3,1x))',1, 2, 3, 4, 5, 6, 7, 8, 9, 10,11,12, 13, 14, 15, 16, 17,18,19,20,21, 22, 23,24
334 !print '(24(i3,1x))',nint(landuse_frac(i,:,j)*100)
335 change_water=change_water+1
336 endif
337 landmask(i,j) = 0
338 ELSE IF ( dominant_index .NE. iswater ) THEN
339 if(landmask(i,j).lt.lwthresh) then
340 !print *,'changing to land at point ',i,j
341 !print '(24(i3,1x))',1, 2, 3, 4, 5, 6, 7, 8, 9, 10,11,12, 13, 14, 15, 16, 17,18,19,20,21, 22, 23,24
342 !print '(24(i3,1x))',nint(landuse_frac(i,:,j)*100)
343 change_land=change_land+1
344 endif
345 landmask(i,j) = 1
346 END IF
347 ivgtyp(i,j) = dominant_index
348 END DO
349 END DO
350
351 ! Compute the dominant SOIL TEXTURE INDEX, TOP.
352
353 iforce = 0
354 DO i = its , MIN(ide-1,ite)
355 DO j = jts , MIN(jde-1,jte)
356 dominant_value = soil_top_cat(i,1,j)
357 dominant_index = 1
358 IF ( landmask(i,j) .GT. lwthresh ) THEN
359 DO l = 2 , num_soil_top_cat
360 IF ( ( l .NE. iswater_soil ) .AND. ( soil_top_cat(i,l,j) .GT. dominant_value ) ) THEN
361 dominant_value = soil_top_cat(i,l,j)
362 dominant_index = l
363 END IF
364 END DO
365 IF ( dominant_value .LT. 0.01 ) THEN
366 iforce = iforce + 1
367 WRITE ( message , FMT = '(A,I4,I4)' ) &
368 'based on landuse, changing soil to land at point ',i,j
369 CALL wrf_debug(1,message)
370 WRITE ( message , FMT = '(16(i3,1x))' ) &
371 1, 2, 3, 4, 5, 6, 7, 8, 9, 10,11,12, 13, 14, 15, 16
372 CALL wrf_debug(1,message)
373 WRITE ( message , FMT = '(16(i3,1x))' ) &
374 nint(soil_top_cat(i,:,j)*100)
375 CALL wrf_debug(1,message)
376 dominant_index = 8
377 END IF
378 ELSE
379 dominant_index = iswater_soil
380 END IF
381 isltyp(i,j) = dominant_index
382 END DO
383 END DO
384
385 if(iforce.ne.0)then
386 WRITE(message,FMT='(A,I4,A,I6)' ) &
387 'forcing artificial silty clay loam at ',iforce,' points, out of ',&
388 (MIN(ide-1,ite)-its+1)*(MIN(jde-1,jte)-jts+1)
389 CALL wrf_debug(0,message)
390 endif
391 print *,'LAND CHANGE = ',change_land
392 print *,'WATER CHANGE = ',change_water
393
394 END SUBROUTINE process_percent_cat_new
395
396 SUBROUTINE process_soil_real ( tsk , tmn , &
397 landmask , sst , &
398 st_input , sm_input , sw_input , st_levels_input , sm_levels_input , sw_levels_input , &
399 zs , dzs , tslb , smois , sh2o , &
400 flag_sst , flag_soilt000, flag_soilm000, &
401 ids , ide , jds , jde , kds , kde , &
402 ims , ime , jms , jme , kms , kme , &
403 its , ite , jts , jte , kts , kte , &
404 sf_surface_physics , num_soil_layers , real_data_init_type , &
405 num_st_levels_input , num_sm_levels_input , num_sw_levels_input , &
406 num_st_levels_alloc , num_sm_levels_alloc , num_sw_levels_alloc )
407
408 IMPLICIT NONE
409
410 INTEGER , INTENT(IN) :: ids , ide , jds , jde , kds , kde , &
411 ims , ime , jms , jme , kms , kme , &
412 its , ite , jts , jte , kts , kte , &
413 sf_surface_physics , num_soil_layers , real_data_init_type , &
414 num_st_levels_input , num_sm_levels_input , num_sw_levels_input , &
415 num_st_levels_alloc , num_sm_levels_alloc , num_sw_levels_alloc
416
417 INTEGER , INTENT(IN) :: flag_sst, flag_soilt000, flag_soilm000
418
419 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(IN) :: landmask , sst
420
421 INTEGER , DIMENSION(1:num_st_levels_input) , INTENT(INOUT) :: st_levels_input
422 INTEGER , DIMENSION(1:num_sm_levels_input) , INTENT(INOUT) :: sm_levels_input
423 INTEGER , DIMENSION(1:num_sw_levels_input) , INTENT(INOUT) :: sw_levels_input
424 REAL , DIMENSION(ims:ime,1:num_st_levels_alloc,jms:jme) , INTENT(INOUT) :: st_input
425 REAL , DIMENSION(ims:ime,1:num_sm_levels_alloc,jms:jme) , INTENT(INOUT) :: sm_input
426 REAL , DIMENSION(ims:ime,1:num_sw_levels_alloc,jms:jme) , INTENT(INOUT) :: sw_input
427
428 REAL, DIMENSION(1:num_soil_layers), INTENT(OUT) :: zs,dzs
429 REAL , DIMENSION(ims:ime,num_soil_layers,jms:jme) , INTENT(OUT) :: tslb , smois , sh2o
430 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(IN) :: tmn
431 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(INOUT) :: tsk
432
433 INTEGER :: i , j , l , dominant_index , num_soil_cat , num_veg_cat
434 REAL :: dominant_value
435
436 ! Initialize the soil depth, and the soil temperature and moisture.
437
438 IF ( ( sf_surface_physics .EQ. 1 ) .AND. ( num_soil_layers .GT. 1 ) ) THEN
439 CALL init_soil_depth_1 ( zs , dzs , num_soil_layers )
440 CALL init_soil_1_real ( tsk , tmn , tslb , zs , dzs , num_soil_layers , real_data_init_type , &
441 landmask , sst , flag_sst , &
442 ids , ide , jds , jde , kds , kde , &
443 ims , ime , jms , jme , kms , kme , &
444 its , ite , jts , jte , kts , kte )
445
446 ELSE IF ( ( sf_surface_physics .EQ. 2 ) .AND. ( num_soil_layers .GT. 1 ) ) THEN
447 CALL init_soil_depth_2 ( zs , dzs , num_soil_layers )
448 CALL init_soil_2_real ( tsk , tmn , smois , sh2o , tslb , &
449 st_input , sm_input , sw_input , landmask , sst , &
450 zs , dzs , &
451 st_levels_input , sm_levels_input , sw_levels_input , &
452 num_soil_layers , num_st_levels_input , num_sm_levels_input , num_sw_levels_input , &
453 num_st_levels_alloc , num_sm_levels_alloc , num_sw_levels_alloc , &
454 flag_sst , flag_soilt000 , flag_soilm000 , &
455 ids , ide , jds , jde , kds , kde , &
456 ims , ime , jms , jme , kms , kme , &
457 its , ite , jts , jte , kts , kte )
458 ! CALL init_soil_old ( tsk , tmn , &
459 ! smois , tslb , zs , dzs , num_soil_layers , &
460 ! st000010_input , st010040_input , st040100_input , st100200_input , &
461 ! st010200_input , &
462 ! sm000010_input , sm010040_input , sm040100_input , sm100200_input , &
463 ! sm010200_input , &
464 ! landmask_input , sst_input , &
465 ! ids , ide , jds , jde , kds , kde , &
466 ! ims , ime , jms , jme , kms , kme , &
467 ! its , ite , jts , jte , kts , kte )
468 ELSE IF ( ( sf_surface_physics .EQ. 3 ) .AND. ( num_soil_layers .GT. 1 ) ) THEN
469 CALL init_soil_depth_3 ( zs , dzs , num_soil_layers )
470 CALL init_soil_3_real ( tsk , tmn , smois , tslb , &
471 st_input , sm_input , landmask , sst , &
472 zs , dzs , &
473 st_levels_input , sm_levels_input , &
474 num_soil_layers , num_st_levels_input , num_sm_levels_input , &
475 num_st_levels_alloc , num_sm_levels_alloc , &
476 flag_sst , flag_soilt000 , flag_soilm000 , &
477 ids , ide , jds , jde , kds , kde , &
478 ims , ime , jms , jme , kms , kme , &
479 its , ite , jts , jte , kts , kte )
480 END IF
481
482 END SUBROUTINE process_soil_real
483
484 SUBROUTINE process_soil_ideal ( xland,xice,vegfra,snow,canwat, &
485 ivgtyp,isltyp,tslb,smois, &
486 tsk,tmn,zs,dzs, &
487 num_soil_layers, &
488 sf_surface_physics , &
489 ids,ide, jds,jde, kds,kde,&
490 ims,ime, jms,jme, kms,kme,&
491 its,ite, jts,jte, kts,kte )
492
493 IMPLICIT NONE
494
495 INTEGER, INTENT(IN) ::ids,ide, jds,jde, kds,kde, &
496 ims,ime, jms,jme, kms,kme, &
497 its,ite, jts,jte, kts,kte
498
499 INTEGER, INTENT(IN) :: num_soil_layers , sf_surface_physics
500
501 REAL, DIMENSION( ims:ime, num_soil_layers, jms:jme ) , INTENT(INOUT) :: smois, tslb
502
503 REAL, DIMENSION(num_soil_layers), INTENT(OUT) :: dzs,zs
504
505 REAL, DIMENSION( ims:ime, jms:jme ) , INTENT(INOUT) :: tsk, tmn
506 REAL, DIMENSION( ims:ime, jms:jme ) , INTENT(OUT) :: xland, snow, canwat, xice, vegfra
507 INTEGER, DIMENSION( ims:ime, jms:jme ) , INTENT(OUT) :: ivgtyp, isltyp
508
509 ! Local variables.
510
511 INTEGER :: itf,jtf
512
513 itf=MIN(ite,ide-1)
514 jtf=MIN(jte,jde-1)
515
516 IF ( ( sf_surface_physics .EQ. 1 ) .AND. ( num_soil_layers .GT. 1 ) ) THEN
517 CALL init_soil_depth_1 ( zs , dzs , num_soil_layers )
518 CALL init_soil_1_ideal(tsk,tmn,tslb,xland, &
519 ivgtyp,zs,dzs,num_soil_layers, &
520 ids,ide, jds,jde, kds,kde, &
521 ims,ime, jms,jme, kms,kme, &
522 its,ite, jts,jte, kts,kte )
523 ELSE IF ( ( sf_surface_physics .EQ. 2 ) .AND. ( num_soil_layers .GT. 1 ) ) THEN
524 CALL init_soil_depth_2 ( zs , dzs , num_soil_layers )
525 CALL init_soil_2_ideal ( xland,xice,vegfra,snow,canwat, &
526 ivgtyp,isltyp,tslb,smois,tmn, &
527 num_soil_layers, &
528 ids,ide, jds,jde, kds,kde, &
529 ims,ime, jms,jme, kms,kme, &
530 its,ite, jts,jte, kts,kte )
531 ELSE IF ( ( sf_surface_physics .EQ. 3 ) .AND. ( num_soil_layers .GT. 1 ) ) THEN
532 CALL init_soil_depth_3 ( zs , dzs , num_soil_layers )
533
534 END IF
535
536 END SUBROUTINE process_soil_ideal
537
538 SUBROUTINE adjust_soil_temp_new ( tmn , sf_surface_physics , &
539 tsk , ter , toposoil , landmask , flag_toposoil , flag_tavgsfc , &
540 st000010 , st010040 , st040100 , st100200 , st010200 , &
541 flag_st000010 , flag_st010040 , flag_st040100 , flag_st100200 , flag_st010200 , &
542 st000007 , st007028 , st028100 , st100255 , &
543 flag_st000007 , flag_st007028 , flag_st028100 , flag_st100255 , &
544 soilt000 , soilt005 , soilt020 , soilt040 , soilt160 , soilt300 , &
545 flag_soilt000 , flag_soilt005 , flag_soilt020 , flag_soilt040 , flag_soilt160 , flag_soilt300 , &
546 ids , ide , jds , jde , kds , kde , &
547 ims , ime , jms , jme , kms , kme , &
548 its , ite , jts , jte , kts , kte )
549
550 IMPLICIT NONE
551
552 INTEGER , INTENT(IN) :: ids , ide , jds , jde , kds , kde , &
553 ims , ime , jms , jme , kms , kme , &
554 its , ite , jts , jte , kts , kte
555
556 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(IN) :: ter , toposoil , landmask
557 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(INOUT) :: tmn , tsk
558 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(INOUT) :: st000010 , st010040 , st040100 , st100200 , st010200
559 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(INOUT) :: st000007 , st007028 , st028100 , st100255
560 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(INOUT) :: soilt000 , soilt005 , soilt020 , soilt040 , soilt160 , soilt300
561
562 INTEGER , INTENT(IN) :: flag_st000010 , flag_st010040 , flag_st040100 , flag_st100200 , flag_st010200
563 INTEGER , INTENT(IN) :: flag_st000007 , flag_st007028 , flag_st028100 , flag_st100255
564 INTEGER , INTENT(IN) :: flag_soilt000 , flag_soilt005 , flag_soilt020 , flag_soilt040 , flag_soilt160 , flag_soilt300
565 INTEGER , INTENT(IN) :: sf_surface_physics , flag_toposoil , flag_tavgsfc
566
567 INTEGER :: i , j
568
569 REAL :: soil_elev_min_val , soil_elev_max_val , soil_elev_min_dif , soil_elev_max_dif
570
571 ! Adjust the annual mean temperature as if it is based on from a sea-level elevation
572 ! if the value used is from the actual annula mean data set. If the input field to
573 ! be used for tmn is one of the first-guess input temp fields, need to do an adjustment
574 ! only on the diff in topo from the model terrain and the first-guess terrain.
575
576 SELECT CASE ( sf_surface_physics )
577
578 CASE (LSMSCHEME)
579 DO j = jts , MIN(jde-1,jte)
580 DO i = its , MIN(ide-1,ite)
581 IF (landmask(i,j) .GT. 0.5 ) THEN
582 tmn(i,j) = tmn(i,j) - 0.0065 * ter(i,j)
583 END IF
584 END DO
585 END DO
586
587 CASE (RUCLSMSCHEME)
588 DO j = jts , MIN(jde-1,jte)
589 DO i = its , MIN(ide-1,ite)
590 IF (landmask(i,j) .GT. 0.5 ) THEN
591 tmn(i,j) = tmn(i,j) - 0.0065 * ter(i,j)
592 END IF
593 END DO
594 END DO
595
596 END SELECT
597
598
599 ! Do we have a soil field with which to modify soil temperatures?
600
601 IF ( flag_toposoil .EQ. 1 ) THEN
602
603 DO j = jts , MIN(jde-1,jte)
604 DO i = its , MIN(ide-1,ite)
605
606 ! Is the toposoil field OK, or is it a subversive soil elevation field. We can tell
607 ! usually by looking at values. Anything less than -1000 m (lower than the Dead Sea) is
608 ! bad. Anything larger than 10 km (taller than Everest) is toast. Also, anything where
609 ! the difference between the soil elevation and the terrain is greater than 3 km means
610 ! that the soil data is either all zeros or that the data are inconsistent. Any of these
611 ! three conditions is grievous enough to induce a WRF fatality. However, if they are at
612 ! a water point, then we can safely ignore them.
613
614 soil_elev_min_val = toposoil(i,j)
615 soil_elev_max_val = toposoil(i,j)
616 soil_elev_min_dif = ter(i,j) - toposoil(i,j)
617 soil_elev_max_dif = ter(i,j) - toposoil(i,j)
618
619 IF ( ( soil_elev_min_val .LT. -1000 ) .AND. ( landmask(i,j) .LT. 0.5 ) ) THEN
620 CYCLE
621 ELSE IF ( ( soil_elev_min_val .LT. -1000 ) .AND. ( landmask(i,j) .GT. 0.5 ) ) THEN
622 !print *,'no soil temperature elevation adjustment, soil height too small = ',toposoil(i,j)
623 cycle
624 ! CALL wrf_error_fatal ( 'TOPOSOIL values have large negative values < -1000 m, unrealistic.' )
625 ENDIF
626
627 IF ( ( soil_elev_max_val .GT. 10000 ) .AND. ( landmask(i,j) .LT. 0.5 ) ) THEN
628 CYCLE
629 ELSE IF ( ( soil_elev_max_val .GT. 10000 ) .AND. ( landmask(i,j) .GT. 0.5 ) ) THEN
630 print *,'no soil temperature elevation adjustment, soil height too high = ',toposoil(i,j)
631 cycle
632 CALL wrf_error_fatal ( 'TOPOSOIL values have large positive values > 10,000 m , unrealistic.' )
633 ENDIF
634
635 IF ( ( ( soil_elev_min_dif .LT. -3000 ) .OR. ( soil_elev_max_dif .GT. 3000 ) ) .AND. &
636 ( landmask(i,j) .LT. 0.5 ) ) THEN
637 CYCLE
638 ELSE IF ( ( ( soil_elev_min_dif .LT. -3000 ) .OR. ( soil_elev_max_dif .GT. 3000 ) ) .AND. &
639 ( landmask(i,j) .GT. 0.5 ) ) THEN
640 print *,'no soil temperature elevation adjustment, diff of soil height and terrain = ',ter(i,j) - toposoil(i,j)
641 cycle
642 CALL wrf_error_fatal ( 'TOPOSOIL difference with terrain elevation differs by more than 3000 m, unrealistic' )
643 ENDIF
644
645 ! For each of the fields that we would like to modify, check to see if it came in from the SI.
646 ! If so, then use a -6.5 K/km lapse rate (based on the elevation diffs). We only adjust when we
647 ! are not at a water point.
648
649 IF (landmask(i,j) .GT. 0.5 ) THEN
650
651 IF ( sf_surface_physics .EQ. SLABSCHEME ) THEN
652 IF ( ( flag_tavgsfc .EQ. 1 ) .OR. &
653 ( flag_st010040 .EQ. 1 ) .OR. &
654 ( flag_st000010 .EQ. 1 ) .OR. &
655 ( flag_soilt020 .EQ. 1 ) .OR. &
656 ( flag_st007028 .EQ. 1 ) ) THEN
657 tmn(i,j) = tmn(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
658 ELSE
659 tmn(i,j) = tmn(i,j) - 0.0065 * ter(i,j)
660 END IF
661 END IF
662
663 tsk(i,j) = tsk(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
664
665 IF ( flag_st000010 .EQ. 1 ) THEN
666 st000010(i,j) = st000010(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
667 END IF
668 IF ( flag_st010040 .EQ. 1 ) THEN
669 st010040(i,j) = st010040(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
670 END IF
671 IF ( flag_st040100 .EQ. 1 ) THEN
672 st040100(i,j) = st040100(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
673 END IF
674 IF ( flag_st100200 .EQ. 1 ) THEN
675 st100200(i,j) = st100200(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
676 END IF
677 IF ( flag_st010200 .EQ. 1 ) THEN
678 st010200(i,j) = st010200(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
679 END IF
680
681 IF ( flag_st000007 .EQ. 1 ) THEN
682 st000007(i,j) = st000007(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
683 END IF
684 IF ( flag_st007028 .EQ. 1 ) THEN
685 st007028(i,j) = st007028(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
686 END IF
687 IF ( flag_st028100 .EQ. 1 ) THEN
688 st028100(i,j) = st028100(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
689 END IF
690 IF ( flag_st100255 .EQ. 1 ) THEN
691 st100255(i,j) = st100255(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
692 END IF
693
694 IF ( flag_soilt000 .EQ. 1 ) THEN
695 soilt000(i,j) = soilt000(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
696 END IF
697 IF ( flag_soilt005 .EQ. 1 ) THEN
698 soilt005(i,j) = soilt005(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
699 END IF
700 IF ( flag_soilt020 .EQ. 1 ) THEN
701 soilt020(i,j) = soilt020(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
702 END IF
703 IF ( flag_soilt040 .EQ. 1 ) THEN
704 soilt040(i,j) = soilt040(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
705 END IF
706 IF ( flag_soilt160 .EQ. 1 ) THEN
707 soilt160(i,j) = soilt160(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
708 END IF
709 IF ( flag_soilt300 .EQ. 1 ) THEN
710 soilt300(i,j) = soilt300(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
711 END IF
712
713 END IF
714 END DO
715 END DO
716
717 END IF
718
719 END SUBROUTINE adjust_soil_temp_new
720
721
722 SUBROUTINE init_soil_depth_1 ( zs , dzs , num_soil_layers )
723
724 IMPLICIT NONE
725
726 INTEGER, INTENT(IN) :: num_soil_layers
727
728 REAL, DIMENSION(1:num_soil_layers), INTENT(OUT) :: zs,dzs
729
730 INTEGER :: l
731
732 ! Define layers (top layer = 0.01 m). Double the thicknesses at each step (dzs values).
733 ! The distance from the ground level to the midpoint of the layer is given by zs.
734
735 ! ------- Ground Level ---------- || || || ||
736 ! || || || || zs(1) = 0.005 m
737 ! -- -- -- -- -- -- -- -- -- || || || \/
738 ! || || ||
739 ! ----------------------------------- || || || \/ dzs(1) = 0.01 m
740 ! || || ||
741 ! || || || zs(2) = 0.02
742 ! -- -- -- -- -- -- -- -- -- || || \/
743 ! || ||
744 ! || ||
745 ! ----------------------------------- || || \/ dzs(2) = 0.02 m
746 ! || ||
747 ! || ||
748 ! || ||
749 ! || || zs(3) = 0.05
750 ! -- -- -- -- -- -- -- -- -- || \/
751 ! ||
752 ! ||
753 ! ||
754 ! ||
755 ! ----------------------------------- \/ dzs(3) = 0.04 m
756
757 IF ( num_soil_layers .NE. 5 ) THEN
758 PRINT '(A)','Usually, the 5-layer diffusion uses 5 layers. Change this in the namelist.'
759 CALL wrf_error_fatal ( '5-layer_diffusion_uses_5_layers' )
760 END IF
761
762 dzs(1)=.01
763 zs(1)=.5*dzs(1)
764
765 DO l=2,num_soil_layers
766 dzs(l)=2*dzs(l-1)
767 zs(l)=zs(l-1)+.5*dzs(l-1)+.5*dzs(l)
768 ENDDO
769
770 END SUBROUTINE init_soil_depth_1
771
772 SUBROUTINE init_soil_depth_2 ( zs , dzs , num_soil_layers )
773
774 IMPLICIT NONE
775
776 INTEGER, INTENT(IN) :: num_soil_layers
777
778 REAL, DIMENSION(1:num_soil_layers), INTENT(OUT) :: zs,dzs
779
780 INTEGER :: l
781
782 dzs = (/ 0.1 , 0.3 , 0.6 , 1.0 /)
783
784 IF ( num_soil_layers .NE. 4 ) THEN
785 PRINT '(A)','Usually, the LSM uses 4 layers. Change this in the namelist.'
786 CALL wrf_error_fatal ( 'LSM_uses_4_layers' )
787 END IF
788
789 zs(1)=.5*dzs(1)
790
791 DO l=2,num_soil_layers
792 zs(l)=zs(l-1)+.5*dzs(l-1)+.5*dzs(l)
793 ENDDO
794
795 END SUBROUTINE init_soil_depth_2
796
797 SUBROUTINE init_soil_depth_3 ( zs , dzs , num_soil_layers )
798
799 IMPLICIT NONE
800
801 INTEGER, INTENT(IN) :: num_soil_layers
802
803 REAL, DIMENSION(1:num_soil_layers), INTENT(OUT) :: zs,dzs
804
805 INTEGER :: l
806
807 CHARACTER (LEN=132) :: message
808
809 ! in RUC LSM ZS - soil levels, and DZS - soil layer thicknesses, not used
810 ! ZS is specified in the namelist: num_soil_layers = 6 or 9.
811 ! Other options with number of levels are possible, but
812 ! WRF users should change consistently the namelist entry with the
813 ! ZS array in this subroutine.
814
815 IF ( num_soil_layers .EQ. 6) THEN
816 zs = (/ 0.00 , 0.05 , 0.20 , 0.40 , 1.60 , 3.00 /)
817 ! dzs = (/ 0.00 , 0.125, 0.175 , 0.70 , 1.30 , 1.40 /)
818 ELSEIF ( num_soil_layers .EQ. 9) THEN
819 zs = (/ 0.00 , 0.05 , 0.20 , 0.40 , 0.60, 1.00, 1.60 , 2.20, 3.00 /)
820 ! dzs = (/ 0.00 , 0.125, 0.175 , 0.70 , 1.30 , 1.40 /)
821 ENDIF
822
823 IF ( num_soil_layers .EQ. 4 .OR. num_soil_layers .EQ. 5 ) THEN
824 write (message, FMT='(A)') 'The RUC LSM uses 6, 9 or more levels. Change this in the namelist.'
825 CALL wrf_error_fatal ( message )
826 END IF
827
828 END SUBROUTINE init_soil_depth_3
829
830 SUBROUTINE init_soil_1_real ( tsk , tmn , tslb , zs , dzs , &
831 num_soil_layers , real_data_init_type , &
832 landmask , sst , flag_sst , &
833 ids , ide , jds , jde , kds , kde , &
834 ims , ime , jms , jme , kms , kme , &
835 its , ite , jts , jte , kts , kte )
836
837 IMPLICIT NONE
838
839 INTEGER , INTENT(IN) :: num_soil_layers , real_data_init_type , &
840 ids , ide , jds , jde , kds , kde , &
841 ims , ime , jms , jme , kms , kme , &
842 its , ite , jts , jte , kts , kte
843
844 INTEGER , INTENT(IN) :: flag_sst
845
846 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(IN) :: landmask , sst
847 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(IN) :: tsk , tmn
848
849 REAL , DIMENSION(num_soil_layers) :: zs , dzs
850
851 REAL , DIMENSION(ims:ime,num_soil_layers,jms:jme) , INTENT(OUT) :: tslb
852
853 INTEGER :: i , j , l
854
855 ! Soil temperature is linearly interpolated between the skin temperature (taken to be at a
856 ! depth of 0.5 cm) and the deep soil, annual temperature (taken to be at a depth of 23 cm).
857 ! The tslb(i,1,j) is the skin temperature, and the tslb(i,num_soil_layers,j) level is the
858 ! annual mean temperature.
859
860 DO j = jts , MIN(jde-1,jte)
861 DO i = its , MIN(ide-1,ite)
862 IF ( landmask(i,j) .GT. 0.5 ) THEN
863 DO l = 1 , num_soil_layers
864 tslb(i,l,j)= ( tsk(i,j) * ( zs(num_soil_layers) - zs(l) ) + &
865 tmn(i,j) * ( zs( l) - zs(1) ) ) / &
866 ( zs(num_soil_layers) - zs(1) )
867 END DO
868 ELSE
869 IF ( ( real_data_init_type .EQ. 1 ) .AND. ( flag_sst .EQ. 1 ) ) THEN
870 DO l = 1 , num_soil_layers
871 tslb(i,l,j)= sst(i,j)
872 END DO
873 ELSE
874 DO l = 1 , num_soil_layers
875 tslb(i,l,j)= tsk(i,j)
876 END DO
877 END IF
878 END IF
879 END DO
880 END DO
881
882 END SUBROUTINE init_soil_1_real
883
884 SUBROUTINE init_soil_1_ideal(tsk,tmn,tslb,xland, &
885 ivgtyp,ZS,DZS,num_soil_layers, &
886 ids,ide, jds,jde, kds,kde, &
887 ims,ime, jms,jme, kms,kme, &
888 its,ite, jts,jte, kts,kte )
889
890 IMPLICIT NONE
891
892 INTEGER, INTENT(IN ) :: ids,ide, jds,jde, kds,kde, &
893 ims,ime, jms,jme, kms,kme, &
894 its,ite, jts,jte, kts,kte
895
896 INTEGER, INTENT(IN ) :: num_soil_layers
897
898 REAL, DIMENSION( ims:ime , 1 , jms:jme ), INTENT(OUT) :: tslb
899 REAL, DIMENSION( ims:ime , jms:jme ), INTENT(OUT) :: xland
900 INTEGER, DIMENSION( ims:ime , jms:jme ), INTENT(OUT) :: ivgtyp
901
902 REAL, DIMENSION(1:), INTENT(IN) :: dzs,zs
903
904 REAL, DIMENSION( ims:ime, jms:jme ) , INTENT(IN) :: tsk, tmn
905
906 ! Lcal variables.
907
908 INTEGER :: l,j,i,itf,jtf
909
910 itf=MIN(ite,ide-1)
911 jtf=MIN(jte,jde-1)
912
913 IF (num_soil_layers.NE.1)THEN
914 DO j=jts,jtf
915 DO l=1,num_soil_layers
916 DO i=its,itf
917 tslb(i,l,j)=( tsk(i,j)*(zs(num_soil_layers)-zs(l)) + tmn(i,j)*(zs(l)-zs(1)) ) / &
918 ( zs(num_soil_layers)-zs(1) )
919 ENDDO
920 ENDDO
921 ENDDO
922 ENDIF
923 DO j=jts,jtf
924 DO i=its,itf
925 xland(i,j) = 2
926 ivgtyp(i,j) = 7
927 ENDDO
928 ENDDO
929
930 END SUBROUTINE init_soil_1_ideal
931
932 SUBROUTINE init_soil_2_real ( tsk , tmn , smois , sh2o , tslb , &
933 st_input , sm_input , sw_input , landmask , sst , &
934 zs , dzs , &
935 st_levels_input , sm_levels_input , sw_levels_input , &
936 num_soil_layers , num_st_levels_input , num_sm_levels_input , num_sw_levels_input , &
937 num_st_levels_alloc , num_sm_levels_alloc , num_sw_levels_alloc , &
938 flag_sst , flag_soilt000 , flag_soilm000 , &
939 ids , ide , jds , jde , kds , kde , &
940 ims , ime , jms , jme , kms , kme , &
941 its , ite , jts , jte , kts , kte )
942
943 IMPLICIT NONE
944
945 INTEGER , INTENT(IN) :: num_soil_layers , &
946 num_st_levels_input , num_sm_levels_input , num_sw_levels_input , &
947 num_st_levels_alloc , num_sm_levels_alloc , num_sw_levels_alloc , &
948 ids , ide , jds , jde , kds , kde , &
949 ims , ime , jms , jme , kms , kme , &
950 its , ite , jts , jte , kts , kte
951
952 INTEGER , INTENT(IN) :: flag_sst, flag_soilt000, flag_soilm000
953
954 INTEGER , DIMENSION(1:num_st_levels_input) , INTENT(INOUT) :: st_levels_input
955 INTEGER , DIMENSION(1:num_sm_levels_input) , INTENT(INOUT) :: sm_levels_input
956 INTEGER , DIMENSION(1:num_sw_levels_input) , INTENT(INOUT) :: sw_levels_input
957
958 REAL , DIMENSION(ims:ime,1:num_st_levels_alloc,jms:jme) , INTENT(INOUT) :: st_input
959 REAL , DIMENSION(ims:ime,1:num_sm_levels_alloc,jms:jme) , INTENT(INOUT) :: sm_input
960 REAL , DIMENSION(ims:ime,1:num_sw_levels_alloc,jms:jme) , INTENT(INOUT) :: sw_input
961 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(IN) :: landmask , sst
962
963 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(IN) :: tmn
964 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(INOUT) :: tsk
965 REAL , DIMENSION(num_soil_layers) :: zs , dzs
966
967 REAL , DIMENSION(ims:ime,num_soil_layers,jms:jme) , INTENT(OUT) :: tslb , smois , sh2o
968
969 REAL , ALLOCATABLE , DIMENSION(:) :: zhave
970
971 INTEGER :: i , j , l , lout , lin , lwant , lhave , num
972 REAL :: temp
973 LOGICAL :: found_levels
974
975 CHARACTER (LEN=132) :: message
976
977 ! Are there any soil temp and moisture levels - ya know, they are mandatory.
978
979 num = num_st_levels_input * num_sm_levels_input
980
981 IF ( num .GE. 1 ) THEN
982
983 ! Ordered levels that we have data for.
984
985 !tgs add option to initialize from RUCLSM
986 IF ( flag_soilt000 .eq. 1 ) THEN
987 write(message, FMT='(A)') ' Assume RUC LSM 6-level input'
988 CALL wrf_message ( message )
989 ALLOCATE ( zhave( MAX(num_st_levels_input,num_sm_levels_input,num_sw_levels_input) ) )
990 ELSE
991 write(message, FMT='(A)') ' Assume Noah LSM input'
992 CALL wrf_message ( message )
993 ALLOCATE ( zhave( MAX(num_st_levels_input,num_sm_levels_input,num_sw_levels_input) +2) )
994 END IF
995
996
997 ! Sort the levels for temperature.
998
999 outert : DO lout = 1 , num_st_levels_input-1
1000 innert : DO lin = lout+1 , num_st_levels_input
1001 IF ( st_levels_input(lout) .GT. st_levels_input(lin) ) THEN
1002 temp = st_levels_input(lout)
1003 st_levels_input(lout) = st_levels_input(lin)
1004 st_levels_input(lin) = NINT(temp)
1005 DO j = jts , MIN(jde-1,jte)
1006 DO i = its , MIN(ide-1,ite)
1007 temp = st_input(i,lout+1,j)
1008 st_input(i,lout+1,j) = st_input(i,lin+1,j)
1009 st_input(i,lin+1,j) = temp
1010 END DO
1011 END DO
1012 END IF
1013 END DO innert
1014 END DO outert
1015 !tgs add IF
1016 IF ( flag_soilt000 .NE. 1 ) THEN
1017 DO j = jts , MIN(jde-1,jte)
1018 DO i = its , MIN(ide-1,ite)
1019 st_input(i,1,j) = tsk(i,j)
1020 st_input(i,num_st_levels_input+2,j) = tmn(i,j)
1021 END DO
1022 END DO
1023 ENDIF
1024
1025 ! Sort the levels for moisture.
1026
1027 outerm: DO lout = 1 , num_sm_levels_input-1
1028 innerm : DO lin = lout+1 , num_sm_levels_input
1029 IF ( sm_levels_input(lout) .GT. sm_levels_input(lin) ) THEN
1030 temp = sm_levels_input(lout)
1031 sm_levels_input(lout) = sm_levels_input(lin)
1032 sm_levels_input(lin) = NINT(temp)
1033 DO j = jts , MIN(jde-1,jte)
1034 DO i = its , MIN(ide-1,ite)
1035 temp = sm_input(i,lout+1,j)
1036 sm_input(i,lout+1,j) = sm_input(i,lin+1,j)
1037 sm_input(i,lin+1,j) = temp
1038 END DO
1039 END DO
1040 END IF
1041 END DO innerm
1042 END DO outerm
1043 !tgs add IF
1044 IF ( flag_soilm000 .NE. 1 ) THEN
1045 DO j = jts , MIN(jde-1,jte)
1046 DO i = its , MIN(ide-1,ite)
1047 sm_input(i,1,j) = sm_input(i,2,j)
1048 sm_input(i,num_sm_levels_input+2,j) = sm_input(i,num_sm_levels_input+1,j)
1049 END DO
1050 END DO
1051 ENDIF
1052
1053 ! Sort the levels for liquid moisture.
1054
1055 outerw: DO lout = 1 , num_sw_levels_input-1
1056 innerw : DO lin = lout+1 , num_sw_levels_input
1057 IF ( sw_levels_input(lout) .GT. sw_levels_input(lin) ) THEN
1058 temp = sw_levels_input(lout)
1059 sw_levels_input(lout) = sw_levels_input(lin)
1060 sw_levels_input(lin) = NINT(temp)
1061 DO j = jts , MIN(jde-1,jte)
1062 DO i = its , MIN(ide-1,ite)
1063 temp = sw_input(i,lout+1,j)
1064 sw_input(i,lout+1,j) = sw_input(i,lin+1,j)
1065 sw_input(i,lin+1,j) = temp
1066 END DO
1067 END DO
1068 END IF
1069 END DO innerw
1070 END DO outerw
1071 IF ( num_sw_levels_input .GT. 1 ) THEN
1072 DO j = jts , MIN(jde-1,jte)
1073 DO i = its , MIN(ide-1,ite)
1074 sw_input(i,1,j) = sw_input(i,2,j)
1075 sw_input(i,num_sw_levels_input+2,j) = sw_input(i,num_sw_levels_input+1,j)
1076 END DO
1077 END DO
1078 END IF
1079
1080 found_levels = .TRUE.
1081
1082 ELSE IF ( ( num .LE. 0 ) .AND. ( start_date .NE. current_date ) ) THEN
1083
1084 found_levels = .FALSE.
1085
1086 ELSE
1087 CALL wrf_error_fatal ( &
1088 'No input soil level data (temperature, moisture or liquid, or all are missing). Required for LSM.' )
1089 END IF
1090
1091 ! Is it OK to continue?
1092
1093 IF ( found_levels ) THEN
1094
1095 !tgs: Here are the levels that we have from the input for temperature.
1096
1097 IF ( flag_soilt000 .EQ. 1 ) THEN
1098 DO l = 1 , num_st_levels_input
1099 zhave(l) = st_levels_input(l) / 100.
1100 END DO
1101
1102 ! Interpolate between the layers we have (zhave) and those that we want (zs).
1103
1104 z_wantt : DO lwant = 1 , num_soil_layers
1105 z_havet : DO lhave = 1 , num_st_levels_input -1
1106 IF ( ( zs(lwant) .GE. zhave(lhave ) ) .AND. &
1107 ( zs(lwant) .LE. zhave(lhave+1) ) ) THEN
1108 DO j = jts , MIN(jde-1,jte)
1109 DO i = its , MIN(ide-1,ite)
1110 tslb(i,lwant,j)= ( st_input(i,lhave,j ) * ( zhave(lhave+1) - zs (lwant) ) + &
1111 st_input(i,lhave+1,j) * ( zs (lwant ) - zhave(lhave) ) ) / &
1112 ( zhave(lhave+1) - zhave(lhave) )
1113 END DO
1114 END DO
1115 EXIT z_havet
1116 END IF
1117 END DO z_havet
1118 END DO z_wantt
1119
1120 ELSE
1121
1122 ! Here are the levels that we have from the input for temperature. The input levels plus
1123 ! two more: the skin temperature at 0 cm, and the annual mean temperature at 300 cm.
1124
1125 zhave(1) = 0.
1126 DO l = 1 , num_st_levels_input
1127 zhave(l+1) = st_levels_input(l) / 100.
1128 END DO
1129 zhave(num_st_levels_input+2) = 300. / 100.
1130
1131 ! Interpolate between the layers we have (zhave) and those that we want (zs).
1132
1133 z_wantt_2: DO lwant = 1 , num_soil_layers
1134 z_havet_2 : DO lhave = 1 , num_st_levels_input +2 -1
1135 IF ( ( zs(lwant) .GE. zhave(lhave ) ) .AND. &
1136 ( zs(lwant) .LE. zhave(lhave+1) ) ) THEN
1137 DO j = jts , MIN(jde-1,jte)
1138 DO i = its , MIN(ide-1,ite)
1139 tslb(i,lwant,j)= ( st_input(i,lhave ,j) * ( zhave(lhave+1) - zs (lwant) ) + &
1140 st_input(i,lhave+1,j) * ( zs (lwant ) - zhave(lhave) ) ) / &
1141 ( zhave(lhave+1) - zhave(lhave) )
1142 END DO
1143 END DO
1144 EXIT z_havet_2
1145 END IF
1146 END DO z_havet_2
1147 END DO z_wantt_2
1148
1149 END IF
1150
1151 !tgs:
1152 IF ( flag_soilm000 .EQ. 1 ) THEN
1153 DO l = 1 , num_sm_levels_input
1154 zhave(l) = sm_levels_input(l) / 100.
1155 END DO
1156
1157 ! Interpolate between the layers we have (zhave) and those that we want (zs).
1158
1159 z_wantm : DO lwant = 1 , num_soil_layers
1160 z_havem : DO lhave = 1 , num_sm_levels_input -1
1161 IF ( ( zs(lwant) .GE. zhave(lhave ) ) .AND. &
1162 ( zs(lwant) .LE. zhave(lhave+1) ) ) THEN
1163 DO j = jts , MIN(jde-1,jte)
1164 DO i = its , MIN(ide-1,ite)
1165 smois(i,lwant,j)= ( sm_input(i,lhave,j ) * ( zhave(lhave+1) - zs (lwant) ) + &
1166 sm_input(i,lhave+1,j) * ( zs (lwant ) - zhave(lhave) ) ) / &
1167 ( zhave(lhave+1) - zhave(lhave) )
1168 END DO
1169 END DO
1170 EXIT z_havem
1171 END IF
1172 END DO z_havem
1173 END DO z_wantm
1174
1175 ELSE
1176 ! Here are the levels that we have from the input for moisture. The input levels plus
1177 ! two more: a value at 0 cm and one at 300 cm. The 0 cm value is taken to be identical
1178 ! to the most shallow layer's value. Similarly, the 300 cm value is taken to be the same
1179 ! as the most deep layer's value.
1180
1181 zhave(1) = 0.
1182 DO l = 1 , num_sm_levels_input
1183 zhave(l+1) = sm_levels_input(l) / 100.
1184 END DO
1185 zhave(num_sm_levels_input+2) = 300. / 100.
1186
1187 ! Interpolate between the layers we have (zhave) and those that we want (zs).
1188
1189 z_wantm_2 : DO lwant = 1 , num_soil_layers
1190 z_havem_2 : DO lhave = 1 , num_sm_levels_input +2 -1
1191 IF ( ( zs(lwant) .GE. zhave(lhave ) ) .AND. &
1192 ( zs(lwant) .LE. zhave(lhave+1) ) ) THEN
1193 DO j = jts , MIN(jde-1,jte)
1194 DO i = its , MIN(ide-1,ite)
1195 smois(i,lwant,j)= ( sm_input(i,lhave ,j) * ( zhave(lhave+1) - zs (lwant) ) + &
1196 sm_input(i,lhave+1,j) * ( zs (lwant ) - zhave(lhave) ) ) / &
1197 ( zhave(lhave+1) - zhave(lhave) )
1198 END DO
1199 END DO
1200 EXIT z_havem_2
1201 END IF
1202 END DO z_havem_2
1203 END DO z_wantm_2
1204 ENDIF
1205
1206 ! Any liquid soil moisture to worry about?
1207
1208 IF ( num_sw_levels_input .GT. 1 ) THEN
1209
1210 zhave(1) = 0.
1211 DO l = 1 , num_sw_levels_input
1212 zhave(l+1) = sw_levels_input(l) / 100.
1213 END DO
1214 zhave(num_sw_levels_input+2) = 300. / 100.
1215
1216 ! Interpolate between the layers we have (zhave) and those that we want (zs).
1217
1218 z_wantw : DO lwant = 1 , num_soil_layers
1219 z_havew : DO lhave = 1 , num_sw_levels_input +2 -1
1220 IF ( ( zs(lwant) .GE. zhave(lhave ) ) .AND. &
1221 ( zs(lwant) .LE. zhave(lhave+1) ) ) THEN
1222 DO j = jts , MIN(jde-1,jte)
1223 DO i = its , MIN(ide-1,ite)
1224 sh2o(i,lwant,j)= ( sw_input(i,lhave ,j) * ( zhave(lhave+1) - zs (lwant) ) + &
1225 sw_input(i,lhave+1,j) * ( zs (lwant ) - zhave(lhave) ) ) / &
1226 ( zhave(lhave+1) - zhave(lhave) )
1227 END DO
1228 END DO
1229 EXIT z_havew
1230 END IF
1231 END DO z_havew
1232 END DO z_wantw
1233
1234 END IF
1235
1236
1237 ! Over water, put in reasonable values for soil temperature and moisture. These won't be
1238 ! used, but they will make a more continuous plot.
1239
1240 IF ( flag_sst .EQ. 1 ) THEN
1241 DO j = jts , MIN(jde-1,jte)
1242 DO i = its , MIN(ide-1,ite)
1243 IF ( landmask(i,j) .LT. 0.5 ) THEN
1244 DO l = 1 , num_soil_layers
1245 tslb(i,l,j)= sst(i,j)
1246 smois(i,l,j)= 1.0
1247 sh2o (i,l,j)= 1.0
1248 END DO
1249 END IF
1250 END DO
1251 END DO
1252 ELSE
1253 DO j = jts , MIN(jde-1,jte)
1254 DO i = its , MIN(ide-1,ite)
1255 IF ( landmask(i,j) .LT. 0.5 ) THEN
1256 DO l = 1 , num_soil_layers
1257 tslb(i,l,j)= tsk(i,j)
1258 smois(i,l,j)= 1.0
1259 sh2o (i,l,j)= 1.0
1260 END DO
1261 END IF
1262 END DO
1263 END DO
1264 END IF
1265
1266 DEALLOCATE (zhave)
1267
1268 END IF
1269
1270 END SUBROUTINE init_soil_2_real
1271
1272 SUBROUTINE init_soil_2_ideal ( xland,xice,vegfra,snow,canwat, &
1273 ivgtyp,isltyp,tslb,smois,tmn, &
1274 num_soil_layers, &
1275 ids,ide, jds,jde, kds,kde, &
1276 ims,ime, jms,jme, kms,kme, &
1277 its,ite, jts,jte, kts,kte )
1278
1279 IMPLICIT NONE
1280
1281 INTEGER, INTENT(IN) ::ids,ide, jds,jde, kds,kde, &
1282 ims,ime, jms,jme, kms,kme, &
1283 its,ite, jts,jte, kts,kte
1284
1285 INTEGER, INTENT(IN) ::num_soil_layers
1286
1287 REAL, DIMENSION( ims:ime, num_soil_layers, jms:jme ) , INTENT(OUT) :: smois, tslb
1288
1289 REAL, DIMENSION( ims:ime, jms:jme ) , INTENT(OUT) :: xland, snow, canwat, xice, vegfra, tmn
1290
1291 INTEGER, DIMENSION( ims:ime, jms:jme ) , INTENT(OUT) :: ivgtyp, isltyp
1292
1293 INTEGER :: icm,jcm,itf,jtf
1294 INTEGER :: i,j,l
1295
1296 itf=min0(ite,ide-1)
1297 jtf=min0(jte,jde-1)
1298
1299 icm = ide/2
1300 jcm = jde/2
1301
1302 DO j=jts,jtf
1303 DO l=1,num_soil_layers
1304 DO i=its,itf
1305
1306 smois(i,1,j)=0.10
1307 smois(i,2,j)=0.10
1308 smois(i,3,j)=0.10
1309 smois(i,4,j)=0.10
1310
1311 tslb(i,1,j)=295.
1312 tslb(i,2,j)=297.
1313 tslb(i,3,j)=293.
1314 tslb(i,4,j)=293.
1315
1316 ENDDO
1317 ENDDO
1318 ENDDO
1319
1320 DO j=jts,jtf
1321 DO i=its,itf
1322 xland(i,j) = 2
1323 tmn(i,j) = 294.
1324 xice(i,j) = 0.
1325 vegfra(i,j) = 0.
1326 snow(i,j) = 0.
1327 canwat(i,j) = 0.
1328 ivgtyp(i,j) = 7
1329 isltyp(i,j) = 8
1330 ENDDO
1331 ENDDO
1332
1333 END SUBROUTINE init_soil_2_ideal
1334
1335 SUBROUTINE init_soil_3_real ( tsk , tmn , smois , tslb , &
1336 st_input , sm_input , landmask, sst, &
1337 zs , dzs , &
1338 st_levels_input , sm_levels_input , &
1339 num_soil_layers , num_st_levels_input , num_sm_levels_input , &
1340 num_st_levels_alloc , num_sm_levels_alloc , &
1341 flag_sst , flag_soilt000 , flag_soilm000 , &
1342 ids , ide , jds , jde , kds , kde , &
1343 ims , ime , jms , jme , kms , kme , &
1344 its , ite , jts , jte , kts , kte )
1345
1346 IMPLICIT NONE
1347
1348 INTEGER , INTENT(IN) :: num_soil_layers , &
1349 num_st_levels_input , num_sm_levels_input , &
1350 num_st_levels_alloc , num_sm_levels_alloc , &
1351 ids , ide , jds , jde , kds , kde , &
1352 ims , ime , jms , jme , kms , kme , &
1353 its , ite , jts , jte , kts , kte
1354
1355 INTEGER , INTENT(IN) :: flag_sst, flag_soilt000, flag_soilm000
1356
1357 INTEGER , DIMENSION(1:num_st_levels_input) , INTENT(INOUT) :: st_levels_input
1358 INTEGER , DIMENSION(1:num_sm_levels_input) , INTENT(INOUT) :: sm_levels_input
1359
1360 REAL , DIMENSION(ims:ime,1:num_st_levels_alloc,jms:jme) , INTENT(INOUT) :: st_input
1361 REAL , DIMENSION(ims:ime,1:num_sm_levels_alloc,jms:jme) , INTENT(INOUT) :: sm_input
1362 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(IN) :: landmask , sst
1363
1364 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(IN) :: tmn
1365 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(INOUT) :: tsk
1366 REAL , DIMENSION(num_soil_layers) :: zs , dzs
1367
1368 REAL , DIMENSION(ims:ime,num_soil_layers,jms:jme) , INTENT(OUT) :: tslb , smois
1369
1370 REAL , ALLOCATABLE , DIMENSION(:) :: zhave
1371
1372 INTEGER :: i , j , l , lout , lin , lwant , lhave, k
1373 REAL :: temp
1374
1375 CHARACTER (LEN=132) :: message
1376
1377 ! Allocate the soil layer array used for interpolating.
1378
1379 IF ( ( num_st_levels_input .LE. 0 ) .OR. &
1380 ( num_sm_levels_input .LE. 0 ) ) THEN
1381 write (message, FMT='(A)')&
1382 'No input soil level data (either temperature or moisture, or both are missing). Required for RUC LSM.'
1383 CALL wrf_error_fatal ( message )
1384 ELSE
1385 IF ( flag_soilt000 .eq. 1 ) THEN
1386 write(message, FMT='(A)') ' Assume RUC LSM 6-level input'
1387 CALL wrf_message ( message )
1388 ALLOCATE ( zhave( MAX(num_st_levels_input,num_sm_levels_input) ) )
1389 ELSE
1390 write(message, FMT='(A)') ' Assume non-RUC LSM input'
1391 CALL wrf_message ( message )
1392 ALLOCATE ( zhave( MAX(num_st_levels_input,num_soil_layers) ) )
1393 END IF
1394 END IF
1395
1396 ! Sort the levels for temperature.
1397
1398 outert : DO lout = 1 , num_st_levels_input-1
1399 innert : DO lin = lout+1 , num_st_levels_input
1400 IF ( st_levels_input(lout) .GT. st_levels_input(lin) ) THEN
1401 temp = st_levels_input(lout)
1402 st_levels_input(lout) = st_levels_input(lin)
1403 st_levels_input(lin) = NINT(temp)
1404 DO j = jts , MIN(jde-1,jte)
1405 DO i = its , MIN(ide-1,ite)
1406 temp = st_input(i,lout,j)
1407 st_input(i,lout,j) = st_input(i,lin,j)
1408 st_input(i,lin,j) = temp
1409 END DO
1410 END DO
1411 END IF
1412 END DO innert
1413 END DO outert
1414
1415 IF ( flag_soilt000 .NE. 1 ) THEN
1416 DO j = jts , MIN(jde-1,jte)
1417 DO i = its , MIN(ide-1,ite)
1418 st_input(i,1,j) = tsk(i,j)
1419 st_input(i,num_st_levels_input+2,j) = tmn(i,j)
1420 END DO
1421 END DO
1422 END IF
1423
1424 ! Sort the levels for moisture.
1425
1426 outerm: DO lout = 1 , num_sm_levels_input-1
1427 innerm : DO lin = lout+1 , num_sm_levels_input
1428 IF ( sm_levels_input(lout) .GT. sm_levels_input(lin) ) THEN
1429 temp = sm_levels_input(lout)
1430 sm_levels_input(lout) = sm_levels_input(lin)
1431 sm_levels_input(lin) = NINT(temp)
1432 DO j = jts , MIN(jde-1,jte)
1433 DO i = its , MIN(ide-1,ite)
1434 temp = sm_input(i,lout,j)
1435 sm_input(i,lout,j) = sm_input(i,lin,j)
1436 sm_input(i,lin,j) = temp
1437 END DO
1438 END DO
1439 END IF
1440 END DO innerm
1441 END DO outerm
1442
1443 IF ( flag_soilm000 .NE. 1 ) THEN
1444 DO j = jts , MIN(jde-1,jte)
1445 DO i = its , MIN(ide-1,ite)
1446 sm_input(i,1,j) = (sm_input(i,2,j)-sm_input(i,3,j))/ &
1447 (st_levels_input(2)-st_levels_input(1))*st_levels_input(1)+ &
1448 sm_input(i,2,j)
1449 ! sm_input(i,1,j) = sm_input(i,2,j)
1450 sm_input(i,num_sm_levels_input+2,j) = sm_input(i,num_sm_levels_input+1,j)
1451 END DO
1452 END DO
1453 END IF
1454
1455 ! Here are the levels that we have from the input for temperature.
1456
1457 IF ( flag_soilt000 .EQ. 1 ) THEN
1458 DO l = 1 , num_st_levels_input
1459 zhave(l) = st_levels_input(l) / 100.
1460 END DO
1461
1462 ! Interpolate between the layers we have (zhave) and those that we want (zs).
1463
1464
1465 z_wantt : DO lwant = 1 , num_soil_layers
1466 z_havet : DO lhave = 1 , num_st_levels_input -1
1467 IF ( ( zs(lwant) .GE. zhave(lhave ) ) .AND. &
1468 ( zs(lwant) .LE. zhave(lhave+1) ) ) THEN
1469 DO j = jts , MIN(jde-1,jte)
1470 DO i = its , MIN(ide-1,ite)
1471 tslb(i,lwant,j)= ( st_input(i,lhave,j ) * ( zhave(lhave+1) - zs (lwant) ) + &
1472 st_input(i,lhave+1,j) * ( zs (lwant ) - zhave(lhave) ) ) / &
1473 ( zhave(lhave+1) - zhave(lhave) )
1474 END DO
1475 END DO
1476 EXIT z_havet
1477 END IF
1478 END DO z_havet
1479 END DO z_wantt
1480
1481 ELSE
1482
1483 zhave(1) = 0.
1484 DO l = 1 , num_st_levels_input
1485 zhave(l+1) = st_levels_input(l) / 100.
1486 END DO
1487 zhave(num_st_levels_input+2) = 300. / 100.
1488
1489 ! Interpolate between the layers we have (zhave) and those that we want (zs).
1490
1491 z_wantt_2 : DO lwant = 1 , num_soil_layers
1492 z_havet_2 : DO lhave = 1 , num_st_levels_input +2
1493 IF ( ( zs(lwant) .GE. zhave(lhave ) ) .AND. &
1494 ( zs(lwant) .LE. zhave(lhave+1) ) ) THEN
1495 DO j = jts , MIN(jde-1,jte)
1496 DO i = its , MIN(ide-1,ite)
1497 tslb(i,lwant,j)= ( st_input(i,lhave,j ) * ( zhave(lhave+1) - zs (lwant) ) + &
1498 st_input(i,lhave+1,j) * ( zs (lwant ) - zhave(lhave) ) ) / &
1499 ( zhave(lhave+1) - zhave(lhave) )
1500 END DO
1501 END DO
1502 EXIT z_havet_2
1503 END IF
1504 END DO z_havet_2
1505 END DO z_wantt_2
1506
1507 END IF
1508
1509 ! Here are the levels that we have from the input for moisture.
1510
1511 IF ( flag_soilm000 .EQ. 1 ) THEN
1512 DO l = 1 , num_sm_levels_input
1513 zhave(l) = sm_levels_input(l) / 100.
1514 END DO
1515
1516 ! Interpolate between the layers we have (zhave) and those that we want (zs).
1517
1518 z_wantm : DO lwant = 1 , num_soil_layers
1519 z_havem : DO lhave = 1 , num_sm_levels_input -1
1520 IF ( ( zs(lwant) .GE. zhave(lhave ) ) .AND. &
1521 ( zs(lwant) .LE. zhave(lhave+1) ) ) THEN
1522 DO j = jts , MIN(jde-1,jte)
1523 DO i = its , MIN(ide-1,ite)
1524 smois(i,lwant,j)= ( sm_input(i,lhave,j ) * ( zhave(lhave+1) - zs (lwant) ) + &
1525 sm_input(i,lhave+1,j) * ( zs (lwant ) - zhave(lhave) ) ) / &
1526 ( zhave(lhave+1) - zhave(lhave) )
1527 END DO
1528 END DO
1529 EXIT z_havem
1530 END IF
1531 END DO z_havem
1532 END DO z_wantm
1533
1534 ELSE
1535
1536 zhave(1) = 0.
1537 DO l = 1 , num_sm_levels_input
1538 zhave(l+1) = sm_levels_input(l) / 100.
1539 END DO
1540 zhave(num_sm_levels_input+2) = 300. / 100.
1541
1542 z_wantm_2 : DO lwant = 1 , num_soil_layers
1543 z_havem_2 : DO lhave = 1 , num_sm_levels_input +2
1544 IF ( ( zs(lwant) .GE. zhave(lhave ) ) .AND. &
1545 ( zs(lwant) .LE. zhave(lhave+1) ) ) THEN
1546 DO j = jts , MIN(jde-1,jte)
1547 DO i = its , MIN(ide-1,ite)
1548 smois(i,lwant,j)= ( sm_input(i,lhave,j ) * ( zhave(lhave+1) - zs (lwant) ) + &
1549 sm_input(i,lhave+1,j) * ( zs (lwant ) - zhave(lhave) ) ) / &
1550 ( zhave(lhave+1) - zhave(lhave) )
1551 END DO
1552 END DO
1553 EXIT z_havem_2
1554 END IF
1555 END DO z_havem_2
1556 END DO z_wantm_2
1557
1558 END IF
1559 ! Over water, put in reasonable values for soil temperature and moisture. These won't be
1560 ! used, but they will make a more continuous plot.
1561
1562 IF ( flag_sst .EQ. 1 ) THEN
1563 DO j = jts , MIN(jde-1,jte)
1564 DO i = its , MIN(ide-1,ite)
1565 IF ( landmask(i,j) .LT. 0.5 ) THEN
1566 DO l = 1 , num_soil_layers
1567 tslb(i,l,j) = sst(i,j)
1568 tsk(i,j) = sst(i,j)
1569 smois(i,l,j)= 1.0
1570 END DO
1571 END IF
1572 END DO
1573 END DO
1574 ELSE
1575 DO j = jts , MIN(jde-1,jte)
1576 DO i = its , MIN(ide-1,ite)
1577 IF ( landmask(i,j) .LT. 0.5 ) THEN
1578 DO l = 1 , num_soil_layers
1579 tslb(i,l,j)= tsk(i,j)
1580 smois(i,l,j)= 1.0
1581 END DO
1582 END IF
1583 END DO
1584 END DO
1585 END IF
1586
1587 DEALLOCATE (zhave)
1588
1589 END SUBROUTINE init_soil_3_real
1590
1591 END MODULE module_soil_pre
1592
1593 #else
1594
1595 MODULE module_soil_pre
1596
1597 USE module_date_time
1598 USE module_state_description
1599
1600 CONTAINS
1601
1602 SUBROUTINE adjust_for_seaice_pre ( xice , landmask , tsk , ivgtyp , vegcat , lu_index , &
1603 xland , landusef , isltyp , soilcat , soilctop , &
1604 soilcbot , tmn , &
1605 seaice_threshold , &
1606 num_veg_cat , num_soil_top_cat , num_soil_bot_cat , &
1607 iswater , isice , &
1608 scheme , &
1609 ids , ide , jds , jde , kds , kde , &
1610 ims , ime , jms , jme , kms , kme , &
1611 its , ite , jts , jte , kts , kte )
1612
1613 IMPLICIT NONE
1614
1615 INTEGER , INTENT(IN) :: ids , ide , jds , jde , kds , kde , &
1616 ims , ime , jms , jme , kms , kme , &
1617 its , ite , jts , jte , kts , kte , &
1618 iswater , isice
1619 INTEGER , INTENT(IN) :: num_veg_cat , num_soil_top_cat , num_soil_bot_cat , scheme
1620
1621 REAL , DIMENSION(ims:ime,1:num_veg_cat,jms:jme) , INTENT(INOUT):: landusef
1622 REAL , DIMENSION(ims:ime,1:num_soil_top_cat,jms:jme) , INTENT(INOUT):: soilctop
1623 REAL , DIMENSION(ims:ime,1:num_soil_bot_cat,jms:jme) , INTENT(INOUT):: soilcbot
1624 INTEGER , DIMENSION(ims:ime,jms:jme), INTENT(OUT) :: isltyp , ivgtyp
1625 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(INOUT) :: landmask , xice , tsk , lu_index , &
1626 vegcat, xland , soilcat , tmn
1627 REAL , INTENT(IN) :: seaice_threshold
1628
1629 INTEGER :: i , j , num_seaice_changes , loop
1630 CHARACTER (LEN=132) :: message
1631
1632 fix_seaice : SELECT CASE ( scheme )
1633
1634 CASE ( SLABSCHEME )
1635 DO j = jts , MIN(jde-1,jte)
1636 DO i = its , MIN(ide-1,ite)
1637 IF ( xice(i,j) .GT. 200.0 ) THEN
1638 xice(i,j) = 0.
1639 num_seaice_changes = num_seaice_changes + 1
1640 END IF
1641 END DO
1642 END DO
1643 IF ( num_seaice_changes .GT. 0 ) THEN
1644 WRITE ( message , FMT='(A,I6)' ) &
1645 'Total pre number of sea ice locations removed (due to FLAG values) = ', &
1646 num_seaice_changes
1647 CALL wrf_debug ( 0 , message )
1648 END IF
1649 num_seaice_changes = 0
1650 DO j = jts , MIN(jde-1,jte)
1651 DO i = its , MIN(ide-1,ite)
1652 IF ( ( xice(i,j) .GE. 0.5 ) .OR. &
1653 ( ( landmask(i,j) .LT. 0.5 ) .AND. ( tsk(i,j) .LT. seaice_threshold ) ) ) THEN
1654 xice(i,j) = 1.
1655 num_seaice_changes = num_seaice_changes + 1
1656 tmn(i,j) = 271.4
1657 vegcat(i,j)=isice
1658 lu_index(i,j)=ivgtyp(i,j)
1659 landmask(i,j)=1.
1660 xland(i,j)=1.
1661 DO loop=1,num_veg_cat
1662 landusef(i,loop,j)=0.
1663 END DO
1664 landusef(i,ivgtyp(i,j),j)=1.
1665
1666 isltyp(i,j) = 16
1667 soilcat(i,j)=isltyp(i,j)
1668 DO loop=1,num_soil_top_cat
1669 soilctop(i,loop,j)=0
1670 END DO
1671 DO loop=1,num_soil_bot_cat
1672 soilcbot(i,loop,j)=0
1673 END DO
1674 soilctop(i,isltyp(i,j),j)=1.
1675 soilcbot(i,isltyp(i,j),j)=1.
1676 END IF
1677 END DO
1678 END DO
1679 IF ( num_seaice_changes .GT. 0 ) THEN
1680 WRITE ( message , FMT='(A,I6)' ) &
1681 'Total pre number of sea ice location changes (water to land) = ', num_seaice_changes
1682 CALL wrf_debug ( 0 , message )
1683 END IF
1684
1685 CASE ( LSMSCHEME )
1686 CASE ( RUCLSMSCHEME )
1687
1688 END SELECT fix_seaice
1689
1690 END SUBROUTINE adjust_for_seaice_pre
1691
1692 SUBROUTINE adjust_for_seaice_post ( xice , landmask , tsk , ivgtyp , vegcat , lu_index , &
1693 xland , landusef , isltyp , soilcat , soilctop , &
1694 soilcbot , tmn , &
1695 tslb , smois , sh2o , &
1696 seaice_threshold , &
1697 num_veg_cat , num_soil_top_cat , num_soil_bot_cat , &
1698 num_soil_layers , &
1699 iswater , isice , &
1700 scheme , &
1701 ids , ide , jds , jde , kds , kde , &
1702 ims , ime , jms , jme , kms , kme , &
1703 its , ite , jts , jte , kts , kte )
1704
1705 IMPLICIT NONE
1706
1707 INTEGER , INTENT(IN) :: ids , ide , jds , jde , kds , kde , &
1708 ims , ime , jms , jme , kms , kme , &
1709 its , ite , jts , jte , kts , kte , &
1710 iswater , isice
1711 INTEGER , INTENT(IN) :: num_veg_cat , num_soil_top_cat , num_soil_bot_cat , scheme
1712 INTEGER , INTENT(IN) :: num_soil_layers
1713
1714 REAL , DIMENSION(ims:ime,1:num_veg_cat,jms:jme) , INTENT(INOUT):: landusef
1715 REAL , DIMENSION(ims:ime,1:num_soil_top_cat,jms:jme) , INTENT(INOUT):: soilctop
1716 REAL , DIMENSION(ims:ime,1:num_soil_bot_cat,jms:jme) , INTENT(INOUT):: soilcbot
1717 REAL , DIMENSION(ims:ime,1:num_soil_layers,jms:jme) , INTENT(INOUT):: tslb , smois , sh2o
1718 INTEGER , DIMENSION(ims:ime,jms:jme), INTENT(OUT) :: isltyp , ivgtyp
1719 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(INOUT) :: landmask , xice , tsk , lu_index , &
1720 vegcat, xland , soilcat , tmn
1721 REAL , INTENT(IN) :: seaice_threshold
1722 REAL :: total_depth , mid_point_depth
1723
1724 INTEGER :: i , j , num_seaice_changes , loop
1725 CHARACTER (LEN=132) :: message
1726
1727 fix_seaice : SELECT CASE ( scheme )
1728
1729 CASE ( SLABSCHEME )
1730
1731 CASE ( LSMSCHEME )
1732 DO j = jts , MIN(jde-1,jte)
1733 DO i = its , MIN(ide-1,ite)
1734 IF ( xice(i,j) .GT. 200.0 ) THEN
1735 xice(i,j) = 0.
1736 num_seaice_changes = num_seaice_changes + 1
1737 END IF
1738 END DO
1739 END DO
1740 IF ( num_seaice_changes .GT. 0 ) THEN
1741 WRITE ( message , FMT='(A,I6)' ) &
1742 'Total post number of sea ice locations removed (due to FLAG values) = ', &
1743 num_seaice_changes
1744 CALL wrf_debug ( 0 , message )
1745 END IF
1746 num_seaice_changes = 0
1747 DO j = jts , MIN(jde-1,jte)
1748 DO i = its , MIN(ide-1,ite)
1749 IF ( ( xice(i,j) .GE. 0.5 ) .OR. &
1750 ( ( landmask(i,j) .LT. 0.5 ) .AND. ( tsk(i,j) .LT. seaice_threshold ) ) ) THEN
1751 xice(i,j) = 1.
1752 num_seaice_changes = num_seaice_changes + 1
1753 tmn(i,j) = 271.16
1754 vegcat(i,j)=isice
1755 lu_index(i,j)=ivgtyp(i,j)
1756 landmask(i,j)=1.
1757 xland(i,j)=1.
1758 DO loop=1,num_veg_cat
1759 landusef(i,loop,j)=0.
1760 END DO
1761 landusef(i,ivgtyp(i,j),j)=1.
1762
1763 isltyp(i,j) = 16
1764 soilcat(i,j)=isltyp(i,j)
1765 DO loop=1,num_soil_top_cat
1766 soilctop(i,loop,j)=0
1767 END DO
1768 DO loop=1,num_soil_bot_cat
1769 soilcbot(i,loop,j)=0
1770 END DO
1771 soilctop(i,isltyp(i,j),j)=1.
1772 soilcbot(i,isltyp(i,j),j)=1.
1773
1774 total_depth = 3. ! ice is 3 m deep, num_soil_layers equispaced layers
1775 DO loop = 1,num_soil_layers
1776 mid_point_depth=(total_depth/num_soil_layers)/2. + &
1777 (loop-1)*(total_depth/num_soil_layers)
1778 tslb(i,loop,j) = ( (total_depth-mid_point_depth)*tsk(i,j) + &
1779 mid_point_depth*tmn(i,j) ) / total_depth
1780 END DO
1781
1782 DO loop=1,num_soil_layers
1783 smois(i,loop,j) = 1.0
1784 sh2o(i,loop,j) = 0.0
1785 END DO
1786 END IF
1787 END DO
1788 END DO
1789 IF ( num_seaice_changes .GT. 0 ) THEN
1790 WRITE ( message , FMT='(A,I6)' ) &
1791 'Total post number of sea ice location changes (water to land) = ', num_seaice_changes
1792 CALL wrf_debug ( 0 , message )
1793 END IF
1794
1795 CASE ( RUCLSMSCHEME )
1796
1797 END SELECT fix_seaice
1798
1799 END SUBROUTINE adjust_for_seaice_post
1800
1801 SUBROUTINE process_percent_cat_new ( landmask , &
1802 landuse_frac , soil_top_cat , soil_bot_cat , &
1803 isltyp , ivgtyp , &
1804 num_veg_cat , num_soil_top_cat , num_soil_bot_cat , &
1805 ids , ide , jds , jde , kds , kde , &
1806 ims , ime , jms , jme , kms , kme , &
1807 its , ite , jts , jte , kts , kte , &
1808 iswater )
1809
1810 IMPLICIT NONE
1811
1812 INTEGER , INTENT(IN) :: ids , ide , jds , jde , kds , kde , &
1813 ims , ime , jms , jme , kms , kme , &
1814 its , ite , jts , jte , kts , kte , &
1815 iswater
1816 INTEGER , INTENT(IN) :: num_veg_cat , num_soil_top_cat , num_soil_bot_cat
1817 REAL , DIMENSION(ims:ime,1:num_veg_cat,jms:jme) , INTENT(INOUT):: landuse_frac
1818 REAL , DIMENSION(ims:ime,1:num_soil_top_cat,jms:jme) , INTENT(IN):: soil_top_cat
1819 REAL , DIMENSION(ims:ime,1:num_soil_bot_cat,jms:jme) , INTENT(IN):: soil_bot_cat
1820 INTEGER , DIMENSION(ims:ime,jms:jme), INTENT(OUT) :: isltyp , ivgtyp
1821 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(INOUT) :: landmask
1822
1823 INTEGER :: i , j , l , ll, dominant_index
1824 REAL :: dominant_value
1825
1826 #ifdef WRF_CHEM
1827 ! REAL :: lwthresh = .99
1828 REAL :: lwthresh = .50
1829 #else
1830 REAL :: lwthresh = .50
1831 #endif
1832
1833 INTEGER , PARAMETER :: iswater_soil = 14
1834 INTEGER :: iforce
1835 CHARACTER (LEN=1024) :: message
1836
1837 ! Sanity check on the 50/50 points
1838
1839 DO j = jts , MIN(jde-1,jte)
1840 DO i = its , MIN(ide-1,ite)
1841 dominant_value = landuse_frac(i,iswater,j)
1842 IF ( dominant_value .EQ. lwthresh ) THEN
1843 DO l = 1 , num_veg_cat
1844 IF ( l .EQ. iswater ) CYCLE
1845 IF ( landuse_frac(i,l,j) .EQ. lwthresh ) THEN
1846 write(message, FMT='(I4,I4,A,I4,A,F12.4)') i,j,' water and category ',l,' both at 50%, landmask is ',landmask(i,j)
1847 call wrf_message ( message )
1848 landuse_frac(i,l,j) = lwthresh - .01
1849 landuse_frac(i,iswater,j) = 1. - (lwthresh + 0.01)
1850 END IF
1851 END DO
1852 END IF
1853 END DO
1854 END DO
1855
1856 ! Compute the dominant VEGETATION INDEX.
1857
1858 DO j = jts , MIN(jde-1,jte)
1859 DO i = its , MIN(ide-1,ite)
1860 dominant_value = landuse_frac(i,1,j)
1861 dominant_index = 1
1862 DO l = 2 , num_veg_cat
1863 IF ( l .EQ. iswater ) THEN
1864 ! wait a bit
1865 ELSE IF ( ( l .NE. iswater ) .AND. ( landuse_frac(i,l,j) .GT. dominant_value ) ) THEN
1866 dominant_value = landuse_frac(i,l,j)
1867 dominant_index = l
1868 END IF
1869 END DO
1870 IF ( landuse_frac(i,iswater,j) .GT. lwthresh ) THEN
1871 dominant_value = landuse_frac(i,iswater,j)
1872 dominant_index = iswater
1873 #if 0
1874 si needs to beef up consistency checks before we can use this part
1875 ELSE IF ( ( landuse_frac(i,iswater,j) .EQ. lwthresh) .AND. &
1876 ( dominant_value .EQ. lwthresh) ) THEN
1877 ! no op
1878 #else
1879 ELSEIF((landuse_frac(i,iswater,j).EQ.lwthresh).AND.(dominant_value.EQ.lwthresh).and.(dominant_index.LT.iswater))THEN
1880 write(message,*)'temporary SI landmask fix'
1881 call wrf_message(trim(message))
1882 ! no op
1883 ELSEIF((landuse_frac(i,iswater,j).EQ.lwthresh).AND.(dominant_value.EQ.lwthresh).and.(dominant_index.GT.iswater))THEN
1884 write(message,*)'temporary SI landmask fix'
1885 call wrf_message(trim(message))
1886 dominant_value = landuse_frac(i,iswater,j)
1887 dominant_index = iswater
1888 #endif
1889 ELSE IF ( ( landuse_frac(i,iswater,j) .EQ. lwthresh ) .AND. &
1890 ( dominant_value .LT. lwthresh ) ) THEN
1891 dominant_value = landuse_frac(i,iswater,j)
1892 dominant_index = iswater
1893 END IF
1894 IF ( dominant_index .EQ. iswater ) THEN
1895 if(landmask(i,j).gt.lwthresh) then
1896 write(message,*) 'changing to water at point ',i,j
1897 call wrf_message(trim(message))
1898 write(message,*) nint(landuse_frac(i,:,j)*100)
1899 call wrf_message(trim(message))
1900 endif
1901 landmask(i,j) = 0
1902 ELSE IF ( dominant_index .NE. iswater ) THEN
1903 if(landmask(i,j).lt.lwthresh) then
1904 write(message,*) 'changing to land at point ',i,j
1905 call wrf_message(trim(message))
1906 write(message,*) nint(landuse_frac(i,:,j)*100)
1907 call wrf_message(trim(message))
1908 endif
1909 landmask(i,j) = 1
1910 END IF
1911 ivgtyp(i,j) = dominant_index
1912 END DO
1913 END DO
1914
1915 ! Compute the dominant SOIL TEXTURE INDEX, TOP.
1916
1917 iforce = 0
1918 DO i = its , MIN(ide-1,ite)
1919 DO j = jts , MIN(jde-1,jte)
1920 dominant_value = soil_top_cat(i,1,j)
1921 dominant_index = 1
1922 IF ( landmask(i,j) .GT. lwthresh ) THEN
1923 DO l = 2 , num_soil_top_cat
1924 IF ( ( l .NE. iswater_soil ) .AND. ( soil_top_cat(i,l,j) .GT. dominant_value ) ) THEN
1925 dominant_value = soil_top_cat(i,l,j)
1926 dominant_index = l
1927 END IF
1928 END DO
1929 IF ( dominant_value .LT. 0.01 ) THEN
1930 iforce = iforce + 1
1931 WRITE ( message , FMT = '(A,I4,I4)' ) &
1932 'based on landuse, changing soil to land at point ',i,j
1933 CALL wrf_debug(1,message)
1934 WRITE ( message , FMT = '(16(i3,1x))' ) &
1935 1, 2, 3, 4, 5, 6, 7, 8, 9, 10,11,12, 13, 14, 15, 16
1936 CALL wrf_debug(1,message)
1937 WRITE ( message , FMT = '(16(i3,1x))' ) &
1938 nint(soil_top_cat(i,:,j)*100)
1939 CALL wrf_debug(1,message)
1940 dominant_index = 8
1941 END IF
1942 ELSE
1943 dominant_index = iswater_soil
1944 END IF
1945 isltyp(i,j) = dominant_index
1946 END DO
1947 END DO
1948
1949 if(iforce.ne.0)then
1950 WRITE(message,FMT='(A,I4,A,I6)' ) &
1951 'forcing artificial silty clay loam at ',iforce,' points, out of ',&
1952 (MIN(ide-1,ite)-its+1)*(MIN(jde-1,jte)-jts+1)
1953 CALL wrf_debug(0,message)
1954 endif
1955
1956 END SUBROUTINE process_percent_cat_new
1957
1958 SUBROUTINE process_soil_real ( tsk , tmn , &
1959 landmask , sst , &
1960 st_input , sm_input , sw_input , st_levels_input , sm_levels_input , sw_levels_input , &
1961 zs , dzs , tslb , smois , sh2o , &
1962 flag_sst , flag_soilt000, flag_soilm000, &
1963 ids , ide , jds , jde , kds , kde , &
1964 ims , ime , jms , jme , kms , kme , &
1965 its , ite , jts , jte , kts , kte , &
1966 sf_surface_physics , num_soil_layers , real_data_init_type , &
1967 num_st_levels_input , num_sm_levels_input , num_sw_levels_input , &
1968 num_st_levels_alloc , num_sm_levels_alloc , num_sw_levels_alloc )
1969
1970 IMPLICIT NONE
1971
1972 INTEGER , INTENT(IN) :: ids , ide , jds , jde , kds , kde , &
1973 ims , ime , jms , jme , kms , kme , &
1974 its , ite , jts , jte , kts , kte , &
1975 sf_surface_physics , num_soil_layers , real_data_init_type , &
1976 num_st_levels_input , num_sm_levels_input , num_sw_levels_input , &
1977 num_st_levels_alloc , num_sm_levels_alloc , num_sw_levels_alloc
1978
1979 INTEGER , INTENT(IN) :: flag_sst, flag_soilt000, flag_soilm000
1980
1981 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(IN) :: landmask , sst
1982
1983 INTEGER , DIMENSION(1:num_st_levels_input) , INTENT(INOUT) :: st_levels_input
1984 INTEGER , DIMENSION(1:num_sm_levels_input) , INTENT(INOUT) :: sm_levels_input
1985 INTEGER , DIMENSION(1:num_sw_levels_input) , INTENT(INOUT) :: sw_levels_input
1986 REAL , DIMENSION(ims:ime,1:num_st_levels_alloc,jms:jme) , INTENT(INOUT) :: st_input
1987 REAL , DIMENSION(ims:ime,1:num_sm_levels_alloc,jms:jme) , INTENT(INOUT) :: sm_input
1988 REAL , DIMENSION(ims:ime,1:num_sw_levels_alloc,jms:jme) , INTENT(INOUT) :: sw_input
1989
1990 REAL, DIMENSION(1:num_soil_layers), INTENT(OUT) :: zs,dzs
1991 REAL , DIMENSION(ims:ime,num_soil_layers,jms:jme) , INTENT(OUT) :: tslb , smois , sh2o
1992 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(IN) :: tmn
1993 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(INOUT) :: tsk
1994
1995 INTEGER :: i , j , l , dominant_index , num_soil_cat , num_veg_cat
1996 REAL :: dominant_value
1997
1998 ! Initialize the soil depth, and the soil temperature and moisture.
1999
2000 IF ( ( sf_surface_physics .EQ. 1 ) .AND. ( num_soil_layers .GT. 1 ) ) THEN
2001 CALL init_soil_depth_1 ( zs , dzs , num_soil_layers )
2002 CALL init_soil_1_real ( tsk , tmn , tslb , zs , dzs , num_soil_layers , real_data_init_type , &
2003 landmask , sst , flag_sst , &
2004 ids , ide , jds , jde , kds , kde , &
2005 ims , ime , jms , jme , kms , kme , &
2006 its , ite , jts , jte , kts , kte )
2007
2008 ELSE IF ( ( sf_surface_physics .EQ. 2 .or. sf_surface_physics .EQ. 99 ) .AND. ( num_soil_layers .GT. 1 ) ) THEN
2009 CALL init_soil_depth_2 ( zs , dzs , num_soil_layers )
2010 CALL init_soil_2_real ( tsk , tmn , smois , sh2o , tslb , &
2011 st_input , sm_input , sw_input , landmask , sst , &
2012 zs , dzs , &
2013 st_levels_input , sm_levels_input , sw_levels_input , &
2014 num_soil_layers , num_st_levels_input , num_sm_levels_input , num_sw_levels_input , &
2015 num_st_levels_alloc , num_sm_levels_alloc , num_sw_levels_alloc , &
2016 flag_sst , flag_soilt000 , flag_soilm000 , &
2017 ids , ide , jds , jde , kds , kde , &
2018 ims , ime , jms , jme , kms , kme , &
2019 its , ite , jts , jte , kts , kte )
2020 ! CALL init_soil_old ( tsk , tmn , &
2021 ! smois , tslb , zs , dzs , num_soil_layers , &
2022 ! st000010_input , st010040_input , st040100_input , st100200_input , &
2023 ! st010200_input , &
2024 ! sm000010_input , sm010040_input , sm040100_input , sm100200_input , &
2025 ! sm010200_input , &
2026 ! landmask_input , sst_input , &
2027 ! ids , ide , jds , jde , kds , kde , &
2028 ! ims , ime , jms , jme , kms , kme , &
2029 ! its , ite , jts , jte , kts , kte )
2030 ELSE IF ( ( sf_surface_physics .EQ. 3 ) .AND. ( num_soil_layers .GT. 1 ) ) THEN
2031 CALL init_soil_depth_3 ( zs , dzs , num_soil_layers )
2032 CALL init_soil_3_real ( tsk , tmn , smois , tslb , &
2033 st_input , sm_input , landmask , sst , &
2034 zs , dzs , &
2035 st_levels_input , sm_levels_input , &
2036 num_soil_layers , num_st_levels_input , num_sm_levels_input , &
2037 num_st_levels_alloc , num_sm_levels_alloc , &
2038 flag_sst , flag_soilt000 , flag_soilm000 , &
2039 ids , ide , jds , jde , kds , kde , &
2040 ims , ime , jms , jme , kms , kme , &
2041 its , ite , jts , jte , kts , kte )
2042 END IF
2043
2044 END SUBROUTINE process_soil_real
2045
2046 SUBROUTINE process_soil_ideal ( xland,xice,vegfra,snow,canwat, &
2047 ivgtyp,isltyp,tslb,smois, &
2048 tsk,tmn,zs,dzs, &
2049 num_soil_layers, &
2050 sf_surface_physics , &
2051 ids,ide, jds,jde, kds,kde,&
2052 ims,ime, jms,jme, kms,kme,&
2053 its,ite, jts,jte, kts,kte )
2054
2055 IMPLICIT NONE
2056
2057 INTEGER, INTENT(IN) ::ids,ide, jds,jde, kds,kde, &
2058 ims,ime, jms,jme, kms,kme, &
2059 its,ite, jts,jte, kts,kte
2060
2061 INTEGER, INTENT(IN) :: num_soil_layers , sf_surface_physics
2062
2063 REAL, DIMENSION( ims:ime, num_soil_layers, jms:jme ) , INTENT(INOUT) :: smois, tslb
2064
2065 REAL, DIMENSION(num_soil_layers), INTENT(OUT) :: dzs,zs
2066
2067 REAL, DIMENSION( ims:ime, jms:jme ) , INTENT(INOUT) :: tsk, tmn
2068 REAL, DIMENSION( ims:ime, jms:jme ) , INTENT(OUT) :: xland, snow, canwat, xice, vegfra
2069 INTEGER, DIMENSION( ims:ime, jms:jme ) , INTENT(OUT) :: ivgtyp, isltyp
2070
2071 ! Local variables.
2072
2073 INTEGER :: itf,jtf
2074
2075 itf=MIN(ite,ide-1)
2076 jtf=MIN(jte,jde-1)
2077
2078 IF ( ( sf_surface_physics .EQ. 1 ) .AND. ( num_soil_layers .GT. 1 ) ) THEN
2079 CALL init_soil_depth_1 ( zs , dzs , num_soil_layers )
2080 CALL init_soil_1_ideal(tsk,tmn,tslb,xland, &
2081 ivgtyp,zs,dzs,num_soil_layers, &
2082 ids,ide, jds,jde, kds,kde, &
2083 ims,ime, jms,jme, kms,kme, &
2084 its,ite, jts,jte, kts,kte )
2085 ELSE IF ( ( sf_surface_physics .EQ. 2 ) .AND. ( num_soil_layers .GT. 1 ) ) THEN
2086 CALL init_soil_depth_2 ( zs , dzs , num_soil_layers )
2087 CALL init_soil_2_ideal ( xland,xice,vegfra,snow,canwat, &
2088 ivgtyp,isltyp,tslb,smois,tmn, &
2089 num_soil_layers, &
2090 ids,ide, jds,jde, kds,kde, &
2091 ims,ime, jms,jme, kms,kme, &
2092 its,ite, jts,jte, kts,kte )
2093 ELSE IF ( ( sf_surface_physics .EQ. 3 ) .AND. ( num_soil_layers .GT. 1 ) ) THEN
2094 CALL init_soil_depth_3 ( zs , dzs , num_soil_layers )
2095
2096 END IF
2097
2098 END SUBROUTINE process_soil_ideal
2099
2100 SUBROUTINE adjust_soil_temp_new ( tmn , sf_surface_physics , &
2101 tsk , ter , toposoil , landmask , flag_toposoil , &
2102 st000010 , st010040 , st040100 , st100200 , st010200 , &
2103 flag_st000010 , flag_st010040 , flag_st040100 , flag_st100200 , flag_st010200 , &
2104 soilt000 , soilt005 , soilt020 , soilt040 , soilt160 , soilt300 , &
2105 flag_soilt000 , flag_soilt005 , flag_soilt020 , flag_soilt040 , flag_soilt160 , flag_soilt300 , &
2106 ids , ide , jds , jde , kds , kde , &
2107 ims , ime , jms , jme , kms , kme , &
2108 its , ite , jts , jte , kts , kte )
2109
2110 IMPLICIT NONE
2111
2112 INTEGER , INTENT(IN) :: ids , ide , jds , jde , kds , kde , &
2113 ims , ime , jms , jme , kms , kme , &
2114 its , ite , jts , jte , kts , kte
2115
2116 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(IN) :: ter , toposoil , landmask
2117 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(INOUT) :: tmn , tsk
2118 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(INOUT) :: st000010 , st010040 , st040100 , st100200 , st010200
2119 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(INOUT) :: soilt000 , soilt005 , soilt020 , soilt040 , soilt160 , soilt300
2120
2121 INTEGER , INTENT(IN) :: flag_st000010 , flag_st010040 , flag_st040100 , flag_st100200 , flag_st010200
2122 INTEGER , INTENT(IN) :: flag_soilt000 , flag_soilt005 , flag_soilt020 , flag_soilt040 , flag_soilt160 , flag_soilt300
2123 INTEGER , INTENT(IN) :: sf_surface_physics , flag_toposoil
2124
2125 INTEGER :: i , j
2126
2127 REAL :: soil_elev_min_val , soil_elev_max_val , soil_elev_min_dif , soil_elev_max_dif
2128
2129 ! Do we have a soil field with which to modify soil temperatures?
2130
2131 IF ( flag_toposoil .EQ. 1 ) THEN
2132
2133 DO j = jts , MIN(jde-1,jte)
2134 DO i = its , MIN(ide-1,ite)
2135
2136 ! Is the toposoil field OK, or is it a subversive soil elevation field. We can tell
2137 ! usually by looking at values. Anything less than -1000 m (lower than the Dead Sea) is
2138 ! bad. Anything larger than 10 km (taller than Everest) is toast. Also, anything where
2139 ! the difference between the soil elevation and the terrain is greater than 3 km means
2140 ! that the soil data is either all zeros or that the data are inconsistent. Any of these
2141 ! three conditions is grievous enough to induce a WRF fatality. However, if they are at
2142 ! a water point, then we can safely ignore them.
2143
2144 soil_elev_min_val = toposoil(i,j)
2145 soil_elev_max_val = toposoil(i,j)
2146 soil_elev_min_dif = ter(i,j) - toposoil(i,j)
2147 soil_elev_max_dif = ter(i,j) - toposoil(i,j)
2148
2149 IF ( ( soil_elev_min_val .LT. -1000 ) .AND. ( landmask(i,j) .LT. 0.5 ) ) THEN
2150 CYCLE
2151 ELSE IF ( ( soil_elev_min_val .LT. -1000 ) .AND. ( landmask(i,j) .GT. 0.5 ) ) THEN
2152 !print *,'no soil temperature elevation adjustment, soil height too small = ',toposoil(i,j)
2153 cycle
2154 ! CALL wrf_error_fatal ( 'TOPOSOIL values have large negative values < -1000 m, unrealistic.' )
2155 ENDIF
2156
2157 IF ( ( soil_elev_max_val .GT. 10000 ) .AND. ( landmask(i,j) .LT. 0.5 ) ) THEN
2158 CYCLE
2159 ELSE IF ( ( soil_elev_max_val .GT. 10000 ) .AND. ( landmask(i,j) .GT. 0.5 ) ) THEN
2160 print *,'no soil temperature elevation adjustment, soil height too high = ',toposoil(i,j)
2161 cycle
2162 CALL wrf_error_fatal ( 'TOPOSOIL values have large positive values > 10,000 m , unrealistic.' )
2163 ENDIF
2164
2165 IF ( ( ( soil_elev_min_dif .LT. -3000 ) .OR. ( soil_elev_max_dif .GT. 3000 ) ) .AND. &
2166 ( landmask(i,j) .LT. 0.5 ) ) THEN
2167 CYCLE
2168 ELSE IF ( ( ( soil_elev_min_dif .LT. -3000 ) .OR. ( soil_elev_max_dif .GT. 3000 ) ) .AND. &
2169 ( landmask(i,j) .GT. 0.5 ) ) THEN
2170 print *,'no soil temperature elevation adjustment, diff of soil height and terrain = ',ter(i,j) - toposoil(i,j)
2171 cycle
2172 CALL wrf_error_fatal ( 'TOPOSOIL difference with terrain elevation differs by more than 3000 m, unrealistic' )
2173 ENDIF
2174
2175 ! For each of the fields that we would like to modify, check to see if it came in from the SI.
2176 ! If so, then use a -6.5 K/km lapse rate (based on the elevation diffs). We only adjust when we
2177 ! are not at a water point.
2178
2179 IF (landmask(i,j) .GT. 0.5 ) THEN
2180
2181 IF ( sf_surface_physics .EQ. 1 ) THEN
2182 tmn(i,j) = tmn(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
2183 END IF
2184
2185 tsk(i,j) = tsk(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
2186
2187 IF ( flag_st000010 .EQ. 1 ) THEN
2188 st000010(i,j) = st000010(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
2189 END IF
2190 IF ( flag_st010040 .EQ. 1 ) THEN
2191 st010040(i,j) = st010040(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
2192 END IF
2193 IF ( flag_st040100 .EQ. 1 ) THEN
2194 st040100(i,j) = st040100(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
2195 END IF
2196 IF ( flag_st100200 .EQ. 1 ) THEN
2197 st100200(i,j) = st100200(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
2198 END IF
2199 IF ( flag_st010200 .EQ. 1 ) THEN
2200 st010200(i,j) = st010200(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
2201 END IF
2202
2203 IF ( flag_soilt000 .EQ. 1 ) THEN
2204 soilt000(i,j) = soilt000(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
2205 END IF
2206 IF ( flag_soilt005 .EQ. 1 ) THEN
2207 soilt005(i,j) = soilt005(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
2208 END IF
2209 IF ( flag_soilt020 .EQ. 1 ) THEN
2210 soilt020(i,j) = soilt020(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
2211 END IF
2212 IF ( flag_soilt040 .EQ. 1 ) THEN
2213 soilt040(i,j) = soilt040(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
2214 END IF
2215 IF ( flag_soilt160 .EQ. 1 ) THEN
2216 soilt160(i,j) = soilt160(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
2217 END IF
2218 IF ( flag_soilt300 .EQ. 1 ) THEN
2219 soilt300(i,j) = soilt300(i,j) - 0.0065 * ( ter(i,j) - toposoil(i,j) )
2220 END IF
2221
2222 END IF
2223 END DO
2224 END DO
2225
2226 END IF
2227
2228 END SUBROUTINE adjust_soil_temp_new
2229
2230
2231 SUBROUTINE init_soil_depth_1 ( zs , dzs , num_soil_layers )
2232
2233 IMPLICIT NONE
2234
2235 INTEGER, INTENT(IN) :: num_soil_layers
2236
2237 REAL, DIMENSION(1:num_soil_layers), INTENT(OUT) :: zs,dzs
2238
2239 INTEGER :: l
2240
2241 CHARACTER (LEN=132) :: message
2242
2243 ! Define layers (top layer = 0.01 m). Double the thicknesses at each step (dzs values).
2244 ! The distance from the ground level to the midpoint of the layer is given by zs.
2245
2246 ! ------- Ground Level ---------- || || || ||
2247 ! || || || || zs(1) = 0.005 m
2248 ! -- -- -- -- -- -- -- -- -- || || || \/
2249 ! || || ||
2250 ! ----------------------------------- || || || \/ dzs(1) = 0.01 m
2251 ! || || ||
2252 ! || || || zs(2) = 0.02
2253 ! -- -- -- -- -- -- -- -- -- || || \/
2254 ! || ||
2255 ! || ||
2256 ! ----------------------------------- || || \/ dzs(2) = 0.02 m
2257 ! || ||
2258 ! || ||
2259 ! || ||
2260 ! || || zs(3) = 0.05
2261 ! -- -- -- -- -- -- -- -- -- || \/
2262 ! ||
2263 ! ||
2264 ! ||
2265 ! ||
2266 ! ----------------------------------- \/ dzs(3) = 0.04 m
2267
2268 IF ( num_soil_layers .NE. 5 ) THEN
2269 write(message,FMT= '(A)') 'Usually, the 5-layer diffusion uses 5 layers. Change this in the namelist.'
2270 CALL wrf_error_fatal ( message )
2271 END IF
2272
2273 dzs(1)=.01
2274 zs(1)=.5*dzs(1)
2275
2276 DO l=2,num_soil_layers
2277 dzs(l)=2*dzs(l-1)
2278 zs(l)=zs(l-1)+.5*dzs(l-1)+.5*dzs(l)
2279 ENDDO
2280
2281 END SUBROUTINE init_soil_depth_1
2282
2283 SUBROUTINE init_soil_depth_2 ( zs , dzs , num_soil_layers )
2284
2285 IMPLICIT NONE
2286
2287 INTEGER, INTENT(IN) :: num_soil_layers
2288
2289 REAL, DIMENSION(1:num_soil_layers), INTENT(OUT) :: zs,dzs
2290
2291 INTEGER :: l
2292
2293 CHARACTER (LEN=132) :: message
2294
2295 dzs = (/ 0.1 , 0.3 , 0.6 , 1.0 /)
2296
2297 IF ( num_soil_layers .NE. 4 ) THEN
2298 write(message,FMT='(A)') 'Usually, the LSM uses 4 layers. Change this in the namelist.'
2299 CALL wrf_error_fatal ( message )
2300 END IF
2301
2302 zs(1)=.5*dzs(1)
2303
2304 DO l=2,num_soil_layers
2305 zs(l)=zs(l-1)+.5*dzs(l-1)+.5*dzs(l)
2306 ENDDO
2307
2308 END SUBROUTINE init_soil_depth_2
2309
2310 SUBROUTINE init_soil_depth_3 ( zs , dzs , num_soil_layers )
2311
2312 IMPLICIT NONE
2313
2314 INTEGER, INTENT(IN) :: num_soil_layers
2315
2316 REAL, DIMENSION(1:num_soil_layers), INTENT(OUT) :: zs,dzs
2317
2318 INTEGER :: l
2319
2320 CHARACTER (LEN=132) :: message
2321
2322 ! in RUC LSM ZS - soil levels, and DZS - soil layer thicknesses, not used
2323 ! ZS is specified in the namelist: num_soil_layers = 6 or 9.
2324 ! Other options with number of levels are possible, but
2325 ! WRF users should change consistently the namelist entry with the
2326 ! ZS array in this subroutine.
2327
2328 IF ( num_soil_layers .EQ. 6) THEN
2329 zs = (/ 0.00 , 0.05 , 0.20 , 0.40 , 1.60 , 3.00 /)
2330 ! dzs = (/ 0.00 , 0.125, 0.175 , 0.70 , 1.30 , 1.40 /)
2331 ELSEIF ( num_soil_layers .EQ. 9) THEN
2332 zs = (/ 0.00 , 0.05 , 0.20 , 0.40 , 0.60, 1.00, 1.60 , 2.20, 3.00 /)
2333 ! dzs = (/ 0.00 , 0.125, 0.175 , 0.70 , 1.30 , 1.40 /)
2334 ENDIF
2335
2336 IF ( num_soil_layers .EQ. 4 .OR. num_soil_layers .EQ. 5 ) THEN
2337 WRITE(message,FMT= '(A)')'Usually, the RUC LSM uses 6, 9 or more levels. Change this in the namelist.'
2338 CALL wrf_error_fatal ( message )
2339 END IF
2340
2341 END SUBROUTINE init_soil_depth_3
2342
2343 SUBROUTINE init_soil_1_real ( tsk , tmn , tslb , zs , dzs , &
2344 num_soil_layers , real_data_init_type , &
2345 landmask , sst , flag_sst , &
2346 ids , ide , jds , jde , kds , kde , &
2347 ims , ime , jms , jme , kms , kme , &
2348 its , ite , jts , jte , kts , kte )
2349
2350 IMPLICIT NONE
2351
2352 INTEGER , INTENT(IN) :: num_soil_layers , real_data_init_type , &
2353 ids , ide , jds , jde , kds , kde , &
2354 ims , ime , jms , jme , kms , kme , &
2355 its , ite , jts , jte , kts , kte
2356
2357 INTEGER , INTENT(IN) :: flag_sst
2358
2359 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(IN) :: landmask , sst
2360 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(IN) :: tsk , tmn
2361
2362 REAL , DIMENSION(num_soil_layers) :: zs , dzs
2363
2364 REAL , DIMENSION(ims:ime,num_soil_layers,jms:jme) , INTENT(OUT) :: tslb
2365
2366 INTEGER :: i , j , l
2367
2368 ! Soil temperature is linearly interpolated between the skin temperature (taken to be at a
2369 ! depth of 0.5 cm) and the deep soil, annual temperature (taken to be at a depth of 23 cm).
2370 ! The tslb(i,1,j) is the skin temperature, and the tslb(i,num_soil_layers,j) level is the
2371 ! annual mean temperature.
2372
2373 DO j = jts , MIN(jde-1,jte)
2374 DO i = its , MIN(ide-1,ite)
2375 IF ( landmask(i,j) .GT. 0.5 ) THEN
2376 DO l = 1 , num_soil_layers
2377 tslb(i,l,j)= ( tsk(i,j) * ( zs(num_soil_layers) - zs(l) ) + &
2378 tmn(i,j) * ( zs( l) - zs(1) ) ) / &
2379 ( zs(num_soil_layers) - zs(1) )
2380 END DO
2381 ELSE
2382 IF ( ( real_data_init_type .EQ. 1 ) .AND. ( flag_sst .EQ. 1 ) ) THEN
2383 DO l = 1 , num_soil_layers
2384 tslb(i,l,j)= sst(i,j)
2385 END DO
2386 ELSE
2387 DO l = 1 , num_soil_layers
2388 tslb(i,l,j)= tsk(i,j)
2389 END DO
2390 END IF
2391 END IF
2392 END DO
2393 END DO
2394
2395 END SUBROUTINE init_soil_1_real
2396
2397 SUBROUTINE init_soil_1_ideal(tsk,tmn,tslb,xland, &
2398 ivgtyp,ZS,DZS,num_soil_layers, &
2399 ids,ide, jds,jde, kds,kde, &
2400 ims,ime, jms,jme, kms,kme, &
2401 its,ite, jts,jte, kts,kte )
2402
2403 IMPLICIT NONE
2404
2405 INTEGER, INTENT(IN ) :: ids,ide, jds,jde, kds,kde, &
2406 ims,ime, jms,jme, kms,kme, &
2407 its,ite, jts,jte, kts,kte
2408
2409 INTEGER, INTENT(IN ) :: num_soil_layers
2410
2411 REAL, DIMENSION( ims:ime , 1 , jms:jme ), INTENT(OUT) :: tslb
2412 REAL, DIMENSION( ims:ime , jms:jme ), INTENT(OUT) :: xland
2413 INTEGER, DIMENSION( ims:ime , jms:jme ), INTENT(OUT) :: ivgtyp
2414
2415 REAL, DIMENSION(1:), INTENT(IN) :: dzs,zs
2416
2417 REAL, DIMENSION( ims:ime, jms:jme ) , INTENT(IN) :: tsk, tmn
2418
2419 ! Lcal variables.
2420
2421 INTEGER :: l,j,i,itf,jtf
2422
2423 itf=MIN(ite,ide-1)
2424 jtf=MIN(jte,jde-1)
2425
2426 IF (num_soil_layers.NE.1)THEN
2427 DO j=jts,jtf
2428 DO l=1,num_soil_layers
2429 DO i=its,itf
2430 tslb(i,l,j)=( tsk(i,j)*(zs(num_soil_layers)-zs(l)) + tmn(i,j)*(zs(l)-zs(1)) ) / &
2431 ( zs(num_soil_layers)-zs(1) )
2432 ENDDO
2433 ENDDO
2434 ENDDO
2435 ENDIF
2436 DO j=jts,jtf
2437 DO i=its,itf
2438 xland(i,j) = 2
2439 ivgtyp(i,j) = 7
2440 ENDDO
2441 ENDDO
2442
2443 END SUBROUTINE init_soil_1_ideal
2444
2445 SUBROUTINE init_soil_2_real ( tsk , tmn , smois , sh2o , tslb , &
2446 st_input , sm_input , sw_input , landmask , sst , &
2447 zs , dzs , &
2448 st_levels_input , sm_levels_input , sw_levels_input , &
2449 num_soil_layers , num_st_levels_input , num_sm_levels_input , num_sw_levels_input , &
2450 num_st_levels_alloc , num_sm_levels_alloc , num_sw_levels_alloc , &
2451 flag_sst , flag_soilt000 , flag_soilm000 , &
2452 ids , ide , jds , jde , kds , kde , &
2453 ims , ime , jms , jme , kms , kme , &
2454 its , ite , jts , jte , kts , kte )
2455
2456 IMPLICIT NONE
2457
2458 INTEGER , INTENT(IN) :: num_soil_layers , &
2459 num_st_levels_input , num_sm_levels_input , num_sw_levels_input , &
2460 num_st_levels_alloc , num_sm_levels_alloc , num_sw_levels_alloc , &
2461 ids , ide , jds , jde , kds , kde , &
2462 ims , ime , jms , jme , kms , kme , &
2463 its , ite , jts , jte , kts , kte
2464
2465 INTEGER , INTENT(IN) :: flag_sst, flag_soilt000, flag_soilm000
2466
2467 INTEGER , DIMENSION(1:num_st_levels_input) , INTENT(INOUT) :: st_levels_input
2468 INTEGER , DIMENSION(1:num_sm_levels_input) , INTENT(INOUT) :: sm_levels_input
2469 INTEGER , DIMENSION(1:num_sw_levels_input) , INTENT(INOUT) :: sw_levels_input
2470
2471 REAL , DIMENSION(ims:ime,1:num_st_levels_alloc,jms:jme) , INTENT(INOUT) :: st_input
2472 REAL , DIMENSION(ims:ime,1:num_sm_levels_alloc,jms:jme) , INTENT(INOUT) :: sm_input
2473 REAL , DIMENSION(ims:ime,1:num_sw_levels_alloc,jms:jme) , INTENT(INOUT) :: sw_input
2474 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(IN) :: landmask , sst
2475
2476 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(IN) :: tmn
2477 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(INOUT) :: tsk
2478 REAL , DIMENSION(num_soil_layers) :: zs , dzs
2479
2480 REAL , DIMENSION(ims:ime,num_soil_layers,jms:jme) , INTENT(OUT) :: tslb , smois , sh2o
2481
2482 REAL , ALLOCATABLE , DIMENSION(:) :: zhave
2483
2484 INTEGER :: i , j , l , lout , lin , lwant , lhave , num
2485 REAL :: temp
2486 LOGICAL :: found_levels
2487
2488 CHARACTER (LEN=132) :: message
2489
2490 ! Are there any soil temp and moisture levels - ya know, they are mandatory.
2491
2492 num = num_st_levels_input * num_sm_levels_input
2493
2494 IF ( num .GE. 1 ) THEN
2495
2496 ! Ordered levels that we have data for.
2497 IF ( flag_soilt000 .eq. 1 ) THEN
2498 write(message, FMT='(A)') ' Assume RUC LSM 6-level input'
2499 CALL wrf_message ( message )
2500 ALLOCATE ( zhave( MAX(num_st_levels_input,num_sm_levels_input,num_sw_levels_input) ) )
2501 ELSE
2502 write(message, FMT='(A)') ' Assume Naoh LSM input'
2503 CALL wrf_message ( message )
2504 ALLOCATE ( zhave( MAX(num_st_levels_input,num_sm_levels_input,num_sw_levels_input) +2) )
2505 END IF
2506
2507
2508 ! Sort the levels for temperature.
2509 !print*,'num_st_levels_input, st_levels_input', num_st_levels_input, st_levels_input
2510 !print*,'num_sm_levels_input,num_sw_levels_input',num_sm_levels_input,num_sw_levels_input
2511
2512 outert : DO lout = 1 , num_st_levels_input-1
2513 innert : DO lin = lout+1 , num_st_levels_input
2514 IF ( st_levels_input(lout) .GT. st_levels_input(lin) ) THEN
2515 temp = st_levels_input(lout)
2516 st_levels_input(lout) = st_levels_input(lin)
2517 st_levels_input(lin) = NINT(temp)
2518 DO j = jts , MIN(jde-1,jte)
2519 DO i = its , MIN(ide-1,ite)
2520 temp = st_input(i,lout+1,j)
2521 st_input(i,lout+1,j) = st_input(i,lin+1,j)
2522 st_input(i,lin+1,j) = temp
2523 END DO
2524 END DO
2525 END IF
2526 END DO innert
2527 END DO outert
2528 !tgs add IF
2529 IF ( flag_soilt000 .NE. 1 ) THEN
2530 DO j = jts , MIN(jde-1,jte)
2531 DO i = its , MIN(ide-1,ite)
2532 st_input(i,1,j) = tsk(i,j)
2533 st_input(i,num_st_levels_input+2,j) = tmn(i,j)
2534 END DO
2535 END DO
2536 ENDIF
2537
2538
2539 #if ( NMM_CORE == 1 )
2540 !new
2541 ! write(0,*) 'st_input(1) in init_soil_2_real'
2542 DO J=MIN(jde-1,jte),JTS, -MIN(jde-1,jte)/15
2543 ! write(0,616) (st_input(I,1,J),I=its , MIN(ide-1,ite),MIN(ide-1,ite)/10)
2544 ENDDO
2545
2546 ! write(0,*) 'st_input(2) in init_soil_2_real'
2547 DO J=MIN(jde-1,jte),JTS, -MIN(jde-1,jte)/15
2548 ! write(0,616) (st_input(I,2,J),I=its , MIN(ide-1,ite),MIN(ide-1,ite)/10)
2549 ENDDO
2550
2551 616 format(20(f4.0,1x))
2552 #endif
2553
2554 ! Sort the levels for moisture.
2555
2556 outerm: DO lout = 1 , num_sm_levels_input-1
2557 innerm : DO lin = lout+1 , num_sm_levels_input
2558 IF ( sm_levels_input(lout) .GT. sm_levels_input(lin) ) THEN
2559 temp = sm_levels_input(lout)
2560 sm_levels_input(lout) = sm_levels_input(lin)
2561 sm_levels_input(lin) = NINT(temp)
2562 DO j = jts , MIN(jde-1,jte)
2563 DO i = its , MIN(ide-1,ite)
2564 temp = sm_input(i,lout+1,j)
2565 sm_input(i,lout+1,j) = sm_input(i,lin+1,j)
2566 sm_input(i,lin+1,j) = temp
2567 END DO
2568 END DO
2569 END IF
2570 END DO innerm
2571 END DO outerm
2572 !tgs add IF
2573 IF ( flag_soilm000 .NE. 1 ) THEN
2574 DO j = jts , MIN(jde-1,jte)
2575 DO i = its , MIN(ide-1,ite)
2576 sm_input(i,1,j) = sm_input(i,2,j)
2577 sm_input(i,num_sm_levels_input+2,j) = sm_input(i,num_sm_levels_input+1,j)
2578 END DO
2579 END DO
2580 ENDIF
2581
2582 ! Sort the levels for liquid moisture.
2583
2584 outerw: DO lout = 1 , num_sw_levels_input-1
2585 innerw : DO lin = lout+1 , num_sw_levels_input
2586 IF ( sw_levels_input(lout) .GT. sw_levels_input(lin) ) THEN
2587 temp = sw_levels_input(lout)
2588 sw_levels_input(lout) = sw_levels_input(lin)
2589 sw_levels_input(lin) = NINT(temp)
2590 DO j = jts , MIN(jde-1,jte)
2591 DO i = its , MIN(ide-1,ite)
2592 temp = sw_input(i,lout+1,j)
2593 sw_input(i,lout+1,j) = sw_input(i,lin+1,j)
2594 sw_input(i,lin+1,j) = temp
2595 END DO
2596 END DO
2597 END IF
2598 END DO innerw
2599 END DO outerw
2600 IF ( num_sw_levels_input .GT. 1 ) THEN
2601 DO j = jts , MIN(jde-1,jte)
2602 DO i = its , MIN(ide-1,ite)
2603 sw_input(i,1,j) = sw_input(i,2,j)
2604 sw_input(i,num_sw_levels_input+2,j) = sw_input(i,num_sw_levels_input+1,j)
2605 END DO
2606 END DO
2607 END IF
2608
2609 found_levels = .TRUE.
2610
2611 ELSE IF ( ( num .LE. 0 ) .AND. ( start_date .NE. current_date ) ) THEN
2612
2613 found_levels = .FALSE.
2614
2615 ELSE
2616 CALL wrf_error_fatal ( &
2617 'No input soil level data (temperature, moisture or liquid, or all are missing). Required for LSM.' )
2618 END IF
2619
2620 ! Is it OK to continue?
2621
2622 IF ( found_levels ) THEN
2623
2624 !tgs add IF
2625 IF ( flag_soilt000 .NE.1) THEN
2626 !tgs initialize from Noah data
2627 ! Here are the levels that we have from the input for temperature. The input levels plus
2628 ! two more: the skin temperature at 0 cm, and the annual mean temperature at 300 cm.
2629
2630 zhave(1) = 0.
2631 DO l = 1 , num_st_levels_input
2632 zhave(l+1) = st_levels_input(l) / 100.
2633 END DO
2634 zhave(num_st_levels_input+2) = 300. / 100.
2635
2636 ! Interpolate between the layers we have (zhave) and those that we want (zs).
2637
2638 z_wantt : DO lwant = 1 , num_soil_layers
2639 z_havet : DO lhave = 1 , num_st_levels_input +2 -1
2640 IF ( ( zs(lwant) .GE. zhave(lhave ) ) .AND. &
2641 ( zs(lwant) .LE. zhave(lhave+1) ) ) THEN
2642 DO j = jts , MIN(jde-1,jte)
2643 DO i = its , MIN(ide-1,ite)
2644 tslb(i,lwant,j)= ( st_input(i,lhave ,j) * ( zhave(lhave+1) - zs (lwant) ) + &
2645 st_input(i,lhave+1,j) * ( zs (lwant ) - zhave(lhave) ) ) / &
2646 ( zhave(lhave+1) - zhave(lhave) )
2647 END DO
2648 END DO
2649 EXIT z_havet
2650 END IF
2651 END DO z_havet
2652 END DO z_wantt
2653
2654 ELSE
2655
2656 !tgs initialize from RUCLSM data
2657 DO l = 1 , num_st_levels_input
2658 zhave(l) = st_levels_input(l) / 100.
2659 END DO
2660
2661 ! Interpolate between the layers we have (zhave) and those that we want (zs).
2662
2663
2664 z_wantt_2 : DO lwant = 1 , num_soil_layers
2665 z_havet_2 : DO lhave = 1 , num_st_levels_input -1
2666 IF ( ( zs(lwant) .GE. zhave(lhave ) ) .AND. &
2667 ( zs(lwant) .LE. zhave(lhave+1) ) ) THEN
2668 DO j = jts , MIN(jde-1,jte)
2669 DO i = its , MIN(ide-1,ite)
2670 tslb(i,lwant,j)= ( st_input(i,lhave,j ) * ( zhave(lhave+1) - zs (lwant) ) + &
2671 st_input(i,lhave+1,j) * ( zs (lwant ) - zhave(lhave) ) ) / &
2672 ( zhave(lhave+1) - zhave(lhave) )
2673 END DO
2674 END DO
2675 EXIT z_havet_2
2676 END IF
2677 END DO z_havet_2
2678 END DO z_wantt_2
2679
2680 ENDIF
2681
2682
2683 IF ( flag_soilm000 .NE. 1 ) THEN
2684 !tgs initialize from Noah
2685 ! Here are the levels that we have from the input for moisture. The input levels plus
2686 ! two more: a value at 0 cm and one at 300 cm. The 0 cm value is taken to be identical
2687 ! to the most shallow layer's value. Similarly, the 300 cm value is taken to be the same
2688 ! as the most deep layer's value.
2689
2690 zhave(1) = 0.
2691 DO l = 1 , num_sm_levels_input
2692 zhave(l+1) = sm_levels_input(l) / 100.
2693 END DO
2694 zhave(num_sm_levels_input+2) = 300. / 100.
2695
2696 ! Interpolate between the layers we have (zhave) and those that we want (zs).
2697
2698 z_wantm : DO lwant = 1 , num_soil_layers
2699 z_havem : DO lhave = 1 , num_sm_levels_input +2 -1
2700 IF ( ( zs(lwant) .GE. zhave(lhave ) ) .AND. &
2701 ( zs(lwant) .LE. zhave(lhave+1) ) ) THEN
2702 DO j = jts , MIN(jde-1,jte)
2703 DO i = its , MIN(ide-1,ite)
2704 smois(i,lwant,j)= ( sm_input(i,lhave ,j) * ( zhave(lhave+1) - zs (lwant) ) + &
2705 sm_input(i,lhave+1,j) * ( zs (lwant ) - zhave(lhave) ) ) / &
2706 ( zhave(lhave+1) - zhave(lhave) )
2707 END DO
2708 END DO
2709 EXIT z_havem
2710 END IF
2711 END DO z_havem
2712 END DO z_wantm
2713
2714 ELSE
2715
2716 !tgs initialize from RUCLSM data
2717 DO l = 1 , num_sm_levels_input
2718 zhave(l) = sm_levels_input(l) / 100.
2719 END DO
2720
2721 ! Interpolate between the layers we have (zhave) and those that we want (zs).
2722
2723 z_wantm_2 : DO lwant = 1 , num_soil_layers
2724 z_havem_2 : DO lhave = 1 , num_sm_levels_input -1
2725 IF ( ( zs(lwant) .GE. zhave(lhave ) ) .AND. &
2726 ( zs(lwant) .LE. zhave(lhave+1) ) ) THEN
2727 DO j = jts , MIN(jde-1,jte)
2728 DO i = its , MIN(ide-1,ite)
2729 smois(i,lwant,j)= ( sm_input(i,lhave,j ) * ( zhave(lhave+1) - zs (lwant) ) + &
2730 sm_input(i,lhave+1,j) * ( zs (lwant ) - zhave(lhave) ) ) / &
2731 ( zhave(lhave+1) - zhave(lhave) )
2732 END DO
2733 END DO
2734 EXIT z_havem_2
2735 END IF
2736 END DO z_havem_2
2737 END DO z_wantm_2
2738
2739 ENDIF
2740 ! Any liquid soil moisture to worry about?
2741
2742 IF ( num_sw_levels_input .GT. 1 ) THEN
2743
2744 zhave(1) = 0.
2745 DO l = 1 , num_sw_levels_input
2746 zhave(l+1) = sw_levels_input(l) / 100.
2747 END DO
2748 zhave(num_sw_levels_input+2) = 300. / 100.
2749
2750 ! Interpolate between the layers we have (zhave) and those that we want (zs).
2751
2752 z_wantw : DO lwant = 1 , num_soil_layers
2753 z_havew : DO lhave = 1 , num_sw_levels_input +2 -1
2754 IF ( ( zs(lwant) .GE. zhave(lhave ) ) .AND. &
2755 ( zs(lwant) .LE. zhave(lhave+1) ) ) THEN
2756 DO j = jts , MIN(jde-1,jte)
2757 DO i = its , MIN(ide-1,ite)
2758 sh2o(i,lwant,j)= ( sw_input(i,lhave ,j) * ( zhave(lhave+1) - zs (lwant) ) + &
2759 sw_input(i,lhave+1,j) * ( zs (lwant ) - zhave(lhave) ) ) / &
2760 ( zhave(lhave+1) - zhave(lhave) )
2761 END DO
2762 END DO
2763 EXIT z_havew
2764 END IF
2765 END DO z_havew
2766 END DO z_wantw
2767
2768 END IF
2769
2770
2771 ! Over water, put in reasonable values for soil temperature and moisture. These won't be
2772 ! used, but they will make a more continuous plot.
2773
2774 IF ( flag_sst .EQ. 1 ) THEN
2775 DO j = jts , MIN(jde-1,jte)
2776 DO i = its , MIN(ide-1,ite)
2777 IF ( landmask(i,j) .LT. 0.5 ) THEN
2778 DO l = 1 , num_soil_layers
2779 tslb(i,l,j)= sst(i,j)
2780 !tgs add line for tsk
2781 tsk(i,j) = sst(i,j)
2782 smois(i,l,j)= 1.0
2783 sh2o (i,l,j)= 1.0
2784 END DO
2785 END IF
2786 END DO
2787 END DO
2788 ELSE
2789 DO j = jts , MIN(jde-1,jte)
2790 DO i = its , MIN(ide-1,ite)
2791 IF ( landmask(i,j) .LT. 0.5 ) THEN
2792 DO l = 1 , num_soil_layers
2793 tslb(i,l,j)= tsk(i,j)
2794 smois(i,l,j)= 1.0
2795 sh2o (i,l,j)= 1.0
2796 END DO
2797 END IF
2798 END DO
2799 END DO
2800 END IF
2801
2802 DEALLOCATE (zhave)
2803
2804 END IF
2805
2806 END SUBROUTINE init_soil_2_real
2807
2808 SUBROUTINE init_soil_2_ideal ( xland,xice,vegfra,snow,canwat, &
2809 ivgtyp,isltyp,tslb,smois,tmn, &
2810 num_soil_layers, &
2811 ids,ide, jds,jde, kds,kde, &
2812 ims,ime, jms,jme, kms,kme, &
2813 its,ite, jts,jte, kts,kte )
2814
2815 IMPLICIT NONE
2816
2817 INTEGER, INTENT(IN) ::ids,ide, jds,jde, kds,kde, &
2818 ims,ime, jms,jme, kms,kme, &
2819 its,ite, jts,jte, kts,kte
2820
2821 INTEGER, INTENT(IN) ::num_soil_layers
2822
2823 REAL, DIMENSION( ims:ime, num_soil_layers, jms:jme ) , INTENT(OUT) :: smois, tslb
2824
2825 REAL, DIMENSION( ims:ime, jms:jme ) , INTENT(OUT) :: xland, snow, canwat, xice, vegfra, tmn
2826
2827 INTEGER, DIMENSION( ims:ime, jms:jme ) , INTENT(OUT) :: ivgtyp, isltyp
2828
2829 INTEGER :: icm,jcm,itf,jtf
2830 INTEGER :: i,j,l
2831
2832 itf=min0(ite,ide-1)
2833 jtf=min0(jte,jde-1)
2834
2835 icm = ide/2
2836 jcm = jde/2
2837
2838 DO j=jts,jtf
2839 DO l=1,num_soil_layers
2840 DO i=its,itf
2841
2842 smois(i,1,j)=0.10
2843 smois(i,2,j)=0.10
2844 smois(i,3,j)=0.10
2845 smois(i,4,j)=0.10
2846
2847 tslb(i,1,j)=295.
2848 tslb(i,2,j)=297.
2849 tslb(i,3,j)=293.
2850 tslb(i,4,j)=293.
2851
2852 ENDDO
2853 ENDDO
2854 ENDDO
2855
2856 DO j=jts,jtf
2857 DO i=its,itf
2858 xland(i,j) = 2
2859 tmn(i,j) = 294.
2860 xice(i,j) = 0.
2861 vegfra(i,j) = 0.
2862 snow(i,j) = 0.
2863 canwat(i,j) = 0.
2864 ivgtyp(i,j) = 7
2865 isltyp(i,j) = 8
2866 ENDDO
2867 ENDDO
2868
2869 END SUBROUTINE init_soil_2_ideal
2870
2871 SUBROUTINE init_soil_3_real ( tsk , tmn , smois , tslb , &
2872 st_input , sm_input , landmask, sst, &
2873 zs , dzs , &
2874 st_levels_input , sm_levels_input , &
2875 num_soil_layers , num_st_levels_input , num_sm_levels_input , &
2876 num_st_levels_alloc , num_sm_levels_alloc , &
2877 flag_sst , flag_soilt000 , flag_soilm000 , &
2878 ids , ide , jds , jde , kds , kde , &
2879 ims , ime , jms , jme , kms , kme , &
2880 its , ite , jts , jte , kts , kte )
2881
2882 IMPLICIT NONE
2883
2884 INTEGER , INTENT(IN) :: num_soil_layers , &
2885 num_st_levels_input , num_sm_levels_input , &
2886 num_st_levels_alloc , num_sm_levels_alloc , &
2887 ids , ide , jds , jde , kds , kde , &
2888 ims , ime , jms , jme , kms , kme , &
2889 its , ite , jts , jte , kts , kte
2890
2891 INTEGER , INTENT(IN) :: flag_sst, flag_soilt000, flag_soilm000
2892
2893 INTEGER , DIMENSION(1:num_st_levels_input) , INTENT(INOUT) :: st_levels_input
2894 INTEGER , DIMENSION(1:num_sm_levels_input) , INTENT(INOUT) :: sm_levels_input
2895
2896 REAL , DIMENSION(ims:ime,1:num_st_levels_alloc,jms:jme) , INTENT(INOUT) :: st_input
2897 REAL , DIMENSION(ims:ime,1:num_sm_levels_alloc,jms:jme) , INTENT(INOUT) :: sm_input
2898 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(IN) :: landmask , sst
2899
2900 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(IN) :: tmn
2901 REAL , DIMENSION(ims:ime,jms:jme) , INTENT(INOUT) :: tsk
2902 REAL , DIMENSION(num_soil_layers) :: zs , dzs
2903
2904 REAL , DIMENSION(ims:ime,num_soil_layers,jms:jme) , INTENT(OUT) :: tslb , smois
2905
2906 REAL , ALLOCATABLE , DIMENSION(:) :: zhave
2907
2908 INTEGER :: i , j , l , lout , lin , lwant , lhave
2909 REAL :: temp
2910
2911 ! Allocate the soil layer array used for interpolating.
2912
2913 IF ( ( num_st_levels_input .LE. 0 ) .OR. &
2914 ( num_sm_levels_input .LE. 0 ) ) THEN
2915 PRINT '(A)','No input soil level data (either temperature or moisture, or both are missing). Required for RUC LSM.'
2916 CALL wrf_error_fatal ( 'no soil data' )
2917 ELSE
2918 IF ( flag_soilt000 .eq. 1 ) THEN
2919 PRINT '(A)',' Assume RUC LSM 6-level input'
2920 ALLOCATE ( zhave( MAX(num_st_levels_input,num_sm_levels_input) ) )
2921 ELSE
2922 PRINT '(A)',' Assume non-RUC LSM input'
2923 ALLOCATE ( zhave( MAX(num_st_levels_input,num_soil_layers) ) )
2924 END IF
2925 END IF
2926
2927 ! Sort the levels for temperature.
2928
2929 outert : DO lout = 1 , num_st_levels_input-1
2930 innert : DO lin = lout+1 , num_st_levels_input
2931 IF ( st_levels_input(lout) .GT. st_levels_input(lin) ) THEN
2932 temp = st_levels_input(lout)
2933 st_levels_input(lout) = st_levels_input(lin)
2934 st_levels_input(lin) = NINT(temp)
2935 DO j = jts , MIN(jde-1,jte)
2936 DO i = its , MIN(ide-1,ite)
2937 temp = st_input(i,lout,j)
2938 st_input(i,lout,j) = st_input(i,lin,j)
2939 st_input(i,lin,j) = temp
2940 END DO
2941 END DO
2942 END IF
2943 END DO innert
2944 END DO outert
2945
2946 IF ( flag_soilt000 .NE. 1 ) THEN
2947 DO j = jts , MIN(jde-1,jte)
2948 DO i = its , MIN(ide-1,ite)
2949 st_input(i,1,j) = tsk(i,j)
2950 st_input(i,num_st_levels_input+2,j) = tmn(i,j)
2951 END DO
2952 END DO
2953 END IF
2954
2955 ! Sort the levels for moisture.
2956
2957 outerm: DO lout = 1 , num_sm_levels_input-1
2958 innerm : DO lin = lout+1 , num_sm_levels_input
2959 IF ( sm_levels_input(lout) .GT. sm_levels_input(lin) ) THEN
2960 temp = sm_levels_input(lout)
2961 sm_levels_input(lout) = sm_levels_input(lin)
2962 sm_levels_input(lin) = NINT(temp)
2963 DO j = jts , MIN(jde-1,jte)
2964 DO i = its , MIN(ide-1,ite)
2965 temp = sm_input(i,lout,j)
2966 sm_input(i,lout,j) = sm_input(i,lin,j)
2967 sm_input(i,lin,j) = temp
2968 END DO
2969 END DO
2970 END IF
2971 END DO innerm
2972 END DO outerm
2973
2974 IF ( flag_soilm000 .NE. 1 ) THEN
2975 DO j = jts , MIN(jde-1,jte)
2976 DO i = its , MIN(ide-1,ite)
2977 sm_input(i,1,j) = sm_input(i,2,j)
2978 sm_input(i,num_sm_levels_input+2,j) = sm_input(i,num_sm_levels_input+1,j)
2979 END DO
2980 END DO
2981 END IF
2982
2983 ! Here are the levels that we have from the input for temperature.
2984
2985 IF ( flag_soilt000 .EQ. 1 ) THEN
2986 DO l = 1 , num_st_levels_input
2987 zhave(l) = st_levels_input(l) / 100.
2988 END DO
2989
2990 ! Interpolate between the layers we have (zhave) and those that we want (zs).
2991
2992 z_wantt : DO lwant = 1 , num_soil_layers
2993 z_havet : DO lhave = 1 , num_st_levels_input -1
2994 IF ( ( zs(lwant) .GE. zhave(lhave ) ) .AND. &
2995 ( zs(lwant) .LE. zhave(lhave+1) ) ) THEN
2996 DO j = jts , MIN(jde-1,jte)
2997 DO i = its , MIN(ide-1,ite)
2998 tslb(i,lwant,j)= ( st_input(i,lhave,j ) * ( zhave(lhave+1) - zs (lwant) ) + &
2999 st_input(i,lhave+1,j) * ( zs (lwant ) - zhave(lhave) ) ) / &
3000 ( zhave(lhave+1) - zhave(lhave) )
3001 END DO
3002 END DO
3003 EXIT z_havet
3004 END IF
3005 END DO z_havet
3006 END DO z_wantt
3007
3008 ELSE
3009
3010 zhave(1) = 0.
3011 DO l = 1 , num_st_levels_input
3012 zhave(l+1) = st_levels_input(l) / 100.
3013 END DO
3014 zhave(num_st_levels_input+2) = 300. / 100.
3015
3016 ! Interpolate between the layers we have (zhave) and those that we want (zs).
3017
3018 z_wantt_2 : DO lwant = 1 , num_soil_layers
3019 z_havet_2 : DO lhave = 1 , num_st_levels_input +2
3020 IF ( ( zs(lwant) .GE. zhave(lhave ) ) .AND. &
3021 ( zs(lwant) .LE. zhave(lhave+1) ) ) THEN
3022 DO j = jts , MIN(jde-1,jte)
3023 DO i = its , MIN(ide-1,ite)
3024 tslb(i,lwant,j)= ( st_input(i,lhave,j ) * ( zhave(lhave+1) - zs (lwant) ) + &
3025 st_input(i,lhave+1,j) * ( zs (lwant ) - zhave(lhave) ) ) / &
3026 ( zhave(lhave+1) - zhave(lhave) )
3027 END DO
3028 END DO
3029 EXIT z_havet_2
3030 END IF
3031 END DO z_havet_2
3032 END DO z_wantt_2
3033
3034 END IF
3035
3036 ! Here are the levels that we have from the input for moisture.
3037
3038 IF ( flag_soilm000 .EQ. 1 ) THEN
3039 DO l = 1 , num_sm_levels_input
3040 zhave(l) = sm_levels_input(l) / 100.
3041 END DO
3042
3043 ! Interpolate between the layers we have (zhave) and those that we want (zs).
3044
3045 z_wantm : DO lwant = 1 , num_soil_layers
3046 z_havem : DO lhave = 1 , num_sm_levels_input -1
3047 IF ( ( zs(lwant) .GE. zhave(lhave ) ) .AND. &
3048 ( zs(lwant) .LE. zhave(lhave+1) ) ) THEN
3049 DO j = jts , MIN(jde-1,jte)
3050 DO i = its , MIN(ide-1,ite)
3051 smois(i,lwant,j)= ( sm_input(i,lhave,j ) * ( zhave(lhave+1) - zs (lwant) ) + &
3052 sm_input(i,lhave+1,j) * ( zs (lwant ) - zhave(lhave) ) ) / &
3053 ( zhave(lhave+1) - zhave(lhave) )
3054 END DO
3055 END DO
3056 EXIT z_havem
3057 END IF
3058 END DO z_havem
3059 END DO z_wantm
3060
3061 ELSE
3062
3063 zhave(1) = 0.
3064 DO l = 1 , num_sm_levels_input
3065 zhave(l+1) = sm_levels_input(l) / 100.
3066 END DO
3067 zhave(num_sm_levels_input+2) = 300. / 100.
3068
3069 z_wantm_2 : DO lwant = 1 , num_soil_layers
3070 z_havem_2 : DO lhave = 1 , num_sm_levels_input +2
3071 IF ( ( zs(lwant) .GE. zhave(lhave ) ) .AND. &
3072 ( zs(lwant) .LE. zhave(lhave+1) ) ) THEN
3073 DO j = jts , MIN(jde-1,jte)
3074 DO i = its , MIN(ide-1,ite)
3075 smois(i,lwant,j)= ( sm_input(i,lhave,j ) * ( zhave(lhave+1) - zs (lwant) ) + &
3076 sm_input(i,lhave+1,j) * ( zs (lwant ) - zhave(lhave) ) ) / &
3077 ( zhave(lhave+1) - zhave(lhave) )
3078 END DO
3079 END DO
3080 EXIT z_havem_2
3081 END IF
3082 END DO z_havem_2
3083 END DO z_wantm_2
3084
3085 END IF
3086 ! Over water, put in reasonable values for soil temperature and moisture. These won't be
3087 ! used, but they will make a more continuous plot.
3088
3089 IF ( flag_sst .EQ. 1 ) THEN
3090 DO j = jts , MIN(jde-1,jte)
3091 DO i = its , MIN(ide-1,ite)
3092 IF ( landmask(i,j) .LT. 0.5 ) THEN
3093 DO l = 1 , num_soil_layers
3094 tslb(i,l,j) = sst(i,j)
3095 tsk(i,j) = sst(i,j)
3096 smois(i,l,j)= 1.0
3097 END DO
3098 END IF
3099 END DO
3100 END DO
3101 ELSE
3102 DO j = jts , MIN(jde-1,jte)
3103 DO i = its , MIN(ide-1,ite)
3104 IF ( landmask(i,j) .LT. 0.5 ) THEN
3105 DO l = 1 , num_soil_layers
3106 tslb(i,l,j)= tsk(i,j)
3107 smois(i,l,j)= 1.0
3108 END DO
3109 END IF
3110 END DO
3111 END DO
3112 END IF
3113
3114 DEALLOCATE (zhave)
3115
3116 END SUBROUTINE init_soil_3_real
3117
3118 END MODULE module_soil_pre
3119
3120 #endif
3121