matplotlib
diff --git a/‎lib/mpl_toolkits/mplot3d/art3d.py
Copy file name to clipboardExpand all lines: lib/mpl_toolkits/mplot3d/art3d.py
+41Lines changed: 41 additions & 0 deletions b/‎lib/mpl_toolkits/mplot3d/art3d.py
Copy file name to clipboardExpand all lines: lib/mpl_toolkits/mplot3d/art3d.py
+41Lines changed: 41 additions & 0 deletions
diff --git a/‎lib/mpl_toolkits/mplot3d/axes3d.py
Copy file name to clipboardExpand all lines: lib/mpl_toolkits/mplot3d/axes3d.py
+17-16Lines changed: 17 additions & 16 deletions b/‎lib/mpl_toolkits/mplot3d/axes3d.py
Copy file name to clipboardExpand all lines: lib/mpl_toolkits/mplot3d/axes3d.py
+17-16Lines changed: 17 additions & 16 deletions
diff --git a/‎lib/mpl_toolkits/mplot3d/tests/test_art3d.py
Copy file name to clipboardExpand all lines: lib/mpl_toolkits/mplot3d/tests/test_art3d.py
+32-1Lines changed: 32 additions & 1 deletion b/‎lib/mpl_toolkits/mplot3d/tests/test_art3d.py
Copy file name to clipboardExpand all lines: lib/mpl_toolkits/mplot3d/tests/test_art3d.py
+32-1Lines changed: 32 additions & 1 deletion
@@ -1185,6 +1185,47 @@ def _zalpha(colors, zs):
     return np.column_stack([rgba[:, :3], rgba[:, 3] * sats])
 
 
+def _all_points_on_plane(xs, ys, zs, atol=1e-8):
+    """
+    Check if all points are on the same plane. Note that NaN values are
+    ignored.
+
+    Parameters
+    ----------
+    xs, ys, zs : array-like
+        The x, y, and z coordinates of the points.
+    atol : float, default: 1e-8
+        The tolerance for the equality check.
+    """
+    xs, ys, zs = np.asarray(xs), np.asarray(ys), np.asarray(zs)
+    points = np.column_stack([xs, ys, zs])
+    points = points[~np.isnan(points).any(axis=1)]
+    # Check for the case where we have less than 3 unique points
+    points = np.unique(points, axis=0)
+    if len(points) <= 3:
+        return True
+    # Calculate the vectors from the first point to all other points
+    vs = (points - points[0])[1:]
+    vs = vs / np.linalg.norm(vs, axis=1)[:, np.newaxis]
+    # Filter out parallel vectors
+    vs = np.unique(vs, axis=0)
+    if len(vs) <= 2:
+        return True
+    # Filter out parallel and antiparallel vectors to the first vector
+    cross_norms = np.linalg.norm(np.cross(vs[0], vs[1:]), axis=1)
+    zero_cross_norms = np.where(np.isclose(cross_norms, 0, atol=atol))[0] + 1
+    vs = np.delete(vs, zero_cross_norms, axis=0)
+    if len(vs) <= 2:
+        return True
+    # Calculate the normal vector from the first three points
+    n = np.cross(vs[0], vs[1])
+    n = n / np.linalg.norm(n)
+    # If the dot product of the normal vector and all other vectors is zero,
+    # all points are on the same plane
+    dots = np.dot(n, vs.transpose())
+    return np.allclose(dots, 0, atol=atol)
+
+
 def _generate_normals(polygons):
     """
     Compute the normals of a list of polygons, one normal per polygon.
 
@@ -1994,9 +1994,8 @@ def fill_between(self, x1, y1, z1, x2, y2, z2, *,
             - 'polygon': The two lines are connected to form a single polygon.
               This is faster and can render more cleanly for simple shapes
               (e.g. for filling between two lines that lie within a plane).
-            - 'auto': If the lines are in a plane parallel to a coordinate axis
-              (one of *x*, *y*, *z* are constant and equal for both lines),
-              'polygon' is used. Otherwise, 'quad' is used.
+            - 'auto': If the points all lie on the same 3D plane, 'polygon' is
+              used. Otherwise, 'quad' is used.
 
         facecolors : list of :mpltype:`color`, default: None
             Colors of each individual patch, or a single color to be used for
@@ -2019,19 +2018,6 @@ def fill_between(self, x1, y1, z1, x2, y2, z2, *,
 
         had_data = self.has_data()
         x1, y1, z1, x2, y2, z2 = cbook._broadcast_with_masks(x1, y1, z1, x2, y2, z2)
-        if mode == 'auto':
-            if ((np.all(x1 == x1[0]) and np.all(x2 == x1[0]))
-                 or (np.all(y1 == y1[0]) and np.all(y2 == y1[0]))
-                 or (np.all(z1 == z1[0]) and np.all(z2 == z1[0]))):
-                mode = 'polygon'
-            else:
-                mode = 'quad'
-
-        if shade is None:
-            if mode == 'quad':
-                shade = True
-            else:
-                shade = False
 
         if facecolors is None:
             facecolors = [self._get_patches_for_fill.get_next_color()]
@@ -2046,6 +2032,21 @@ def fill_between(self, x1, y1, z1, x2, y2, z2, *,
                                  f"size ({x1.size})")
         where = where & ~np.isnan(x1)  # NaNs were broadcast in _broadcast_with_masks
 
+        if mode == 'auto':
+            if art3d._all_points_on_plane(np.concatenate((x1[where], x2[where])),
+                                          np.concatenate((y1[where], y2[where])),
+                                          np.concatenate((z1[where], z2[where])),
+                                          atol=1e-12):
+                mode = 'polygon'
+            else:
+                mode = 'quad'
+
+        if shade is None:
+            if mode == 'quad':
+                shade = True
+            else:
+                shade = False
+
         polys = []
         for idx0, idx1 in cbook.contiguous_regions(where):
             x1i = x1[idx0:idx1]
 
@@ -3,7 +3,7 @@
 import matplotlib.pyplot as plt
 
 from matplotlib.backend_bases import MouseEvent
-from mpl_toolkits.mplot3d.art3d import Line3DCollection
+from mpl_toolkits.mplot3d.art3d import Line3DCollection, _all_points_on_plane
 
 
 def test_scatter_3d_projection_conservation():
@@ -54,3 +54,34 @@ def test_zordered_error():
     ax.add_collection(Line3DCollection(lc))
     ax.scatter(*pc, visible=False)
     plt.draw()
+
+
+def test_all_points_on_plane():
+    # Non-coplanar points
+    points = np.array([[0, 0, 0], [1, 0, 0], [0, 1, 0], [0, 0, 1]])
+    assert not _all_points_on_plane(*points.T)
+
+    # Duplicate points
+    points = np.array([[0, 0, 0], [1, 0, 0], [0, 1, 0], [0, 0, 0]])
+    assert _all_points_on_plane(*points.T)
+
+    # NaN values
+    points = np.array([[0, 0, 0], [1, 0, 0], [0, 1, 0], [0, 0, np.nan]])
+    assert _all_points_on_plane(*points.T)
+
+    # Less than 3 unique points
+    points = np.array([[0, 0, 0], [0, 0, 0], [0, 0, 0]])
+    assert _all_points_on_plane(*points.T)
+
+    # All points lie on a line
+    points = np.array([[0, 0, 0], [0, 1, 0], [0, 2, 0], [0, 3, 0]])
+    assert _all_points_on_plane(*points.T)
+
+    # All points lie on two lines, with antiparallel vectors
+    points = np.array([[-2, 2, 0], [-1, 1, 0], [1, -1, 0],
+                       [0, 0, 0], [2, 0, 0], [1, 0, 0]])
+    assert _all_points_on_plane(*points.T)
+
+    # All points lie on a plane
+    points = np.array([[0, 0, 0], [0, 1, 0], [1, 0, 0], [1, 1, 0], [1, 2, 0]])
+    assert _all_points_on_plane(*points.T)