Change global_ocean ISC test cases to use BedMachine

testing_and_setup/compass/ocean/global_ocean/CUSP12/spin_up/config_driver.xml

@@ -1,4 +1,4 @@
-<driver_script name="run_test.py">
+<driver_script name="run.py">
 	<case name="spin_up1">
 		<step executable="./run.py" quiet="true" pre_message=" * Running spin_up1" post_message=" - Complete"/>
 	</case>

testing_and_setup/compass/ocean/global_ocean/CUSP8/spin_up/config_driver.xml

@@ -1,4 +1,4 @@
-<driver_script name="run_test.py">
+<driver_script name="run.py">
 	<case name="spin_up1">
 		<step executable="./run.py" quiet="true" pre_message=" * Running spin_up1" post_message=" - Complete"/>
 	</case>

testing_and_setup/compass/ocean/global_ocean/EC60to30/spin_up/config_driver.xml

@@ -1,4 +1,4 @@
-<driver_script name="run_test.py">
+<driver_script name="run.py">
 	<case name="prep_spin_up1">
 		<step executable="./run.py" quiet="true" pre_message=" * Running prep_spin_up1" post_message=" - Complete"/>
 	</case>

testing_and_setup/compass/ocean/global_ocean/EC60to30wISC/init/config_initial_state.xml

@@ -0,0 +1,59 @@
+<?xml version="1.0"?>
+<config case="initial_state">
+	<get_file dest_path="initial_condition_database" file_name="PotentialTemperature.01.filled.60levels.PHC.151106.nc">
+		<mirror protocol="wget" url="https://web.lcrc.anl.gov/public/e3sm/mpas_standalonedata/mpas-ocean/initial_condition_database/"/>
+	</get_file>
+
+	<get_file dest_path="initial_condition_database" file_name="Salinity.01.filled.60levels.PHC.151106.nc">
+		<mirror protocol="wget" url="https://web.lcrc.anl.gov/public/e3sm/mpas_standalonedata/mpas-ocean/initial_condition_database/"/>
+	</get_file>
+
+
+	<get_file dest_path="initial_condition_database" file_name="windStress.ncep_1958-2000avg.interp3600x2431.151106.nc">
+		<mirror protocol="wget" url="https://web.lcrc.anl.gov/public/e3sm/mpas_standalonedata/mpas-ocean/initial_condition_database/"/>
+	</get_file>
+
+	<get_file dest_path="bathymetry_database" file_name="BedMachineAntarctica_and_GEBCO_2019_0.05_degree.200128.nc">
+		<mirror protocol="wget" url="https://web.lcrc.anl.gov/public/e3sm/mpas_standalonedata/mpas-ocean/bathymetry_database/"/>
+	</get_file>
+
+	<get_file dest_path="initial_condition_database" file_name="chlorophyllA_monthly_averages_1deg.151201.nc">
+		<mirror protocol="wget" url="https://web.lcrc.anl.gov/public/e3sm/mpas_standalonedata/mpas-ocean/initial_condition_database/"/>
+	</get_file>
+
+
+	<add_link source="../culled_mesh/culled_mesh.nc" dest="mesh.nc"/>
+	<add_link source="../culled_mesh/culled_graph.info" dest="graph.info"/>
+	<add_link source="../culled_mesh/critical_passages_mask_final.nc" dest="critical_passages.nc"/>
+    <add_link source_path="initial_condition_database" source="Salinity.01.filled.60levels.PHC.151106.nc" dest="layer_depth.nc"/>
+	<add_link source_path="initial_condition_database" source="PotentialTemperature.01.filled.60levels.PHC.151106.nc" dest="temperature.nc"/>
+	<add_link source_path="initial_condition_database" source="Salinity.01.filled.60levels.PHC.151106.nc" dest="salinity.nc"/>
+	<add_link source_path="initial_condition_database" source="windStress.ncep_1958-2000avg.interp3600x2431.151106.nc" dest="wind_stress.nc"/>
+	<add_link source_path="bathymetry_database" source="BedMachineAntarctica_and_GEBCO_2019_0.05_degree.200128.nc" dest="topography.nc"/>
+	<add_link source_path="initial_condition_database" source="chlorophyllA_monthly_averages_1deg.151201.nc" dest="swData.nc"/>
+
+	<namelist name="namelist.ocean" mode="init">
+		<template file="template_initial_state.xml" path_base="script_configuration_dir"/>
+		<template file="template_critical_passages.xml" path_base="script_core_dir" path="global_ocean"/>
+		<template file="template_init_with_land_ice.xml" path_base="script_configuration_dir"/>
+		<option name="config_iterative_init_variable">'landIcePressure_from_top_density'</option>
+		<option name="config_global_ocean_depth_file">'layer_depth.nc'</option>
+		<option name="config_global_ocean_minimum_depth">10</option>
+		<option name="config_global_ocean_deepen_critical_passages">.false.</option>
+        	<option name="config_global_ocean_depth_conversion_factor">0.01</option>
+        	<option name="config_global_ocean_tracer_depth_conversion_factor">0.01</option>
+	</namelist>
+
+	<streams name="streams.ocean" keep="immutable" mode="init">
+		<template file="template_initial_state.xml" path_base="script_configuration_dir"/>
+		<template file="template_critical_passages.xml" path_base="script_core_dir" path="global_ocean"/>
+	</streams>
+
+	<run_script name="run.py">
+		<step executable="gpmetis">
+			<argument flag="graph.info">4</argument>
+		</step>
+		<model_run procs="4" threads="1" namelist="namelist.ocean" streams="streams.ocean"/>
+    </run_script>
+
+</config>

testing_and_setup/compass/ocean/global_ocean/EC60to30wISC/spin_up/config_driver.xml

@@ -1,4 +1,4 @@
-<driver_script name="run_test.py">
+<driver_script name="run.py">
 	<case name="prep_spin_up1">
 		<step executable="./run.py" quiet="true" pre_message=" * Running prep_spin_up1" post_message=" - Complete"/>
 	</case>

testing_and_setup/compass/ocean/global_ocean/QU60/spin_up/config_driver.xml

@@ -1,4 +1,4 @@
-<driver_script name="run_test.py">
+<driver_script name="run.py">
 	<case name="prep_spin_up1">
 		<step executable="./run.py" quiet="true" pre_message=" * Running prep_spin_up1" post_message=" - Complete"/>
 	</case>

testing_and_setup/compass/ocean/global_ocean/SO60to10wISC/init/config_base_mesh.xml

@@ -0,0 +1,13 @@
+<?xml version="1.0"?>
+<config case="base_mesh">
+
+	<add_link source_path="mpas_model" source="testing_and_setup/compass/ocean/jigsaw_to_MPAS" dest="jigsaw_to_MPAS"/>
+	<add_link source_path="script_test_dir" source="." dest="define_base_mesh"/>
+	<add_link source_path="script_test_dir" source="high_res_region.geojson" dest="high_res_region.geojson"/>
+
+	<run_script name="run.py">
+		<step executable="python">
+			<argument flag="-m">jigsaw_to_MPAS.build_mesh</argument>
+		</step>
+	</run_script>
+</config>

testing_and_setup/compass/ocean/global_ocean/SO60to10wISC/init/define_base_mesh.py

@@ -1,39 +1,38 @@
-# /usr/bin/env python
-"""
-% Create cell width array for this mesh on a regular latitude-longitude grid.
-% Outputs:
-%    cellWidth - m x n array, entries are desired cell width in km
-%    lat - latitude, vector of length m, with entries between -90 and 90, degrees
-%    lon - longitude, vector of length n, with entries between -180 and 180, degrees
-"""
 import numpy as np
 import jigsaw_to_MPAS.mesh_definition_tools as mdt
+from jigsaw_to_MPAS.coastal_tools import signed_distance_from_geojson, \
+    compute_cell_width
+from geometric_features import read_feature_collection
+import xarray
 
 # Uncomment to plot the cell size distribution.
-#import matplotlib
-#matplotlib.use('Agg')
-#import matplotlib.pyplot as plt
+# import matplotlib
+# matplotlib.use('Agg')
+# import matplotlib.pyplot as plt
 
 
 def cellWidthVsLatLon():
-    lat = np.arange(-90, 90.01, 0.1)
-    # Note that longitude step is 10 degrees, but should only be used if mesh does not 
-    # vary with longitude. Otherwise, set to 0.1 degrees.
-    lon = np.arange(-180, 180.01, 10.0)
-
-    # define uniform distributions
-    cellWidth10 = 10.0 * np.ones(lat.size)
-    cellWidth30 = 30.0 * np.ones(lat.size)
-
-    # Southern transition
-    latTransition = -48.0
-    latWidthTransition = 10.0
-    cellWidthSouth = mdt.mergeCellWidthVsLat(
-        lat,
-        cellWidth10,
-        cellWidth30,
-        latTransition,
-        latWidthTransition)
+    """
+    Create cell width array for this mesh on a regular latitude-longitude grid.
+    Returns
+    -------
+       cellWidth : numpy.ndarray
+            m x n array, entries are desired cell width in km
+       lat : numpy.ndarray
+            latitude, vector of length m, with entries between -90 and 90,
+            degrees
+       lon : numpy.ndarray
+            longitude, vector of length n, with entries between -180 and 180,
+            degrees
+    """
+    dlon = 0.1
+    dlat = dlon
+    nlon = int(360./dlon) + 1
+    nlat = int(180./dlat) + 1
+    lon = np.linspace(-180., 180., nlon)
+    lat = np.linspace(-90., 90., nlat)
+
+    cellWidthSouth = 30. * np.ones((len(lat)))
 
     # Transition at Equator
     cellWidthNorth = mdt.EC_CellWidthVsLat(lat)
@@ -46,17 +45,44 @@ def cellWidthVsLatLon():
         latTransition,
         latWidthTransition)
 
+    _, cellWidth = np.meshgrid(lon, cellWidthVsLat)
+
+    # now, add the high-res region
+    fc = read_feature_collection('high_res_region.geojson')
+
+    signed_distance = signed_distance_from_geojson(fc, lon, lat,
+                                                   max_length=0.25)
+
+    da = xarray.DataArray(signed_distance,
+                          dims=['y', 'x'],
+                          coords={'y': lat, 'x': lon},
+                          name='signed_distance')
+    cw_filename = 'signed_distance.nc'
+    da.to_netcdf(cw_filename)
+
+    # multiply by 5 because transition_width gets multiplied by 0.2 in
+    # compute_cell_width
+    # Equivalent to 10 degrees latitude
+    trans_width = 5*1100e3
+    # The last term compensates for the offset in compute_cell_width.
+    # The middle of the transition is ~2.5 degrees (300 km) south of the
+    # region boundary to best match previous transition at 48 S. (The mean lat
+    # of the boundary is 45.5 S.)
+    trans_start = -300e3 - 0.5 * trans_width
+    dx_min = 10.
+
+    cellWidth = compute_cell_width(signed_distance, cellWidth, lon,
+                                   lat, dx_min, trans_start, trans_width,
+                                   restrict_box={'include': [], 'exclude': []})
+
     # Uncomment to plot the cell size distribution.
-    #plt.plot(lat,cellWidthVsLat)
-    #plt.grid(True)
-    #plt.xlabel('latitude')
-    #plt.ylabel('grid cell size')
-    #plt.title('SO60to10wISC, transition at 55S')
-    #plt.savefig('cellWidthVsLat.pdf')
-    #plt.savefig('cellWidthVsLat.png')
-
-    cellWidth = np.ones((lat.size, lon.size))
-    for i in range(lon.size):
-        cellWidth[:, i] = cellWidthVsLat
+    # Lon, Lat = np.meshgrid(lon, lat)
+    # ax = plt.subplot(111)
+    # plt.pcolormesh(Lon, Lat, cellWidth)
+    # plt.colorbar()
+    # ax.set_aspect('equal')
+    # ax.autoscale(tight=True)
+    # plt.tight_layout()
+    # plt.savefig('cellWidthVsLat.png', dpi=200)
 
     return cellWidth, lon, lat

testing_and_setup/compass/ocean/global_ocean/SO60to10wISC/init/high_res_region.geojson

@@ -0,0 +1,73 @@
+{
+  "type": "FeatureCollection",
+  "features": [
+    {
+      "type": "Feature",
+      "properties": {
+        "name": "SO60to10 high res region",
+        "component": "ocean",
+        "object": "region",
+        "author": "Xylar Asay-Davis"
+      },
+      "geometry": {
+        "type": "Polygon",
+        "coordinates": [
+          [
+            [
+              -75.5859375,
+              -48.92249926375823
+            ],
+            [
+              -96.50390625,
+              -54.67383096593114
+            ],
+            [
+              -124.45312499999999,
+              -53.95608553098789
+            ],
+            [
+              -180,
+              -53
+            ],
+            [
+              -180,
+              -90
+            ],
+            [
+              180,
+              -90
+            ],
+            [
+              180,
+              -53
+            ],
+            [
+              168.3984375,
+              -51.17934297928927
+            ],
+            [
+              121.640625,
+              -45.82879925192133
+            ],
+            [
+              22.148437499999996,
+              -37.99616267972812
+            ],
+            [
+              -61.17187499999999,
+              -34.30714385628803
+            ],
+            [
+              -68.90625,
+              -40.97989806962013
+            ],
+            [
+              -75.5859375,
+              -48.92249926375823
+            ]
+          ]
+        ]
+      }
+    }
+  ]
+}

testing_and_setup/compass/ocean/global_ocean/SO60to10wISC/spin_up/config_driver.xml

@@ -1,4 +1,4 @@
-<driver_script name="run_test.py">
+<driver_script name="run.py">
 	<case name="prep_spin_up1">
 		<step executable="./run.py" quiet="true" pre_message=" * Running prep_spin_up1" post_message=" - Complete"/>
 	</case>

testing_and_setup/compass/ocean/global_ocean/SOQU60to15/spin_up/config_driver.xml

@@ -1,4 +1,4 @@
-<driver_script name="run_test.py">
+<driver_script name="run.py">
 	<case name="spin_up1">
 		<step executable="./run.py" quiet="true" pre_message=" * Running spin_up1" post_message=" - Complete"/>
 	</case>

testing_and_setup/compass/ocean/global_ocean/config_files_ISC/config_initial_state.xml

@@ -16,8 +16,8 @@
 		<mirror protocol="wget" url="https://web.lcrc.anl.gov/public/e3sm/mpas_standalonedata/mpas-ocean/initial_condition_database/"/>
 	</get_file>
 
-	<get_file dest_path="initial_condition_database" file_name="landIceTopo_0.1deg_smooth_0.1deg.nc">
-		<mirror protocol="wget" url="https://web.lcrc.anl.gov/public/e3sm/mpas_standalonedata/mpas-ocean/initial_condition_database/"/>
+	<get_file dest_path="bathymetry_database" file_name="BedMachineAntarctica_and_GEBCO_2019_0.05_degree.200128.nc">
+		<mirror protocol="wget" url="https://web.lcrc.anl.gov/public/e3sm/mpas_standalonedata/mpas-ocean/bathymetry_database/"/>
 	</get_file>
 
 	<get_file dest_path="initial_condition_database" file_name="chlorophyllA_monthly_averages_1deg.151201.nc">
@@ -33,7 +33,7 @@
 	<add_link source_path="initial_condition_database" source="PTemp.Jan_p3.filled.mpas100levs.160127.nc" dest="temperature.nc"/>
 	<add_link source_path="initial_condition_database" source="Salt.Jan_p3.noBlackCaspian.filled.mpas100levs.160127.nc" dest="salinity.nc"/>
 	<add_link source_path="initial_condition_database" source="windStress.ncep_1958-2000avg.interp3600x2431.151106.nc" dest="wind_stress.nc"/>
-	<add_link source_path="initial_condition_database" source="landIceTopo_0.1deg_smooth_0.1deg.nc" dest="topography.nc"/>
+	<add_link source_path="bathymetry_database" source="BedMachineAntarctica_and_GEBCO_2019_0.05_degree.200128.nc" dest="topography.nc"/>
 	<add_link source_path="initial_condition_database" source="chlorophyllA_monthly_averages_1deg.151201.nc" dest="swData.nc"/>
 
 	<namelist name="namelist.ocean" mode="init">

testing_and_setup/compass/ocean/global_ocean/config_files_ISC/template_initial_state.xml

@@ -0,0 +1,11 @@
+<template>
+	<streams>
+		<stream name="land_ice_mask_init">
+			<attribute name="type">input</attribute>
+			<attribute name="filename_template">land_ice_mask.nc</attribute>
+			<add_contents>
+				<member name="landIceMask" type="var"/>
+			</add_contents>
+		</stream>
+	</streams>
+</template>

testing_and_setup/compass/ocean/global_ocean/cull_mesh.py

@@ -63,11 +63,10 @@
 # Add "land" coverage from either the full ice sheet or just the grounded
 # part
 if options.with_cavities:
-    fcAntarcticLand = gf.read(componentName='bedmap2', objectType='region',
+    fcAntarcticLand = gf.read(componentName='bedmachine', objectType='region',
                               featureNames=['AntarcticGroundedIceCoverage'])
-
 else:
-    fcAntarcticLand = gf.read(componentName='bedmap2', objectType='region',
+    fcAntarcticLand = gf.read(componentName='bedmachine', objectType='region',
                               featureNames=['AntarcticIceCoverage'])
 
 fcLandCoverage.merge(fcAntarcticLand)
@@ -145,8 +144,9 @@
                  format=netcdfFormat)
 
 if options.with_cavities:
-    fcAntarcticIce = gf.read(componentName='bedmap2', objectType='region',
+    fcAntarcticIce = gf.read(componentName='bedmachine', objectType='region',
                              featureNames=['AntarcticIceCoverage'])
+    fcAntarcticIce.to_geojson('ice_coverage.geojson')
     dsMask = conversion.mask(dsCulledMesh, fcMask=fcAntarcticIce)
     landIceMask = dsMask.regionCellMasks.isel(nRegions=0)
     dsLandIceMask = xarray.Dataset()