chgres_cube  1.9.0
•All Data Structures Namespaces Files Functions Variables Pages
atmosphere.F90
Go to the documentation of this file.
1 
4 
19  module atmosphere
20 
21  use esmf
22 
23  use atmosphere_target_data, only : lev_target, levp1_target, nvcoord_target, &
24  vcoord_target, delp_target_grid, &
25  dzdt_target_grid, ps_target_grid, &
26  temp_target_grid, tracers_target_grid, &
27  u_s_target_grid, v_s_target_grid, &
28  u_w_target_grid, v_w_target_grid, &
29  zh_target_grid, qnwfa_climo_target_grid, &
30  qnifa_climo_target_grid
31 
32  use atm_input_data, only : lev_input, &
33  levp1_input, &
34  tracers_input_grid, &
35  dzdt_input_grid, &
36  ps_input_grid, &
37  wind_input_grid, &
38  temp_input_grid, &
39  pres_input_grid, &
40  terrain_input_grid, &
41  read_input_atm_data, &
42  cleanup_input_atm_data
43 
44  use model_grid, only : target_grid, &
45  latitude_s_target_grid, &
46  longitude_s_target_grid, &
47  latitude_w_target_grid, &
48  longitude_w_target_grid, &
49  terrain_target_grid
50 
51  use program_setup, only : vcoord_file_target_grid, &
52  wam_cold_start, &
53  cycle_year, cycle_mon, &
54  cycle_day, cycle_hour, &
55  regional, &
56  tracers, num_tracers, &
57  num_tracers_input, &
58  atm_weight_file, &
59  use_thomp_mp_climo
60 
61  use thompson_mp_climo_data, only : read_thomp_mp_climo_data, &
62  cleanup_thomp_mp_climo_input_data, &
63  qnifa_climo_input_grid, &
64  qnwfa_climo_input_grid, &
65  thomp_pres_climo_input_grid, &
66  lev_thomp_mp_climo
67 
68  use write_data, only : write_fv3_atm_header_netcdf, &
69  write_fv3_atm_bndy_data_netcdf, &
70  write_fv3_atm_data_netcdf
71 
72  use utilities, only : error_handler
73  implicit none
74 
75  private
76 
77  type(esmf_field) :: dzdt_b4adj_target_grid
78  type(esmf_field), allocatable :: tracers_b4adj_target_grid(:)
79  type(esmf_field) :: ps_b4adj_target_grid
80  type(esmf_field) :: pres_target_grid
81  type(esmf_field) :: pres_b4adj_target_grid
82  type(esmf_field) :: temp_b4adj_target_grid
83  type(esmf_field) :: terrain_interp_to_target_grid
84  type(esmf_field) :: wind_target_grid
85  type(esmf_field) :: wind_b4adj_target_grid
86  type(esmf_field) :: wind_s_target_grid
87  type(esmf_field) :: wind_w_target_grid
88 
89 ! Fields associated with thompson microphysics climatological tracers.
90 
91  type(esmf_field) :: qnifa_climo_b4adj_target_grid
93  type(esmf_field) :: qnwfa_climo_b4adj_target_grid
95  type(esmf_field) :: thomp_pres_climo_b4adj_target_grid
97 
98  public :: atmosphere_driver
99  public :: read_vcoord_info
100 
101  contains
102 
107  subroutine atmosphere_driver(localpet)
108 
109  use mpi
110 
111  implicit none
112 
113  integer, intent(in) :: localpet
114 
115  integer :: isrctermprocessing
116  integer :: rc, n
117 
118  type(esmf_regridmethod_flag) :: method
119  type(esmf_routehandle) :: regrid_bl
120 
121  real(esmf_kind_r8), parameter :: p0=101325.0
122  real(esmf_kind_r8), parameter :: rd = 287.058
123  real(esmf_kind_r8), parameter :: grav = 9.81
124  real(esmf_kind_r8), parameter :: lapse = -6.5e-03
125 
126  real(esmf_kind_r8), parameter :: exponent = rd*lapse/grav
127  real(esmf_kind_r8), parameter :: one_over_exponent = 1.0 / exponent
128 
129  real(esmf_kind_r8), pointer :: psptr(:,:), tempptr(:,:,:)
130 
131 !-----------------------------------------------------------------------------------
132 ! Read atmospheric fields on the input grid.
133 !-----------------------------------------------------------------------------------
134 
135  call read_input_atm_data(localpet)
136 
137 !-----------------------------------------------------------------------------------
138 ! Read vertical coordinate info for target grid.
139 !-----------------------------------------------------------------------------------
140 
141  call read_vcoord_info
142 
143 !-----------------------------------------------------------------------------------
144 ! Create target grid field objects to hold data before vertical adjustment.
145 !-----------------------------------------------------------------------------------
146 
147  call create_atm_b4adj_esmf_fields
148 
149 !-----------------------------------------------------------------------------------
150 ! Horizontally interpolate. If specified, use weights from file.
151 !-----------------------------------------------------------------------------------
152 
153  isrctermprocessing = 1
154 
155  if (trim(atm_weight_file) /= "NULL") then
156 
157  print*,"- CALL FieldSMMStore FOR ATMOSPHERIC FIELDS."
158 
159  call esmf_fieldsmmstore(temp_input_grid, &
160  temp_b4adj_target_grid, &
161  atm_weight_file, &
162  routehandle=regrid_bl, &
163  srctermprocessing=isrctermprocessing, rc=rc)
164  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
165  call error_handler("IN FieldSMMStore", rc)
166 
167  else
168 
169  print*,"- CALL FieldRegridStore FOR ATMOSPHERIC FIELDS."
170 
171  method=esmf_regridmethod_bilinear
172 
173  call esmf_fieldregridstore(temp_input_grid, &
174  temp_b4adj_target_grid, &
175  polemethod=esmf_polemethod_allavg, &
176  srctermprocessing=isrctermprocessing, &
177  routehandle=regrid_bl, &
178  regridmethod=method, rc=rc)
179  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
180  call error_handler("IN FieldRegridStore", rc)
181 
182  endif
183 
184  print*,"- CALL Field_Regrid FOR TEMPERATURE."
185  call esmf_fieldregrid(temp_input_grid, &
186  temp_b4adj_target_grid, &
187  routehandle=regrid_bl, &
188  termorderflag=esmf_termorder_srcseq, &
189  rc=rc)
190  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
191  call error_handler("IN FieldRegrid", rc)
192 
193  print*,"- CALL Field_Regrid FOR PRESSURE."
194  call esmf_fieldregrid(pres_input_grid, &
195  pres_b4adj_target_grid, &
196  routehandle=regrid_bl, &
197  termorderflag=esmf_termorder_srcseq, &
198  rc=rc)
199  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
200  call error_handler("IN FieldRegrid", rc)
201 
202  do n = 1, num_tracers_input
203  print*,"- CALL Field_Regrid FOR TRACER ", trim(tracers(n))
204  call esmf_fieldregrid(tracers_input_grid(n), &
205  tracers_b4adj_target_grid(n), &
206  routehandle=regrid_bl, &
207  termorderflag=esmf_termorder_srcseq, rc=rc)
208  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
209  call error_handler("IN FieldRegrid", rc)
210 
211  enddo
212 
213  print*,"- CALL Field_Regrid FOR VERTICAL VELOCITY."
214  call esmf_fieldregrid(dzdt_input_grid, &
215  dzdt_b4adj_target_grid, &
216  routehandle=regrid_bl, &
217  termorderflag=esmf_termorder_srcseq, rc=rc)
218  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
219  call error_handler("IN FieldRegrid", rc)
220 
221  nullify(tempptr)
222  print*,"- CALL FieldGet FOR INPUT GRID VERTICAL VEL."
223  call esmf_fieldget(dzdt_input_grid, &
224  farrayptr=tempptr, rc=rc)
225  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
226  call error_handler("IN FieldGet", rc)
227 
228  print*, "MIN MAX W INPUT = ", minval(tempptr), maxval(tempptr)
229 
230  nullify(tempptr)
231  print*,"- CALL FieldGet FOR VERTICAL VEL B4ADJ."
232  call esmf_fieldget(dzdt_b4adj_target_grid, &
233  farrayptr=tempptr, rc=rc)
234  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
235  call error_handler("IN FieldGet", rc)
236 
237  print*, "MIN MAX W B4ADJ = ", minval(tempptr), maxval(tempptr)
238 
239  nullify(psptr)
240  print*,"- CALL FieldGet FOR INPUT SURFACE PRESSURE."
241  call esmf_fieldget(ps_input_grid, &
242  farrayptr=psptr, rc=rc)
243  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
244  call error_handler("IN FieldGet", rc)
245 
246 !------------------------------------------------------------------------------------
247 ! Assume standard lapse rate when interpolating pressure (per Phil Pegion).
248 !------------------------------------------------------------------------------------
249 
250  psptr = (psptr/p0)**exponent
251 
252  print*,"- CALL Field_Regrid FOR SURFACE PRESSURE."
253  call esmf_fieldregrid(ps_input_grid, &
254  ps_b4adj_target_grid, &
255  routehandle=regrid_bl, &
256  termorderflag=esmf_termorder_srcseq, &
257  rc=rc)
258  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
259  call error_handler("IN FieldRegrid", rc)
260 
261  nullify(psptr)
262  print*,"- CALL FieldGet FOR INPUT SURFACE PRESSURE B4ADJ."
263  call esmf_fieldget(ps_b4adj_target_grid, &
264  farrayptr=psptr, rc=rc)
265  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
266  call error_handler("IN FieldGet", rc)
267 
268  psptr = p0 * psptr**one_over_exponent
269 
270  print*,"- CALL Field_Regrid FOR TERRAIN."
271  call esmf_fieldregrid(terrain_input_grid, &
272  terrain_interp_to_target_grid, &
273  routehandle=regrid_bl, &
274  termorderflag=esmf_termorder_srcseq, &
275  rc=rc)
276  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
277  call error_handler("IN FieldRegrid", rc)
278 
279  print*,"- CALL Field_Regrid FOR 3-D WIND."
280  call esmf_fieldregrid(wind_input_grid, &
281  wind_b4adj_target_grid, &
282  routehandle=regrid_bl, &
283  termorderflag=esmf_termorder_srcseq, rc=rc)
284  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
285  call error_handler("IN FieldRegrid", rc)
286 
287  print*,"- CALL FieldRegridRelease."
288  call esmf_fieldregridrelease(routehandle=regrid_bl, rc=rc)
289  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
290  call error_handler("IN FieldRegridRelease", rc)
291 
292 !-----------------------------------------------------------------------------------
293 ! Deallocate input fields.
294 !-----------------------------------------------------------------------------------
295 
296  call cleanup_input_atm_data
297 
298 !-----------------------------------------------------------------------------------
299 ! Create target grid field objects to hold data after vertical interpolation.
300 !-----------------------------------------------------------------------------------
301 
302  call create_atm_esmf_fields
303 
304 !-----------------------------------------------------------------------------------
305 ! Adjust surface pressure for terrain differences.
306 !-----------------------------------------------------------------------------------
307 
308  call newps(localpet)
309 
310 !-----------------------------------------------------------------------------------
311 ! Compute 3-d pressure based on adjusted surface pressure.
312 !-----------------------------------------------------------------------------------
313 
314  call newpr1(localpet)
315 
316 !-----------------------------------------------------------------------------------
317 ! Vertically interpolate.
318 !-----------------------------------------------------------------------------------
319 
320  call vintg
321 
322  if( wam_cold_start ) then
323  call vintg_wam (cycle_year,cycle_mon,cycle_day,cycle_hour)
324  endif
325 
326 !-----------------------------------------------------------------------------------
327 ! Compute height.
328 !-----------------------------------------------------------------------------------
329 
330  call compute_zh
331 
332 !-----------------------------------------------------------------------------------
333 ! Free up memory.
334 !-----------------------------------------------------------------------------------
335 
336  call cleanup_target_atm_b4adj_data
337 
338 !-----------------------------------------------------------------------------------
339 ! Interpolate winds to 'd' grid.
340 !-----------------------------------------------------------------------------------
341 
342  isrctermprocessing = 1
343  method=esmf_regridmethod_bilinear
344 
345  print*,"- CALL FieldRegridStore FOR 3D-WIND WEST EDGE."
346  call esmf_fieldregridstore(wind_target_grid, &
347  wind_w_target_grid, &
348  polemethod=esmf_polemethod_allavg, &
349  srctermprocessing=isrctermprocessing, &
350  routehandle=regrid_bl, &
351  extrapmethod=esmf_extrapmethod_nearest_stod, &
352  regridmethod=method, rc=rc)
353  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
354  call error_handler("IN FieldRegridStore", rc)
355 
356  print*,"- CALL Field_Regrid FOR 3-D WIND WEST EDGE."
357  call esmf_fieldregrid(wind_target_grid, &
358  wind_w_target_grid, &
359  routehandle=regrid_bl, &
360  termorderflag=esmf_termorder_srcseq, rc=rc)
361  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
362  call error_handler("IN FieldRegrid", rc)
363 
364  print*,"- CALL FieldRegridRelease."
365  call esmf_fieldregridrelease(routehandle=regrid_bl, rc=rc)
366  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
367  call error_handler("IN FieldRegridRelease", rc)
368 
369  isrctermprocessing = 1
370  method=esmf_regridmethod_bilinear
371 
372  print*,"- CALL FieldRegridStore FOR 3D-WIND SOUTH EDGE."
373  call esmf_fieldregridstore(wind_target_grid, &
374  wind_s_target_grid, &
375  polemethod=esmf_polemethod_allavg, &
376  srctermprocessing=isrctermprocessing, &
377  routehandle=regrid_bl, &
378  extrapmethod=esmf_extrapmethod_nearest_stod, &
379  regridmethod=method, rc=rc)
380  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
381  call error_handler("IN FieldRegridStore", rc)
382 
383  print*,"- CALL Field_Regrid FOR 3-D WIND SOUTH EDGE."
384  call esmf_fieldregrid(wind_target_grid, &
385  wind_s_target_grid, &
386  routehandle=regrid_bl, &
387  termorderflag=esmf_termorder_srcseq, rc=rc)
388  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
389  call error_handler("IN FieldRegrid", rc)
390 
391  print*,"- CALL FieldRegridRelease."
392  call esmf_fieldregridrelease(routehandle=regrid_bl, rc=rc)
393  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
394  call error_handler("IN FieldRegridRelease", rc)
395 
396 !-----------------------------------------------------------------------------------
397 ! Convert from 3-d to 2-d cartesian winds.
398 !-----------------------------------------------------------------------------------
399 
400  call convert_winds
401 
402 !-----------------------------------------------------------------------------------
403 ! If selected, process thompson microphysics climatological fields.
404 !-----------------------------------------------------------------------------------
405 
406  if (use_thomp_mp_climo) then
407  call read_thomp_mp_climo_data
408  call horiz_interp_thomp_mp_climo
409  call vintg_thomp_mp_climo
410  endif
411 
412 !-----------------------------------------------------------------------------------
413 ! Write target data to file.
414 !-----------------------------------------------------------------------------------
415 
416  call write_fv3_atm_header_netcdf(localpet)
417  if (regional <= 1) call write_fv3_atm_data_netcdf(localpet)
418  if (regional >= 1) call write_fv3_atm_bndy_data_netcdf(localpet)
419 
420 !-----------------------------------------------------------------------------------
421 ! Free up memory.
422 !-----------------------------------------------------------------------------------
423 
424  call cleanup_all_target_atm_data
425 
426  end subroutine atmosphere_driver
427 
433  subroutine create_atm_b4adj_esmf_fields
434 
435  implicit none
436 
437  integer :: rc, n
438 
439  allocate(tracers_b4adj_target_grid(num_tracers_input))
440 
441  do n = 1, num_tracers_input
442  print*,"- CALL FieldCreate FOR TARGET GRID TRACER BEFORE ADJUSTMENT ", trim(tracers(n))
443  tracers_b4adj_target_grid(n) = esmf_fieldcreate(target_grid, &
444  typekind=esmf_typekind_r8, &
445  staggerloc=esmf_staggerloc_center, &
446  ungriddedlbound=(/1/), &
447  ungriddedubound=(/lev_input/), rc=rc)
448  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
449  call error_handler("IN FieldCreate", rc)
450  enddo
451 
452  print*,"- CALL FieldCreate FOR TARGET GRID TEMPERATURE BEFORE ADJUSTMENT."
453  temp_b4adj_target_grid = esmf_fieldcreate(target_grid, &
454  typekind=esmf_typekind_r8, &
455  staggerloc=esmf_staggerloc_center, &
456  ungriddedlbound=(/1/), &
457  ungriddedubound=(/lev_input/), rc=rc)
458  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
459  call error_handler("IN FieldCreate", rc)
460 
461  print*,"- CALL FieldCreate FOR TARGET GRID PRESSURE BEFORE ADJUSTMENT."
462  pres_b4adj_target_grid = esmf_fieldcreate(target_grid, &
463  typekind=esmf_typekind_r8, &
464  staggerloc=esmf_staggerloc_center, &
465  ungriddedlbound=(/1/), &
466  ungriddedubound=(/lev_input/), rc=rc)
467  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
468  call error_handler("IN FieldCreate", rc)
469 
470  print*,"- CALL FieldCreate FOR TARGET GRID VERTICAL VELOCITY BEFORE ADJUSTMENT."
471  dzdt_b4adj_target_grid = esmf_fieldcreate(target_grid, &
472  typekind=esmf_typekind_r8, &
473  staggerloc=esmf_staggerloc_center, &
474  ungriddedlbound=(/1/), &
475  ungriddedubound=(/lev_input/), rc=rc)
476  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
477  call error_handler("IN FieldCreate", rc)
478 
479  print*,"- CALL FieldCreate FOR TARGET GRID UNSTAGGERED WINDS BEFORE ADJUSTMENT."
480  wind_b4adj_target_grid = esmf_fieldcreate(target_grid, &
481  typekind=esmf_typekind_r8, &
482  staggerloc=esmf_staggerloc_center, &
483  ungriddedlbound=(/1,1/), &
484  ungriddedubound=(/lev_input,3/), rc=rc)
485  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
486  call error_handler("IN FieldCreate", rc)
487 
488  print*,"- CALL FieldCreate FOR TARGET TERRAIN."
489  terrain_interp_to_target_grid = esmf_fieldcreate(target_grid, &
490  typekind=esmf_typekind_r8, &
491  staggerloc=esmf_staggerloc_center, rc=rc)
492  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
493  call error_handler("IN FieldCreate", rc)
494 
495  print*,"- CALL FieldCreate FOR TARGET SURFACE PRESSURE BEFORE ADJUSTMENT."
496  ps_b4adj_target_grid = esmf_fieldcreate(target_grid, &
497  typekind=esmf_typekind_r8, &
498  staggerloc=esmf_staggerloc_center, rc=rc)
499  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
500  call error_handler("IN FieldCreate", rc)
501 
502  end subroutine create_atm_b4adj_esmf_fields
503 
507  subroutine create_atm_esmf_fields
508 
509  implicit none
510 
511  integer :: rc, n
512 
513  allocate(tracers_target_grid(num_tracers))
514 
515  do n = 1, num_tracers
516  print*,"- CALL FieldCreate FOR TARGET GRID TRACERS ", trim(tracers(n))
517  tracers_target_grid(n) = esmf_fieldcreate(target_grid, &
518  typekind=esmf_typekind_r8, &
519  staggerloc=esmf_staggerloc_center, &
520  ungriddedlbound=(/1/), &
521  ungriddedubound=(/lev_target/), rc=rc)
522  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
523  call error_handler("IN FieldCreate", rc)
524  enddo
525 
526  print*,"- CALL FieldCreate FOR TARGET GRID TEMPERATURE."
527  temp_target_grid = esmf_fieldcreate(target_grid, &
528  typekind=esmf_typekind_r8, &
529  staggerloc=esmf_staggerloc_center, &
530  ungriddedlbound=(/1/), &
531  ungriddedubound=(/lev_target/), rc=rc)
532  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
533  call error_handler("IN FieldCreate", rc)
534 
535  print*,"- CALL FieldCreate FOR TARGET GRID PRESSURE."
536  pres_target_grid = esmf_fieldcreate(target_grid, &
537  typekind=esmf_typekind_r8, &
538  staggerloc=esmf_staggerloc_center, &
539  ungriddedlbound=(/1/), &
540  ungriddedubound=(/lev_target/), rc=rc)
541  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
542  call error_handler("IN FieldCreate", rc)
543 
544  print*,"- CALL FieldCreate FOR TARGET GRID VERTICAL VELOCITY."
545  dzdt_target_grid = esmf_fieldcreate(target_grid, &
546  typekind=esmf_typekind_r8, &
547  staggerloc=esmf_staggerloc_center, &
548  ungriddedlbound=(/1/), &
549  ungriddedubound=(/lev_target/), rc=rc)
550  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
551  call error_handler("IN FieldCreate", rc)
552 
553  print*,"- CALL FieldCreate FOR TARGET GRID DELP."
554  delp_target_grid = esmf_fieldcreate(target_grid, &
555  typekind=esmf_typekind_r8, &
556  staggerloc=esmf_staggerloc_center, &
557  ungriddedlbound=(/1/), &
558  ungriddedubound=(/lev_target/), rc=rc)
559  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
560  call error_handler("IN FieldCreate", rc)
561 
562  print*,"- CALL FieldCreate FOR TARGET HEIGHT."
563  zh_target_grid = esmf_fieldcreate(target_grid, &
564  typekind=esmf_typekind_r8, &
565  staggerloc=esmf_staggerloc_center, &
566  ungriddedlbound=(/1/), &
567  ungriddedubound=(/levp1_target/), rc=rc)
568  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
569  call error_handler("IN FieldCreate", rc)
570 
571  print*,"- CALL FieldCreate FOR TARGET UNSTAGGERED 3D-WIND."
572  wind_target_grid = esmf_fieldcreate(target_grid, &
573  typekind=esmf_typekind_r8, &
574  staggerloc=esmf_staggerloc_center, &
575  ungriddedlbound=(/1,1/), &
576  ungriddedubound=(/lev_target,3/), rc=rc)
577  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
578  call error_handler("IN FieldCreate", rc)
579 
580  print*,"- CALL FieldCreate FOR TARGET U_S."
581  u_s_target_grid = esmf_fieldcreate(target_grid, &
582  typekind=esmf_typekind_r8, &
583  staggerloc=esmf_staggerloc_edge2, &
584  ungriddedlbound=(/1/), &
585  ungriddedubound=(/lev_target/), rc=rc)
586  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
587  call error_handler("IN FieldCreate", rc)
588 
589  print*,"- CALL FieldCreate FOR TARGET V_S."
590  v_s_target_grid = esmf_fieldcreate(target_grid, &
591  typekind=esmf_typekind_r8, &
592  staggerloc=esmf_staggerloc_edge2, &
593  ungriddedlbound=(/1/), &
594  ungriddedubound=(/lev_target/), rc=rc)
595  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
596  call error_handler("IN FieldCreate", rc)
597 
598  print*,"- CALL FieldCreate FOR TARGET 3D-WIND_S."
599  wind_s_target_grid = esmf_fieldcreate(target_grid, &
600  typekind=esmf_typekind_r8, &
601  staggerloc=esmf_staggerloc_edge2, &
602  ungriddedlbound=(/1,1/), &
603  ungriddedubound=(/lev_target,3/), rc=rc)
604  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
605  call error_handler("IN FieldCreate", rc)
606 
607  print*,"- CALL FieldCreate FOR TARGET U_W."
608  u_w_target_grid = esmf_fieldcreate(target_grid, &
609  typekind=esmf_typekind_r8, &
610  staggerloc=esmf_staggerloc_edge1, &
611  ungriddedlbound=(/1/), &
612  ungriddedubound=(/lev_target/), rc=rc)
613  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
614  call error_handler("IN FieldCreate", rc)
615 
616  print*,"- CALL FieldCreate FOR TARGET V_W."
617  v_w_target_grid = esmf_fieldcreate(target_grid, &
618  typekind=esmf_typekind_r8, &
619  staggerloc=esmf_staggerloc_edge1, &
620  ungriddedlbound=(/1/), &
621  ungriddedubound=(/lev_target/), rc=rc)
622  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
623  call error_handler("IN FieldCreate", rc)
624 
625  print*,"- CALL FieldCreate FOR TARGET 3D-WIND_W."
626  wind_w_target_grid = esmf_fieldcreate(target_grid, &
627  typekind=esmf_typekind_r8, &
628  staggerloc=esmf_staggerloc_edge1, &
629  ungriddedlbound=(/1,1/), &
630  ungriddedubound=(/lev_target,3/), rc=rc)
631  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
632  call error_handler("IN FieldCreate", rc)
633 
634  print*,"- CALL FieldCreate FOR TARGET SURFACE PRESSURE."
635  ps_target_grid = esmf_fieldcreate(target_grid, &
636  typekind=esmf_typekind_r8, &
637  staggerloc=esmf_staggerloc_center, rc=rc)
638  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
639  call error_handler("IN FieldCreate", rc)
640 
641  end subroutine create_atm_esmf_fields
642 
646  subroutine convert_winds
647 
648  implicit none
649 
650  integer :: clb(4), cub(4)
651  integer :: i, j, k, rc
652 
653  real(esmf_kind_r8), pointer :: latptr(:,:)
654  real(esmf_kind_r8), pointer :: lonptr(:,:)
655  real(esmf_kind_r8), pointer :: uptr(:,:,:)
656  real(esmf_kind_r8), pointer :: vptr(:,:,:)
657  real(esmf_kind_r8), pointer :: windptr(:,:,:,:)
658  real(esmf_kind_r8) :: latrad, lonrad
659 
660 !-----------------------------------------------------------------------------------
661 ! Convert from 3-d cartesian to 2-cartesian winds
662 !-----------------------------------------------------------------------------------
663 
664  print*,'- CONVERT WINDS.'
665 
666  print*,"- CALL FieldGet FOR 3-D WIND_S."
667  call esmf_fieldget(wind_s_target_grid, &
668  computationallbound=clb, &
669  computationalubound=cub, &
670  farrayptr=windptr, rc=rc)
671  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
672  call error_handler("IN FieldGet", rc)
673 
674  print*,"- CALL FieldGet FOR U_S."
675  call esmf_fieldget(u_s_target_grid, &
676  farrayptr=uptr, rc=rc)
677  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
678  call error_handler("IN FieldGet", rc)
679 
680  print*,"- CALL FieldGet FOR V_S."
681  call esmf_fieldget(v_s_target_grid, &
682  farrayptr=vptr, rc=rc)
683  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
684  call error_handler("IN FieldGet", rc)
685 
686  print*,"- CALL FieldGet FOR LATITUDE_S."
687  call esmf_fieldget(latitude_s_target_grid, &
688  farrayptr=latptr, rc=rc)
689  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
690  call error_handler("IN FieldGet", rc)
691 
692  print*,"- CALL FieldGet FOR LONGITUDE_S."
693  call esmf_fieldget(longitude_s_target_grid, &
694  farrayptr=lonptr, rc=rc)
695  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
696  call error_handler("IN FieldGet", rc)
697 
698  do i = clb(1), cub(1)
699  do j = clb(2), cub(2)
700  latrad = latptr(i,j) * acos(-1.) / 180.0
701  lonrad = lonptr(i,j) * acos(-1.) / 180.0
702  do k = clb(3), cub(3)
703  uptr(i,j,k) = windptr(i,j,k,1) * cos(lonrad) + windptr(i,j,k,2) * sin(lonrad)
704  vptr(i,j,k) = -windptr(i,j,k,1) * sin(latrad) * sin(lonrad) + &
705  windptr(i,j,k,2) * sin(latrad) * cos(lonrad) + &
706  windptr(i,j,k,3) * cos(latrad)
707  enddo
708  enddo
709  enddo
710 
711 
712  print*,"- CALL FieldGet FOR 3-D WIND_W."
713  call esmf_fieldget(wind_w_target_grid, &
714  computationallbound=clb, &
715  computationalubound=cub, &
716  farrayptr=windptr, rc=rc)
717  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
718  call error_handler("IN FieldGet", rc)
719 
720  print*,"- CALL FieldGet FOR U_W."
721  call esmf_fieldget(u_w_target_grid, &
722  farrayptr=uptr, rc=rc)
723  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
724  call error_handler("IN FieldGet", rc)
725 
726  print*,"- CALL FieldGet FOR V_W."
727  call esmf_fieldget(v_w_target_grid, &
728  farrayptr=vptr, rc=rc)
729  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
730  call error_handler("IN FieldGet", rc)
731 
732  print*,"- CALL FieldGet FOR LATITUDE_W."
733  call esmf_fieldget(latitude_w_target_grid, &
734  farrayptr=latptr, rc=rc)
735  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
736  call error_handler("IN FieldGet", rc)
737 
738  print*,"- CALL FieldGet FOR LONGITUDE_W."
739  call esmf_fieldget(longitude_w_target_grid, &
740  farrayptr=lonptr, rc=rc)
741  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
742  call error_handler("IN FieldGet", rc)
743 
744  do i = clb(1), cub(1)
745  do j = clb(2), cub(2)
746  latrad = latptr(i,j) * acos(-1.) / 180.0
747  lonrad = lonptr(i,j) * acos(-1.) / 180.0
748  do k = clb(3), cub(3)
749  uptr(i,j,k) = windptr(i,j,k,1) * cos(lonrad) + windptr(i,j,k,2) * sin(lonrad)
750  vptr(i,j,k) = -windptr(i,j,k,1) * sin(latrad) * sin(lonrad) + &
751  windptr(i,j,k,2) * sin(latrad) * cos(lonrad) + &
752  windptr(i,j,k,3) * cos(latrad)
753  enddo
754  enddo
755  enddo
756 
757  end subroutine convert_winds
758 
790  subroutine newpr1(localpet)
791  implicit none
792 
793  integer, intent(in) :: localpet
794 
795  integer :: idsl, idvc, rc
796  integer :: i, j, k, clb(3), cub(3)
797 
798  real(esmf_kind_r8), parameter :: rd=287.05
799  real(esmf_kind_r8), parameter :: cp=1004.6
800  real(esmf_kind_r8), parameter :: rocp=rd/cp
801  real(esmf_kind_r8), parameter :: rocp1=rocp+1
802  real(esmf_kind_r8), parameter :: rocpr=1/rocp
803 
804  real(esmf_kind_r8), pointer :: delp_ptr(:,:,:)
805  real(esmf_kind_r8), pointer :: pptr(:,:,:) ! adjusted 3-d p.
806  real(esmf_kind_r8), pointer :: psptr(:,:) ! adjusted surface p.
807  real(esmf_kind_r8) :: ak, bk
808  real(esmf_kind_r8), allocatable :: pi(:,:,:)
809 
810  print*,"COMPUTE 3-D PRESSURE FROM ADJUSTED SURFACE PRESSURE."
811 
812  idvc = 2 ! hard wire for now.
813  idsl = 2 ! hard wire for now.
814 
815  print*,"- CALL FieldGet FOR 3-D PRES."
816  call esmf_fieldget(pres_target_grid, &
817  computationallbound=clb, &
818  computationalubound=cub, &
819  farrayptr=pptr, rc=rc)
820  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
821  call error_handler("IN FieldGet", rc)
822 
823  print*,"- CALL FieldGet FOR DELP."
824  call esmf_fieldget(delp_target_grid, &
825  computationallbound=clb, &
826  computationalubound=cub, &
827  farrayptr=delp_ptr, rc=rc)
828  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
829  call error_handler("IN FieldGet", rc)
830 
831  print*,"- CALL FieldGet FOR SURFACE PRESSURE AFTER ADJUSTMENT"
832  call esmf_fieldget(ps_target_grid, &
833  farrayptr=psptr, rc=rc)
834  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
835  call error_handler("IN FieldGet", rc)
836 
837  allocate(pi(clb(1):cub(1),clb(2):cub(2),1:levp1_target))
838 
839  if(idvc.eq.2) then
840  do k=1,levp1_target
841  ak = vcoord_target(k,1)
842  bk = vcoord_target(k,2)
843  do i= clb(1), cub(1)
844  do j= clb(2), cub(2)
845  pi(i,j,k) = ak + bk*psptr(i,j)
846  enddo
847  enddo
848  enddo
849  do k=1,lev_target
850  do i= clb(1), cub(1)
851  do j= clb(2), cub(2)
852  delp_ptr(i,j,k) = pi(i,j,k) - pi(i,j,k+1)
853  enddo
854  enddo
855  enddo
856  else
857  call error_handler("PROGRAM ONLY WORKS WITH IDVC 2", 1)
858  endif
859 
860  if(idsl.eq.2) then
861  do k=1,lev_target
862  do i= clb(1), cub(1)
863  do j= clb(2), cub(2)
864  pptr(i,j,k) = (pi(i,j,k)+pi(i,j,k+1))/2.0
865  enddo
866  enddo
867  enddo
868  else
869  do k=1,lev_target
870  do i= clb(1), cub(1)
871  do j= clb(2), cub(2)
872  pptr(i,j,k) = ((pi(i,j,k)**rocp1-pi(i,j,k+1)**rocp1)/ &
873  (rocp1*(pi(i,j,k)-pi(i,j,k+1))))**rocpr
874  enddo
875  enddo
876  enddo
877  endif
878 
879  deallocate(pi)
880 
881  if (localpet == 0) then
882  print*,'new pres ',pptr(clb(1),clb(2),:)
883  print*,'delp ',delp_ptr(clb(1),clb(2),:)
884  endif
885 
886  end subroutine newpr1
887 
901  subroutine newps(localpet)
902 
903  implicit none
904 
905  integer, intent(in) :: localpet
906  integer :: i, j, k, ii
907  integer :: clb(3), cub(3), ls, rc
908 
909  real(esmf_kind_r8), pointer :: pptr(:,:,:)
910  real(esmf_kind_r8), pointer :: psptr(:,:)
911  real(esmf_kind_r8), pointer :: psnewptr(:,:) ! adjusted surface p.
912  real(esmf_kind_r8), pointer :: tptr(:,:,:)
913  real(esmf_kind_r8), pointer :: qptr(:,:,:)
914  real(esmf_kind_r8), pointer :: zsptr(:,:)
915  real(esmf_kind_r8), pointer :: zsnewptr(:,:)
916  real(esmf_kind_r8), allocatable :: zu(:,:)
917  real(esmf_kind_r8), parameter :: beta=-6.5e-3
918  real(esmf_kind_r8), parameter :: epsilon=1.e-9
919  real(esmf_kind_r8), parameter :: g=9.80665
920  real(esmf_kind_r8), parameter :: rd=287.05
921  real(esmf_kind_r8), parameter :: rv=461.50
922  real(esmf_kind_r8), parameter :: gor=g/rd
923  real(esmf_kind_r8), parameter :: fv=rv/rd-1.
924  real(esmf_kind_r8) :: ftv, fgam, apu, fz0
925  real(esmf_kind_r8) :: atvu, atv, fz1, fp0
926  real(esmf_kind_r8) :: apd, azd, agam, azu
927  real(esmf_kind_r8) :: atvd, fp1, gamma, pu
928  real(esmf_kind_r8) :: tvu, pd, tvd
929  real(esmf_kind_r8) :: at, aq, ap, az
930 
931  ftv(at,aq)=at*(1+fv*aq)
932  fgam(apu,atvu,apd,atvd)=-gor*log(atvd/atvu)/log(apd/apu)
933  fz0(ap,atv,azd,apd)=azd+atv/gor*log(apd/ap)
934  fz1(ap,atv,azd,apd,agam)=azd-atv/agam*((apd/ap)**(-agam/gor)-1)
935  fp0(az,azu,apu,atvu)=apu*exp(-gor/atvu*(az-azu))
936  fp1(az,azu,apu,atvu,agam)=apu*(1+agam/atvu*(az-azu))**(-gor/agam)
937 
938  print*,"- ADJUST SURFACE PRESSURE FOR NEW TERRAIN."
939 
940  print*,"- CALL FieldGet FOR 3-D PRES."
941  call esmf_fieldget(pres_b4adj_target_grid, &
942  computationallbound=clb, &
943  computationalubound=cub, &
944  farrayptr=pptr, rc=rc)
945  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
946  call error_handler("IN FieldGet", rc)
947 
948  if(localpet==0) then
949  print*,'old pres ',pptr(clb(1),clb(2),:)
950  endif
951 
952  print*,"- CALL FieldGet FOR TEMPERATURE"
953  call esmf_fieldget(temp_b4adj_target_grid, &
954  farrayptr=tptr, rc=rc)
955  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
956  call error_handler("IN FieldGet", rc)
957 
958 ! Find specific humidity in the array of tracer fields.
959 
960  do ii = 1, num_tracers
961  if (trim(tracers(ii)) == "sphum") exit
962  enddo
963 
964  print*,"- CALL FieldGet FOR SPECIFIC HUMIDITY"
965  call esmf_fieldget(tracers_b4adj_target_grid(ii), &
966  farrayptr=qptr, rc=rc)
967  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
968  call error_handler("IN FieldGet", rc)
969 
970  print*,"- CALL FieldGet FOR SURFACE PRESSURE BEFORE ADJUSTMENT"
971  call esmf_fieldget(ps_b4adj_target_grid, &
972  farrayptr=psptr, rc=rc)
973  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
974  call error_handler("IN FieldGet", rc)
975 
976  print*,"- CALL FieldGet FOR SURFACE PRESSURE AFTER ADJUSTMENT"
977  call esmf_fieldget(ps_target_grid, &
978  farrayptr=psnewptr, rc=rc)
979  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
980  call error_handler("IN FieldGet", rc)
981 
982  print*,"- CALL FieldGet FOR OLD TERRAIN"
983  call esmf_fieldget(terrain_interp_to_target_grid, &
984  farrayptr=zsptr, rc=rc)
985  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
986  call error_handler("IN FieldGet", rc)
987 
988  print*,"- CALL FieldGet FOR NEW TERRAIN"
989  call esmf_fieldget(terrain_target_grid, &
990  farrayptr=zsnewptr, rc=rc)
991  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
992  call error_handler("IN FieldGet", rc)
993 
994  allocate(zu(clb(1):cub(1),clb(2):cub(2)))
995 
996 !-----------------------------------------------------------------------------------
997 ! Note, this routine was adapted from the spectral GFS which labeled the lowest
998 ! model layer as '1'.
999 !-----------------------------------------------------------------------------------
1000 
1001 !-----------------------------------------------------------------------------------
1002 ! Compute surface pressure below the original ground.
1003 !-----------------------------------------------------------------------------------
1004 
1005  ls=0
1006  k=1
1007  gamma=beta
1008  do i=clb(1), cub(1)
1009  do j=clb(2), cub(2)
1010  pu=pptr(i,j,k)
1011  tvu=ftv(tptr(i,j,k),qptr(i,j,k))
1012  zu(i,j)=fz1(pu,tvu,zsptr(i,j),psptr(i,j),gamma)
1013  if(zsnewptr(i,j).le.zu(i,j)) then
1014  pu=pptr(i,j,k)
1015  tvu=ftv(tptr(i,j,k),qptr(i,j,k))
1016  if(abs(gamma).gt.epsilon) then
1017  psnewptr(i,j)=fp1(zsnewptr(i,j),zu(i,j),pu,tvu,gamma)
1018  else
1019  psnewptr(i,j)=fp0(zsnewptr(i,j),zu(i,j),pu,tvu)
1020  endif
1021  else
1022  psnewptr(i,j)=0
1023  ls=ls+1
1024  endif
1025  enddo
1026  enddo
1027 
1028 !-----------------------------------------------------------------------------------
1029 ! Compute surface pressure above the original ground.
1030 !-----------------------------------------------------------------------------------
1031 
1032  do k=2,cub(3)
1033  if(ls.gt.0) then
1034  do i=clb(1),cub(1)
1035  do j=clb(2),cub(2)
1036  if(psnewptr(i,j).eq.0) then
1037  pu=pptr(i,j,k)
1038  tvu=ftv(tptr(i,j,k),qptr(i,j,k))
1039  pd=pptr(i,j,k-1)
1040  tvd=ftv(tptr(i,j,k-1),qptr(i,j,k-1))
1041  gamma=fgam(pu,tvu,pd,tvd)
1042  if(abs(gamma).gt.epsilon) then
1043  zu(i,j)=fz1(pu,tvu,zu(i,j),pd,gamma)
1044  else
1045  zu(i,j)=fz0(pu,tvu,zu(i,j),pd)
1046  endif
1047  if(zsnewptr(i,j).le.zu(i,j)) then
1048  if(abs(gamma).gt.epsilon) then
1049  psnewptr(i,j)=fp1(zsnewptr(i,j),zu(i,j),pu,tvu,gamma)
1050  else
1051  psnewptr(i,j)=fp0(zsnewptr(i,j),zu(i,j),pu,tvu)
1052  endif
1053  ls=ls-1
1054  endif
1055  endif
1056  enddo
1057  enddo
1058  endif
1059  enddo
1060 
1061 !-----------------------------------------------------------------------------------
1062 ! Compute surface pressure over the top.
1063 !-----------------------------------------------------------------------------------
1064 
1065 
1066  if(ls.gt.0) then
1067  k=cub(3)
1068  gamma=0
1069  do i=clb(1),cub(1)
1070  do j=clb(2),cub(2)
1071  if(psnewptr(i,j).eq.0) then
1072  pu=pptr(i,j,k)
1073  tvu=ftv(tptr(i,j,k),qptr(i,j,k))
1074  psnewptr(i,j)=fp0(zsnewptr(i,j),zu(i,j),pu,tvu)
1075  endif
1076  enddo
1077  enddo
1078  endif
1079 
1080  deallocate(zu)
1081 
1082  if (localpet == 0) then
1083 ! do i=clb(1),cub(1)
1084 ! do j=clb(2),cub(2)
1085  do i=clb(1),clb(1)
1086  do j=clb(2),clb(2)
1087  print*,'sfcp adjust ',(zsnewptr(i,j)-zsptr(i,j)), psptr(i,j),psnewptr(i,j)
1088  enddo
1089  enddo
1090  endif
1091 
1092  end subroutine newps
1093 
1098  subroutine read_vcoord_info
1099  implicit none
1100 
1101  integer :: istat, n, k
1102 
1103  print*
1104  print*,"OPEN VERTICAL COORD FILE: ", trim(vcoord_file_target_grid)
1105  open(14, file=trim(vcoord_file_target_grid), form='formatted', iostat=istat)
1106  if (istat /= 0) then
1107  call error_handler("OPENING VERTICAL COORD FILE", istat)
1108  endif
1109 
1110  read(14, *, iostat=istat) nvcoord_target, lev_target
1111  if (istat /= 0) then
1112  call error_handler("READING VERTICAL COORD FILE", istat)
1113  endif
1114 
1115  levp1_target = lev_target + 1
1116 
1117  allocate(vcoord_target(levp1_target, nvcoord_target))
1118  read(14, *, iostat=istat) ((vcoord_target(n,k), k=1,nvcoord_target), n=1,levp1_target)
1119  if (istat /= 0) then
1120  call error_handler("READING VERTICAL COORD FILE", istat)
1121  endif
1122 
1123  print*
1124 
1125  close(14)
1126 
1127  end subroutine read_vcoord_info
1128 
1133  subroutine horiz_interp_thomp_mp_climo
1134 
1135  implicit none
1136 
1137  integer :: isrctermprocessing, rc
1138 
1139  type(esmf_regridmethod_flag) :: method
1140  type(esmf_routehandle) :: regrid_bl
1141 
1142  isrctermprocessing=1
1143 
1144  print*,"- CALL FieldCreate FOR TARGET GRID THOMP CLIMO QNIFA BEFORE ADJUSTMENT."
1145  qnifa_climo_b4adj_target_grid = esmf_fieldcreate(target_grid, &
1146  typekind=esmf_typekind_r8, &
1147  staggerloc=esmf_staggerloc_center, &
1148  ungriddedlbound=(/1/), &
1149  ungriddedubound=(/lev_thomp_mp_climo/), rc=rc)
1150  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1151  call error_handler("IN FieldCreate", rc)
1152 
1153  print*,"- CALL FieldCreate FOR TARGET GRID THOMP CLIMO QNWFA BEFORE ADJUSTMENT."
1154  qnwfa_climo_b4adj_target_grid = esmf_fieldcreate(target_grid, &
1155  typekind=esmf_typekind_r8, &
1156  staggerloc=esmf_staggerloc_center, &
1157  ungriddedlbound=(/1/), &
1158  ungriddedubound=(/lev_thomp_mp_climo/), rc=rc)
1159  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1160  call error_handler("IN FieldCreate", rc)
1161 
1162  print*,"- CALL FieldCreate FOR TARGET GRID THOMP CLIMO PRESSURE BEFORE ADJUSTMENT."
1163  thomp_pres_climo_b4adj_target_grid = esmf_fieldcreate(target_grid, &
1164  typekind=esmf_typekind_r8, &
1165  staggerloc=esmf_staggerloc_center, &
1166  ungriddedlbound=(/1/), &
1167  ungriddedubound=(/lev_thomp_mp_climo/), rc=rc)
1168  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1169  call error_handler("IN FieldCreate", rc)
1170 
1171  print*,"- CALL FieldCreate FOR TARGET GRID THOMP CLIMO QNIFA."
1172  qnifa_climo_target_grid = esmf_fieldcreate(target_grid, &
1173  typekind=esmf_typekind_r8, &
1174  staggerloc=esmf_staggerloc_center, &
1175  ungriddedlbound=(/1/), &
1176  ungriddedubound=(/lev_target/), rc=rc)
1177  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1178  call error_handler("IN FieldCreate", rc)
1179 
1180  print*,"- CALL FieldCreate FOR TARGET GRID THOMP CLIMO QNWFA."
1181  qnwfa_climo_target_grid = esmf_fieldcreate(target_grid, &
1182  typekind=esmf_typekind_r8, &
1183  staggerloc=esmf_staggerloc_center, &
1184  ungriddedlbound=(/1/), &
1185  ungriddedubound=(/lev_target/), rc=rc)
1186  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1187  call error_handler("IN FieldCreate", rc)
1188 
1189  print*,"- CALL FieldRegridStore FOR THOMPSON CLIMO FIELDS."
1190 
1191  method=esmf_regridmethod_bilinear
1192 
1193  call esmf_fieldregridstore(qnifa_climo_input_grid, &
1194  qnifa_climo_b4adj_target_grid, &
1195  polemethod=esmf_polemethod_allavg, &
1196  srctermprocessing=isrctermprocessing, &
1197  routehandle=regrid_bl, &
1198  regridmethod=method, rc=rc)
1199  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1200  call error_handler("IN FieldRegridStore", rc)
1201 
1202  print*,"- CALL Field_Regrid FOR THOMP CLIMO QNIFA."
1203  call esmf_fieldregrid(qnifa_climo_input_grid, &
1204  qnifa_climo_b4adj_target_grid, &
1205  routehandle=regrid_bl, &
1206  termorderflag=esmf_termorder_srcseq, rc=rc)
1207  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1208  call error_handler("IN FieldRegrid", rc)
1209 
1210  print*,"- CALL Field_Regrid FOR THOMP CLIMO QNWFA."
1211  call esmf_fieldregrid(qnwfa_climo_input_grid, &
1212  qnwfa_climo_b4adj_target_grid, &
1213  routehandle=regrid_bl, &
1214  termorderflag=esmf_termorder_srcseq, rc=rc)
1215  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1216  call error_handler("IN FieldRegrid", rc)
1217 
1218  print*,"- CALL Field_Regrid FOR THOMP PRESSURE."
1219  call esmf_fieldregrid(thomp_pres_climo_input_grid, &
1220  thomp_pres_climo_b4adj_target_grid, &
1221  routehandle=regrid_bl, &
1222  termorderflag=esmf_termorder_srcseq, rc=rc)
1223  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1224  call error_handler("IN FieldRegrid", rc)
1225 
1226  print*,"- CALL FieldRegridRelease."
1227  call esmf_fieldregridrelease(routehandle=regrid_bl, rc=rc)
1228  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1229  call error_handler("IN FieldRegridRelease", rc)
1230 
1231 !-----------------------------------------------------------------------------------
1232 ! Free up input data memory.
1233 !-----------------------------------------------------------------------------------
1234 
1235  call cleanup_thomp_mp_climo_input_data
1236 
1237  end subroutine horiz_interp_thomp_mp_climo
1238 
1245  SUBROUTINE vintg_thomp_mp_climo
1246 
1247  implicit none
1248 
1249  INTEGER :: CLB(3), CUB(3), RC
1250  INTEGER :: IM, KM1, KM2, NT
1251  INTEGER :: I, J, K
1252 
1253  REAL(ESMF_KIND_R8), ALLOCATABLE :: Z1(:,:,:), Z2(:,:,:)
1254  REAL(ESMF_KIND_R8), ALLOCATABLE :: C1(:,:,:,:),C2(:,:,:,:)
1255 
1256  REAL(ESMF_KIND_R8), POINTER :: QNIFA1PTR(:,:,:) ! input
1257  REAL(ESMF_KIND_R8), POINTER :: QNIFA2PTR(:,:,:) ! target
1258  REAL(ESMF_KIND_R8), POINTER :: QNWFA1PTR(:,:,:) ! input
1259  REAL(ESMF_KIND_R8), POINTER :: QNWFA2PTR(:,:,:) ! target
1260  REAL(ESMF_KIND_R8), POINTER :: P1PTR(:,:,:) ! input pressure
1261  REAL(ESMF_KIND_R8), POINTER :: P2PTR(:,:,:) ! target pressure
1262 
1263  print*,"- VERTICALY INTERPOLATE THOMP MP CLIMO TRACERS."
1264 
1265  print*,"- CALL FieldGet FOR 3-D THOMP PRES."
1266  call esmf_fieldget(thomp_pres_climo_b4adj_target_grid, &
1267  computationallbound=clb, &
1268  computationalubound=cub, &
1269  farrayptr=p1ptr, rc=rc)
1270  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1271  call error_handler("IN FieldGet", rc)
1272 
1273 ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1274 ! The '1'/'2' arrays hold fields before/after interpolation.
1275 ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1276 
1277  nt= 2 ! number of thomp tracers
1278 
1279  ALLOCATE(z1(clb(1):cub(1),clb(2):cub(2),lev_thomp_mp_climo))
1280  ALLOCATE(z2(clb(1):cub(1),clb(2):cub(2),lev_target))
1281  ALLOCATE(c1(clb(1):cub(1),clb(2):cub(2),lev_thomp_mp_climo,nt))
1282  ALLOCATE(c2(clb(1):cub(1),clb(2):cub(2),lev_target,nt))
1283 
1284  z1 = -log(p1ptr)
1285 
1286  print*,"- CALL FieldGet FOR 3-D ADJUSTED PRESS"
1287  call esmf_fieldget(pres_target_grid, &
1288  farrayptr=p2ptr, rc=rc)
1289  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1290  call error_handler("IN FieldGet", rc)
1291 
1292  z2 = -log(p2ptr)
1293 
1294 !print*,'pres check 1 ', p1ptr(clb(1),clb(2),:)
1295 !print*,'pres check 2 ', p2ptr(clb(1),clb(2),:)
1296 
1297  print*,"- CALL FieldGet FOR qnifa before vertical adjustment."
1298  call esmf_fieldget(qnifa_climo_b4adj_target_grid, &
1299  farrayptr=qnifa1ptr, rc=rc)
1300  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1301  call error_handler("IN FieldGet", rc)
1302 
1303  c1(:,:,:,1) = qnifa1ptr(:,:,:)
1304 
1305  print*,"- CALL FieldGet FOR qnwfa before vertical adjustment."
1306  call esmf_fieldget(qnwfa_climo_b4adj_target_grid, &
1307  farrayptr=qnwfa1ptr, rc=rc)
1308  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1309  call error_handler("IN FieldGet", rc)
1310 
1311  c1(:,:,:,2) = qnwfa1ptr(:,:,:)
1312 
1313 ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1314 ! PERFORM LAGRANGIAN ONE-DIMENSIONAL INTERPOLATION
1315 ! THAT IS 4TH-ORDER IN INTERIOR, 2ND-ORDER IN OUTSIDE INTERVALS
1316 ! AND 1ST-ORDER FOR EXTRAPOLATION.
1317 ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1318 
1319  im = (cub(1)-clb(1)+1) * (cub(2)-clb(2)+1)
1320  km1= lev_thomp_mp_climo
1321  km2= lev_target
1322 
1323  CALL terp3(im,1,1,1,1,nt,(im*km1),(im*km2), &
1324  km1,im,im,z1,c1,km2,im,im,z2,c2)
1325 
1326  print*,"- CALL FieldGet FOR ADJUSTED climo qnifa."
1327  call esmf_fieldget(qnifa_climo_target_grid, &
1328  farrayptr=qnifa2ptr, rc=rc)
1329  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1330  call error_handler("IN FieldGet", rc)
1331 
1332  print*,"- CALL FieldGet FOR ADJUSTED climo qnwfa."
1333  call esmf_fieldget(qnwfa_climo_target_grid, &
1334  farrayptr=qnwfa2ptr, rc=rc)
1335  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1336  call error_handler("IN FieldGet", rc)
1337 
1338  DO k=1,lev_target
1339  DO i=clb(1),cub(1)
1340  DO j=clb(2),cub(2)
1341  qnifa2ptr(i,j,k) = c2(i,j,k,1)
1342  qnwfa2ptr(i,j,k) = c2(i,j,k,2)
1343  ENDDO
1344  ENDDO
1345  ENDDO
1346 
1347  DEALLOCATE (z1, z2, c1, c2)
1348 
1349  call esmf_fielddestroy(qnifa_climo_b4adj_target_grid, rc=rc)
1350  call esmf_fielddestroy(qnwfa_climo_b4adj_target_grid, rc=rc)
1351  call esmf_fielddestroy(thomp_pres_climo_b4adj_target_grid, rc=rc)
1352 
1353  END SUBROUTINE vintg_thomp_mp_climo
1354 
1355 
1368  SUBROUTINE vintg_wam (YEAR,MONTH,DAY,HOUR)
1369 
1370  IMPLICIT NONE
1371 
1372  include 'mpif.h'
1373 
1374  INTEGER, INTENT(IN) :: YEAR,MONTH,DAY,HOUR
1375 
1376  REAL(ESMF_KIND_R8), PARAMETER :: AMO = 15.9994 ! molecular weight of o
1377  REAL(ESMF_KIND_R8), PARAMETER :: AMO2 = 31.999 !molecular weight of o2
1378  REAL(ESMF_KIND_R8), PARAMETER :: AMN2 = 28.013 !molecular weight of n2
1379 
1380  REAL(ESMF_KIND_R8) :: COE,WFUN(10),DEGLAT,HOLD
1381  REAL(ESMF_KIND_R8) :: SUMMASS,QVMASS,O3MASS
1382  INTEGER :: I, J, K, II, CLB(3), CUB(3), RC, KREF
1383  INTEGER :: IDAT(8),JDOW,JDAY,ICDAY
1384 
1385  REAL(ESMF_KIND_R8), ALLOCATABLE :: TEMP(:),ON(:),O2N(:),N2N(:),PRMB(:)
1386 
1387  REAL(ESMF_KIND_R8), POINTER :: LATPTR(:,:) ! output latitude
1388  REAL(ESMF_KIND_R8), POINTER :: P1PTR(:,:,:) ! input pressure
1389  REAL(ESMF_KIND_R8), POINTER :: P2PTR(:,:,:) ! output pressure
1390  REAL(ESMF_KIND_R8), POINTER :: DZDT2PTR(:,:,:) ! output vvel
1391  REAL(ESMF_KIND_R8), POINTER :: T2PTR(:,:,:) ! output temperature
1392  REAL(ESMF_KIND_R8), POINTER :: Q2PTR(:,:,:) ! output tracer
1393  REAL(ESMF_KIND_R8), POINTER :: QVPTR(:,:,:) ! output tracer
1394  REAL(ESMF_KIND_R8), POINTER :: QOPTR(:,:,:) ! output tracer
1395  REAL(ESMF_KIND_R8), POINTER :: O2PTR(:,:,:) ! output tracer
1396  REAL(ESMF_KIND_R8), POINTER :: O3PTR(:,:,:) ! output tracer
1397  REAL(ESMF_KIND_R8), POINTER :: WIND2PTR(:,:,:,:) ! output wind (x,y,z components)
1398 
1399 ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1400 
1401  print*,"VINTG_WAM:- VERTICALY EXTEND FIELDS FOR WAM COLD START."
1402 
1403 ! prepare date
1404  idat = 0
1405  jdow = 0
1406  jday = 0
1407  icday = 0
1408  idat(1)=year
1409  idat(2)=month
1410  idat(3)=day
1411  idat(5)=hour
1412  CALL w3doxdat(idat,jdow,icday,jday)
1413  print *,"VINTG_WAM: WAM START DATE FOR ICDAY=",icday
1414 
1415 ! prepare weighting function
1416  DO k=1,10
1417  wfun(k) = (k-1.0) / 9.0
1418  ENDDO
1419 
1420  ALLOCATE(temp(lev_target))
1421  ALLOCATE(prmb(lev_target))
1422  ALLOCATE( on(lev_target))
1423  ALLOCATE( o2n(lev_target))
1424  ALLOCATE( n2n(lev_target))
1425 
1426 ! p1 (pascal)
1427  print*,"VINTG_WAM:- CALL FieldGet FOR 3-D PRES."
1428  call esmf_fieldget(pres_b4adj_target_grid, &
1429  computationallbound=clb, &
1430  computationalubound=cub, &
1431  farrayptr=p1ptr, rc=rc)
1432  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1433  call error_handler("IN FieldGet", rc)
1434 !print*,"VINTG_WAM: p1ptr ",(p1ptr(1,1,k),k=1,LEV_INPUT)
1435 
1436 ! p2 (pascal)
1437  print*,"VINTG_WAM:- CALL FieldGet FOR 3-D ADJUSTED PRESS"
1438  call esmf_fieldget(pres_target_grid, &
1439  farrayptr=p2ptr, rc=rc)
1440  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1441  call error_handler("IN FieldGet", rc)
1442 !print*,"VINTG_WAM: p2ptr ",(p2ptr(1,1,k),k=1,LEV_TARGET)
1443 
1444 ! latitude in degree
1445  print*,"VINTG_WAM - CALL FieldGet FOR LATITUDE_S."
1446  call esmf_fieldget(latitude_s_target_grid, &
1447  farrayptr=latptr, rc=rc)
1448  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1449  call error_handler("IN FieldGet", rc)
1450 !print*,"VINTG_WAM: latptr ",(latptr(1,j),j=clb(2),cub(2))
1451 
1452 ! temp
1453  print*,"VINTG_WAM:- CALL FieldGet FOR 3-D ADJUSTED TEMP."
1454  call esmf_fieldget(temp_target_grid, &
1455  farrayptr=t2ptr, rc=rc)
1456  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1457  call error_handler("IN FieldGet", rc)
1458 
1459 ! dzdt
1460  print*,"VINTG_WAM:- CALL FieldGet FOR ADJUSTED VERTICAL VELOCITY."
1461  call esmf_fieldget(dzdt_target_grid, &
1462  farrayptr=dzdt2ptr, rc=rc)
1463  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1464  call error_handler("IN FieldGet", rc)
1465 
1466 ! wind
1467  print*,"VINTG_WAM:- CALL FieldGet FOR 3-D ADJUSTED WIND."
1468  call esmf_fieldget(wind_target_grid, &
1469  farrayptr=wind2ptr, rc=rc)
1470  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1471  call error_handler("IN FieldGet", rc)
1472 
1473 !
1474 ! determine vertical blending point and modified extrapolation values
1475 !
1476  DO i=clb(1),cub(1)
1477  DO j=clb(2),cub(2)
1478 
1479  DO k=1,lev_target
1480  IF(p2ptr(i,j,k).le.p1ptr(i,j,lev_input)) THEN
1481  kref =k-1
1482 !x print*,'VINTG_WAM: KREF P1 P2 ',KREF,P1PTR(I,J,LEV_INPUT),P2PTR(I,J,K)
1483  GO TO 11
1484  ENDIF
1485  ENDDO
1486  11 CONTINUE
1487 !
1488  DO k=kref,lev_target
1489  coe = p2ptr(i,j,k) / p2ptr(i,j,kref)
1490  wind2ptr(i,j,k,1) = coe*wind2ptr(i,j,k,1)
1491  wind2ptr(i,j,k,2) = coe*wind2ptr(i,j,k,2)
1492  wind2ptr(i,j,k,3) = coe*wind2ptr(i,j,k,3)
1493  dzdt2ptr(i,j,k) = coe*dzdt2ptr(i,j,k)
1494  ENDDO
1495 
1496  ENDDO
1497  ENDDO
1498 
1499 !
1500 ! point necessary tracers
1501 !
1502  DO ii = 1, num_tracers
1503 
1504  print*,"VINTG_WAM:- CALL FieldGet FOR 3-D TRACER ", trim(tracers(ii))
1505  call esmf_fieldget(tracers_target_grid(ii), &
1506  farrayptr=q2ptr, rc=rc)
1507  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1508  call error_handler("IN FieldGet", rc)
1509 
1510  DO j=clb(2),cub(2)
1511  DO i=clb(1),cub(1)
1512  DO k=1,lev_target
1513  IF(p2ptr(i,j,k).le.p1ptr(i,j,lev_input)) THEN
1514  kref =k-1
1515  GO TO 22
1516  ENDIF
1517  ENDDO
1518  22 CONTINUE
1519 !
1520  DO k=kref,lev_target
1521  coe = min(1.0, p2ptr(i,j,k) / p2ptr(i,j,kref) )
1522  q2ptr(i,j,k) = coe * q2ptr(i,j,k)
1523  ENDDO
1524  ENDDO
1525  ENDDO
1526 
1527  IF (trim(tracers(ii)) == "sphum") qvptr => q2ptr
1528  IF (trim(tracers(ii)) == "spo" ) qoptr => q2ptr
1529  IF (trim(tracers(ii)) == "spo2" ) o2ptr => q2ptr
1530  IF (trim(tracers(ii)) == "spo3" ) o3ptr => q2ptr
1531 
1532  ENDDO
1533 
1534 !
1535 ! obtained wam gases distribution and temperature profile
1536 !
1537  DO i=clb(1),cub(1)
1538  DO j=clb(2),cub(2)
1539 !
1540  deglat = latptr(i,j)
1541  DO k=1,lev_target
1542  prmb(k) = p2ptr(i,j,k) * 0.01
1543  ENDDO
1544  CALL gettemp(icday,1,deglat,1,prmb,lev_target,temp,on,o2n,n2n)
1545 !
1546  DO k=1,lev_target
1547  summass = on(k)*amo+o2n(k)*amo2+n2n(k)*amn2
1548  qvmass = summass*qvptr(i,j,k)/(1.-qvptr(i,j,k))
1549  summass = summass+qvmass
1550  o3mass = summass*o3ptr(i,j,k)
1551  summass = summass+o3mass
1552  hold = 1.0 / summass
1553  qoptr(i,j,k) = on(k)*amo *hold
1554  o2ptr(i,j,k) = o2n(k)*amo2*hold
1555  o3ptr(i,j,k) = o3mass * hold
1556  qvptr(i,j,k) = qvmass * hold
1557  ENDDO
1558 !
1559  DO k=1,lev_target
1560  IF(p2ptr(i,j,k).le.p1ptr(i,j,lev_input)) THEN
1561  kref =k-1
1562  GO TO 33
1563  ENDIF
1564  ENDDO
1565  33 CONTINUE
1566 !
1567  DO k=kref,lev_target
1568  t2ptr(i,j,k) = temp(k)
1569  ENDDO
1570  DO k=kref-10,kref-1
1571  t2ptr(i,j,k) = wfun(k-kref+11) * temp(k) + &
1572  (1.- wfun(k-kref+11)) * t2ptr(i,j,k)
1573  ENDDO
1574  ENDDO
1575  ENDDO
1576 
1577  DEALLOCATE (temp, prmb, on, o2n, n2n)
1578 
1579  END SUBROUTINE vintg_wam
1580 
1594  SUBROUTINE vintg
1595  use mpi
1596 
1597  IMPLICIT NONE
1598 
1599  REAL(ESMF_KIND_R8), PARAMETER :: DLTDZ=-6.5e-3*287.05/9.80665
1600  REAL(ESMF_KIND_R8), PARAMETER :: DLPVDRT=-2.5e6/461.50
1601  REAL(ESMF_KIND_R8), PARAMETER :: ONE = 1.0_esmf_kind_r8
1602 
1603  INTEGER :: I, J, K, CLB(3), CUB(3), RC
1604  INTEGER :: IM, KM1, KM2, NT, II
1605 
1606  REAL(ESMF_KIND_R8) :: DZ
1607  REAL(ESMF_KIND_R8), ALLOCATABLE :: Z1(:,:,:), Z2(:,:,:)
1608  REAL(ESMF_KIND_R8), ALLOCATABLE :: C1(:,:,:,:),C2(:,:,:,:)
1609 
1610  REAL(ESMF_KIND_R8), POINTER :: P1PTR(:,:,:) ! input pressure
1611  REAL(ESMF_KIND_R8), POINTER :: P2PTR(:,:,:) ! output pressure
1612  REAL(ESMF_KIND_R8), POINTER :: DZDT1PTR(:,:,:) ! input vvel
1613  REAL(ESMF_KIND_R8), POINTER :: DZDT2PTR(:,:,:) ! output vvel
1614  REAL(ESMF_KIND_R8), POINTER :: T1PTR(:,:,:) ! input temperature
1615  REAL(ESMF_KIND_R8), POINTER :: T2PTR(:,:,:) ! output temperature
1616  REAL(ESMF_KIND_R8), POINTER :: Q1PTR(:,:,:) ! input tracer
1617  REAL(ESMF_KIND_R8), POINTER :: Q2PTR(:,:,:) ! output tracer
1618  REAL(ESMF_KIND_R8), POINTER :: WIND1PTR(:,:,:,:) ! input wind (x,y,z components)
1619  REAL(ESMF_KIND_R8), POINTER :: WIND2PTR(:,:,:,:) ! input wind (x,y,z components)
1620 
1621 ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1622 ! COMPUTE LOG PRESSURE INTERPOLATING COORDINATE
1623 ! AND COPY INPUT WIND, TEMPERATURE, HUMIDITY AND OTHER TRACERS
1624 ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1625 
1626  print*,"- VERTICALY INTERPOLATE FIELDS."
1627 
1628  print*,"- CALL FieldGet FOR 3-D PRES."
1629  call esmf_fieldget(pres_b4adj_target_grid, &
1630  computationallbound=clb, &
1631  computationalubound=cub, &
1632  farrayptr=p1ptr, rc=rc)
1633  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1634  call error_handler("IN FieldGet", rc)
1635 
1636 ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1637 ! The '1'/'2' arrays hold fields before/after interpolation.
1638 ! Note the 'z' component of the horizontal wind will be treated as a
1639 ! tracer. So add one extra third dimension to these 3-d arrays.
1640 ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1641 
1642  ALLOCATE(z1(clb(1):cub(1),clb(2):cub(2),lev_input))
1643  ALLOCATE(z2(clb(1):cub(1),clb(2):cub(2),lev_target))
1644  ALLOCATE(c1(clb(1):cub(1),clb(2):cub(2),lev_input,num_tracers_input+5))
1645  ALLOCATE(c2(clb(1):cub(1),clb(2):cub(2),lev_target,num_tracers_input+5))
1646 
1647  z1 = -log(p1ptr)
1648 
1649  print*,"- CALL FieldGet FOR 3-D ADJUSTED PRESS"
1650  call esmf_fieldget(pres_target_grid, &
1651  farrayptr=p2ptr, rc=rc)
1652  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1653  call error_handler("IN FieldGet", rc)
1654 
1655  z2 = -log(p2ptr)
1656 
1657  print*,"- CALL FieldGet FOR 3-D WIND."
1658  call esmf_fieldget(wind_b4adj_target_grid, &
1659  farrayptr=wind1ptr, rc=rc)
1660  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1661  call error_handler("IN FieldGet", rc)
1662 
1663  c1(:,:,:,1) = wind1ptr(:,:,:,1)
1664  c1(:,:,:,2) = wind1ptr(:,:,:,2)
1665  c1(:,:,:,3) = wind1ptr(:,:,:,3)
1666 
1667  print*,"- CALL FieldGet FOR VERTICAL VELOCITY."
1668  call esmf_fieldget(dzdt_b4adj_target_grid, &
1669  farrayptr=dzdt1ptr, rc=rc)
1670  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1671  call error_handler("IN FieldGet", rc)
1672 
1673  c1(:,:,:,4) = dzdt1ptr(:,:,:)
1674  print*,"MIN MAX W TARGETB4 IN VINTG = ", minval(dzdt1ptr(:,:,:)), maxval(dzdt1ptr(:,:,:))
1675 
1676  print*,"- CALL FieldGet FOR 3-D TEMP."
1677  call esmf_fieldget(temp_b4adj_target_grid, &
1678  farrayptr=t1ptr, rc=rc)
1679  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1680  call error_handler("IN FieldGet", rc)
1681 
1682  c1(:,:,:,5) = t1ptr(:,:,:)
1683 
1684  DO i = 1, num_tracers_input
1685 
1686  print*,"- CALL FieldGet FOR 3-D TRACERS ", trim(tracers(i))
1687  call esmf_fieldget(tracers_b4adj_target_grid(i), &
1688  farrayptr=q1ptr, rc=rc)
1689  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1690  call error_handler("IN FieldGet", rc)
1691 
1692  c1(:,:,:,5+i) = q1ptr(:,:,:)
1693 
1694  ENDDO
1695 
1696 ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1697 ! PERFORM LAGRANGIAN ONE-DIMENSIONAL INTERPOLATION
1698 ! THAT IS 4TH-ORDER IN INTERIOR, 2ND-ORDER IN OUTSIDE INTERVALS
1699 ! AND 1ST-ORDER FOR EXTRAPOLATION.
1700 ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1701 
1702  im = (cub(1)-clb(1)+1) * (cub(2)-clb(2)+1)
1703  km1= lev_input
1704  km2= lev_target
1705  nt= num_tracers_input + 1 ! treat 'z' wind as tracer.
1706 
1707  CALL terp3(im,1,1,1,1,4+nt,(im*km1),(im*km2), &
1708  km1,im,im,z1,c1,km2,im,im,z2,c2)
1709 ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1710 ! COPY OUTPUT WIND, TEMPERATURE, HUMIDITY AND OTHER TRACERS
1711 ! EXCEPT BELOW THE INPUT DOMAIN, LET TEMPERATURE INCREASE WITH A FIXED
1712 ! LAPSE RATE AND LET THE RELATIVE HUMIDITY REMAIN CONSTANT.
1713 ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1714 
1715  print*,"- CALL FieldGet FOR 3-D ADJUSTED TEMP."
1716  call esmf_fieldget(temp_target_grid, &
1717  farrayptr=t2ptr, rc=rc)
1718  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1719  call error_handler("IN FieldGet", rc)
1720 
1721  print*,"- CALL FieldGet FOR ADJUSTED VERTICAL VELOCITY."
1722  call esmf_fieldget(dzdt_target_grid, &
1723  farrayptr=dzdt2ptr, rc=rc)
1724  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1725  call error_handler("IN FieldGet", rc)
1726 
1727  print*,"- CALL FieldGet FOR 3-D ADJUSTED WIND."
1728  call esmf_fieldget(wind_target_grid, &
1729  farrayptr=wind2ptr, rc=rc)
1730  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1731  call error_handler("IN FieldGet", rc)
1732 
1733  DO k=1,lev_target
1734  DO i=clb(1),cub(1)
1735  DO j=clb(2),cub(2)
1736  wind2ptr(i,j,k,1)=c2(i,j,k,1)
1737  wind2ptr(i,j,k,2)=c2(i,j,k,2)
1738  wind2ptr(i,j,k,3)=c2(i,j,k,3)
1739  dzdt2ptr(i,j,k)=c2(i,j,k,4)
1740  dz=z2(i,j,k)-z1(i,j,1)
1741  IF(dz.GE.0) THEN
1742  t2ptr(i,j,k)=c2(i,j,k,5)
1743  ELSE
1744  t2ptr(i,j,k)=c1(i,j,1,5)*exp(dltdz*dz)
1745  ENDIF
1746  ENDDO
1747  ENDDO
1748  ENDDO
1749 
1750  DO ii = 1, num_tracers_input
1751 
1752  print*,"- CALL FieldGet FOR 3-D TRACER ", trim(tracers(ii))
1753  call esmf_fieldget(tracers_target_grid(ii), &
1754  farrayptr=q2ptr, rc=rc)
1755  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
1756  call error_handler("IN FieldGet", rc)
1757 
1758  IF (trim(tracers(ii)) == "sphum") THEN ! specific humidity
1759 
1760  DO k=1,lev_target
1761  DO i=clb(1),cub(1)
1762  DO j=clb(2),cub(2)
1763  dz=z2(i,j,k)-z1(i,j,1)
1764  IF(dz.GE.0) THEN
1765  q2ptr(i,j,k) = c2(i,j,k,5+ii)
1766  ELSE
1767  q2ptr(i,j,k) = c1(i,j,1,5+ii)*exp(dlpvdrt*(one/t2ptr(i,j,k)-one/t1ptr(i,j,1))-dz)
1768  ENDIF
1769  ENDDO
1770  ENDDO
1771  ENDDO
1772 
1773  ELSE ! all other tracers
1774 
1775  DO k=1,lev_target
1776  DO i=clb(1),cub(1)
1777  DO j=clb(2),cub(2)
1778  q2ptr(i,j,k) = c2(i,j,k,5+ii)
1779  ENDDO
1780  ENDDO
1781  ENDDO
1782 
1783  ENDIF
1784 
1785  ENDDO
1786 
1787  DEALLOCATE (z1, z2, c1, c2)
1788 
1789  END SUBROUTINE vintg
1790 
1827  SUBROUTINE terp3(IM,IXZ1,IXQ1,IXZ2,IXQ2,NM,NXQ1,NXQ2, &
1828  KM1,KXZ1,KXQ1,Z1,Q1,KM2,KXZ2,KXQ2,Z2,Q2)
1829  IMPLICIT NONE
1830  INTEGER IM,IXZ1,IXQ1,IXZ2,IXQ2,NM,NXQ1,NXQ2
1831  INTEGER KM1,KXZ1,KXQ1,KM2,KXZ2,KXQ2
1832  INTEGER I,K1,K2,N
1833  INTEGER K1S(IM,KM2)
1834  REAL(ESMF_KIND_R8), PARAMETER :: ONE = 1.0_esmf_kind_r8
1835  REAL(ESMF_KIND_R8) :: Z1(1+(IM-1)*IXZ1+(KM1-1)*KXZ1)
1836  REAL(ESMF_KIND_R8) :: Q1(1+(IM-1)*IXQ1+(KM1-1)*KXQ1+(NM-1)*NXQ1)
1837  REAL(ESMF_KIND_R8) :: Z2(1+(IM-1)*IXZ2+(KM2-1)*KXZ2)
1838  REAL(ESMF_KIND_R8) :: Q2(1+(IM-1)*IXQ2+(KM2-1)*KXQ2+(NM-1)*NXQ2)
1839 ! REAL(ESMF_KIND_R8) :: J2(1+(IM-1)*IXQ2+(KM2-1)*KXQ2+(NM-1)*NXQ2)
1840  REAL(ESMF_KIND_R8) :: FFA(IM),FFB(IM),FFC(IM),FFD(IM)
1841  REAL(ESMF_KIND_R8) :: GGA(IM),GGB(IM),GGC(IM),GGD(IM)
1842  REAL(ESMF_KIND_R8) :: Z1A,Z1B,Z1C,Z1D,Q1A,Q1B,Q1C,Q1D,Z2S,Q2S
1843 ! REAL(ESMF_KIND_R8) :: J2S
1844 
1845 ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1846 ! FIND THE SURROUNDING INPUT INTERVAL FOR EACH OUTPUT POINT.
1847  CALL rsearch(im,km1,ixz1,kxz1,z1,km2,ixz2,kxz2,z2,1,im,k1s)
1848 
1849 ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1850 ! GENERALLY INTERPOLATE CUBICALLY WITH MONOTONIC CONSTRAINT
1851 ! FROM TWO NEAREST INPUT POINTS ON EITHER SIDE OF THE OUTPUT POINT,
1852 ! BUT WITHIN THE TWO EDGE INTERVALS INTERPOLATE LINEARLY.
1853 ! KEEP THE OUTPUT FIELDS CONSTANT OUTSIDE THE INPUT DOMAIN.
1854 
1855 !$OMP PARALLEL DO DEFAULT(PRIVATE) SHARED(IM,IXZ1,IXQ1,IXZ2), &
1856 !$OMP& SHARED(IXQ2,NM,NXQ1,NXQ2,KM1,KXZ1,KXQ1,Z1,Q1,KM2,KXZ2), &
1857 !$OMP& SHARED(KXQ2,Z2,Q2,K1S)
1858  DO k2=1,km2
1859  DO i=1,im
1860  k1=k1s(i,k2)
1861  IF(k1.EQ.1.OR.k1.EQ.km1-1) THEN
1862  z2s=z2(1+(i-1)*ixz2+(k2-1)*kxz2)
1863  z1a=z1(1+(i-1)*ixz1+(k1-1)*kxz1)
1864  z1b=z1(1+(i-1)*ixz1+(k1+0)*kxz1)
1865  ffa(i)=(z2s-z1b)/(z1a-z1b)
1866  ffb(i)=(z2s-z1a)/(z1b-z1a)
1867  gga(i)=one/(z1a-z1b)
1868  ggb(i)=one/(z1b-z1a)
1869  ELSEIF(k1.GT.1.AND.k1.LT.km1-1) THEN
1870  z2s=z2(1+(i-1)*ixz2+(k2-1)*kxz2)
1871  z1a=z1(1+(i-1)*ixz1+(k1-2)*kxz1)
1872  z1b=z1(1+(i-1)*ixz1+(k1-1)*kxz1)
1873  z1c=z1(1+(i-1)*ixz1+(k1+0)*kxz1)
1874  z1d=z1(1+(i-1)*ixz1+(k1+1)*kxz1)
1875  ffa(i)=(z2s-z1b)/(z1a-z1b)* &
1876  (z2s-z1c)/(z1a-z1c)* &
1877  (z2s-z1d)/(z1a-z1d)
1878  ffb(i)=(z2s-z1a)/(z1b-z1a)* &
1879  (z2s-z1c)/(z1b-z1c)* &
1880  (z2s-z1d)/(z1b-z1d)
1881  ffc(i)=(z2s-z1a)/(z1c-z1a)* &
1882  (z2s-z1b)/(z1c-z1b)* &
1883  (z2s-z1d)/(z1c-z1d)
1884  ffd(i)=(z2s-z1a)/(z1d-z1a)* &
1885  (z2s-z1b)/(z1d-z1b)* &
1886  (z2s-z1c)/(z1d-z1c)
1887  gga(i)= one/(z1a-z1b)* &
1888  (z2s-z1c)/(z1a-z1c)* &
1889  (z2s-z1d)/(z1a-z1d)+ &
1890  (z2s-z1b)/(z1a-z1b)* &
1891  one/(z1a-z1c)* &
1892  (z2s-z1d)/(z1a-z1d)+ &
1893  (z2s-z1b)/(z1a-z1b)* &
1894  (z2s-z1c)/(z1a-z1c)* &
1895  one/(z1a-z1d)
1896  ggb(i)= one/(z1b-z1a)* &
1897  (z2s-z1c)/(z1b-z1c)* &
1898  (z2s-z1d)/(z1b-z1d)+ &
1899  (z2s-z1a)/(z1b-z1a)* &
1900  one/(z1b-z1c)* &
1901  (z2s-z1d)/(z1b-z1d)+ &
1902  (z2s-z1a)/(z1b-z1a)* &
1903  (z2s-z1c)/(z1b-z1c)* &
1904  one/(z1b-z1d)
1905  ggc(i)= one/(z1c-z1a)* &
1906  (z2s-z1b)/(z1c-z1b)* &
1907  (z2s-z1d)/(z1c-z1d)+ &
1908  (z2s-z1a)/(z1c-z1a)* &
1909  one/(z1c-z1b)* &
1910  (z2s-z1d)/(z1c-z1d)+ &
1911  (z2s-z1a)/(z1c-z1a)* &
1912  (z2s-z1b)/(z1c-z1b)* &
1913  one/(z1c-z1d)
1914  ggd(i)= one/(z1d-z1a)* &
1915  (z2s-z1b)/(z1d-z1b)* &
1916  (z2s-z1c)/(z1d-z1c)+ &
1917  (z2s-z1a)/(z1d-z1a)* &
1918  one/(z1d-z1b)* &
1919  (z2s-z1c)/(z1d-z1c)+ &
1920  (z2s-z1a)/(z1d-z1a)* &
1921  (z2s-z1b)/(z1d-z1b)* &
1922  one/(z1d-z1c)
1923  ENDIF
1924  ENDDO
1925 
1926 ! INTERPOLATE.
1927  DO n=1,nm
1928  DO i=1,im
1929  k1=k1s(i,k2)
1930  IF(k1.EQ.0) THEN
1931  q2s=q1(1+(i-1)*ixq1+(n-1)*nxq1)
1932 ! J2S=0
1933  ELSEIF(k1.EQ.km1) THEN
1934  q2s=q1(1+(i-1)*ixq1+(km1-1)*kxq1+(n-1)*nxq1)
1935 ! J2S=0
1936  ELSEIF(k1.EQ.1.OR.k1.EQ.km1-1) THEN
1937  q1a=q1(1+(i-1)*ixq1+(k1-1)*kxq1+(n-1)*nxq1)
1938  q1b=q1(1+(i-1)*ixq1+(k1+0)*kxq1+(n-1)*nxq1)
1939  q2s=ffa(i)*q1a+ffb(i)*q1b
1940 ! J2S=GGA(I)*Q1A+GGB(I)*Q1B
1941  ELSE
1942  q1a=q1(1+(i-1)*ixq1+(k1-2)*kxq1+(n-1)*nxq1)
1943  q1b=q1(1+(i-1)*ixq1+(k1-1)*kxq1+(n-1)*nxq1)
1944  q1c=q1(1+(i-1)*ixq1+(k1+0)*kxq1+(n-1)*nxq1)
1945  q1d=q1(1+(i-1)*ixq1+(k1+1)*kxq1+(n-1)*nxq1)
1946  q2s=ffa(i)*q1a+ffb(i)*q1b+ffc(i)*q1c+ffd(i)*q1d
1947 ! J2S=GGA(I)*Q1A+GGB(I)*Q1B+GGC(I)*Q1C+GGD(I)*Q1D
1948  IF(q2s.LT.min(q1b,q1c)) THEN
1949  q2s=min(q1b,q1c)
1950 ! J2S=0
1951  ELSEIF(q2s.GT.max(q1b,q1c)) THEN
1952  q2s=max(q1b,q1c)
1953 ! J2S=0
1954  ENDIF
1955  ENDIF
1956  q2(1+(i-1)*ixq2+(k2-1)*kxq2+(n-1)*nxq2)=q2s
1957 ! J2(1+(I-1)*IXQ2+(K2-1)*KXQ2+(N-1)*NXQ2)=J2S
1958  ENDDO
1959  ENDDO
1960  ENDDO
1961 !$OMP END PARALLEL DO
1962 
1963  END SUBROUTINE terp3
1964 
2021  SUBROUTINE rsearch(IM,KM1,IXZ1,KXZ1,Z1,KM2,IXZ2,KXZ2,Z2,IXL2,KXL2,L2)
2022  IMPLICIT NONE
2023 
2024  INTEGER,INTENT(IN) :: IM,KM1,IXZ1,KXZ1,KM2,IXZ2,KXZ2,IXL2,KXL2
2025  INTEGER,INTENT(OUT) :: L2(1+(IM-1)*IXL2+(KM2-1)*KXL2)
2026 
2027  REAL(ESMF_KIND_R8),INTENT(IN) :: Z1(1+(IM-1)*IXZ1+(KM1-1)*KXZ1)
2028  REAL(ESMF_KIND_R8),INTENT(IN) :: Z2(1+(IM-1)*IXZ2+(KM2-1)*KXZ2)
2029 
2030  INTEGER :: I,K2,L
2031 
2032  REAL(ESMF_KIND_R8) :: Z
2033 
2034 
2035 ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
2036 ! FIND THE SURROUNDING INPUT INTERVAL FOR EACH OUTPUT POINT.
2037  DO i=1,im
2038  IF (z1(1+(i-1)*ixz1).LE.z1(1+(i-1)*ixz1+(km1-1)*kxz1)) THEN
2039 ! INPUT COORDINATE IS MONOTONICALLY ASCENDING.
2040  DO k2=1,km2
2041  z=z2(1+(i-1)*ixz2+(k2-1)*kxz2)
2042  l=0
2043  DO
2044  IF(z.LT.z1(1+(i-1)*ixz1+l*kxz1)) EXIT
2045  l=l+1
2046  IF(l.EQ.km1) EXIT
2047  ENDDO
2048  l2(1+(i-1)*ixl2+(k2-1)*kxl2)=l
2049  ENDDO
2050  ELSE
2051 ! INPUT COORDINATE IS MONOTONICALLY DESCENDING.
2052  DO k2=1,km2
2053  z=z2(1+(i-1)*ixz2+(k2-1)*kxz2)
2054  l=0
2055  DO
2056  IF(z.GT.z1(1+(i-1)*ixz1+l*kxz1)) EXIT
2057  l=l+1
2058  IF(l.EQ.km1) EXIT
2059  ENDDO
2060  l2(1+(i-1)*ixl2+(k2-1)*kxl2)=l
2061  ENDDO
2062  ENDIF
2063  ENDDO
2064 
2065  END SUBROUTINE rsearch
2066 
2069  subroutine compute_zh
2070 
2071  implicit none
2072 
2073  integer :: i,ii, j,k, rc, clb(2), cub(2)
2074 
2075  real(esmf_kind_r8), allocatable :: pe0(:), pn0(:)
2076  real(esmf_kind_r8), pointer :: psptr(:,:)
2077  real(esmf_kind_r8), pointer :: zhsfcptr(:,:)
2078  real(esmf_kind_r8), pointer :: zhptr(:,:,:)
2079  real(esmf_kind_r8), pointer :: tptr(:,:,:)
2080  real(esmf_kind_r8), pointer :: qptr(:,:,:)
2081  real(esmf_kind_r8) :: ak, bk, zvir, grd
2082  real(esmf_kind_r8), parameter :: grav = 9.80665
2083  real(esmf_kind_r8), parameter :: rdgas = 287.05
2084  real(esmf_kind_r8), parameter :: rvgas = 461.50
2085 
2086  print*,"- COMPUTE HEIGHT"
2087 
2088  print*,"- CALL FieldGet FOR SURFACE PRESSURE"
2089  call esmf_fieldget(ps_target_grid, &
2090  computationallbound=clb, &
2091  computationalubound=cub, &
2092  farrayptr=psptr, rc=rc)
2093  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
2094  call error_handler("IN FieldGet", rc)
2095 
2096  print*,"- CALL FieldGet FOR TERRAIN HEIGHT"
2097  call esmf_fieldget(terrain_target_grid, &
2098  farrayptr=zhsfcptr, rc=rc)
2099  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
2100  call error_handler("IN FieldGet", rc)
2101 
2102  print*,"- CALL FieldGet FOR HEIGHT"
2103  call esmf_fieldget(zh_target_grid, &
2104  farrayptr=zhptr, rc=rc)
2105  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
2106  call error_handler("IN FieldGet", rc)
2107 
2108  print*,"- CALL FieldGet FOR TEMPERATURE"
2109  call esmf_fieldget(temp_target_grid, &
2110  farrayptr=tptr, rc=rc)
2111  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
2112  call error_handler("IN FieldGet", rc)
2113 
2114  do ii = 1, num_tracers
2115  if (trim(tracers(ii)) == "sphum") exit
2116  enddo
2117 
2118  print*,"- CALL FieldGet FOR SPECIFIC HUMIDITY"
2119  call esmf_fieldget(tracers_target_grid(ii), &
2120  farrayptr=qptr, rc=rc)
2121  if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
2122  call error_handler("IN FieldGet", rc)
2123 
2124  grd = grav/rdgas
2125  zvir = rvgas/rdgas - 1.0_esmf_kind_r8
2126 
2127  allocate(pe0(levp1_target))
2128  allocate(pn0(levp1_target))
2129 
2130  do j = clb(2), cub(2)
2131  do i = clb(1), cub(1)
2132 
2133  do k = 1, levp1_target
2134  ak = vcoord_target(k,1)
2135  ak = max(ak, 1.e-9)
2136  bk = vcoord_target(k,2)
2137 
2138  pe0(k) = ak + bk*psptr(i,j)
2139  pn0(k) = log(pe0(k))
2140  enddo
2141 
2142  zhptr(i,j,1) = zhsfcptr(i,j)
2143 
2144  do k = 2, levp1_target
2145  zhptr(i,j,k) = zhptr(i,j,k-1)+tptr(i,j,k-1)*(1.+zvir*qptr(i,j,k-1))* &
2146  (pn0(k-1)-pn0(k))/grd
2147  enddo
2148 
2149  enddo
2150  enddo
2151 
2152  deallocate(pe0, pn0)
2153 
2154  end subroutine compute_zh
2155 
2159  subroutine cleanup_target_atm_b4adj_data
2160 
2161  implicit none
2162 
2163  integer :: i, rc
2164 
2165  print*,"- DESTROY TARGET GRID ATMOSPHERIC BEFORE ADJUSTMENT FIELDS."
2166 
2167  call esmf_fielddestroy(wind_b4adj_target_grid, rc=rc)
2168  call esmf_fielddestroy(dzdt_b4adj_target_grid, rc=rc)
2169  call esmf_fielddestroy(ps_b4adj_target_grid, rc=rc)
2170  call esmf_fielddestroy(pres_b4adj_target_grid, rc=rc)
2171  call esmf_fielddestroy(temp_b4adj_target_grid, rc=rc)
2172  call esmf_fielddestroy(terrain_interp_to_target_grid, rc=rc)
2173 
2174  do i = 1, num_tracers_input
2175  call esmf_fielddestroy(tracers_b4adj_target_grid(i), rc=rc)
2176  enddo
2177 
2178  deallocate(tracers_b4adj_target_grid)
2179 
2180  end subroutine cleanup_target_atm_b4adj_data
2181 
2184  subroutine cleanup_all_target_atm_data
2185 
2186  use atmosphere_target_data, only : cleanup_atmosphere_target_data
2187 
2188  implicit none
2189 
2190  integer :: rc
2191 
2192  print*,"- DESTROY LOCAL TARGET GRID ATMOSPHERIC FIELDS."
2193 
2194  call esmf_fielddestroy(wind_target_grid, rc=rc)
2195  call esmf_fielddestroy(wind_s_target_grid, rc=rc)
2196  call esmf_fielddestroy(wind_w_target_grid, rc=rc)
2197  call esmf_fielddestroy(pres_target_grid, rc=rc)
2198 
2199  call cleanup_atmosphere_target_data
2200 
2201  end subroutine cleanup_all_target_atm_data
2202 
2203  end module atmosphere
atmosphere::compute_zh
subroutine compute_zh
Compute vertical level height.
Definition: atmosphere.F90:2070
atmosphere::newps
subroutine newps(localpet)
Computes adjusted surface pressure given a new terrain height.
Definition: atmosphere.F90:902
program_setup
This module contains code to read the setup namelist file, handle the varmap file for GRIB2 data...
Definition: program_setup.F90:9
program_setup::num_tracers_input
integer, public num_tracers_input
Number of atmospheric tracers in input file.
Definition: program_setup.F90:61
atmosphere_target_data::v_w_target_grid
type(esmf_field), public v_w_target_grid
V-wind, &#39;west&#39; edge of grid cell.
Definition: atmosphere_target_data.F90:31
thompson_mp_climo_data::lev_thomp_mp_climo
integer, public lev_thomp_mp_climo
number of vert lvls of Thompson climo data
Definition: thompson_mp_climo_data.F90:25
model_grid::longitude_w_target_grid
type(esmf_field), public longitude_w_target_grid
longitude of &#39;west&#39; edge of grid box, target grid
Definition: model_grid.F90:90
atmosphere::cleanup_target_atm_b4adj_data
subroutine cleanup_target_atm_b4adj_data
Cleanup atmospheric field (before adjustment) objects.
Definition: atmosphere.F90:2160
atmosphere::create_atm_esmf_fields
subroutine create_atm_esmf_fields
Create target grid field objects.
Definition: atmosphere.F90:508
atm_input_data
Read atmospheric data on the input grid.
Definition: atm_input_data.F90:15
atmosphere::wind_target_grid
type(esmf_field) wind_target_grid
3-d wind, grid box center
Definition: atmosphere.F90:84
program_setup::cycle_mon
integer, public cycle_mon
Cycle month.
Definition: program_setup.F90:85
atmosphere::temp_b4adj_target_grid
type(esmf_field) temp_b4adj_target_grid
temp before vert adj
Definition: atmosphere.F90:82
program_setup::regional
integer, public regional
For regional target grids.
Definition: program_setup.F90:88
thompson_mp_climo_data::qnifa_climo_input_grid
type(esmf_field), public qnifa_climo_input_grid
number concentration of ice friendly nuclei.
Definition: thompson_mp_climo_data.F90:31
atm_input_data::dzdt_input_grid
type(esmf_field), public dzdt_input_grid
vert velocity
Definition: atm_input_data.F90:48
program_setup::tracers
character(len=20), dimension(max_tracers), public tracers
Name of each atmos tracer to be processed.
Definition: program_setup.F90:66
atm_input_data::read_input_atm_data
subroutine, public read_input_atm_data(localpet)
Read input grid atmospheric data driver.
Definition: atm_input_data.F90:76
atmosphere_target_data::vcoord_target
real(esmf_kind_r8), dimension(:,:), allocatable, public vcoord_target
Vertical coordinate.
Definition: atmosphere_target_data.F90:21
atmosphere::atmosphere_driver
subroutine, public atmosphere_driver(localpet)
Driver routine to process for atmospheric fields.
Definition: atmosphere.F90:108
thompson_mp_climo_data::qnwfa_climo_input_grid
type(esmf_field), public qnwfa_climo_input_grid
number concentration of water friendly nuclei.
Definition: thompson_mp_climo_data.F90:34
program_setup::cycle_year
integer, public cycle_year
Cycle year.
Definition: program_setup.F90:84
thompson_mp_climo_data::read_thomp_mp_climo_data
subroutine, public read_thomp_mp_climo_data
Read Thompson climatological MP data file and time interpolate data to current cycle time...
Definition: thompson_mp_climo_data.F90:51
atmosphere::read_vcoord_info
subroutine, public read_vcoord_info
Reads model vertical coordinate definition file (as specified by namelist variable vcoord_file_target...
Definition: atmosphere.F90:1099
atm_input_data::tracers_input_grid
type(esmf_field), dimension(:), allocatable, public tracers_input_grid
tracers
Definition: atm_input_data.F90:58
model_grid::longitude_s_target_grid
type(esmf_field), public longitude_s_target_grid
longitude of &#39;south&#39; edge of grid box, target grid
Definition: model_grid.F90:87
model_grid
Sets up the ESMF grid objects for the input data grid and target FV3 grid.
Definition: model_grid.F90:9
atmosphere_target_data::qnifa_climo_target_grid
type(esmf_field), public qnifa_climo_target_grid
Number concentration of ice friendly aerosols.
Definition: atmosphere_target_data.F90:33
atm_input_data::pres_input_grid
type(esmf_field), public pres_input_grid
3-d pressure
Definition: atm_input_data.F90:50
atm_input_data::levp1_input
integer, public levp1_input
number of atmos layer interfaces
Definition: atm_input_data.F90:61
program_setup::cycle_day
integer, public cycle_day
Cycle day.
Definition: program_setup.F90:86
model_grid::latitude_s_target_grid
type(esmf_field), public latitude_s_target_grid
latitude of &#39;south&#39; edge of grid box, target grid
Definition: model_grid.F90:79
atmosphere_target_data::zh_target_grid
type(esmf_field), public zh_target_grid
3-d height.
Definition: atmosphere_target_data.F90:32
atmosphere::pres_b4adj_target_grid
type(esmf_field) pres_b4adj_target_grid
3-d pres before terrain adj
Definition: atmosphere.F90:81
model_grid::target_grid
type(esmf_grid), public target_grid
target grid esmf grid object.
Definition: model_grid.F90:54
model_grid::latitude_w_target_grid
type(esmf_field), public latitude_w_target_grid
latitude of &#39;west&#39; edge of grid box, target grid
Definition: model_grid.F90:82
atmosphere_target_data::delp_target_grid
type(esmf_field), public delp_target_grid
Pressure thickness.
Definition: atmosphere_target_data.F90:23
atmosphere_target_data::lev_target
integer, public lev_target
Number of vertical levels.
Definition: atmosphere_target_data.F90:17
thompson_mp_climo_data
Module to read the Thompson climatological MP data file and set up the associated esmf field and grid...
Definition: thompson_mp_climo_data.F90:9
atmosphere::cleanup_all_target_atm_data
subroutine cleanup_all_target_atm_data
Cleanup target grid atmospheric field objects.
Definition: atmosphere.F90:2185
program_setup::wam_cold_start
logical, public wam_cold_start
When true, cold start for whole atmosphere model.
Definition: program_setup.F90:97
atmosphere::horiz_interp_thomp_mp_climo
subroutine horiz_interp_thomp_mp_climo
Horizontally interpolate thompson microphysics data to the target model grid.
Definition: atmosphere.F90:1134
program_setup::vcoord_file_target_grid
character(len=500), public vcoord_file_target_grid
Vertical coordinate definition file.
Definition: program_setup.F90:39
atm_input_data::wind_input_grid
type(esmf_field), public wind_input_grid
3-component wind
Definition: atm_input_data.F90:57
atmosphere_target_data
Module to hold variables and ESMF fields associated with the target grid atmospheric data...
Definition: atmosphere_target_data.F90:9
atmosphere::newpr1
subroutine newpr1(localpet)
Computes 3-D pressure given an adjusted surface pressure.
Definition: atmosphere.F90:791
atmosphere::convert_winds
subroutine convert_winds
Convert 3-d component winds to u and v.
Definition: atmosphere.F90:647
atmosphere::pres_target_grid
type(esmf_field) pres_target_grid
3-d pressure
Definition: atmosphere.F90:80
atmosphere::terrain_interp_to_target_grid
type(esmf_field) terrain_interp_to_target_grid
Input grid terrain interpolated to target grid.
Definition: atmosphere.F90:83
atmosphere::wind_b4adj_target_grid
type(esmf_field) wind_b4adj_target_grid
3-d wind before vert adj
Definition: atmosphere.F90:85
atm_input_data::terrain_input_grid
type(esmf_field), public terrain_input_grid
terrain height
Definition: atm_input_data.F90:52
atm_input_data::ps_input_grid
type(esmf_field), public ps_input_grid
surface pressure
Definition: atm_input_data.F90:51
program_setup::atm_weight_file
character(len=500), public atm_weight_file
File containing pre-computed weights to horizontally interpolate atmospheric fields.
Definition: program_setup.F90:42
atmosphere_target_data::qnwfa_climo_target_grid
type(esmf_field), public qnwfa_climo_target_grid
Number concentration of water friendly aerosols.
Definition: atmosphere_target_data.F90:35
atmosphere_target_data::v_s_target_grid
type(esmf_field), public v_s_target_grid
V-wind, &#39;south&#39; edge of grid cell.
Definition: atmosphere_target_data.F90:29
atmosphere::thomp_pres_climo_b4adj_target_grid
type(esmf_field) thomp_pres_climo_b4adj_target_grid
pressure of each level on target grid
Definition: atmosphere.F90:95
atmosphere_target_data::nvcoord_target
integer, public nvcoord_target
Number of vertical coordinate variables.
Definition: atmosphere_target_data.F90:19
atmosphere_target_data::u_s_target_grid
type(esmf_field), public u_s_target_grid
U-wind, &#39;south&#39; edge of grid cell.
Definition: atmosphere_target_data.F90:28
atmosphere_target_data::ps_target_grid
type(esmf_field), public ps_target_grid
Surface pressure.
Definition: atmosphere_target_data.F90:25
atmosphere::tracers_b4adj_target_grid
type(esmf_field), dimension(:), allocatable tracers_b4adj_target_grid
tracers before vert adj
Definition: atmosphere.F90:78
atmosphere_target_data::u_w_target_grid
type(esmf_field), public u_w_target_grid
U-wind, &#39;west&#39; edge of grid cell.
Definition: atmosphere_target_data.F90:30
gettemp
subroutine gettemp(iday, nday, xlat, nlat, pr, np, temp, n_o, n_o2, n_n2)
Entry routine to get WAM needed temperature and composition profiles.
Definition: wam_climo_data.f90:968
atm_input_data::cleanup_input_atm_data
subroutine, public cleanup_input_atm_data
Free up memory associated with atm data.
Definition: atm_input_data.F90:3243
thompson_mp_climo_data::thomp_pres_climo_input_grid
type(esmf_field), public thomp_pres_climo_input_grid
3-d pressure of the Thompson climo data points
Definition: thompson_mp_climo_data.F90:37
atmosphere
Process atmospheric fields.
Definition: atmosphere.F90:19
atmosphere::wind_w_target_grid
type(esmf_field) wind_w_target_grid
3-d wind, &#39;west&#39; edge
Definition: atmosphere.F90:87
thompson_mp_climo_data::cleanup_thomp_mp_climo_input_data
subroutine, public cleanup_thomp_mp_climo_input_data
Free up memory associated with this module.
Definition: thompson_mp_climo_data.F90:329
program_setup::num_tracers
integer, public num_tracers
Number of atmospheric tracers to be processed.
Definition: program_setup.F90:60
atm_input_data::temp_input_grid
type(esmf_field), public temp_input_grid
temperature
Definition: atm_input_data.F90:53
model_grid::terrain_target_grid
type(esmf_field), public terrain_target_grid
terrain height target grid
Definition: model_grid.F90:96
atmosphere::ps_b4adj_target_grid
type(esmf_field) ps_b4adj_target_grid
sfc pres before terrain adj
Definition: atmosphere.F90:79
atmosphere_target_data::temp_target_grid
type(esmf_field), public temp_target_grid
Temperautre.
Definition: atmosphere_target_data.F90:26
atmosphere_target_data::levp1_target
integer, public levp1_target
Number of vertical levels plus 1.
Definition: atmosphere_target_data.F90:18
program_setup::use_thomp_mp_climo
logical, public use_thomp_mp_climo
When true, read and process Thompson MP climatological tracers.
Definition: program_setup.F90:125
atmosphere_target_data::tracers_target_grid
type(esmf_field), dimension(:), allocatable, public tracers_target_grid
Tracers.
Definition: atmosphere_target_data.F90:27
program_setup::cycle_hour
integer, public cycle_hour
Cycle hour.
Definition: program_setup.F90:87
atm_input_data::lev_input
integer, public lev_input
number of atmospheric layers
Definition: atm_input_data.F90:60
atmosphere::wind_s_target_grid
type(esmf_field) wind_s_target_grid
3-d wind, &#39;south&#39; edge
Definition: atmosphere.F90:86
atmosphere_target_data::dzdt_target_grid
type(esmf_field), public dzdt_target_grid
Vertical velocity.
Definition: atmosphere_target_data.F90:24
atmosphere::qnifa_climo_b4adj_target_grid
type(esmf_field) qnifa_climo_b4adj_target_grid
number concentration of ice friendly aerosols before vert adj
Definition: atmosphere.F90:91
atmosphere_target_data::cleanup_atmosphere_target_data
subroutine, public cleanup_atmosphere_target_data
Free up memory for fields and variables in this module.
Definition: atmosphere_target_data.F90:46
atmosphere::create_atm_b4adj_esmf_fields
subroutine create_atm_b4adj_esmf_fields
Create target grid field objects to hold data before vertical interpolation.
Definition: atmosphere.F90:434
atmosphere::dzdt_b4adj_target_grid
type(esmf_field) dzdt_b4adj_target_grid
vertical vel before vert adj
Definition: atmosphere.F90:77
atmosphere::qnwfa_climo_b4adj_target_grid
type(esmf_field) qnwfa_climo_b4adj_target_grid
number concentration of water friendly aerosols before vert adj
Definition: atmosphere.F90:93