Context Navigation

cgutils.F90 @ 9411952

FPv9.3.2

Last change on this file since 9411952 was 9411952, checked in by Don Morton <Don.Morton@…>, 3 years ago
Synchronised repo with the current CTBTO master. The sole change is the addition of the flexpart_code/checkgrib/ dir
Property mode set to `100644`
File size: 70.8 KB

Line
1	MODULE cgutils
2
3	IMPLICIT NONE
4
5	PRIVATE
6
7	! All publicly accessible functions and variables need to be
8	! declared as such here. All public names are prefixed with
9	! "gributils_"
10	PUBLIC :: cgutils_check_ecmwf_grib, &
11	& cgutils_check_ncep_grib2, &
12	& cgutils_check_fv3_grib2
13
14	! These codes depict the origin/model of the GRIB file
15	INTEGER, PARAMETER :: GRIB_CENTRE_NCEP = 7, &
16	& GRIB_CENTRE_ECMWF = 98
17
18
19	! Return codes
20	INTEGER, PARAMETER :: RETURN_CODE_NORMAL = 0, &
21	& RETURN_CODE_GRIBFILE_ERROR = 1, &
22	& RETURN_CODE_UTILITY_ERROR = 2
23
24	! String lengths
25	INTEGER, PARAMETER :: GRIB_ATTRIBUTE_STR_LEN = 64
26
27	! Global (within this module) variables
28	! All global variables are prefixed with "global_"
29
30	! 3D variables of interest, T/F flag for each level present
31	!
32	! Variables names correspond more to "Flexpart" names than
33	! actual GRIB names
34	!
35	! ECMWF uses - u, v, t, q, etadot
36	!
37	! NCEP uses - u, v, t, rh, w
38	LOGICAL, DIMENSION(:), ALLOCATABLE :: global_levels_u, &
39	& global_levels_v, &
40	& global_levels_t, &
41	& global_levels_q, &
42	& global_levels_etadot, &
43	& global_levels_w, &
44	& global_levels_rh
45
46
47	! 2D variables of interest, T/F flag for each
48	!
49	! Variables names correspond more to "Flexpart" names than
50	! actual GRIB names
51	!
52	! ECMWF uses - ps, sd, msl, tcc, u10, v10, t2, td2, lsprec,
53	! convprec, shf, sr, ewss, nsss, oro, excessoro, lsm
54	!
55	! NCEP uses - ps, sd, msl, u10, v10, t2, td2, tcc, lsprec, convprec,
56	! oro, lsm, hmix, rh2, tsig1, usig1, vsig1
57	LOGICAL :: global_2dvar_ps, &
58	& global_2dvar_sd, &
59	& global_2dvar_msl, &
60	& global_2dvar_tcc, &
61	& global_2dvar_u10, &
62	& global_2dvar_v10, &
63	& global_2dvar_t2, &
64	& global_2dvar_td2, &
65	& global_2dvar_lsprec, &
66	& global_2dvar_convprec, &
67	& global_2dvar_shf, &
68	& global_2dvar_sr, &
69	& global_2dvar_ewss, &
70	& global_2dvar_nsss, &
71	& global_2dvar_oro, &
72	& global_2dvar_excessoro, &
73	& global_2dvar_lsm, &
74	& global_2dvar_hmix, &
75	& global_2dvar_rh2, &
76	& global_2dvar_tsig1, &
77	& global_2dvar_usig1, &
78	& global_2dvar_vsig1
79
80
81	CONTAINS
82
83	SUBROUTINE allocate_and_init_flags(expected_num_levels)
84	IMPLICIT NONE
85
86	! This is used by ECMWF, NCEP and FV3 routines. There are
87	! several variables used by one, that isn't used by the
88	! other, but I think the memory load is negligible and
89	! prefer to take this opportunity to present a little bit
90	! of clean, reusable code.
91
92	INTEGER, INTENT(IN) :: expected_num_levels
93
94	! 3D variables (NCEP, FV3 and ECMWF)
95	ALLOCATE(global_levels_u(expected_num_levels))
96	global_levels_u = .FALSE.
97
98	ALLOCATE(global_levels_v(expected_num_levels))
99	global_levels_v = .FALSE.
100
101	ALLOCATE(global_levels_t(expected_num_levels))
102	global_levels_t = .FALSE.
103
104	ALLOCATE(global_levels_q(expected_num_levels))
105	global_levels_q = .FALSE.
106
107	ALLOCATE(global_levels_etadot(expected_num_levels))
108	global_levels_etadot = .FALSE.
109
110	ALLOCATE(global_levels_w(expected_num_levels))
111	global_levels_w = .FALSE.
112
113	ALLOCATE(global_levels_rh(expected_num_levels))
114	global_levels_rh = .FALSE.
115
116
117
118	! 2D variables (both NCEP, ECMWF and FV3)
119	global_2dvar_ps = .FALSE.
120	global_2dvar_sd = .FALSE.
121	global_2dvar_msl = .FALSE.
122	global_2dvar_tcc = .FALSE.
123	global_2dvar_u10 = .FALSE.
124	global_2dvar_v10 = .FALSE.
125	global_2dvar_t2 = .FALSE.
126	global_2dvar_td2 = .FALSE.
127	global_2dvar_lsprec = .FALSE.
128	global_2dvar_convprec = .FALSE.
129	global_2dvar_shf = .FALSE.
130	global_2dvar_sr = .FALSE.
131	global_2dvar_ewss = .FALSE.
132	global_2dvar_nsss = .FALSE.
133	global_2dvar_oro = .FALSE.
134	global_2dvar_excessoro = .FALSE.
135	global_2dvar_lsm = .FALSE.
136	global_2dvar_hmix = .FALSE.
137	global_2dvar_rh2 = .FALSE.
138	global_2dvar_tsig1 = .FALSE.
139	global_2dvar_usig1 = .FALSE.
140	global_2dvar_vsig1 = .FALSE.
141
142	END SUBROUTINE allocate_and_init_flags
143
144
145
146	SUBROUTINE delete_levels_flags
147	IMPLICIT NONE
148	DEALLOCATE(global_levels_u)
149	DEALLOCATE(global_levels_v)
150	DEALLOCATE(global_levels_t)
151	DEALLOCATE(global_levels_q)
152	DEALLOCATE(global_levels_etadot)
153	DEALLOCATE(global_levels_w)
154	DEALLOCATE(global_levels_rh)
155	END SUBROUTINE delete_levels_flags
156
157
158
159	LOGICAL FUNCTION ecmwf_all_2dvars_present()
160	IMPLICIT NONE
161
162	! Checks logical flags for each 2d variable, testing for
163	! all true.
164	!
165	! For sake of "situational awareness," this routine will
166	! go through the logical flags of each 2d array of interest,
167	! printing error message for each one with a problem, then
168	! returning the error code after all that.
169
170	ecmwf_all_2dvars_present = .TRUE.
171
172	! surface levels
173	ecmwf_all_2dvars_present = &
174	& (ecmwf_all_2dvars_present .AND. global_2dvar_ps)
175	IF ( .NOT. global_2dvar_ps) &
176	& WRITE(6,*) 'ERROR: sp: not found'
177
178	ecmwf_all_2dvars_present = &
179	& (ecmwf_all_2dvars_present .AND. global_2dvar_sd)
180	IF ( .NOT. global_2dvar_sd) &
181	& WRITE(6,*) 'ERROR: sd (or sde): not found'
182
183	ecmwf_all_2dvars_present = &
184	& (ecmwf_all_2dvars_present .AND. global_2dvar_tcc)
185	IF ( .NOT. global_2dvar_tcc) &
186	& WRITE(6,*) 'ERROR: tcc: not found'
187
188	ecmwf_all_2dvars_present = &
189	& (ecmwf_all_2dvars_present .AND. global_2dvar_lsprec)
190	IF ( .NOT. global_2dvar_lsprec) &
191	& WRITE(6,*) 'ERROR: lsp: not found'
192
193	ecmwf_all_2dvars_present = &
194	& (ecmwf_all_2dvars_present .AND. global_2dvar_convprec)
195	IF ( .NOT. global_2dvar_convprec) &
196	& WRITE(6,*) 'ERROR: acpcp: not found'
197
198	ecmwf_all_2dvars_present = &
199	& (ecmwf_all_2dvars_present .AND. global_2dvar_shf)
200	IF ( .NOT. global_2dvar_shf) &
201	& WRITE(6,*) 'ERROR: sshf: not found'
202
203	ecmwf_all_2dvars_present = &
204	& (ecmwf_all_2dvars_present .AND. global_2dvar_sr)
205	IF ( .NOT. global_2dvar_sr) &
206	& WRITE(6,*) 'ERROR: ssr: not found'
207
208	ecmwf_all_2dvars_present = &
209	& (ecmwf_all_2dvars_present .AND. global_2dvar_ewss)
210	IF ( .NOT. global_2dvar_ewss) &
211	& WRITE(6,*) 'ERROR: ewss: not found'
212
213	ecmwf_all_2dvars_present = &
214	& (ecmwf_all_2dvars_present .AND. global_2dvar_nsss)
215	IF ( .NOT. global_2dvar_nsss) &
216	& WRITE(6,*) 'ERROR: nsss: not found'
217
218	ecmwf_all_2dvars_present = &
219	& (ecmwf_all_2dvars_present .AND. global_2dvar_oro)
220	IF ( .NOT. global_2dvar_oro) &
221	& WRITE(6,*) 'ERROR: z: not found'
222
223	ecmwf_all_2dvars_present = &
224	& (ecmwf_all_2dvars_present .AND. global_2dvar_excessoro)
225	IF ( .NOT. global_2dvar_excessoro) &
226	& WRITE(6,*) 'ERROR: sdor: not found'
227
228	ecmwf_all_2dvars_present = &
229	& (ecmwf_all_2dvars_present .AND. global_2dvar_lsm)
230	IF ( .NOT. global_2dvar_lsm) &
231	& WRITE(6,*) 'ERROR: lsm: not found'
232
233
234	! heightAboveGround levels
235	ecmwf_all_2dvars_present = &
236	& (ecmwf_all_2dvars_present .AND. global_2dvar_u10)
237	IF ( .NOT. global_2dvar_u10) &
238	& WRITE(6,*) 'ERROR: 10u: not found'
239
240	ecmwf_all_2dvars_present = &
241	& (ecmwf_all_2dvars_present .AND. global_2dvar_v10)
242	IF ( .NOT. global_2dvar_v10) &
243	& WRITE(6,*) 'ERROR: 10v: not found'
244
245	ecmwf_all_2dvars_present = &
246	& (ecmwf_all_2dvars_present .AND. global_2dvar_t2)
247	IF ( .NOT. global_2dvar_t2) &
248	& WRITE(6,*) 'ERROR: 2t: not found'
249
250	ecmwf_all_2dvars_present = &
251	& (ecmwf_all_2dvars_present .AND. global_2dvar_td2)
252	IF ( .NOT. global_2dvar_td2) &
253	& WRITE(6,*) 'ERROR: 2d: not found'
254
255
256	! meanSea levels
257	ecmwf_all_2dvars_present = &
258	& (ecmwf_all_2dvars_present .AND. global_2dvar_msl)
259	IF ( .NOT. global_2dvar_msl) &
260	& WRITE(6,*) 'ERROR: msl: not found'
261
262
263	RETURN
264
265
266
267	END FUNCTION ecmwf_all_2dvars_present
268
269
270	LOGICAL FUNCTION ncep_all_2dvars_present()
271	IMPLICIT NONE
272
273	! Checks logical flags for each 2d variable, testing for
274	! all true.
275	!
276	! For sake of "situational awareness," this routine will
277	! go through the logical flags of each 2d array of interest,
278	! printing error message for each one with a problem, then
279	! returning the error code after all that.
280
281	ncep_all_2dvars_present = .TRUE.
282
283	! surface levels
284	ncep_all_2dvars_present = &
285	& (ncep_all_2dvars_present .AND. global_2dvar_ps)
286	IF ( .NOT. global_2dvar_ps) &
287	& WRITE(6,*) 'ERROR: sp: not found'
288
289	ncep_all_2dvars_present = &
290	& (ncep_all_2dvars_present .AND. global_2dvar_sd)
291	IF ( .NOT. global_2dvar_sd) &
292	& WRITE(6,*) 'ERROR: sdwe: not found'
293
294
295
296
297	!!! 2019-01-20 - These two variables are currently not part of what I find in
298	!!! NCEP GRIB files at CTBTO. I think they used to be available. So, I'm
299	!!! commenting out their testing right now.
300	! ncep_all_2dvars_present = &
301	!& (ncep_all_2dvars_present .AND. global_2dvar_lsprec)
302	! IF ( .NOT. global_2dvar_lsprec) &
303	!& WRITE(6,*) 'ERROR: tp: not found'
304	!
305	! ncep_all_2dvars_present = &
306	!& (ncep_all_2dvars_present .AND. global_2dvar_convprec)
307	! IF ( .NOT. global_2dvar_convprec) &
308	!& WRITE(6,*) 'ERROR: acpcp: not found'
309
310	ncep_all_2dvars_present = &
311	& (ncep_all_2dvars_present .AND. global_2dvar_oro)
312	IF ( .NOT. global_2dvar_oro) &
313	& WRITE(6,*) 'ERROR: orog: not found'
314
315	ncep_all_2dvars_present = &
316	& (ncep_all_2dvars_present .AND. global_2dvar_lsm)
317	IF ( .NOT. global_2dvar_lsm) &
318	& WRITE(6,*) 'ERROR: lsm: not found'
319
320	ncep_all_2dvars_present = &
321	& (ncep_all_2dvars_present .AND. global_2dvar_hmix)
322	IF ( .NOT. global_2dvar_hmix) &
323	& WRITE(6,*) 'ERROR: hpbl: not found'
324
325	! heightAboveGround levels
326	ncep_all_2dvars_present = &
327	& (ncep_all_2dvars_present .AND. global_2dvar_u10)
328	IF ( .NOT. global_2dvar_u10) &
329	& WRITE(6,*) 'ERROR: 10u: not found'
330
331	ncep_all_2dvars_present = &
332	& (ncep_all_2dvars_present .AND. global_2dvar_v10)
333	IF ( .NOT. global_2dvar_v10) &
334	& WRITE(6,*) 'ERROR: 10v: not found'
335
336	ncep_all_2dvars_present = &
337	& (ncep_all_2dvars_present .AND. global_2dvar_t2)
338	IF ( .NOT. global_2dvar_t2) &
339	& WRITE(6,*) 'ERROR: 2t: not found'
340
341	ncep_all_2dvars_present = &
342	& (ncep_all_2dvars_present .AND. global_2dvar_td2)
343	IF ( .NOT. global_2dvar_td2) &
344	& WRITE(6,*) 'ERROR: 2d: not found'
345
346	ncep_all_2dvars_present = &
347	& (ncep_all_2dvars_present .AND. global_2dvar_rh2)
348	IF ( .NOT. global_2dvar_rh2) &
349	& WRITE(6,*) 'ERROR: 2m r: not found'
350
351
352	! meanSea levels
353	ncep_all_2dvars_present = &
354	& (ncep_all_2dvars_present .AND. global_2dvar_msl)
355	IF ( .NOT. global_2dvar_msl) &
356	& WRITE(6,*) 'ERROR: prmsl: not found'
357
358	! sigma levels
359	ncep_all_2dvars_present = &
360	& (ncep_all_2dvars_present .AND. global_2dvar_tsig1)
361	IF ( .NOT. global_2dvar_tsig1) &
362	& WRITE(6,*) 'ERROR: sigma t: not found'
363
364	ncep_all_2dvars_present = &
365	& (ncep_all_2dvars_present .AND. global_2dvar_usig1)
366	IF ( .NOT. global_2dvar_usig1) &
367	& WRITE(6,*) 'ERROR: sigma u: not found'
368
369	ncep_all_2dvars_present = &
370	& (ncep_all_2dvars_present .AND. global_2dvar_vsig1)
371	IF ( .NOT. global_2dvar_vsig1) &
372	& WRITE(6,*) 'ERROR: sigma v: not found'
373
374
375
376
377
378
379	!!! 2019-01-20 - These variable is currently not part of what I find in
380	!!! NCEP GRIB files at CTBTO. I think it used to be available. So, I'm
381	!!! commenting out the testing right now.
382	! Low cloud layer levels
383	! ncep_all_2dvars_present = &
384	!& (ncep_all_2dvars_present .AND. global_2dvar_tcc)
385	! IF ( .NOT. global_2dvar_tcc) &
386	!& WRITE(6,*) 'ERROR: tcc: not found'
387
388
389	RETURN
390
391
392
393	END FUNCTION ncep_all_2dvars_present
394
395
396	LOGICAL FUNCTION fv3_all_2dvars_present()
397	IMPLICIT NONE
398
399	! Checks logical flags for each 2d variable, testing for
400	! all true.
401	!
402	! For sake of "situational awareness," this routine will
403	! go through the logical flags of each 2d array of interest,
404	! printing error message for each one with a problem, then
405	! returning the error code after all that.
406
407	fv3_all_2dvars_present = .TRUE.
408
409	! surface levels
410	fv3_all_2dvars_present = &
411	& (fv3_all_2dvars_present .AND. global_2dvar_ps)
412	IF ( .NOT. global_2dvar_ps) &
413	& WRITE(6,*) 'ERROR: sp: not found'
414
415	fv3_all_2dvars_present = &
416	& (fv3_all_2dvars_present .AND. global_2dvar_sd)
417	IF ( .NOT. global_2dvar_sd) &
418	& WRITE(6,*) 'ERROR: sdwe: not found'
419
420
421
422
423	!!! 2019-01-20 - These two variables are currently not part of what I find in
424	!!! NCEP GRIB files at CTBTO. I think they used to be available. So, I'm
425	!!! commenting out their testing right now.
426	! ncep_all_2dvars_present = &
427	!& (ncep_all_2dvars_present .AND. global_2dvar_lsprec)
428	! IF ( .NOT. global_2dvar_lsprec) &
429	!& WRITE(6,*) 'ERROR: tp: not found'
430	!
431	! ncep_all_2dvars_present = &
432	!& (ncep_all_2dvars_present .AND. global_2dvar_convprec)
433	! IF ( .NOT. global_2dvar_convprec) &
434	!& WRITE(6,*) 'ERROR: acpcp: not found'
435
436	fv3_all_2dvars_present = &
437	& (fv3_all_2dvars_present .AND. global_2dvar_oro)
438	IF ( .NOT. global_2dvar_oro) &
439	& WRITE(6,*) 'ERROR: orog: not found'
440
441	fv3_all_2dvars_present = &
442	& (fv3_all_2dvars_present .AND. global_2dvar_lsm)
443	IF ( .NOT. global_2dvar_lsm) &
444	& WRITE(6,*) 'ERROR: lsm: not found'
445
446	fv3_all_2dvars_present = &
447	& (fv3_all_2dvars_present .AND. global_2dvar_hmix)
448	IF ( .NOT. global_2dvar_hmix) &
449	& WRITE(6,*) 'ERROR: hpbl: not found'
450
451	! heightAboveGround levels
452	fv3_all_2dvars_present = &
453	& (fv3_all_2dvars_present .AND. global_2dvar_u10)
454	IF ( .NOT. global_2dvar_u10) &
455	& WRITE(6,*) 'ERROR: 10u: not found'
456
457	fv3_all_2dvars_present = &
458	& (fv3_all_2dvars_present .AND. global_2dvar_v10)
459	IF ( .NOT. global_2dvar_v10) &
460	& WRITE(6,*) 'ERROR: 10v: not found'
461
462	fv3_all_2dvars_present = &
463	& (fv3_all_2dvars_present .AND. global_2dvar_t2)
464	IF ( .NOT. global_2dvar_t2) &
465	& WRITE(6,*) 'ERROR: 2t: not found'
466
467	fv3_all_2dvars_present = &
468	& (fv3_all_2dvars_present .AND. global_2dvar_td2)
469	IF ( .NOT. global_2dvar_td2) &
470	& WRITE(6,*) 'ERROR: 2d: not found'
471
472	fv3_all_2dvars_present = &
473	& (fv3_all_2dvars_present .AND. global_2dvar_rh2)
474	IF ( .NOT. global_2dvar_rh2) &
475	& WRITE(6,*) 'ERROR: 2m r: not found'
476
477
478	! meanSea levels
479	fv3_all_2dvars_present = &
480	& (fv3_all_2dvars_present .AND. global_2dvar_msl)
481	IF ( .NOT. global_2dvar_msl) &
482	& WRITE(6,*) 'ERROR: prmsl: not found'
483
484	! sigma levels
485	fv3_all_2dvars_present = &
486	& (fv3_all_2dvars_present .AND. global_2dvar_tsig1)
487	IF ( .NOT. global_2dvar_tsig1) &
488	& WRITE(6,*) 'ERROR: sigma t: not found'
489
490	fv3_all_2dvars_present = &
491	& (fv3_all_2dvars_present .AND. global_2dvar_usig1)
492	IF ( .NOT. global_2dvar_usig1) &
493	& WRITE(6,*) 'ERROR: sigma u: not found'
494
495	fv3_all_2dvars_present = &
496	& (fv3_all_2dvars_present .AND. global_2dvar_vsig1)
497	IF ( .NOT. global_2dvar_vsig1) &
498	& WRITE(6,*) 'ERROR: sigma v: not found'
499
500
501
502
503
504
505	!!! 2019-01-20 - These variable is currently not part of what I find in
506	!!! NCEP GRIB files at CTBTO. I think it used to be available. So, I'm
507	!!! commenting out the testing right now.
508	! Low cloud layer levels
509	! ncep_all_2dvars_present = &
510	!& (ncep_all_2dvars_present .AND. global_2dvar_tcc)
511	! IF ( .NOT. global_2dvar_tcc) &
512	!& WRITE(6,*) 'ERROR: tcc: not found'
513
514
515	RETURN
516
517
518
519	END FUNCTION fv3_all_2dvars_present
520
521
522	LOGICAL FUNCTION ecmwf_all_levels_present()
523	IMPLICIT NONE
524
525	! Checks logical flags for each 3d variable, testing for
526	! all true.
527	!
528	! For sake of "situational awareness," this routine will
529	! go through the logical flags of each 3d array of interest,
530	! printing error message for each one with a problem, then
531	! returning the error code after all that.
532
533	ecmwf_all_levels_present = .TRUE.
534
535	ecmwf_all_levels_present = &
536	& (ecmwf_all_levels_present .AND. ALL(global_levels_u))
537	IF ( .NOT. ALL(global_levels_u) ) &
538	& WRITE(6,*) 'ERROR: u: not all levels found'
539
540	ecmwf_all_levels_present = &
541	& (ecmwf_all_levels_present .AND. ALL(global_levels_v))
542	IF ( .NOT. ALL(global_levels_v) ) &
543	& WRITE(6,*) 'ERROR: v: not all levels found'
544
545	ecmwf_all_levels_present = &
546	& (ecmwf_all_levels_present .AND. ALL(global_levels_t))
547	IF ( .NOT. ALL(global_levels_t) ) &
548	& WRITE(6,*) 'ERROR: t: not all levels found'
549
550	ecmwf_all_levels_present = &
551	& (ecmwf_all_levels_present .AND. ALL(global_levels_q))
552	IF ( .NOT. ALL(global_levels_q) ) &
553	& WRITE(6,*) 'ERROR: q: not all levels found'
554
555	ecmwf_all_levels_present = &
556	& (ecmwf_all_levels_present .AND. ALL(global_levels_etadot))
557	IF ( .NOT. ALL(global_levels_etadot) ) &
558	& WRITE(6,*) 'ERROR: etadot: not all levels found'
559
560
561	END FUNCTION ecmwf_all_levels_present
562
563
564
565	LOGICAL FUNCTION ncep_all_levels_present()
566	IMPLICIT NONE
567
568	! Checks logical flags for each 3d variable, testing for
569	! all true.
570	!
571	! For sake of "situational awareness," this routine will
572	! go through the logical flags of each 3d array of interest,
573	! printing error message for each one with a problem, then
574	! returning the error code after all that.
575
576	ncep_all_levels_present = .TRUE.
577
578	ncep_all_levels_present = &
579	& (ncep_all_levels_present .AND. ALL(global_levels_u))
580	IF ( .NOT. ALL(global_levels_u) ) &
581	& WRITE(6,*) 'ERROR: u: not all levels found'
582
583	ncep_all_levels_present = &
584	& (ncep_all_levels_present .AND. ALL(global_levels_v))
585	IF ( .NOT. ALL(global_levels_v) ) &
586	& WRITE(6,*) 'ERROR: v: not all levels found'
587
588	ncep_all_levels_present = &
589	& (ncep_all_levels_present .AND. ALL(global_levels_t))
590	IF ( .NOT. ALL(global_levels_t) ) &
591	& WRITE(6,*) 'ERROR: t: not all levels found'
592
593	ncep_all_levels_present = &
594	& (ncep_all_levels_present .AND. ALL(global_levels_rh))
595	IF ( .NOT. ALL(global_levels_rh) ) &
596	& WRITE(6,*) 'ERROR: r: not all levels found'
597
598	! We can't check w like the others, because it only seems to be included for messages
599	! with pressure levels 100 mb and higher. So, this is checked elsewhere
600
601	! ncep_all_levels_present = &
602	!& (ncep_all_levels_present .AND. ALL(global_levels_w))
603	! IF ( .NOT. ALL(global_levels_w) ) &
604	!& WRITE(6,*) 'ERROR: w: not all levels found'
605
606	RETURN
607
608	END FUNCTION ncep_all_levels_present
609
610
611	LOGICAL FUNCTION fv3_all_levels_present()
612	IMPLICIT NONE
613
614	! Checks logical flags for each 3d variable, testing for
615	! all true.
616	!
617	! For sake of "situational awareness," this routine will
618	! go through the logical flags of each 3d array of interest,
619	! printing error message for each one with a problem, then
620	! returning the error code after all that.
621
622	fv3_all_levels_present = .TRUE.
623
624	fv3_all_levels_present = &
625	& (fv3_all_levels_present .AND. ALL(global_levels_u))
626	IF ( .NOT. ALL(global_levels_u) ) &
627	& WRITE(6,*) 'ERROR: u: not all levels found'
628
629	fv3_all_levels_present = &
630	& (fv3_all_levels_present .AND. ALL(global_levels_v))
631	IF ( .NOT. ALL(global_levels_v) ) &
632	& WRITE(6,*) 'ERROR: v: not all levels found'
633
634	fv3_all_levels_present = &
635	& (fv3_all_levels_present .AND. ALL(global_levels_t))
636	IF ( .NOT. ALL(global_levels_t) ) &
637	& WRITE(6,*) 'ERROR: t: not all levels found'
638
639	fv3_all_levels_present = &
640	& (fv3_all_levels_present .AND. ALL(global_levels_rh))
641	IF ( .NOT. ALL(global_levels_rh) ) &
642	& WRITE(6,*) 'ERROR: r: not all levels found'
643
644	! We can't check w like the others, because it only seems to be included for messages
645	! with pressure levels 100 mb and higher. So, this is checked elsewhere
646
647	! ncep_all_levels_present = &
648	!& (ncep_all_levels_present .AND. ALL(global_levels_w))
649	! IF ( .NOT. ALL(global_levels_w) ) &
650	!& WRITE(6,*) 'ERROR: w: not all levels found'
651
652	RETURN
653
654	END FUNCTION fv3_all_levels_present
655
656
657	LOGICAL FUNCTION ncep_all_expected_w_levels_present(expected_num_levels, millibars)
658	IMPLICIT NONE
659
660
661	! It seems like the w messages are only available at 100mb and below, so I can't check
662	! this variable easily like I can the other 3d variables. I need to explictly go through
663	! each value in the millibars array and, if the value is 100 or greaters, use its index to
664	! determine if the corresponding value in the flag array is .TRUE. or .FALSE.
665
666
667	INTEGER, INTENT(IN) :: expected_num_levels
668	INTEGER, INTENT(IN), DIMENSION(expected_num_levels) :: millibars
669
670
671	INTEGER :: i, plevel
672	LOGICAL :: all_present
673
674	all_present = .TRUE.
675	DO i=1,expected_num_levels
676	plevel = millibars(i)
677	IF (plevel .GE. 100) all_present = (all_present .AND. global_levels_w(i))
678	END DO
679
680	ncep_all_expected_w_levels_present = all_present
681
682	RETURN
683
684	END FUNCTION ncep_all_expected_w_levels_present
685
686	LOGICAL FUNCTION fv3_all_expected_w_levels_present(expected_num_levels, millibars)
687	IMPLICIT NONE
688
689
690	! It seems like the w messages are only available at 100mb and below, so I can't check
691	! this variable easily like I can the other 3d variables. I need to explictly go through
692	! each value in the millibars array and, if the value is 100 or greaters, use its index to
693	! determine if the corresponding value in the flag array is .TRUE. or .FALSE.
694
695
696	INTEGER, INTENT(IN) :: expected_num_levels
697	INTEGER, INTENT(IN), DIMENSION(expected_num_levels) :: millibars
698
699
700	INTEGER :: i, plevel
701	LOGICAL :: all_present
702
703	all_present = .TRUE.
704	DO i=1,expected_num_levels
705	plevel = millibars(i)
706	IF (plevel .GE. 100) all_present = (all_present .AND. global_levels_w(i))
707	END DO
708
709	fv3_all_expected_w_levels_present = all_present
710
711	RETURN
712
713	END FUNCTION fv3_all_expected_w_levels_present
714
715
716
717	LOGICAL FUNCTION message_is_good(igrib, check_anl, check_fcst, check_grib2)
718	IMPLICIT NONE
719
720	! This function is intended to drive several checks on a
721	! grib message and, as requirements change, this function
722	! may become more complex.
723	!
724	! For now, in January 2019, it checks for message being
725	! grib version 2, and, depending on the settings of
726	! check_anl and check_fcst, checks whether the message
727	! is of the specified type
728	!
729	! August 2019 - the checking for GRIB2 message is now conditional
730	!
731	! It is assumed that igrib refers to an already acccessed
732	! and active grib message
733
734	INTEGER, INTENT(IN) :: igrib
735	LOGICAL, INTENT(IN) :: check_anl, check_fcst, check_grib2
736
737	LOGICAL :: fcst_anl_good
738
739	message_is_good = .TRUE.
740
741
742	IF (check_grib2) THEN
743	IF (.NOT. message_is_grib2(igrib)) THEN
744	WRITE(6,*) 'Grib message not GRIB2'
745	message_is_good = .FALSE.
746	RETURN
747	ENDIF
748	ENDIF
749
750	CALL anl_fcst_check(igrib, check_anl, check_fcst, fcst_anl_good)
751	IF (.NOT. fcst_anl_good) THEN
752	WRITE(6,*) 'Failed fcst/anl check'
753	message_is_good = .FALSE.
754	RETURN
755	ENDIF
756
757	END FUNCTION message_is_good
758
759
760
761	LOGICAL FUNCTION message_is_forecast(igrib)
762
763	USE grib_api
764
765	IMPLICIT NONE
766
767	INTEGER, INTENT(IN) :: igrib
768	! igrib is an id for a valid, already-fetched, grib message
769
770	CHARACTER(LEN=8) :: data_type ! Looking for 'fc' or 'an'
771
772	message_is_forecast = .TRUE.
773	CALL grib_get(igrib, 'dataType', data_type)
774	IF (TRIM(data_type) .NE. 'fc') message_is_forecast = .FALSE.
775
776
777	RETURN
778
779	END FUNCTION message_is_forecast
780
781	LOGICAL FUNCTION message_is_analysis(igrib)
782
783	USE grib_api
784
785	IMPLICIT NONE
786
787	INTEGER, INTENT(IN) :: igrib
788	! igrib is an id for a valid, already-fetched, grib message
789
790	CHARACTER(LEN=8) :: data_type ! Looking for 'fc' or 'an'
791
792	message_is_analysis = .TRUE.
793	CALL grib_get(igrib, 'dataType', data_type)
794	IF (data_type .NE. 'an') message_is_analysis = .FALSE.
795
796
797	RETURN
798
799	END FUNCTION message_is_analysis
800
801
802	SUBROUTINE anl_fcst_check(igrib, check_for_anl, check_for_fcst, result)
803
804	USE grib_api
805
806	IMPLICIT NONE
807
808	! If either (but not both) the check_for_anl or check_for_fcst
809	! flag is true, then check that message is good for that
810
811	INTEGER, INTENT(IN) :: igrib
812	! igrib is an id for a valid, already-fetched, grib message
813	LOGICAL, INTENT(IN) :: check_for_anl, check_for_fcst
814	LOGICAL, INTENT(OUT) :: result
815
816	! Only one of the two checks can be true
817	IF (check_for_anl .AND. check_for_fcst) THEN
818	WRITE(0,*) 'Only one of fcst and anl check can be requested'
819	CALL EXIT(RETURN_CODE_UTILITY_ERROR)
820	ENDIF
821
822	! It's possible and valid for this subroutine to be called where
823	! both flags are False. That's OK, just set result to True. No
824	! request was made, so the test passes
825	IF (.NOT.(check_for_anl .OR. check_for_fcst)) THEN
826	result = .TRUE.
827	RETURN
828	ENDIF
829
830	! OK, if we're here, exactly one of the checks is true
831	IF (check_for_anl) THEN
832	result = .TRUE.
833	IF (.NOT. message_is_analysis(igrib)) THEN
834	result = .FALSE.
835	RETURN
836	ENDIF
837	ENDIF
838
839	IF (check_for_fcst) THEN
840	result = .TRUE.
841	IF (.NOT. message_is_forecast(igrib)) THEN
842	result = .FALSE.
843	RETURN
844	ENDIF
845	ENDIF
846
847
848	END SUBROUTINE anl_fcst_check
849
850
851
852
853
854	LOGICAL FUNCTION message_is_grib2(igrib)
855
856	USE grib_api
857
858	IMPLICIT NONE
859
860	INTEGER, INTENT(IN) :: igrib
861	! igrib is an id for a valid, already-fetched, grib message
862
863	INTEGER :: edition_number
864
865	message_is_grib2 = .TRUE.
866	CALL grib_get(igrib, 'editionNumber', edition_number)
867	IF (edition_number .NE. 2) message_is_grib2 = .FALSE.
868
869
870	RETURN
871	END FUNCTION message_is_grib2
872
873
874
875	SUBROUTINE build_millibars_array(ifile, expected_num_levels, &
876	& millibars, return_code)
877
878	USE grib_api
879
880	IMPLICIT NONE
881
882	INTEGER, INTENT(IN) :: ifile, expected_num_levels
883	INTEGER, INTENT(OUT), DIMENSION(expected_num_levels) :: millibars
884	INTEGER, INTENT(OUT) :: return_code
885
886
887
888	! This routine is for the old GFS (pre-FV3) GRIB files, where it is assumed
889	! that all isobaric-level variables will have the same number of levels.
890	!
891	! Goes through each message in the already-open GRIBFILE
892	! indicated by "ifile" and creates a sorted array of all of the
893	! pressures found in isobaric level types. When we're done,
894	! we have an inventory of all of the pressure levels expected by
895	! 3d isobaric-level variables.
896	!
897	! Note that the returned millibars array is NOT sorted. I can't
898	! think of a reason right now why that would be necessary.
899	!
900	! If, in creating the inventory, we end up with fewer or more
901	! levels than expected by "expected_num_levels" we return with
902	! an error
903
904	INTEGER :: level_count ! Keep track of number of new levels found
905	LOGICAL :: end_of_file
906	INTEGER :: igrib, iret
907
908	CHARACTER(LEN=GRIB_ATTRIBUTE_STR_LEN) :: gribmsg_leveltype, &
909	& gribmsg_shortname
910	INTEGER :: gribmsg_level
911
912
913
914	return_code = RETURN_CODE_NORMAL
915	level_count = 0
916
917	millibars = -999
918
919	!PRINT *, 'millibars: ', millibars
920
921	! Iterate through the entire file, keeping track of what we've found
922	end_of_file = .FALSE.
923	DO WHILE (.NOT. end_of_file)
924	CALL grib_new_from_file(ifile, igrib, iret)
925
926
927	IF (iret .eq. GRIB_END_OF_FILE) THEN
928	end_of_file = .TRUE.
929	ELSE
930	CALL grib_get(igrib, 'shortName', gribmsg_shortname)
931	CALL grib_get(igrib, 'typeOfLevel', gribmsg_leveltype)
932	CALL grib_get(igrib, 'level', gribmsg_level)
933
934	! We're only interested in the isobaric level messages
935	IF ( TRIM(gribmsg_leveltype) == 'isobaricInhPa' ) THEN
936
937
938	CALL grib_get(igrib, 'level', gribmsg_level)
939	!PRINT *, 'level: ', gribmsg_level
940
941
942	! If this level is not already in the array, update
943	! the counter (make sure we haven't exceeded expected
944	! number of levels), then add to array.
945	IF (myfindloc(expected_num_levels, &
946	& millibars, VALUE=gribmsg_level) .EQ. 0) THEN
947	level_count = level_count + 1
948	IF (level_count .LE. expected_num_levels) THEN
949	millibars(level_count) = gribmsg_level
950	ELSE
951	WRITE(6,"('ERROR: Found more than ', I3, ' levels')") &
952	& expected_num_levels
953	return_code = RETURN_CODE_GRIBFILE_ERROR
954	RETURN
955	ENDIF
956
957	ENDIF
958
959
960	ENDIF
961	END IF
962
963
964	END DO ! End of main loop for reading each GRIB message
965
966	! Make sure we found the number if isobaric levels that we expected
967	IF (level_count .LT. expected_num_levels) THEN
968	WRITE(6,"('ERROR: Only found ', I3, ' levels')") level_count
969	return_code = RETURN_CODE_GRIBFILE_ERROR
970	RETURN
971	ENDIF
972
973	! OK, if we made it here, we should have a full millibars array
974	! containing all of the isobaric levels we found. NOTE THAT I
975	! AM NOT BOTHERING TO SORT THIS - I DON'T SEE A REASON FOR THAT
976	! SO FAR.
977
978	RETURN
979
980
981	END SUBROUTINE build_millibars_array
982
983
984
985
986
987
988	SUBROUTINE build_fv3_millibars_array(ifile, expected_num_levels, &
989	& millibars, return_code)
990
991	USE grib_api
992
993	IMPLICIT NONE
994
995	INTEGER, INTENT(IN) :: ifile, expected_num_levels
996	INTEGER, INTENT(OUT), DIMENSION(expected_num_levels) :: millibars
997	INTEGER, INTENT(OUT) :: return_code
998
999	! Goes through each message in the already-open GRIBFILE
1000	! indicated by "ifile" and creates a sorted array of all of the
1001	! pressures found in isobaric level types. When we're done,
1002	! we have an inventory of all of the pressure levels expected by
1003	! 3d isobaric-level variables.
1004	!
1005	! Note that the returned millibars array is NOT sorted. I can't
1006	! think of a reason right now why that would be necessary.
1007	!
1008	! If, in creating the inventory, we end up with fewer or more
1009	! levels than expected by "expected_num_levels" we return with
1010	! an error
1011	!
1012	! MODIFICATION - DJM - 2019-07-27 - New NCEP FV3 files have additional
1013	! pressure levels at 15 mb and 40 mb, but the variables at those
1014	! levels are not exactly those at the more standard levels. The
1015	! logic of this routine previously added these levels to the list of
1016	! millibars, but then other routines would fail because they wouldn't
1017	! find expected variables at these levels. So, I have modified this
1018	! routine to simple ignore pressure levels of 15 and 40 millibars.
1019	! This isn't necessarily the most robust way to do things, but it is
1020	! what it is right now.
1021
1022	INTEGER :: level_count ! Keep track of number of new levels found
1023	LOGICAL :: end_of_file
1024	INTEGER :: igrib, iret
1025
1026	CHARACTER(LEN=GRIB_ATTRIBUTE_STR_LEN) :: gribmsg_leveltype, &
1027	& gribmsg_shortname
1028	INTEGER :: gribmsg_level
1029
1030
1031
1032	return_code = RETURN_CODE_NORMAL
1033	level_count = 0
1034
1035	millibars = -999
1036
1037	!PRINT *, 'millibars: ', millibars
1038
1039	! Iterate through the entire file, keeping track of what we've found
1040	end_of_file = .FALSE.
1041	DO WHILE (.NOT. end_of_file)
1042	CALL grib_new_from_file(ifile, igrib, iret)
1043
1044
1045	IF (iret .eq. GRIB_END_OF_FILE) THEN
1046	end_of_file = .TRUE.
1047	ELSE
1048	CALL grib_get(igrib, 'shortName', gribmsg_shortname)
1049	CALL grib_get(igrib, 'typeOfLevel', gribmsg_leveltype)
1050	CALL grib_get(igrib, 'level', gribmsg_level)
1051
1052	! We're only interested in the isobaric level messages
1053	IF ( TRIM(gribmsg_leveltype) == 'isobaricInhPa' ) THEN
1054
1055
1056	CALL grib_get(igrib, 'level', gribmsg_level)
1057	!PRINT *, 'level: ', gribmsg_level
1058
1059	IF ( (gribmsg_level .EQ. 15) .OR. (gribmsg_level .EQ. 40) ) THEN
1060	! Bypass this message - we are ignoring 15 and 40 mb
1061	! pressure levels. See comments above
1062	CONTINUE
1063	ELSE
1064
1065	! If this level is not already in the array, update
1066	! the counter (make sure we haven't exceeded expected
1067	! number of levels), then add to array.
1068	IF (myfindloc(expected_num_levels, &
1069	& millibars, VALUE=gribmsg_level) .EQ. 0) THEN
1070	level_count = level_count + 1
1071	IF (level_count .LE. expected_num_levels) THEN
1072	millibars(level_count) = gribmsg_level
1073	ELSE
1074	WRITE(6,"('ERROR: Found more than ', I3, ' levels')") &
1075	& expected_num_levels
1076	return_code = RETURN_CODE_GRIBFILE_ERROR
1077	RETURN
1078	ENDIF
1079
1080	ENDIF
1081	ENDIF
1082
1083
1084
1085
1086	ENDIF
1087	END IF
1088
1089
1090	END DO ! End of main loop for reading each GRIB message
1091
1092	! Make sure we found the number of isobaric levels that we expected
1093	IF (level_count .LT. expected_num_levels) THEN
1094	WRITE(6,"('ERROR: Only found ', I3, ' levels')") level_count
1095	return_code = RETURN_CODE_GRIBFILE_ERROR
1096	RETURN
1097	ENDIF
1098
1099	! OK, if we made it here, we should have a full millibars array
1100	! containing all of the isobaric levels we found. NOTE THAT I
1101	! AM NOT BOTHERING TO SORT THIS - I DON'T SEE A REASON FOR THAT
1102	! SO FAR.
1103
1104	RETURN
1105
1106
1107	END SUBROUTINE build_fv3_millibars_array
1108
1109
1110
1111
1112
1113
1114
1115	INTEGER FUNCTION myfindloc(n, array, value)
1116
1117
1118	INTEGER, INTENT(IN) :: n
1119	INTEGER, INTENT(IN), DIMENSION(n) :: array
1120	INTEGER, INTENT(IN) :: value
1121
1122	! My own implementation of FINDLOC, which is part of the
1123	! Fortran 2008 standard
1124
1125	! Returns index of first instance of "value" in the array. Returns
1126	! zero if not found. Assumes single-dimension INTEGER array
1127
1128	! Default value
1129
1130	INTEGER :: i
1131	LOGICAL :: found
1132
1133	found = .FALSE.
1134	myfindloc = 0
1135	DO i=1,SIZE(array)
1136	IF (array(i) .EQ. value) THEN
1137	myfindloc = i
1138	EXIT
1139	ENDIF
1140	END DO
1141
1142	RETURN
1143
1144
1145
1146	END FUNCTION myfindloc
1147
1148
1149
1150
1151
1152	SUBROUTINE write_err(shortname, level)
1153
1154	IMPLICIT NONE
1155
1156	CHARACTER(LEN=GRIB_ATTRIBUTE_STR_LEN), INTENT(IN) :: shortname
1157	INTEGER, INTENT(IN) :: level
1158
1159	WRITE(6,*) 'ERROR: Bad gribmsg_shortname: ', TRIM(shortname),&
1160	& ', gribmsg_level: ', level
1161
1162	END SUBROUTINE write_err
1163
1164
1165	INTEGER FUNCTION cgutils_check_ecmwf_grib(filepath, expected_num_levels, &
1166	& check_fcst, check_anl, check_grib2) &
1167	& RESULT(return_code)
1168
1169	USE grib_api
1170
1171	IMPLICIT NONE
1172
1173	CHARACTER(LEN=256), INTENT(IN) :: filepath
1174	INTEGER, INTENT(IN) :: expected_num_levels
1175	LOGICAL, INTENT(IN) :: check_fcst, check_anl, check_grib2
1176
1177	! return_code: possible values defined in module header, RETURN_CODE_*
1178
1179	INTEGER :: ifile, igrib, iret, grib_centre
1180
1181	CHARACTER(LEN=GRIB_ATTRIBUTE_STR_LEN) :: gribmsg_leveltype, &
1182	& gribmsg_shortname
1183	INTEGER :: gribmsg_level
1184
1185	LOGICAL :: end_of_file
1186	!!LOGICAL :: fcst_anl_good
1187
1188
1189	! Start with return_code 0. Anything wrong will change it
1190	return_code = RETURN_CODE_NORMAL
1191
1192	CALL grib_open_file(ifile, filepath, 'r', iret)
1193	IF (iret == 0) THEN
1194	! Use first record to detect centre, which is assumed constant
1195	! amongst all messages
1196	CALL grib_new_from_file(ifile, igrib, iret)
1197	CALL grib_get(igrib, 'centre', grib_centre)
1198	CALL grib_close_file(ifile)
1199	ELSE
1200	WRITE(0,*) 'ERROR: problem opening GRIB file: '//TRIM(filepath)
1201	return_code = RETURN_CODE_UTILITY_ERROR
1202	RETURN
1203	END IF
1204
1205
1206	IF (grib_centre .NE. GRIB_CENTRE_ECMWF) THEN
1207	WRITE(6,*) 'ERROR: Unexpected grib_centre: ', grib_centre
1208	return_code = RETURN_CODE_GRIBFILE_ERROR
1209	RETURN
1210	ENDIF
1211
1212
1213	CALL allocate_and_init_flags(expected_num_levels)
1214
1215	! Iterate through the entire file, keeping track of what we've found
1216	end_of_file = .FALSE.
1217	CALL grib_open_file(ifile, filepath, 'r', iret)
1218
1219	DO WHILE (.NOT. end_of_file)
1220	CALL grib_new_from_file(ifile, igrib, iret)
1221	IF (iret .eq. GRIB_END_OF_FILE) THEN
1222	end_of_file = .TRUE.
1223	ELSE
1224
1225
1226	CALL grib_get(igrib, 'shortName', gribmsg_shortname)
1227	CALL grib_get(igrib, 'typeOfLevel', gribmsg_leveltype)
1228	CALL grib_get(igrib, 'level', gribmsg_level)
1229
1230	! Check general attributes of this message. Note that
1231	! all messages will be checked, even those we might not
1232	! care about. As of 2019-01-17, the checks are for GRIB2
1233	! and for fcst/anl if specified by "check_anl" or "check_fcst"
1234	IF (.NOT. message_is_good(igrib, &
1235	& check_anl, check_fcst, check_grib2)) THEN
1236	CALL write_err(gribmsg_shortname, gribmsg_level)
1237	return_code = RETURN_CODE_GRIBFILE_ERROR
1238	RETURN
1239	ENDIF
1240
1241	! Filter for 3D variables of interest
1242	IF ( TRIM(gribmsg_leveltype) == 'hybrid' ) THEN
1243
1244
1245	CALL grib_get(igrib, 'level', gribmsg_level)
1246
1247	! Check that there aren't more levels than expected
1248	IF (gribmsg_level .GT. expected_num_levels) THEN
1249	WRITE(6,*) 'ERROR: gribmsg_shortname: ', &
1250	& TRIM(gribmsg_shortname), &
1251	& ', gribmsg_level: ', gribmsg_level, &
1252	& ', expected max levels: ', &
1253	& expected_num_levels
1254	return_code = RETURN_CODE_GRIBFILE_ERROR
1255	RETURN
1256	ENDIF
1257
1258
1259	SELECT CASE(gribmsg_shortname)
1260	CASE('u')
1261	global_levels_u(gribmsg_level) = .TRUE.
1262	CASE('v')
1263	global_levels_v(gribmsg_level) = .TRUE.
1264	CASE('t')
1265	global_levels_t(gribmsg_level) = .TRUE.
1266	CASE('q')
1267	global_levels_q(gribmsg_level) = .TRUE.
1268	CASE('etadot')
1269	global_levels_etadot(gribmsg_level) = .TRUE.
1270	END SELECT
1271
1272
1273	!!!!!! ! End of section for 3d variables of interest
1274
1275	! Filter for 2d variables of interest with level type
1276	! of "surface"
1277	ELSE IF ( TRIM(gribmsg_leveltype) == 'surface' ) THEN
1278
1279
1280	SELECT CASE(gribmsg_shortname)
1281	CASE('sp')
1282	global_2dvar_ps = .TRUE.
1283
1284	! For snow depth, I've been seeing both 'sd' and
1285	! 'sde' used in the CTBTO ECMWF files
1286	CASE('sd')
1287	global_2dvar_sd = .TRUE.
1288	CASE('sde')
1289	global_2dvar_sd = .TRUE.
1290
1291	CASE('tcc')
1292	global_2dvar_tcc = .TRUE.
1293	CASE('lsp')
1294	global_2dvar_lsprec = .TRUE.
1295	CASE('acpcp')
1296	global_2dvar_convprec = .TRUE.
1297	CASE('sshf')
1298	global_2dvar_shf = .TRUE.
1299	CASE('ssr')
1300	global_2dvar_sr = .TRUE.
1301	CASE('ewss')
1302	global_2dvar_ewss = .TRUE.
1303	CASE('nsss')
1304	global_2dvar_nsss = .TRUE.
1305	CASE('z')
1306	global_2dvar_oro = .TRUE.
1307	CASE('sdor')
1308	global_2dvar_excessoro = .TRUE.
1309	CASE('lsm')
1310	global_2dvar_lsm = .TRUE.
1311
1312
1313	! GRIB 1 variables found in the EFYYMMDDHH files
1314	! Note that these are also found in GRIB2 files
1315	! under different level types, except for cp, which
1316	! is acpcp in GRIB2. If there is ever a mixture
1317	! (which I think is unlikely) in a GRIB file, then
1318	! the world will likely end and this won't work
1319	! correctly.
1320	CASE('cp')
1321	global_2dvar_convprec = .TRUE.
1322	CASE('10u')
1323	global_2dvar_u10 = .TRUE.
1324	CASE('10v')
1325	global_2dvar_v10 = .TRUE.
1326	CASE('2t')
1327	global_2dvar_t2 = .TRUE.
1328	CASE('2d')
1329	global_2dvar_td2 = .TRUE.
1330	CASE('msl')
1331	global_2dvar_msl = .TRUE.
1332
1333
1334
1335	END SELECT
1336
1337
1338	!!!!!!!! End of section for 2d "surface" variables of interest
1339
1340
1341	ELSE IF ( TRIM(gribmsg_leveltype) == 'heightAboveGround' ) THEN
1342
1343	SELECT CASE(gribmsg_shortname)
1344	CASE('10u')
1345	global_2dvar_u10 = .TRUE.
1346	CASE('10v')
1347	global_2dvar_v10 = .TRUE.
1348	CASE('2t')
1349	global_2dvar_t2 = .TRUE.
1350	CASE('2d')
1351	global_2dvar_td2 = .TRUE.
1352	END SELECT
1353
1354
1355	!!!!!!!! End of section for 2d "heightAboveGround"
1356	!!!!!!!! variables of interest
1357
1358	ELSE IF ( TRIM(gribmsg_leveltype) == 'meanSea' ) THEN
1359
1360	SELECT CASE(gribmsg_shortname)
1361	CASE('msl')
1362	global_2dvar_msl = .TRUE.
1363	END SELECT
1364
1365
1366	!!!!!!!! End of section for 2d "meanSea"
1367	!!!!!!!! variables of interest
1368
1369
1370	END IF ! Selection of message types (3d, surface, etc.)
1371
1372
1373	END IF
1374
1375
1376
1377
1378	END DO ! End of main loop for reading each GRIB message
1379	CALL grib_close_file(ifile)
1380
1381	! This function will have the side effect of printing any
1382	! problems to stdout. If false, we force exit
1383	IF (ecmwf_all_levels_present()) THEN
1384	return_code = RETURN_CODE_NORMAL
1385	ELSE
1386	return_code = RETURN_CODE_GRIBFILE_ERROR
1387	CALL delete_levels_flags
1388	RETURN
1389	ENDIF
1390
1391	IF (ecmwf_all_2dvars_present()) THEN
1392	return_code = RETURN_CODE_NORMAL
1393	ELSE
1394	return_code = RETURN_CODE_GRIBFILE_ERROR
1395	RETURN
1396	ENDIF
1397
1398
1399	RETURN
1400
1401	END FUNCTION cgutils_check_ecmwf_grib
1402
1403
1404
1405	INTEGER FUNCTION cgutils_check_ncep_grib2(filepath, expected_num_levels, &
1406	& check_fcst, check_anl, check_grib2) &
1407	& RESULT(return_code)
1408
1409	USE grib_api
1410
1411	IMPLICIT NONE
1412
1413	CHARACTER(LEN=256), INTENT(IN) :: filepath
1414	INTEGER, INTENT(IN) :: expected_num_levels
1415	LOGICAL, INTENT(IN) :: check_fcst, check_anl, check_grib2
1416
1417	! return_code: possible values defined in module header, RETURN_CODE_*
1418
1419	INTEGER :: ifile, igrib, iret, grib_centre
1420
1421	CHARACTER(LEN=GRIB_ATTRIBUTE_STR_LEN) :: gribmsg_leveltype, &
1422	& gribmsg_shortname
1423	INTEGER :: gribmsg_level, plevel_idx
1424
1425	LOGICAL :: end_of_file
1426
1427	INTEGER, ALLOCATABLE, DIMENSION(:) :: millibars
1428
1429
1430	! Start with return_code 0. Anything wrong will change it
1431	return_code = RETURN_CODE_NORMAL
1432
1433	CALL grib_open_file(ifile, filepath, 'r', iret)
1434	IF (iret == 0) THEN
1435	! Use first record to detect centre, which is assumed constant
1436	! amongst all messages
1437	CALL grib_new_from_file(ifile, igrib, iret)
1438	CALL grib_get(igrib, 'centre', grib_centre)
1439	CALL grib_close_file(ifile)
1440	ELSE
1441	WRITE(0,*) 'ERROR: problem opening GRIB file: '//TRIM(filepath)
1442	return_code = RETURN_CODE_UTILITY_ERROR
1443	RETURN
1444	END IF
1445
1446
1447	IF (grib_centre .NE. GRIB_CENTRE_NCEP) THEN
1448	WRITE(6,*) 'ERROR: Unexpected grib_centre: ', grib_centre
1449	return_code = RETURN_CODE_GRIBFILE_ERROR
1450	RETURN
1451	ENDIF
1452
1453
1454
1455	! Pass through the file once in order to build an inventory of
1456	! isobaric pressure levels so that later we can make sure that
1457	! each 3D variable has each of these pressure levels. While
1458	! we're doing the inventory we can check that the number of
1459	! isobaric levels found is indeed what we expect.
1460	CALL grib_open_file(ifile, filepath, 'r', iret)
1461	IF (iret == 0) THEN
1462	ALLOCATE(millibars(expected_num_levels))
1463	CALL build_millibars_array(ifile, expected_num_levels, &
1464	& millibars, return_code)
1465	CALL grib_close_file(ifile)
1466	ELSE
1467	WRITE(0,*) 'ERROR: problem opening GRIB file: '//TRIM(filepath)
1468	return_code = RETURN_CODE_UTILITY_ERROR
1469	RETURN
1470	END IF
1471
1472	! Check that build_millibars_array() had no errors
1473	! We will assume that any error messages were written from the
1474	! subroutine
1475	IF (return_code .NE. RETURN_CODE_NORMAL) THEN
1476	return_code = RETURN_CODE_GRIBFILE_ERROR
1477	RETURN
1478	ENDIF
1479
1480	! At this point, millibars should be an array of all the pressure
1481	! levels (isobaricInhPa) that we found.
1482	!PRINT *, 'millibars: ', millibars
1483
1484	CALL allocate_and_init_flags(expected_num_levels)
1485
1486
1487	! Now, make our primary checking pass through the GRIB file.
1488
1489	! For 3D isobaric-level variables, as we read a level, we will
1490	! find its index in millibars array, and use that to index the
1491	! correct flag element in our logical array for each 3D variable
1492
1493	! Iterate through the entire file, keeping track of what we've found
1494	end_of_file = .FALSE.
1495	CALL grib_open_file(ifile, filepath, 'r', iret)
1496
1497	DO WHILE (.NOT. end_of_file)
1498	CALL grib_new_from_file(ifile, igrib, iret)
1499	IF (iret .eq. GRIB_END_OF_FILE) THEN
1500	end_of_file = .TRUE.
1501	ELSE
1502
1503
1504	CALL grib_get(igrib, 'shortName', gribmsg_shortname)
1505	CALL grib_get(igrib, 'typeOfLevel', gribmsg_leveltype)
1506	CALL grib_get(igrib, 'level', gribmsg_level)
1507
1508	! Check general attributes of this message. Note that
1509	! all messages will be checked, even those we might not
1510	! care about. As of 2019-01-17, the checks are for GRIB2
1511	! and for fcst/anl if specified by "check_anl" or "check_fcst"
1512	IF (.NOT. message_is_good(igrib, &
1513	& check_anl, check_fcst, check_grib2)) THEN
1514	CALL write_err(gribmsg_shortname, gribmsg_level)
1515	return_code = RETURN_CODE_GRIBFILE_ERROR
1516	RETURN
1517	ENDIF
1518
1519
1520
1521
1522
1523
1524
1525	!PRINT *, 'shortname, leveltype, level, plevel_idx: ', &
1526	!& gribmsg_shortname, gribmsg_leveltype, &
1527	!& gribmsg_level, plevel_idx
1528
1529
1530
1531
1532
1533
1534	! Filter for 3D variables of interest
1535	IF ( TRIM(gribmsg_leveltype) == 'isobaricInhPa' ) THEN
1536
1537
1538	CALL grib_get(igrib, 'level', gribmsg_level)
1539
1540	! Unsorted array "millibars," created above, tells us what index we
1541	! should use to access the retrieved level flag as we maintain an
1542	! inventory of what we found, and what we didn't
1543	! Use my function myfindloc() to find the index for a specified pressure
1544	! level.
1545
1546	plevel_idx = myfindloc(expected_num_levels, millibars, gribmsg_level)
1547	IF (plevel_idx .LE. 0) THEN
1548	CALL write_err(gribmsg_shortname, gribmsg_level)
1549	return_code = RETURN_CODE_GRIBFILE_ERROR
1550	RETURN
1551	ENDIF
1552
1553
1554	!PRINT *, 'shortname, leveltype, level, plevel_idx: ', &
1555	!& gribmsg_shortname, gribmsg_leveltype, &
1556	!& gribmsg_level, plevel_idx
1557
1558	SELECT CASE(gribmsg_shortname)
1559	CASE('u')
1560	global_levels_u(plevel_idx) = .TRUE.
1561	CASE('v')
1562	global_levels_v(plevel_idx) = .TRUE.
1563	CASE('t')
1564	global_levels_t(plevel_idx) = .TRUE.
1565	CASE('r')
1566	global_levels_rh(plevel_idx) = .TRUE.
1567	CASE('w')
1568	global_levels_w(plevel_idx) = .TRUE.
1569	END SELECT
1570
1571	!!!!!! ! End of section for 3d variables of interest
1572
1573
1574
1575
1576
1577	! Filter for 2d variables of interest with level type
1578	! of "surface"
1579	ELSE IF ( TRIM(gribmsg_leveltype) == 'surface' ) THEN
1580
1581
1582	SELECT CASE(gribmsg_shortname)
1583	CASE('sp')
1584	global_2dvar_ps = .TRUE.
1585	CASE('sdwe')
1586	global_2dvar_sd = .TRUE.
1587	CASE('tp')
1588	global_2dvar_lsprec = .TRUE.
1589	CASE('acpcp')
1590	global_2dvar_convprec = .TRUE.
1591	CASE('orog')
1592	global_2dvar_oro = .TRUE.
1593	CASE('lsm')
1594	global_2dvar_lsm = .TRUE.
1595	CASE('hpbl')
1596	global_2dvar_hmix = .TRUE.
1597	END SELECT
1598
1599
1600	!!!!!!!! End of section for 2d "surface" variables of interest
1601
1602
1603	ELSE IF ( TRIM(gribmsg_leveltype) == 'heightAboveGround' ) THEN
1604
1605	SELECT CASE(gribmsg_shortname)
1606	CASE('10u')
1607	global_2dvar_u10 = .TRUE.
1608	CASE('10v')
1609	global_2dvar_v10 = .TRUE.
1610	CASE('2t')
1611	global_2dvar_t2 = .TRUE.
1612	CASE('2d')
1613	global_2dvar_td2 = .TRUE.
1614
1615	! 2m RH is different depending on version of gribapi - in 1.14 it's just r, so
1616	! one needs to also make sure that level is 2, but in later gribapi it seems to be
1617	! 2r
1618	CASE('r')
1619	IF (gribmsg_level .EQ. 2) global_2dvar_rh2 = .TRUE.
1620	CASE('2r')
1621	global_2dvar_rh2 = .TRUE.
1622	END SELECT
1623
1624	!!!!!!!! End of section for 2d "heightAboveGround"
1625	!!!!!!!! variables of interest
1626
1627	ELSE IF ( TRIM(gribmsg_leveltype) == 'meanSea' ) THEN
1628
1629	SELECT CASE(gribmsg_shortname)
1630	CASE('prmsl')
1631	global_2dvar_msl = .TRUE.
1632	END SELECT
1633
1634
1635	!!!!!!!! End of section for 2d "meanSea"
1636	!!!!!!!! variables of interest
1637
1638	ELSE IF ( TRIM(gribmsg_leveltype) == 'Low cloud layer' ) THEN
1639
1640	SELECT CASE(gribmsg_shortname)
1641	CASE('tcc')
1642	global_2dvar_tcc = .TRUE.
1643	END SELECT
1644
1645
1646	!!!!!!!! End of section for 2d "Low cloud layer"
1647	!!!!!!!! variables of interest
1648
1649	ELSE IF ( TRIM(gribmsg_leveltype) == 'sigma' ) THEN
1650
1651	SELECT CASE(gribmsg_shortname)
1652	CASE('t')
1653	global_2dvar_tsig1 = .TRUE.
1654	CASE('u')
1655	global_2dvar_usig1 = .TRUE.
1656	CASE('v')
1657	global_2dvar_vsig1 = .TRUE.
1658	END SELECT
1659
1660
1661	!!!!!!!! End of section for 2d "sigma"
1662	!!!!!!!! variables of interest
1663
1664
1665
1666
1667
1668
1669
1670	END IF ! Selection of message types (3d, surface, etc.)
1671
1672
1673	END IF
1674
1675
1676
1677
1678	END DO ! End of main loop for reading each GRIB message
1679	CALL grib_close_file(ifile)
1680
1681
1682
1683
1684
1685	! Check that all the 3d levels we flagged in the above loop resulted
1686	! in all flags being true
1687
1688	! This function will have the side effect of printing any
1689	! problems to stdout. If false, we force exit
1690	IF (ncep_all_levels_present()) THEN
1691	return_code = RETURN_CODE_NORMAL
1692	ELSE
1693	return_code = RETURN_CODE_GRIBFILE_ERROR
1694	CALL delete_levels_flags
1695	RETURN
1696	ENDIF
1697
1698	! The w (vertical wind) variable has to be handled differently because it
1699	! is only available in the NCEP files for 100 mb and higher pressure levels
1700	IF (ncep_all_expected_w_levels_present(expected_num_levels, &
1701	& millibars)) THEN
1702	return_code = RETURN_CODE_NORMAL
1703	ELSE
1704	return_code = RETURN_CODE_GRIBFILE_ERROR
1705	WRITE(6,*) 'ERROR: Not all expected w levels are present'
1706	CALL delete_levels_flags
1707	RETURN
1708	ENDIF
1709
1710
1711
1712	IF (ncep_all_2dvars_present()) THEN
1713	return_code = RETURN_CODE_NORMAL
1714	ELSE
1715	return_code = RETURN_CODE_GRIBFILE_ERROR
1716	RETURN
1717	ENDIF
1718
1719
1720	! Be sure to deallocate our allocatable arrays
1721	DEALLOCATE(millibars)
1722	CALL delete_levels_flags
1723
1724
1725
1726
1727	RETURN
1728
1729	END FUNCTION cgutils_check_ncep_grib2
1730
1731
1732	INTEGER FUNCTION cgutils_check_fv3_grib2(filepath, expected_num_levels, &
1733	& check_fcst, check_anl, check_grib2) &
1734	& RESULT(return_code)
1735
1736	USE grib_api
1737
1738	IMPLICIT NONE
1739
1740	CHARACTER(LEN=256), INTENT(IN) :: filepath
1741	INTEGER, INTENT(IN) :: expected_num_levels
1742	LOGICAL, INTENT(IN) :: check_fcst, check_anl, check_grib2
1743
1744	! return_code: possible values defined in module header, RETURN_CODE_*
1745
1746	INTEGER :: ifile, igrib, iret, grib_centre
1747
1748	CHARACTER(LEN=GRIB_ATTRIBUTE_STR_LEN) :: gribmsg_leveltype, &
1749	& gribmsg_shortname
1750	INTEGER :: gribmsg_level, plevel_idx
1751
1752	LOGICAL :: end_of_file
1753
1754	INTEGER, ALLOCATABLE, DIMENSION(:) :: millibars
1755
1756
1757	! Start with return_code 0. Anything wrong will change it
1758	return_code = RETURN_CODE_NORMAL
1759
1760	CALL grib_open_file(ifile, filepath, 'r', iret)
1761	IF (iret == 0) THEN
1762	! Use first record to detect centre, which is assumed constant
1763	! amongst all messages
1764	CALL grib_new_from_file(ifile, igrib, iret)
1765	CALL grib_get(igrib, 'centre', grib_centre)
1766	CALL grib_close_file(ifile)
1767	ELSE
1768	WRITE(0,*) 'ERROR: problem opening GRIB file: '//TRIM(filepath)
1769	return_code = RETURN_CODE_UTILITY_ERROR
1770	RETURN
1771	END IF
1772
1773
1774	IF (grib_centre .NE. GRIB_CENTRE_NCEP) THEN
1775	WRITE(6,*) 'ERROR: Unexpected grib_centre: ', grib_centre
1776	return_code = RETURN_CODE_GRIBFILE_ERROR
1777	RETURN
1778	ENDIF
1779
1780
1781
1782	! Pass through the file once in order to build an inventory of
1783	! isobaric pressure levels so that later we can make sure that
1784	! each 3D variable has each of these pressure levels. While
1785	! we're doing the inventory we can check that the number of
1786	! isobaric levels found is indeed what we expect.
1787	CALL grib_open_file(ifile, filepath, 'r', iret)
1788	IF (iret == 0) THEN
1789	ALLOCATE(millibars(expected_num_levels))
1790	CALL build_fv3_millibars_array(ifile, expected_num_levels, &
1791	& millibars, return_code)
1792	CALL grib_close_file(ifile)
1793	ELSE
1794	WRITE(0,*) 'ERROR: problem opening GRIB file: '//TRIM(filepath)
1795	return_code = RETURN_CODE_UTILITY_ERROR
1796	RETURN
1797	END IF
1798
1799	! Check that build_millibars_array() had no errors
1800	! We will assume that any error messages were written from the
1801	! subroutine
1802	IF (return_code .NE. RETURN_CODE_NORMAL) THEN
1803	return_code = RETURN_CODE_GRIBFILE_ERROR
1804	RETURN
1805	ENDIF
1806
1807	! At this point, millibars should be an array of all the pressure
1808	! levels (isobaricInhPa) that we found.
1809	!PRINT *, 'millibars: ', millibars
1810
1811	CALL allocate_and_init_flags(expected_num_levels)
1812
1813
1814	! Now, make our primary checking pass through the GRIB file.
1815
1816	! For 3D isobaric-level variables, as we read a level, we will
1817	! find its index in millibars array, and use that to index the
1818	! correct flag element in our logical array for each 3D variable
1819
1820	! Iterate through the entire file, keeping track of what we've found
1821	end_of_file = .FALSE.
1822	CALL grib_open_file(ifile, filepath, 'r', iret)
1823
1824	DO WHILE (.NOT. end_of_file)
1825	CALL grib_new_from_file(ifile, igrib, iret)
1826	IF (iret .eq. GRIB_END_OF_FILE) THEN
1827	end_of_file = .TRUE.
1828	ELSE
1829
1830
1831	CALL grib_get(igrib, 'shortName', gribmsg_shortname)
1832	CALL grib_get(igrib, 'typeOfLevel', gribmsg_leveltype)
1833	CALL grib_get(igrib, 'level', gribmsg_level)
1834
1835	! Check general attributes of this message. Note that
1836	! all messages will be checked, even those we might not
1837	! care about. As of 2019-01-17, the checks are for GRIB2
1838	! and for fcst/anl if specified by "check_anl" or "check_fcst"
1839	IF (.NOT. message_is_good(igrib, &
1840	& check_anl, check_fcst, check_grib2)) THEN
1841	CALL write_err(gribmsg_shortname, gribmsg_level)
1842	return_code = RETURN_CODE_GRIBFILE_ERROR
1843	RETURN
1844	ENDIF
1845
1846
1847	!PRINT *, 'shortname, leveltype, level, plevel_idx: ', &
1848	!& gribmsg_shortname, gribmsg_leveltype, &
1849	!& gribmsg_level, plevel_idx
1850
1851
1852	! Filter for 3D variables of interest
1853	IF ( TRIM(gribmsg_leveltype) == 'isobaricInhPa' ) THEN
1854
1855
1856	CALL grib_get(igrib, 'level', gribmsg_level)
1857
1858	! Unsorted array "millibars," created above, tells us what index we
1859	! should use to access the retrieved level flag as we maintain an
1860	! inventory of what we found, and what we didn't
1861	! Use my function myfindloc() to find the index for a specified pressure
1862	! level.
1863
1864
1865
1866	! This is annoyingly "special case" code. Because of new FV3 files, "partial"
1867	! levels are present in the GRIB files. Other routines compute the
1868	! list of "valid" pressure levels, but as we go through and process
1869	! each message, as we're doing here, if it's on one of the invalid
1870	! pressure levels then we simply want to ignore it. If it's a valid
1871	! pressure level, and it's a variable we're interested in, we mark it
1872	! as present on that pressure level
1873
1874	plevel_idx = myfindloc(expected_num_levels, millibars, gribmsg_level)
1875	IF (plevel_idx .LE. 0) THEN
1876	!CALL write_err(gribmsg_shortname, gribmsg_level)
1877	!return_code = RETURN_CODE_GRIBFILE_ERROR
1878	!RETURN
1879	CONTINUE
1880	ELSE
1881
1882
1883	!PRINT *, 'shortname, leveltype, level, plevel_idx: ', &
1884	!& gribmsg_shortname, gribmsg_leveltype, &
1885	!& gribmsg_level, plevel_idx
1886
1887	SELECT CASE(gribmsg_shortname)
1888	CASE('u')
1889	global_levels_u(plevel_idx) = .TRUE.
1890	CASE('v')
1891	global_levels_v(plevel_idx) = .TRUE.
1892	CASE('t')
1893	global_levels_t(plevel_idx) = .TRUE.
1894	CASE('r')
1895	global_levels_rh(plevel_idx) = .TRUE.
1896	CASE('w')
1897	global_levels_w(plevel_idx) = .TRUE.
1898	END SELECT
1899	ENDIF
1900
1901	!!!!!! ! End of section for 3d variables of interest
1902
1903
1904
1905
1906
1907	! Filter for 2d variables of interest with level type
1908	! of "surface"
1909	ELSE IF ( TRIM(gribmsg_leveltype) == 'surface' ) THEN
1910
1911
1912	SELECT CASE(gribmsg_shortname)
1913	CASE('sp')
1914	global_2dvar_ps = .TRUE.
1915	CASE('sdwe')
1916	global_2dvar_sd = .TRUE.
1917	CASE('tp')
1918	global_2dvar_lsprec = .TRUE.
1919	CASE('acpcp')
1920	global_2dvar_convprec = .TRUE.
1921	CASE('orog')
1922	global_2dvar_oro = .TRUE.
1923	CASE('lsm')
1924	global_2dvar_lsm = .TRUE.
1925	CASE('hpbl')
1926	global_2dvar_hmix = .TRUE.
1927	END SELECT
1928
1929
1930	!!!!!!!! End of section for 2d "surface" variables of interest
1931
1932
1933	ELSE IF ( TRIM(gribmsg_leveltype) == 'heightAboveGround' ) THEN
1934
1935	SELECT CASE(gribmsg_shortname)
1936	CASE('10u')
1937	global_2dvar_u10 = .TRUE.
1938	CASE('10v')
1939	global_2dvar_v10 = .TRUE.
1940	CASE('2t')
1941	global_2dvar_t2 = .TRUE.
1942	CASE('2d')
1943	global_2dvar_td2 = .TRUE.
1944
1945	! 2m RH is different depending on version of gribapi - in 1.14 it's just r, so
1946	! one needs to also make sure that level is 2, but in later gribapi it seems to be
1947	! 2r
1948	CASE('r')
1949	IF (gribmsg_level .EQ. 2) global_2dvar_rh2 = .TRUE.
1950	CASE('2r')
1951	global_2dvar_rh2 = .TRUE.
1952	END SELECT
1953
1954	!!!!!!!! End of section for 2d "heightAboveGround"
1955	!!!!!!!! variables of interest
1956
1957	ELSE IF ( TRIM(gribmsg_leveltype) == 'meanSea' ) THEN
1958
1959	SELECT CASE(gribmsg_shortname)
1960	CASE('prmsl')
1961	global_2dvar_msl = .TRUE.
1962	END SELECT
1963
1964
1965	!!!!!!!! End of section for 2d "meanSea"
1966	!!!!!!!! variables of interest
1967
1968	ELSE IF ( TRIM(gribmsg_leveltype) == 'Low cloud layer' ) THEN
1969
1970	SELECT CASE(gribmsg_shortname)
1971	CASE('tcc')
1972	global_2dvar_tcc = .TRUE.
1973	END SELECT
1974
1975
1976	!!!!!!!! End of section for 2d "Low cloud layer"
1977	!!!!!!!! variables of interest
1978
1979	ELSE IF ( TRIM(gribmsg_leveltype) == 'sigma' ) THEN
1980
1981	SELECT CASE(gribmsg_shortname)
1982	CASE('t')
1983	global_2dvar_tsig1 = .TRUE.
1984	CASE('u')
1985	global_2dvar_usig1 = .TRUE.
1986	CASE('v')
1987	global_2dvar_vsig1 = .TRUE.
1988	END SELECT
1989
1990
1991	!!!!!!!! End of section for 2d "sigma"
1992	!!!!!!!! variables of interest
1993
1994
1995	END IF ! Selection of message types (3d, surface, etc.)
1996
1997
1998	END IF
1999
2000
2001
2002
2003	END DO ! End of main loop for reading each GRIB message
2004	CALL grib_close_file(ifile)
2005
2006
2007
2008
2009
2010	! Check that all the 3d levels we flagged in the above loop resulted
2011	! in all flags being true
2012
2013	! This function will have the side effect of printing any
2014	! problems to stdout. If false, we force exit
2015	IF (ncep_all_levels_present()) THEN
2016	return_code = RETURN_CODE_NORMAL
2017	ELSE
2018	return_code = RETURN_CODE_GRIBFILE_ERROR
2019	CALL delete_levels_flags
2020	RETURN
2021	ENDIF
2022
2023	! The w (vertical wind) variable has to be handled differently because it
2024	! is only available in the NCEP files for 100 mb and higher pressure levels
2025	IF (fv3_all_expected_w_levels_present(expected_num_levels, &
2026	& millibars)) THEN
2027	return_code = RETURN_CODE_NORMAL
2028	ELSE
2029	return_code = RETURN_CODE_GRIBFILE_ERROR
2030	WRITE(6,*) 'ERROR: Not all expected w levels are present'
2031	CALL delete_levels_flags
2032	RETURN
2033	ENDIF
2034
2035
2036
2037	IF (fv3_all_2dvars_present()) THEN
2038	return_code = RETURN_CODE_NORMAL
2039	ELSE
2040	return_code = RETURN_CODE_GRIBFILE_ERROR
2041	RETURN
2042	ENDIF
2043
2044
2045	! Be sure to deallocate our allocatable arrays
2046	DEALLOCATE(millibars)
2047	CALL delete_levels_flags
2048
2049
2050
2051
2052	RETURN
2053
2054	END FUNCTION cgutils_check_fv3_grib2
2055
2056
2057
2058
2059
2060	END MODULE cgutils
2061

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

Download in other formats:

Original Format