Skip to content

Real time stageiv precip

We have a few applications that would benefit from having real time StageIV precipitation back up and running.

To get there, we need a pipeline that can download data on a schedule, open it in xarray, and append it to a zarr store.

[edit: 6/10/24] 2016-present NetCDF formatted data can be retrieved from https://water.noaa.gov/resources/downloads/precip/stageIV/ -- I have not done any investigation to see if it is the same as the grib source below. It does appear that the legacy grib data are available there -- summary of all that is available here: https://water.noaa.gov/about/precipitation-data-access

Realtime data for the last two weeks can be retrieved from: https://nomads.ncep.noaa.gov/pub/data/nccf/com/pcpanl/prod/

We will also need to download a bit of the archives from here: https://data.eol.ucar.edu/dataset/21.093 to back fill what's missing from the store that is already available on WMA object storage.

This paper is the initial reference for this work: https://www.usgs.gov/publications/new-service-interface-river-forecasting-center-derived-quantitative-precipitation

It appears that cfgrib and xarray will get us where we need to go -- I'll do a bit of hacking to make sure we aren't barking up the wrong tree.

Yeah, this looks promising.

A little validation that we can reproduce the data in the existing archive would be a good step to make sure things look ok.

import eccodes
import cfgrib
import xarray

xarray.open_dataset("st4_conus.2024060500.01h.grb2", engine='cfgrib')
Here be details of the data

image

Note that cfgrib does not have the correct grib table (metadata lookup for the internal grib codes) for this data. NetCDF-Java ToolsUI is a good utility to review these data more accurately. Here is the raw ncdump from toolsui:

Here be more details of the data
netcdf C:/Users/dblodgett/Downloads/st4_conus.2024060500.01h.grb2 {
  dimensions:
    x = 1121;
    y = 881;
    reftime = 1;
    time = 1;
  variables:
    int PolarStereographic_Projection;
      :grid_mapping_name = "stereographic";
      :longitude_of_projection_origin = 255.0; // double
      :latitude_of_projection_origin = 90.0; // double
      :scale_factor_at_projection_origin = 0.9330127018922193; // double
      :earth_radius = 6371229.0; // double

    float x(x=1121);
      :standard_name = "projection_x_coordinate";
      :units = "km";

    float y(y=881);
      :standard_name = "projection_y_coordinate";
      :units = "km";

    double reftime;
      :units = "Hour since 2024-06-04T23:00:00Z";
      :standard_name = "forecast_reference_time";
      :long_name = "GRIB reference time";
      :calendar = "proleptic_gregorian";

    double time(time=1);
      :units = "Hour since 2024-06-04T23:00:00Z";
      :standard_name = "time";
      :long_name = "GRIB forecast or observation time";
      :calendar = "proleptic_gregorian";
      :bounds = "time_bounds";

    double time_bounds(time=1, 2);
      :units = "Hour since 2024-06-04T23:00:00Z";
      :long_name = "bounds for time";

    float Total_precipitation_surface_1_Hour_Accumulation(time=1, y=881, x=1121);
      :long_name = "Total precipitation (1_Hour Accumulation) @ Ground or water surface";
      :units = "kg.m-2";
      :abbreviation = "APCP";
      :missing_value = NaNf; // float
      :grid_mapping = "PolarStereographic_Projection";
      :coordinates = "reftime time y x ";
      :Grib_Statistical_Interval_Type = "Accumulation";
      :Grib_Variable_Id = "VAR_0-1-8_L1_I1_Hour_S1";
      :Grib2_Parameter = 0, 1, 8; // int
      :Grib2_Parameter_Discipline = "Meteorological products";
      :Grib2_Parameter_Category = "Moisture";
      :Grib2_Parameter_Name = "Total precipitation";
      :Grib2_Level_Type = 1; // int
      :Grib2_Level_Desc = "Ground or water surface";
      :Grib2_Generating_Process_Type = "Forecast";
      :Grib2_Statistical_Process_Type = "Accumulation";

  // global attributes:
  :Originating_or_generating_Center = "US National Weather Service, National Centres for Environmental Prediction (NCEP)";
  :Originating_or_generating_Subcenter = "Environmental Modeling Center";
  :GRIB_table_version = "2,1";
  :Type_of_generating_process = "Forecast";
  :Analysis_or_forecast_generating_process_identifier_defined_by_originating_centre = "River Forecast Center Quantitative Precipitation estimate mosaic generated by NCEP";
  :file_format = "GRIB-2";
  :Conventions = "CF-1.6";
  :history = "Read using CDM IOSP GribCollection v3";
  :featureType = "GRID";
}

Here is the zarr collection we are working with:

Here be details of the data

> url <- "https://usgs.osn.mghpcc.org/mdmf/gdp/stageiv_combined.zarr/"
> rnz::nzdump(url)
zarr {
dimensions:
time = 191948 ;
y = 881 ;
x = 1121 ;
variables:
	<f4 Total_precipitation_surface_1_Hour_Accumulation(time, y, x) ;
		Total_precipitation_surface_1_Hour_Accumulation:cell_methods = time: sum (interval: 1 hr) ;
		Total_precipitation_surface_1_Hour_Accumulation:coordinates = time lat lon ;
		Total_precipitation_surface_1_Hour_Accumulation:description = Total precipitation ;
		Total_precipitation_surface_1_Hour_Accumulation:grid_mapping = crs: lat lon ;
		Total_precipitation_surface_1_Hour_Accumulation:long_name = Total precipitation (1_Hour Accumulation) @ Ground or water surface ;
		Total_precipitation_surface_1_Hour_Accumulation:missing_value =  ;
		Total_precipitation_surface_1_Hour_Accumulation:standard_name = precipitation_amount ;
		Total_precipitation_surface_1_Hour_Accumulation:units = kg m^-2 ;
	<i4 crs() ;
		crs:crs_wkt = GEOGCRS["WGS 84",ENSEMBLE["World Geodetic System 1984 ensemble",MEMBER["World Geodetic System 1984 (Transit)"],MEMBER["World Geodetic System 1984 (G730)"],MEMBER["World Geodetic System 1984 (G873)"],MEMBER["World Geodetic System 1984 (G1150)"],MEMBER["World Geodetic System 1984 (G1674)"],MEMBER["World Geodetic System 1984 (G1762)"],MEMBER["World Geodetic System 1984 (G2139)"],ELLIPSOID["WGS 84",6378137,298.257223563,LENGTHUNIT["metre",1]],ENSEMBLEACCURACY[2.0]],PRIMEM["Greenwich",0,ANGLEUNIT["degree",0.0174532925199433]],CS[ellipsoidal,2],AXIS["geodetic latitude (Lat)",north,ORDER[1],ANGLEUNIT["degree",0.0174532925199433]],AXIS["geodetic longitude (Lon)",east,ORDER[2],ANGLEUNIT["degree",0.0174532925199433]],USAGE[SCOPE["Horizontal component of 3D system."],AREA["World."],BBOX[-90,-180,90,180]],ID["EPSG",4326]] ;
		crs:geographic_crs_name = WGS 84 ;
		crs:grid_mapping_name = latitude_longitude ;
		crs:horizontal_datum_name = World Geodetic System 1984 ensemble ;
		crs:inverse_flattening = 298.257223563 ;
		crs:longitude_of_prime_meridian = 0 ;
		crs:prime_meridian_name = Greenwich ;
		crs:reference_ellipsoid_name = WGS 84 ;
		crs:semi_major_axis = 6378137 ;
		crs:semi_minor_axis = 6356752.31424518 ;
	<f4 lat(x, y) ;
		lat:units = degrees_north ;
	<f4 lon(x, y) ;
		lon:units = degrees_east ;
	<f8 time(time) ;
		time:bounds = time_bounds ;
		time:calendar = proleptic_gregorian ;
		time:long_name = GRIB forecast or observation time ;
		time:standard_name = time ;
		time:units = hours since 2001-12-31T23:00:00+00:00 ;
	<f8 time_bounds(time, time_bounds_1) ;
		time_bounds:calendar = proleptic_gregorian ;
		time_bounds:long_name = bounds for time ;
		time_bounds:units = hours since 2001-12-31T23:00:00+00:00 ;

// global attributes:
		:Conventions = CF-1.4 ;
		:Generating_process_or_model = River Forecast Center Quantitative Precipitation estimate mosaic generated by NCEP ;
		:Metadata_Conventions = Unidata Dataset Discovery v1.0 ;
		:Originating_or_generating_Center = US National Weather Service, National Centres for Environmental Prediction (NCEP) ;
		:Originating_or_generating_Subcenter = Environmental Modeling Center ;
		:_CoordSysBuilder = ucar.nc2.dataset.conv.CF1Convention ;
		:acknowledgement = Data originated at the National Center for Environmental Prediction: http://www.emc.ncep.noaa.gov/mmb/ylin/pcpanl/stage4/ ;
		:cdm_data_type = Grid ;
		:contributor_name = Center for Integrated Data Analytics ;
		:contributor_role = Aggregation and Serving via THREDDS ;
		:date_created = ongoing ;
		:featureType = GRID ;
		:geospatial_lat_max = 53 ;
		:geospatial_lat_min = 24 ;
		:geospatial_lat_resolution = 4 ;
		:geospatial_lat_units = km ;
		:geospatial_lon_max = -66 ;
		:geospatial_lon_min = -125 ;
		:geospatial_lon_resolution = 4 ;
		:geospatial_lon_units = km ;
		:geospatial_vertical_units = km ;
		:history = Aggregated from NCEP Stage IV Analysis by dblodgett@usgs.gov 2016-01-01 ;
		:id =  stageiv ;
		:institution = National Center for Environmental Prediction ;
		:keywords = Atmosphere > Precipitation > Precipitation Amount   Atmosphere > Precipitation > Rain ;
		:keywords_vocabulary = GCMD Science Keywords ;
		:license = Freely available ;
		:naming_authority = cida.usgs.gov ;
		:processing_level = Aggregated .grib archive from the NCEP. The original Hydrologic Rainfall Analysis Projection grid cell locations redefined by unprojecting grid cell locations using spherical earth projection algorithms.    The cell lat/lon locations should be interpreted as being on a modern geotetic datum. This transformation is needed because the original grid cell location's lat/lon values were on a geodetic datum but were assumed to be on a spherical datum.    See: https://pubs.usgs.gov/fs/2013/3035/pdf/fs2013-3035.pdf for more information. Note that this publication reverences the National Precipitation Verification Unit.   The data represented here are aggregated from the same source as the National Precipitation Verification Unit used prior to being shut down. Data processing notes can be found here: http://www.emc.ncep.noaa.gov/mmb/ylin/pcpanl/ ;
		:project = Various ;
		:summary = Mosaicked into a national product at NCEP, from the regional hourly multi-sensor (radar+gauges) precipitation analyses produced by the 12 River Forecast Centers over CONUS. Some manual QC done at the RFCs.  The Stage II/IV job is run at 33min past the top of each hour. Hourly Stage IV is re-made hourly (if there is new input after valid time for the next 23 hours, then again at 1/3/5/7 days after valid time.  Source: http://www.emc.ncep.noaa.gov/mmb/ylin/pcpanl/stage4/ Note: This dataset is not truely national. Some RFCs do not produce hourly precip estimates and others may not be present in the dataset. Please review the NCEP documentation for more information. ;
		:time_coverage_end = present ;
		:time_coverage_start = 2002-01-01T00:00:00Z ;
		:title = United States Stage IV Quantitative Precipitation Archive ;
}

Note that the coordinate variables are corrected in the zarr store -- we only need to append to time, time_bounds, and Total_precipitation_surface_1_Hour_Accumulation.

Edited by Blodgett, David L.