forked from hschwaiger-usgs/volcano-ash3d-metreader
-
Notifications
You must be signed in to change notification settings - Fork 0
/
MetReader.F90
4638 lines (4039 loc) · 214 KB
/
MetReader.F90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
module MetReader
integer, parameter,private :: sp = 4 ! single precision
integer, parameter,private :: dp = 8 ! double precision
integer, parameter :: MR_MAXVARS = 50 ! Maximum number of variables in fixed arrays
real(kind=dp), parameter :: MR_EPS_SMALL = 1.0e-7_dp ! Small number
real(kind=dp), parameter :: MR_EPS_TINY = 1.0e-12_dp ! Very small number
real(kind=sp), parameter :: RAD_EARTH_MET = 6371.229_sp ! Radius of Earth in km
real(kind=sp), parameter :: DEG2RAD_MET = 1.7453292519943295e-2_sp
real(kind=sp), parameter :: MR_MIN_DZ = 1.0e-4_sp ! used in reassiging z for low-level negative GPH
integer,public :: MR_global_essential = 6
integer,public :: MR_global_production = 6
integer,public :: MR_global_debug = 6
integer,public :: MR_global_info = 6
integer,public :: MR_global_log = 9
integer,public :: MR_global_error = 0
logical :: MR_useCompGrid = .true. ! Reset this to .false. if you only need the Met grid
logical :: MR_useCompTime = .true. ! Reset this to .false. if you only need the time of the file
integer,public :: MR_iwind ! MR_IWIND specifies the type of wind input to the model:
! MR_IWIND=1 if a 1-D wind sounding is use,
! =2 if a 3-D grid is read from a ASCII file.
! =3 if a single, multistep 3-D file is used
! =4 if multiple 3-D NetCDF files are used
! =5 if multiple file with multiple steps are used
integer,public :: MR_iwindformat ! MR_iwindformat specifies the format of the met data
! 0 Custom format based on template
! 1 ASCII profile
! 2 Radiosonde data
! 3 NARR3D 32 km North America files (32 km) :: ds608.0
! 4 NAM Regional North America 221 (32 km)
! 5 NAM AK 216 (45 km)
! 6 NAM Regional 104 (90 km)
! 7 NAM CONUS 212 (40 km)
! 8 NAM CONUS 218 (12 km)
! 9 NAM CONUS 227 (5.08 km)
! 10 NAM AK 242 (11.25 km)
! 11 NAM HI 196 (2.5 km)
! 12 NAM AK 198 (5.953 km)
! 13 NAM AK 91 (2.976 km)
! 14 NAM CONUS 1227 (3.0 km)
! 20 GFS 0.5
! 21 GFS 1.0
! 22 GFS 0.25
! 23 NCEP / DOE reanalysis 2.5 degree files :: ds091.0
! 24 NASA-MERRA-2 reanalysis 0.625/0.5 degree files
! 25 NCEP/NCAR reanalysis 2.5 degree files :: ds090.0 iwind=4 or 5
! 26 JRA-55 :: ds628.0 iwind=5
! 27 NOAA-CIRES reanalysis 2.0 degree files :: ds131.2 iwind=5
! 28 ECMWF Interim Reanalysis (ERA-Interim) :: ds627.0 requires catted grib files
! 29 ECMWF ERA5 :: ds630.0 iwind=5
! 30 ECMWF 20-Century (ERA-20C) :: ds626.0 iwind=5
! 32 Air Force Weather Agency subcenter = 0
! 33 CCSM3.0 Community Atmosphere Model (CAM)
! 40 NASA-GEOS Cp
! 41 NASA-GEOS Np
! 50 WRF - output
logical,public :: MR_Use_RDA = .true.
integer,public :: MR_iGridCode ! MR_iGridCode specifies the NCEP grid described in:
! http://www.nco.ncep.noaa.gov/pmb/docs/on388/tableb.html
integer,public :: MR_idataFormat ! Specifies the data model used
! =1 ASCII
! =2 netcdf
! =3 grib1 or grib2
! These variables describe the full list of windfiles read
integer ,public :: MR_iwindfiles ! number of files provided
#ifdef USEPOINTERS
! character(len=130), allocatable,dimension(:) ,public :: fc_windfilename
character(len=130), pointer ,dimension(:) ,public :: MR_windfiles ! name of file
#else
character(len=130), allocatable,dimension(:) ,public :: MR_windfiles ! name of file
#endif
character(len=130) ,public :: MR_iw5_root
character(len=42) ,public :: MR_iw5_prefix
character(len=24) ,public :: MR_iw5_suffix1 ! e.g. 1912060100_1912063021.nc
character(len=24) ,public :: MR_iw5_suffix2
real(kind=dp) , allocatable,dimension(:) ,public :: MR_windfile_starthour ! start hour of the file
integer , allocatable,dimension(:) ,public :: MR_windfiles_nt_fullmet ! number of steps in files
real(kind=dp) , allocatable,dimension(:,:),public :: MR_windfile_stephour ! offset hours of step
character(len=80) ,public :: MR_iwf_template ! name of the template file
logical , allocatable,dimension(:) ,public :: MR_windfiles_Have_GRIB_index
character(len=130), allocatable,dimension(:) ,public :: MR_windfiles_GRIB_index ! name of grib index file
! These variables are a list of the same data above, but specific to the simulation
! duration
integer ,public :: MR_MetSteps_Total
integer ,public :: MR_iMetStep_Now
character(len=130), allocatable,dimension(:),public :: MR_MetStep_File
integer , allocatable,dimension(:),public :: MR_MetStep_findex
integer , allocatable,dimension(:),public :: MR_MetStep_tindex
real(kind=dp) , allocatable,dimension(:),public :: MR_MetStep_Hour_since_baseyear
real(kind=dp) , allocatable,dimension(:),public :: MR_MetStep_Interval
integer , allocatable,dimension(:),public :: MR_MetStep_year
integer , allocatable,dimension(:),public :: MR_MetStep_month
integer , allocatable,dimension(:),public :: MR_MetStep_day
integer , allocatable,dimension(:),public :: MR_MetStep_DOY
real(kind=dp) , allocatable,dimension(:),public :: MR_MetStep_Hour_Of_Day
logical :: MR_Reannalysis = .false.
integer ,dimension(:), allocatable,public :: MR_iwind5_year
! Native grid of Met file using Pressure as vertical coordinate
#ifdef USEPOINTERS
integer ,dimension(:,:) ,pointer, public :: MR_dum2d_met_int => null() ! Used for categorical variables
real(kind=sp),dimension(:,:) ,pointer, public :: MR_dum2d_met => null() ! Used for surface variables
real(kind=sp),dimension(:,:,:),pointer, public :: MR_dum3d_metP => null()
real(kind=sp),dimension(:,:,:),pointer, public :: MR_dum3d2_metP => null()
real(kind=sp),dimension(:,:,:),pointer, public :: MR_geoH_metP_last => null() ! These are needed for compH interpolation
real(kind=sp),dimension(:,:,:),pointer, public :: MR_geoH_metP_next => null()
real(kind=sp),dimension(:,:,:),pointer, public :: MR_vx_metP_last => null() !These might need to be stored to avoid a
real(kind=sp),dimension(:,:,:),pointer, public :: MR_vx_metP_next => null() !second reading
real(kind=sp),dimension(:,:,:),pointer, public :: MR_vy_metP_last => null()
real(kind=sp),dimension(:,:,:),pointer, public :: MR_vy_metP_next => null()
#else
integer ,dimension(:,:) ,allocatable, public :: MR_dum2d_met_int ! Used for categorical variables
real(kind=sp),dimension(:,:) ,allocatable, public :: MR_dum2d_met ! Used for surface variables
real(kind=sp),dimension(:,:,:),allocatable, public :: MR_dum3d_metP
real(kind=sp),dimension(:,:,:),allocatable, public :: MR_dum3d2_metP
real(kind=sp),dimension(:,:,:),allocatable, public :: MR_geoH_metP_last ! These are needed for compH interpolation
real(kind=sp),dimension(:,:,:),allocatable, public :: MR_geoH_metP_next
real(kind=sp),dimension(:,:,:),allocatable, public :: MR_vx_metP_last !These might need to be stored to avoid a
real(kind=sp),dimension(:,:,:),allocatable, public :: MR_vx_metP_next !second reading
real(kind=sp),dimension(:,:,:),allocatable, public :: MR_vy_metP_last
real(kind=sp),dimension(:,:,:),allocatable, public :: MR_vy_metP_next
#endif
logical, public :: MR_Save_Velocities = .false.
real(kind=sp), public :: MR_Max_geoH_metP_predicted ! Highest expected height in Met file, based on pressure
real(kind=sp) :: Max_geoH_metP_last
real(kind=sp) :: Max_geoH_metP_next
real(kind=sp) :: Max_geoH_metP
real(kind=sp) :: Suppl_H
!Parameter iHeightHandler specifies what to do if the maximum height
!of the simulation region exceeds the maximum height in the wind files.
!If iHeightHandler = 1, stop the program if the plume height exceeds mesoscale height
! 2, wind direction at levels above the highest node
! equal that of the highest node. Temperatures in the
! upper nodes don't change between 11 and 20 km; above
! 20 km they increase by 2 C/km, as in the Standard
! atmosphere. A warning is written to the log file.
integer :: MR_iHeightHandler = 2
! The following are used for sonde data
integer :: MR_nSnd_Locs = 1 ! Number of Sonde locations
integer :: MR_Snd_nt_fullmet = 1 ! Number of times at the Sonde locations
integer :: MR_Snd_nvars = 5 ! Number of Sonde variables (P,H,U,V,T,+user-specified)
logical :: Snd_Have_PT = .false.
logical :: Snd_Have_Coord = .false. ! If the 1-d data have the optional projection params
! then it will be used, otherwise vel will be relative
! to comp grid.
real(kind=sp),dimension(:,:,:,:),allocatable, public :: MR_SndVars_metP ! (MR_nSnd_Locs,MR_Snd_nt_fullmet,MR_Snd_nvars,300)
integer ,dimension(:,:) ,allocatable, public :: MR_Snd_np_fullmet ! Number of pressure values for each location/time
integer ,dimension(:) ,allocatable, public :: MR_SndVarsID ! Lists which vars are in which columns of MR_SndVars_metP
real(kind=sp),dimension(:,:,:) ,allocatable, public :: MR_Snd2Comp_tri_map_wgt ! weights of nearby sondes for every comp point
integer ,dimension(:,:,:) ,allocatable, public :: MR_Snd2Comp_tri_map_idx ! sonde index of weights
! Native grid of Met file using Height as vertical coordinate
! (resampled onto z-gridpoints of computational grid)
#ifdef USEPOINTERS
real(kind=sp),dimension(:,:,:),pointer, public :: MR_dum3d_metH => null()
#else
real(kind=sp),dimension(:,:,:),allocatable, public :: MR_dum3d_metH
#endif
real(kind=sp),dimension(:,:,:),allocatable, public :: MR_u_ER_metP ! For the cases where Met is proj and comp
real(kind=sp),dimension(:,:,:),allocatable, public :: MR_v_ER_metP ! different we need to rotate so these
! store Earth-Relative velocities on MetP
! Computations grid
#ifdef USEPOINTERS
integer ,dimension(:,:) ,pointer, public :: MR_dum2d_comp_int => null() ! Used for categorical variables
real(kind=sp),dimension(:,:) ,pointer, public :: MR_dum2d_comp => null() ! Used for surface variables
real(kind=sp),dimension(:,:,:),pointer, public :: MR_dum3d_compH => null()
real(kind=sp),dimension(:,:,:),pointer, public :: MR_dum3d_compH_2 => null() ! Used only when a vector field
! rotation is needed
#else
integer ,dimension(:,:) ,allocatable, public :: MR_dum2d_comp_int ! Used for categorical variables
real(kind=sp),dimension(:,:) ,allocatable, public :: MR_dum2d_comp ! Used for surface variables
real(kind=sp),dimension(:,:,:),allocatable, public :: MR_dum3d_compH
real(kind=sp),dimension(:,:,:),allocatable, public :: MR_dum3d_compH_2 ! Used only when a vector field
! rotation is needed
#endif
integer :: nt_fullmet ! length of t of file
integer :: nx_fullmet ! length of x or lon of full met grid
integer :: ny_fullmet ! length of y or lat of full met grid
integer :: np_fullmet ! length of pressure of full met grid for H,U,V
integer :: np_fullmet_pad = 1 ! We might need to pad the top of the pressure grid.
integer :: neta_fullmet ! Only used by WRF
real(kind=sp) :: Pressure_Conv_Fac = 100.0_sp ! factor for converting from Met file units to Pa
#ifdef USEPOINTERS
real(kind=sp),dimension(:), pointer :: x_fullmet_sp => null() ! x-coordinates of full met grid
real(kind=sp),dimension(:), pointer :: y_fullmet_sp => null() ! y-coordinates of full met grid
real(kind=sp),dimension(:), pointer :: p_fullmet_sp => null() ! z-coordinates of full met grid for H
real(kind=sp),dimension(:,:), pointer :: levs_fullmet_sp => null() ! This hold each of the numbered level coordinates:
! i.e. isobaric, isobaric1, isobaric2, but also
! height_above_ground, depth_below_surface_layer, etc.
! p_fullmet_sp,p_fullmet_[Vz,RH]sp are copies of one
! of the slices
integer,dimension(:), pointer :: nlevs_fullmet ! length of z coordinate
integer,dimension(:), pointer :: levs_code ! code indicating how to map to the GPH grid
! 0 = no mapping (not a pressure coordinate)
! 1 = one-to-one mapping (U,V)
! 2 = upper truncation (missing upper levels)
! 3 = interpolation (missing mid-levels)
! 4 = more levels than GPH grid
#else
real(kind=sp),dimension(:), allocatable :: x_fullmet_sp ! x-coordinates of full met grid
real(kind=sp),dimension(:), allocatable :: y_fullmet_sp ! y-coordinates of full met grid
real(kind=sp),dimension(:), allocatable :: p_fullmet_sp ! z-coordinates of full met grid for H
real(kind=sp),dimension(:,:),allocatable :: levs_fullmet_sp ! This hold each of the numbered level coordinates:
! i.e. isobaric, isobaric1, isobaric2, but also
! height_above_ground, depth_below_surface_layer, etc.
! p_fullmet_sp,p_fullmet_[Vz,RH]sp are copies of one
! of the slices
integer,dimension(:), allocatable :: nlevs_fullmet
integer,dimension(:), allocatable :: levs_code ! code indicating how to map to the GPH grid
! 0 = no mapping (not a pressure coordinate)
! 1 = one-to-one mapping (U,V)
! 2 = upper truncation (missing upper levels)
! 3 = interpolation (missing mid-levels)
! 4 = more levels than GPH grid
#endif
logical :: IsLatLon_MetGrid
logical :: IsGlobal_MetGrid = .false. ! Not all Lon/Lat grids are periodic
logical :: IsLatLon_CompGrid
logical :: IsPeriodic_CompGrid = .false.
logical :: UseFullMetGrid = .false. ! This is the special case where the comp grid
! equals the Met grid
integer :: nx_submet ! length of x or lon of sub-grid
integer :: ny_submet ! length of y or lat of sub-grid
#ifdef USEPOINTERS
real(kind=sp),dimension(:), pointer:: x_submet_sp => null() ! x-coordinates of met sub-grid
real(kind=sp),dimension(:), pointer:: y_submet_sp => null() ! y-coordinates of met sub-grid
real(kind=sp),dimension(:), pointer:: z_approx => null() ! zpproximate altidue from STD Atmos and press (in km)
#else
real(kind=sp),dimension(:), allocatable :: x_submet_sp ! x-coordinates of met sub-grid
real(kind=sp),dimension(:), allocatable :: y_submet_sp ! y-coordinates of met sub-grid
real(kind=sp),dimension(:), allocatable :: z_approx ! zpproximate altidue from STD Atmos and press (in km)
#endif
logical :: x_inverted = .false.
logical :: y_inverted = .false. ! Some LatLon grids start at the North Pole and increment down
logical :: z_inverted = .false. ! Some grids give top pressure first
logical :: y_pad_North = .false. ! Some computational grids will require values above the top met point
logical :: y_pad_South = .false. !
logical :: isGridRelative = .true. ! Most windfiles, whether Lat/Lon or projected, give
! velocities relative to the grid of the file. Some (NARR)
! give velocities relative to earth coordinates and need to
! be rotated
real(kind=dp),dimension(:,:) ,allocatable :: theta_Met ! Earth to grid
real(kind=dp),dimension(:,:) ,allocatable :: theta_Comp ! Earth to grid
! Met copies of projection variables, used for proj call on Met Grid
integer :: Met_iprojflag
real(kind=8) :: Met_Re
real(kind=8) :: Met_k0
real(kind=8) :: Met_phi0 != 90.0_ip ! latitude of projection point
real(kind=8) :: Met_lam0 != -135.0_ip ! longitude of projection point
real(kind=8) :: Met_lam1,Met_phi1
real(kind=8) :: Met_lam2,Met_phi2
integer :: Comp_iprojflag
real(kind=8) :: Comp_Re
real(kind=8) :: Comp_k0
real(kind=8) :: Comp_phi0 != 90.0_ip ! latitude of projection point
real(kind=8) :: Comp_lam0 != -135.0_ip ! longitude of projection point
real(kind=8) :: Comp_lam1,Comp_phi1
real(kind=8) :: Comp_lam2,Comp_phi2
integer :: Map_Case
! Some geometry terms
#ifdef USEPOINTERS
real(kind=sp),dimension(:,:) ,pointer :: rdphi_MetP_sp => null()
real(kind=sp),dimension(:,:,:) ,pointer :: rdlambda_MetP_sp => null()
real(kind=sp),dimension(:) ,pointer :: MR_dx_met, MR_dx_submet => null()
real(kind=sp),dimension(:) ,pointer :: MR_dy_met, MR_dy_submet => null()
#else
real(kind=sp),dimension(:,:) ,allocatable,private :: rdphi_MetP_sp
real(kind=sp),dimension(:,:,:) ,allocatable,private :: rdlambda_MetP_sp
real(kind=sp),dimension(:) ,allocatable :: MR_dx_met, MR_dx_submet
real(kind=sp),dimension(:) ,allocatable :: MR_dy_met, MR_dy_submet
#endif
logical :: IsRegular_MetGrid ! True if the grid-spacing is uniform
real(kind=sp) :: dx_met_const, dy_met_const
! There are some computational grid variables we might need, so make local copies
integer :: MR_BaseYear = 1900 ! This should be reset in calling program
logical :: MR_useLeap = .true. ! This too
integer :: MR_Comp_StartYear
integer :: MR_Comp_StartMonth
integer :: MR_Comp_StartDay
real(kind=dp) :: MR_Comp_StartHour = 0.0_dp ! Note that these must be double-precision to
real(kind=dp) :: MR_Comp_Time_in_hours = 0.0_dp ! be passed correctly to HoursSince
integer :: nx_comp
integer :: ny_comp
integer :: nz_comp
real(kind=sp) :: dx_comp,dy_comp
real(kind=sp) :: MaxZ_comp_sp
#ifdef USEPOINTERS
real(kind=sp),dimension(:), pointer :: x_comp_sp => null() ! x-coordinates of computational grid
real(kind=sp),dimension(:), pointer :: y_comp_sp => null() ! y-coordinates of computational grid
real(kind=sp),dimension(:), pointer :: z_comp_sp => null() ! z-coordinates of computational grid
real(kind=sp),dimension(:,:), pointer :: CompPoint_X_on_Met_sp => null() ! x-coord (on Met grid) of comp point
real(kind=sp),dimension(:,:), pointer :: CompPoint_Y_on_Met_sp => null() ! y-coord (on Met grid) of comp point
integer ,dimension(:,:,:),pointer :: CompPoint_on_subMet_idx => null() ! index on met sub-grid of comp point
real(kind=sp),dimension(:,:,:),pointer :: bilin_map_wgt => null()
#else
real(kind=sp),dimension(:), allocatable :: x_comp_sp ! x-coordinates of computational grid
real(kind=sp),dimension(:), allocatable :: y_comp_sp ! y-coordinates of computational grid
real(kind=sp),dimension(:), allocatable :: z_comp_sp ! z-coordinates of computational grid
real(kind=sp),dimension(:,:), allocatable :: CompPoint_X_on_Met_sp ! x-coord (on Met grid) of comp point
real(kind=sp),dimension(:,:), allocatable :: CompPoint_Y_on_Met_sp ! y-coord (on Met grid) of comp point
integer ,dimension(:,:,:),allocatable :: CompPoint_on_subMet_idx ! index on met sub-grid of comp point
real(kind=sp),dimension(:,:,:),allocatable :: bilin_map_wgt
#endif
! Here are a few variables needed for sigma-altitude coordinates
logical :: MR_use_SigmaAlt = .false.
integer :: MR_ZScaling_ID = 0 ! = 0 for no scaling (i.e. s = z)
! = 1 for altitude shifting (s=z=zsurf)
! = 2 for sigma-altitude (s=(z-surf)/(top-surf))
real(kind=sp) :: MR_ztop
#ifdef USEPOINTERS
real(kind=sp),dimension(:), pointer :: s_comp_sp => null() ! s-coordinates (scaled z) of computational grid
real(kind=sp),dimension(:,:), pointer :: MR_zsurf => null() ! surface elevation in km
real(kind=sp),dimension(:,:), pointer :: MR_jacob => null() ! Jacobian of trans. = MR_ztop-MR_zsurf
#else
real(kind=sp),dimension(:), allocatable :: s_comp_sp ! s-coordinates (scaled z) of computational grid
real(kind=sp),dimension(:,:), allocatable :: MR_zsurf ! surface elevation in km
real(kind=sp),dimension(:,:), allocatable :: MR_jacob ! Jacobian of trans. = MR_ztop-MR_zsurf
#endif
logical :: FoundFillVAttr = .false.
real(kind=sp) :: fill_value_sp = -9999.0_sp ! Initialize to the commonly used fill value
character(len=30),dimension(9) :: Met_dim_names ! name of dimension
logical ,dimension(9) :: Met_dim_IsAvailable
real(kind=sp) ,dimension(9) :: Met_dim_fac = 1.0_sp
! Here is the list of variables that can be read. Each iwindformat will
! have just a sub-set availible with specific names. For now, allocate
! space for 50 variable names
! Mechanical / State variables
! 1 = Geopotential Height
! 2 = Vx
! 3 = Vy
! 4 = Vz
! 5 = Temperature
! Surface
! 10 = Planetary Boundary Layer Height
! 11 = U @ 10m
! 12 = V @ 10m
! 13 = Friction velocity
! 15 = Snow cover
! 16 = Soil moisture
! 17 = Surface Roughness
! 18 = Wind gust speed
! 19 = surface temperature
! Atmospheric Structure
! 20 = pressure at lower cloud base
! 21 = pressure at lower cloud top
! 22 = temperature at lower cloud top
! 23 = Total Cloud cover
! 24 = Cloud cover (low)
! 25 = Cloud cover (convective)
! Moisture
! 30 = Rel. Hum
! 31 = QV (specific humidity)
! 32 = QL (liquid)
! 33 = QI (ice)
! Precipitation
! 40 = Categorical rain
! 41 = Categorical snow
! 42 = Categorical frozen rain
! 43 = Categorical ice
! 44 = Precipitation rate large-scale (liquid)
! 45 = Precipitation rate convective (liquid)
! 46 = Precipitation rate large-scale (ice)
! 47 = Precipitation rate convective (ice)
logical ,dimension(MR_MAXVARS) :: Met_var_IsAvailable ! true if iwf contains the var
character(len=71),dimension(MR_MAXVARS) :: Met_var_NC_names ! name in the file
character(len=71),dimension(MR_MAXVARS) :: Met_var_GRIB_names ! name in the file
character(len=5) ,dimension(MR_MAXVARS) :: Met_var_WMO_names ! WMO version of the name
integer ,dimension(MR_MAXVARS) :: Met_var_ndim !
integer ,dimension(MR_MAXVARS) :: Met_var_zdim_idx ! The index of this coordinate (used in Met_var_nlevs)
integer ,dimension(MR_MAXVARS) :: Met_var_zdim_ncid ! The dimID of the dimension in the nc file
integer :: nlev_coords_detected = 0
integer ,dimension(MR_MAXVARS,4) :: Met_var_GRIB2_DPcPnSt ! Grib2 files have variables identified by
! discpln,param_cat,param_num,surf_class
integer ,dimension(MR_MAXVARS) :: Met_var_GRIB1_Param ! indicatorOfParameter
integer ,dimension(MR_MAXVARS) :: Met_var_GRIB1_Table ! table2Version
character(len=3) ,dimension(MR_MAXVARS) :: Met_var_GRIB1_St ! level type (pl, src, 116 etc)
integer :: MR_GRIB_Version = 0
!logical ,dimension(MR_MAXVARS) :: Met_var_IsFloat ! true if kind=4 otherwise false
real(kind=sp) ,dimension(MR_MAXVARS) :: Met_var_conversion_factor
integer ,dimension(MR_MAXVARS) :: Met_var_nlevs
! Variables needed by netcdf reader
real(kind=sp) :: iwf25_scale_facs(MR_MAXVARS)
real(kind=sp) :: iwf25_offsets(MR_MAXVARS)
real(kind=sp) :: x_in_iwf25_sp(192)
real(kind=sp) :: y_in_iwf25_sp(94)
! Here is the mapping for bilinear weighting coeffiecients (amap) and
! indices (imap) from the 1.875-deg 2d grid to the 2.5-deg
#ifdef USEPOINTERS
real(kind=sp),dimension(:,:,:) ,pointer :: amap_iwf25 => null()
integer ,dimension(:,:,:) ,pointer :: imap_iwf25 => null()
#else
real(kind=sp),dimension(:,:,:) ,allocatable :: amap_iwf25
integer ,dimension(:,:,:) ,allocatable :: imap_iwf25
#endif
integer(kind=sp),dimension(:,:,:),allocatable :: tmpsurf2d_short
integer :: istart, iend
integer :: jstart, jend
integer :: ilhalf_fm_l, ilhalf_fm_r
integer :: irhalf_fm_l, irhalf_fm_r
integer :: ilhalf_nx, irhalf_nx
logical :: wrapgrid
real(kind=sp) ,dimension(:) ,allocatable :: temp1d_sp
real(kind=sp) ,dimension(:,:,:) ,allocatable :: temp2d_sp
real(kind=sp) ,dimension(:,:,:,:) ,allocatable :: temp3d_sp
integer(kind=sp),dimension(:,:,:) ,allocatable :: temp2d_int
integer(kind=sp),dimension(:,:,:) ,allocatable :: temp2d_short
integer(kind=sp),dimension(:,:,:,:) ,allocatable :: temp3d_short
real(kind=4),dimension(:,:),allocatable :: Met_Proj_lat
real(kind=4),dimension(:,:),allocatable :: Met_Proj_lon
! Status variables for error-checking
logical :: Check_prereq_conditions = .true.
logical :: CALLED_MR_Allocate_FullMetFileList = .false.
logical :: CALLED_MR_Read_Met_DimVars = .false.
logical :: CALLED_MR_Set_CompProjection = .false.
logical :: CALLED_MR_Initialize_Met_Grids = .false.
logical :: CALLED_MR_Set_Met_Times = .false.
contains
!##############################################################################
!
! MR_Reset_Memory
!
! This subroutine reinitializes MetReader by deallocating all MR variables.
! This is useful if a program needs to use multiple types of wind files.
!
!##############################################################################
subroutine MR_Reset_Memory
implicit none
write(MR_global_production,*)"-------------------------------------------------------"
write(MR_global_production,*)"-------- Resetting all MetReader Memory ---------------"
write(MR_global_production,*)"-------------------------------------------------------"
if(allocated(MR_windfile_starthour ))deallocate(MR_windfile_starthour)
if(allocated(MR_windfiles_nt_fullmet ))deallocate(MR_windfiles_nt_fullmet)
if(allocated(MR_windfile_stephour ))deallocate(MR_windfile_stephour)
if(allocated(MR_windfiles_Have_GRIB_index ))deallocate(MR_windfiles_Have_GRIB_index)
if(allocated(MR_windfiles_GRIB_index ))deallocate(MR_windfiles_GRIB_index)
if(allocated(MR_MetStep_File ))deallocate(MR_MetStep_File)
if(allocated(MR_MetStep_findex ))deallocate(MR_MetStep_findex)
if(allocated(MR_MetStep_tindex ))deallocate(MR_MetStep_tindex)
if(allocated(MR_MetStep_Hour_since_baseyear))deallocate(MR_MetStep_Hour_since_baseyear)
if(allocated(MR_MetStep_Interval ))deallocate(MR_MetStep_Interval)
if(allocated(MR_MetStep_year ))deallocate(MR_MetStep_year)
if(allocated(MR_MetStep_month ))deallocate(MR_MetStep_month)
if(allocated(MR_MetStep_day ))deallocate(MR_MetStep_day)
if(allocated(MR_MetStep_DOY ))deallocate(MR_MetStep_DOY)
if(allocated(MR_MetStep_Hour_Of_Day ))deallocate(MR_MetStep_Hour_Of_Day)
if(allocated(MR_iwind5_year ))deallocate(MR_iwind5_year)
#ifdef USEPOINTERS
if(associated(MR_windfiles ))deallocate(MR_windfiles)
if(associated(MR_dum2d_met_int ))deallocate(MR_dum2d_met_int)
if(associated(MR_dum2d_met ))deallocate(MR_dum2d_met)
if(associated(MR_dum3d_metP ))deallocate(MR_dum3d_metP)
if(associated(MR_dum3d2_metP ))deallocate(MR_dum3d2_metP)
if(associated(MR_geoH_metP_last ))deallocate(MR_geoH_metP_last)
if(associated(MR_geoH_metP_next ))deallocate(MR_geoH_metP_next)
if(associated(MR_vx_metP_last ))deallocate(MR_vx_metP_last)
if(associated(MR_vx_metP_next ))deallocate(MR_vx_metP_next)
if(associated(MR_vy_metP_last ))deallocate(MR_vy_metP_last)
if(associated(MR_vy_metP_next ))deallocate(MR_vy_metP_next)
if(associated(MR_dum3d_metH ))deallocate(MR_dum3d_metH)
if(associated(MR_dum2d_comp_int ))deallocate(MR_dum2d_comp_int)
if(associated(MR_dum2d_comp ))deallocate(MR_dum2d_comp)
if(associated(MR_dum3d_compH ))deallocate(MR_dum3d_compH)
if(associated(MR_dum3d_compH_2 ))deallocate(MR_dum3d_compH_2)
if(associated(x_fullmet_sp ))deallocate(x_fullmet_sp)
if(associated(y_fullmet_sp ))deallocate(y_fullmet_sp)
if(associated(p_fullmet_sp ))deallocate(p_fullmet_sp)
if(associated(x_submet_sp ))deallocate(x_submet_sp)
if(associated(y_submet_sp ))deallocate(y_submet_sp)
if(associated(levs_fullmet_sp ))deallocate(levs_fullmet_sp)
if(associated(nlevs_fullmet ))deallocate(nlevs_fullmet)
if(associated(levs_code ))deallocate(levs_code)
if(associated(z_approx ))deallocate(z_approx)
if(associated(rdphi_MetP_sp ))deallocate(rdphi_MetP_sp)
if(associated(rdlambda_MetP_sp ))deallocate(rdlambda_MetP_sp)
if(associated(MR_dx_met ))deallocate(MR_dx_met)
if(associated(MR_dx_submet ))deallocate(MR_dx_submet)
if(associated(MR_dy_met ))deallocate(MR_dy_met)
if(associated(MR_dy_submet ))deallocate(MR_dy_submet)
if(associated(x_comp_sp ))deallocate(x_comp_sp)
if(associated(y_comp_sp ))deallocate(y_comp_sp)
if(associated(z_comp_sp ))deallocate(z_comp_sp)
if(associated(CompPoint_X_on_Met_sp ))deallocate(CompPoint_X_on_Met_sp)
if(associated(CompPoint_Y_on_Met_sp ))deallocate(CompPoint_Y_on_Met_sp)
if(associated(CompPoint_on_subMet_idx ))deallocate(CompPoint_on_subMet_idx)
if(associated(bilin_map_wgt ))deallocate(bilin_map_wgt)
if(associated(amap_iwf25 ))deallocate(amap_iwf25)
if(associated(imap_iwf25 ))deallocate(imap_iwf25)
if(associated(s_comp_sp ))deallocate(s_comp_sp)
if(associated(MR_zsurf ))deallocate(MR_zsurf)
if(associated(MR_jacob ))deallocate(MR_jacob)
#else
if(allocated(MR_windfiles ))deallocate(MR_windfiles)
if(allocated(MR_dum2d_met_int ))deallocate(MR_dum2d_met_int)
if(allocated(MR_dum2d_met ))deallocate(MR_dum2d_met)
if(allocated(MR_dum3d_metP ))deallocate(MR_dum3d_metP)
if(allocated(MR_dum3d2_metP ))deallocate(MR_dum3d2_metP)
if(allocated(MR_geoH_metP_last ))deallocate(MR_geoH_metP_last)
if(allocated(MR_geoH_metP_next ))deallocate(MR_geoH_metP_next)
if(allocated(MR_vx_metP_last ))deallocate(MR_vx_metP_last)
if(allocated(MR_vx_metP_next ))deallocate(MR_vx_metP_next)
if(allocated(MR_vy_metP_last ))deallocate(MR_vy_metP_last)
if(allocated(MR_vy_metP_next ))deallocate(MR_vy_metP_next)
if(allocated(MR_dum3d_metH ))deallocate(MR_dum3d_metH)
if(allocated(MR_dum2d_comp_int ))deallocate(MR_dum2d_comp_int)
if(allocated(MR_dum2d_comp ))deallocate(MR_dum2d_comp)
if(allocated(MR_dum3d_compH ))deallocate(MR_dum3d_compH)
if(allocated(MR_dum3d_compH_2 ))deallocate(MR_dum3d_compH_2)
if(allocated(x_fullmet_sp ))deallocate(x_fullmet_sp)
if(allocated(y_fullmet_sp ))deallocate(y_fullmet_sp)
if(allocated(p_fullmet_sp ))deallocate(p_fullmet_sp)
if(allocated(x_submet_sp ))deallocate(x_submet_sp)
if(allocated(y_submet_sp ))deallocate(y_submet_sp)
if(allocated(levs_fullmet_sp ))deallocate(levs_fullmet_sp)
if(allocated(nlevs_fullmet ))deallocate(nlevs_fullmet)
if(allocated(levs_code ))deallocate(levs_code)
if(allocated(z_approx ))deallocate(z_approx)
if(allocated(rdphi_MetP_sp ))deallocate(rdphi_MetP_sp)
if(allocated(rdlambda_MetP_sp ))deallocate(rdlambda_MetP_sp)
if(allocated(MR_dx_met ))deallocate(MR_dx_met)
if(allocated(MR_dx_submet ))deallocate(MR_dx_submet)
if(allocated(MR_dy_met ))deallocate(MR_dy_met)
if(allocated(MR_dy_submet ))deallocate(MR_dy_submet)
if(allocated(x_comp_sp ))deallocate(x_comp_sp)
if(allocated(y_comp_sp ))deallocate(y_comp_sp)
if(allocated(z_comp_sp ))deallocate(z_comp_sp)
if(allocated(CompPoint_X_on_Met_sp ))deallocate(CompPoint_X_on_Met_sp)
if(allocated(CompPoint_Y_on_Met_sp ))deallocate(CompPoint_Y_on_Met_sp)
if(allocated(CompPoint_on_subMet_idx ))deallocate(CompPoint_on_subMet_idx)
if(allocated(bilin_map_wgt ))deallocate(bilin_map_wgt)
if(allocated(amap_iwf25 ))deallocate(amap_iwf25)
if(allocated(imap_iwf25 ))deallocate(imap_iwf25)
if(allocated(s_comp_sp ))deallocate(s_comp_sp)
if(allocated(MR_zsurf ))deallocate(MR_zsurf)
if(allocated(MR_jacob ))deallocate(MR_jacob)
#endif
if(allocated(MR_u_ER_metP ))deallocate(MR_u_ER_metP)
if(allocated(MR_v_ER_metP ))deallocate(MR_v_ER_metP)
if(allocated(theta_Met ))deallocate(theta_Met)
if(allocated(theta_Comp ))deallocate(theta_Comp)
if(allocated(tmpsurf2d_short ))deallocate(tmpsurf2d_short)
if(allocated(temp1d_sp ))deallocate(temp1d_sp)
if(allocated(temp2d_sp ))deallocate(temp2d_sp)
if(allocated(temp3d_sp ))deallocate(temp3d_sp)
if(allocated(temp2d_int ))deallocate(temp2d_int)
if(allocated(temp2d_short ))deallocate(temp2d_short)
if(allocated(temp3d_short ))deallocate(temp3d_short)
if(allocated(Met_Proj_lat ))deallocate(Met_Proj_lat)
if(allocated(Met_Proj_lon ))deallocate(Met_Proj_lon)
nlev_coords_detected = 0
end subroutine MR_Reset_Memory
!##############################################################################
!
! MR_Allocate_FullMetFileList
!
! This subroutine allocates the list of windfiles and does some
! error-checking based on iwind and iwindfiles.
!
! From the calling program, this is called once the information about
! the NWP files is available (e.g. the iwind, iwindformat, grid ID, data
! format, and number of windfiles.
!
! Variables allocated:
! MR_windfiles(MR_iwindfiles)
! MR_windfiles_nt_fullmet(MR_iwindfiles)
! and possibly:
! MR_windfiles_Have_GRIB_index(MR_iwindfiles)
! MR_windfiles_GRIB_index(MR_iwindfiles)
!
! The next step is for the calling program to populate MR_windfiles
!
!##############################################################################
subroutine MR_Allocate_FullMetFileList(iw,iwf,igrid,idf,iwfiles)
implicit none
integer,intent(in) :: iw
integer,intent(in) :: iwf
integer,intent(in) :: igrid
integer,intent(in) :: idf
integer,intent(in) :: iwfiles
integer :: i
MR_iwind = iw
MR_iwindformat = iwf
MR_iGridCode = igrid
MR_idataFormat = idf
MR_iwindfiles = iwfiles
write(MR_global_production,*)"--------------------------------------------------------------------------------"
write(MR_global_production,*)"---------- MR_Allocate_FullMetFileList ----------"
write(MR_global_production,*)"--------------------------------------------------------------------------------"
if ((MR_iwind.ne.1).and.(MR_iwind.ne.2).and. &
(MR_iwind.ne.3).and.(MR_iwind.ne.4).and. &
(MR_iwind.ne.5)) then
write(MR_global_error,*)'MR_iwind must be between 1 and 5. Program stopped'
write(MR_global_error,*)' MR_iwind = ',MR_iwind
write(MR_global_error,*)' MR_IWIND OPTIONS:'
write(MR_global_error,*)' MR_iwind = 1 read from a 1-D wind sounding'
write(MR_global_error,*)' 2 read from 3D gridded ASCII files'
write(MR_global_error,*)' 3 read from a single, multistep file'
write(MR_global_error,*)' 4 read from multiple files'
write(MR_global_error,*)' 5 read variables from separate files'
stop 1
endif
! Initialize the dimension and variable arrays. Select slots in these arrays will be
! overwritten from the calls in the case block below
Met_dim_IsAvailable=.false.
Met_var_IsAvailable=.false.
Met_var_IsAvailable(1:MR_MAXVARS) = .false.
Met_var_zdim_idx(1:MR_MAXVARS) = 0
Met_var_zdim_ncid(1:MR_MAXVARS) = 0
Met_var_GRIB2_DPcPnSt(1:MR_MAXVARS,1:4) = 0
Met_var_GRIB1_Param(1:MR_MAXVARS) = 0
Met_var_GRIB1_St(1:MR_MAXVARS) = ""
Met_var_conversion_factor(1:MR_MAXVARS) = 1.0_sp
Met_var_nlevs(1:MR_MAXVARS) = 0
isGridRelative = .true.
!--------------------------------
! Mechanical / State variables
!--------------------------------
! Note: All default names are that assigned by netcdf-java in the grib-to-nc conversion
! Command used :: java -Xmx2048m -classpath ~/ncj/netcdfAll-4.5.jar ucar.nc2.dataset.NetcdfDataset \
! -in ${GribFile} -out ${NCFile} -IsLargeFile
!
! Alternatively, one could use ncl_convert2nc
! Geopotential Height (m^2/s^2)
Met_var_NC_names(1) = "Geopotential_height_isobaric"
Met_var_GRIB_names(1) = "gh"
Met_var_WMO_names(1) = "HGT"
Met_var_GRIB2_DPcPnSt(1,1:4) = (/0, 3, 5, 100/)
Met_var_GRIB1_Param(1) = 7
Met_var_GRIB1_St(1) = "pl"
Met_var_ndim(1) = 4
! Velocity component in x (or E) direction (m/s)
Met_var_NC_names(2) = "u-component_of_wind_isobaric"
Met_var_GRIB_names(2) = "u"
Met_var_WMO_names(2) = "UGRD"
Met_var_GRIB2_DPcPnSt(2,1:4) = (/0, 2, 2, 100/)
Met_var_GRIB1_Param(2) = 33
Met_var_GRIB1_St(2) = "pl"
Met_var_ndim(2) = 4
! Velocity component in y (or N) direction (m/s)
Met_var_NC_names(3) = "v-component_of_wind_isobaric"
Met_var_GRIB_names(3) = "v"
Met_var_WMO_names(3) = "VGRD"
Met_var_GRIB2_DPcPnSt(3,1:4) = (/0, 2, 3, 100/)
Met_var_GRIB1_Param(3) = 34
Met_var_GRIB1_St(3) = "pl"
Met_var_ndim(3) = 4
! Velocity component in z direction (Pa/s)
Met_var_NC_names(4) = "Vertical_velocity_pressure_isobaric"
Met_var_GRIB_names(4) = "w"
Met_var_WMO_names(4) = "VVEL"
Met_var_GRIB2_DPcPnSt(4,1:4) = (/0, 2, 8, 100/)
Met_var_GRIB1_Param(4) = 39
Met_var_GRIB1_St(4) = "pl"
Met_var_ndim(4) = 4
! Temperature (K)
Met_var_NC_names(5) = "Temperature_isobaric"
Met_var_GRIB_names(5) = "t"
Met_var_WMO_names(5) = "TMP"
Met_var_GRIB2_DPcPnSt(5,1:4) = (/0, 0, 0, 100/)
Met_var_GRIB1_Param(5) = 11
Met_var_GRIB1_St(5) = "pl"
Met_var_ndim(5) = 4
!--------------------------------
! Surface
!--------------------------------
! Height of planetary boundary layer (m)
Met_var_NC_names(10) = "Planetary_Boundary_Layer_Height_surface"
Met_var_GRIB_names(10) = "hpbl"
Met_var_WMO_names(10) = "HPBL"
Met_var_GRIB2_DPcPnSt(10,1:4) = (/0, 3, 196, 1/)
Met_var_GRIB1_Param(10) = 221
Met_var_GRIB1_St(10) = "1"
Met_var_ndim(10) = 3
! Velocity component in x (or E) direction at 10 m above ground surface (m/s)
Met_var_NC_names(11) = "u-component_of_wind_height_above_ground"
Met_var_GRIB_names(11) = "u"
Met_var_WMO_names(11) = "UGRD"
Met_var_GRIB2_DPcPnSt(11,1:4) = (/0, 2, 2, 103/)
Met_var_GRIB1_Param(11) = 33
Met_var_GRIB1_St(11) = "105"
Met_var_ndim(11) = 4
! Velocity component in y (or N) direction at 10 m above ground surface (m/s)
Met_var_NC_names(12) = "v-component_of_wind_height_above_ground"
Met_var_GRIB_names(12) = "v"
Met_var_WMO_names(12) = "VGRD"
Met_var_GRIB2_DPcPnSt(12,1:4) = (/0, 2, 3, 103/)
Met_var_GRIB1_Param(12) = 34
Met_var_GRIB1_St(12) = "105"
Met_var_ndim(12) = 4
! Friction velocity (m/s)
Met_var_NC_names(13) = "Frictional_Velocity_surface"
Met_var_GRIB_names(13) = "fricv"
Met_var_WMO_names(13) = "FRICV"
Met_var_GRIB2_DPcPnSt(13,1:4) = (/0, 2, 197, 1/)
Met_var_GRIB1_Param(13) = 253
Met_var_GRIB1_St(13) = "1"
Met_var_ndim(13) = 3
! Snow depth (m)
Met_var_NC_names(15) = "Snow_depth_surface"
Met_var_GRIB_names(15) = "sd"
Met_var_WMO_names(15) = "SNOD"
Met_var_GRIB2_DPcPnSt(15,1:4) = (/0, 1, 11, 1/)
Met_var_GRIB1_Param(15) = 66
Met_var_GRIB1_St(15) = "1"
Met_var_ndim(15) = 3
! Soil Moisture (fraction)
Met_var_NC_names(16) = "Volumetric_Soil_Moisture_Content_depth_below_surface_layer"
Met_var_GRIB_names(16) = "soilw"
Met_var_WMO_names(16) = "SOILW"
Met_var_GRIB2_DPcPnSt(16,1:4) = (/2, 0, 192, 106/)
Met_var_GRIB1_Param(16) = 0
Met_var_GRIB1_St(16) = ""
Met_var_ndim(16) = 4
! Surface roughness (m)
Met_var_NC_names(17) = "Surface_roughness_surface"
Met_var_GRIB_names(17) = "sr"
Met_var_WMO_names(17) = "SFCR"
Met_var_GRIB2_DPcPnSt(17,1:4) = (/2, 0, 1, 1/)
Met_var_GRIB1_Param(17) = 0
Met_var_GRIB1_St(17) = ""
Met_var_ndim(17) = 3
! Wind gust speed (m/s)
Met_var_NC_names(18) = "Wind_speed_gust_surface"
Met_var_GRIB_names(18) = "gust"
Met_var_WMO_names(18) = "GUST"
Met_var_GRIB2_DPcPnSt(18,1:4) = (/0, 2, 22, 1/)
Met_var_GRIB1_Param(18) = 0
Met_var_GRIB1_St(18) = ""
Met_var_ndim(18) = 3
! Surface temperature (K)
Met_var_NC_names(19) = "Temperature_surface"
Met_var_WMO_names(19) = ""
Met_var_GRIB2_DPcPnSt(19,1:4) = (/0, 0, 0, 0/)
Met_var_GRIB1_Param(19) = 11
Met_var_GRIB1_St(19) = "1"
Met_var_ndim(19) = 3
!--------------------------------
! Atmospheric Structure
!--------------------------------
! Pressure at base of lower cloud level (Pa)
Met_var_NC_names(20) = "Pressure_cloud_base"
Met_var_GRIB_names(20) = "pres"
Met_var_WMO_names(20) = "PRES"
Met_var_GRIB2_DPcPnSt(20,1:4) = (/0, 3, 0, 2/)
Met_var_GRIB1_Param(20) = 1
Met_var_GRIB1_St(20) = "2"
Met_var_ndim(20) = 3
! Pressure at top of lower cloud level (Pa)
Met_var_NC_names(21) = "Pressure_cloud_tops"
Met_var_GRIB_names(21) = "pres"
Met_var_WMO_names(21) = "PRES"
Met_var_GRIB2_DPcPnSt(21,1:4) = (/0, 3, 0, 3/)
Met_var_GRIB1_Param(21) = 1
Met_var_GRIB1_St(21) = "3"
Met_var_ndim(21) = 3
! Temperature at the top of lower cloud level (K)
Met_var_NC_names(22) = "Temperature_cloud_tops"
Met_var_GRIB_names(22) = "t"
Met_var_WMO_names(22) = "TMP"
Met_var_GRIB2_DPcPnSt(22,1:4) = (/0, 0, 0, 3/)
Met_var_GRIB1_Param(22) = 11
Met_var_GRIB1_St(22) = "3"
Met_var_ndim(22) = 3
! Total cloud cover (%)
Met_var_NC_names(23) = "Total_cloud_cover_entire_atmosphere"
Met_var_GRIB_names(23) = "tcc"
Met_var_WMO_names(23) = "TCDC"
Met_var_GRIB2_DPcPnSt(23,1:4) = (/0, 6, 1, 200/)
Met_var_GRIB1_Param(23) = 71
Met_var_GRIB1_St(23) = "200"
Met_var_ndim(23) = 3
! Cloud cover of lower cloud level (%)
Met_var_NC_names(24) = "Low_cloud_cover_low_cloud"
Met_var_GRIB_names(24) = "lcc"
Met_var_WMO_names(24) = "LCDC"
Met_var_GRIB2_DPcPnSt(24,1:4) = (/0, 6, 3, 214/)
Met_var_GRIB1_Param(24) = 0
Met_var_GRIB1_St(24) = ""
Met_var_ndim(24) = 3
!--------------------------------
! Moisture
!--------------------------------
! Relative humidity (%)
Met_var_NC_names(30) = "Relative_humidity_isobaric"
Met_var_GRIB_names(30) = "r"
Met_var_WMO_names(30) = "RH"
Met_var_GRIB2_DPcPnSt(30,1:4) = (/0, 1, 1, 100/)
Met_var_GRIB1_Param(30) = 52
Met_var_GRIB1_St(30) = "pl"
Met_var_ndim(30) = 4
! Specific humidity (kg/kg)
Met_var_NC_names(31) = "Specific_humidity_isobaric"
Met_var_GRIB_names(31) = "q"
Met_var_WMO_names(31) = "SPFH"
Met_var_GRIB2_DPcPnSt(31,1:4) = (/0, 1, 0, 100/)
Met_var_GRIB1_Param(31) = 51
Met_var_GRIB1_St(31) = "pl"
Met_var_ndim(31) = 4
! Cloud water mixing ratio (kg/kg)
Met_var_NC_names(32) = "Cloud_mixing_ratio_isobaric"
Met_var_GRIB_names(32) = "clwmr"
Met_var_WMO_names(32) = "CLWMR"
Met_var_GRIB2_DPcPnSt(32,1:4) = (/0, 1, 22, 100/)
Met_var_ndim(32) = 4
! Snow mixing ratio (kg/kg)
Met_var_NC_names(33) = "Snow_mixing_ratio_isobaric"
Met_var_GRIB_names(33) = "snmr"
Met_var_WMO_names(33) = "SNMR"
Met_var_GRIB2_DPcPnSt(33,1:4) = (/0, 1, 25, 100/)
Met_var_GRIB1_Param(33) = 0
Met_var_GRIB1_St(33) = ""
Met_var_ndim(33) = 4
!--------------------------------
! Precipitation
!--------------------------------
! Categorical rain at ground surface (0/1 no/yes)
Met_var_NC_names(40) = "Categorical_Rain_surface"
Met_var_GRIB_names(40) = "crain"
Met_var_WMO_names(40) = "CRAIN"
Met_var_GRIB2_DPcPnSt(40,1:4) = (/0, 1, 192, 1/)
Met_var_GRIB1_Param(40) = 140
Met_var_GRIB1_St(40) = "1"
Met_var_ndim(40) = 3
! Categorical snow at ground surface (0/1 no/yes)
Met_var_NC_names(41) = "Categorical_Snow_surface"
Met_var_GRIB_names(41) = "csnow"
Met_var_WMO_names(41) = "CSNOW"
Met_var_GRIB2_DPcPnSt(41,1:4) = (/0, 1, 195, 1/)
Met_var_GRIB1_Param(41) = 143
Met_var_GRIB1_St(41) = "1"
Met_var_ndim(41) = 3
! Categorical freezing rain at ground surface (0/1 no/yes)
Met_var_NC_names(42) = "Categorical_Freezing_Rain_surface"
Met_var_GRIB_names(42) = "cfrzr"
Met_var_WMO_names(42) = "CFRZR"
Met_var_GRIB2_DPcPnSt(42,1:4) = (/0, 1, 193, 1/)
Met_var_GRIB1_Param(42) = 141
Met_var_GRIB1_St(42) = "1"
Met_var_ndim(42) = 3
! Categorical ice pellets at ground surface (0/1 no/yes)
Met_var_NC_names(43) = "Categorical_Ice_Pellets_surface"
Met_var_GRIB_names(43) = "cicep"
Met_var_WMO_names(43) = "CICEP"
Met_var_GRIB2_DPcPnSt(43,1:4) = (/0, 1, 194, 1/)
Met_var_GRIB1_Param(43) = 142
Met_var_GRIB1_St(43) = "1"
Met_var_ndim(43) = 3
! Precipitation rate at surface (kg/m2/s)
Met_var_NC_names(44) = "Precipitation_rate_surface"
Met_var_GRIB_names(44) = "prate"
Met_var_WMO_names(44) = "PRATE"
Met_var_GRIB2_DPcPnSt(44,1:4) = (/0, 1, 7, 1/)
Met_var_GRIB1_Param(44) = 59
Met_var_GRIB1_St(44) = "pl"
Met_var_ndim(44) = 3
! Convective liquid precipitation rate at surface (kg/m2/s)
Met_var_NC_names(45) = "Precip.rate convective (liquid)"
Met_var_WMO_names(45) = "CPRAT"
Met_var_GRIB1_Param(45) = 0
Met_var_GRIB1_St(45) = ""
Met_var_ndim(45) = 3
! Large-scale precipitation rate at surface (kg/m2/s)
Met_var_NC_names(46) = "Precip.rate large-scale (ice)"
Met_var_WMO_names(46) = ""
Met_var_GRIB1_Param(46) = 0
Met_var_GRIB1_St(46) = ""
Met_var_ndim(46) = 3
! Convective frozen precipitation rate at surface for (kg/m2/s)
Met_var_NC_names(47) = "Precip.rate convective (ice)"
Met_var_WMO_names(47) = ""
Met_var_GRIB1_Param(47) = 0
Met_var_GRIB1_St(47) = ""
Met_var_ndim(47) = 3
if (MR_iwindformat.eq.0) then
! Custom format based on template
MR_Reannalysis = .false.
write(MR_global_info,*)" NWP format to be used = ",MR_iwindformat,&
"Custom format based on template"
! This expects that MR_iwf_template has been filled by the calling program
call MR_Read_Met_Template
elseif (MR_iwindformat.eq.1) then
! ASCII profile
MR_Reannalysis = .false.
write(MR_global_info,*)" NWP format to be used = ",MR_iwindformat,&
"ASCII profile"
elseif (MR_iwindformat.eq.2) then
! Radiosonde data
MR_Reannalysis = .false.
write(MR_global_info,*)" NWP format to be used = ",MR_iwindformat,&
"Radiosonde data"
elseif (MR_iwindformat.eq.3) then
! NARR3D NAM221 32 km North America files
! https://rda.ucar.edu/datasets/ds608.0/
! merged_AWIP32.2018062000(.nc)
! Note that winds are "earth-relative" and must be rotated!
! See http://www.emc.ncep.noaa.gov/mmb/rreanl/faq.html#eta-winds
write(MR_global_info,*)" NWP format to be used = ",MR_iwindformat,&
"NARR3D NAM221 32 km North America files with Re=6367.47"
MR_iGridCode = 1221 ! This is almost NAM221, but uses a diff Re
call MR_Set_Met_NCEPGeoGrid(MR_iGridCode)
MR_Reannalysis = .true.
Met_dim_IsAvailable(1)=.true.; Met_dim_names(1) = "reftime"
Met_dim_IsAvailable(2)=.true.; Met_dim_names(2) = "isobaric2"
Met_dim_IsAvailable(3)=.true.; Met_dim_names(3) = "y"
Met_dim_IsAvailable(4)=.true.; Met_dim_names(4) = "x"
! Mechanical / State variables
Met_var_IsAvailable(1)=.true. ! Geopotential Height
Met_var_IsAvailable(2)=.true.; Met_var_NC_names(2)="u_wind_isobaric"
Met_var_IsAvailable(3)=.true.; Met_var_NC_names(3)="v_wind_isobaric"
Met_var_IsAvailable(4)=.true.; Met_var_NC_names(4)="Pressure_vertical_velocity_isobaric"
Met_var_IsAvailable(5)=.true.; Met_var_NC_names(5)="Temp_isobaric"
! Surface
Met_var_IsAvailable(10)=.true.; Met_var_NC_names(10)="Planetary_boundary_layer_height_surface"
Met_var_IsAvailable(11)=.true.; Met_var_NC_names(11)="u_wind_height_above_ground"
Met_var_IsAvailable(12)=.true.; Met_var_NC_names(12)="v_wind_height_above_ground"
Met_var_IsAvailable(13)=.true.; Met_var_NC_names(13)="Surface_friction_velocity_surface"
Met_var_IsAvailable(15)=.true.
! Atmospheric Structure
Met_var_IsAvailable(20)=.true.
Met_var_IsAvailable(21)=.true.
! Moisture
Met_var_IsAvailable(31)=.true.
Met_var_IsAvailable(32)=.true.; Met_var_NC_names(32)="Cloud_water_isobaric"
Met_var_IsAvailable(33)=.true.; Met_var_NC_names(33)="Ice_mixing_ratio_isobaric"
! Precipitation
Met_var_IsAvailable(40)=.true.; Met_var_NC_names(40)="Categorical_rain_yes1_no0_surface"
Met_var_IsAvailable(41)=.true.; Met_var_NC_names(41)="Categorical_snow_yes1_no0_surface"
Met_var_IsAvailable(42)=.true.; Met_var_NC_names(42)="Categorical_freezing_rain_yes1_no0_surface"
Met_var_IsAvailable(43)=.true.; Met_var_NC_names(43)="Categorical_ice_pellets_yes1_no0_surface"
Met_var_IsAvailable(44)=.true.
fill_value_sp = -9999.0_sp
elseif (MR_iwindformat.eq.4) then
! NAM Regional North America 221 32 km North America files
! http://motherlode.ucar.edu/native/conduit/data/nccf/com/nam/prod/
! nam.t00z.awip3200.tm00.grib2
write(MR_global_info,*)" NWP format to be used = ",MR_iwindformat,&
"NAM221 32 km North America files with Re=6371.229"
MR_iGridCode = 221
call MR_Set_Met_NCEPGeoGrid(MR_iGridCode)