diff --git a/source/blender/blenlib/BLI_math_geom.h b/source/blender/blenlib/BLI_math_geom.h
index d056c42e019..e24f3c3bde8 100644
--- a/source/blender/blenlib/BLI_math_geom.h
+++ b/source/blender/blenlib/BLI_math_geom.h
@@ -1164,8 +1164,8 @@ void box_minmax_bounds_m4(float min[3], float max[3], float boundbox[2][3], floa
 /** \name Mapping
  * \{ */
 
-void map_to_tube(float *r_u, float *r_v, float x, float y, float z);
-void map_to_sphere(float *r_u, float *r_v, float x, float y, float z);
+bool map_to_tube(float *r_u, float *r_v, float x, float y, float z);
+bool map_to_sphere(float *r_u, float *r_v, float x, float y, float z);
 void map_to_plane_v2_v3v3(float r_co[2], const float co[3], const float no[3]);
 void map_to_plane_axis_angle_v2_v3v3fl(float r_co[2],
                                        const float co[3],
diff --git a/source/blender/blenlib/intern/math_geom.c b/source/blender/blenlib/intern/math_geom.c
index 08152976f7d..415e21cd272 100644
--- a/source/blender/blenlib/intern/math_geom.c
+++ b/source/blender/blenlib/intern/math_geom.c
@@ -4914,39 +4914,59 @@ void box_minmax_bounds_m4(float min[3], float max[3], float boundbox[2][3], floa
 
 /********************************** Mapping **********************************/
 
-void map_to_tube(float *r_u, float *r_v, const float x, const float y, const float z)
+static float snap_coordinate(float u)
 {
-  float len;
-
-  *r_v = (z + 1.0f) / 2.0f;
-
-  len = sqrtf(x * x + y * y);
-  if (len > 0.0f) {
-    *r_u = (1.0f - (atan2f(x / len, y / len) / (float)M_PI)) / 2.0f;
+  /* Adjust a coordinate value `u` to obtain a value inside the (closed) unit interval.
+   *   i.e. 0.0 <= snap_coordinate(u) <= 1.0.
+   * Due to round-off errors, it is possible that the value of `u` may be close to the boundary of
+   * the unit interval, but not exactly on it. In order to handle these cases, `snap_coordinate`
+   * checks whether `u` is within `epsilon` of the boundary, and if so, it snaps the return value
+   * to the boundary. */
+  if (u < 0.0f) {
+    u += 1.0f; /* Get back into the unit interval. */
   }
-  else {
-    *r_v = *r_u = 0.0f; /* to avoid un-initialized variables */
+  BLI_assert(0.0f <= u);
+  BLI_assert(u <= 1.0f);
+  const float epsilon = 0.25f / 65536.0f; /* i.e. Quarter of a texel on a 65536 x 65536 texture. */
+  if (u < epsilon) {
+    return 0.0f; /* `u` is close to 0, just return 0. */
   }
+  if (1.0f - epsilon < u) {
+    return 1.0f; /* `u` is close to 1, just return 1. */
+  }
+  return u;
 }
 
-void map_to_sphere(float *r_u, float *r_v, const float x, const float y, const float z)
+bool map_to_tube(float *r_u, float *r_v, const float x, const float y, const float z)
 {
-  float len;
-
-  len = sqrtf(x * x + y * y + z * z);
-  if (len > 0.0f) {
-    if (UNLIKELY(x == 0.0f && y == 0.0f)) {
-      *r_u = 0.0f; /* Otherwise domain error. */
-    }
-    else {
-      *r_u = (1.0f - atan2f(x, y) / (float)M_PI) / 2.0f;
-    }
-
-    *r_v = 1.0f - saacos(z / len) / (float)M_PI;
+  bool regular = true;
+  if (x * x + y * y < 1e-6f * 1e-6f) {
+    regular = false; /* We're too close to the cylinder's axis. */
+    *r_u = 0.5f;
   }
   else {
-    *r_v = *r_u = 0.0f; /* to avoid un-initialized variables */
+    /* The "Regular" case, just compute the coordinate. */
+    *r_u = snap_coordinate(atan2f(x, -y) / (float)(2.0f * M_PI));
   }
+  *r_v = (z + 1.0f) / 2.0f;
+  return regular;
+}
+
+bool map_to_sphere(float *r_u, float *r_v, const float x, const float y, const float z)
+{
+  bool regular = true;
+  const float epsilon = 0.25f / 65536.0f; /* i.e. Quarter of a texel on a 65536 x 65536 texture. */
+  const float len_xy = sqrtf(x * x + y * y);
+  if (len_xy <= fabsf(z) * epsilon) {
+    regular = false; /* We're on the line that runs through the north and south poles. */
+    *r_u = 0.5f;
+  }
+  else {
+    /* The "Regular" case, just compute the coordinate. */
+    *r_u = snap_coordinate(atan2f(x, -y) / (float)(2.0f * M_PI));
+  }
+  *r_v = snap_coordinate(atan2f(len_xy, -z) / (float)M_PI);
+  return regular;
 }
 
 void map_to_plane_v2_v3v3(float r_co[2], const float co[3], const float no[3])
diff --git a/source/blender/editors/uvedit/uvedit_unwrap_ops.c b/source/blender/editors/uvedit/uvedit_unwrap_ops.c
index 18e920c55f4..0f8830a4521 100644
--- a/source/blender/editors/uvedit/uvedit_unwrap_ops.c
+++ b/source/blender/editors/uvedit/uvedit_unwrap_ops.c
@@ -1318,6 +1318,9 @@ void ED_uvedit_live_unwrap_end(short cancel)
 #define POLAR_ZX 0
 #define POLAR_ZY 1
 
+#define PINCH 0
+#define FAN 1
+
 static void uv_map_transform_calc_bounds(BMEditMesh *em, float r_min[3], float r_max[3])
 {
   BMFace *efa;
@@ -1457,50 +1460,37 @@ static void uv_map_rotation_matrix_ex(float result[4][4],
   mul_m4_series(result, rotup, rotside, viewmatrix, rotobj);
 }
 
-static void uv_map_rotation_matrix(float result[4][4],
-                                   RegionView3D *rv3d,
-                                   Object *ob,
-                                   float upangledeg,
-                                   float sideangledeg,
-                                   float radius)
+static void uv_map_transform(bContext *C, wmOperator *op, float rotmat[3][3])
 {
-  const float offset[4] = {0};
-  uv_map_rotation_matrix_ex(result, rv3d, ob, upangledeg, sideangledeg, radius, offset);
-}
-
-static void uv_map_transform(bContext *C, wmOperator *op, float rotmat[4][4])
-{
-  /* context checks are messy here, making it work in both 3d view and uv editor */
   Object *obedit = CTX_data_edit_object(C);
   RegionView3D *rv3d = CTX_wm_region_view3d(C);
-  /* common operator properties */
-  int align = RNA_enum_get(op->ptr, "align");
-  int direction = RNA_enum_get(op->ptr, "direction");
-  float radius = RNA_struct_find_property(op->ptr, "radius") ? RNA_float_get(op->ptr, "radius") :
-                                                               1.0f;
-  float upangledeg, sideangledeg;
 
-  if (direction == VIEW_ON_EQUATOR) {
-    upangledeg = 90.0f;
-    sideangledeg = 0.0f;
-  }
-  else {
-    upangledeg = 0.0f;
-    if (align == POLAR_ZY) {
-      sideangledeg = 0.0f;
-    }
-    else {
-      sideangledeg = 90.0f;
-    }
-  }
+  const int align = RNA_enum_get(op->ptr, "align");
+  const int direction = RNA_enum_get(op->ptr, "direction");
+  const float radius = RNA_struct_find_property(op->ptr, "radius") ?
+                           RNA_float_get(op->ptr, "radius") :
+                           1.0f;
 
-  /* be compatible to the "old" sphere/cylinder mode */
+  /* Be compatible to the "old" sphere/cylinder mode. */
   if (direction == ALIGN_TO_OBJECT) {
-    unit_m4(rotmat);
-  }
-  else {
-    uv_map_rotation_matrix(rotmat, rv3d, obedit, upangledeg, sideangledeg, radius);
+    unit_m3(rotmat);
+
+    if (align == POLAR_ZY) {
+      rotmat[0][0] = 0.0f;
+      rotmat[0][1] = 1.0f;
+      rotmat[1][0] = -1.0f;
+      rotmat[1][1] = 0.0f;
+    }
+    return;
   }
+
+  const float up_angle_deg = (direction == VIEW_ON_EQUATOR) ? 90.0f : 0.0f;
+  const float side_angle_deg = (align == POLAR_ZY) == (direction == VIEW_ON_EQUATOR) ? 90.0f :
+                                                                                       0.0f;
+  const float offset[4] = {0};
+  float rotmat4[4][4];
+  uv_map_rotation_matrix_ex(rotmat4, rv3d, obedit, up_angle_deg, side_angle_deg, radius, offset);
+  copy_m3_m4(rotmat, rotmat4);
 }
 
 static void uv_transform_properties(wmOperatorType *ot, int radius)
@@ -1521,6 +1511,12 @@ static void uv_transform_properties(wmOperatorType *ot, int radius)
       {0, NULL, 0, NULL, NULL},
   };
 
+  static const EnumPropertyItem pole_items[] = {
+      {PINCH, "PINCH", 0, "Pinch", "UVs are pinched at the poles"},
+      {FAN, "FAN", 0, "Fan", "UVs are fanned at the poles"},
+      {0, NULL, 0, NULL, NULL},
+  };
+
   RNA_def_enum(ot->srna,
                "direction",
                direction_items,
@@ -1530,9 +1526,10 @@ static void uv_transform_properties(wmOperatorType *ot, int radius)
   RNA_def_enum(ot->srna,
                "align",
                align_items,
-               VIEW_ON_EQUATOR,
+               POLAR_ZX,
                "Align",
                "How to determine rotation around the pole");
+  RNA_def_enum(ot->srna, "pole", pole_items, PINCH, "Pole", "How to handle faces at the poles");
   if (radius) {
     RNA_def_float(ot->srna,
                   "radius",
@@ -2730,55 +2727,120 @@ void UV_OT_reset(wmOperatorType *ot)
 /** \name Sphere UV Project Operator
  * \{ */
 
-static void uv_sphere_project(float target[2],
-                              const float source[3],
-                              const float center[3],
-                              const float rotmat[4][4])
+static void uv_map_mirror(BMFace *efa,
+                          const bool *regular,
+                          const bool fan,
+                          const int cd_loop_uv_offset)
 {
-  float pv[3];
+  /* A heuristic to improve alignment of faces near the seam.
+   * In simple terms, we're looking for faces which span more
+   * than 0.5 units in the *u* coordinate.
+   * If we find such a face, we try and improve the unwrapping
+   * by adding (1.0, 0.0) onto some of the face's UVs.
 
-  sub_v3_v3v3(pv, source, center);
-  mul_m4_v3(rotmat, pv);
+   * Note that this is only a heuristic. The property we're
+   * attempting to maintain is that the winding of the face
+   * in UV space corresponds with the handedness of the face
+   * in 3D space w.r.t to the unwrapping. Even for triangles,
+   * that property is somewhat complicated to evaluate.
+   */
 
-  map_to_sphere(&target[0], &target[1], pv[0], pv[1], pv[2]);
-
-  /* split line is always zero */
-  if (target[0] >= 1.0f) {
-    target[0] -= 1.0f;
-  }
-}
-
-static void uv_map_mirror(BMEditMesh *em, BMFace *efa)
-{
+  float right_u = -1.0e30f;
   BMLoop *l;
   BMIter liter;
-  MLoopUV *luv;
   float **uvs = BLI_array_alloca(uvs, efa->len);
-  float dx;
-  int i, mi;
-
-  const int cd_loop_uv_offset = CustomData_get_offset(&em->bm->ldata, CD_MLOOPUV);
-
-  BM_ITER_ELEM_INDEX (l, &liter, efa, BM_LOOPS_OF_FACE, i) {
-    luv = BM_ELEM_CD_GET_VOID_P(l, cd_loop_uv_offset);
-    uvs[i] = luv->uv;
+  int j;
+  BM_ITER_ELEM_INDEX (l, &liter, efa, BM_LOOPS_OF_FACE, j) {
+    MLoopUV *luv = BM_ELEM_CD_GET_VOID_P(l, cd_loop_uv_offset);
+    uvs[j] = luv->uv;
+    if (luv->uv[0] >= 1.0f) {
+      luv->uv[0] -= 1.0f;
+    }
+    right_u = max_ff(right_u, luv->uv[0]);
   }
 
-  mi = 0;
-  for (i = 1; i < efa->len; i++) {
-    if (uvs[i][0] > uvs[mi][0]) {
-      mi = i;
+  float left_u = 1.0e30f;
+  BM_ITER_ELEM (l, &liter, efa, BM_LOOPS_OF_FACE) {
+    MLoopUV *luv = BM_ELEM_CD_GET_VOID_P(l, cd_loop_uv_offset);
+    if (right_u <= luv->uv[0] + 0.5f) {
+      left_u = min_ff(left_u, luv->uv[0]);
     }
   }
 
-  for (i = 0; i < efa->len; i++) {
-    if (i != mi) {
-      dx = uvs[mi][0] - uvs[i][0];
-      if (dx > 0.5f) {
-        uvs[i][0] += 1.0f;
+  BM_ITER_ELEM (l, &liter, efa, BM_LOOPS_OF_FACE) {
+    MLoopUV *luv = BM_ELEM_CD_GET_VOID_P(l, cd_loop_uv_offset);
+    if (luv->uv[0] + 0.5f < right_u) {
+      if (2 * luv->uv[0] + 1.0f < left_u + right_u) {
+        luv->uv[0] += 1.0f;
       }
     }
   }
+  if (!fan) {
+    return;
+  }
+
+  /* Another heuristic, this time, we attempt to "fan"
+   * the UVs of faces which pass through one of the poles
+   * of the unwrapping. */
+
+  /* Need to recompute min and max. */
+  float minmax_u[2] = {1.0e30f, -1.0e30f};
+  int pole_count = 0;
+  for (int i = 0; i < efa->len; i++) {
+    if (regular[i]) {
+      minmax_u[0] = min_ff(minmax_u[0], uvs[i][0]);
+      minmax_u[1] = max_ff(minmax_u[1], uvs[i][0]);
+    }
+    else {
+      pole_count++;
+    }
+  }
+  if (pole_count == 0 || pole_count == efa->len) {
+    return;
+  }
+  for (int i = 0; i < efa->len; i++) {
+    if (regular[i]) {
+      continue;
+    }
+    float u = 0.0f;
+    float sum = 0.0f;
+    const int i_plus = (i + 1) % efa->len;
+    const int i_minus = (i + efa->len - 1) % efa->len;
+    if (regular[i_plus]) {
+      u += uvs[i_plus][0];
+      sum += 1.0f;
+    }
+    if (regular[i_minus]) {
+      u += uvs[i_minus][0];
+      sum += 1.0f;
+    }
+    if (sum == 0) {
+      u += minmax_u[0] + minmax_u[1];
+      sum += 2.0f;
+    }
+    uvs[i][0] = u / sum;
+  }
+}
+
+static void uv_sphere_project(BMFace *efa,
+                              const float center[3],
+                              const float rotmat[3][3],
+                              const bool fan,
+                              const int cd_loop_uv_offset)
+{
+  bool *regular = BLI_array_alloca(regular, efa->len);
+  int i;
+  BMLoop *l;
+  BMIter iter;
+  BM_ITER_ELEM_INDEX (l, &iter, efa, BM_LOOPS_OF_FACE, i) {
+    MLoopUV *luv = BM_ELEM_CD_GET_VOID_P(l, cd_loop_uv_offset);
+    float pv[3];
+    sub_v3_v3v3(pv, l->v->co, center);
+    mul_m3_v3(rotmat, pv);
+    regular[i] = map_to_sphere(&luv->uv[0], &luv->uv[1], pv[0], pv[1], pv[2]);
+  }
+
+  uv_map_mirror(efa, regular, fan, cd_loop_uv_offset);
 }
 
 static int sphere_project_exec(bContext *C, wmOperator *op)
@@ -2800,9 +2862,7 @@ static int sphere_project_exec(bContext *C, wmOperator *op)
     Object *obedit = objects[ob_index];
     BMEditMesh *em = BKE_editmesh_from_object(obedit);
     BMFace *efa;
-    BMLoop *l;
-    BMIter iter, liter;
-    MLoopUV *luv;
+    BMIter iter;
 
     if (em->bm->totfacesel == 0) {
       continue;
@@ -2814,11 +2874,13 @@ static int sphere_project_exec(bContext *C, wmOperator *op)
     }
 
     const int cd_loop_uv_offset = CustomData_get_offset(&em->bm->ldata, CD_MLOOPUV);
-    float center[3], rotmat[4][4];
+    float center[3], rotmat[3][3];
 
     uv_map_transform(C, op, rotmat);
     uv_map_transform_center(scene, v3d, obedit, em, center, NULL);
 
+    const bool fan = RNA_enum_get(op->ptr, "pole");
+
     BM_ITER_MESH (efa, &iter, em->bm, BM_FACES_OF_MESH) {
       if (!BM_elem_flag_test(efa, BM_ELEM_SELECT)) {
         continue;
@@ -2831,13 +2893,7 @@ static int sphere_project_exec(bContext *C, wmOperator *op)
         }
       }
 
-      BM_ITER_ELEM (l, &liter, efa, BM_LOOPS_OF_FACE) {
-        luv = BM_ELEM_CD_GET_VOID_P(l, cd_loop_uv_offset);
-
-        uv_sphere_project(luv->uv, l->v->co, center, rotmat);
-      }
-
-      uv_map_mirror(em, efa);
+      uv_sphere_project(efa, center, rotmat, fan, cd_loop_uv_offset);
     }
 
     const bool per_face_aspect = true;
@@ -2875,22 +2931,25 @@ void UV_OT_sphere_project(wmOperatorType *ot)
 /** \name Cylinder UV Project Operator
  * \{ */
 
-static void uv_cylinder_project(float target[2],
-                                const float source[3],
+static void uv_cylinder_project(BMFace *efa,
                                 const float center[3],
-                                const float rotmat[4][4])
+                                const float rotmat[3][3],
+                                const bool fan,
+                                const int cd_loop_uv_offset)
 {
-  float pv[3];
-
-  sub_v3_v3v3(pv, source, center);
-  mul_m4_v3(rotmat, pv);
-
-  map_to_tube(&target[0], &target[1], pv[0], pv[1], pv[2]);
-
-  /* split line is always zero */
-  if (target[0] >= 1.0f) {
-    target[0] -= 1.0f;
+  bool *regular = BLI_array_alloca(regular, efa->len);
+  int i;
+  BMLoop *l;
+  BMIter iter;
+  BM_ITER_ELEM_INDEX (l, &iter, efa, BM_LOOPS_OF_FACE, i) {
+    MLoopUV *luv = BM_ELEM_CD_GET_VOID_P(l, cd_loop_uv_offset);
+    float pv[3];
+    sub_v3_v3v3(pv, l->v->co, center);
+    mul_m3_v3(rotmat, pv);
+    regular[i] = map_to_tube(&luv->uv[0], &luv->uv[1], pv[0], pv[1], pv[2]);
   }
+
+  uv_map_mirror(efa, regular, fan, cd_loop_uv_offset);
 }
 
 static int cylinder_project_exec(bContext *C, wmOperator *op)
@@ -2912,9 +2971,7 @@ static int cylinder_project_exec(bContext *C, wmOperator *op)
     Object *obedit = objects[ob_index];
     BMEditMesh *em = BKE_editmesh_from_object(obedit);
     BMFace *efa;
-    BMLoop *l;
-    BMIter iter, liter;
-    MLoopUV *luv;
+    BMIter iter;
 
     if (em->bm->totfacesel == 0) {
       continue;
@@ -2926,11 +2983,13 @@ static int cylinder_project_exec(bContext *C, wmOperator *op)
     }
 
     const int cd_loop_uv_offset = CustomData_get_offset(&em->bm->ldata, CD_MLOOPUV);
-    float center[3], rotmat[4][4];
+    float center[3], rotmat[3][3];
 
     uv_map_transform(C, op, rotmat);
     uv_map_transform_center(scene, v3d, obedit, em, center, NULL);
 
+    const bool fan = RNA_enum_get(op->ptr, "pole");
+
     BM_ITER_MESH (efa, &iter, em->bm, BM_FACES_OF_MESH) {
       if (!BM_elem_flag_test(efa, BM_ELEM_SELECT)) {
         continue;
@@ -2941,12 +3000,7 @@ static int cylinder_project_exec(bContext *C, wmOperator *op)
         continue;
       }
 
-      BM_ITER_ELEM (l, &liter, efa, BM_LOOPS_OF_FACE) {
-        luv = BM_ELEM_CD_GET_VOID_P(l, cd_loop_uv_offset);
-        uv_cylinder_project(luv->uv, l->v->co, center, rotmat);
-      }
-
-      uv_map_mirror(em, efa);
+      uv_cylinder_project(efa, center, rotmat, fan, cd_loop_uv_offset);
     }
 
     const bool per_face_aspect = true;