fixing cmake

src/CMakeLists.txt

@@ -5,7 +5,7 @@
 
 # These are the fortran source files.
 set(fortran_src addfield.F90 addgrid.F90 addlocal.F90 cmplxpack.F90 compack.F90
-comunpack.f drstemplates.F90 g2_gbytesc.F90 g2grids.F90 gb_info.F90
+comunpack.F90 drstemplates.F90 g2_gbytesc.F90 g2grids.F90 gb_info.F90
 getdim.F90 getfield.F90 getg2i.F90 getg2ir.F90 getgb2.F90 getgb2l.F90 getgb2p.F90
 getgb2r.F90 getgb2rp.F90 getgb2s.F90 getidx.F90 getlocal.F90 getpoly.F90
 gettemplates.F90 gf_free.F90 gf_getfld.F90 gf_unpack1.F90 gf_unpack2.F90

src/addfield.F90

@@ -412,8 +412,8 @@ subroutine addfield(cgrib,lcgrib,ipdsnum,ipdstmpl,ipdstmplen, &
   ! If specifying a previously defined bit-map, make sure
   ! one already exists in the current GRIB message.
   if ((ibmap.eq.254).and.(.not.isprevbmap)) then
-     print *,'addfield: Requested previously defined bitmap, ',
-     & ' but one does not exist in the current GRIB message.'
+     print *,'addfield: Requested previously defined bitmap, ', &
+          ' but one does not exist in the current GRIB message.'
      ierr=8
      return
   endif

src/comunpack.F90

@@ -0,0 +1,319 @@
+!> @file
+!> @brief Unpack a data field that was packed using a complex packing
+!> algorithm as defined in the GRIB2 documention.
+!> @author Stephen Gilbert @date 2000-06-21
+
+!> Unpack a data field that was packed using a complex packing
+!> algorithm as defined in the GRIB2 documention.
+!>
+!> This subroutine Supports GRIB2 complex packing templates with or
+!> without spatial differences: It supports GRIB2 complex packing
+!> templates with or without spatial differences: Data Representation
+!> Tables [5.2]
+!> (https://www.nco.ncep.noaa.gov/pmb/docs/grib2/grib2_doc/grib2_table5-2.shtml)
+!> and
+!> [5.3](https://www.nco.ncep.noaa.gov/pmb/docs/grib2/grib2_doc/grib2_table5-3.shtml)).
+!>
+!> @param[in] cpack The packed data field (character*1 array).
+!> @param[in] len Length of packed field cpack.
+!> @param[in] lensec Length of section 7 (used for error checking).
+!> @param[in] idrsnum Data Representation Template number 5.N. Must
+!> equal 2 or 3.
+!> @param[in] idrstmpl The array of values for Data Representation
+!> Template 5.2 or 5.3.
+!> @param[in] ndpts The number of data values to unpack.
+!> @param[out] fld Contains the unpacked data values.
+!> @param[out] ier Error return:
+!> - 0 = No error.
+!> - 1 = Problem - inconsistent group lengths of widths.
+!>
+!> @author Stephen Gilbert @date 2000-06-21
+subroutine comunpack(cpack,len,lensec,idrsnum,idrstmpl,ndpts, &
+     fld,ier)
+
+  character(len=1),intent(in) :: cpack(len)
+  integer,intent(in) :: ndpts,len
+  integer,intent(in) :: idrstmpl(*)
+  real,intent(out) :: fld(ndpts)
+
+  integer,allocatable :: ifld(:),ifldmiss(:)
+  integer(4) :: ieee
+  integer,allocatable :: gref(:),gwidth(:),glen(:)
+  real :: ref,bscale,dscale,rmiss1,rmiss2
+  !      real :: fldo(6045)
+  integer :: totBit, totLen
+
+  ier=0
+  !print *,'IDRSTMPL: ',(idrstmpl(j),j=1,16)
+  ieee = idrstmpl(1)
+  call rdieee(ieee,ref,1)
+  bscale = 2.0**real(idrstmpl(2))
+  dscale = 10.0**real(-idrstmpl(3))
+  nbitsgref = idrstmpl(4)
+  itype = idrstmpl(5)
+  ngroups = idrstmpl(10)
+  nbitsgwidth = idrstmpl(12)
+  nbitsglen = idrstmpl(16)
+  if (idrsnum.eq.3) then
+     nbitsd=idrstmpl(18)*8
+  endif
+
+  !   Constant field
+
+  if (ngroups.eq.0) then
+     do j=1,ndpts
+        fld(j)=ref
+     enddo
+     return
+  endif
+
+  iofst=0
+  allocate(ifld(ndpts),stat=is)
+  !print *,'ALLOC ifld: ',is,ndpts
+  allocate(gref(ngroups),stat=is)
+  !print *,'ALLOC gref: ',is
+  allocate(gwidth(ngroups),stat=is)
+  !print *,'ALLOC gwidth: ',is
+  !
+  !  Get missing values, if supplied
+  !
+  if ( idrstmpl(7).eq.1 ) then
+     if (itype.eq.0) then
+        call rdieee(idrstmpl(8),rmiss1,1)
+     else
+        rmiss1=real(idrstmpl(8))
+     endif
+  elseif ( idrstmpl(7).eq.2 ) then
+     if (itype.eq.0) then
+        call rdieee(idrstmpl(8),rmiss1,1)
+        call rdieee(idrstmpl(9),rmiss2,1)
+     else
+        rmiss1=real(idrstmpl(8))
+        rmiss2=real(idrstmpl(9))
+     endif
+  endif
+  !print *,'RMISSs: ',rmiss1,rmiss2,ref
+  ! 
+  !  Extract Spatial differencing values, if using DRS Template 5.3
+  !
+  if (idrsnum.eq.3) then
+     if (nbitsd.ne.0) then
+        call g2_gbytec(cpack,ival1,iofst,nbitsd)
+        iofst=iofst+nbitsd
+        if (idrstmpl(17).eq.2) then
+           call g2_gbytec(cpack,ival2,iofst,nbitsd)
+           iofst=iofst+nbitsd
+        endif
+        call g2_gbytec(cpack,isign,iofst,1)
+        iofst=iofst+1
+        call g2_gbytec(cpack,minsd,iofst,nbitsd-1)
+        iofst=iofst+nbitsd-1
+        if (isign.eq.1) minsd=-minsd
+     else
+        ival1=0
+        ival2=0
+        minsd=0
+     endif
+     !print *,'SDu ',ival1,ival2,minsd,nbitsd
+  endif
+  !
+  !  Extract Each Group's reference value
+  !
+  !print *,'SAG1: ',nbitsgref,ngroups,iofst
+  if (nbitsgref.ne.0) then
+     call g2_gbytesc(cpack,gref,iofst,nbitsgref,0,ngroups)
+     itemp=nbitsgref*ngroups
+     iofst=iofst+(itemp)
+     if (mod(itemp,8).ne.0) iofst=iofst+(8-mod(itemp,8))
+  else
+     gref(1:ngroups)=0
+  endif
+  !write(78,*)'GREFs: ',(gref(j),j=1,ngroups)
+  !
+  !  Extract Each Group's bit width
+  !
+  !print *,'SAG2: ',nbitsgwidth,ngroups,iofst,idrstmpl(11)
+  if (nbitsgwidth.ne.0) then
+     call g2_gbytesc(cpack,gwidth,iofst,nbitsgwidth,0,ngroups)
+     itemp=nbitsgwidth*ngroups
+     iofst=iofst+(itemp)
+     if (mod(itemp,8).ne.0) iofst=iofst+(8-mod(itemp,8))
+  else
+     gwidth(1:ngroups)=0
+  endif
+  do j=1,ngroups
+     gwidth(j)=gwidth(j)+idrstmpl(11)
+  enddo
+  !write(78,*)'GWIDTHs: ',(gwidth(j),j=1,ngroups)
+  !
+  !  Extract Each Group's length (number of values in each group)
+  !
+  allocate(glen(ngroups),stat=is)
+  !print *,'ALLOC glen: ',is
+  !print *,'SAG3: ',nbitsglen,ngroups,iofst,idrstmpl(14),idrstmpl(13)
+  if (nbitsglen.ne.0) then
+     call g2_gbytesc(cpack,glen,iofst,nbitsglen,0,ngroups)
+     itemp=nbitsglen*ngroups
+     iofst=iofst+(itemp)
+     if (mod(itemp,8).ne.0) iofst=iofst+(8-mod(itemp,8))
+  else
+     glen(1:ngroups)=0
+  endif
+  do j=1,ngroups
+     glen(j)=(glen(j)*idrstmpl(14))+idrstmpl(13)
+  enddo
+  glen(ngroups)=idrstmpl(15)
+  !write(78,*)'GLENs: ',(glen(j),j=1,ngroups)
+  !print *,'GLENsum: ',sum(glen)
+  !
+  !  Test to see if the group widths and lengths are consistent with number of
+  !  values, and length of section 7.
+  !
+  totBit = 0
+  totLen = 0
+  do j=1,ngroups
+     totBit = totBit + (gwidth(j)*glen(j));
+     totLen = totLen + glen(j);
+  enddo
+  if (totLen .NE. ndpts) then
+     ier=1
+     return
+  endif
+  if ( (totBit/8) .GT. lensec) then
+     ier=1
+     return
+  endif
+  !
+  !  For each group, unpack data values
+  !
+  if ( idrstmpl(7).eq.0 ) then        ! no missing values
+     n=1
+     do j=1,ngroups
+        !write(78,*)'NGP ',j,gwidth(j),glen(j),gref(j)
+        if (gwidth(j).ne.0) then
+           call g2_gbytesc(cpack,ifld(n),iofst,gwidth(j),0,glen(j))
+           do k=1,glen(j)
+              ifld(n)=ifld(n)+gref(j)
+              n=n+1
+           enddo
+        else
+           ifld(n:n+glen(j)-1)=gref(j)
+           n=n+glen(j)
+        endif
+        iofst=iofst+(gwidth(j)*glen(j))
+     enddo
+  elseif ( idrstmpl(7).eq.1.OR.idrstmpl(7).eq.2 ) then
+     ! missing values included
+     allocate(ifldmiss(ndpts))
+     !ifldmiss=0
+     n=1
+     non=1
+     do j=1,ngroups
+        !print *,'SAGNGP ',j,gwidth(j),glen(j),gref(j)
+        if (gwidth(j).ne.0) then
+           msng1=(2**gwidth(j))-1
+           msng2=msng1-1
+           call g2_gbytesc(cpack,ifld(n),iofst,gwidth(j),0,glen(j))
+           iofst=iofst+(gwidth(j)*glen(j))
+           do k=1,glen(j)
+              if (ifld(n).eq.msng1) then
+                 ifldmiss(n)=1
+              elseif (idrstmpl(7).eq.2.AND.ifld(n).eq.msng2) then
+                 ifldmiss(n)=2
+              else
+                 ifldmiss(n)=0
+                 ifld(non)=ifld(n)+gref(j)
+                 non=non+1
+              endif
+              n=n+1
+           enddo
+        else
+           msng1=(2**nbitsgref)-1
+           msng2=msng1-1
+           if (gref(j).eq.msng1) then
+              ifldmiss(n:n+glen(j)-1)=1
+              !ifld(n:n+glen(j)-1)=0
+           elseif (idrstmpl(7).eq.2.AND.gref(j).eq.msng2) then
+              ifldmiss(n:n+glen(j)-1)=2
+              !ifld(n:n+glen(j)-1)=0
+           else
+              ifldmiss(n:n+glen(j)-1)=0
+              ifld(non:non+glen(j)-1)=gref(j)
+              non=non+glen(j)
+           endif
+           n=n+glen(j)
+        endif
+     enddo
+  endif
+  !write(78,*)'IFLDs: ',(ifld(j),j=1,ndpts)
+
+  if ( allocated(gref) ) deallocate(gref)
+  if ( allocated(gwidth) ) deallocate(gwidth)
+  if ( allocated(glen) ) deallocate(glen)
+  !
+  !  If using spatial differences, add overall min value, and
+  !  sum up recursively
+  !
+  if (idrsnum.eq.3) then         ! spatial differencing
+     if (idrstmpl(17).eq.1) then      ! first order
+        ifld(1)=ival1
+        if ( idrstmpl(7).eq.0 ) then        ! no missing values
+           itemp=ndpts
+        else
+           itemp=non-1
+        endif
+        do n=2,itemp
+           ifld(n)=ifld(n)+minsd
+           ifld(n)=ifld(n)+ifld(n-1)
+        enddo
+     elseif (idrstmpl(17).eq.2) then    ! second order
+        ifld(1)=ival1
+        ifld(2)=ival2
+        if ( idrstmpl(7).eq.0 ) then        ! no missing values
+           itemp=ndpts
+        else
+           itemp=non-1
+        endif
+        do n=3,itemp
+           ifld(n)=ifld(n)+minsd
+           ifld(n)=ifld(n)+(2*ifld(n-1))-ifld(n-2)
+        enddo
+     endif
+     !write(78,*)'IFLDs: ',(ifld(j),j=1,ndpts)
+  endif
+  !
+  !  Scale data back to original form
+  !
+  !print *,'SAGT: ',ref,bscale,dscale
+  if ( idrstmpl(7).eq.0 ) then        ! no missing values
+     do n=1,ndpts
+        fld(n)=((real(ifld(n))*bscale)+ref)*dscale
+        !write(78,*)'SAG ',n,fld(n),ifld(n),bscale,ref,1./dscale
+     enddo
+  elseif ( idrstmpl(7).eq.1.OR.idrstmpl(7).eq.2 ) then
+     ! missing values included
+     non=1
+     do n=1,ndpts
+        if ( ifldmiss(n).eq.0 ) then
+           fld(n)=((real(ifld(non))*bscale)+ref)*dscale
+           !print *,'SAG ',n,fld(n),ifld(non),bscale,ref,dscale
+           non=non+1
+        elseif ( ifldmiss(n).eq.1 ) then
+           fld(n)=rmiss1
+        elseif ( ifldmiss(n).eq.2 ) then
+           fld(n)=rmiss2
+        endif
+     enddo
+     if ( allocated(ifldmiss) ) deallocate(ifldmiss)
+  endif
+
+  if ( allocated(ifld) ) deallocate(ifld)
+
+  !open(10,form='unformatted',recl=24180,access='direct') 
+  !read(10,rec=1) (fldo(k),k=1,6045)
+  !do i =1,6045
+  !  print *,i,fldo(i),fld(i),fldo(i)-fld(i)
+  !enddo
+
+  return
+end subroutine comunpack