RE: Shift in URMA QPE06 compared to RFC QPE06 All - Home - RTMA/URMA Discussion Group

I don't know how AWIPS works behind the scenes.

The two plots I sent yesterday had different background map resolutions; also the reply was sent via Thunderbird mail, so if you were not viewing them directly on VLAB interface you might have seen one plot showing up twice.  Re-posting the two PNG files (same background map resolution) directly on VLAB interface for clarity:

1) pcpurma_g188.2018061312.06h.grb2.png: the URMA file (on the NW grid) correctly plotted

2) pcpurma_g188.2018061312.06h.grb2_off_by_half_dxdy.png: the same URMA file incorrectly plotted, assuming that each grid box's coordinate is defined at the "lower-left" corner of the grid box, rather than at the center of the box.

Ying

Ying Lin
NCEP/EMC/Verification, Post-processing and Product Generation Branch

I don't know if any of this will help clarify anything, but NCRFC was involved in various grid shift problems in A2 mpe for a year or two.  The HRAP lat/lon of the xmrg created by mpe represents the southwest corner of the HRAP box.  HRAP is a subset of MDR which is a subset of the LFM grid, and that's how that was defined.  However, the lat/lon point encoded in the grib1 represents the *center* of the grib1 grid box.  There should be a token in your Localization perspective HydroApps->Hydro->Apps_defaults->SITE(yoursiteid) which is:

grib_shift  : 0.5

When mpe runs gribit on the xmrg to create the grib1 file, it will shift the xmrg lat/lon by 0.5 cells to the N and E (which is the center).  The grib1 file is then transmitted on the SBN, ingested by NCEP, etc.

NCRFC constructed our own HRAP overlay shapefile by looping the OHD Fortran code to create the lat/lon coordinates of the southwest corners of all our HRAP boxes.  Our GIS focal then constructed the shape file in the HRAP projection, by connecting the lat/lons.  It's often instructive to load a grib into an independent third-party application, such as NCEI's WCT viewer, which I believe allows one to also load shapefiles.  One can then compare the grib with the HRAP overlay and see if there are improper shifts.  This helps narrow down if it's a local grib encoding problem, AWIPS viewer problem, etc.

-Larry Ellis
NCRFC HAS

[The reply an hour ago was done though Thunderbird email - nearly unreadable.  Here's another try using VLab interface directly.}

Larry,

Thank you very much for weighing in.  I was not aware of this 0.5 cell shift done at RFC - have been assuming that QPEs from RFCs (GRIB1 files sent via SBN) use the lower-left corner lat/lon as co-ord for grid box and do a 0.5 cell shift here.  When did this shifting at the RFCs start (starting with A2 MPE?), and has it been used by all RFCs in creating their QPEs?  And, does the activation of the grib_shift alters the (lat,lon) of the (1,1) point in the grib file (is that the lower-left corner of the (1,1) point, or the center of the (1,1) point?  In the case of NCRFC this is currently
Lat1 37.870000 Long1 -105.515000.

On the subject of AWIPS plotting of precip URMA: I asked our local GRIB2 expert earlier today whether there is a parameter in GRIB2 meta data that can be used to indicate where the lat/lon for a grid box is actually valid at (e.g. lower-left or center of box), thinking that that might be used by AWIPS graphics.  He told me he's not aware of such an option.

Ying
 

[I've been trying to move this thread to the 'QPE' sub-forum, so far no luck yet.  Requested help from VLab Admin - if they can resolve this, the thread will be moved to the 'QPE' sub forum from greater RFC visibility, in case some RFCs have subscribed to the sub-forum but not the main RTMA/URMA Feedback forum.]

Hi Larry,  in re the 1/2 cell shift in A2 MPE you mentioned yesterday: has the shifting been done at NCRFC for some time?  How far back does it go?  I dug up some QPEs we received from NCRFC in Jan 2006, more than 12 years ago, and the (1,1) point has the same lat/lon as it does now.  Wouldn't any grid shifting affect the lat/lon of the (1,1) point, since the grid boxes themselves shouldn't move, it's only how the lat/lons are defined at?  Or did the 1/2 cell shifting in MPE done at NCRFC happen earlier than Jan 2006? 

With help from VLab Admin, this thread is now in the QPE sub-forum of the RTMA/URMA feedback forum, so that VLab notifications for this topic can reach people who have subscribed to the QPE sub-forum but not the main RTMA/URMA forum (main forum subscribers do receive sub-forum posts).

Ying

A history of QPE download formats is given at:
https://water.weather.gov/precip/download.php

Note what's said under "Data Formats Pre-March 22nd, 2017 - Derived from Stage III QPE".  For the shapefile format, it says the latitude and longitude are for the HRAP grid "point".  That's a massive conceptual error...HRAP is not a point, it's an *area*, and the data value represents *mean depth* over the HRAP area, which when combined with the HRAP area denotes the *volume* of precipitation.  For other meteorological parameters, such as pressure, temperature, wind, etc., the lat/lon of a value encoded in shapefile or grib represents an infinitely small sample *point* on a continuous surface.  

For the ESRI shapefile format, GIS software wanted the lat/lon encoded to be the center of the HRAP box, but before fall 2007 it was incorrectly the southwest corner of the HRAP box (same lat/lon as the xmrg format).  They went back and recomputed the shapefile lat/lon coordinates to be the center of the HRAP box.  I'm guessing the 2007 date was probably also about the time they started to look at the encoding of other formats, such as grib1 encoding?  I don't remember dates for certain, but I am certain the grib_shift token was introduced in the AWIPS-1 era, and it's in our AWIPS-1 /awips/hydroapps/.Apps_defaults_site file.  

(Note, in the current download metadata, they mention August 9, 2007, but in 2010, the date was October 16, 2007.  Reference the internet Wayback machine, https://web.archive.org/web/20100513032047/http://water.weather.gov:80/precip/download.php)

I don't remember if we ever sent the native xmrg format via the SBN.  The current AWIPS-2 mpe directory structure has provision to send xmrg's via the SBN, and the directory is being populated, but those data do not appear to be transmitted.

Not meaning to throw Paul Tilles under the bus, but he and/or his group have been the longtime programmers for the gribit software, especially in the AWIPS-1 world.  He might still have some history on changes made to that program, at least in the AWIPS-1 era.

- Larry

I kept forgetting that replying via email (at least through Thunderbird interface) results in an awful-looking posting through VLab.  Here's the post directly through VLAB:

hanks Larry.

I went back further in digging for old QPEs.  In 2001/06/20 (17 years ago), QPEs from NCRFC did have different lat/lon for the (1,1) point:
    Lat1 37.869000 Long1 -104.390000 for 2001/06/20
vs.
     Lat1 37.870000 Long1 -105.515000 (current)
the longitude difference  way more than half a grid cell difference though.  Same array dimension of (450 x 350).   Was there a general shifting  of the NCRFC array westward between 2001 and 2016?

The NWRFC QPE's grids appeared to have stayed the same in the past 17 years.  OHRFC's are different:

    Lat1 36.340000 Long1 -89.896000 in 2001
vs.
   36.188000 Long1 -91.347000
(again, rather large difference in longitudes)

ABRFC's grid stayed the same.

We did start getting GRIB1 files from the RFCs as far back as 2001.   If there is a 1/2 grid cell shift, it would be reflected in the lat/lon of the (1,1) point, wouldn't it?  Does that mean that at present, some RFCs are using the lower-left corner definition, some are using the center-of-grid-box definition?

RFCs:  can you confirm whether the QPEs (GRIB1 files) you sent to NCEP via SBN have the lat/lon co-ordinate defined at the center of the grid box, or at the lower-left corner?  I.e.  the (lat1,lon1) in the grib header - are they the lat/lon at the center of the (1,1) point, or at the lower-left corner of the (1,1) point?  If you know there was shift of 1/2 of a grid cell, and remember when it happened, that would be even better.

Thank you.

Ying

Yes, in regard to MSR's 1,1, grid...when our local HRAP mask was created in the early 1990's, we didn't have any good maps of the Souris River basin in Canada, which we have to model (it flows south into North Dakota then back north into Saskatchewan).  Sometime in the early 2000's, after delineating Souris River basins in IHABBS, we shifted our whole grid west to cover the Souris.  Despite a lot of coordination with OH, the Corps of Engineers and others, that caused problems for a lot of external users.

What's a good Linux utility to see the metadata from a grib1?  I use xxd to look a the first couple lines of the header.  I guess I'd look at a grib from around 2006 compared to today's, there should be a 0.5 degree shift to the N and E.

Hello Ying:

I have been following the post and what Larry mentioned is as much and more than I know about the grib cell centers. We have been using two sets of software to generate the grib files that we disseminate via AWIPS. Off the top of my head, I am not sure whether they have different WMO headers, but when I am back in the office can provide examples of both if you have the means to compare them next week. One set, is converted from xmrg to grib1 using the old gribits program which has been around and un-changed for a long time. The other is via the AWIPS2 MPE software directly. During the early period of the AWIPS2 migration and up until about January of 2017 we had a couple of DRs that had to be addressed with this grid-cell shift in the AWIPS2 MPE software. As far as I can remember though I thought the main the issue was within the MPE display and do not know whether any changes were made to the export of the data in the process of working on the former issue. The last DR associated with this was AWIPS2 Redmine DR 19589.

David Welch

LMRFC

Thanks Larry and David.

A good GRIB1 inventory/decoding utility is 'wgrib': http://www.cpc.ncep.noaa.gov/products/wesley/wgrib.html.  It should work on a Linux.  Here's a "verbose" inventory of one of yesterday's NCRFC QPE:

g14a2 $ wgrib -V QPE.157.2018062015.01h
rec 1:21:date 2018062014 APCP kpds5=61 kpds6=1 kpds7=0 levels=(0,0) grid=255 sfc 0-1hr acc: bitmap: 0 undef
  APCP=Total precipitation [kg/m^2]
  timerange 4 P1 0 P2 1 TimeU 1  nx 450 ny 350 GDS grid 5 num_in_ave 0 missing 0
  center 9 subcenter 157 process 161 Table 2 scan: WE:SN winds(grid)
  polar stereo: Lat1 37.870000 Long1 -105.515000 Orient -105.000000
     north pole (450 x 350) Dx 4762 Dy 4762 scan 64 mode 8
  min/max data 0 20.5  num bits 8  BDS_Ref 0  DecScale 0 BinScale -3

Today I dug out some QPE's from 14 years ago:

tfe09:wgrib -V QPE.009.157 (file from 2004/6/20, at the time each RFC's QPE was piled up into a single big file for the day):
rec 1:21:date 2004062000 APCP kpds5=61 kpds6=1 kpds7=0 levels=(0,0) grid=255 sfc
 0-1hr acc: bitmap: 0 undef
  APCP=Total precipitation [kg/m^2]
  timerange 4 P1 0 P2 1 TimeU 1  nx 450 ny 350 GDS grid 5 num_in_ave 0 missing 0
  center 9 subcenter 157 process 161 Table 2 scan: WE:SN winds(grid)
  polar stereo: Lat1 37.870000 Long1 -105.515000 Orient -105.000000
     north pole (450 x 350) Dx 4762 Dy 4762 scan 64 mode 8
  min/max data 0 2.875  num bits 5  BDS_Ref 0  DecScale 0 BinScale

Note the grid in 2004 appeared to be exactly the same as it is now.  If the 1/2 cell shift in "A2 MPE" happened in 2006, it seems it would have caused a shift in the Lat1/Long1 (lat/lon of the (1,1) grid box)) in the GRIB file.

Also looked at the LMRFC QPEs:

g14a2 $ wgrib -V QPE.154.2018062015.01h

Undefined parameter table (center 9-154 table 128), using NCEP-opn
rec 1:21:date 2018062014 O3TOT kpds5=237 kpds6=1 kpds7=0 levels=(0,0) grid=255 sfc 0-1hr acc: bitmap: 0 undef
  O3TOT=Total ozone [kg/m^2]
  timerange 4 P1 0 P2 1 TimeU 1  nx 419 ny 419 GDS grid 5 num_in_ave 0 missing 0
  center 9 subcenter 154 process 152 Table 128 scan: WE:SN winds(grid)
  polar stereo: Lat1 28.686000 Long1 -98.791000 Orient -105.000000
     north pole (419 x 419) Dx 4762 Dy 4762 scan 64 mode 8
  min/max data 0 39.25  num bits 9  BDS_Ref 0  DecScale 0 BinScale -3

vs.

tfe09:wgrib -V QPE.009.154  (from 2004/06/20)
rec 1:21:date 2004062001 APCP kpds5=61 kpds6=1 kpds7=0 levels=(0,0) grid=255 sfc
 0-1hr acc: bitmap: 0 undef
  APCP=Total precipitation [kg/m^2]
  timerange 4 P1 0 P2 1 TimeU 1  nx 419 ny 419 GDS grid 5 num_in_ave 0 missing 0
  center 9 subcenter 154 process 161 Table 2 scan: WE:SN winds(grid)
  polar stereo: Lat1 28.686000 Long1 -98.792000 Orient -105.000000
     north pole (419 x 419) Dx 4762 Dy 4762 scan 64 mode 8
  min/max data 0 32.75  num bits 9  BDS_Ref 0  DecScale 0 BinScale -3

There is a 0.001 deg difference in Long1, which is about 97m - less than 1/20 of the distance of half of a grid cell shift.  It seems more likely the result of a change in round-off errors when there is, say, a change in the grib encoding program.

Ying

Okay everyone, there appears to definitely be a problem.  Last night I saved an hour of mpe, immediately copied the resultant grib to my work PC and uploaded to Google drive.  I also did the same with our local HRAP shapefile.  I just now loaded these in qgis at home, and the cells in the grib are definitely offset from the HRAP grid.  Visually, it doesn't look like exactly half a grid, but more like 3/4 a grid.  However, the HRAP shapefile I grabbed was reprojected from native HRAP projection to Unprojected WGS84 for use in AWIPS, and I'm not sure I trust that transformation.  But it looks like it's time to go back and examine everything computational in A2mpe.

It'd be interesting grab something from around 2012, we have some old A1 gribit executables from that time frame.  I do know when I would tail the A1mpe log file, every time we saved an hour, it would it explicitly print out it was shifting the grib lat/lon by 0.5 grids to the north and east.  We lost that helpful information with A2.

Larry

Second time around...original message disappeared.  Testing in NCEI's wct shows current A2mpe saved grib's are shifted to the southwest by half a grid.

When wgrib -V <filename> says: "Lat1 37.870000 Long1 -105.515000" what is that supposed to represent?  What it actually represents is the southwest corner of our local HRAP 1,1 (southwest corner of national HRAP grid 390,380).

Attached are a couple simple Fortran programs from the dmip project.  Compile with gfortran. 

In wct, I queried the lat/lon near the center of the northeast quarter of a positive grib cell.  I entered that into "latlon_to_hrap".  Truncating those results to the whole number, I entered those HRAP coordinates into "get_hrap_latlon".  That gave the lat/lon of the southwest corner of the HRAP box.  Mousing to those coordinates in WCT should have placed me at the southwest corner of that grib cell, if the grib were properly aligned with the HRAP grid.  What I got was basically the centroid of the grib cell.

One would expect that result if the grib encoding software were simply giving it the southwest corner of the HRAP and not computing the centroid of the HRAP (i.e. shift the HRAP lat/lon by 0.5 cells to the North and East.

We import these grib's into CHPS, however I can verify these are plotting in the correct location in the CHPS spatial display.  CHPS has an import routine where the grid geometry is specified, and it's probable the shift was corrected when the import configuration was done.

Larry

Hi Larry,

Thank you very much for your perseverance in using VLAB messaging (I lost ~40 minutes of work last week on a disappeared message too, now I compose a message on on a local text editor first) and for testing/verifying the QPE grib file in WCT and CHPS displays.  "wgrib -V" output of "Lat1 37.870000 Long1 -105.515000" is the (lat,lon) of the (1,1) point of the NCRFC's GRIB1 QPE file.  That's why I keep thinking that RFCs still defines each grid box by the lat/lon of the lower-left corner of the box.  If A2mpe saved the GRIB file using "lat/lon defined at center of the grid box" , how is such a shift manifested in the GRIB file, if not in the GRIB header's grid definition?  

Here's a google doc summary of this issue: https://docs.google.com/document/d/1ZaYswOSPCjwQWCmj2-oF_PmOJEvcPiRGMgCxvWuwpC4/edit?usp=sharing. 

I went to some pains to look at the (Lat1,Lon1) of all 12 RFC's QPE over the years to see if there is a change that corresponds to the A2mpe 1/2 cell shift (not seeing it yet).  You can edit the doc (or comment, if you prefer). 

I sent an inquiry to another GRIB expert asking if there is an indicator in GRIB1 that shows where the lat/lon of a grid box is defined - last week I learned there isn't such a thing in GRIB2, and it seemed unlikely it could be in the much-simpler GRIB1.  We'll see what the expert says.  

Ying

I don't know if this muddies the water any.  In our CHPS/FEWS, we import the mpe grib1 files for our QPE.  In the qpe Import definition in our ModuleConfigFiles, we specify the gridStartPoint is the SW corner, which is why our QPE imports correctly.  However, checking the other import xml files, we have gridStartPoint set to SW for *all* the grib's we import, including those from NOHRSC, MRMS, RTMA, forecast grib's from GFE, etc., and also netCDF's from the National Water Model.

Larry

Hi Larry,

I suspect that the GRIB1 RFC QPEs have been using the "southwest/lower-left" corner of grid box lat/lon definition all along (anyone knows otherwise?) but then what was the "A2mpe" 1/2-cell shift about, and if a shift did happen, how was it reflected in the resulting GRIB file, if not in the grid spec?.  What you said late last Friday in re NCEI WCT plot might provide a clue.  I'm rather swamped now but when the flood recedes a bit I'd like to ask the WCT developers (anyone has a contact for WCT?) how they "know" exactly how a grid box is defined, in their plots of RFC QPEs. If you don't mind I'll quote what you said ("Testing in NCEI's wct shows current A2mpe saved grib's are shifted to the southwest by half a grid") when contacting WCT. 

Ying

Today I heard from the expert on whether GRIB1 has a way to specify whether lat/lon of a grid box (e.g. lat1/lon1 for the (1,1) point) is defined at the center or the lower-left corner of the grid box.  The answer is no.

Ying

A ticket was filed with NCF (#220084.  Possible AWIPS2 problem?) on 22 Feb to investigate this issue.  This is the document included with the ticket: https://docs.google.com/document/d/1ihMv-03ZiwLZOlGU472TBMTXfhjvOGh8fUS7l_uw4c4/edit?usp=sharing. ​​​​​

If you would like to add to the document but do not yet have write permission, please drop me a line and I'll add your name to the list of people who can edit it. 

Ying