Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

ECMWF GRIB Decode Success with new AEC/CCSDS Library? #461

Closed
JKrobNESDIS opened this issue Dec 1, 2023 · 28 comments
Closed

ECMWF GRIB Decode Success with new AEC/CCSDS Library? #461

JKrobNESDIS opened this issue Dec 1, 2023 · 28 comments
Assignees
Labels
question Further information is requested

Comments

@JKrobNESDIS
Copy link

JKrobNESDIS commented Dec 1, 2023

All,

With the release of the new AEC/CCSDS Decode Library, I have been trying to decode the public ECMWF GRIB data...with no success.

https://data.ecmwf.int/forecasts/

Function m_split in 'decode.c' keeps responding with "return M_ERROR - 5"

I'm calling g2c_aecunpackf from the fortran gf_unpack.f & I can see all my arguments are passed intact so....

If someone could test out the ECMWF GRIB from the pure C routines to confirm/deny whether it works correctly would be most helpful.

Thanks in advance,
Jeff Krob
NOAA/NESDIS

@JKrobNESDIS JKrobNESDIS changed the title ECMWF GRIB Decode Success with new AEC Library? ECMWF GRIB Decode Success with new AEC/CCSDS Library? Dec 1, 2023
@EricEngle-NOAA
Copy link
Contributor

Hi @JKrobNESDIS, can you pass along or point me to the ECMWF GRIB2 files you are working with? Even if just 1 message from the file?

@EricEngle-NOAA
Copy link
Contributor

@JKrobNESDIS I quickly read earlier and did not notice the link. I will download a file and take a look.

@EricEngle-NOAA
Copy link
Contributor

I have been looking at this over the last couple of days. I can confirm that this is a bug in the implementation of the AEC compression in the g2c library. I have identified the bug and am working a fix immediately.

@edhartnett would we be able to issue a v1.8.1 patch release?

@edwardhartnett
Copy link
Contributor

We can do a new release if you need it.

@JKrobNESDIS
Copy link
Author

JKrobNESDIS commented Dec 15, 2023 via email

@EricEngle-NOAA
Copy link
Contributor

Hi @JKrobNESDIS. Would you be willing to share your code? Or least the portion that does the ECMWF GRIB2 reading/unpacking using g2c?

@JKrobNESDIS
Copy link
Author

JKrobNESDIS commented Dec 18, 2023 via email

@EricEngle-NOAA
Copy link
Contributor

Feel free to share with me via NOAA email or through the NOAA Google Drive.

edwardhartnett pushed a commit that referenced this issue Dec 20, 2023
* Update for aecunpack.c

This commit fixes a bug in aecunpack.c where the bytes from the decoded
AEC buffer stream are not being properly copied into the ifld array. The
AEC stream is byte aligned.

This commit references #461

* Update aecunpack.c

Commenting out ifld1 since the code logic that uses it is commented out.

This commit references #461

---------

Co-authored-by: Eric Engle <EricEngle-NOAA@users.noreply.github.com>
@edwardhartnett
Copy link
Contributor

Should this issue be closed?

@JKrobNESDIS
Copy link
Author

JKrobNESDIS commented Jun 17, 2024 via email

@JKrobNESDIS
Copy link
Author

JKrobNESDIS commented Jun 17, 2024 via email

@EricEngle-NOAA
Copy link
Contributor

EricEngle-NOAA commented Jun 17, 2024

I will try to devote some time to this in the next couple of weeks. We are spinning up development on NBM v5.0.

@JKrobNESDIS
Copy link
Author

JKrobNESDIS commented Jun 17, 2024 via email

@edwardhartnett
Copy link
Contributor

Any progress on this issue? I'm getting ready to do the next release...

@EricEngle-NOAA
Copy link
Contributor

No progress on this -- have other works tasks in the way right now.

@EricEngle-NOAA
Copy link
Contributor

@JKrobNESDIS - Can you provide code here where your application interfaces with g2c? If I remember correctly, your application is Fortran?

@JKrobNESDIS
Copy link
Author

JKrobNESDIS commented Jul 10, 2024 via email

@EricEngle-NOAA
Copy link
Contributor

Yes. Probably best to just put code text here.

@JKrobNESDIS
Copy link
Author

JKrobNESDIS commented Jul 10, 2024 via email

@JKrobNESDIS
Copy link
Author

Eric - here is my version of gf_unpack7.f that works fine unpacking all other NCEP GRIB

  subroutine gf_unpack7(cgrib,lcgrib,iofst,igdsnum,igdstmpl,
 &		  	idrsnum,idrstmpl,ndpts,fld,ierr)

interface
    function g2c_aecunpackf(cgrib,lensec,idrstmpl,ndpts,fld)  
 &                     bind(c, name="g2c_aecunpackf")
      use iso_c_binding
      character(kind = c_char), intent(in) :: cgrib(*)
      integer(c_int), value, intent(in) :: lensec
  integer(c_int), intent(in) :: idrstmpl(*)
  integer(c_int), value, intent(in) :: ndpts
      real,intent(out) :: fld(*)
  integer(c_int) :: g2c_aecunpackf
    end function g2c_aecunpackf
end interface


  elseif (idrsnum.eq.42) then
	WRITE (99,*) '*****GF_UNPACK7 lensec-5, ndpts	= ',
 &                 lensec-5, ndpts
  WRITE(99,*)'******* CCSDS/AEC UNCOMPRESS ******'
	
ret = g2c_aecunpackf(cgrib(ipos), lensec-5, idrstmpl, ndpts, fld)

	OPEN (99,FILE='NGRB2PCG32.OUT',ACCESS='APPEND',
 &            STATUS='OLD',IOSTAT=JERR)
	WRITE(99,*)'*****AEC END*******',ret

Thanks,
Jeff

@JKrobNESDIS
Copy link
Author

JKrobNESDIS commented Jul 24, 2024 via email

@EricEngle-NOAA
Copy link
Contributor

EricEngle-NOAA commented Jul 24, 2024

@JKrobNESDIS - I have been looking into this over the past few days. When I use grib2io (Python), which links to g2c through a Cython extension module, I do not have any decode/unpack issues when reading from AEC encoded ECMWF GRIB2 files. Here an screenshot example of this in a Jupyter notebook

Screenshot 2024-07-24 at 7 54 37 AM

I am thinking that the Fortran interface to g2_aecunpackf might not be configured correctly. Looking at the function definition from grib2.h

int g2c_aecunpackf(unsigned char *cpack, size_t len, int *idrstmpl, size_t ndpts, float *fld);

Try the following Fortran interface:

interface
  function g2c_aecunpackf(cpack, lensec, idrstmpl, ndpts, fld) bind(C, name="g2c_aecunpackf")
    use, intrinsic :: iso_c_binding
    implicit none
    ! Arguments
    character(kind=c_char), intent(in), dimension(*) :: cpack   ! Equivalent to unsigned char *
    integer(kind=c_size_t), value :: lensec
    integer(kind=c_int), dimension(*), intent(inout) :: idrstmpl
    integer(kind=c_size_t), value :: ndpts
    real(kind=c_float), dimension(*), intent(inout) :: fld
    ! Return type
    integer(kind=c_int) :: g2c_aecunpackf
  end function g2c_aecunpackf
end interface

which only differs from your version by making lecsec and ndpts type c_size_t and I think on most systems, the c size_t type is unsigned long.

@JKrobNESDIS
Copy link
Author

JKrobNESDIS commented Jul 24, 2024 via email

@EricEngle-NOAA
Copy link
Contributor

@JKrobNESDIS - Sorry that is has taken me a few days to reply here. Last week I was able to get an example Fortran program working related to this issue. Here is the example program.

program test
   use grib_mod
   implicit none

   type(gribfield) :: gfld

   character(len=:), allocatable :: filename

   integer :: j, k, ifile, iret
   integer :: jdisc, jpdtn, jgdtn
   integer :: firstval, lastval
   real :: fldmax, fldmin
   logical :: unpack = .true.

   integer, dimension(200) :: jids, jpdt, jgdt

   j = 0
   k = 0
   ifile = 10
   iret = 0
   filename = "./20231201000000-24h-oper-fc.grib2"

   ! Open GRIB2 file
   call baopenr(ifile, filename, iret)
   print *, "iret from baopenr = ", iret

   ! Set GRIB2 field identification values to search for
   jdisc = -1
   jids(:) = -9999
   jpdtn = -1
   jpdt(:) = -9999
   jgdtn = -1
   jgdt(:) = -9999

   do

      ! Get field from file
      call getgb2(ifile, 0, j, jdisc, jids, jpdtn, jpdt, jgdtn, jgdt, unpack, k, gfld, iret)
      print *, "iret from getgb2 = ", iret
      if (iret .ne. 0) exit
      print *, "j, k = ", j, k

      print *, "data representation template number = ", gfld%idrtnum

      ! Process field ...
      if (unpack) then
         firstval = gfld%fld(1)
         lastval = gfld%fld(gfld%ndpts)
         fldmax = maxval(gfld%fld)
         fldmin = minval(gfld%fld)
         print *, "min, max of data = ", fldmin, fldmax
      end if

      ! Free memory when done with field
      call gf_free(gfld)
      print *, "Free gfld..."

      ! Set j to k in order to advance to next message
      j = k

   end do

   call gf_finalize()
   print *, "Free g2 library memory..."

end program test

The program is using a modified local version of gf_unpack7 which contains an interface to g2c's aecunpackf() defined as:

interface
   function g2c_aecunpackf(cgrib,lensec,idrstmpl,ndpts,fld) bind(c, name="g2c_aecunpackf")
      use iso_c_binding
      character(kind=c_char), intent(in), dimension(*) :: cgrib   ! Equivalent to unsigned char *
      integer(c_int), value, intent(in) :: lensec
      integer(kind=c_int), dimension(*), intent(inout) :: idrstmpl
      integer(kind=c_size_t), value :: ndpts
      real(kind=c_float), dimension(*), intent(inout) :: fld
      integer(kind=c_int) :: g2c_aecunpackf
   end function g2c_aecunpackf
end interface

and is called in gf_unpack7 as like this

ret = g2c_aecunpackf(cgrib(ipos), lensec-5, idrstmpl, int(ndpts,kind=8), fld)

Note that I did not need any modified local copies of any of the AEC C source files from g2c. This program runs successfully on macOS and Linux (Fedora 40).

I am wondering about the following in your situation:

  1. How was g2c and g2 built on your Windows system? MS compilers or Intel?
  2. What Fortran compiler are you using for the program?

@JKrobNESDIS
Copy link
Author

JKrobNESDIS commented Aug 6, 2024 via email

@EricEngle-NOAA
Copy link
Contributor

No problem. There could be an incompatibility between g2c (built with MSVC++) and g2 (Intel Fortran) and your program (Intel Fortran) what does not show itself until program execution.

@edwardhartnett
Copy link
Contributor

@EricEngle-NOAA are you submitting those g2c changes as a PR on g2c?

@EricEngle-NOAA
Copy link
Contributor

@edwardhartnett - As of now there are no changes to g2c needed. The AEC compression is working correctly from within the g2c library. Changes are definitely needed in NCEPLIBS-g2, but we can perhaps discuss that in an issue for g2, likely in NOAA-EMC/NCEPLIBS-g2#705 or NOAA-EMC/NCEPLIBS-g2#458. IMO, this specific issue can be closed.

@github-project-automation github-project-automation bot moved this from To do to Done in NCEPLIBS-g2c-2.0.0 Aug 7, 2024
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
question Further information is requested
Projects
No open projects
Development

No branches or pull requests

3 participants