Closed kevmcgrath closed 1 year ago
The min()
and max()
methods for NumPy arrays do not ignore nans -- similar to np.amin()
and np.amax()
. To ignore the nans, use np.nanmin()
and np.nanmax()
and you will see the correct min and max values.
Ahhh, thanks for pointing me to those NumPy commands. These are now showing the proper values.
I have additional information that will clarify the issue that I'm seeing. While I'm able to read all URMA precip GRIBs with grib2io and the values look good, the problem is that packing and writing to a new GRIB results in a file with all 0 values for some specific input files.
Note that URMA precip GRIBs get updated multiple times after their initial release. I'm finding that grib2io is able to properly pack and write to GRIB when the input data is from the last 24 hours. It's the URMA precip GRIBs that have valid times older than 24 hours that result in grib2io writing out all 0 values.
One possible smoking gun is that initially, the URMA precip GRIBs use this packing:
Section 5: dataRepresentationTemplateNumber = 2 - Grid Point Data - Complex Packing (see Template 5.2)
...but after about a day, it changes to this:
Section 5: dataRepresentationTemplateNumber = 3 - Grid Point Data - Complex Packing and Spatial Differencing (see Template 5.3)
Perhaps there are other differences that don't play nicely with grib2io. Note that pygrib isn't impacted by these differences.
Below I'll include messages for "good" and "bad" URMA precip GRIBs.
"Good" GRIB (grib2io is able to read the data, pack it, and write out the true values to a new GRIB) - AT THE TIME OF THIS WRITING: /lfs/h1/ops/prod/com/urma/v2.10/pcpurma.20230426/pcpurma_wexp.2023042612.06h.grb2
Section 0: discipline = 0 - Meteorological Products Section 1: originatingCenter = 7 - US National Weather Service - NCEP (WMC) Section 1: originatingSubCenter = 4 - Environmental Modeling Center Section 1: masterTableInfo = 2 - Version Implemented on 4 November 2003 Section 1: localTableInfo = 1 - Number of local table version used. Section 1: significanceOfReferenceTime = 0 - Analysis Section 1: year = 2023 Section 1: month = 4 Section 1: day = 26 Section 1: hour = 6 Section 1: minute = 0 Section 1: second = 0 Section 1: refDate = 2023-04-26 06:00:00 Section 1: productionStatus = 0 - Operational Products Section 1: typeOfData = 0 - Analysis Products Section 3: interpretationOfListOfNumbers = 0 - There is no appended list Section 3: gridDefinitionTemplateNumber = 30 - Lambert Conformal (Can be Secant, Tangent, Conical, or Bipolar) Section 3: shapeOfEarth = 1 - Earth assumed spherical with radius specified (in m) by data producer Section 3: earthRadius = 6371200.0 Section 3: earthMajorAxis = None Section 3: earthMinorAxis = None Section 3: resolutionAndComponentFlags = [0, 0, 1, 1, 1, 0, 0, 0] Section 3: ny = 1597 Section 3: nx = 2345 Section 3: scanModeFlags = [0, 1, 0, 0, 0, 0, 0, 0] Section 3: gridOrientation = 265.0 Section 3: gridlengthXDirection = 2539.703 Section 3: gridlengthYDirection = 2539.703 Section 3: latitudeFirstGridpoint = 19.228976 Section 3: latitudeSouthernPole = -90.0 Section 3: latitudeTrueScale = 25.0 Section 3: longitudeFirstGridpoint = 233.723448 Section 3: longitudeSouthernPole = 0.0 Section 3: projParameters = {'a': 6371200.0, 'b': 6371200.0, 'proj': 'lcc', 'lat_1': 25.0, 'lat_2': 25.0, 'lat_0': 25.0, 'lon_0': 265.0} Section 3: projectionCenterFlag = 0 Section 3: standardLatitude1 = 25.0 Section 3: standardLatitude2 = 25.0 Section 4: productDefinitionTemplateNumber = 8 - Average, accumulation, extreme values or other statistically processed values at a horizontal level or in a horizontal layer in a continuous or non-continuous time interval. (see Template 4.8) Section 4: fullName = Total Precipitation Section 4: units = kg m-2 Section 4: shortName = APCP Section 4: leadTime = 6:00:00 Section 4: unitOfFirstFixedSurface = unknown Section 4: valueOfFirstFixedSurface = 0.0 Section 4: unitOfSecondFixedSurface = None Section 4: valueOfSecondFixedSurface = 0.0 Section 4: validDate = 2023-04-26 12:00:00 Section 4: duration = 6:00:00 Section 4: level = surface Section 4: backgroundGeneratingProcessIdentifier = 0 Section 4: dayOfEndOfTimePeriod = 26 Section 4: forecastTime = 0 Section 4: generatingProcess = 118 - UnRestricted Mesoscale Analysis (URMA) Section 4: hourOfEndOfTimePeriod = 12 Section 4: hoursAfterDataCutoff = 0 Section 4: minuteOfEndOfTimePeriod = 0 Section 4: minutesAfterDataCutoff = 0 Section 4: monthOfEndOfTimePeriod = 4 Section 4: numberOfMissingValues = 255 Section 4: numberOfTimeRanges = 1 Section 4: parameterCategory = 1 Section 4: parameterNumber = 8 Section 4: scaleFactorOfFirstFixedSurface = 0 Section 4: scaleFactorOfSecondFixedSurface = 0 Section 4: scaledValueOfFirstFixedSurface = 0 Section 4: scaledValueOfSecondFixedSurface = 0 Section 4: secondOfEndOfTimePeriod = 0 Section 4: statisticalProcess = 1 - Accumulation Section 4: timeIncrementOfSuccessiveFields = 0 Section 4: timeRangeOfStatisticalProcess = 6 Section 4: typeOfFirstFixedSurface = 1 - ['Ground or Water Surface', 'unknown'] Section 4: typeOfGeneratingProcess = 2 - Forecast Section 4: typeOfSecondFixedSurface = 255 - ['Missing', 'unknown'] Section 4: typeOfTimeIncrementOfStatisticalProcess = 2 - Successive times processed have same start time of forecast, forecast time is incremented. Section 4: unitOfTimeRange = 1 - Hour Section 4: unitOfTimeRangeOfStatisticalProcess = 1 - Hour Section 4: unitOfTimeRangeOfSuccessiveFields = 255 - Missing Section 4: yearOfEndOfTimePeriod = 2023 Section 5: dataRepresentationTemplateNumber = 2 - Grid Point Data - Complex Packing (see Template 5.2) Section 5: numberOfPackedValues = 1644513 Section 5: typeOfValues = 0 - Floating Point Section 5: binScaleFactor = 4 Section 5: decScaleFactor = 4 Section 5: groupLengthIncrement = 1 Section 5: groupSplittingMethod = 1 - General Group Splitting Section 5: lengthOfLastGroup = 257 Section 5: nBitsGroupWidth = 4 Section 5: nBitsPacking = 15 Section 5: nBitsScaledGroupLength = 10 Section 5: nGroups = 15292 Section 5: priMissingValue = None Section 5: refGroupLength = 1 Section 5: refGroupWidth = 0 Section 5: refValue = 0.0 Section 5: secMissingValue = None Section 5: typeOfMissingValueManagement = 0 - No explicit missing values included within the data values Section 6: bitMapFlag = 0 - A bit map applies to this product and is specified in this section.
"Bad" GRIB (grib2io is able to read the data, pack it, but writing the data out to GRIB results in all 0 values): /lfs/h1/ops/prod/com/urma/v2.10/pcpurma.20230425/pcpurma_wexp.2023042512.06h.grb2
Section 0: discipline = 0 - Meteorological Products Section 1: originatingCenter = 7 - US National Weather Service - NCEP (WMC) Section 1: originatingSubCenter = 4 - Environmental Modeling Center Section 1: masterTableInfo = 1 - Version Implemented on 7 November 2001 Section 1: localTableInfo = 0 - Local tables not used. Only table entries and templates from the current master table are valid. Section 1: significanceOfReferenceTime = 0 - Analysis Section 1: year = 2023 Section 1: month = 4 Section 1: day = 25 Section 1: hour = 6 Section 1: minute = 0 Section 1: second = 0 Section 1: refDate = 2023-04-25 06:00:00 Section 1: productionStatus = 0 - Operational Products Section 1: typeOfData = 0 - Analysis Products Section 3: interpretationOfListOfNumbers = 0 - There is no appended list Section 3: gridDefinitionTemplateNumber = 30 - Lambert Conformal (Can be Secant, Tangent, Conical, or Bipolar) Section 3: shapeOfEarth = 1 - Earth assumed spherical with radius specified (in m) by data producer Section 3: earthRadius = 6371200.0 Section 3: earthMajorAxis = None Section 3: earthMinorAxis = None Section 3: resolutionAndComponentFlags = [0, 0, 1, 1, 1, 0, 0, 0] Section 3: ny = 1597 Section 3: nx = 2345 Section 3: scanModeFlags = [0, 1, 0, 0, 0, 0, 0, 0] Section 3: gridOrientation = 265.0 Section 3: gridlengthXDirection = 2539.703 Section 3: gridlengthYDirection = 2539.703 Section 3: latitudeFirstGridpoint = 19.228976 Section 3: latitudeSouthernPole = -90.0 Section 3: latitudeTrueScale = 25.0 Section 3: longitudeFirstGridpoint = 233.723448 Section 3: longitudeSouthernPole = 0.0 Section 3: projParameters = {'a': 6371200.0, 'b': 6371200.0, 'proj': 'lcc', 'lat_1': 25.0, 'lat_2': 25.0, 'lat_0': 25.0, 'lon_0': 265.0} Section 3: projectionCenterFlag = 0 Section 3: standardLatitude1 = 25.0 Section 3: standardLatitude2 = 25.0 Section 4: productDefinitionTemplateNumber = 8 - Average, accumulation, extreme values or other statistically processed values at a horizontal level or in a horizontal layer in a continuous or non-continuous time interval. (see Template 4.8) Section 4: fullName = Total Precipitation Section 4: units = kg m-2 Section 4: shortName = APCP Section 4: leadTime = 6:00:00 Section 4: unitOfFirstFixedSurface = unknown Section 4: valueOfFirstFixedSurface = 2550.0 Section 4: unitOfSecondFixedSurface = None Section 4: valueOfSecondFixedSurface = 2550.0 Section 4: validDate = 2023-04-25 12:00:00 Section 4: duration = 6:00:00 Section 4: level = surface Section 4: backgroundGeneratingProcessIdentifier = 0 Section 4: dayOfEndOfTimePeriod = 25 Section 4: forecastTime = 0 Section 4: generatingProcess = 118 - UnRestricted Mesoscale Analysis (URMA) Section 4: hourOfEndOfTimePeriod = 12 Section 4: hoursAfterDataCutoff = 0 Section 4: minuteOfEndOfTimePeriod = 0 Section 4: minutesAfterDataCutoff = 0 Section 4: monthOfEndOfTimePeriod = 4 Section 4: numberOfMissingValues = 0 Section 4: numberOfTimeRanges = 1 Section 4: parameterCategory = 1 Section 4: parameterNumber = 8 Section 4: scaleFactorOfFirstFixedSurface = -1 Section 4: scaleFactorOfSecondFixedSurface = -127 Section 4: scaledValueOfFirstFixedSurface = 255 Section 4: scaledValueOfSecondFixedSurface = -2147483647 Section 4: secondOfEndOfTimePeriod = 0 Section 4: statisticalProcess = 1 - Accumulation Section 4: timeIncrementOfSuccessiveFields = 0 Section 4: timeRangeOfStatisticalProcess = 6 Section 4: typeOfFirstFixedSurface = 1 - ['Ground or Water Surface', 'unknown'] Section 4: typeOfGeneratingProcess = 2 - Forecast Section 4: typeOfSecondFixedSurface = 255 - ['Missing', 'unknown'] Section 4: typeOfTimeIncrementOfStatisticalProcess = 2 - Successive times processed have same start time of forecast, forecast time is incremented. Section 4: unitOfTimeRange = 1 - Hour Section 4: unitOfTimeRangeOfStatisticalProcess = 1 - Hour Section 4: unitOfTimeRangeOfSuccessiveFields = 255 - Missing Section 4: yearOfEndOfTimePeriod = 2023 Section 5: dataRepresentationTemplateNumber = 3 - Grid Point Data - Complex Packing and Spatial Differencing (see Template 5.3) Section 5: numberOfPackedValues = 3744965 Section 5: typeOfValues = 0 - Floating Point Section 5: binScaleFactor = -5 Section 5: decScaleFactor = 0 Section 5: groupLengthIncrement = 1 Section 5: groupSplittingMethod = 1 - General Group Splitting Section 5: lengthOfLastGroup = 167 Section 5: nBitsGroupWidth = 4 Section 5: nBitsPacking = 11 Section 5: nBitsScaledGroupLength = 10 Section 5: nBytesSpatialDifference = 2 Section 5: nGroups = 43341 Section 5: priMissingValue = 9.999000260554009e+20 Section 5: refGroupLength = 1 Section 5: refGroupWidth = 0 Section 5: refValue = 0.0 Section 5: secMissingValue = nan Section 5: spatialDifferenceOrder = 1 - First-Order Spatial Differencing Section 5: typeOfMissingValueManagement = 1 - Primary missing values included within the data values Section 6: bitMapFlag = 255 - A bit map does not apply to this product.
What is your process for writing out the "bad" GRIB message?
Here are examples of reading/writing/reading both a "good" and a "bad" URMA precip GRIB.
"Good" GRIB (at the time of this writing): /lfs/h1/ops/prod/com/urma/v2.10/pcpurma.20230427/pcpurma_wexp.2023042706.06h.grb2
---------------------------- * Examine source GRIB $ ll /lfs/h1/ops/prod/com/urma/v2.10/pcpurma.20230427/pcpurma_wexp.2023042706.06h.grb2 -rw-rw-r-- 1 ops.prod prod 881599 Apr 27 12:02 /lfs/h1/ops/prod/com/urma/v2.10/pcpurma.20230427/pcpurma_wexp.2023042706.06h.grb2 $ gdalinfo /lfs/h1/ops/prod/com/urma/v2.10/pcpurma.20230427/pcpurma_wexp.2023042706.06h.grb2 -mm | grep Min Computed Min/Max=0.000,114.899 $ degrib /lfs/h1/ops/prod/com/urma/v2.10/pcpurma.20230427/pcpurma_wexp.2023042706.06h.grb2 -C -stdout | grep Packing Info | Packing that was used | 2 (Grid point data - complex packing) ---------------------------- * Open source GRIB with grib2io and write to new GRIB >>> import grib2io >>> import numpy as np >>> sourceGrib='/lfs/h1/ops/prod/com/urma/v2.10/pcpurma.20230427/pcpurma_wexp.2023042706.06h.grb2' >>> grbs = grib2io.open(sourceGrib) >>> grbs mode = rb name = /lfs/h1/ops/prod/com/urma/v2.10/pcpurma.20230427/pcpurma_wexp.2023042706.06h.grb2 messages = 1 current_message = 0 size = 881599 closed = False variables = ('APCP',) levels = ('surface',) >>> message = grbs[0] >>> message Section 0: discipline = 0 - Meteorological Products Section 1: originatingCenter = 7 - US National Weather Service - NCEP (WMC) Section 1: originatingSubCenter = 4 - Environmental Modeling Center Section 1: masterTableInfo = 2 - Version Implemented on 4 November 2003 Section 1: localTableInfo = 1 - Number of local table version used. Section 1: significanceOfReferenceTime = 0 - Analysis Section 1: year = 2023 Section 1: month = 4 Section 1: day = 27 Section 1: hour = 0 Section 1: minute = 0 Section 1: second = 0 Section 1: refDate = 2023-04-27 00:00:00 Section 1: productionStatus = 0 - Operational Products Section 1: typeOfData = 0 - Analysis Products Section 3: interpretationOfListOfNumbers = 0 - There is no appended list Section 3: gridDefinitionTemplateNumber = 30 - Lambert Conformal (Can be Secant, Tangent, Conical, or Bipolar) Section 3: shapeOfEarth = 1 - Earth assumed spherical with radius specified (in m) by data producer Section 3: earthRadius = 6371200.0 Section 3: earthMajorAxis = None Section 3: earthMinorAxis = None Section 3: resolutionAndComponentFlags = [0, 0, 1, 1, 1, 0, 0, 0] Section 3: ny = 1597 Section 3: nx = 2345 Section 3: scanModeFlags = [0, 1, 0, 0, 0, 0, 0, 0] Section 3: gridOrientation = 265.0 Section 3: gridlengthXDirection = 2539.703 Section 3: gridlengthYDirection = 2539.703 Section 3: latitudeFirstGridpoint = 19.228976 Section 3: latitudeSouthernPole = -90.0 Section 3: latitudeTrueScale = 25.0 Section 3: longitudeFirstGridpoint = 233.723448 Section 3: longitudeSouthernPole = 0.0 Section 3: projParameters = {'a': 6371200.0, 'b': 6371200.0, 'proj': 'lcc', 'lat_1': 25.0, 'lat_2': 25.0, 'lat_0': 25.0, 'lon_0': 265.0} Section 3: projectionCenterFlag = 0 Section 3: standardLatitude1 = 25.0 Section 3: standardLatitude2 = 25.0 Section 4: productDefinitionTemplateNumber = 8 - Average, accumulation, extreme values or other statistically processed values at a horizontal level or in a horizontal layer in a continuous or non-continuous time interval. (see Template 4.8) Section 4: fullName = Total Precipitation Section 4: units = kg m-2 Section 4: shortName = APCP Section 4: leadTime = 6:00:00 Section 4: unitOfFirstFixedSurface = unknown Section 4: valueOfFirstFixedSurface = 0.0 Section 4: unitOfSecondFixedSurface = None Section 4: valueOfSecondFixedSurface = 0.0 Section 4: validDate = 2023-04-27 06:00:00 Section 4: duration = 6:00:00 Section 4: level = surface Section 4: backgroundGeneratingProcessIdentifier = 0 Section 4: dayOfEndOfTimePeriod = 27 Section 4: forecastTime = 0 Section 4: generatingProcess = 118 - UnRestricted Mesoscale Analysis (URMA) Section 4: hourOfEndOfTimePeriod = 6 Section 4: hoursAfterDataCutoff = 0 Section 4: minuteOfEndOfTimePeriod = 0 Section 4: minutesAfterDataCutoff = 0 Section 4: monthOfEndOfTimePeriod = 4 Section 4: numberOfMissingValues = 255 Section 4: numberOfTimeRanges = 1 Section 4: parameterCategory = 1 Section 4: parameterNumber = 8 Section 4: scaleFactorOfFirstFixedSurface = 0 Section 4: scaleFactorOfSecondFixedSurface = 0 Section 4: scaledValueOfFirstFixedSurface = 0 Section 4: scaledValueOfSecondFixedSurface = 0 Section 4: secondOfEndOfTimePeriod = 0 Section 4: statisticalProcess = 1 - Accumulation Section 4: timeIncrementOfSuccessiveFields = 0 Section 4: timeRangeOfStatisticalProcess = 6 Section 4: typeOfFirstFixedSurface = 1 - ['Ground or Water Surface', 'unknown'] Section 4: typeOfGeneratingProcess = 2 - Forecast Section 4: typeOfSecondFixedSurface = 255 - ['Missing', 'unknown'] Section 4: typeOfTimeIncrementOfStatisticalProcess = 2 - Successive times processed have same start time of forecast, forecast time is incremented. Section 4: unitOfTimeRange = 1 - Hour Section 4: unitOfTimeRangeOfStatisticalProcess = 1 - Hour Section 4: unitOfTimeRangeOfSuccessiveFields = 255 - Missing Section 4: yearOfEndOfTimePeriod = 2023 Section 5: dataRepresentationTemplateNumber = 2 - Grid Point Data - Complex Packing (see Template 5.2) Section 5: numberOfPackedValues = 1383199 Section 5: typeOfValues = 0 - Floating Point Section 5: binScaleFactor = 5 Section 5: decScaleFactor = 4 Section 5: groupLengthIncrement = 1 Section 5: groupSplittingMethod = 1 - General Group Splitting Section 5: lengthOfLastGroup = 212 Section 5: nBitsGroupWidth = 4 Section 5: nBitsPacking = 15 Section 5: nBitsScaledGroupLength = 10 Section 5: nGroups = 18182 Section 5: priMissingValue = None Section 5: refGroupLength = 1 Section 5: refGroupWidth = 0 Section 5: refValue = 0.0 Section 5: secMissingValue = None Section 5: typeOfMissingValueManagement = 0 - No explicit missing values included within the data values Section 6: bitMapFlag = 0 - A bit map applies to this product and is specified in this section. >>> np.nanmin(message.data) 0.0 >>> np.nanmax(message.data) 114.8992 >>> message.pack() >>> grbout = grib2io.open('/u/mdl.wsup/scratch/grib1.out', mode='w') >>> grbout.write(message) >>> grbout.close() --------------------------- * Examine new GRIB from the CLI (looks good) $ ll /u/mdl.wsup/scratch/grib1.out -rw-r--r-- 1 mdl.wsup mdl 895549 Apr 27 12:59 /u/mdl.wsup/scratch/grib1.out $ gdalinfo /u/mdl.wsup/scratch/grib1.out -mm|grep Min Computed Min/Max=0.000,114.899 $ degrib /u/mdl.wsup/scratch/grib1.out -C -stdout | grep Packing Info | Packing that was used | 2 (Grid point data - complex packing) --------------------------- * Examine new GRIB with grib2io (looks good) >>> sourceGrib='/u/mdl.wsup/scratch/grib1.out' >>> grbs = grib2io.open(sourceGrib) >>> grbs mode = rb name = /u/mdl.wsup/scratch/grib1.out messages = 1 current_message = 0 size = 895549 closed = False variables = ('APCP',) levels = ('surface',) >>> grbs[0] Section 0: discipline = 0 - Meteorological Products Section 1: originatingCenter = 7 - US National Weather Service - NCEP (WMC) Section 1: originatingSubCenter = 4 - Environmental Modeling Center Section 1: masterTableInfo = 2 - Version Implemented on 4 November 2003 Section 1: localTableInfo = 1 - Number of local table version used. Section 1: significanceOfReferenceTime = 0 - Analysis Section 1: year = 2023 Section 1: month = 4 Section 1: day = 27 Section 1: hour = 0 Section 1: minute = 0 Section 1: second = 0 Section 1: refDate = 2023-04-27 00:00:00 Section 1: productionStatus = 0 - Operational Products Section 1: typeOfData = 0 - Analysis Products Section 3: interpretationOfListOfNumbers = 0 - There is no appended list Section 3: gridDefinitionTemplateNumber = 30 - Lambert Conformal (Can be Secant, Tangent, Conical, or Bipolar) Section 3: shapeOfEarth = 1 - Earth assumed spherical with radius specified (in m) by data producer Section 3: earthRadius = 6371200.0 Section 3: earthMajorAxis = None Section 3: earthMinorAxis = None Section 3: resolutionAndComponentFlags = [0, 0, 1, 1, 1, 0, 0, 0] Section 3: ny = 1597 Section 3: nx = 2345 Section 3: scanModeFlags = [0, 1, 0, 0, 0, 0, 0, 0] Section 3: gridOrientation = 265.0 Section 3: gridlengthXDirection = 2539.703 Section 3: gridlengthYDirection = 2539.703 Section 3: latitudeFirstGridpoint = 19.228976 Section 3: latitudeSouthernPole = -90.0 Section 3: latitudeTrueScale = 25.0 Section 3: longitudeFirstGridpoint = 233.723448 Section 3: longitudeSouthernPole = 0.0 Section 3: projParameters = {'a': 6371200.0, 'b': 6371200.0, 'proj': 'lcc', 'lat_1': 25.0, 'lat_2': 25.0, 'lat_0': 25.0, 'lon_0': 265.0} Section 3: projectionCenterFlag = 0 Section 3: standardLatitude1 = 25.0 Section 3: standardLatitude2 = 25.0 Section 4: productDefinitionTemplateNumber = 8 - Average, accumulation, extreme values or other statistically processed values at a horizontal level or in a horizontal layer in a continuous or non-continuous time interval. (see Template 4.8) Section 4: fullName = Total Precipitation Section 4: units = kg m-2 Section 4: shortName = APCP Section 4: leadTime = 6:00:00 Section 4: unitOfFirstFixedSurface = unknown Section 4: valueOfFirstFixedSurface = 0.0 Section 4: unitOfSecondFixedSurface = None Section 4: valueOfSecondFixedSurface = 0.0 Section 4: validDate = 2023-04-27 06:00:00 Section 4: duration = 6:00:00 Section 4: level = surface Section 4: backgroundGeneratingProcessIdentifier = 0 Section 4: dayOfEndOfTimePeriod = 27 Section 4: forecastTime = 0 Section 4: generatingProcess = 118 - UnRestricted Mesoscale Analysis (URMA) Section 4: hourOfEndOfTimePeriod = 6 Section 4: hoursAfterDataCutoff = 0 Section 4: minuteOfEndOfTimePeriod = 0 Section 4: minutesAfterDataCutoff = 0 Section 4: monthOfEndOfTimePeriod = 4 Section 4: numberOfMissingValues = 255 Section 4: numberOfTimeRanges = 1 Section 4: parameterCategory = 1 Section 4: parameterNumber = 8 Section 4: scaleFactorOfFirstFixedSurface = 0 Section 4: scaleFactorOfSecondFixedSurface = 0 Section 4: scaledValueOfFirstFixedSurface = 0 Section 4: scaledValueOfSecondFixedSurface = 0 Section 4: secondOfEndOfTimePeriod = 0 Section 4: statisticalProcess = 1 - Accumulation Section 4: timeIncrementOfSuccessiveFields = 0 Section 4: timeRangeOfStatisticalProcess = 6 Section 4: typeOfFirstFixedSurface = 1 - ['Ground or Water Surface', 'unknown'] Section 4: typeOfGeneratingProcess = 2 - Forecast Section 4: typeOfSecondFixedSurface = 255 - ['Missing', 'unknown'] Section 4: typeOfTimeIncrementOfStatisticalProcess = 2 - Successive times processed have same start time of forecast, forecast time is incremented. Section 4: unitOfTimeRange = 1 - Hour Section 4: unitOfTimeRangeOfStatisticalProcess = 1 - Hour Section 4: unitOfTimeRangeOfSuccessiveFields = 255 - Missing Section 4: yearOfEndOfTimePeriod = 2023 Section 5: dataRepresentationTemplateNumber = 2 - Grid Point Data - Complex Packing (see Template 5.2) Section 5: numberOfPackedValues = 1383199 Section 5: typeOfValues = 0 - Floating Point Section 5: binScaleFactor = 5 Section 5: decScaleFactor = 4 Section 5: groupLengthIncrement = 1 Section 5: groupSplittingMethod = 1 - General Group Splitting Section 5: lengthOfLastGroup = 212 Section 5: nBitsGroupWidth = 4 Section 5: nBitsPacking = 15 Section 5: nBitsScaledGroupLength = 9 Section 5: nGroups = 32622 Section 5: priMissingValue = None Section 5: refGroupLength = 1 Section 5: refGroupWidth = 0 Section 5: refValue = 0.0 Section 5: secMissingValue = None Section 5: typeOfMissingValueManagement = 0 - No explicit missing values included within the data values Section 6: bitMapFlag = 0 - A bit map applies to this product and is specified in this section. >>> np.nanmin(grbs[0].data) 0.0 >>> np.nanmax(grbs[0].data) 114.8992
"Bad" GRIB (grib2io is able to read the data, pack it, but writing the data out to GRIB results in a very small file containing values of 0): /lfs/h1/ops/prod/com/urma/v2.10/pcpurma.20230426/pcpurma_wexp.2023042600.06h.grb2
---------------------------- * Examine source GRIB $ ll /lfs/h1/ops/prod/com/urma/v2.10/pcpurma.20230426/pcpurma_wexp.2023042600.06h.grb2 -rw-rw-r-- 1 ops.prod prod 541324 Apr 27 03:04 /lfs/h1/ops/prod/com/urma/v2.10/pcpurma.20230426/pcpurma_wexp.2023042600.06h.grb2 $ gdalinfo /lfs/h1/ops/prod/com/urma/v2.10/pcpurma.20230426/pcpurma_wexp.2023042600.06h.grb2 -mm|grep Min Computed Min/Max=0.000,81.188 $ degrib /lfs/h1/ops/prod/com/urma/v2.10/pcpurma.20230426/pcpurma_wexp.2023042600.06h.grb2 -C -stdout | grep Packing Info | Packing that was used | 3 (Grid point data - complex packing and spatial differencing) ---------------------------- * Open source GRIB with grib2io and write to new GRIB >>> sourceGrib='/lfs/h1/ops/prod/com/urma/v2.10/pcpurma.20230426/pcpurma_wexp.2023042600.06h.grb2' >>> grbs = grib2io.open(sourceGrib) >>> grbs mode = rb name = /lfs/h1/ops/prod/com/urma/v2.10/pcpurma.20230426/pcpurma_wexp.2023042600.06h.grb2 messages = 1 current_message = 0 size = 541324 closed = False variables = ('APCP',) levels = ('surface',) >>> message = grbs[0] >>> message Section 0: discipline = 0 - Meteorological Products Section 1: originatingCenter = 7 - US National Weather Service - NCEP (WMC) Section 1: originatingSubCenter = 4 - Environmental Modeling Center Section 1: masterTableInfo = 1 - Version Implemented on 7 November 2001 Section 1: localTableInfo = 0 - Local tables not used. Only table entries and templates from the current master table are valid. Section 1: significanceOfReferenceTime = 0 - Analysis Section 1: year = 2023 Section 1: month = 4 Section 1: day = 25 Section 1: hour = 18 Section 1: minute = 0 Section 1: second = 0 Section 1: refDate = 2023-04-25 18:00:00 Section 1: productionStatus = 0 - Operational Products Section 1: typeOfData = 0 - Analysis Products Section 3: interpretationOfListOfNumbers = 0 - There is no appended list Section 3: gridDefinitionTemplateNumber = 30 - Lambert Conformal (Can be Secant, Tangent, Conical, or Bipolar) Section 3: shapeOfEarth = 1 - Earth assumed spherical with radius specified (in m) by data producer Section 3: earthRadius = 6371200.0 Section 3: earthMajorAxis = None Section 3: earthMinorAxis = None Section 3: resolutionAndComponentFlags = [0, 0, 1, 1, 1, 0, 0, 0] Section 3: ny = 1597 Section 3: nx = 2345 Section 3: scanModeFlags = [0, 1, 0, 0, 0, 0, 0, 0] Section 3: gridOrientation = 265.0 Section 3: gridlengthXDirection = 2539.703 Section 3: gridlengthYDirection = 2539.703 Section 3: latitudeFirstGridpoint = 19.228976 Section 3: latitudeSouthernPole = -90.0 Section 3: latitudeTrueScale = 25.0 Section 3: longitudeFirstGridpoint = 233.723448 Section 3: longitudeSouthernPole = 0.0 Section 3: projParameters = {'a': 6371200.0, 'b': 6371200.0, 'proj': 'lcc', 'lat_1': 25.0, 'lat_2': 25.0, 'lat_0': 25.0, 'lon_0': 265.0} Section 3: projectionCenterFlag = 0 Section 3: standardLatitude1 = 25.0 Section 3: standardLatitude2 = 25.0 Section 4: productDefinitionTemplateNumber = 8 - Average, accumulation, extreme values or other statistically processed values at a horizontal level or in a horizontal layer in a continuous or non-continuous time interval. (see Template 4.8) Section 4: fullName = Total Precipitation Section 4: units = kg m-2 Section 4: shortName = APCP Section 4: leadTime = 6:00:00 Section 4: unitOfFirstFixedSurface = unknown Section 4: valueOfFirstFixedSurface = 2550.0 Section 4: unitOfSecondFixedSurface = None Section 4: valueOfSecondFixedSurface = 2550.0 Section 4: validDate = 2023-04-26 00:00:00 Section 4: duration = 6:00:00 Section 4: level = surface Section 4: backgroundGeneratingProcessIdentifier = 0 Section 4: dayOfEndOfTimePeriod = 26 Section 4: forecastTime = 0 Section 4: generatingProcess = 118 - UnRestricted Mesoscale Analysis (URMA) Section 4: hourOfEndOfTimePeriod = 0 Section 4: hoursAfterDataCutoff = 0 Section 4: minuteOfEndOfTimePeriod = 0 Section 4: minutesAfterDataCutoff = 0 Section 4: monthOfEndOfTimePeriod = 4 Section 4: numberOfMissingValues = 0 Section 4: numberOfTimeRanges = 1 Section 4: parameterCategory = 1 Section 4: parameterNumber = 8 Section 4: scaleFactorOfFirstFixedSurface = -1 Section 4: scaleFactorOfSecondFixedSurface = -127 Section 4: scaledValueOfFirstFixedSurface = 255 Section 4: scaledValueOfSecondFixedSurface = -2147483647 Section 4: secondOfEndOfTimePeriod = 0 Section 4: statisticalProcess = 1 - Accumulation Section 4: timeIncrementOfSuccessiveFields = 0 Section 4: timeRangeOfStatisticalProcess = 6 Section 4: typeOfFirstFixedSurface = 1 - ['Ground or Water Surface', 'unknown'] Section 4: typeOfGeneratingProcess = 2 - Forecast Section 4: typeOfSecondFixedSurface = 255 - ['Missing', 'unknown'] Section 4: typeOfTimeIncrementOfStatisticalProcess = 2 - Successive times processed have same start time of forecast, forecast time is incremented. Section 4: unitOfTimeRange = 1 - Hour Section 4: unitOfTimeRangeOfStatisticalProcess = 1 - Hour Section 4: unitOfTimeRangeOfSuccessiveFields = 255 - Missing Section 4: yearOfEndOfTimePeriod = 2023 Section 5: dataRepresentationTemplateNumber = 3 - Grid Point Data - Complex Packing and Spatial Differencing (see Template 5.3) Section 5: numberOfPackedValues = 3744965 Section 5: typeOfValues = 0 - Floating Point Section 5: binScaleFactor = -5 Section 5: decScaleFactor = 0 Section 5: groupLengthIncrement = 1 Section 5: groupSplittingMethod = 1 - General Group Splitting Section 5: lengthOfLastGroup = 29959729 Section 5: nBitsGroupWidth = 0 Section 5: nBitsPacking = 1 Section 5: nBitsScaledGroupLength = 0 Section 5: nBytesSpatialDifference = 1 Section 5: nGroups = 0 Section 5: priMissingValue = 9.999000260554009e+20 Section 5: refGroupLength = 0 Section 5: refGroupWidth = 0 Section 5: refValue = 0.0 Section 5: secMissingValue = nan Section 5: spatialDifferenceOrder = 1 - First-Order Spatial Differencing Section 5: typeOfMissingValueManagement = 1 - Primary missing values included within the data values Section 6: bitMapFlag = 255 - A bit map does not apply to this product. >>> np.nanmin(message.data) 0.0 >>> np.nanmax(message.data) 81.1875 >>> message.pack() >>> grbout = grib2io.open('/u/mdl.wsup/scratch/grib2.out',mode='w') >>> grbout.write(message) >>> grbout.close() --------------------------- * Examine new GRIB from the CLI - File size is very small and max value is 0 $ ll /u/mdl.wsup/scratch/grib2.out -rw-r--r-- 1 mdl.wsup mdl 242 Apr 27 13:12 /u/mdl.wsup/scratch/grib2.out $ gdalinfo /u/mdl.wsup/scratch/grib2.out -mm|grep Min Computed Min/Max=0.000,0.000 $ degrib /u/mdl.wsup/scratch/grib2.out -C -stdout | grep Packing Info | Packing that was used | 3 (Grid point data - complex packing and spatial differencing) --------------------------- * Examine new GRIB with grib2io - Max value is 0 >>> sourceGrib='/u/mdl.wsup/scratch/grib2.out' >>> grbs = grib2io.open(sourceGrib) >>> grbs mode = rb name = /u/mdl.wsup/scratch/grib2.out messages = 1 current_message = 0 size = 242 closed = False variables = ('APCP',) levels = ('surface',) >>> grbs[0] Section 0: discipline = 0 - Meteorological Products Section 1: originatingCenter = 7 - US National Weather Service - NCEP (WMC) Section 1: originatingSubCenter = 4 - Environmental Modeling Center Section 1: masterTableInfo = 1 - Version Implemented on 7 November 2001 Section 1: localTableInfo = 0 - Local tables not used. Only table entries and templates from the current master table are valid. Section 1: significanceOfReferenceTime = 0 - Analysis Section 1: year = 2023 Section 1: month = 4 Section 1: day = 25 Section 1: hour = 18 Section 1: minute = 0 Section 1: second = 0 Section 1: refDate = 2023-04-25 18:00:00 Section 1: productionStatus = 0 - Operational Products Section 1: typeOfData = 0 - Analysis Products Section 3: interpretationOfListOfNumbers = 0 - There is no appended list Section 3: gridDefinitionTemplateNumber = 30 - Lambert Conformal (Can be Secant, Tangent, Conical, or Bipolar) Section 3: shapeOfEarth = 1 - Earth assumed spherical with radius specified (in m) by data producer Section 3: earthRadius = 6371200.0 Section 3: earthMajorAxis = None Section 3: earthMinorAxis = None Section 3: resolutionAndComponentFlags = [0, 0, 1, 1, 1, 0, 0, 0] Section 3: ny = 1597 Section 3: nx = 2345 Section 3: scanModeFlags = [0, 1, 0, 0, 0, 0, 0, 0] Section 3: gridOrientation = 265.0 Section 3: gridlengthXDirection = 2539.703 Section 3: gridlengthYDirection = 2539.703 Section 3: latitudeFirstGridpoint = 19.228976 Section 3: latitudeSouthernPole = -90.0 Section 3: latitudeTrueScale = 25.0 Section 3: longitudeFirstGridpoint = 233.723448 Section 3: longitudeSouthernPole = 0.0 Section 3: projParameters = {'a': 6371200.0, 'b': 6371200.0, 'proj': 'lcc', 'lat_1': 25.0, 'lat_2': 25.0, 'lat_0': 25.0, 'lon_0': 265.0} Section 3: projectionCenterFlag = 0 Section 3: standardLatitude1 = 25.0 Section 3: standardLatitude2 = 25.0 Section 4: productDefinitionTemplateNumber = 8 - Average, accumulation, extreme values or other statistically processed values at a horizontal level or in a horizontal layer in a continuous or non-continuous time interval. (see Template 4.8) Section 4: fullName = Total Precipitation Section 4: units = kg m-2 Section 4: shortName = APCP Section 4: leadTime = 6:00:00 Section 4: unitOfFirstFixedSurface = unknown Section 4: valueOfFirstFixedSurface = 2550.0 Section 4: unitOfSecondFixedSurface = None Section 4: valueOfSecondFixedSurface = 2550.0 Section 4: validDate = 2023-04-26 00:00:00 Section 4: duration = 6:00:00 Section 4: level = surface Section 4: backgroundGeneratingProcessIdentifier = 0 Section 4: dayOfEndOfTimePeriod = 26 Section 4: forecastTime = 0 Section 4: generatingProcess = 118 - UnRestricted Mesoscale Analysis (URMA) Section 4: hourOfEndOfTimePeriod = 0 Section 4: hoursAfterDataCutoff = 0 Section 4: minuteOfEndOfTimePeriod = 0 Section 4: minutesAfterDataCutoff = 0 Section 4: monthOfEndOfTimePeriod = 4 Section 4: numberOfMissingValues = 0 Section 4: numberOfTimeRanges = 1 Section 4: parameterCategory = 1 Section 4: parameterNumber = 8 Section 4: scaleFactorOfFirstFixedSurface = -1 Section 4: scaleFactorOfSecondFixedSurface = -127 Section 4: scaledValueOfFirstFixedSurface = 255 Section 4: scaledValueOfSecondFixedSurface = -2147483647 Section 4: secondOfEndOfTimePeriod = 0 Section 4: statisticalProcess = 1 - Accumulation Section 4: timeIncrementOfSuccessiveFields = 0 Section 4: timeRangeOfStatisticalProcess = 6 Section 4: typeOfFirstFixedSurface = 1 - ['Ground or Water Surface', 'unknown'] Section 4: typeOfGeneratingProcess = 2 - Forecast Section 4: typeOfSecondFixedSurface = 255 - ['Missing', 'unknown'] Section 4: typeOfTimeIncrementOfStatisticalProcess = 2 - Successive times processed have same start time of forecast, forecast time is incremented. Section 4: unitOfTimeRange = 1 - Hour Section 4: unitOfTimeRangeOfStatisticalProcess = 1 - Hour Section 4: unitOfTimeRangeOfSuccessiveFields = 255 - Missing Section 4: yearOfEndOfTimePeriod = 2023 Section 5: dataRepresentationTemplateNumber = 3 - Grid Point Data - Complex Packing and Spatial Differencing (see Template 5.3) Section 5: numberOfPackedValues = 3744965 Section 5: typeOfValues = 0 - Floating Point Section 5: binScaleFactor = -5 Section 5: decScaleFactor = 0 Section 5: groupLengthIncrement = 1 Section 5: groupSplittingMethod = 1 - General Group Splitting Section 5: lengthOfLastGroup = 29959729 Section 5: nBitsGroupWidth = 0 Section 5: nBitsPacking = 1 Section 5: nBitsScaledGroupLength = 0 Section 5: nBytesSpatialDifference = 1 Section 5: nGroups = 0 Section 5: priMissingValue = 9.999000260554009e+20 Section 5: refGroupLength = 0 Section 5: refGroupWidth = 0 Section 5: refValue = 0.0 Section 5: secMissingValue = nan Section 5: spatialDifferenceOrder = 1 - First-Order Spatial Differencing Section 5: typeOfMissingValueManagement = 1 - Primary missing values included within the data values Section 6: bitMapFlag = 255 - A bit map does not apply to this product. >>> np.nanmin(grbs[0].data) 0.0 >>> np.nanmax(grbs[0].data) 0.0
Got this working...see above commit references.
Note that the output GRIB2 message is larger than the input, despite being the same data and using the same packing scheme. wgrib2 does not use the g2clib library and instead has it own packing/unpacking routines.
Excellent! The WSUP team appreciates you addressing this and have taken note about output file sizes. Larger files full of good data are better than smaller sizes containing bad data, but that's just me.
When I read URMA QPE01 and QPE06 GRIBs using the latest version of grib2io on WCOSS2 (2.0.0b2), all of the values are all 'nan'.
Acquire QPE06 data:
Set up environment:
Python:
I can see via gdalinfo and QGIS that there are valid values in the GRIB file: