Polygons for marsh migration under sea-level rise in Chesapeake Bay

This process step and subsequent process steps were performed by the same person, Zafer Defne, in ArcGIS Pro (ver.3.3.2) unless otherwise stated. For complex operations, names of specific tools used are given in CAPITAL letters (any critical parameters used are given in parentheses, separated by a semicolon, immediately after the tool name). The input and output file names are provided in [square brackets] when necessary. Units for length and area calculations are meters (m) and square meters (m2) unless otherwise stated.
This processing step defines the extent of migration area.
a) Download the NOAA Sea level rise wetland impacts and migration datasets for the Chesapeake Bay domain [md_marshmigration_2016.zip, va_marshmigration_2016.zip, nc_marshmigration_2016.zip]. These datasets present the potential distribution of each wetland type based on their elevation and frequency of inundation under different SLR rates.
b) From each scenario select classes of Palustrine Emergent Wetland (15), Brackish/Transition Wetland (17), Estuarine Wetland (18) by EXTRACT_BY_ATTRIBUTES(VALUE = 15 Or VALUE = 17 Or VALUE = 18).
c) Mosaic all four rasters from previous step into to a new raster for Maryland (MD) while setting the coordinate system to NAD1983 (2011) UTM Zone 18N. MOSAIC_TO_NEW_RASTER(Pixel Type="1bit"; Spatial Reference=NAD_1983_2011_UTM_Zone_18N; Mosaic Operator="Maximum"). Export only the values greater than 1 as a new raster and convert to polygon. RASTER_TO_POLYGON(Field = "VALUE").
Repeat the same steps for Virginia (VA) and North Carolina (NC) domains. And merge the results in to a single feature dataset [mm_MD_VA_NC_merge_UTM18_2011.shp].

This processing step delineates the migration area based on the USGS National Hydrography Dataset Plus High Resolution (NHD Plus HR) catchments map. First, the catchments intersecting the marsh units are selected as immediate catchments. Then, if the inland boundary of the migration layer spans across multiple catchments, those catchments are merged to the immediate catchments by proximity. Catchments areas beyond the inland boundary of migration layer are clipped to the boundary. The resulting product is the marsh migration polygons which is the migration zone delineated with catchments.
a) Download Chesapeake Bay marsh units dataset [mu_UVVR_CB.shp] and NHD Plus HR catchments datasets for the domain. [mu_UVVR_CB.shp]; NHDPlusCatchment feature datasets [NHDPlusCatchment_0206.shp], [NHDPlusCatchment_0207.shp], [NHDPlusCatchment_0208.shp] and [NHDPlusCatchment_0301.shp].
b) Project all source feature datasets to the same coordinate system with the previous processing step, NAD1983 (2011) UTM Zone 18N. Resulting projected data sets are: Marsh units [mu_UVVR_CB_UTM18_2011.shp]; NHDPlusCatchment feature datasets [NHDPlusCatchment_0206_UTM18_2011.shp], [NHDPlusCatchment_0207_UTM18_2011.shp], [NHDPlusCatchment_0208_UTM18_2011.shp] and [NHDPlusCatchment_0301_UTM18_2011.shp]
c) Merge all catchment feature datasets to create the total catchment area [HDPlusCatchment_merge_0206_0207_0208_0301.shp].
d) Select parts of migration area that intersect with marsh units. SELECT_LAYER_BY_LOCATION(Input features=mm_MD_VA_NC_merge_UTM18_2011.shp; Selecting features=mu_UVVR_CB_UTM18_2011.shp; Relationship=Intersect; Search distance =1 meter). Export features as [mmp_by_mu.shp].
e) SELECT_LAYER_BY_LOCATION(Input features=HDPlusCatchment_merge_0206_0207_0208_0301.shp; Selecting features= mmp_by_mu.shp, Relationship=Intersect). Export features as [cathments_selected.shp].
f) Intersect the two feature datasets to identify the immediate catchment polygons. PAIRWISE_INTERSECT(Input features=catchments_selected.shp, mmp_by_mu.shp ). Export features as [mmp_iniside_cat.shp].
g) Remove the immediate catchments from the overall marsh migration area so that the catchments in the migration area but beyond the immediate ones can be processed. PAIRWISE_ERASE(Input features=mm_MD_VA_NC_merge_UTM18_2011.shp; Erase features=mmp_inside_cat.shp; Output=mm_outside_cat.shp). MULTIPART_TO_SINGLEPART(Input features=mm_outside_cat.shp; Output features=mm_outside_cat_single.shp) and SELECT_LAYER_BY_LOCATION(Input features=mm_outside_cat_single.shp; Selecting features=mu_UVVR_CB_UTM18_2011.shp; Relationship=Intersect; Search distance =1 meter). Export features as [mmp_outside_cat.shp].
h) Union the immediate catchments with the remote catchments to create a union set of catchments. UNION(Input features=mmp_iniside_cat.shp, mmp_outside_cat.shp; Gaps allowed=Yes, Output features=mmp_all.shp).
i) Prepare and split the union of catchments to immediate catchments and remote catchments so that the remote catchments can be merged by proximity. Specifically, ERASE marsh units from the union of catchments and MULTIPART_TO_SINGLEPART to obtain [mmp_all_single.shp]. Create two feature datasets from this dataset: For immediate catchments SELECT_LAYER_BY_LOCATION(Input features=mmp_all_single.shp; Selecting features=mu_UVVR_CB_UTM18_2011.shp; Relationship=Intersect; Search distance=0.001 meter) and export as [input_mmu.shp]. For remote catchments SELECT_LAYER_BY_LOCATION(Input features=mmp_all_single.shp; Selecting features=mu_UVVR_CB_UTM18_2011.shp; Relationship=Intersect; Search distance=0.001 meter, Invert Spatial Relationship=Yes) and export as [near_mmu.shp].
j) Run the HYDUNITLOOP tool to merge catchments. At each iteration step, the script uses NEAR tool to find remote catchments from the [input_mmu.shp] within 1 meter of an immediate catchment [near_mmu.shp] and merges them to it using UNION and DISSOLVE tools. The catchments merged with an immediate catchment are removed from the remote catchments dataset. The script iterates until there is no change in the number of remote and immediate units. The remaining catchments are those that are more than 1 meter away from any immediate ones and they are discarded. The resulting final dataset is [final_mmu.shp]. HYDUNITLOOP(Input features=input_mmu.shp; Near features=near_mmu.shp; Output features=final_mmu.shp].
k) Because of the resolution difference between the datasets, processing sometimes creates sliver polygons around the edges. To clean up the final output, remove detached marsh migration polygons with area smaller than a threshold value. Specifically, first PAIRWISE_DISSOLVE(Input features=final_mmu.shp; Output feature class=final_mmu_diss.shp) and CALCULATE_GEOMETRY_ATTRIBUTES(Input features=final_mmu_diss.shp; Geometry attributes field=ASQM; Property=Area (geodesic); Area unit=Square meters; Coordinate system=NAD_1983_2011_Contiguous_USA_Albers). Then, select polygons with areas greater than 300 square meter to export [fmmu_gt_300.shp] and use them to select from [final_mmu.shp] dataset by SELECT_BY_LOCATION(Input feature=final_mmu.shp; Relationship=Within; Selecting features=fmmu_gt_300.shp; Search distance=0 meters). Finally, export the selected features as [marsh_migration_polygons.shp]. l) Keep only the field ASQM and project to Web Mercator projection to obtain [marsh_mig_poly_CB.shp].