diff --git a/gallery/.buildinfo b/gallery/.buildinfo index 80859b295d..1da86e5851 100644 --- a/gallery/.buildinfo +++ b/gallery/.buildinfo @@ -1,4 +1,4 @@ # Sphinx build info version 1 # This file hashes the configuration used when building these files. When it is not found, a full rebuild will be done. -config: 3cf162944646a407121badf1683ef09d +config: 40297df7aaecaef04ecfc45d4de625c9 tags: 645f666f9bcd5a90fca523b33c5a78b7 diff --git a/gallery/_sources/notebooks.rst.txt b/gallery/_sources/notebooks.rst.txt index 2423188e0b..1c85dde44b 100644 --- a/gallery/_sources/notebooks.rst.txt +++ b/gallery/_sources/notebooks.rst.txt @@ -91,9 +91,13 @@ Examples using GeoClaw if they are below sea level (but protected by dikes). * `MakeFlagregionsCoast <_static/apps/notebooks/geoclaw/MakeFlagregionsCoast.html>`_ - illustrating making a ruled rectangle for use as a flagregion using the + illustrates making a ruled rectangle for use as a flagregion using the marching front algorithm. +* `IslandBuffering <_static/apps/notebooks/geoclaw/IslandBuffering.html>`_ + illustrates how to make a ruled rectangle surrounding an island with + a buffer zone that extends out some distance that is independent + of water depth. .. _notebooks_tsunami-examples: diff --git a/gallery/_static/apps/notebooks/geoclaw/IslandBuffering.html b/gallery/_static/apps/notebooks/geoclaw/IslandBuffering.html new file mode 100644 index 0000000000..e77d4d0275 --- /dev/null +++ b/gallery/_static/apps/notebooks/geoclaw/IslandBuffering.html @@ -0,0 +1,14423 @@ + + + + +IslandBuffering + + + + + + + + + + + + + + + + + + + + + + +
+
+ +
+
+
+

Flexible buffering on islands for AMR flagregions

+
+
+
+
+
+
+

This Jupyter notebook can be found in collection of Clawpack apps as the file $CLAW/apps/notebooks/geoclaw/IslandBuffering.ipynb. +To run this notebook, install Clawpack, and clone the apps repository. +A static view of this and other notebooks can be found in the Clawpack Gallery of Jupyter notebooks.

+ +
+
+
+
+
+
+

The MakeFlagRegionsCoast.ipynb and MarchingFront.ipynb notebooks illustrate tools to create AMR flagregions by choosing points along shore that satisfy elevation constraints.

+

Here, we will demonstrate another tool developed specifically for island regions. Since the topography can vary on each side of an island (for instance, barrier islands along the coast), we will explore a method to create a buffer zone surrounding the island, independent of water depth, and show how to make a Ruled Rectangle for AMR flagregions.

+

The procedure primarily involves taking the convex hull of the island region and expanding it on all sides; this will become our buffer zone. We demonstrate this in the first example. It is important to note that the expansion process requires the shape to be convex.

+

For certain shapes of islands, the convex hull will cover a region larger than we would like. In the second example, we explore a method to create a buffer zone that better approximates the island's shape. In particular, we partition the island, find the convex hull of each partition, and use the same buffering process on each convex hull.

+

To run the notebook, you will need the Python package scikit-image. Installation instructions here: https://scikit-image.org/docs/dev/install.html

+ +
+
+
+
+
+
+

Version

Written by: Mirah Shi (github @mirahshi), 2020-05-31
+Runs with Clawpack v5.8.0

+ +
+
+
+
+
+
+

Examples

+
+
+
+
+
+
+

First import some modules and set up plotting tools.

+ +
+
+
+
+
+
In [1]:
+
+
+
%matplotlib inline
+from pylab import *
+import numpy as np
+import matplotlib.pyplot as plt
+
+import netCDF4
+
+from clawpack.visclaw import plottools
+from clawpack.amrclaw import region_tools
+from clawpack.geoclaw import topotools
+from clawpack.geoclaw.util import haversine
+
+from skimage.morphology import convex_hull, convex_hull_object
+from skimage.measure import label
+
+from scipy.spatial import ConvexHull
+
+ +
+
+
+ +
+
+
+
In [2]:
+
+
+
norm = mpl.colors.Normalize(vmin=0.,vmax=1.)
+
+ +
+
+
+ +
+
+
+
+

Sample island topography

+
+
+
+
+
+
+

We crop a 1/3-arc-second DEM of Key West, available from the NCEI THREDDS server.

+ +
+
+
+
+
+
In [3]:
+
+
+
path = 'https://www.ngdc.noaa.gov/thredds/dodsC/regional/key_west_13_mhw_2011.nc'
+
+extent = [-82.17, -81.28, 24.51, 24.86]
+topo = topotools.read_netcdf(path, extent=extent)
+
+topo.plot()
+
+ +
+
+
+ +
+
+ + +
+ +
Out[3]:
+ + + + +
+
<matplotlib.axes._subplots.AxesSubplot at 0x11b70c9b0>
+
+ +
+ +
+ +
+ + + + +
+ +
+ +
+ +
+
+ +
+
+
+
+

For simplicity, we first convert the topography to binary, setting all points above sea level to be 1 and all points sea level or below to be 0. As a result, buffering will not depend on elevation.

+ +
+
+
+
+
+
In [4]:
+
+
+
Ztopo = np.where(topo.Z <= 0, 0, 1)
+
+ +
+
+
+ +
+
+
+
In [5]:
+
+
+
figure(figsize=(8,4))
+plottools.pcolorcells(topo.X, topo.Y, Ztopo, cmap='Greys', norm=norm)
+gca().set_aspect(1./cos(48*pi/180.))
+
+ +
+
+
+ +
+
+ + +
+ +
+ + + + +
+ +
+ +
+ +
+
+ +
+
+
+
+

Note: ideally Ztopo will only consist of island(s). If not, you can use a masked array to mask out coastal land as 0.

+ +
+
+
+
+
+
+

Creating a buffer region around the convex hull

+
+
+
+
+
+
+

We take the convex hull of all points above land. grid_points_in_poly tests whether points in the rectangular array lie within the convex hull.

+ +
+
+
+
+
+
In [6]:
+
+
+
coords = np.transpose(np.nonzero(Ztopo))
+
+hull = ConvexHull(coords)
+vertices = hull.points[hull.vertices]
+
+mask = convex_hull.grid_points_in_poly(Ztopo.shape, vertices)
+
+ +
+
+
+ +
+
+
+
In [7]:
+
+
+
figure(figsize=(8,4))
+plottools.pcolorcells(topo.X, topo.Y, mask, cmap='Greys', norm=norm)
+gca().set_aspect(1./cos(48*pi/180.))
+title('Convex hull')
+
+ +
+
+
+ +
+
+ + +
+ +
Out[7]:
+ + + + +
+
Text(0.5, 1.0, 'Convex hull')
+
+ +
+ +
+ +
+ + + + +
+ +
+ +
+ +
+
+ +
+
+
+
+

The convex hull can be specified as a Ruled Rectangle using region_tools.ruledrectangle_covering_selected_points.

+

Or, we can slightly extend the refinement region with a buffer function that expands the boundary of the polygon by a specified distance or fraction of degree.

+

For each side of the convex hull, we take the normalized perpendicular vector. Scaling this vector allows us to find lines parallel to the original edges, now expanded outward by some distance or fraction of degree in longitude/latitude. The intersection of these lines will give us the vertices of the buffered convex hull.

+

Note: vertices, in the order given, must define a convex polygon. Otherwise, pushing the edges outward may not result in a well-defined polygon.

+ +
+
+
+
+
+
In [8]:
+
+
+
def buffer(Ztopo, vertices, scale, method='meters', plot=False):
+    
+    '''
+    Given vertices of a convex polygon and a scale by which to buffer, 
+    returns expanded convex hull.
+    
+    vertices is an (N, 2) array of (row, column) indices.
+    
+    scale determines the width added to each side of the convex polygon.
+    
+    If method=='meters', then the polygon will be expanded by the specified 
+    meters on each side. If method=='degrees', then the polygon will be 
+    expanded by the specified fraction of degree on each side.
+        
+    Returns an array of the same shape as Ztopo, indicating points inside 
+    the buffered region as 1 and 0 otherwise.     
+    '''
+    
+    vertices[:,0] = [topo.Y[:,0][int(y)] for y in vertices[:,0]]
+    vertices[:,1] = [topo.X[0][int(x)] for x in vertices[:,1]]
+
+    for i in range(len(vertices) - 1):
+        if (vertices[i+1] == vertices[i]).all():
+            np.delete(vertices, i, axis=0)
+            
+    vertices = np.vstack((vertices, vertices[0]))
+
+    v = [np.subtract(vertices[i+1], vertices[i]) for i in range(len(vertices)-1)]
+    v = [[v[i][0], v[i][1] / abs(cos(radians(vertices[i][0])))] for i in range(len(vertices)-1)]
+    
+    if method == 'meters':
+        d = [haversine(vertices[i][1], vertices[i][0], vertices[i+1][1], vertices[i+1][0]) for i in range(len(vertices)-1)]
+        u = [v[i] / d[i] for i in range(len(v))] 
+    elif method == 'degrees':
+        u = [v[i] / np.linalg.norm(v[i]) for i in range(len(v))]
+
+    # rotate unit vectors (270 degrees if clockwise, 90 degrees if counterclockwise) 
+    # to find perpendicular vectors
+    n0 = [v[0][1], -v[0][0]]
+    if np.dot(n0, v[1]) < 0:
+        p = [[u[i][1], -1 * u[i][0]] for i in range(len(u))]
+    if np.dot(n0, v[1]) > 0:
+        p = [[-1 * u[i][1], u[i][0]] for i in range(len(u))]
+
+    p_scaled = np.array(p) * scale
+    
+    # points at the end of the scaled vector
+    points1 = [vertices[i] + p_scaled[i] for i in range(len(p_scaled))]
+    points2 = [vertices[i+1] + p_scaled[i] for i in range(len(p_scaled))]
+
+    for i in range(len(points1)):
+        if (points2[i][0] - points1[i][0]) == 0:
+            # bump point up by a small value if slope is 0
+            points1[i] += 10 ** -8
+
+    m = [(points2[i][1] - points1[i][1]) / (points2[i][0] - points1[i][0]) for i in range(len(points1))]
+
+    intersect_pts = []
+    for i in range(len(m)-1):
+        x = (points1[i+1][1] - points1[i][1] + m[i]*points1[i][0] - m[i+1]*points1[i+1][0]) / (m[i] - m[i+1])
+        y = m[i] * (x - points1[i][0]) + points1[i][1]
+        intersect_pts.append((x,y))
+
+    x = (points1[-1][1] - points1[0][1] + m[0]*points1[0][0] - m[-1]*points1[-1][0]) / (m[0] - m[-1])
+    y = m[0] * (x - points1[0][0]) + points1[0][1]
+    intersect_pts.append((x, y))
+
+    x_space = np.linspace(0, 120, 500)
+    line_eqns = [m[i] * (x_space - points1[i][0]) for i in range(len(m))]
+
+    intersect_pts.append(intersect_pts[0])
+
+    x, y = zip(*vertices)
+    buffered_x, buffered_y = zip(*intersect_pts)
+    
+    buffered_vertices = []
+    for lat, long in intersect_pts:
+        row = np.absolute(topo.Y[:,0] - lat).argmin()
+        col = np.absolute(topo.X[0] - long).argmin()
+        buffered_vertices.append((row, col))
+    
+    buffered_mask = convex_hull.grid_points_in_poly(Ztopo.shape, buffered_vertices)
+    
+    if plot == True:
+        
+        figure(figsize=(8,4))
+        fill(buffered_y, buffered_x, color='lightblue', label='buffer region')
+        fill(y, x, color='black', label='convex hull')
+        title('Buffer region around mask')
+        legend(loc='upper left')
+                
+    return buffered_mask
+
+ +
+
+
+ +
+
+
+
+

Plot the result of expanding our convex hull:

+ +
+
+
+
+
+
In [9]:
+
+
+
buffered_mask = buffer(topo.Z, vertices, scale=5000, method='meters', plot=True)
+
+ +
+
+
+ +
+
+ + +
+ +
+ + + + +
+ +
+ +
+ +
+
+ +
+
+
+
+

Creating a Ruled Rectangle

+
+
+
+
+
+
+

Now we can define a Ruled Rectangle around our buffered convex hull. See the flagregions documentation for writing out Ruled Rectangles as a data file to feed into GeoClaw.

+ +
+
+
+
+
+
In [10]:
+
+
+
# extend plotting region
+extend = [-.2, .2, -.15, .15]
+extended_region = [x + y for x, y in zip(extent, extend)]
+
+figure(figsize=(12,5))
+rr = region_tools.ruledrectangle_covering_selected_points(topo.X, topo.Y, buffered_mask, 
+                                                          ixy='y', method=0,
+                                                          padding=0, verbose=True)
+xv,yv = rr.vertices()
+plottools.pcolorcells(topo.X, topo.Y, mask, cmap='Greys', norm=norm)
+axis(extended_region)
+gca().set_aspect(1./cos(48*pi/180.))
+plot(xv, yv, 'r')
+title("With ixy = 'y'")
+
+ +
+
+
+ +
+
+ + +
+ +
+ + +
+
Extending rectangles to cover grid cells
+RuledRectangle covers 25018169 grid points
+
+
+
+ +
+ +
Out[10]:
+ + + + +
+
Text(0.5, 1.0, "With ixy = 'y'")
+
+ +
+ +
+ +
+ + + + +
+ +
+ +
+ +
+
+ +
+
+
+
+

The original convex hull is shown in black, and the Ruled Rectangle covers the buffered region. A few edges of the Ruled Rectangle are not extended as in the buffer because of the topography's extent. We can fix this by manually extending the Ruled Rectangle slightly to the left, right, and bottom:

+ +
+
+
+
+
+
In [11]:
+
+
+
# extend left side
+rr.lower = [lower - 0.04 if lower == min(rr.lower) else lower for lower in rr.lower]
+
+# extend bottom
+rr.s = [s - 0.01 if s == min(rr.s) else s for s in rr.s]
+
+# extend right side
+rr.upper = [upper + 0.04 if upper == max(rr.upper) else upper for upper in rr.upper]
+
+figure(figsize=(12,5))
+xv,yv = rr.vertices()
+plottools.pcolorcells(topo.X, topo.Y, mask, cmap='Greys', norm=norm)
+axis(extended_region)
+gca().set_aspect(1./cos(48*pi/180.))
+plot(xv, yv, 'r')
+title("With ixy = 'y'")
+
+ +
+
+
+ +
+
+ + +
+ +
Out[11]:
+ + + + +
+
Text(0.5, 1.0, "With ixy = 'y'")
+
+ +
+ +
+ +
+ + + + +
+ +
+ +
+ +
+
+ +
+
+
+
+

Using better convex hull approximations

+
+
+
+
+
+
+

In some cases, the convex hull around an island will span a larger region than we would like. In this example, we will repeat the buffering procedure above but use multiple convex hulls to better approximate the shape of the island.

+ +
+
+
+
+
+
+

We crop a 1/3-arc-second DEM of Cape Hatteras, available from the NCEI THREDDS server.

+ +
+
+
+
+
+
In [12]:
+
+
+
path = 'https://www.ngdc.noaa.gov/thredds/dodsC/regional/cape_hatteras_13_mhw_2006.nc'
+
+extent = [-75.85, -75.3, 35.0, 35.55]
+topo = topotools.read_netcdf(path, extent=extent)
+
+topo.plot()
+
+ +
+
+
+ +
+
+ + +
+ +
Out[12]:
+ + + + +
+
<matplotlib.axes._subplots.AxesSubplot at 0x123741978>
+
+ +
+ +
+ +
+ + + + +
+ +
+ +
+ +
+
+ +
+
+
+
+

Like in the previous example, we convert the topography into binary and find the convex hull:

+ +
+
+
+
+
+
In [13]:
+
+
+
Ztopo = np.where(topo.Z <= 0, 0, 1)
+
+ +
+
+
+ +
+
+
+
In [14]:
+
+
+
figure(figsize=(12,6))
+plottools.pcolorcells(topo.X, topo.Y, Ztopo, cmap='Greys', norm=norm)
+gca().set_aspect(1./cos(48*pi/180.))
+
+ +
+
+
+ +
+
+ + +
+ +
+ + + + +
+ +
+ +
+ +
+
+ +
+
+
+
In [15]:
+
+
+
coords = np.transpose(np.nonzero(Ztopo))
+
+hull = ConvexHull(coords)
+vertices = hull.points[hull.vertices]
+
+x,y = zip(*vertices)
+
+mask = convex_hull.grid_points_in_poly(Ztopo.shape, vertices)
+
+ +
+
+
+ +
+
+
+
In [16]:
+
+
+
figure(figsize=(12,6))
+plottools.pcolorcells(topo.X, topo.Y, mask, cmap='Greys', norm=norm)
+gca().set_aspect(1./cos(48*pi/180.))
+title('Convex hull (non-optimal)')
+
+ +
+
+
+ +
+
+ + +
+ +
Out[16]:
+ + + + +
+
Text(0.5, 1.0, 'Convex hull (non-optimal)')
+
+ +
+ +
+ +
+ + + + +
+ +
+ +
+ +
+
+ +
+
+
+
+

Creating a partitioned convex hull

+
+
+
+
+
+
+

The convex hull above is a fairly inefficient approximation of the island region. Our goal is to narrow the area of refinement to better fit the shape of the island.

+

To do this, we will recursively split the region in half and create seperate convex hulls for each half. Once the convexity measure of all partitions rises above a certain threshold, we return the union of convex hulls.

+ +
+
+
+
+
+
+

First, we define convexity as

+$$\text{convexity }=\frac{\text{island area}}{\text{convex hull area}}$$ +
+
+
+
+
+
In [17]:
+
+
+
def measure_convexity(Ztopo, mask):   
+    '''Return convexity measure as defined above'''
+    
+    land_area = np.count_nonzero(Ztopo)
+    chull_area = np.count_nonzero(mask)
+    
+    # if partition does not contain land, compute convexity measure to be 1 
+    if land_area == 0:
+        return 1
+    else:
+        return land_area / chull_area
+
+ +
+
+
+ +
+
+
+
+

If there are multiple islands in our region of interest, we will treat them as separate objects. In addition to taking the convex hull of each partition, then, we will also take the convex hull of each connected component. To do this, we label and create a hull object for each component:

+ +
+
+
+
+
+
In [18]:
+
+
+
def chull_objects(Ztopo):
+    '''Return convex hull of entire topo and list of convex hull objects for each connected component'''
+    
+    chull = np.zeros(Ztopo.shape, dtype=bool)
+    obj = []
+    
+    labeled_topo = label(Ztopo, connectivity=2, background=0)
+        
+    for i in range(1, labeled_topo.max() + 1):
+        coords = np.transpose(np.nonzero(labeled_topo == i))
+        
+        # ignore small dots
+        if len(coords) < 3:
+            continue
+        
+        # if points are on same line, group with next labeled object
+        if (all([c[0] == coords[0][0] for c in coords]) or all([c[1] == coords[0][1] for c in coords])):
+            labeled_topo = np.where(labeled_topo == i, i+1, labeled_topo)
+        else:        
+            hull = ConvexHull(coords)
+            obj.append(hull)
+        
+            vertices = hull.points[hull.vertices]        
+            mask = convex_hull.grid_points_in_poly(labeled_topo.shape, vertices)
+            chull = np.logical_or(chull, mask)
+    
+    return chull, obj
+
+ +
+
+
+ +
+
+
+
+

Now, we create the recursive function. In each iteration, the function splits the binary topograpy (or a section of the topography) in half and computes the convexity measure. If the convexity falls below our threshold in any half, we constrain Ztopo to that half by redefining the starting and ending indices. The function repeatedly partitions the topography, until the convexity of all partitions reach the treshold, while adding the hull objects to a list.

+

Note that a lower threshold will give a coarser approximation, while a higher threshold will give a finer approximation.

+

Also note that the function splits the topography vertically. If you would like to split the topography horizontally, you can simply change the slicing.

+ +
+
+
+
+
+
In [19]:
+
+
+
def partitioned_chull(Ztopo, threshold, start=0, convexity=0, end=0):
+    '''Finds each convex object in a partitioned convex hull'''
+    
+    global objects
+    
+    if convexity >= threshold:
+        return convex_hull_object(Ztopo)
+    
+    mid_col = int(np.floor((start + end) / 2))
+    
+    Ztopo[:, 0:start] = 0
+    Ztopo[:, end:len(Ztopo)] = 0
+    
+    topo_left = Ztopo.copy()
+    topo_left[:, mid_col:] = 0
+    topo_right = Ztopo.copy()
+    topo_right[:, 0:mid_col] = 0
+    
+    chull_left, obj_left = chull_objects(topo_left)
+    chull_right, obj_right = chull_objects(topo_right)
+    
+    conv_left = measure_convexity(topo_left, chull_left)
+    conv_right = measure_convexity(topo_right, chull_right)
+    
+    print('Convexity of left half: ', conv_left)
+    print('Convexity of right half: ', conv_right) 
+    
+    if conv_left >= threshold:
+        objects.extend(obj_left)
+        print('Left half satisfied')
+    if conv_right >= threshold:
+        objects.extend(obj_right)
+        print('Right half satisfied')    
+
+    if conv_left < threshold or conv_right < threshold: 
+        
+        chull_left = partitioned_chull(Ztopo, threshold=threshold, convexity=conv_left, start=start, end=mid_col)
+        chull_right = partitioned_chull(Ztopo, threshold=threshold, convexity=conv_right, start=mid_col+1, end=end)
+
+        return np.logical_or(chull_left, chull_right)
+    
+    else:
+
+        return np.logical_or(chull_left, chull_right)
+
+ +
+
+
+ +
+
+
+
+

Let's try it on our Cape Hatteras topography:

+ +
+
+
+
+
+
In [20]:
+
+
+
objects = []
+chull = partitioned_chull(Ztopo, threshold=0.4, end=len(Ztopo))
+
+ +
+
+
+ +
+
+ + +
+ +
+ + +
+
Convexity of left half:  0.46072521666290794
+Convexity of right half:  0.3253215514315942
+Left half satisfied
+Convexity of left half:  0.32538890812711435
+Convexity of right half:  1
+Right half satisfied
+Convexity of left half:  0.489490434318509
+Convexity of right half:  0.5094851236964442
+Left half satisfied
+Right half satisfied
+
+
+
+ +
+
+ +
+
+
+
+

We can combine all of our convex hull objects into a single mask by taking the logical OR:

+ +
+
+
+
+
+
In [21]:
+
+
+
mask = np.zeros(Ztopo.shape, dtype=bool)
+
+for obj in objects:
+    
+    vertices = obj.points[obj.vertices]    
+    
+    obj_mask = convex_hull.grid_points_in_poly(Ztopo.shape, vertices)
+    mask = np.logical_or(mask, obj_mask)
+
+ +
+
+
+ +
+
+
+
In [22]:
+
+
+
figure(figsize=(12,6))
+plottools.pcolorcells(topo.X, topo.Y, mask, cmap='Greys', norm=norm)
+gca().set_aspect(1./cos(48*pi/180.))
+title('Partitioned convex hull')
+
+ +
+
+
+ +
+
+ + +
+ +
Out[22]:
+ + + + +
+
Text(0.5, 1.0, 'Partitioned convex hull')
+
+ +
+ +
+ +
+ + + + +
+ +
+ +
+ +
+
+ +
+
+
+
+

Creating ruled rectangles around the buffered region

+
+
+
+
+
+
+

Now we can apply the buffering function on each convex hull object. We can then either define Ruled Rectangles for each convex hull object or define a Ruled Rectangle for the entire mask. The first option may be useful if island points are fairly scattered in the topography. The second option may be simpler if the points are more concentrated.

+ +
+
+
+
+
+
+

For this example, we will skip some small objects (fewer than 1000 points). You might want to experiment with which objects to include. You will also likely want to toggle the scale to change the buffering amount.

+ +
+
+
+
+
+
In [23]:
+
+
+
rrs = []
+buffered_mask = np.zeros(Ztopo.shape, dtype=bool)
+
+for obj in objects:
+    
+    vertices = obj.points[obj.vertices]    
+    
+    obj_mask = convex_hull.grid_points_in_poly(Ztopo.shape, vertices)
+    area = np.count_nonzero(obj_mask)   
+    
+    if area < 1000:
+        continue
+    
+    buffered_obj = buffer(topo.Z, vertices, scale=0.02, method='degrees', plot=False)
+    
+    # for option 1
+    rr = region_tools.ruledrectangle_covering_selected_points(topo.X, topo.Y, buffered_obj, 
+                                                         ixy='x', method=0,
+                                                         padding=0, verbose=True)
+    rrs.append(rr)  
+    
+    # for option 2
+    buffered_mask = np.logical_or(buffered_mask, buffered_obj)
+
+ +
+
+
+ +
+
+ + +
+ +
+ + +
+
Extending rectangles to cover grid cells
+RuledRectangle covers 565123 grid points
+Extending rectangles to cover grid cells
+RuledRectangle covers 1635229 grid points
+Extending rectangles to cover grid cells
+RuledRectangle covers 219413 grid points
+Extending rectangles to cover grid cells
+RuledRectangle covers 178325 grid points
+Extending rectangles to cover grid cells
+RuledRectangle covers 1572834 grid points
+Extending rectangles to cover grid cells
+RuledRectangle covers 1574951 grid points
+
+
+
+ +
+
+ +
+
+
+
+

Plot Ruled Rectangles covering each object:

+ +
+
+
+
+
+
In [24]:
+
+
+
figure(figsize=(12,5))
+axis(extent) 
+title("With ixy = 'x'")
+plottools.pcolorcells(topo.X, topo.Y, mask, cmap='Greys', norm=norm) 
+
+for rr in rrs:
+    
+    xv,yv = rr.vertices()
+    
+    plt.gca().set_aspect(1./cos(48*pi/180.))
+    plt.plot(xv, yv, 'r')    
+
+ +
+
+
+ +
+
+ + +
+ +
+ + + + +
+ +
+ +
+ +
+
+ +
+
+
+
+

Create and plot Ruled Rectangle covering the entire buffered region:

+ +
+
+
+
+
+
In [25]:
+
+
+
figure(figsize=(12,5))
+subplot(121)
+rr2 = region_tools.ruledrectangle_covering_selected_points(topo.X, topo.Y, buffered_mask, 
+                                                          ixy='x', method=0,
+                                                          padding=0, verbose=False)
+xv,yv = rr2.vertices()
+plottools.pcolorcells(topo.X, topo.Y, mask, cmap='Greys', norm=norm)
+axis(extent) 
+gca().set_aspect(1./cos(48*pi/180.))
+plot(xv, yv, 'r')
+title("With ixy = 'x'")
+
+subplot(122)
+rr2 = region_tools.ruledrectangle_covering_selected_points(topo.X, topo.Y, buffered_mask, 
+                                                          ixy='y', method=0,
+                                                          padding=0, verbose=False)
+xv,yv = rr2.vertices()
+plottools.pcolorcells(topo.X, topo.Y, mask, cmap='Greys', norm=norm)
+axis(extent) 
+gca().set_aspect(1./cos(48*pi/180.))
+plot(xv, yv, 'r')
+title("With ixy = 'y'")
+
+ +
+
+
+ +
+
+ + +
+ +
Out[25]:
+ + + + +
+
Text(0.5, 1.0, "With ixy = 'y'")
+
+ +
+ +
+ +
+ + + + +
+ +
+ +
+ +
+
+ +
+
+
+ + + + + + diff --git a/gallery/_static/documentation_options.js b/gallery/_static/documentation_options.js index d9d9985639..7ea22573dc 100644 --- a/gallery/_static/documentation_options.js +++ b/gallery/_static/documentation_options.js @@ -1,6 +1,6 @@ var DOCUMENTATION_OPTIONS = { URL_ROOT: document.getElementById("documentation_options").getAttribute('data-url_root'), - VERSION: '5.7.x', + VERSION: '5.8.x', LANGUAGE: 'None', COLLAPSE_INDEX: false, BUILDER: 'html', diff --git a/gallery/gallery/amrclaw/examples/acoustics_1d_heterogeneous/README.html b/gallery/gallery/amrclaw/examples/acoustics_1d_heterogeneous/README.html index c61da0b610..8a7027f457 100644 --- a/gallery/gallery/amrclaw/examples/acoustics_1d_heterogeneous/README.html +++ b/gallery/gallery/amrclaw/examples/acoustics_1d_heterogeneous/README.html @@ -4,7 +4,7 @@ - Acoustics 1D – Heterogeneous example — Clawpack 5.7.x documentation + Acoustics 1D – Heterogeneous example — Clawpack 5.8.x documentation @@ -63,7 +63,7 @@

Navigation

  • modules |
    • -
  • Clawpack 5.7.x documentation »
    • +
  • Clawpack 5.8.x documentation »
    • @@ -92,7 +92,7 @@

    Navigation

    -

    Version 5.7.x

    +

    Version 5.8.x

    Related Topics