chgres_cube 1.14.0
Loading...
Searching...
No Matches
thompson_mp_climo_data.F90
Go to the documentation of this file.
1
4
9 module thompson_mp_climo_data
10
11 use esmf
12 use netcdf
13 use program_setup, only : cycle_mon, cycle_day, cycle_hour, &
14 thomp_mp_climo_file
15 use utilities, only : error_handler, netcdf_err
16
17 implicit none
18
19 private
20
21 integer :: i_thomp_mp_climo
22
23 integer :: j_thomp_mp_climo
24
25 integer, public :: lev_thomp_mp_climo
26
27
28 type(esmf_grid) :: thomp_mp_climo_grid
29
30
31 type(esmf_field), public :: qnifa_climo_input_grid
32
34 type(esmf_field), public :: qnwfa_climo_input_grid
35
37 type(esmf_field), public :: thomp_pres_climo_input_grid
38
40
41 public :: read_thomp_mp_climo_data
42 public :: cleanup_thomp_mp_climo_input_data
43
44 contains
45
50 subroutine read_thomp_mp_climo_data
51
52 implicit none
53
54 integer :: error, ncid, rc, clb(2), cub(2)
55 integer :: i, j, localpet, npets, id_var
56 integer :: jda(8), jdow, jdoy, jday, id_dim
57 integer :: mm, mmm, mmp, mon1, mon2, idum(3)
58
59 real(esmf_kind_r8), allocatable :: dummy3d(:,:,:)
60 real(esmf_kind_r8), allocatable :: dummy3d_mon1(:,:,:)
61 real(esmf_kind_r8), allocatable :: dummy3d_mon2(:,:,:)
62 real(esmf_kind_r8), pointer :: lat_ptr(:,:), lon_ptr(:,:)
63 real(esmf_kind_r8), allocatable :: lons(:), lats(:)
64 real :: rjday, dayhf(13), wei1m, wei2m
65
66 type(esmf_vm) :: vm
67
68 type(esmf_polekind_flag) :: polekindflag(2)
69
70 data dayhf/ 15.5, 45.0, 74.5,105.0,135.5,166.0, &
71 196.5,227.5,258.0,288.5,319.0,349.5,380.5/
72
73!-----------------------------------------------------------------------------------
74! Open the file and read the grid dimensions and latitude/longitude.
75!-----------------------------------------------------------------------------------
76
77 print*,"- CALL VMGetGlobal"
78 call esmf_vmgetglobal(vm, rc=rc)
79 if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
80 call error_handler("IN VMGetGlobal", rc)
81
82 print*,"- CALL VMGet"
83 call esmf_vmget(vm, localpet=localpet, petcount=npets, rc=rc)
84 if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
85 call error_handler("IN VMGet", rc)
86
87 if (localpet == 0) then
88 print*,"- READ THOMP_MP_CLIMO_FILE: ", trim(thomp_mp_climo_file)
89 error=nf90_open(trim(thomp_mp_climo_file),nf90_nowrite,ncid)
90 call netcdf_err(error, 'opening: '//trim(thomp_mp_climo_file) )
91
92 error=nf90_inq_dimid(ncid, 'lat', id_dim)
93 call netcdf_err(error, 'reading lat id')
94 error=nf90_inquire_dimension(ncid,id_dim,len=idum(1))
95 call netcdf_err(error, 'reading lat')
96
97 error=nf90_inq_dimid(ncid, 'lon', id_dim)
98 call netcdf_err(error, 'reading lon id')
99 error=nf90_inquire_dimension(ncid,id_dim,len=idum(2))
100 call netcdf_err(error, 'reading lon')
101
102 error=nf90_inq_dimid(ncid, 'plev', id_dim)
103 call netcdf_err(error, 'reading plev id')
104 error=nf90_inquire_dimension(ncid,id_dim,len=idum(3))
105 call netcdf_err(error, 'reading plev')
106 endif
107
108 call esmf_vmbroadcast(vm, idum, 3, 0, rc=rc)
109 if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
110 call error_handler("IN VMBroadcast", rc)
111
112 j_thomp_mp_climo = idum(1)
113 i_thomp_mp_climo = idum(2)
114 lev_thomp_mp_climo = idum(3)
115
116 allocate(lons(i_thomp_mp_climo))
117 allocate(lats(j_thomp_mp_climo))
118
119 if (localpet == 0) then
120 error=nf90_inq_varid(ncid, 'lon', id_var)
121 call netcdf_err(error, 'reading lon field id' )
122 error=nf90_get_var(ncid, id_var, lons)
123 call netcdf_err(error, 'reading grid longitude' )
124 error=nf90_inq_varid(ncid, 'lat', id_var)
125 call netcdf_err(error, 'reading lat field id' )
126 error=nf90_get_var(ncid, id_var, lats)
127 call netcdf_err(error, 'reading grid latitude' )
128 endif
129
130 call esmf_vmbroadcast(vm, lons, i_thomp_mp_climo, 0, rc=rc)
131 if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
132 call error_handler("IN VMGet", rc)
133
134 call esmf_vmbroadcast(vm, lats, j_thomp_mp_climo, 0, rc=rc)
135 if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
136 call error_handler("IN VMGet", rc)
137
138!-----------------------------------------------------------------------------------
139! Now that we have the grid information, create the esmf grid object.
140!-----------------------------------------------------------------------------------
141
142 polekindflag(1:2) = esmf_polekind_monopole
143
144 print*,"- CALL GridCreate1PeriDim FOR THOMP MP CLIMO GRID."
145 thomp_mp_climo_grid = esmf_gridcreate1peridim(minindex=(/1,1/), &
146 maxindex=(/i_thomp_mp_climo,j_thomp_mp_climo/), &
147 polekindflag=polekindflag, &
148 periodicdim=1, &
149 poledim=2, &
150 coordsys=esmf_coordsys_sph_deg, &
151 regdecomp=(/1,npets/), &
152 indexflag=esmf_index_global, rc=rc)
153 if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
154 call error_handler("IN GridCreate1PeriDim", rc)
155
156 print*,"- CALL GridAddCoord FOR THOMP MP CLIMO GRID."
157 call esmf_gridaddcoord(thomp_mp_climo_grid, &
158 staggerloc=esmf_staggerloc_center, rc=rc)
159 if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
160 call error_handler("IN GridAddCoord", rc)
161
162!-----------------------------------------------------------------------------------
163! Set the grid object lat/lon.
164!-----------------------------------------------------------------------------------
165
166 print*,"- CALL GridGetCoord FOR INPUT GRID X-COORD."
167 nullify(lon_ptr)
168 call esmf_gridgetcoord(thomp_mp_climo_grid, &
169 staggerloc=esmf_staggerloc_center, &
170 coorddim=1, &
171 farrayptr=lon_ptr, rc=rc)
172 if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
173 call error_handler("IN GridGetCoord", rc)
174
175 print*,"- CALL GridGetCoord FOR INPUT GRID Y-COORD."
176 nullify(lat_ptr)
177 call esmf_gridgetcoord(thomp_mp_climo_grid, &
178 staggerloc=esmf_staggerloc_center, &
179 coorddim=2, &
180 computationallbound=clb, &
181 computationalubound=cub, &
182 farrayptr=lat_ptr, rc=rc)
183 if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
184 call error_handler("IN GridGetCoord", rc)
185
186 do i = clb(1), cub(1)
187 lon_ptr(i,:) = lons(i)
188 enddo
189
190 do j = clb(2), cub(2)
191 lat_ptr(:,j) = lats(j)
192 enddo
193
194!-----------------------------------------------------------------------------------
195! Create esmf fields for the two tracers and 3-d pressure.
196!-----------------------------------------------------------------------------------
197
198 print*,"- CALL FieldCreate FOR QNIFA INPUT CLIMO."
199 qnifa_climo_input_grid = esmf_fieldcreate(thomp_mp_climo_grid, &
200 typekind=esmf_typekind_r8, &
201 staggerloc=esmf_staggerloc_center, &
202 ungriddedlbound=(/1/), &
203 ungriddedubound=(/lev_thomp_mp_climo/), rc=rc)
204 if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
205 call error_handler("IN FieldCreate", rc)
206
207 print*,"- CALL FieldCreate FOR QNWFA INPUT CLIMO."
208 qnwfa_climo_input_grid = esmf_fieldcreate(thomp_mp_climo_grid, &
209 typekind=esmf_typekind_r8, &
210 staggerloc=esmf_staggerloc_center, &
211 ungriddedlbound=(/1/), &
212 ungriddedubound=(/lev_thomp_mp_climo/), rc=rc)
213 if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
214 call error_handler("IN FieldCreate", rc)
215
216 print*,"- CALL FieldCreate FOR THOMP PRESS CLIMO."
217 thomp_pres_climo_input_grid = esmf_fieldcreate(thomp_mp_climo_grid, &
218 typekind=esmf_typekind_r8, &
219 staggerloc=esmf_staggerloc_center, &
220 ungriddedlbound=(/1/), &
221 ungriddedubound=(/lev_thomp_mp_climo/), rc=rc)
222 if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
223 call error_handler("IN FieldCreate", rc)
224
225!-----------------------------------------------------------------------------------
226! Data are monthly and valid at the 15th of the month. Compute time interpolation
227! weights for the current cycle.
228!-----------------------------------------------------------------------------------
229
230 jda=0
231 jda(1) = 2007
232 if (cycle_mon == 2 .and. cycle_day == 29) then ! leap year
233 jda(2) = 3
234 jda(3) = 1
235 else
236 jda(2) = cycle_mon
237 jda(3) = cycle_day
238 endif
239
240 jda(5) = cycle_hour
241 jdow = 0
242 jdoy = 0
243 jday = 0
244 call w3doxdat(jda,jdow,jdoy,jday)
245 rjday = float(jdoy) + float(jda(5)) / 24.
246 if(rjday < dayhf(1)) rjday = rjday + 365.
247
248 do mm=1,12
249 mmm = mm
250 mmp = mm + 1
251 if(rjday >= dayhf(mmm) .and. rjday < dayhf(mmp)) then
252 mon1 = mmm
253 mon2 = mmp
254 exit
255 endif
256 enddo
257
258 wei1m = (dayhf(mon2)-rjday)/(dayhf(mon2)-dayhf(mon1))
259 wei2m = 1.0 - wei1m
260
261 if (mon2==13) mon2=1
262
263 print*,"- BOUNDING MONTHS AND INTERPOLATION WEIGHTS: ", mon1, wei1m, mon2, wei2m
264
265!-----------------------------------------------------------------------------------
266! Read tracers and 3-d pressure for each bounding month. Then linearly
267! interpolate in time.
268!-----------------------------------------------------------------------------------
269
270 if (localpet == 0) then
271 allocate(dummy3d(i_thomp_mp_climo, j_thomp_mp_climo, lev_thomp_mp_climo))
272 dummy3d = 0.0
273 allocate(dummy3d_mon1(i_thomp_mp_climo, j_thomp_mp_climo, lev_thomp_mp_climo))
274 dummy3d_mon1 = 0.0
275 allocate(dummy3d_mon2(i_thomp_mp_climo, j_thomp_mp_climo, lev_thomp_mp_climo))
276 dummy3d_mon2 = 0.0
277 else
278 allocate(dummy3d(0,0,0))
279 allocate(dummy3d_mon1(0,0,0))
280 allocate(dummy3d_mon2(0,0,0))
281 endif
282
283 if (localpet == 0) then
284 print*,"- READ QNIFA FOR BOUNDING MONTH 1"
285 error=nf90_inq_varid(ncid, 'nifa', id_var)
286 call netcdf_err(error, 'reading nifa field id' )
287 error=nf90_get_var(ncid, id_var, dummy3d_mon1, start=(/1,1,1,mon1/), &
288 count=(/i_thomp_mp_climo,j_thomp_mp_climo,lev_thomp_mp_climo,1/) )
289 call netcdf_err(error, 'reading nifa month1 field' )
290 print*,"- READ QNIFA FOR BOUNDING MONTH 2"
291 error=nf90_get_var(ncid, id_var, dummy3d_mon2, start=(/1,1,1,mon2/), &
292 count=(/i_thomp_mp_climo,j_thomp_mp_climo,lev_thomp_mp_climo,1/) )
293 call netcdf_err(error, 'reading nifa month2 field' )
294 dummy3d(:,:,:) = wei1m * dummy3d_mon1 + wei2m * dummy3d_mon2
295 endif
296
297 print*,"- CALL FieldScatter FOR qnifa input climo."
298 call esmf_fieldscatter(qnifa_climo_input_grid, dummy3d, rootpet=0, rc=rc)
299 if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
300 call error_handler("IN FieldScatter", rc)
301
302 if (localpet == 0) then
303 print*,"- READ QNWFA FOR BOUNDING MONTH 1"
304 error=nf90_inq_varid(ncid, 'nwfa', id_var)
305 call netcdf_err(error, 'reading nwfa field id' )
306 error=nf90_get_var(ncid, id_var, dummy3d_mon1, start=(/1,1,1,mon1/), &
307 count=(/i_thomp_mp_climo,j_thomp_mp_climo,lev_thomp_mp_climo,1/) )
308 call netcdf_err(error, 'reading nwfa month1 field' )
309 print*,"- READ QNWFA FOR BOUNDING MONTH 2"
310 error=nf90_get_var(ncid, id_var, dummy3d_mon2, start=(/1,1,1,mon2/), &
311 count=(/i_thomp_mp_climo,j_thomp_mp_climo,lev_thomp_mp_climo,1/) )
312 call netcdf_err(error, 'reading nwfa month2 field' )
313 dummy3d(:,:,:) = wei1m * dummy3d_mon1 + wei2m * dummy3d_mon2
314 endif
315
316 print*,"- CALL FieldScatter FOR qnwfa input climo."
317 call esmf_fieldscatter(qnwfa_climo_input_grid, dummy3d, rootpet=0, rc=rc)
318 if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
319 call error_handler("IN FieldScatter", rc)
320
321 if (localpet == 0) then
322 print*,"- READ PRESSURE FOR BOUNDING MONTH 1"
323 error=nf90_inq_varid(ncid, 'prs', id_var)
324 call netcdf_err(error, 'reading prs field id' )
325 error=nf90_get_var(ncid, id_var, dummy3d_mon1, start=(/1,1,1,mon1/), &
326 count=(/i_thomp_mp_climo,j_thomp_mp_climo,lev_thomp_mp_climo,1/) )
327 call netcdf_err(error, 'reading prs month1 field' )
328 print*,"- READ PRESSURE FOR BOUNDING MONTH 2"
329 error=nf90_get_var(ncid, id_var, dummy3d_mon2, start=(/1,1,1,mon2/), &
330 count=(/i_thomp_mp_climo,j_thomp_mp_climo,lev_thomp_mp_climo,1/) )
331 call netcdf_err(error, 'reading prs month2 field' )
332 dummy3d(:,:,:) = wei1m * dummy3d_mon1 + wei2m * dummy3d_mon2
333 error=nf90_close(ncid)
334 endif
335
336 print*,"- CALL FieldScatter FOR thomp press."
337 call esmf_fieldscatter(thomp_pres_climo_input_grid, dummy3d, rootpet=0, rc=rc)
338 if(esmf_logfounderror(rctocheck=rc,msg=esmf_logerr_passthru,line=__line__,file=__file__)) &
339 call error_handler("IN FieldScatter", rc)
340
341 deallocate(lons, lats, dummy3d, dummy3d_mon1, dummy3d_mon2)
342
343 end subroutine read_thomp_mp_climo_data
344
348 subroutine cleanup_thomp_mp_climo_input_data
349
350 implicit none
351
352 integer :: rc
353
354 call esmf_griddestroy(thomp_mp_climo_grid, rc=rc)
355 call esmf_fielddestroy(thomp_pres_climo_input_grid, rc=rc)
356 call esmf_fielddestroy(qnifa_climo_input_grid, rc=rc)
357 call esmf_fielddestroy(qnwfa_climo_input_grid, rc=rc)
358
359 end subroutine cleanup_thomp_mp_climo_input_data
360
361 end module thompson_mp_climo_data
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::thomp_mp_climo_file
character(len=500), public thomp_mp_climo_file
Path/name to the Thompson MP climatology file.
Definition program_setup.F90:44
program_setup::cycle_mon
integer, public cycle_mon
Cycle month.
Definition program_setup.F90:92
program_setup::cycle_day
integer, public cycle_day
Cycle day.
Definition program_setup.F90:93
program_setup::cycle_hour
integer, public cycle_hour
Cycle hour.
Definition program_setup.F90:94
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::j_thomp_mp_climo
integer j_thomp_mp_climo
j-dimension of Thompson climo data
Definition thompson_mp_climo_data.F90:23
thompson_mp_climo_data::thomp_mp_climo_grid
type(esmf_grid) thomp_mp_climo_grid
esmf grid object for Thompson data grid
Definition thompson_mp_climo_data.F90:28
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
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:349
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
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
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
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
thompson_mp_climo_data::i_thomp_mp_climo
integer i_thomp_mp_climo
i-dimension of Thompson climo data
Definition thompson_mp_climo_data.F90:21