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