search

USACE / cumulus

Cumulus project issue tracking and project planning
MIT License
3 stars 2 forks source link

[BUG REPORT]: NULL values in WPC grid results #316

Closed adamscarberry closed 1 year ago

adamscarberry commented 2 years ago

Report from Ashley in SAW

image

Direct link of download: https://cumulus-api.corps.cloud/cumulus/download/dss7/download_a3bbd4f2-128b-45fb-8157-9fa4ce2480b8.dss

image

adamscarberry commented 2 years ago

Further analysis for WPC QPF, but on a different date/time.

Confirmation of which tif generated the DSS record below (packager payload snippet):

{
    "key": "cumulus/products/wpc-qpf-2p5km/p06m_2022110112f018.tif",
    "bucket": "castle-data-develop",
    "dss_unit": "MM",
    "dss_cpart": "PRECIP",
    "dss_dpart": "02NOV2022:0000",
    "dss_epart": "02NOV2022:0600",
    "dss_fpart": "WPC-QPF-2.5KM",
    "dss_datatype": "PER-CUM"
  }

DSS v7 File results on DSS-Vue 3.3.8 Beta

Max Data Value: -2.0E38
Min Data Value: 2.0E38
Mean Data Value: NaN
Number Of Ranges: 3

image

Grib Stats

> gdalinfo p06m_2022110112f018.grb -stats

Driver: GRIB/GRIdded Binary (.grb, .grb2)
Files: p06m_2022110112f018.grb
Size is 2145, 1597
Coordinate System is:
PROJCRS["unnamed",
    BASEGEOGCRS["Coordinate System imported from GRIB file",
        DATUM["unnamed",
            ELLIPSOID["Sphere",6371200,0,
                LENGTHUNIT["metre",1,
                    ID["EPSG",9001]]]],
        PRIMEM["Greenwich",0,
            ANGLEUNIT["degree",0.0174532925199433,
                ID["EPSG",9122]]]],
    CONVERSION["Lambert Conic Conformal (2SP)",
        METHOD["Lambert Conic Conformal (2SP)",
            ID["EPSG",9802]],
        PARAMETER["Latitude of false origin",25,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8821]],
        PARAMETER["Longitude of false origin",-95,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8822]],
        PARAMETER["Latitude of 1st standard parallel",25,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8823]],
        PARAMETER["Latitude of 2nd standard parallel",25,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8824]],
        PARAMETER["Easting at false origin",0,
            LENGTHUNIT["metre",1],
            ID["EPSG",8826]],
        PARAMETER["Northing at false origin",0,
            LENGTHUNIT["metre",1],
            ID["EPSG",8827]]],
    CS[Cartesian,2],
        AXIS["(E)",east,
            ORDER[1],
            LENGTHUNIT["Metre",1]],
        AXIS["(N)",north,
            ORDER[2],
            LENGTHUNIT["Metre",1]]]
Data axis to CRS axis mapping: 1,2
Origin = (-2764474.350731992628425,3790844.774292394984514)
Pixel Size = (2539.702999999999975,-2539.702000000000226)
Corner Coordinates:
Upper Left  (-2764474.351, 3790844.774) (132d 0'57.11"W, 54d12'45.58"N)
Lower Left  (-2764474.351, -265059.320) (121d33'48.72"W, 20d10'43.04"N)
Upper Right ( 2683188.584, 3790844.774) ( 59d 1'17.21"W, 54d22'46.87"N)
Lower Right ( 2683188.584, -265059.320) ( 69d11'54.73"W, 20d19' 6.40"N)
Center      (  -40642.883, 1762892.727) ( 95d27'46.89"W, 40d38'51.51"N)
Band 1 Block=2145x1 Type=Float64, ColorInterp=Undefined
  Description = 0--1[-] SFC="Ground or water surface"
  Minimum=0.000, Maximum=30.500, Mean=0.435, StdDev=1.584
  Metadata:
    GRIB_COMMENT=06 hr Total precipitation [kg/(m^2)]
    GRIB_DISCIPLINE=0(Meteorological)
    GRIB_ELEMENT=APCP06
    GRIB_FORECAST_SECONDS=43200
    GRIB_IDS=CENTER=7(US-NCEP) SUBCENTER=5(Hydrometeorological Prediction Center) MASTER_TABLE=2 LOCAL_TABLE=1 SIGNF_REF_TIME=1(Start_of_Forecast) REF_TIME=2022-11-01T12:00:00Z PROD_STATUS=0(Operational) TYPE=1(Forecast)
    GRIB_PDS_PDTN=8
    GRIB_PDS_TEMPLATE_ASSEMBLED_VALUES=1 8 2 0 180 255 255 1 12 1 0 0 1 0 -1 2022 11 2 6 0 0 1 0 1 255 1 6 1 0
    GRIB_PDS_TEMPLATE_NUMBERS=1 8 2 0 180 0 255 255 1 0 0 0 12 1 0 0 0 0 0 1 0 128 0 0 1 7 230 11 2 6 0 0 1 0 0 0 0 1 255 1 0 0 0 6 1 0 0 0 0
    GRIB_REF_TIME=1667304000
    GRIB_SHORT_NAME=0--1-SFC
    GRIB_UNIT=[kg/(m^2)]
    GRIB_VALID_TIME=1667368800
    STATISTICS_MAXIMUM=30.5
    STATISTICS_MEAN=0.43528866770884
    STATISTICS_MINIMUM=0
    STATISTICS_STDDEV=1.5841776450257
    STATISTICS_VALID_PERCENT=100

Tif Stats

> gdalinfo p06m_2022110112f018.tif -stats
Driver: GTiff/GeoTIFF
Files: p06m_2022110112f018.tif
Size is 2145, 1597
Coordinate System is:
PROJCRS["unnamed",
    BASEGEOGCRS["Coordinate System imported from GRIB file",
        DATUM["unnamed",
            ELLIPSOID["Sphere",6371200,0,
                LENGTHUNIT["metre",1,
                    ID["EPSG",9001]]]],
        PRIMEM["Greenwich",0,
            ANGLEUNIT["degree",0.0174532925199433,
                ID["EPSG",9122]]]],
    CONVERSION["Lambert Conic Conformal (2SP)",
        METHOD["Lambert Conic Conformal (2SP)",
            ID["EPSG",9802]],
        PARAMETER["Latitude of false origin",25,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8821]],
        PARAMETER["Longitude of false origin",-95,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8822]],
        PARAMETER["Latitude of 1st standard parallel",25,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8823]],
        PARAMETER["Latitude of 2nd standard parallel",25,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8824]],
        PARAMETER["Easting at false origin",0,
            LENGTHUNIT["metre",1],
            ID["EPSG",8826]],
        PARAMETER["Northing at false origin",0,
            LENGTHUNIT["metre",1],
            ID["EPSG",8827]]],
    CS[Cartesian,2],
        AXIS["easting",east,
            ORDER[1],
            LENGTHUNIT["metre",1,
                ID["EPSG",9001]]],
        AXIS["northing",north,
            ORDER[2],
            LENGTHUNIT["metre",1,
                ID["EPSG",9001]]]]
Data axis to CRS axis mapping: 1,2
Origin = (-2764474.350731992628425,3790844.774292394984514)
Pixel Size = (2539.702999999999975,-2539.702000000000226)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  COMPRESSION=LZW
  INTERLEAVE=BAND
  LAYOUT=COG
Corner Coordinates:
Upper Left  (-2764474.351, 3790844.774) (132d 0'57.11"W, 54d12'45.58"N)
Lower Left  (-2764474.351, -265059.320) (121d33'48.72"W, 20d10'43.04"N)
Upper Right ( 2683188.584, 3790844.774) ( 59d 1'17.21"W, 54d22'46.87"N)
Lower Right ( 2683188.584, -265059.320) ( 69d11'54.73"W, 20d19' 6.40"N)
Center      (  -40642.883, 1762892.727) ( 95d27'46.89"W, 40d38'51.51"N)
Band 1 Block=512x512 Type=Float64, ColorInterp=Gray
  Description = 0--1[-] SFC="Ground or water surface"
  Minimum=0.000, Maximum=30.500, Mean=0.435, StdDev=1.584
  Overviews: 1072x798, 536x399, 268x199
  Metadata:
    GRIB_COMMENT=06 hr Total precipitation [kg/(m^2)]
    GRIB_DISCIPLINE=0(Meteorological)
    GRIB_ELEMENT=APCP06
    GRIB_FORECAST_SECONDS=43200
    GRIB_IDS=CENTER=7(US-NCEP) SUBCENTER=5(Hydrometeorological Prediction Center) MASTER_TABLE=2 LOCAL_TABLE=1 SIGNF_REF_TIME=1(Start_of_Forecast) REF_TIME=2022-11-01T12:00:00Z PROD_STATUS=0(Operational) TYPE=1(Forecast)
    GRIB_PDS_PDTN=8
    GRIB_PDS_TEMPLATE_ASSEMBLED_VALUES=1 8 2 0 180 255 255 1 12 1 0 0 1 0 -1 2022 11 2 6 0 0 1 0 1 255 1 6 1 0
    GRIB_PDS_TEMPLATE_NUMBERS=1 8 2 0 180 0 255 255 1 0 0 0 12 1 0 0 0 0 0 1 0 128 0 0 1 7 230 11 2 6 0 0 1 0 0 0 0 1 255 1 0 0 0 6 1 0 0 0 0
    GRIB_REF_TIME=1667304000
    GRIB_SHORT_NAME=0--1-SFC
    GRIB_UNIT=[kg/(m^2)]
    GRIB_VALID_TIME=1667368800
    STATISTICS_MAXIMUM=30.5
    STATISTICS_MEAN=0.43528866770883
    STATISTICS_MINIMUM=0
    STATISTICS_STDDEV=1.5841776450257
    STATISTICS_VALID_PERCENT=100
adamscarberry commented 2 years ago

Keep in mind both the grib and tif file above are CONUS (not clipped). The clipping occurs during the DSS packaging process.

adamscarberry commented 2 years ago

Here are the stats (abbreviated) for the clipped (Cape Fear River) tif for the same hour as above.

 gdalinfo p06m_2022110112f018.tif -stats

...

Metadata:
    GRIB_COMMENT=06 hr Total precipitation [kg/(m^2)]
    GRIB_DISCIPLINE=0(Meteorological)
    GRIB_ELEMENT=APCP06
    GRIB_FORECAST_SECONDS=43200
    GRIB_IDS=CENTER=7(US-NCEP) SUBCENTER=5(Hydrometeorological Prediction Center) MASTER_TABLE=2 LOCAL_TABLE=1 SIGNF_REF_TIME=1(Start_of_Forecast) REF_TIME=2022-11-01T12:00:00Z PROD_STATUS=0(Operational) TYPE=1(Forecast)
    GRIB_PDS_PDTN=8
    GRIB_PDS_TEMPLATE_ASSEMBLED_VALUES=1 8 2 0 180 255 255 1 12 1 0 0 1 0 -1 2022 11 2 6 0 0 1 0 1 255 1 6 1 0
    GRIB_PDS_TEMPLATE_NUMBERS=1 8 2 0 180 0 255 255 1 0 0 0 12 1 0 0 0 0 0 1 0 128 0 0 1 7 230 11 2 6 0 0 1 0 0 0 0 1 255 1 0 0 0 6 1 0 0 0 0
    GRIB_REF_TIME=1667304000
    GRIB_SHORT_NAME=0--1-SFC
    GRIB_UNIT=[kg/(m^2)]
    GRIB_VALID_TIME=1667368800
    STATISTICS_MAXIMUM=0
    STATISTICS_MEAN=0
    STATISTICS_MINIMUM=0
    STATISTICS_STDDEV=0
    STATISTICS_VALID_PERCENT=100
ktarbet commented 2 years ago

@katfeingold and I found where the Infinity is coming from in GridInfoFlat.java

public boolean unflatten (int array[]) {
//....
 // deal with special cases of min, max, and mean values for all-null cases
    if (_minDataValue > _maxDataValue){
        _maxDataValue = RMAConst.UNDEF_FLOAT;    
        _minDataValue = RMAConst.UNDEF_FLOAT;
        _meanDataValue = RMAConst.UNDEF_FLOAT;      
    }
// where: public static final float UNDEF_FLOAT = Float.NEGATIVE_INFINITY;

the Infinity is intentional, but we could change this if other programs don't rely on that behavior?

FHanbali commented 2 years ago

@adamscarberry can you share with me please the source grb files for the forecast timewindow 28July2022 throuhg 31July at least?

@ktarbet I believe the issue is not actually about the statistics reported in the Gridded DSS Record's header. Rather it's about the array of values for grid cells, they seem to be coming out as -infiniti during conversion from original values of 0.

jeffsuperglide commented 2 years ago

@katfeingold and I found where the Infinity is coming from in GridInfoFlat.java

public boolean unflatten (int array[]) {
//....
 // deal with special cases of min, max, and mean values for all-null cases
    if (_minDataValue > _maxDataValue){
      _maxDataValue = RMAConst.UNDEF_FLOAT;    
      _minDataValue = RMAConst.UNDEF_FLOAT;
      _meanDataValue = RMAConst.UNDEF_FLOAT;      
    }
// where: public static final float UNDEF_FLOAT = Float.NEGATIVE_INFINITY;

the Infinity is intentional, but we could change this if other programs don't rely on that behavior?

@ktarbet I thought we found this about a month ago? I remember some DSS6 code that computes stats again before this evaluation.

GDAL is used in tiffdss to get the data array. Stats are computed in tiffdss with concepts taken from DSS6/7 code. @ktarbet is correct that it's about the grid data array and not stats computed from that array. Evaluate the grib, GeoTiff, and the clipped GeoTiff only pulling out the data array.

adamscarberry commented 2 years ago

This is a sample data set (based on my stats output above) that I can get you quickly:

Archive.zip

ktarbet commented 2 years ago

@jeffsuperglide Yes! we recently had https://github.com/USACE/cumulus/issues/231 and https://github.com/USACE/cumulus/issues/240

More details on what we saw yesterday in the debugger:

All the data in the grid-data-array was listed as missing. Then the statistics were re-computed during the DSS7->DSS6 conversion which creates the Infinity values for (min,max,mean), using the java code block above.

We need to change either:

In a previous iteration, of this issue, we made a change to DssVue to handle this case better.

@FHanbali Was there was a recent change to MetVue to handle this case?

FHanbali commented 2 years ago

Thanks Adam. Geotiff sample shows all zero (valid) values. Either Geotiff conversion to DSS7 or DSS7 to DSS6 conversion is causing null values to be stored in final DSS product.

jeffsuperglide commented 2 years ago

This is a sample data set (based on my stats output above) that I can get you quickly:

Archive.zip

Two Tiff files are provided in the archive, conus_p06m_2022110112f018.tif and clipped_p06m_2022110112f018.tif. conus has a no data value (nullValue) of 0 and the clipped no data value is -9999. A no data value of 0 produces the following because all values == no data value.

Max Data Value: -2.0E38
Min Data Value: 2.0E38
Mean Data Value: NaN

Having -9999 (or != 0) results in zeros for stats. The packager uses gdal.Warp() without defining the no data value; therefore, the no data value is taken from the source Tiff. Packager writer dss7 could check for that no data value and set to some non-zero value.

ktarbet commented 2 years ago

:: draft in progress table of workflow, input welcome

Step file description range of data no data value program
1 raw grib file 6 hour precipitation (mm) 0-30.5 not defined raw input
2 tIf - unclipped 6 hour precipitation 0-24.9 not defined python/gdal
3 tIf - clipped 6 hour precipitation 0-0 not defined python/gdal
3.5 in-memory grid 6 hour precipitation 0-0 0 python/gdal
4 DSS v7 6 hour precipitation all missing -3.4028235E38 tif2dss
5 DSS v6 6 hour precipitation all missing -3.4028235E38 java/heclib
6 MetVue grid processing ? ? java/heclib
7 HMS hydrology ? ? java/heclib
jeffsuperglide commented 1 year ago

@FHanbali, @ktarbet, @adamscarberry, and @brettpalmberg

Tested locally, pushed to develop, and tested on develop-cumulus.corps.cloud. Made a couple requests for Savannah River Basin, WPC QPF. Below is example screen shot of what zero precip looks like.

image
adamscarberry commented 1 year ago

@jeffsuperglide and @brettpalmberg We ready for push to stable?

jeffsuperglide commented 1 year ago

@jeffsuperglide and @brettpalmberg We ready for push to stable?

I say yes.

FHanbali commented 1 year ago

Let me know please if https://cumulus.corps.cloud/ is ready to have folks try out the fix.

adamscarberry commented 1 year ago

@FHanbali @jeffsuperglide's patch has been applied to stable. Ready for final verification.

brettpalmberg commented 1 year ago

Kind ping 😃. This work is complete and awaiting final verification. Hoping we can officially close this issue.

@FHanbali @jeffsuperglide @adamscarberry @ktarbet

FHanbali commented 1 year ago

Sorry, I should have been more clear with my feedback. Tested this a couple of weeks ago and worked as expected. Also got verification from field users. All good, thank you again.

brettpalmberg commented 1 year ago

Excellent! No worries and good teamwork all. Closing this issue now.