Show bathymetry plot in notebook + some other fixes

demos/access_om2-forced.ipynb

@@ -607,9 +607,10 @@
     "## Directory where ocean model cut-outs go before processing\n",
     "tmp_dir = f\"{gdata}/{expt_name}\"\n",
     "\n",
-    "for i in [run_dir, tmp_dir, input_dir]:\n",
-    "    if not os.path.exists(i):\n",
-    "        os.makedirs(str(i))"
+    "## if directories don't exist, create them\n",
+    "for path in (run_dir, tmp_dir, input_dir):\n",
+    "    if not os.path.exists(path):\n",
+    "        os.makedirs(path)"
    ]
   },
   {
@@ -737,7 +738,7 @@
    "source": [
     "## Step 4: Set up bathymetry\n",
     "\n",
-    "Similarly to ocean forcing, we point our 'bathymetry' method at the location of the file of choice, and pass it a dictionary mapping variable names. This time we don't need to preprocess the topography since it's just a 2D field and easier to deal with. Afterwards you can run `expt.topog` and have a look at your domain. "
+    "Similarly to ocean forcing, we point the `setup_bathymetry` method at the location of the file of choice. In this case, we don't need to preprocess the bathymetry since it's just a 2D field and easier to deal with. Afterwards, the bathymetry is stored in `expt.bathymetry`."
    ]
   },
   {
@@ -812,7 +813,7 @@
     }
    ],
    "source": [
-    "expt.topog.depth.plot()"
+    "expt.bathymetry.depth.plot()"
    ]
   },
   {
@@ -823,7 +824,7 @@
     "\n",
     "This cuts out and interpolates the initial condition as well as all boundaries (unless you don't pass it boundaries).\n",
     "\n",
-    "The dictionary maps the MOM6 variable names to what they're called in your ocean input file. Notice how the horizontal dimensions are `x` and `y` in the GLORYS reanalysis example, vs `xh`, `yh`, `xq`, `yq`. This is because ACCESS-OM2-01 is on a `B` grid, so we need to differentiate between `q` and `t` points. \n",
+    "The dictionary maps the MOM6 variable names to what they're called in your ocean input file. Notice how the horizontal dimensions are `x` and `y` in the GLORYS reanalysis example, vs `xh`, `yh`, `xq`, `yq`. This is because ACCESS-OM2-01 is on a `B` grid, so we need to differentiate between the `q` and `t` points. \n",
     "\n",
     "If one of your segments is land, you can delete its string from the 'boundaries' list. You'll need to update `MOM_input` to reflect this though so it knows how many segments to look for, and their orientations. "
    ]
@@ -1055,7 +1056,7 @@
     "\n",
     "Hopefully your model is running. If not, the first thing you should do is reduce the timestep. You can do this by adding `#override DT=XXXX` to your `MOM_override` file. \n",
     "\n",
-    "If there's strange behaviour on your boundaries, you could play around with the `nudging timescale` (an example is already included in the `MOM_override` file). Sometimes, if your boundary has a lot going on (like all of the eddies spinning off the ACC), it can be hard to avoid these edge effects. This is because the chaotic, submesoscale structures developed within the regional domain won't match those at the boundary. "
+    "If there's strange behaviour on your boundaries, you could play around with the `nudging timescale` (an example is already included in the `MOM_override` file). Sometimes, if your boundary has a lot going on (like all of the eddies spinning off the western boundary currents or off the Antarctic Circumpolar current), it can be hard to avoid these edge effects. This is because the chaotic, submesoscale structures developed within the regional domain won't match those at the boundary. "
    ]
   }
  ],

regional_mom6/regional_mom6.py

@@ -296,7 +296,7 @@ def rectangular_hgrid(λ, φ):
         np.diff(λ), dλ * np.ones(np.size(λ) - 1)
     ), "provided array of longitudes must be uniformly spaced"
 
-    # dx = R * cos(φ) * np.deg2rad(dλ) / 2
+    # dx = R * cos(np.deg2rad(φ)) * np.deg2rad(dλ) / 2
     # Note: division by 2 because we're on the supergrid
     dx = np.broadcast_to(
         R * np.cos(np.deg2rad(φ)) * np.deg2rad(dλ) / 2,
@@ -364,7 +364,7 @@ class experiment:
 
     Methods in this class generate the various input files needed for a MOM6
     experiment forced with open boundary conditions (OBCs). The code is agnostic
-    to the user's choice of boundary forcing, bathymetry and surface forcing;
+    to the user's choice of boundary forcing, bathymetry, and surface forcing;
     users need to prescribe what variables are all called via mapping dictionaries
     from MOM6 variable/coordinate name to the name in the input dataset.
 
@@ -825,7 +825,8 @@ def initial_condition(
         eta_out.yh.attrs = ic_raw_tracers.lat.attrs
         eta_out.attrs = ic_raw_eta.attrs
 
-        # if temp units are K, convert to C
+        ## if min(temp) > 100 then assume that units must be degrees K
+        ## (otherwise we can't be on Earth) and convert to degrees C
         if np.min(tracers_out["temp"].isel({"zl": 0})) > 100:
             tracers_out["temp"] -= 273.15
 
@@ -898,8 +899,8 @@ def rectangular_boundary(
                 ``'north'``, or ``'south'``.
             segment_number (int): Number the segments according to how they'll be specified in
                 the ``MOM_input``.
-            arakawa_grid (Optional[str]): Arakawa grid staggering of input; either ``'A'``, ``'B'``,
-                or ``'C'``.
+            arakawa_grid (Optional[str]): Arakawa grid staggering type of the boundary forcing.
+                Either ``'A'`` (default), ``'B'``, or ``'C'``.
         """
 
         print("Processing {} boundary...".format(orientation), end="")
@@ -1339,7 +1340,7 @@ def FRE_tools(self, layout=None):
             "Running GFDL's FRE Tools. The following information is all printed by the FRE tools themselves"
         )
         if not (self.mom_input_dir / "bathymetry.nc").exists():
-            print("No bathymetry file! Need to run .setup_bathymetry method first")
+            print("No bathymetry file! Need to run setup_bathymetry method first")
             return
 
         for p in self.mom_input_dir.glob("mask_table*"):
@@ -1755,7 +1756,7 @@ def __init__(
 
     def rectangular_brushcut(self):
         """
-        Cut out and interpolate tracers. This method assumes that the boundary
+        Cut out and interpolate tracers. ``rectangular_brushcut`` assumes that the boundary
         is a simple Northern, Southern, Eastern, or Western boundary.
         """
         if self.orientation == "north":
@@ -1949,7 +1950,7 @@ def rectangular_brushcut(self):
             },
         }
 
-        ### Generate our dz variable. This needs to be in layer thicknesses
+        ### Generate the dz variable; need needs to be in layer thicknesses
         dz = segment_out[self.z].diff(self.z)
         dz.name = "dz"
         dz = xr.concat([dz, dz[-1]], dim=self.z)