search

CEMPD / VERDI

This is the repo for the VERDI project, written in java.
GNU General Public License v3.0
17 stars 13 forks source link

VERDI shapefile export of GRIDCRO2D and GRIDDOT2D variables from 2016_12SE1 mcip output #225

Closed lizadams closed 3 years ago

lizadams commented 3 years ago

Describe the bug This bug is related to #223, but demonstrates it using different input datasets and Harker's Island in NC as geographic reference point.

To Reproduce Steps to reproduce the behavior:

  1. Load datasets into verdi verdi.sh -f $cwd/GRIDCRO2D_160702.nc -f $cwd/GRIDDOT2D_160702.nc -s "DLUSE[1]" -g tile -s "LATU[2]" -g tile

  2. Zoom into the region over Harker's island for both plots Click on Controls > Set Row and Column Ranges > Select Rows 59 .. 60, Select Columns 89 .. 89.

  3. Use the Externalize this view into a Floating Window to place the plots side by side.

    Screen Shot 2021-03-04 at 11 06 01 AM

  4. Use File > Export as Image/GIS > Select the File Format: Shapefile Specify the filenames: LATU_GRIDDOT2D_160702_zoomed_in_Rows59_60_Col_89_89_island_near_Hatterous DLUSE_GRIDCRO2D_160702_Zoomed_in_ROW59_60_COL_89_89

  5. Load files into QGIS Layer > Add Vector Layer > Categorize the variable and zoom into layer

QGIS_LATU_GRIDCRO2D_and_LATU_GRIDDOT2D

A zoomed out image in QGIS showing all of the outer banks, and the full extent of both exported grid cells.

QGIS_showing_LATU_GRIDDOT2D_shifted_to_West

Expected behavior The Grid Cell boundary between rows 59 and 60 for the LATU_GRIDDOT2D file should be directly over Harker's Island according to the VERDI Tile Plot Display. The exported shapefile shows this boundary shifted to the West and perhaps a bit South as well.

Using QGIS, examining the Properties > Information, the following information including the Extent for each of the files was obtained.

Screen Shot 2021-03-04 at 12 32 43 PM Screen Shot 2021-03-04 at 12 32 29 PM
lizadams commented 3 years ago

I also plotted the LANDMASK variable from the WRF output used to generate the MCIP output. I believe this also shows that the exported shapefile is shifted to the south west. Used VERDI > Control > Show Lat/Lon Then put my cursor at the lower left corner of the grid cell to record the lat/lon values displayed by VERDI. I then compared that with the lower left corner of the grid cell lat/lon values displayed in QGIS. The exported shapefile is shifted to the west of where it should be according to the VERDI values of the llc.

Screen Shot 2021-03-04 at 4 10 49 PM Screen Shot 2021-03-04 at 4 37 29 PM

It would be easier to compare the LWMASK and LANDMASK with the geographical features, if the grid resolution were finer. These files are 12km grid resolution, not fine enough for correlating the data with the islands.

lizadams commented 3 years ago

In the following documentation they mention that the corner_lats and corner_lons global attributes in a WPS-generated GEOGRID file must be known: https://ral.ucar.edu/sites/default/files/public/WRFHydro_GIS_Preprocessor_v5.pdf

"Compatibility of the pre-processing tools with subsetted domainsMany programs such as NCO and NCL will allow a user to subset a 2D netCDF file. However, only files that are generated using the WPS geogrid.exe program may be used as input to the WRF-Hydro ArcGIS Pre-processing tools. This is because the georeferencing of GEOGRID files requires the corner latitude and longitude to be known. This is stored in the ‘corner_lats’ and ‘corner_lons’ global attributes in a WPS-generated GEOGRID file. If a user wishes to subset a WRF domain, they must also update the values in these attributes in order for the domain to be properly georeferenced."

The ncdump GRIDCRO2D_04OTC2_160930.nc file contains the LAT and LON values a the cell centers. LON:var_desc = "longitude at cell centers

The following links describes the GEOGRID file, and lat/lon values for the MASS, U, V, and unstaggered variables. https://www2.mmm.ucar.edu/wrf/users/docs/user_guide_V3/users_guide_chap3.htm#_WPS_Output_Fields

yadongxuEPA commented 3 years ago

Retested with VERDI_2.0_linux64_20210427.tar.gz, the zoomed in tile plots in VERDI is shown as below: Issue_225_VERDI_testing_04272021_in_VERDI The corresponding exported shapefiles in QGIS is shown as below: issue_225_VERDI_testing_04272021_in_QGIS_1 issue_225_VERDI_testing_04272021_in_QGIS_2 We can see that the grid cell boundary between rows 59 and 60 for the LATU_GRIDDOT2D file is directly over Harker's Island, which is consistent with the location shown in the VERDI Tile Plot Display.

The projection information in the exported shapefiles is changed into: PROJCS["unknown", GEOGCS["GCS_unknown", DATUM["D_unknown", SPHEROID["unknown", 6370000.0, 0.0]], PRIMEM["Greenwich", 0.0], UNIT["degree", 0.017453292519943295], AXIS["Longitude", EAST], AXIS["Latitude", NORTH]], PROJECTION["Lambert_Conformal_Conic"], PARAMETER["central_meridian", -97.0], PARAMETER["latitude_of_origin", 40.0], PARAMETER["standard_parallel_1", 45.0], PARAMETER["false_easting", 0.0], PARAMETER["false_northing", 0.0], PARAMETER["scale_factor", 1.0], PARAMETER["standard_parallel_2", 33.0], UNIT["m", 1.0], AXIS["x", EAST], AXIS["y", NORTH]]

yadongxuEPA commented 3 years ago

Retested with VERDI_2.1_linux64_20210906.tar.gz, confirmed that the upgraded openJDK has not caused any problem on this issue. issue_225_VERDI_testing_20210906_in_QGIS_1