Scheduled Downtime
On Friday 21 April 2023 @ 5pm MT, this website will be down for maintenance and expected to return online the morning of 24 April 2023 at the latest

Issue with Custom Land Use Data (BIL) Not Being Applied in geogrid.exe

aba0327

New member
I am trying to use a custom EGIS Land Use dataset for WRF-UCM by converting it into a binary (BIL) file. However, after running geogrid.exe, the generated geo_em.d03.nc file does not reflect my custom dataset, and LU_INDEX remains at 17 across the entire domain.

Steps I Followed:

  1. Converted an EGIS shapefile to raster using Python (rasterio).
  2. Saved it as a BIL filewith EHdr format (8-bit categorical).
    • Python:
      ds["MODIS_CODE"] = ds["MODIS_CODE"].astype(np.uint8)
      ds = ds.to_crs("EPSG:4326")
      
      resolution = 0.00083333
      
      bounds = ds.total_bounds
      left, bottom, right, top = bounds
      
      width = int((right - left) / resolution)
      height = int((top - bottom) / resolution)
      
      raster_array = np.full((height, width), 0, dtype=np.uint8)
      
      shapes = [(geom, value) for geom, value in zip(ds.geometry, ds["MODIS_CODE"])]
      
      rasterized = rasterize(
          shapes,
          out_shape=(height, width),
          transform=from_origin(left, top, resolution, resolution),
          fill=0,
          dtype=np.uint8
      )
      
      binary_output = "./EGIS_LANDUSE.2023.v8/EGIS_LANDUSE.bil"
      
      with rasterio.open(
          binary_output,
          "w",
          driver="EHdr",
          width=width,
          height=height,
          count=1,
          dtype=np.uint8,
          crs="EPSG:4326",
          transform=from_origin(left, top, resolution, resolution),
      ) as dst:
          dst.write(rasterized, 1)
  3. Configured GEOGRID.TBLto reference my dataset (egis_landuse).
    • INI:
      name=LANDUSEF
              priority = 2
              dest_type = categorical
              z_dim_name = land_cat
              dominant = LU_INDEX
              landmask_water =  egis_landuse:17,21      # Calculate a landmask from this field
              interp_option =  egis_landuse:nearest_neighbor
              rel_path =   egis_landuse:egis_landuse/
  4. Ran geogrid.exe, and the log confirms that my dataset is being processed.
  5. Checked geo_em.d03.nc, but my land use data is not applied.

Verification Checks & Findings

  1. Endian Check:
    • Used od -An -t d1 00001-01890.00001-01200 | head -n 20
    • Values appear correctly stored as uint8 integers.
  2. Coordinate System Check:
    • ncdump -h geo_em.d03.nc | grep -i 'MAP_PROJ' → MAP_PROJ = 1 (Lambert Conformal Projection).
    • ncdump -h geo_em.d03.nc | grep -i 'DX\|DY' → DX = 0.00083333, DY = 0.00083333.
  3. BIL File Header & Data Integrity:
    • BIL file size matches expected dimensions.
    • hexdump -C shows diverse land use values.
    • row_order = top_bottom was confirmed in the metadata.

Problem & Observations:

  • geogrid.log confirms that egis_landuse is being processed.
  • However, in geo_em.d03.nc, LU_INDEX remains 17 across all points instead of reflecting EGIS data.
  • When using geog_res_path = egis_landuse+default, LU_INDEX is set to a single code (not reflecting EGIS values).
  • When using geog_res_path = default+egis_landuse, only the default dataset is applied.

Questions:

  1. Is using rasterio instead of GDAL causing an issue with the BIL file format?
  2. Does WRF require a specific Endian format for BIL files? How can I ensure my file is formatted correctly?
  3. Should I reproject my dataset to match WRF’s Lambert Projection instead of using EPSG:4326?
  4. How can I verify that geogrid.exe is correctly reading my binary file and not defaulting to MODIS?
  5. Are there specific settings in GEOGRID.TBL that must be changed to force egis_landuse to overwrite default land use data?
 
Apologies for the delay in response, and thank you for your patience. Are you still experiencing this issue? If so,

1) Have you plotted the landuse data from the input EGIS data to verify the landuse data is accurate in the input file?
2) Can you attach your namelist.wps file, geogrid.log, and the modified GEOGRID.TBL so we can take a look?
 
Dear kwerner,

Thank you for your response. I appreciate your help.

I am still experiencing the issue and have not yet been able to resolve it.

  1. I have converted the shapefile to a binary file and visualized the output. I have attached the visualization, where it appears that the land use codes have been correctly assigned. The regions with code 0 correspond to areas where no data was present in the original shapefile.1741842972753.png
  2. I have attached the requested files: namelist.wps, geogrid.log, and GEOGRID.TBL (land use-related sections only).
    • namelist.wps
    • geogrid.log
    • Modified GEOGRID.TBL: I have included only the land use-related sections. Other parts of the file remain unchanged.
      • INI:
        ===============================
        name=LANDUSEF
                priority = 2
                dest_type = categorical
                z_dim_name = land_cat
                dominant = LU_INDEX
                landmask_water =  egis_landuse:17,21      # Calculate a landmask from this field
                interp_option =  egis_landuse:nearest_neighbor
                rel_path =   egis_landuse:egis_landuse/
        name=LANDUSEF
                priority = 1
                dest_type = categorical
                z_dim_name = land_cat
                dominant = LU_INDEX
                landmask_water =    nlcd2006_9s:17         # Calculate a landmask from this field
                landmask_water =   nlcd2006_30s:17         # Calculate a landmask from this field
                landmask_water =    nlcd2011_9s:17         # Calculate a landmask from this field
                landmask_water =       nlcd2006:17         # Calculate a landmask from this field
                landmask_water =       ssib_10m:16         # Calculate a landmask from this field
                landmask_water =        ssib_5m:16         # Calculate a landmask from this field
                landmask_water =      modis_15s:17         # Calculate a landmask from this field
                landmask_water = modis_15s_lake:17,21      # Calculate a landmask from this field
                landmask_water =      modis_30s:17         # Calculate a landmask from this field
                landmask_water = modis_30s_lake:17,21      # Calculate a landmask from this field
                landmask_water =       usgs_30s:16         # Calculate a landmask from this field
                landmask_water =     usgs_lakes:16,28      # Calculate a landmask from this field
                landmask_water =    modis_lakes:17,21      # Calculate a landmask from this field
                landmask_water =     modis_2010:17         # Calculate a landmask from this field
                landmask_water = modis_2010_lake:17,21     # Calculate a landmask from this field
                landmask_water =     modis_2050:17         # Calculate a landmask from this field
                landmask_water = modis_2050_lake:17,21     # Calculate a landmask from this field
                landmask_water =        usgs_2m:16         # Calculate a landmask from this field
                landmask_water =        usgs_5m:16         # Calculate a landmask from this field
                landmask_water =       usgs_10m:16         # Calculate a landmask from this field
                landmask_water =         lowres:17,21      # Calculate a landmask from this field
                landmask_water =        default:17,21      # Calculate a landmask from this field
                interp_option =    nlcd2006_9s:average_gcell(0.0)
                interp_option =   nlcd2006_30s:average_gcell(0.0)
                interp_option =    nlcd2011_9s:average_gcell(0.0) 
                interp_option =       nlcd2006:nearest_neighbor
                interp_option =       ssib_10m:four_pt
                interp_option =        ssib_5m:four_pt
                interp_option =      modis_15s:nearest_neighbor
                interp_option = modis_15s_lake:nearest_neighbor
                interp_option =      modis_30s:nearest_neighbor
                interp_option = modis_30s_lake:nearest_neighbor
                interp_option =       usgs_30s:nearest_neighbor
                interp_option =     usgs_lakes:nearest_neighbor
                interp_option =    modis_lakes:nearest_neighbor
                interp_option =     modis_2010:nearest_neighbor
                interp_option = modis_2010_lake:nearest_neighbor
                interp_option =     modis_2050:nearest_neighbor
                interp_option = modis_2050_lake:nearest_neighbor
                interp_option =        usgs_2m:four_pt
                interp_option =        usgs_5m:four_pt
                interp_option =       usgs_10m:four_pt
                interp_option =         lowres:average_gcell(4.0)+four_pt
                interp_option =        default:nearest_neighbor
                rel_path =    nlcd2006_9s:nlcd2006_ll_9s/
                rel_path =   nlcd2006_30s:nlcd2006_ll_30s/
                rel_path =    nlcd2011_9s:nlcd2011_ll_9s/
                rel_path =       nlcd2006:nlcd2006_ll_30s/
                rel_path =       ssib_10m:ssib_landuse_10m/
                rel_path =        ssib_5m:ssib_landuse_5m/
                rel_path =      modis_15s:modis_landuse_20class_15s/
                rel_path = modis_15s_lake:modis_landuse_20class_15s_with_lakes/
                rel_path =      modis_30s:modis_landuse_20class_30s/
                rel_path = modis_30s_lake:modis_landuse_20class_30s_with_lakes/
                rel_path =       usgs_30s:landuse_30s/
                rel_path =     usgs_lakes:landuse_30s_with_lakes/
                rel_path =    modis_lakes:modis_landuse_21class_30s/
                rel_path =     modis_2010:slucm/modis_landuse_20class_15s_2010s/
                rel_path = modis_2010_lake:slucm/modis_landuse_20class_15s_2010s_with_lakes/
                rel_path =     modis_2050:slucm/modis_landuse_20class_15s_2050s/
                rel_path = modis_2050_lake:slucm/modis_landuse_20class_15s_2050s_with_lakes/
                rel_path =        usgs_2m:landuse_2m/
                rel_path =        usgs_5m:landuse_5m/
                rel_path =       usgs_10m:landuse_10m/
                rel_path =         lowres:modis_landuse_20class_5m_with_lakes/
                rel_path =        default:modis_landuse_20class_30s_with_lakes/
        ===============================
    • index: I am attaching the index file in case it might be useful for debugging.
      • type=categorical
        category_min=1
        category_max=21
        projection=regular_ll
        dx=0.00083333
        dy=0.00083333
        known_x=1.0
        known_y=1.0
        known_lat=34.0245836581763
        known_lon=126.142083298846
        wordsize=1
        row_order=top_bottom
        tile_x=1110
        tile_y=1000
        tile_z=1
        missing_value=0
        units="category"
        description="EGIS Landuse"
        mminlu="MODIFIED_IGBP_MODIS_NOAH"
        iswater=17
        islake=21
        isice=15
        isurban=13
Please let me know if you need any additional information. Thank you again for your support.

Best regards,
aba0327
 

Attachments

  • geogrid.log
    86.8 KB · Views: 1
  • namelist.wps
    950 bytes · Views: 1
Thanks for sending those. I'm not positive if this is causing the issue, but I believe the format on the GEOGRID.TBL is important, meaning that you will need to add the divider between the two LANDUSEF entries in your table. Currently, it's listed as:

Code:
===============================
name=LANDUSEF
        priority = 2
        dest_type = categorical
        z_dim_name = land_cat
        dominant = LU_INDEX
        landmask_water =  egis_landuse:17,21      # Calculate a landmask from this field
        interp_option =  egis_landuse:nearest_neighbor
        rel_path =   egis_landuse:egis_landuse/
name=LANDUSEF
        priority = 1
        dest_type = categorical
        z_dim_name = land_cat
        dominant = LU_INDEX
        landmask_water =    nlcd2006_9s:17         # Calculate a landmask from this field
        landmask_water =   nlcd2006_30s:17         #
        ...

Try adding the divider, as is shown below:

Code:
===============================
name=LANDUSEF
        priority = 2
        dest_type = categorical
        z_dim_name = land_cat
        dominant = LU_INDEX
        landmask_water =  egis_landuse:17,21      # Calculate a landmask from this field
        interp_option =  egis_landuse:nearest_neighbor
        rel_path =   egis_landuse:egis_landuse/
===============================
name=LANDUSEF
        priority = 1
        dest_type = categorical
        z_dim_name = land_cat
        dominant = LU_INDEX
        landmask_water =    nlcd2006_9s:17         # Calculate a landmask from this field
        landmask_water =   nlcd2006_30s:17         #
        ...

and then rerun geogrid and see if there is any difference. If not, will you package up your EGIS data and share it with me so I can test it out? If it's too large to attach to this thread, see the forum home page for instructions on sharing large files. Thanks!
 
Top