86 integer isec,hrs_obs_cutoff,min_obs_cutoff
87 integer sec_intvl,stat_miss_val,time_inc_betwn_succ_fld
88 integer perturb_num,num_ens_fcst,prob_num,tot_num_prob
89 character*80 type_of_time_inc,stat_unit_time_key_succ
90 logical*1,
allocatable :: bmap(:)
92 integer,
allocatable :: mg(:)
94 integer,
parameter :: max_bytes=1000*1300000
95 integer,
parameter :: max_numbit=16
96 integer,
parameter :: lugi=650
97 character*255 fl_nametbl,fl_gdss3
98 logical :: first_grbtbl
100 public num_pset,pset,nrecout,gribit2,grib_info_init,first_grbtbl,grib_info_finalize,read_grib2_head,read_grib2_sngle
101 real(8),
EXTERNAL :: timef
108 subroutine grib_info_init()
115 character(len=80) outfile
116 character(len=10) envvar
123 if(pset%grid_num==0) &
125 if(trim(pset%sub_center)==
'') &
126 pset%sub_center=
"ncep_emc"
127 if(trim(pset%version_no)==
'') &
128 pset%version_no=
"v2003"
129 if(trim(pset%local_table_vers_no)==
'') &
130 pset%local_table_vers_no=
"local_table_no"
131 if(trim(pset%sigreftime)==
'') &
132 pset%sigreftime=
"fcst"
133 if(trim(pset%prod_status)==
'') &
134 pset%prod_status=
"oper_test"
135 if(trim(pset%data_type)==
'') &
136 pset%data_type=
"fcst"
137 if(trim(pset%orig_center)==
'') &
138 pset%orig_center=
"nws_ncep"
139 if(trim(pset%time_range_unit)==
'') &
140 pset%time_range_unit=
"hour"
141 if(trim(pset%gen_proc_type)==
'') &
142 pset%gen_proc_type=
"fcst"
143 if(trim(pset%gen_proc)==
'') &
144 pset%gen_proc=
"gfs_avn"
145 if(trim(pset%packing_method)==
'') &
146 pset%packing_method=
"jpeg"
147 if(trim(pset%field_datatype)==
'') &
148 pset%field_datatype=
"flting_pnt"
149 if(trim(pset%comprs_type)==
'') &
150 pset%comprs_type=
"lossless"
151 if(trim(pset%type_ens_fcst)==
'') &
152 pset%type_ens_fcst=
"pos_pert_fcst"
153 if(trim(pset%type_derived_fcst)==
'') &
154 pset%type_derived_fcst=
"unweighted_mean_all_mem"
163 type_of_time_inc=
'same_start_time_fcst_fcst_time_inc'
164 stat_unit_time_key_succ=
'missing'
165 time_inc_betwn_succ_fld=0
171 if(first_grbtbl)
then
172 fl_nametbl=
'params_grib2_tbl_new'
173 call open_and_read_4dot2( fl_nametbl, ierr )
174 if ( ierr /= 0 )
then
175 print*,
'Couldnt open table file - return code was ',ierr
183 end subroutine grib_info_init
187 subroutine grib_info_finalize
193 call close_4dot2(ierr)
195 end subroutine grib_info_finalize
199 subroutine gribit2(post_fname)
200 use ctlblk_mod,
only : im,jm,im_jm,num_procs,me,ista,iend,jsta,jend,ifhr,sdat,ihrst,imin, &
201 mpi_comm_comp,ntlfld,fld_info,datapd,icnt,idsp
207 character(255),
intent(in) :: post_fname
210 integer i,j,k,n,nm,nprm,nlvl,fldlvl1,fldlvl2,cstart,cgrblen,ierr
212 integer fh, clength,lunout
213 integer idisc,icatg,iparm,itblinfo,ntrange,leng_time_range_stat
214 integer,
allocatable :: nfld_pe(:),snfld_pe(:),enfld_pe(:)
215 integer(4),
allocatable :: isdsp(:),iscnt(:),ircnt(:),irdsp(:)
216 integer status(mpi_status_size)
217 integer(kind=MPI_OFFSET_KIND) idisp
218 integer,
allocatable :: ista_pe(:),iend_pe(:)
219 integer,
allocatable :: jsta_pe(:),jend_pe(:)
220 integer,
allocatable :: grbmsglen(:)
221 real,
allocatable :: datafld(:,:)
222 real,
allocatable :: datafldtmp(:)
223 logical,
parameter :: debugprint = .false.
225 character(1),
dimension(:),
allocatable :: cgrib
226 real :: level_unit_conversion
231 allocate(cgrib(max_bytes))
239 nmod=mod(ntlfld,num_procs)
241 allocate(snfld_pe(num_procs),enfld_pe(num_procs),nfld_pe(num_procs))
244 snfld_pe(n)=nfpe*(n-1)+1
245 enfld_pe(n)=snfld_pe(n)+nfpe-1
248 snfld_pe(n)=nfpe*nmod+nf*(n-1-nmod)+1
249 enfld_pe(n)=snfld_pe(n)+nf-1
258 allocate(ista_pe(num_procs),iend_pe(num_procs))
259 call mpi_allgather(ista,1,mpi_integer,ista_pe,1, &
260 mpi_integer,mpi_comm_comp,ierr)
261 call mpi_allgather(iend,1,mpi_integer,iend_pe,1, &
262 mpi_integer,mpi_comm_comp,ierr)
264 allocate(jsta_pe(num_procs),jend_pe(num_procs))
265 call mpi_allgather(jsta,1,mpi_integer,jsta_pe,1, &
266 mpi_integer,mpi_comm_comp,ierr)
267 call mpi_allgather(jend,1,mpi_integer,jend_pe,1, &
268 mpi_integer,mpi_comm_comp,ierr)
275 allocate(bmap(im_jm))
281 if(num_procs==1)
then
284 allocate(datafld(im_jm,ntlfld) )
286 datafld=reshape(datapd,(/im_jm,ntlfld/))
298 call baopenw(lunout,trim(post_fname),ierr)
299 print*,
'write_grib2: opened ',lunout, &
300 'for grib2 data ',trim(post_fname), &
301 'return code is ',ierr
304 nprm=fld_info(i)%ifld
306 fldlvl1=fld_info(i+snfld_pe(me+1)-1)%lvl1
307 fldlvl2=fld_info(i+snfld_pe(me+1)-1)%lvl2
308 if(trim(pset%param(nprm)%table_info)==
'NCEP')
then
317 call search_for_4dot2_entry( &
318 pset%param(nprm)%pname, &
320 idisc, icatg, iparm, ierr)
322 write(6,
'(3(A,I4),A,A)')
' discipline ',idisc, &
323 ' category ',icatg, &
324 ' parameter ',iparm, &
325 ' for var ',trim(pset%param(nprm)%pname)
326 if(index(pset%param(nprm)%shortname,
'IFI_FLIGHT_LEVEL')>0)
then
327 level_unit_conversion=0.3048
329 level_unit_conversion=1
331 call gengrb2msg(idisc,icatg, iparm,nprm,nlvl,fldlvl1,fldlvl2, &
332 fld_info(i)%ntrange,fld_info(i)%tinvstat,datafld(:,i), &
333 cgrib,clength,level_unit_conversion)
335 call wryte(lunout, clength, cgrib)
337 print *,
'WRONG, could not find ',trim(pset%param(nprm)%pname), &
338 " in WMO and NCEP table!, ierr=", ierr
343 call baclose(lunout,ierr)
344 print *,
'finish one grib file'
351 allocate(iscnt(num_procs),isdsp(num_procs))
352 allocate(ircnt(num_procs),irdsp(num_procs))
355 iscnt(n)=(jend_pe(me+1)-jsta_pe(me+1)+1)*(iend_pe(me+1)-ista_pe(me+1)+1)*nfld_pe(n)
356 if(n<num_procs)isdsp(n+1)=isdsp(n)+iscnt(n)
361 ircnt(n)=(jend_pe(n)-jsta_pe(n)+1)*(iend_pe(n)-ista_pe(n)+1)*nfld_pe(me+1)
362 if(n<num_procs)irdsp(n+1)=irdsp(n)+ircnt(n)
366 allocate(datafldtmp(im_jm*nfld_pe(me+1)) )
367 allocate(datafld(im_jm,nfld_pe(me+1)) )
369 call mpi_alltoallv(datapd,iscnt,isdsp,mpi_real, &
370 datafldtmp,ircnt,irdsp,mpi_real,mpi_comm_comp,ierr)
377 do j=jsta_pe(n),jend_pe(n)
378 do i=ista_pe(n),iend_pe(n)
380 datafld((j-1)*im+i,k)=datafldtmp(nm)
386 deallocate(datafldtmp)
396 nprm=fld_info(i+snfld_pe(me+1)-1)%ifld
397 nlvl=fld_info(i+snfld_pe(me+1)-1)%lvl
398 fldlvl1=fld_info(i+snfld_pe(me+1)-1)%lvl1
399 fldlvl2=fld_info(i+snfld_pe(me+1)-1)%lvl2
400 ntrange=fld_info(i+snfld_pe(me+1)-1)%ntrange
401 leng_time_range_stat=fld_info(i+snfld_pe(me+1)-1)%tinvstat
402 if(trim(pset%param(nprm)%table_info)==
'NCEP')
then
411 call search_for_4dot2_entry( &
412 pset%param(nprm)%pname, &
414 idisc, icatg, iparm, ierr)
417 write(6,
'(3(A,I4),A,A)')
' discipline ',idisc, &
418 ' category ',icatg, &
419 ' parameter ',iparm, &
420 ' for var ',trim(pset%param(nprm)%pname)
425 if(index(pset%param(nprm)%shortname,
'IFI_FLIGHT_LEVEL')>0)
then
426 level_unit_conversion=0.3048
428 level_unit_conversion=1
430 call gengrb2msg(idisc,icatg, iparm,nprm,nlvl,fldlvl1,fldlvl2,ntrange, &
431 leng_time_range_stat,datafld(:,i),cgrib(cstart),clength, &
432 level_unit_conversion)
433 cstart=cstart+clength
437 print *,
'WRONG, could not find ',trim(pset%param(nprm)%pname), &
438 " in WMO and NCEP table!"
451 call mpi_barrier(mpi_comm_comp,ierr)
454 call mpi_file_open(mpi_comm_comp,trim(post_fname), &
455 mpi_mode_create+mpi_mode_wronly,mpi_info_null,fh,ierr)
459 allocate(grbmsglen(num_procs))
460 call mpi_allgather(cgrblen,1,mpi_integer,grbmsglen,1,mpi_integer, &
467 idisp=idisp+grbmsglen(n)
470 call mpi_file_write_at(fh,idisp,cgrib,cgrblen,mpi_character,status,ierr)
472 call mpi_file_close(fh,ierr)
476 deallocate(grbmsglen)
477 deallocate(iscnt,isdsp,ircnt,irdsp)
481 deallocate(datafld,bmap,mg)
482 deallocate(nfld_pe,snfld_pe,enfld_pe,jsta_pe,jend_pe)
485 end subroutine gribit2
490 subroutine gengrb2msg(idisc,icatg, iparm,nprm,nlvl,fldlvl1,fldlvl2,ntrange,tinvstat, &
491 datafld1,cgrib,lengrib,level_unit_conversion)
495 use ctlblk_mod,
only : im,jm,im_jm,ifhr,idat,sdat,ihrst,ifmin,imin,fld_info,spval, &
497 use gridspec_mod,
only: maptype
498 use grib2_all_tables_module,
only: g2sec0,g2sec1, &
499 g2sec4_temp0,g2sec4_temp8,g2sec4_temp9,g2sec4_temp44, &
500 g2sec4_temp46,g2sec4_temp48,g2sec5_temp0,g2sec5_temp2, &
501 g2sec5_temp3,g2sec5_temp40,get_g2_sec5packingmethod, &
506 integer,
intent(in) :: idisc,icatg, iparm,nprm,fldlvl1,fldlvl2,ntrange
507 integer,
intent(inout) :: nlvl,tinvstat
508 real,
dimension(:),
intent(in) :: datafld1
509 character(1),
intent(inout) :: cgrib(max_bytes)
510 integer,
intent(inout) :: lengrib
511 real,
intent(in) :: level_unit_conversion
513 integer,
parameter :: igdsmaxlen=200
515 integer,
parameter :: ipdstmplenmax=100
516 integer,
parameter :: ipdstmp4_0len=15
517 integer,
parameter :: ipdstmp4_1len=18
518 integer,
parameter :: ipdstmp4_8len=29
519 integer,
parameter :: ipdstmp4_9len=36
520 integer,
parameter :: ipdstmp4_11len=32
521 integer,
parameter :: ipdstmp4_12len=31
522 integer,
parameter :: ipdstmp4_44len=21
523 integer,
parameter :: ipdstmp4_46len=35
524 integer,
parameter :: ipdstmp4_48len=26
525 integer,
parameter :: ipdstmp4_49len=29
527 integer,
parameter :: idrstmplenmax=50
528 integer,
parameter :: idrstmp5_0len=5
529 integer,
parameter :: idrstmp5_2len=16
530 integer,
parameter :: idrstmp5_3len=18
531 integer,
parameter :: idrstmp5_40len=7
537 integer ipdstmpl(ipdstmplenmax)
539 integer idrstmpl(idrstmplenmax)
542 integer igdtlen,igdtn
544 integer ideflist(100)
545 integer idrsnum,numcoord,ipdsnum
546 integer scaled_val_fixed_sfc2,scale_fct_fixed_sfc1
547 integer scaled_val_fixed_sfc1,scale_fct_fixed_sfc2
548 character(80) fixed_sfc2_type
549 integer idec_scl,ibin_scl,ibmap,inumbits
550 character(80) prob_type
552 integer igdstmpl(igdsmaxlen)
553 integer lat1,lon1,lat2,lon2,lad,ds1
557 integer ierr,ifhrorig,ihr_start
558 integer gefs1,gefs2,gefs3,gefs_status
559 character(len=4) cdum
560 integer perturb_num,num_ens_fcst,e1_type
566 if(trim(pset%gen_proc)==
'gefs')
then
567 call getenv(
'e1',cdum)
568 read(cdum,
'(I4)',iostat=gefs_status)gefs1
571 if(gefs_status /= 0) print *, &
572 "GEFS Run: Could not read e1 envir. var, User needs to set in script"
574 call getenv(
'e2',cdum)
575 read(cdum,
'(I4)',iostat=gefs_status)gefs2
578 if(gefs_status /= 0) print *, &
579 "GEFS Run: Could not read e2 envir. var, User needs to set in script"
583 call getenv(
'e3',cdum)
584 read(cdum,
'(I4)',iostat=gefs_status)gefs3
587 if(gefs_status /= 0) print *, &
588 "GEFS Run: Could not read e3 envir. var, User needs to set in script"
590 print*,
'GEFS env var ',e1_type,perturb_num,num_ens_fcst
593 print *,
"Processing for GEFS and default setting is tmpl4_1 and tmpl4_11"
594 if (trim(pset%param(nprm)%pdstmpl)==
'tmpl4_0')
then
595 pset%param(nprm)%pdstmpl=
'tmpl4_1'
596 elseif (trim(pset%param(nprm)%pdstmpl)==
'tmpl4_8')
then
597 pset%param(nprm)%pdstmpl=
'tmpl4_11'
598 elseif (trim(pset%param(nprm)%pdstmpl)==
'tmpl4_48')
then
599 pset%param(nprm)%pdstmpl=
'tmpl4_49'
606 call g2sec0(idisc,listsec0)
629 call g2sec1(pset%orig_center,pset%sub_center, &
630 pset%version_no,pset%local_table_vers_no,&
631 pset%sigreftime,nint(sdat(3)),nint(sdat(1)),nint(sdat(2)),ihrst,imin, &
632 isec,pset%prod_status,pset%data_type,listsec1)
635 if(trim(pset%gen_proc)==
'gefs')
then
638 if(e1_type==1.or.e1_type==2)
then
640 elseif(e1_type==3.or.e1_type==4)
then
643 print *,
"After g2sec1 call we need to set listsec1(2) = ",listsec1(2)
644 print *,
"After g2sec1 call we need to set listsec1(13) = ",listsec1(13)
649 call gribcreate(cgrib,max_bytes,listsec0,listsec1,ierr)
658 ldfgrd=(maptype==203.and.(trim(pset%param(nprm)%pname)==
'ugrd'.or. &
659 trim(pset%param(nprm)%pname)==
'vgrd'))
660 call getgds(ldfgrd,igdsmaxlen,igdtlen,igds,igdstmpl)
664 call addgrid(cgrib,max_bytes,igds,igdstmpl,igdtlen,ideflist,idefnum,ierr)
702 if(fldlvl1==0.and.fldlvl2==0)
then
704 if(
size(pset%param(nprm)%level)>1.and.
size(pset%param(nprm)%level)>=nlvl)
then
705 scaled_val_fixed_sfc1=nint(pset%param(nprm)%level(nlvl))
706 else if(
size(pset%param(nprm)%level)==1)
then
707 scaled_val_fixed_sfc1=nint(pset%param(nprm)%level(1))
709 scaled_val_fixed_sfc1=0
711 if(
size(pset%param(nprm)%scale_fact_fixed_sfc1)>1.and. &
712 size(pset%param(nprm)%scale_fact_fixed_sfc1)>=nlvl)
then
713 scale_fct_fixed_sfc1=pset%param(nprm)%scale_fact_fixed_sfc1(nlvl)
714 else if(
size(pset%param(nprm)%scale_fact_fixed_sfc1)==1)
then
715 scale_fct_fixed_sfc1=pset%param(nprm)%scale_fact_fixed_sfc1(1)
717 scale_fct_fixed_sfc1=0
722 scaled_val_fixed_sfc1=fldlvl1
723 scale_fct_fixed_sfc1=0
724 scaled_val_fixed_sfc2=fldlvl2
725 scale_fct_fixed_sfc2=0
728 fixed_sfc2_type=pset%param(nprm)%fixed_sfc2_type
729 if(
size(pset%param(nprm)%level2)>1.and.
size(pset%param(nprm)%level2)<nlvl)
then
733 if (
associated(pset%param(nprm)%level2))
then
734 if(
size(pset%param(nprm)%level2)>1.and.
size(pset%param(nprm)%level2)>=nlvl)
then
735 scaled_val_fixed_sfc2=nint(pset%param(nprm)%level2(nlvl))
736 else if(
size(pset%param(nprm)%level2)==1)
then
737 scaled_val_fixed_sfc2=nint(pset%param(nprm)%level2(1))
739 scaled_val_fixed_sfc2=0
742 scaled_val_fixed_sfc2=0
745 if(
size(pset%param(nprm)%scale_fact_fixed_sfc2)>1 .and. &
746 size(pset%param(nprm)%scale_fact_fixed_sfc2)>=nlvl)
then
747 scale_fct_fixed_sfc2=pset%param(nprm)%scale_fact_fixed_sfc2(nlvl)
748 else if(
size(pset%param(nprm)%scale_fact_fixed_sfc2)==1)
then
749 scale_fct_fixed_sfc2=pset%param(nprm)%scale_fact_fixed_sfc2(1)
751 scale_fct_fixed_sfc2=0
757 if(len_trim(fixed_sfc2_type) == 0)
then
760 fixed_sfc2_type =
'missing'
761 pset%param(nprm)%fixed_sfc2_type =
'missing'
764 if(abs(level_unit_conversion-1)>1e-4)
then
767 scaled_val_fixed_sfc1=nint(scaled_val_fixed_sfc1*real(level_unit_conversion,kind=8))
768 scaled_val_fixed_sfc2=nint(scaled_val_fixed_sfc2*real(level_unit_conversion,kind=8))
772 ihr_start = ifhr-tinvstat
773 if((modelname==
'RAPR'.and.vtimeunits==
'FMIN').or.(modelname==
'FV3R'.and.pset%time_range_unit==
"minute"))
then
775 ifhr = ifhr*60 + ifmin
777 tinvstat = tinvstat*60 + ifmin
779 ihr_start = max(0,ifhr-tinvstat)
782 pset%time_range_unit=
"minute"
784 ifhr = ifhr*60 + ifmin
785 ihr_start = max(0,ifhr-tinvstat*60)
794 if(trim(pset%param(nprm)%pdstmpl)==
'tmpl4_0')
then
796 ipdstmpllen=ipdstmp4_0len
797 call g2sec4_temp0(icatg,iparm,pset%gen_proc_type, &
798 pset%gen_proc,hrs_obs_cutoff,min_obs_cutoff, &
799 pset%time_range_unit,ifhr, &
800 pset%param(nprm)%fixed_sfc1_type, &
801 scale_fct_fixed_sfc1, &
802 scaled_val_fixed_sfc1, &
804 scale_fct_fixed_sfc2, &
805 scaled_val_fixed_sfc2, &
806 ipdstmpl(1:ipdstmpllen))
808 elseif(trim(pset%param(nprm)%pdstmpl)==
'tmpl4_1')
then
810 ipdstmpllen=ipdstmp4_1len
811 call g2sec4_temp1(icatg,iparm,pset%gen_proc_type, &
812 pset%gen_proc,hrs_obs_cutoff,min_obs_cutoff, &
813 pset%time_range_unit,ifhr, &
814 pset%param(nprm)%fixed_sfc1_type, &
815 scale_fct_fixed_sfc1, &
816 scaled_val_fixed_sfc1, &
818 scale_fct_fixed_sfc2, &
819 scaled_val_fixed_sfc2, &
820 pset%type_ens_fcst,perturb_num,num_ens_fcst, &
821 ipdstmpl(1:ipdstmpllen))
824 elseif(trim(pset%param(nprm)%pdstmpl)==
'tmpl4_8')
then
827 ipdstmpllen=ipdstmp4_8len
828 call g2sec4_temp8(icatg,iparm,pset%gen_proc_type, &
829 pset%gen_proc,hrs_obs_cutoff,min_obs_cutoff, &
830 pset%time_range_unit,ihr_start, &
831 pset%param(nprm)%fixed_sfc1_type, &
832 scale_fct_fixed_sfc1, &
833 scaled_val_fixed_sfc1, &
834 pset%param(nprm)%fixed_sfc2_type, &
835 scale_fct_fixed_sfc2, &
836 scaled_val_fixed_sfc2, &
837 idat(3),idat(1),idat(2),idat(4),idat(5), &
838 sec_intvl,ntrange,stat_miss_val, &
839 pset%param(nprm)%stats_proc,type_of_time_inc, &
840 pset%time_range_unit, tinvstat, &
841 stat_unit_time_key_succ,time_inc_betwn_succ_fld, &
842 ipdstmpl(1:ipdstmpllen))
845 elseif(trim(pset%param(nprm)%pdstmpl)==
'tmpl4_9')
then
848 ipdstmpllen=ipdstmp4_9len
849 call g2sec4_temp9(icatg,iparm,pset%gen_proc_type, &
850 pset%gen_proc,hrs_obs_cutoff,min_obs_cutoff, &
851 pset%time_range_unit,ihr_start, &
852 pset%param(nprm)%fixed_sfc1_type, &
853 scale_fct_fixed_sfc1, &
854 scaled_val_fixed_sfc1, &
855 pset%param(nprm)%fixed_sfc2_type, &
856 scale_fct_fixed_sfc2, &
857 scaled_val_fixed_sfc2, &
858 prob_num,tot_num_prob, &
859 pset%param(nprm)%prob_type, &
860 pset%param(nprm)%scale_fact_lower_limit, &
861 pset%param(nprm)%scale_val_lower_limit, &
862 pset%param(nprm)%scale_fact_upper_limit, &
863 pset%param(nprm)%scale_val_upper_limit, &
864 idat(3),idat(1),idat(2),idat(4),idat(5), &
865 sec_intvl,ntrange,stat_miss_val, &
866 pset%param(nprm)%stats_proc,type_of_time_inc, &
867 pset%time_range_unit, tinvstat, &
868 stat_unit_time_key_succ,time_inc_betwn_succ_fld, &
869 ipdstmpl(1:ipdstmpllen))
872 elseif(trim(pset%param(nprm)%pdstmpl)==
'tmpl4_11')
then
874 ipdstmpllen=ipdstmp4_11len
875 call g2sec4_temp11(icatg,iparm,pset%gen_proc_type, &
876 pset%gen_proc,hrs_obs_cutoff,min_obs_cutoff, &
877 pset%time_range_unit,ifhr-tinvstat, &
878 pset%param(nprm)%fixed_sfc1_type, &
879 scale_fct_fixed_sfc1, &
880 scaled_val_fixed_sfc1, &
881 pset%param(nprm)%fixed_sfc2_type, &
882 scale_fct_fixed_sfc2, &
883 scaled_val_fixed_sfc2, &
884 pset%type_ens_fcst,perturb_num,num_ens_fcst, &
885 idat(3),idat(1),idat(2),idat(4),idat(5), &
886 sec_intvl,ntrange,stat_miss_val, &
887 pset%param(nprm)%stats_proc,type_of_time_inc, &
888 pset%time_range_unit, tinvstat, &
889 stat_unit_time_key_succ,time_inc_betwn_succ_fld, &
890 ipdstmpl(1:ipdstmpllen))
893 elseif(trim(pset%param(nprm)%pdstmpl)==
'tmpl4_12')
then
895 ipdstmpllen=ipdstmp4_12len
896 call g2sec4_temp12(icatg,iparm,pset%gen_proc_type, &
897 pset%gen_proc,hrs_obs_cutoff,min_obs_cutoff, &
898 pset%time_range_unit,ifhr-tinvstat, &
899 pset%param(nprm)%fixed_sfc1_type, &
900 scale_fct_fixed_sfc1, &
901 scaled_val_fixed_sfc1, &
903 scale_fct_fixed_sfc2, &
904 scaled_val_fixed_sfc2, &
905 pset%type_derived_fcst,num_ens_fcst, &
906 idat(3),idat(1),idat(2),idat(4),idat(5), &
907 sec_intvl,ntrange,stat_miss_val, &
908 pset%param(nprm)%stats_proc,type_of_time_inc, &
909 pset%time_range_unit, tinvstat, &
910 stat_unit_time_key_succ,time_inc_betwn_succ_fld, &
911 ipdstmpl(1:ipdstmpllen))
914 elseif(trim(pset%param(nprm)%pdstmpl)==
'tmpl4_44')
then
917 ipdstmpllen=ipdstmp4_44len
918 call g2sec4_temp44(icatg,iparm,pset%param(nprm)%aerosol_type, &
919 pset%param(nprm)%typ_intvl_size, &
920 pset%param(nprm)%scale_fact_1st_size, &
921 pset%param(nprm)%scale_val_1st_size, &
922 pset%param(nprm)%scale_fact_2nd_size, &
923 pset%param(nprm)%scale_val_2nd_size, &
924 pset%gen_proc_type, &
925 pset%gen_proc,hrs_obs_cutoff,min_obs_cutoff, &
926 pset%time_range_unit,ifhr, &
927 pset%param(nprm)%fixed_sfc1_type, &
928 scale_fct_fixed_sfc1, &
929 scaled_val_fixed_sfc1, &
930 pset%param(nprm)%fixed_sfc2_type, &
931 scale_fct_fixed_sfc2, &
932 scaled_val_fixed_sfc2, &
933 ipdstmpl(1:ipdstmpllen))
936 elseif(trim(pset%param(nprm)%pdstmpl)==
'tmpl4_46')
then
939 ipdstmpllen=ipdstmp4_46len
940 call g2sec4_temp46(icatg,iparm,pset%param(nprm)%aerosol_type, &
941 pset%param(nprm)%typ_intvl_size, &
942 pset%param(nprm)%scale_fact_1st_size, &
943 pset%param(nprm)%scale_val_1st_size, &
944 pset%param(nprm)%scale_fact_2nd_size, &
945 pset%param(nprm)%scale_val_2nd_size, &
946 pset%gen_proc_type, &
947 pset%gen_proc,hrs_obs_cutoff,min_obs_cutoff, &
948 pset%time_range_unit,ifhr, &
949 pset%param(nprm)%fixed_sfc1_type, &
950 scale_fct_fixed_sfc1, &
951 scaled_val_fixed_sfc1, &
952 pset%param(nprm)%fixed_sfc2_type, &
953 scale_fct_fixed_sfc2, &
954 scaled_val_fixed_sfc2, &
955 idat(3),idat(1),idat(2),idat(4),idat(5), &
956 sec_intvl,ntrange,stat_miss_val, &
957 pset%param(nprm)%stats_proc,type_of_time_inc, &
958 pset%time_range_unit, tinvstat, &
959 stat_unit_time_key_succ,time_inc_betwn_succ_fld, &
960 ipdstmpl(1:ipdstmpllen))
964 elseif(trim(pset%param(nprm)%pdstmpl)==
'tmpl4_48')
then
967 ipdstmpllen=ipdstmp4_48len
968 call g2sec4_temp48(icatg,iparm,pset%param(nprm)%aerosol_type, &
969 pset%param(nprm)%typ_intvl_size, &
970 pset%param(nprm)%scale_fact_1st_size, &
971 pset%param(nprm)%scale_val_1st_size, &
972 pset%param(nprm)%scale_fact_2nd_size, &
973 pset%param(nprm)%scale_val_2nd_size, &
974 pset%param(nprm)%typ_intvl_wvlen, &
975 pset%param(nprm)%scale_fact_1st_wvlen, &
976 pset%param(nprm)%scale_val_1st_wvlen, &
977 pset%param(nprm)%scale_fact_2nd_wvlen, &
978 pset%param(nprm)%scale_val_2nd_wvlen, &
979 pset%gen_proc_type, &
980 pset%gen_proc,hrs_obs_cutoff,min_obs_cutoff, &
981 pset%time_range_unit,ifhr, &
982 pset%param(nprm)%fixed_sfc1_type, &
983 scale_fct_fixed_sfc1, &
984 scaled_val_fixed_sfc1, &
985 pset%param(nprm)%fixed_sfc2_type, &
986 scale_fct_fixed_sfc2, &
987 scaled_val_fixed_sfc2, &
988 ipdstmpl(1:ipdstmpllen))
992 elseif(trim(pset%param(nprm)%pdstmpl)==
'tmpl4_49')
then
995 ipdstmpllen=ipdstmp4_49len
996 call g2sec4_temp49(icatg,iparm,pset%param(nprm)%aerosol_type, &
997 pset%param(nprm)%typ_intvl_size, &
998 pset%param(nprm)%scale_fact_1st_size, &
999 pset%param(nprm)%scale_val_1st_size, &
1000 pset%param(nprm)%scale_fact_2nd_size, &
1001 pset%param(nprm)%scale_val_2nd_size, &
1002 pset%param(nprm)%typ_intvl_wvlen, &
1003 pset%param(nprm)%scale_fact_1st_wvlen, &
1004 pset%param(nprm)%scale_val_1st_wvlen, &
1005 pset%param(nprm)%scale_fact_2nd_wvlen, &
1006 pset%param(nprm)%scale_val_2nd_wvlen, &
1007 pset%gen_proc_type, &
1008 pset%gen_proc,hrs_obs_cutoff,min_obs_cutoff, &
1009 pset%time_range_unit,ifhr, &
1010 pset%param(nprm)%fixed_sfc1_type, &
1011 scale_fct_fixed_sfc1, &
1012 scaled_val_fixed_sfc1, &
1013 pset%param(nprm)%fixed_sfc2_type, &
1014 scale_fct_fixed_sfc2, &
1015 scaled_val_fixed_sfc2, &
1016 pset%type_ens_fcst,perturb_num,num_ens_fcst, &
1017 ipdstmpl(1:ipdstmpllen))
1023 if((modelname==
'RAPR'.or.modelname==
'FV3R').and.vtimeunits==
'FMIN')
then
1036 call get_g2_sec5packingmethod(pset%packing_method,idrsnum,ierr)
1037 if(maxval(datafld1)==minval(datafld1))
then
1041 if(modelname==
'RAPR')
then
1042 if((abs(maxval(datafld1)-minval(datafld1)) < 1.1) .and. (datafld1(1) > 500.0))
then
1047 if(trim(pset%param(nprm)%shortname)==
'UGRD_ON_SPEC_HGT_LVL_ABOVE_GRND_10m'.or.&
1048 trim(pset%param(nprm)%shortname)==
'VGRD_ON_SPEC_HGT_LVL_ABOVE_GRND_10m'.or.&
1049 trim(pset%param(nprm)%shortname)==
'VWSH_ON_SPEC_HGT_LVL_ABOVE_GRND'.or.&
1050 trim(pset%param(nprm)%shortname)==
'VUCSH_ON_SPEC_HGT_LVL_ABOVE_GRND_0-1km'.or.&
1051 trim(pset%param(nprm)%shortname)==
'VVCSH_ON_SPEC_HGT_LVL_ABOVE_GRND_0-1km'.or.&
1052 trim(pset%param(nprm)%shortname)==
'VUCSH_ON_SPEC_HGT_LVL_ABOVE_GRND_0-6km'.or.&
1053 trim(pset%param(nprm)%shortname)==
'VVCSH_ON_SPEC_HGT_LVL_ABOVE_GRND_0-6km')
then
1066 if(any(datafld1>=spval))
then
1068 where(abs(datafld1)>=spval)bmap=.false.
1071 if(
size(pset%param(nprm)%level)==
size(pset%param(nprm)%scale))
then
1072 if(
size(pset%param(nprm)%level)==1) nlvl=1
1074 fldscl=nint(pset%param(nprm)%scale(nlvl))
1078 else if (
size(pset%param(nprm)%scale)==1)
then
1079 fldscl=nint(pset%param(nprm)%scale(1))
1083 if(trim(pset%param(nprm)%pname)==
'APCP' .or. &
1084 trim(pset%param(nprm)%pname)==
'ACPCP' .or. &
1085 trim(pset%param(nprm)%pname)==
'NCPCP')mxbit=22
1086 if(mxbit>16)print*,
'increased MXBIT for ', pset%param(nprm)%pname,mxbit
1087 call g2getbits(mxbit,ibmap,fldscl,
size(datafld1),bmap,datafld1,ibin_scl,idec_scl,inumbits)
1089 if( idrsnum==40 )
then
1090 idrstmplen=idrstmp5_40len
1091 call g2sec5_temp40(idec_scl,ibin_scl,inumbits,pset%comprs_type,idrstmpl(1:idrstmplen))
1092 elseif( idrsnum==2 )
then
1093 idrstmplen=idrstmp5_2len
1094 call g2sec5_temp2(idec_scl,ibin_scl,idrstmpl(1:idrstmplen))
1095 elseif( idrsnum==3 )
then
1096 idrstmplen=idrstmp5_3len
1097 call g2sec5_temp3(idec_scl,ibin_scl,pset%order_of_sptdiff,idrstmpl(1:idrstmplen))
1098 elseif( idrsnum==0 )
then
1099 idrstmplen=idrstmp5_0len
1100 call g2sec5_temp0(idec_scl,ibin_scl,inumbits,idrstmpl(1:idrstmplen))
1112 call addfield(cgrib,max_bytes,ipdsnum,ipdstmpl(1:ipdstmpllen), &
1113 ipdstmpllen,coordlist,numcoord,idrsnum,idrstmpl, &
1114 idrstmplen ,datafld1,im_jm,ibmap,bmap,ierr)
1121 call gribend(cgrib,max_bytes,lengrib,ierr)
1124 end subroutine gengrb2msg
1139 subroutine read_grib2_head(filenameG2,nx,ny,nz,rlonmin,rlatmax,rdx,rdy)
1143 character*256,
intent(in) :: filenameg2
1144 integer,
intent(out) :: nx,ny,nz
1145 real,
intent(out) :: rlonmin,rlatmax
1146 real*8,
intent(out) :: rdx,rdy
1149 type(gribfield) :: gfld
1150 logical :: expand=.true.
1152 character(len=1),
allocatable,
dimension(:) :: cgrib
1153 integer,
parameter :: msk1=32000
1154 integer :: lskip, lgrib,iseek
1156 integer :: icount , lengrib
1157 integer :: listsec0(3)
1158 integer :: listsec1(13)
1159 integer year, month, day, hour, minute, second, fcst
1160 integer :: numfields,numlocal,maxlocal,ierr
1161 integer :: grib_edition
1164 real :: scale_factor,scale_factor2
1167 integer :: nn,n,j,iret
1168 real :: fldmax,fldmin,sum
1183 call baopenr(ifile,trim(filenameg2),iret)
1187 call skgb(ifile,iseek,msk1,lskip,lgrib)
1189 if (lgrib.eq.0)
then
1193 if (lgrib.gt.currlen)
then
1194 if (
allocated(cgrib))
deallocate(cgrib)
1195 allocate(cgrib(lgrib))
1199 call baread(ifile,lskip,lgrib,lengrib,cgrib)
1200 if(lgrib.ne.lengrib)
then
1201 write(*,*)
'ERROR, read_grib2 lgrib ne lengrib', &
1208 call gb_info(cgrib,lengrib,listsec0,listsec1, &
1209 numfields,numlocal,maxlocal,ierr)
1211 write(6,*)
'Error querying GRIB2 message',ierr
1215 grib_edition=listsec0(2)
1216 if (grib_edition.ne.2)
then
1224 call gf_getfld(cgrib,lengrib,n,.false.,expand,gfld,ierr)
1225 year =gfld%idsect(6)
1226 month =gfld%idsect(7)
1228 hour =gfld%idsect(9)
1229 minute=gfld%idsect(10)
1230 second=gfld%idsect(11)
1231 write(*,*)
'year,month,day,hour,minute,second='
1232 write(*,*) year,month,day,hour,minute,second
1233 write(*,*)
'source center =',gfld%idsect(1)
1234 write(*,*)
'Indicator of model =',gfld%ipdtmpl(5)
1235 write(*,*)
'observation level (m)=',gfld%ipdtmpl(12)
1236 write(*,*)
'map projection=',gfld%igdtnum
1237 if (gfld%igdtnum.eq.0)
then
1239 nx = gfld%igdtmpl(8)
1240 ny = gfld%igdtmpl(9)
1242 rdx = gfld%igdtmpl(17)/scale_factor
1243 rdy = gfld%igdtmpl(18)/scale_factor
1244 rlatmax = gfld%igdtmpl(12)/scale_factor
1245 rlonmin = gfld%igdtmpl(13)/scale_factor
1249 else if (gfld%igdtnum.eq.1)
then
1250 nx = gfld%igdtmpl(8)
1251 ny = gfld%igdtmpl(9)
1253 rdx = gfld%igdtmpl(17)/scale_factor
1254 rdy = gfld%igdtmpl(18)/scale_factor
1255 rlatmax = gfld%igdtmpl(12)/scale_factor
1256 rlonmin = gfld%igdtmpl(13)/scale_factor
1260 else if (gfld%igdtnum.eq.30)
then
1261 nx = gfld%igdtmpl(8)
1262 ny = gfld%igdtmpl(9)
1264 rdx = gfld%igdtmpl(15)/scale_factor2
1265 rdy = gfld%igdtmpl(16)/scale_factor2
1266 rlatmax = gfld%igdtmpl(10)/scale_factor
1267 rlonmin = gfld%igdtmpl(11)/scale_factor
1272 write(*,*)
'unknown projection'
1278 CALL baclose(ifile,ierr)
1280 if (
allocated(cgrib))
deallocate(cgrib)
1283 end subroutine read_grib2_head
1291 subroutine read_grib2_sngle(filenameG2,ntot,height,var)
1295 character*256,
intent(in) :: filenameg2
1296 integer,
intent(in) :: ntot
1297 real,
intent(out) :: var(ntot)
1298 integer,
intent(out) :: height
1301 type(gribfield) :: gfld
1302 logical :: expand=.true.
1304 character(len=1),
allocatable,
dimension(:) :: cgrib
1305 integer,
parameter :: msk1=32000
1306 integer :: lskip, lgrib,iseek
1308 integer :: icount , lengrib
1309 integer :: listsec0(3)
1310 integer :: listsec1(13)
1311 integer year, month, day, hour, minute, second, fcst
1312 integer :: numfields,numlocal,maxlocal,ierr
1313 integer :: grib_edition
1316 real :: dx,dy,lat1,lon1,rtnum, nlat
1317 real :: ref_value,bin_scale_fac,dec_scale_fac,bit_number,field_type
1319 real :: scale_factor,scale_factor2
1322 integer :: nn,n,j,iret
1323 real :: fldmax,fldmin,sum
1338 call baopenr(ifile,trim(filenameg2),iret)
1342 call skgb(ifile,iseek,msk1,lskip,lgrib)
1344 if (lgrib.eq.0)
then
1348 if (lgrib.gt.currlen)
then
1349 if (
allocated(cgrib))
deallocate(cgrib)
1350 allocate(cgrib(lgrib))
1354 call baread(ifile,lskip,lgrib,lengrib,cgrib)
1355 if(lgrib.ne.lengrib)
then
1356 write(*,*)
'ERROR, read_grib2 lgrib ne lengrib', &
1364 call gb_info(cgrib,lengrib,listsec0,listsec1, &
1365 numfields,numlocal,maxlocal,ierr)
1367 write(6,*)
'Error querying GRIB2 message',ierr
1371 grib_edition=listsec0(2)
1372 if (grib_edition.ne.2)
then
1380 call gf_getfld(cgrib,lengrib,n,.false.,expand,gfld,ierr)
1381 year =gfld%idsect(6)
1382 month =gfld%idsect(7)
1384 hour =gfld%idsect(9)
1385 minute=gfld%idsect(10)
1386 second=gfld%idsect(11)
1393 height=gfld%ipdtmpl(12)
1394 if (gfld%igdtnum.eq.0)
then
1396 nx = gfld%igdtmpl(8)
1397 ny = gfld%igdtmpl(9)
1398 dx = gfld%igdtmpl(17)/scale_factor
1399 dy = gfld%igdtmpl(18)/scale_factor
1400 lat1 = gfld%igdtmpl(12)/scale_factor
1401 lon1 = gfld%igdtmpl(13)/scale_factor
1405 else if (gfld%igdtnum.eq.1)
then
1406 nx = gfld%igdtmpl(8)
1407 ny = gfld%igdtmpl(9)
1408 dx = gfld%igdtmpl(17)/scale_factor
1409 dy = gfld%igdtmpl(18)/scale_factor
1410 lat1 = gfld%igdtmpl(12)/scale_factor
1411 lon1 = gfld%igdtmpl(13)/scale_factor
1415 else if (gfld%igdtnum.eq.30)
then
1416 nx = gfld%igdtmpl(8)
1417 ny = gfld%igdtmpl(9)
1418 dx = gfld%igdtmpl(15)/scale_factor2
1419 dy = gfld%igdtmpl(16)/scale_factor2
1420 lat1 = gfld%igdtmpl(10)/scale_factor
1421 lon1 = gfld%igdtmpl(11)/scale_factor
1426 rtnum = gfld%idrtnum
1428 ref_value = gfld%idrtmpl(1)
1429 bin_scale_fac = gfld%idrtmpl(2)
1430 dec_scale_fac = gfld%idrtmpl(3)
1431 bit_number = gfld%idrtmpl(4)
1432 field_type = gfld%idrtmpl(5)
1433 bit_map = gfld%ibmap
1440 else if (gfld%igdtnum.eq.40)
then
1441 nx = gfld%igdtmpl(8)
1442 ny = gfld%igdtmpl(9)
1443 lat1 = gfld%igdtmpl(12)/scale_factor
1444 lon1 = gfld%igdtmpl(13)/scale_factor
1445 dx = gfld%igdtmpl(17)/scale_factor
1446 nlat = gfld%igdtmpl(18)
1447 write(*,*) gfld%igdtnum, nx, ny, lat1, lon1, dx, nlat
1449 write(*,*)
'unknown projection'
1454 num_fields:
do n = 1, numfields
1456 call gf_getfld(cgrib,lengrib,n,.true.,expand,gfld,ierr)
1458 write(*,*)
' ERROR extracting field gf_getfld = ',ierr
1467 if(ntot .ne. gfld%ngrdpts)
then
1468 write(*,*)
'Error, wrong dimension ',ntot, gfld%ngrdpts
1480 CALL baclose(ifile,ierr)
1481 if (
allocated(cgrib))
deallocate(cgrib)
1485 end subroutine read_grib2_sngle
1490 subroutine g2sec3tmpl40(nx,nY,lat1,lon1,lat2,lon2,lad,ds1,len3,igds,ifield3)
1493 integer(4),
intent(inout) :: nx,ny,lat1,lon1,lat2,lon2,lad,ds1,len3
1494 integer(4),
intent(inout) :: ifield3(len3),igds(5)
1530 end subroutine g2sec3tmpl40
1534 subroutine g2getbits(MXBIT,ibm,scl,len,bmap,g,ibs,ids,nbits)
1555 INTEGER,
INTENT(IN) :: mxbit,ibm,len
1556 LOGICAL*1,
INTENT(IN) :: bmap(len)
1557 REAL,
INTENT(IN) :: scl,g(len)
1558 INTEGER,
INTENT(OUT) :: ibs,ids,nbits
1563 real,
PARAMETER :: alog2=0.69314718056,hpeps=0.500001
1566 INTEGER :: i,i1,icnt,ipo,le,irange
1567 REAL :: ground,gmin,gmax,s,rmin,rmax,range,rr,rng2,po,rln2
1569 DATA rln2/0.69314718/
1578 gmax = max(gmax,g(i))
1579 gmin = min(gmin,g(i))
1594 gmax = max(gmax,g(i))
1595 gmin = min(gmin,g(i))
1612 IF ( range <= 1.e-30 )
THEN
1617 IF ( scl == 0.0 )
THEN
1620 ELSE IF ( scl > 0.0 )
THEN
1621 ipo = int(alog10( range ))
1623 if(ipo<0.and.ipo+scl<-20)
then
1624 print *,
'for small range,ipo=',ipo,
'ipo+scl=',ipo+scl,
'scl=',scl
1629 IF ( range < 1.00 ) ipo = ipo - 1
1630 po = float(ipo) - scl + 1.
1632 rr = range * 10. ** ( -po )
1633 nbits = int( alog( rr ) / rln2 ) + 1
1636 rng2 = range * 2. ** (-ibs)
1637 nbits = int( alog( rng2 ) / rln2 ) + 1
1643 IF(abs(gmin) >= 1.)
THEN
1644 ids = -int(alog10(abs(gmin)))
1645 ELSE IF (abs(gmin) < 1.0.AND.abs(gmin) > 0.0)
THEN
1646 ids = -int(alog10(abs(gmin)))+1
1651 nbits = min(nbits,mxbit)
1654 IF ( scl > 0.0 )
THEN
1661 gmax = max(gmax,g(i)*s)
1662 gmin = min(gmin,g(i)*s)
1678 gmax = max(gmax,g(i)*s)
1679 gmin = min(gmin,g(i)*s)
1689 if(gmax == gmin)
then
1692 ibs = nint(alog(range/(2.**nbits-0.5))/alog2+hpeps)
1700 END subroutine g2getbits
1704 subroutine getgds(ldfgrd,len3,ifield3len,igds,ifield3)
1706 use ctlblk_mod,
only : im,jm,gdsdegr,modelname
1707 use gridspec_mod,
only: dxval,dyval,cenlat,cenlon,latstart,lonstart,latlast, &
1708 & lonlast,maptype,standlon,latstartv,cenlatv,lonstartv, &
1709 cenlonv,truelat1,truelat2,latstart_r,lonstart_r, &
1714 logical,
intent(in) :: ldfgrd
1715 integer(4),
intent(in) :: len3
1716 integer(4),
intent(out) :: ifield3len
1717 integer(4),
intent(inout) :: ifield3(len3),igds(5)
1728 IF(maptype == 1)
THEN
1736 ifield3(10) = latstart
1737 ifield3(11) = lonstart
1740 if (modelname ==
'FV3R')
then
1744 ifield3(13) = truelat1
1745 ifield3(14) = standlon
1754 ifield3(19) = truelat1
1755 ifield3(20) = truelat2
1758 ELSE IF(maptype == 2)
THEN
1766 ifield3(10) = latstart
1767 ifield3(11) = lonstart
1769 ifield3(13) = truelat1
1770 ifield3(14) = standlon
1777 ELSE IF(maptype==3)
THEN
1785 ifield3(10) = latstart
1786 ifield3(11) = lonstart
1788 ifield3(13) = truelat1
1789 ifield3(14) = latlast
1790 ifield3(15) = lonlast
1797 ELSE IF(maptype == 203)
THEN
1806 if(.not.ldfgrd)
then
1809 ifield3(11) = 45000000
1811 ifield3(12) = latstart
1812 ifield3(13) = lonstart
1814 ifield3(15) = cenlat
1815 ifield3(16) = cenlon
1821 ELSE IF(maptype == 205 .OR. maptype == 6)
THEN
1830 if(.not.ldfgrd)
then
1833 ifield3(11) = 45000000
1835 ifield3(12) = latstart
1836 ifield3(13) = lonstart
1838 ifield3(15) = cenlat
1839 ifield3(16) = cenlon
1843 ifield3(20) = latlast
1844 ifield3(21) = lonlast
1847 ELSE IF(maptype == 207)
THEN
1856 if(.not.ldfgrd)
then
1859 ifield3(11) = 45000000
1861 ifield3(12) = latstart_r
1862 ifield3(13) = lonstart_r
1865 if(modelname==
'FV3R')
then
1869 ifield3(15) = latlast_r
1870 ifield3(16) = lonlast_r
1874 ifield3(20) = cenlat-90.*gdsdegr
1875 ifield3(21) = cenlon
1881 ELSE IF(maptype == 4 )
THEN
1891 ifield3(12) = latstart
1892 ifield3(13) = lonstart
1894 ifield3(15) = latlast
1895 ifield3(16) = lonlast
1896 ifield3(17) = nint(360./(im)*1000000.)
1897 ifield3(18) = nint(jm/2.0)
1898 if( latstart < latlast )
then
1905 ELSE IF(maptype == 0 )
THEN
1915 ifield3(12) = latstart
1916 ifield3(13) = lonstart
1918 ifield3(15) = latlast
1919 ifield3(16) = lonlast
1930 if( latstart < latlast )
then
1939 end subroutine getgds
1943 end module grib2_module
Parameters that are used to read in Perl XML processed flat file and handle parameter marshalling for...
Parameters that are used to read in Perl XML processed flat file and handle parameter marshalling for...