Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Normal
Revision Log

module_soil_pre.F @ 3567

Last change on this file since 3567 was 2759, checked in by aslmd, 2 years ago
adding unmodified code from WRFV3.0.1.1, expurged from useless data +1M size
File size: 144.4 KB

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

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format