Welcome, Guest
Username: Password: Remember me
General topics

TOPIC: Domain size from grib parameters

Domain size from grib parameters 9 years 2 days ago #394

I would like to get model domain parameters from grib file directly to my python script. Perfect output would be the one that "gl" gives, like:
Input geometry:

Data comes from: 96 99 99
Projection : 3
Nlon,Nlat,Nlev : 39 39 65
South,West : 58.020000000000003 9.4260000000000002
North,East : 60.209217032207675 18.875895662695324
Dlon,Dlat : 11000.000000000000 11000.000000000000
Projlon,Projlat(1,2) : -10.000000000000000 59.200000000000003 59.200000000000003
but at least at the moment I prefer to use the grib_api directly from python and with grib_api I'm only able to get:
South, West

So, I guess that coordinates North and East are actually not inside the grib file and are computed in gl? Is this correct? Or I just do not know the right parameters for grib_api. For example, for South and West I can use key names like "latitudeOfFirstGridPointInDegrees", but I cannot find anything for the upper right corner.
Also, if there is a simple proj command which can give me the right answer, it would be ok. Any comments appreciated.

Best regards
Andres Luhamaa

Re:Domain size from grib parameters 9 years 2 days ago #395

  • Ulf Andrae
  • Ulf Andrae's Avatar
  • Administrator
  • Posts: 288
  • Thank you received: 32

You are right that the corner information is not always in the grib file, thus gl calculates it. In gl_grib_api we extract the following geometrical properties in fill_geo.f90:
   IF ( scani_pos == 1 ) THEN
       CALL grib_get(igr,'longitudeOfFirstGridPointInDegrees',geo%west )
       CALL grib_get(igr,'longitudeOfLastGridPointInDegrees' ,geo%east,ierr)
       IF ( ierr /= 0 ) geo%east = 0.0_jprb
       CALL grib_get(igr,'longitudeOfFirstGridPointInDegrees',geo%east )
       CALL grib_get(igr,'longitudeOfLastGridPointInDegrees' ,geo%west )
    IF ( scanj_pos == 1 ) THEN
       CALL grib_get(igr,'latitudeOfFirstGridPointInDegrees' ,geo%south)
       CALL grib_get(igr,'latitudeOfLastGridPointInDegrees'  ,geo%north,ierr)
       IF ( ierr /= 0 ) geo%north = 0.0_jprb
       CALL grib_get(igr,'latitudeOfFirstGridPointInDegrees' ,geo%north)
       CALL grib_get(igr,'latitudeOfLastGridPointInDegrees'  ,geo%south)

As you can see we try to extract all corners, but they are not always returned. In the case they are not available we calculate the corners in check_geo.f90

I hope this helps.


Re:Domain size from grib parameters 9 years 2 days ago #396

I think it is clear now, thanks a lot.
Time to create page: 0.069 seconds