Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Blame
Revision Log

module_fddaobs_driver.F @ 400

Last change on this file since 400 was 1, checked in by lfita, 11 years ago

-- --- Opening of the WRF+LMDZ coupling repository --- -- -

WRF: version v3.3
LMDZ: version v1818

More details in:

File size: 29.6 KB

Line
1	!WRF:MODEL_LAYER:PHYSICS
2	MODULE module_fddaobs_driver
3
4	! This obs-nudging FDDA module (RTFDDA) is developed by the
5	! NCAR/RAL/NSAP (National Security Application Programs), under the
6	! sponsorship of ATEC (Army Test and Evaluation Commands). ATEC is
7	! acknowledged for releasing this capability for WRF community
8	! research applications.
9	!
10	! The NCAR/RAL RTFDDA module was adapted, and significantly modified
11	! from the obs-nudging module in the standard MM5V3.1 which was originally
12	! developed by PSU (Stauffer and Seaman, 1994).
13	!
14	! Yubao Liu (NCAR/RAL): lead developer of the RTFDDA module
15	! Al Bourgeois (NCAR/RAL): lead engineer implementing RTFDDA into WRF-ARW
16	! Nov. 2006
17	!
18	! References:
19	!
20	! Liu, Y., A. Bourgeois, T. Warner, S. Swerdlin and J. Hacker, 2005: An
21	! implementation of obs-nudging-based FDDA into WRF for supporting
22	! ATEC test operations. 2005 WRF user workshop. Paper 10.7.
23	!
24	! Liu, Y., A. Bourgeois, T. Warner, S. Swerdlin and W. Yu, 2006: An update
25	! on "obs-nudging"-based FDDA for WRF-ARW: Verification using OSSE
26	! and performance of real-time forecasts. 2006 WRF user workshop. Paper 4.7.
27
28	!
29	! Stauffer, D.R., and N.L. Seaman, 1994: Multi-scale four-dimensional data
30	! assimilation. J. Appl. Meteor., 33, 416-434.
31	!
32	! http://www.rap.ucar.edu/projects/armyrange/references.html
33	!
34
35	CONTAINS
36
37	!-----------------------------------------------------------------------
38	SUBROUTINE fddaobs_driver( inest, domid, parid, restart, &
39	config_flags, &
40	nudge_opt, iprt_errob, iprt_nudob, &
41	fdasta, fdaend, &
42	nudge_wind, nudge_temp, nudge_mois, &
43	nudge_pstr, &
44	coef_wind, coef_temp, coef_mois, &
45	coef_pstr, rinxy, rinsig, &
46	npfi, ionf, &
47	obs_prt_max, obs_prt_freq, idynin, dtramp, &
48	parent_grid_ratio, maxdom, itimestep, &
49	xtime, &
50	dt, gmt, julday, &
51	#if ( EM_CORE == 1 )
52	fdob, &
53	#endif
54	max_obs, nobs_ndg_vars, &
55	nobs_err_flds, nstat, varobs, errf, dx, &
56	KPBL, HT, mut, muu, muv, &
57	msftx, msfty, msfux, msfuy, msfvx, msfvy, p_phy, t_tendf, t0, &
58	ub, vb, tb, qvb, pbase, ptop, pp, phb, ph, &
59	uratx, vratx, tratx, ru_tendf, rv_tendf, &
60	moist_tend, savwt, &
61	regime, pblh, z_at_w, &
62	z, &
63	ids,ide, jds,jde, kds,kde, & ! domain dims
64	ims,ime, jms,jme, kms,kme, & ! memory dims
65	its,ite, jts,jte, kts,kte ) ! tile dims
66
67	!-----------------------------------------------------------------------
68	USE module_domain
69	USE module_bc
70	USE module_model_constants, ONLY : g, rcp
71	USE module_fddaobs_rtfdda
72
73	! This driver calls subroutines for fdda obs-nudging and
74	! returns computed tendencies
75
76	!
77	!-----------------------------------------------------------------------
78	IMPLICIT NONE
79	!-----------------------------------------------------------------------
80	! taken from MM5 code - 03 Feb 2004.
81	!-----------------------------------------------------------------------
82
83	!=======================================================================
84	! Definitions
85	!-----------
86	!-- KPBL vertical layer index for PBL top
87	!-- HT terrain height (m)
88	!-- p_phy pressure (Pa)
89	!-- t_tendf temperature tendency
90
91	INTEGER, intent(in) :: ids,ide, jds,jde, kds,kde ! domain dims.
92	INTEGER, intent(in) :: ims,ime, jms,jme, kms,kme ! memory dims.
93	INTEGER, intent(in) :: its,ite, jts,jte, kts,kte ! tile dims.
94
95	INTEGER, intent(in) :: inest
96	INTEGER, intent(in) :: maxdom
97	INTEGER, intent(in) :: domid(maxdom) ! Domain IDs
98	INTEGER, intent(in) :: parid(maxdom) ! Parent domain IDs
99	LOGICAL, intent(in) :: restart
100	TYPE(grid_config_rec_type), INTENT(IN ) :: config_flags
101	INTEGER, intent(in) :: itimestep
102	INTEGER, intent(in) :: nudge_opt
103	LOGICAL, intent(in) :: iprt_errob
104	LOGICAL, intent(in) :: iprt_nudob
105	REAL, intent(in) :: fdasta
106	REAL, intent(in) :: fdaend
107	INTEGER, intent(in) :: nudge_wind
108	INTEGER, intent(in) :: nudge_temp
109	INTEGER, intent(in) :: nudge_mois
110	INTEGER, intent(in) :: nudge_pstr
111	REAL, intent(in) :: coef_wind
112	REAL, intent(in) :: coef_temp
113	REAL, intent(in) :: coef_mois
114	REAL, intent(in) :: coef_pstr
115	REAL, intent(inout) :: rinxy
116	REAL, intent(inout) :: rinsig
117	INTEGER, intent(in) :: npfi
118	INTEGER, intent(in) :: ionf
119	INTEGER, intent(in) :: obs_prt_max ! max number of obs in printout
120	INTEGER, intent(in) :: obs_prt_freq ! frequency (in obs index) printout
121	INTEGER, intent(in) :: idynin
122	REAL, intent(inout) :: dtramp
123	INTEGER, intent(in) :: parent_grid_ratio
124	REAL, intent(in) :: xtime ! forecast time in minutes
125	REAL, intent(in) :: dt
126	REAL, intent(in) :: gmt
127	INTEGER, intent(in) :: julday
128	INTEGER, intent(in) :: max_obs ! max number of observations
129	INTEGER, intent(in) :: nobs_ndg_vars
130	INTEGER, intent(in) :: nobs_err_flds
131	INTEGER, intent(in) :: nstat
132	REAL, intent(inout) :: varobs(nobs_ndg_vars, max_obs)
133	REAL, intent(inout) :: errf(nobs_err_flds, max_obs)
134	REAL, intent(in) :: dx ! this-domain grid cell-size (m)
135	INTEGER, INTENT(IN) :: kpbl( ims:ime, jms:jme ) ! vertical layer index for PBL top
136	REAL, INTENT(IN) :: ht( ims:ime, jms:jme )
137	REAL, INTENT(IN) :: mut( ims:ime , jms:jme ) ! Air mass on t-grid
138	REAL, INTENT(IN) :: muu( ims:ime , jms:jme ) ! Air mass on u-grid
139	REAL, INTENT(IN) :: muv( ims:ime , jms:jme ) ! Air mass on v-grid
140	REAL, INTENT(IN) :: msftx( ims:ime , jms:jme ) ! Map scale on t-grid
141	REAL, INTENT(IN) :: msfty( ims:ime , jms:jme ) ! Map scale on t-grid
142	REAL, INTENT(IN) :: msfux( ims:ime , jms:jme ) ! Map scale on u-grid
143	REAL, INTENT(IN) :: msfuy( ims:ime , jms:jme ) ! Map scale on u-grid
144	REAL, INTENT(IN) :: msfvx( ims:ime , jms:jme ) ! Map scale on v-grid
145	REAL, INTENT(IN) :: msfvy( ims:ime , jms:jme ) ! Map scale on v-grid
146
147	REAL, INTENT(IN) :: p_phy( ims:ime, kms:kme, jms:jme )
148	REAL, INTENT(INOUT) :: t_tendf( ims:ime, kms:kme, jms:jme )
149	REAL, INTENT(IN) :: t0
150	REAL, INTENT(INOUT) :: savwt( nobs_ndg_vars, ims:ime, kms:kme, jms:jme )
151	REAL, INTENT(INOUT) :: regime( ims:ime, jms:jme )
152	REAL, INTENT(IN) :: pblh( ims:ime, jms:jme )
153	REAL, INTENT(IN) :: z_at_w( ims:ime, kms:kme, jms:jme ) ! Model ht(m) asl, f-levs
154	REAL, INTENT(IN) :: z( ims:ime, kms:kme, jms:jme ) ! Model ht(m) asl, h-levs
155
156	#if ( EM_CORE == 1 )
157	TYPE(fdob_type), intent(inout) :: fdob
158	#endif
159
160	REAL, INTENT(IN) :: ub( ims:ime, kms:kme, jms:jme )
161	REAL, INTENT(IN) :: vb( ims:ime, kms:kme, jms:jme )
162	REAL, INTENT(IN) :: tb( ims:ime, kms:kme, jms:jme )
163	REAL, INTENT(IN) :: qvb( ims:ime, kms:kme, jms:jme )
164	REAL, INTENT(IN) :: pbase( ims:ime, kms:kme, jms:jme ) ! Base press. (Pa)
165	REAL, INTENT(IN) :: ptop
166	REAL, INTENT(IN) :: pp( ims:ime, kms:kme, jms:jme ) ! Press. pert. (Pa)
167	REAL, INTENT(IN) :: phb( ims:ime, kms:kme, jms:jme ) ! Base geopotential
168	REAL, INTENT(IN) :: ph( ims:ime, kms:kme, jms:jme ) ! Perturbation geopotential
169	REAL, INTENT(IN) :: uratx( ims:ime, jms:jme ) ! On mass points
170	REAL, INTENT(IN) :: vratx( ims:ime, jms:jme ) ! "
171	REAL, INTENT(IN) :: tratx( ims:ime, jms:jme ) ! "
172	REAL, INTENT(INOUT) :: ru_tendf( ims:ime, kms:kme, jms:jme )
173	REAL, INTENT(INOUT) :: rv_tendf( ims:ime, kms:kme, jms:jme )
174	REAL, INTENT(INOUT) :: moist_tend( ims:ime, kms:kme, jms:jme )
175
176	! Local variables
177	logical :: nudge_flag ! Flag for doing nudging
178	integer :: KTAU ! Forecast timestep
179	real :: dtmin ! dt in minutes
180	integer :: i, j, k ! Loop counters.
181	integer :: idom ! Loop counter.
182	integer :: nsta ! Number of observation stations
183	integer :: infr ! Frequency for obs input and error calc
184	integer :: idarst ! Flag for calling sub errob on restart
185	real :: dtr ! Abs value of dtramp (for dynamic init)
186	real :: tconst ! Reciprocal of dtr
187	real :: vih_uv(its:ite,jts:jte,2) ! Vert infl heights abv grd for LML obs (wind)
188	real :: vih_t (its:ite,jts:jte,2) ! Vert infl heights abv grd for LML obs (temp)
189	real :: vih_q (its:ite,jts:jte,2) ! Vert infl heights abv grd for LML obs (mois)
190	integer :: vik_uv(its:ite,jts:jte,2) ! Vert infl k-levels for LML obs (wind)
191	integer :: vik_t (its:ite,jts:jte,2) ! Vert infl k-levels for LML obs (temp)
192	integer :: vik_q (its:ite,jts:jte,2) ! Vert infl k-levels for LML obs (mois)
193	real :: z_at_p( kms:kme ) ! Height at p levels
194	#ifdef RAL
195	real :: HTIJ(ids:ide, jds:jde) = 0. ! Terrain ht on global grid
196	#endif
197	character(len=200) :: msg ! Argument to wrf_message
198
199	#if ( EM_CORE == 1 )
200	nudge_flag = (nudge_opt .eq. 1)
201
202	if (.not. nudge_flag) return
203
204	!----------------------------------------------------------------------
205	! ************* BEGIN FDDA SETUP SECTION *************
206
207	! Calculate forecast time.
208	dtmin = dt/60.
209	ktau = itimestep - 1 !ktau corresponds to xtime
210
211	! Set NSTA to zero on startup, or else retrieve value from last pass.
212	IF(ktau.EQ.fdob%ktaur) THEN
213	if (iprt_nudob) then
214	write(msg,'(a,i2,a)') 'OBS NUDGING is requested on a total of ', &
215	fdob%domain_tot,' domain(s).'
216	call wrf_message(msg)
217	endif
218	nsta=0.
219	ELSE
220	nsta=fdob%nstat
221	ENDIF
222
223	infr = ionf(parent_grid_ratio*fdob%levidn(inest))
224	nsta=fdob%nstat
225	idarst = 0
226	IF(restart .AND. ktau.EQ.fdob%ktaur) idarst=1
227
228	CALL wrf_debug(100,'in PSU FDDA scheme')
229
230	! Make sure regime array is set over entire grid
231	! (ajb: Copied code from fddagd)
232	IF( config_flags%bl_pbl_physics /= 1 &
233	.AND. config_flags%bl_pbl_physics /= 5 &
234	.AND. config_flags%bl_pbl_physics /= 6 &
235	.AND. config_flags%bl_pbl_physics /= 7 &
236	.AND. config_flags%bl_pbl_physics /= 99 ) THEN
237	DO j = jts, jte
238	DO i = its, ite
239	IF( pblh(i,j) > z_at_w(i,2,j)-ht(i,j) ) THEN
240	regime(i,j) = 4.0
241	ELSE
242	regime(i,j) = 1.0
243	ENDIF
244	ENDDO
245	ENDDO
246	ENDIF
247
248	! Compute VIF heights for each grid column (used for LML obs)
249	if(fdob%sfc_scheme_vert.EQ.0) then
250	if(nudge_wind.EQ.1 .AND. NSTA.GT.0) then
251	CALL compute_VIH( fdob%vif_uv, fdob%vif_max, &
252	fdob%vif_fullmin, fdob%vif_rampmin, &
253	regime, pblh, &
254	ht, z, vih_uv, &
255	ids,ide, jds,jde, kds,kde, &
256	ims,ime, jms,jme, kms,kme, &
257	its,ite, jts,jte, kts,kte )
258	endif
259	if(nudge_temp.EQ.1 .AND. NSTA.GT.0) then
260	CALL compute_VIH( fdob%vif_t, fdob%vif_max, &
261	fdob%vif_fullmin, fdob%vif_rampmin, &
262	regime, pblh, &
263	ht, z, vih_t, &
264	ids,ide, jds,jde, kds,kde, &
265	ims,ime, jms,jme, kms,kme, &
266	its,ite, jts,jte, kts,kte )
267	endif
268	if(nudge_mois.EQ.1 .AND. NSTA.GT.0) then
269	CALL compute_VIH( fdob%vif_q, fdob%vif_max, &
270	fdob%vif_fullmin, fdob%vif_rampmin, &
271	regime, pblh, &
272	ht, z, vih_q, &
273	ids,ide, jds,jde, kds,kde, &
274	ims,ime, jms,jme, kms,kme, &
275	its,ite, jts,jte, kts,kte )
276	endif
277	endif
278	!******************* END AJB MOVE TO SLAB *************************
279
280	! COMPUTE ERROR BETWEEN OBSERVATIONS and MODEL
281	IF( nsta.GT.0 ) THEN
282	IF( MOD(ktau,infr).EQ.0 .OR. idarst.EQ.1) THEN
283
284	CALL errob(inest, ub, vb, tb, t0, qvb, pbase, pp, rcp, &
285	z, &
286	uratx, vratx, tratx, kpbl, &
287	nobs_ndg_vars, nobs_err_flds, max_obs, maxdom, &
288	fdob%levidn, parid, fdob%nstat, fdob%nstaw, &
289	nudge_wind, nudge_temp, nudge_mois, nudge_pstr, &
290	fdob%timeob, fdob%rio, fdob%rjo, fdob%rko, &
291	varobs, errf, ktau, xtime, &
292	parent_grid_ratio, npfi, &
293	obs_prt_max, obs_prt_freq, iprt_errob, &
294	fdob%obsprt, fdob%stnidprt, &
295	fdob%latprt, fdob%lonprt, &
296	fdob%mlatprt, fdob%mlonprt, &
297	ids,ide, jds,jde, kds,kde, &
298	ims,ime, jms,jme, kms,kme, &
299	its,ite, jts,jte, kts,kte)
300	ENDIF
301	ENDIF
302
303	fdob%tfaci=1.0
304	IF(idynin.EQ.1.AND.nudge_opt.EQ.1) THEN
305	dtr=ABS(dtramp)
306	tconst=1./dtr
307	! FDAEND(IN) IS THE TIME IN MINUTES TO END THE DYNAMIC INITIALIZATION CY
308	IF(xtime.LT.fdaend-dtr)THEN
309	fdob%tfaci=1.
310	ELSEIF(xtime.GE.fdaend-dtr.AND.xtime.LE.fdaend) THEN
311	fdob%tfaci=(fdaend-xtime)*tconst
312	ELSE
313	fdob%tfaci=0.0
314	ENDIF
315	IF(ktau.EQ.fdob%ktaur.OR.MOD(ktau,10).EQ.0) THEN
316	IF (iprt_nudob) &
317	PRINT*,' DYNINOBS: IN,KTAU,XTIME,FDAEND,DTRAMP,DTR,TCONST', &
318	',TFACI: ',INEST,KTAU,XTIME,FDAEND,DTRAMP,DTR,TCONST, &
319	fdob%TFACI
320	ENDIF
321	ENDIF
322
323	#ifdef RAL
324	! MEIXU: collect terrain array HT into a global array HTIJ
325	CALL loc2glob(HT, HTIJ, "2D", "REAL", &
326	ids,ide, jds,jde, kds,kde, &
327	ims,ime, jms,jme, kms,kme )
328	! MEIXU end
329	#endif
330	!
331	! ************* END FDDA SETUP SECTION *************
332	!----------------------------------------------------------------------
333
334	!----------------------------------------------------------------------
335	! ************* BEGIN NUDGING SECTION *************
336
337	DO J = jts, jte
338	!
339	! IF NUDGING SURFACE WINDS IN THE BOUNDARY LAYER, IF IWINDS(INEST+2)=1
340	! USE A SIMILARITY CORRECTION BASED ON ROUGHNESS TO APPLY 10M
341	! WIND TO THE SURFACE LAYER (K=KL) AT 40M. TO DO THIS WE MUST
342	! STORE ROUGHNESS AND REGIME FOR EACH J SLICE AFTER THE CALL TO
343	! HIRPBL FOR LATER USE IN BLNUDGD.
344	!
345	!--- OBS NUDGING FOR TEMP AND MOISTURE
346	!
347	NSTA=NSTAT
348	IF(J .GT. 2 .and. J .LT. jde-1) THEN
349	IF(nudge_temp.EQ.1 .AND. NSTA.GT.0) &
350	THEN
351	! write(6,*) 'calling nudob: IVAR=3, J = ',j
352	CALL nudob(J, 3, t_tendf(ims,kms,j), &
353	inest, restart, ktau, fdob%ktaur, xtime, &
354	mut(ims,j), msftx(ims,j), msfty(ims,j), &
355	nobs_ndg_vars, nobs_err_flds, max_obs, maxdom, &
356	npfi, ionf, rinxy, fdob%window, &
357	fdob%nudge_t_pbl, &
358	fdob%sfcfact, fdob%sfcfacr, &
359	fdob%levidn, &
360	parid, nstat, &
361	fdob%rinfmn, fdob%rinfmx, fdob%pfree, &
362	fdob%dcon, fdob%tfaci, &
363	fdob%sfc_scheme_horiz, fdob%sfc_scheme_vert, &
364	fdob%max_sndng_gap, &
365	fdob%lev_in_ob, fdob%plfo, fdob%nlevs_ob, &
366	parent_grid_ratio, dx, dtmin, fdob%rio, fdob%rjo, &
367	fdob%rko, fdob%timeob, varobs, errf, &
368	pbase(ims,kms,j), ptop, pp(ims,kms,j), &
369	nudge_wind, nudge_temp, nudge_mois, &
370	coef_wind, coef_temp, coef_mois, &
371	savwt(1,ims,kms,j), kpbl(ims,j), 0, &
372	vih_t(its,j,1), vih_t(its,j,2), ht(ims,j), &
373	z(ims,kms,j), &
374	iprt_nudob, &
375	ids,ide, jds,jde, kds,kde, & ! domain dims
376	ims,ime, jms,jme, kms,kme, & ! memory dims
377	its,ite, jts,jte, kts,kte ) ! tile dims
378	! write(6,*) 'return from nudob: IVAR=3, J = ',j
379	ENDIF
380
381	IF(nudge_mois.EQ.1 .AND. NSTA.GT.0) &
382	THEN
383	! write(6,*) 'calling nudob: IVAR=4, J = ',j
384	CALL nudob(J, 4, moist_tend(ims,kms,j), &
385	inest, restart, ktau, fdob%ktaur, xtime, &
386	mut(ims,j), msftx(ims,j), msfty(ims,j), &
387	nobs_ndg_vars, nobs_err_flds, max_obs, maxdom, &
388	npfi, ionf, rinxy, fdob%window, &
389	fdob%nudge_q_pbl, &
390	fdob%sfcfact, fdob%sfcfacr, &
391	fdob%levidn, &
392	parid, nstat, &
393	fdob%rinfmn, fdob%rinfmx, fdob%pfree, &
394	fdob%dcon, fdob%tfaci, &
395	fdob%sfc_scheme_horiz, fdob%sfc_scheme_vert, &
396	fdob%max_sndng_gap, &
397	fdob%lev_in_ob, fdob%plfo, fdob%nlevs_ob, &
398	parent_grid_ratio, dx, dtmin, fdob%rio, fdob%rjo, &
399	fdob%rko, fdob%timeob, varobs, errf, &
400	pbase(ims,kms,j), ptop, pp(ims,kms,j), &
401	nudge_wind, nudge_temp, nudge_mois, &
402	coef_wind, coef_temp, coef_mois, &
403	savwt(1,ims,kms,j), kpbl(ims,j), 0, &
404	vih_q(its,j,1), vih_q(its,j,2), ht(ims,j), &
405	z(ims,kms,j), &
406	iprt_nudob, &
407	ids,ide, jds,jde, kds,kde, & ! domain dims
408	ims,ime, jms,jme, kms,kme, & ! memory dims
409	its,ite, jts,jte, kts,kte ) ! tile dims
410	! write(6,*) 'return from nudob: IVAR=4, J = ',j
411	ENDIF
412	ENDIF
413
414	IF(nudge_wind.EQ.1 .AND. NSTA.GT.0) &
415	THEN
416	! write(6,*) 'calling nudob: IVAR=1, J = ',j
417	CALL nudob(J, 1, ru_tendf(ims,kms,j), &
418	inest, restart, ktau, fdob%ktaur, xtime, &
419	muu(ims,j), msfux(ims,j), msfuy(ims,j), &
420	nobs_ndg_vars, nobs_err_flds, max_obs, maxdom, &
421	npfi, ionf, rinxy, fdob%window, &
422	fdob%nudge_uv_pbl, &
423	fdob%sfcfact, fdob%sfcfacr, &
424	fdob%levidn, &
425	parid, nstat, &
426	fdob%rinfmn, fdob%rinfmx, fdob%pfree, &
427	fdob%dcon, fdob%tfaci, &
428	fdob%sfc_scheme_horiz, fdob%sfc_scheme_vert, &
429	fdob%max_sndng_gap, &
430	fdob%lev_in_ob, fdob%plfo, fdob%nlevs_ob, &
431	parent_grid_ratio, dx, dtmin, fdob%rio, fdob%rjo, &
432	fdob%rko, fdob%timeob, varobs, errf, &
433	pbase(ims,kms,j), ptop, pp(ims,kms,j), &
434	nudge_wind, nudge_temp, nudge_mois, &
435	coef_wind, coef_temp, coef_mois, &
436	savwt(1,ims,kms,j), kpbl(ims,j), 0, &
437	vih_uv(its,j,1), vih_uv(its,j,2), ht(ims,j), &
438	z(ims,kms,j), &
439	iprt_nudob, &
440	ids,ide, jds,jde, kds,kde, & ! domain dims
441	ims,ime, jms,jme, kms,kme, & ! memory dims
442	its,ite, jts,jte, kts,kte ) ! tile dims
443	! write(6,*) 'return from nudob: IVAR=1, J = ',j
444
445	! write(6,*) 'calling nudob: IVAR=2, J = ',j
446	CALL nudob(J, 2, rv_tendf(ims,kms,j), &
447	inest, restart, ktau, fdob%ktaur, xtime, &
448	muv(ims,j), msfvx(ims,j), msfvy(ims,j), &
449	nobs_ndg_vars, nobs_err_flds, max_obs, maxdom, &
450	npfi, ionf, rinxy, fdob%window, &
451	fdob%nudge_uv_pbl, &
452	fdob%sfcfact, fdob%sfcfacr, &
453	fdob%levidn, &
454	parid, nstat, &
455	fdob%rinfmn, fdob%rinfmx, fdob%pfree, &
456	fdob%dcon, fdob%tfaci, &
457	fdob%sfc_scheme_horiz, fdob%sfc_scheme_vert, &
458	fdob%max_sndng_gap, &
459	fdob%lev_in_ob, fdob%plfo, fdob%nlevs_ob, &
460	parent_grid_ratio, dx, dtmin, fdob%rio, fdob%rjo, &
461	fdob%rko, fdob%timeob, varobs, errf, &
462	pbase(ims,kms,j), ptop, pp(ims,kms,j), &
463	nudge_wind, nudge_temp, nudge_mois, &
464	coef_wind, coef_temp, coef_mois, &
465	savwt(1,ims,kms,j), kpbl(ims,j), 0, &
466	vih_uv(its,j,1), vih_uv(its,j,2), ht(ims,j), &
467	z(ims,kms,j), &
468	iprt_nudob, &
469	ids,ide, jds,jde, kds,kde, & ! domain dims
470	ims,ime, jms,jme, kms,kme, & ! memory dims
471	its,ite, jts,jte, kts,kte ) ! tile dims
472	! write(6,*) 'return from nudob: IVAR=2, J = ',j
473	ENDIF
474	ENDDO
475	!
476	! --- END OF 4DDA
477	!
478	RETURN
479	#endif
480	END SUBROUTINE fddaobs_driver
481
482	SUBROUTINE compute_VIH(vif, hmax, fullmin, rampmin, &
483	regime, pblh, terrh, z, vih, &
484	ids,ide, jds,jde, kds,kde, & ! domain dims
485	ims,ime, jms,jme, kms,kme, &
486	its,ite, jts,jte, kts,kte)
487
488	USE module_fddaobs_rtfdda
489
490	IMPLICIT NONE
491	!*******************************************************************************
492	!*** COMPUTE HEIGHTS FOR SURFACE OBS VERTICAL INFLUENCE FUNCTION ***
493	!*******************************************************************************
494
495	REAL, INTENT(IN) :: vif(6) ! Vert infl params for regimes
496	REAL, INTENT(IN) :: hmax ! Max height to apply nudging
497	REAL, INTENT(IN) :: fullmin ! Min height of full nudging
498	REAL, INTENT(IN) :: rampmin ! Min height to ramp full-to-0
499	REAL, INTENT(IN) :: regime(ims:ime,jms:jme) ! Stability regime
500	REAL, INTENT(IN) :: pblh(ims:ime,jms:jme) ! PBL height (m)
501	REAL, INTENT(IN) :: terrh(ims:ime,jms:jme) ! Terrain ht (m)
502	REAL, INTENT(IN) :: z(ims:ime,kms:kme,jms:jme) ! Ht (m) above sl (half levs)
503	REAL, INTENT(OUT) :: vih(its:ite,jts:jte,2) ! Vt infl hts abv grd for LML obs
504	! INTEGER, INTENT(OUT) :: vik(its:ite,jts:jte,2) ! Vert infl k-levels for LML obs
505	INTEGER, INTENT(IN) :: ids,ide, jds,jde, kds,kde ! domain dims.
506	INTEGER, INTENT(IN) :: ims,ime, jms,jme, kms,kme ! memory dims.
507	INTEGER, INTENT(IN) :: its,ite, jts,jte, kts,kte ! tile dims.
508
509	! Local variables
510	real :: fullr(its:ite) ! Height up to full vertical weighting
511	real :: rampr(its:ite) ! Height of ramp-down to zero weighting
512	character(len=200) :: msg ! Argument to wrf_message
513	integer:: i, j ! Loop counters
514
515	integer k ! ajb test
516
517	! Do J-slabs
518	do j = jts, jte
519
520	! Set fullr and rampr values according to regimes
521	do i = its, ite
522
523	if(regime(i,j).eq.1.0) then ! REGIME 1
524	fullr(i) = vif(1)
525	rampr(i) = vif(2)
526	elseif(regime(i,j).eq.2.0) then ! REGIME 2
527	fullr(i) = vif(3)
528	rampr(i) = vif(4)
529	elseif(regime(i,j).eq.3.0 .or. regime(i,j).eq.4.0) then ! REGIME 4
530	fullr(i) = vif(5)
531	rampr(i) = vif(6)
532	else
533	write(msg,'(a,f5.1,2(a,i4))') 'Unknown regime type ', regime(i,j), &
534	' at grid coordinate i = ',i,' j = ',j
535	call wrf_message(msg)
536	call wrf_error_fatal ( 'fddaobs_driver: compute_VIH STOP' )
537
538	endif
539
540	enddo
541
542	! Get vert infl heights for LML obs, from fullr, rampr, and pblh
543	CALL get_vif_hts_slab(fullr, rampr, pblh(ims,j), &
544	hmax, fullmin, rampmin, &
545	vih(its,j,1), vih(its,j,2), &
546	ims,ime, its,ite)
547	enddo
548	END SUBROUTINE compute_VIH
549
550	SUBROUTINE get_vif_hts_slab(fullr, rampr, pblh, hmax, fullmin, rampmin, &
551	ht1, ht2, ims,ime, its,ite)
552	! Compute VIF heights
553
554	IMPLICIT NONE
555
556	REAL, INTENT(IN) :: fullr(its:ite) ! Height up to full vertical weighting
557	REAL, INTENT(IN) :: rampr(its:ite) ! Height of ramp-down to zero weighting
558	REAL, INTENT(IN) :: pblh(ims:ime) ! PBL height (m)
559	REAL, INTENT(IN) :: hmax ! Max height to apply nudging
560	REAL, INTENT(IN) :: fullmin ! Min height of full nudging
561	REAL, INTENT(IN) :: rampmin ! Min height to ramp full-to-0
562	REAL, INTENT(OUT) :: ht1(its:ite) ! Vert infl fcn height 1
563	REAL, INTENT(OUT) :: ht2(its:ite) ! Vert infl fcn height 2
564	INTEGER, INTENT(IN) :: ims,ime ! Memory dims.
565	INTEGER, INTENT(IN) :: its,ite ! Tile dims.
566
567	! Local variables
568	integer :: i
569
570	do i = its, ite
571
572	! Determine lower height (below which the VIF=1 for full weighting)
573	if(fullr(i).ge.0.0) then ! fullr is height from ground
574	ht1(i) = fullr(i)
575	else ! fullr is relative to pbl (-5000 bias)
576	ht1(i) = pblh(i) - (fullr(i)+5000.)
577	endif
578
579	! Height ht1 can be no smaller than fullmin
580	ht1(i) = max(fullmin,ht1(i))
581
582	! Determine upper height (to which the VIF ramps down to zero weighting)
583	! NOTE: Height of ramp-down (ht2-ht1) can be no smaller than rampmin
584
585	if(rampr(i).ge.0.0) then
586	! rampr is height from ht1
587	ht2(i) = ht1(i) + max(rampmin,rampr(i))
588	else
589	! rampr is relative to pbl (-5000 bias)
590	ht2(i) = max( ht1(i)+rampmin, pblh(i)-(rampr(i)+5000.) )
591	endif
592
593	! Apply hmax
594	ht1(i) = min(ht1(i), hmax-rampmin)
595	ht2(i) = min(ht2(i), hmax)
596	enddo
597	END SUBROUTINE get_vif_hts_slab
598
599	SUBROUTINE get_vik_slab( h, hlevs, ht, vik, ims,ime, kms,kme, its,ite, kts,kte )
600
601	! Compute VIK values from heights on j-slab
602
603	IMPLICIT NONE
604
605	REAL, INTENT(IN) :: h(its:ite) ! height (m) above ground, on j-slab
606	REAL, INTENT(IN) :: hlevs(ims:ime,kms:kme) ! hgt (m) abv grd at modl levs (slab)
607	REAL, INTENT(IN) :: ht(ims:ime) ! terrain height (m) (slab)
608	INTEGER, INTENT(OUT):: vik(its:ite) ! vert infl k levels (slab)
609	INTEGER, INTENT(IN) :: ims,ime, kms,kme ! memory dims
610	INTEGER, INTENT(IN) :: its,ite, kts,kte ! tile dims
611
612	! Local variables
613	integer :: i
614	integer :: k
615	real :: ht_ag(kts:kte)
616
617	do i = its, ite
618
619	! Get column of height-above-ground values for this i coord
620	do k = kts,kte
621	ht_ag(k) = hlevs(i,k) - ht(i)
622	enddo
623	! Get k levels that correspond to height values
624	vik(i) = ht_to_k( h(i), ht_ag, kts,kte )
625	enddo
626	END SUBROUTINE get_vik_slab
627
628	INTEGER FUNCTION ht_to_k( h, hlevs, kts,kte )
629	IMPLICIT NONE
630
631	REAL, INTENT(IN) :: h ! height value (m)
632	REAL, INTENT(IN) :: hlevs(kts:kte) ! model height levels
633	INTEGER, INTENT(IN) :: kts,kte ! tile dims
634
635	! Local variables
636	INTEGER :: k ! loop counter
637	INTEGER :: klo ! lower k bound
638
639	KLEVS: do k = kts, kte
640	klo = k-1
641	if(h .le. hlevs(k)) then
642	EXIT KLEVS
643	endif
644	enddo KLEVS
645	klo = max0(1,klo)
646	ht_to_k = min0(kte,klo)
647	RETURN
648	END FUNCTION ht_to_k
649
650	END MODULE module_fddaobs_driver

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

Download in other formats:

Original Format