diff --git a/Marlin/src/config/examples/Micromake/C1/README.md b/Marlin/src/config/examples/Micromake/C1/README.md
index 5e5f0afafaa60013c30b9dda839dc7300aba07ec..0111f6f00a0ee83f2447a52fbe7c7cee148006d7 100644
--- a/Marlin/src/config/examples/Micromake/C1/README.md
+++ b/Marlin/src/config/examples/Micromake/C1/README.md
@@ -13,4 +13,3 @@ Configuration files for Micromake C1 with…
   - 128 STEPS configured with jumper on the motherboard (all open for 128 Steps).
   - Capacitive Probe (Adjust offsets at your convenience)
   - French language with no accents for Japanese LCD.
-  
diff --git a/Marlin/src/config/examples/delta/FLSUN/auto_calibrate/Configuration.h b/Marlin/src/config/examples/delta/FLSUN/auto_calibrate/Configuration.h
index 57bb3c3fb3006fdc8ce25c04d53de7d25e65b8fc..17cff0ef91f3b88a0e0d430680368a6774129758 100644
--- a/Marlin/src/config/examples/delta/FLSUN/auto_calibrate/Configuration.h
+++ b/Marlin/src/config/examples/delta/FLSUN/auto_calibrate/Configuration.h
@@ -496,6 +496,12 @@
   #if ENABLED(DELTA_AUTO_CALIBRATION)
     // set the default number of probe points : n*n (1 -> 7)
     #define DELTA_CALIBRATION_DEFAULT_POINTS 4
+
+    // Enable and set these values based on results of 'G33 A1'
+    //#define H_FACTOR 1.01
+    //#define R_FACTOR 2.61
+    //#define A_FACTOR 0.87
+
   #endif
 
   #if ENABLED(DELTA_AUTO_CALIBRATION) || ENABLED(DELTA_CALIBRATION_MENU)
diff --git a/Marlin/src/config/examples/delta/FLSUN/kossel_mini/Configuration.h b/Marlin/src/config/examples/delta/FLSUN/kossel_mini/Configuration.h
index 65a34da45f6c0ff62f9f6adfacf62fb5807d8ce2..49ac6ede2afae5ce753a530027619ad13f8dc644 100644
--- a/Marlin/src/config/examples/delta/FLSUN/kossel_mini/Configuration.h
+++ b/Marlin/src/config/examples/delta/FLSUN/kossel_mini/Configuration.h
@@ -496,6 +496,12 @@
   #if ENABLED(DELTA_AUTO_CALIBRATION)
     // set the default number of probe points : n*n (1 -> 7)
     #define DELTA_CALIBRATION_DEFAULT_POINTS 4
+
+    // Enable and set these values based on results of 'G33 A1'
+    //#define H_FACTOR 1.01
+    //#define R_FACTOR 2.61
+    //#define A_FACTOR 0.87
+
   #endif
 
   #if ENABLED(DELTA_AUTO_CALIBRATION) || ENABLED(DELTA_CALIBRATION_MENU)
diff --git a/Marlin/src/config/examples/delta/generic/Configuration.h b/Marlin/src/config/examples/delta/generic/Configuration.h
index ebf8a73901307530ea97865fd32c1e1290dbe65b..2542f3f37067119ef7adfbcade271f9fbfa95010 100644
--- a/Marlin/src/config/examples/delta/generic/Configuration.h
+++ b/Marlin/src/config/examples/delta/generic/Configuration.h
@@ -486,6 +486,12 @@
   #if ENABLED(DELTA_AUTO_CALIBRATION)
     // set the default number of probe points : n*n (1 -> 7)
     #define DELTA_CALIBRATION_DEFAULT_POINTS 4
+
+    // Enable and set these values based on results of 'G33 A1'
+    //#define H_FACTOR 1.01
+    //#define R_FACTOR 2.61
+    //#define A_FACTOR 0.87
+
   #endif
 
   #if ENABLED(DELTA_AUTO_CALIBRATION) || ENABLED(DELTA_CALIBRATION_MENU)
diff --git a/Marlin/src/config/examples/delta/kossel_mini/Configuration.h b/Marlin/src/config/examples/delta/kossel_mini/Configuration.h
index 764ec46021192d0c588e66490527f2125f320549..94a74254d8163736a6e8188385453d840e9e0266 100644
--- a/Marlin/src/config/examples/delta/kossel_mini/Configuration.h
+++ b/Marlin/src/config/examples/delta/kossel_mini/Configuration.h
@@ -486,6 +486,12 @@
   #if ENABLED(DELTA_AUTO_CALIBRATION)
     // set the default number of probe points : n*n (1 -> 7)
     #define DELTA_CALIBRATION_DEFAULT_POINTS 4
+
+    // Enable and set these values based on results of 'G33 A1'
+    //#define H_FACTOR 1.01
+    //#define R_FACTOR 2.61
+    //#define A_FACTOR 0.87
+
   #endif
 
   #if ENABLED(DELTA_AUTO_CALIBRATION) || ENABLED(DELTA_CALIBRATION_MENU)
diff --git a/Marlin/src/config/examples/delta/kossel_pro/Configuration.h b/Marlin/src/config/examples/delta/kossel_pro/Configuration.h
index c1f1f56c0602616b8d10d506b9ab4b24cfe2d4ad..dce7a25b0209d84eb8f4179c91ad0178e4a3ab34 100644
--- a/Marlin/src/config/examples/delta/kossel_pro/Configuration.h
+++ b/Marlin/src/config/examples/delta/kossel_pro/Configuration.h
@@ -472,6 +472,12 @@
   #if ENABLED(DELTA_AUTO_CALIBRATION)
     // set the default number of probe points : n*n (1 -> 7)
     #define DELTA_CALIBRATION_DEFAULT_POINTS 4
+
+    // Enable and set these values based on results of 'G33 A1'
+    //#define H_FACTOR 1.01
+    //#define R_FACTOR 2.61
+    //#define A_FACTOR 0.87
+
   #endif
 
   #if ENABLED(DELTA_AUTO_CALIBRATION) || ENABLED(DELTA_CALIBRATION_MENU)
diff --git a/Marlin/src/config/examples/delta/kossel_xl/Configuration.h b/Marlin/src/config/examples/delta/kossel_xl/Configuration.h
index 24355c424c418c628e335c98e0f8b594e796e08a..ba73a95e59c2adb3236d409065aba841892fc419 100644
--- a/Marlin/src/config/examples/delta/kossel_xl/Configuration.h
+++ b/Marlin/src/config/examples/delta/kossel_xl/Configuration.h
@@ -490,6 +490,12 @@
   #if ENABLED(DELTA_AUTO_CALIBRATION)
     // set the default number of probe points : n*n (1 -> 7)
     #define DELTA_CALIBRATION_DEFAULT_POINTS 4
+
+    // Enable and set these values based on results of 'G33 A1'
+    //#define H_FACTOR 1.01
+    //#define R_FACTOR 2.61
+    //#define A_FACTOR 0.87
+
   #endif
 
   #if ENABLED(DELTA_AUTO_CALIBRATION) || ENABLED(DELTA_CALIBRATION_MENU)
diff --git a/Marlin/src/gcode/calibrate/G33.cpp b/Marlin/src/gcode/calibrate/G33.cpp
index 5746c3c4dedff7325cbec13f427f727d62428dc0..abb8eed4b275dff04db9a6875f28ad4a7382648f 100644
--- a/Marlin/src/gcode/calibrate/G33.cpp
+++ b/Marlin/src/gcode/calibrate/G33.cpp
@@ -37,35 +37,6 @@
   #include "../../feature/bedlevel/bedlevel.h"
 #endif
 
-/**
- * G33 - Delta '1-4-7-point' Auto-Calibration
- *       Calibrate height, endstops, delta radius, and tower angles.
- *
- * Parameters:
- *
- *   Pn  Number of probe points:
- *
- *      P0     No probe. Normalize only.
- *      P1     Probe center and set height only.
- *      P2     Probe center and towers. Set height, endstops, and delta radius.
- *      P3     Probe all positions: center, towers and opposite towers. Set all.
- *      P4-P7  Probe all positions at different locations and average them.
- *
- *   T0  Don't calibrate tower angle corrections
- *
- *   Cn.nn Calibration precision; when omitted calibrates to maximum precision
- *
- *   Fn  Force to run at least n iterations and takes the best result
- *
- *   Vn  Verbose level:
- *
- *      V0  Dry-run mode. Report settings and probe results. No calibration.
- *      V1  Report settings
- *      V2  Report settings and probe results
- *
- *   E   Engage the probe for each point
- */
-
 static void print_signed_float(const char * const prefix, const float &f) {
   SERIAL_PROTOCOLPGM("  ");
   serialprintPGM(prefix);
@@ -77,21 +48,55 @@ static void print_signed_float(const char * const prefix, const float &f) {
 static void print_G33_settings(const bool end_stops, const bool tower_angles) {
   SERIAL_PROTOCOLPAIR(".Height:", DELTA_HEIGHT + home_offset[Z_AXIS]);
   if (end_stops) {
-    print_signed_float(PSTR("  Ex"), delta_endstop_adj[A_AXIS]);
+    print_signed_float(PSTR("Ex"), delta_endstop_adj[A_AXIS]);
     print_signed_float(PSTR("Ey"), delta_endstop_adj[B_AXIS]);
     print_signed_float(PSTR("Ez"), delta_endstop_adj[C_AXIS]);
-    SERIAL_PROTOCOLPAIR("    Radius:", delta_radius);
   }
-  SERIAL_EOL();
+  if (end_stops && tower_angles) {
+    SERIAL_PROTOCOLPAIR("  Radius:", delta_radius);
+    SERIAL_EOL();
+    SERIAL_CHAR('.');
+    SERIAL_PROTOCOL_SP(13);
+  }
   if (tower_angles) {
-    SERIAL_PROTOCOLPGM(".Tower angle :  ");
     print_signed_float(PSTR("Tx"), delta_tower_angle_trim[A_AXIS]);
     print_signed_float(PSTR("Ty"), delta_tower_angle_trim[B_AXIS]);
     print_signed_float(PSTR("Tz"), delta_tower_angle_trim[C_AXIS]);
+  }
+  if ((!end_stops && tower_angles) || (end_stops && !tower_angles)) { // XOR
+    SERIAL_PROTOCOLPAIR("  Radius:", delta_radius);
+  }
+  SERIAL_EOL();
+}
+
+static void print_G33_results(const float z_at_pt[13], const bool tower_points, const bool opposite_points) {
+  SERIAL_PROTOCOLPGM(".    ");
+  print_signed_float(PSTR("c"), z_at_pt[0]);
+  if (tower_points) {
+    print_signed_float(PSTR(" x"), z_at_pt[1]);
+    print_signed_float(PSTR(" y"), z_at_pt[5]);
+    print_signed_float(PSTR(" z"), z_at_pt[9]);
+  }
+  if (tower_points && opposite_points) {
     SERIAL_EOL();
+    SERIAL_CHAR('.');
+    SERIAL_PROTOCOL_SP(13);
   }
+  if (opposite_points) {
+    print_signed_float(PSTR("yz"), z_at_pt[7]);
+    print_signed_float(PSTR("zx"), z_at_pt[11]);
+    print_signed_float(PSTR("xy"), z_at_pt[3]);
+  }
+  SERIAL_EOL();
 }
 
+/**
+ * After G33:
+ *  - Move to the print ceiling (DELTA_HOME_TO_SAFE_ZONE only)
+ *  - Stow the probe
+ *  - Restore endstops state
+ *  - Select the old tool, if needed
+ */
 static void G33_cleanup(
   #if HOTENDS > 1
     const uint8_t old_tool_index
@@ -107,6 +112,216 @@ static void G33_cleanup(
   #endif
 }
 
+static float probe_G33_points(float z_at_pt[13], const int8_t probe_points, const bool towers_set, const bool stow_after_each) {
+  const bool _0p_calibration      = probe_points == 0,
+             _1p_calibration      = probe_points == 1,
+             _4p_calibration      = probe_points == 2,
+             _4p_opposite_points  = _4p_calibration && !towers_set,
+             _7p_calibration      = probe_points >= 3 || probe_points == 0,
+             _7p_half_circle      = probe_points == 3,
+             _7p_double_circle    = probe_points == 5,
+             _7p_triple_circle    = probe_points == 6,
+             _7p_quadruple_circle = probe_points == 7,
+             _7p_intermed_points  = probe_points >= 4,
+             _7p_multi_circle     = probe_points >= 5;
+
+  #if DISABLED(PROBE_MANUALLY)
+    const float dx = (X_PROBE_OFFSET_FROM_EXTRUDER),
+                dy = (Y_PROBE_OFFSET_FROM_EXTRUDER);
+  #endif
+
+  for (uint8_t i = 0; i < COUNT(z_at_pt); i++) z_at_pt[i] = 0.0;
+
+  if (!_0p_calibration) {
+
+    if (!_7p_half_circle && !_7p_triple_circle) { // probe the center
+      #if ENABLED(PROBE_MANUALLY)
+        z_at_pt[0] += lcd_probe_pt(0, 0);
+      #else
+        z_at_pt[0] += probe_pt(dx, dy, stow_after_each, 1, false);
+      #endif
+    }
+
+    if (_7p_calibration) { // probe extra center points
+      for (int8_t axis = _7p_multi_circle ? COUNT(z_at_pt) - 2 : COUNT(z_at_pt) - 4; axis > 0; axis -= _7p_multi_circle ? 2 : 4) {
+        const float a = RADIANS(180 + 30 * axis), r = delta_calibration_radius * 0.1;
+        #if ENABLED(PROBE_MANUALLY)
+          z_at_pt[0] += lcd_probe_pt(cos(a) * r, sin(a) * r);
+        #else
+          z_at_pt[0] += probe_pt(cos(a) * r + dx, sin(a) * r + dy, stow_after_each, 1);
+        #endif
+      }
+      z_at_pt[0] /= float(_7p_double_circle ? 7 : probe_points);
+    }
+
+    if (!_1p_calibration) {  // probe the radius
+      bool zig_zag = true;
+      const uint8_t start = _4p_opposite_points ? 3 : 1,
+                    step = _4p_calibration ? 4 : _7p_half_circle ? 2 : 1;
+      for (uint8_t axis = start; axis < COUNT(z_at_pt); axis += step) {
+        const float zigadd = (zig_zag ? 0.5 : 0.0),
+                    offset_circles = _7p_quadruple_circle ? zigadd + 1.0 :
+                                     _7p_triple_circle    ? zigadd + 0.5 :
+                                     _7p_double_circle    ? zigadd : 0;
+        for (float circles = -offset_circles ; circles <= offset_circles; circles++) {
+          const float a = RADIANS(180 + 30 * axis),
+                      r = delta_calibration_radius * (1 + circles * (zig_zag ? 0.1 : -0.1));
+          #if ENABLED(PROBE_MANUALLY)
+            z_at_pt[axis] += lcd_probe_pt(cos(a) * r, sin(a) * r);
+          #else
+            z_at_pt[axis] += probe_pt(cos(a) * r + dx, sin(a) * r + dy, stow_after_each, 1);
+          #endif
+        }
+        zig_zag = !zig_zag;
+        z_at_pt[axis] /= (2 * offset_circles + 1);
+      }
+    }
+
+    if (_7p_intermed_points) // average intermediates to tower and opposites
+      for (uint8_t axis = 1; axis < COUNT(z_at_pt); axis += 2)
+        z_at_pt[axis] = (z_at_pt[axis] + (z_at_pt[axis + 1] + z_at_pt[(axis + 10) % 12 + 1]) / 2.0) / 2.0;
+
+    float S1 = z_at_pt[0],
+          S2 = sq(z_at_pt[0]);
+    int16_t N = 1;
+    if (!_1p_calibration) // std dev from zero plane
+      for (uint8_t axis = (_4p_opposite_points ? 3 : 1); axis < COUNT(z_at_pt); axis += (_4p_calibration ? 4 : 2)) {
+        S1 += z_at_pt[axis];
+        S2 += sq(z_at_pt[axis]);
+        N++;
+      }
+    return round(SQRT(S2 / N) * 1000.0) / 1000.0 + 0.00001;
+  }
+
+  return 0.00001;
+}
+
+#if DISABLED(PROBE_MANUALLY)
+
+  static void G33_auto_tune() {
+    float z_at_pt[13]      = { 0.0 },
+          z_at_pt_base[13] = { 0.0 },
+          z_temp, h_fac = 0.0, r_fac = 0.0, a_fac = 0.0, norm = 0.8;
+
+    #define ZP(N,I) ((N) * z_at_pt[I])
+    #define Z06(I)  ZP(6, I)
+    #define Z03(I)  ZP(3, I)
+    #define Z02(I)  ZP(2, I)
+    #define Z01(I)  ZP(1, I)
+    #define Z32(I)  ZP(3/2, I)
+
+    SERIAL_PROTOCOLPGM("AUTO TUNE baseline");
+    SERIAL_EOL();
+    probe_G33_points(z_at_pt_base, 3, true, false);
+    print_G33_results(z_at_pt_base, true, true);
+
+    LOOP_XYZ(axis) {
+      delta_endstop_adj[axis] -= 1.0;
+
+      endstops.enable(true);
+      if (!home_delta()) return;
+      endstops.not_homing();
+
+      SERIAL_PROTOCOLPGM("Tuning E");
+      SERIAL_CHAR(tolower(axis_codes[axis]));
+      SERIAL_EOL();
+
+      probe_G33_points(z_at_pt, 3, true, false);
+      for (int8_t i = 0; i < COUNT(z_at_pt); i++) z_at_pt[i] -= z_at_pt_base[i];
+      print_G33_results(z_at_pt, true, true);
+      delta_endstop_adj[axis] += 1.0;
+      switch (axis) {
+        case A_AXIS :
+          h_fac += 4.0 / (Z03(0) +Z01(1)                         +Z32(11) +Z32(3)); // Offset by X-tower end-stop
+          break;
+        case B_AXIS :
+          h_fac += 4.0 / (Z03(0)         +Z01(5)         +Z32(7)          +Z32(3)); // Offset by Y-tower end-stop
+          break;
+        case C_AXIS :
+          h_fac += 4.0 / (Z03(0)                 +Z01(9) +Z32(7) +Z32(11)        ); // Offset by Z-tower end-stop
+          break;
+      }
+    }
+    h_fac /= 3.0;
+    h_fac *= norm; // Normalize to 1.02 for Kossel mini
+
+    for (int8_t zig_zag = -1; zig_zag < 2; zig_zag += 2) {
+      delta_radius += 1.0 * zig_zag;
+      recalc_delta_settings(delta_radius, delta_diagonal_rod, delta_tower_angle_trim);
+
+      endstops.enable(true);
+      if (!home_delta()) return;
+      endstops.not_homing();
+
+      SERIAL_PROTOCOLPGM("Tuning R");
+      SERIAL_PROTOCOL(zig_zag == -1 ? "-" : "+");
+      SERIAL_EOL();
+      probe_G33_points(z_at_pt, 3, true, false);
+      for (int8_t i = 0; i < COUNT(z_at_pt); i++) z_at_pt[i] -= z_at_pt_base[i];
+      print_G33_results(z_at_pt, true, true);
+      delta_radius -= 1.0 * zig_zag;
+      recalc_delta_settings(delta_radius, delta_diagonal_rod, delta_tower_angle_trim);
+      r_fac -= zig_zag * 6.0 / (Z03(1) + Z03(5) + Z03(9) + Z03(7) + Z03(11) + Z03(3)); // Offset by delta radius
+    }
+    r_fac /= 2.0;
+    r_fac *= 3 * norm; // Normalize to 2.25 for Kossel mini
+
+    LOOP_XYZ(axis) {
+      delta_tower_angle_trim[axis] += 1.0;
+      delta_endstop_adj[(axis + 1) % 3] -= 1.0 / 4.5;
+      delta_endstop_adj[(axis + 2) % 3] += 1.0 / 4.5;
+      z_temp = MAX3(delta_endstop_adj[A_AXIS], delta_endstop_adj[B_AXIS], delta_endstop_adj[C_AXIS]);
+      home_offset[Z_AXIS] -= z_temp;
+      LOOP_XYZ(axis) delta_endstop_adj[axis] -= z_temp;
+      recalc_delta_settings(delta_radius, delta_diagonal_rod, delta_tower_angle_trim);
+
+      endstops.enable(true);
+      if (!home_delta()) return;
+      endstops.not_homing();
+
+      SERIAL_PROTOCOLPGM("Tuning T");
+      SERIAL_CHAR(tolower(axis_codes[axis]));
+      SERIAL_EOL();
+
+      probe_G33_points(z_at_pt, 3, true, false);
+      for (int8_t i = 0; i < COUNT(z_at_pt); i++) z_at_pt[i] -= z_at_pt_base[i];
+      print_G33_results(z_at_pt, true, true);
+
+      delta_tower_angle_trim[axis] -= 1.0;
+      delta_endstop_adj[(axis+1) % 3] += 1.0/4.5;
+      delta_endstop_adj[(axis+2) % 3] -= 1.0/4.5;
+      z_temp = MAX3(delta_endstop_adj[A_AXIS], delta_endstop_adj[B_AXIS], delta_endstop_adj[C_AXIS]);
+      home_offset[Z_AXIS] -= z_temp;
+      LOOP_XYZ(axis) delta_endstop_adj[axis] -= z_temp;
+      recalc_delta_settings(delta_radius, delta_diagonal_rod, delta_tower_angle_trim);
+      switch (axis) {
+        case A_AXIS :
+        a_fac += 4.0 / (       Z06(5) -Z06(9)         +Z06(11) -Z06(3)); // Offset by alpha tower angle
+        break;
+        case B_AXIS :
+        a_fac += 4.0 / (-Z06(1)       +Z06(9) -Z06(7)          +Z06(3)); // Offset by beta tower angle
+        break;
+        case C_AXIS :
+        a_fac += 4.0 / (Z06(1) -Z06(5)        +Z06(7) -Z06(11)        ); // Offset by gamma tower angle
+        break;
+      }
+    }
+    a_fac /= 3.0;
+    a_fac *= norm; // Normalize to 0.83 for Kossel mini
+
+    endstops.enable(true);
+    if (!home_delta()) return;
+    endstops.not_homing();
+    print_signed_float(PSTR( "H_FACTOR: "), h_fac);
+    print_signed_float(PSTR(" R_FACTOR: "), r_fac);
+    print_signed_float(PSTR(" A_FACTOR: "), a_fac);
+    SERIAL_EOL();
+    SERIAL_PROTOCOLPGM("Copy these values to Configuration.h");
+    SERIAL_EOL();
+  }
+
+#endif // !PROBE_MANUALLY
+
 /**
  * G33 - Delta '1-4-7-point' Auto-Calibration
  *       Calibrate height, endstops, delta radius, and tower angles.
@@ -114,21 +329,21 @@ static void G33_cleanup(
  * Parameters:
  *
  *   Pn  Number of probe points:
- *
  *      P0     No probe. Normalize only.
  *      P1     Probe center and set height only.
- *      P2     Probe center and towers. Set height, endstops, and delta radius.
+ *      P2     Probe center and towers. Set height, endstops and delta radius.
  *      P3     Probe all positions: center, towers and opposite towers. Set all.
  *      P4-P7  Probe all positions at different locations and average them.
  *
- *   T0  Don't calibrate tower angle corrections
+ *   T   Don't calibrate tower angle corrections
  *
- *   Cn.nn Calibration precision; when omitted calibrates to maximum precision
+ *   Cn.nn  Calibration precision; when omitted calibrates to maximum precision
  *
  *   Fn  Force to run at least n iterations and takes the best result
  *
- *   Vn  Verbose level:
+ *   A   Auto tune calibartion factors (set in Configuration.h)
  *
+ *   Vn  Verbose level:
  *      V0  Dry-run mode. Report settings and probe results. No calibration.
  *      V1  Report settings
  *      V2  Report settings and probe results
@@ -162,26 +377,24 @@ void GcodeSuite::G33() {
   }
 
   const bool towers_set           = !parser.boolval('T'),
+             auto_tune            = parser.boolval('A'),
              stow_after_each      = parser.boolval('E'),
              _0p_calibration      = probe_points == 0,
              _1p_calibration      = probe_points == 1,
              _4p_calibration      = probe_points == 2,
-             _4p_towers_points    = _4p_calibration && towers_set,
-             _4p_opposite_points  = _4p_calibration && !towers_set,
-             _7p_calibration      = probe_points >= 3 || _0p_calibration,
-             _7p_half_circle      = probe_points == 3,
+             _tower_results       = (_4p_calibration && towers_set)
+                                    || probe_points >= 3 || probe_points == 0,
+             _opposite_results    = (_4p_calibration && !towers_set)
+                                    || probe_points >= 3 || probe_points == 0,
+             _endstop_results     = probe_points != 1,
+             _angle_results       = (probe_points >= 3 || probe_points == 0) && towers_set,
              _7p_double_circle    = probe_points == 5,
              _7p_triple_circle    = probe_points == 6,
-             _7p_quadruple_circle = probe_points == 7,
-             _7p_multi_circle     = _7p_double_circle || _7p_triple_circle || _7p_quadruple_circle,
-             _7p_intermed_points  = _7p_calibration && !_7p_half_circle;
+             _7p_quadruple_circle = probe_points == 7;
   const static char save_message[] PROGMEM = "Save with M500 and/or copy to Configuration.h";
-  const float dx = (X_PROBE_OFFSET_FROM_EXTRUDER),
-              dy = (Y_PROBE_OFFSET_FROM_EXTRUDER);
   int8_t iterations = 0;
   float test_precision,
         zero_std_dev = (verbose_level ? 999.0 : 0.0), // 0.0 in dry-run mode : forced end
-        zero_std_dev_old = zero_std_dev,
         zero_std_dev_min = zero_std_dev,
         e_old[ABC] = {
           delta_endstop_adj[A_AXIS],
@@ -196,12 +409,14 @@ void GcodeSuite::G33() {
           delta_tower_angle_trim[C_AXIS]
         };
 
+  SERIAL_PROTOCOLLNPGM("G33 Auto Calibrate");
+
   if (!_1p_calibration && !_0p_calibration) {  // test if the outer radius is reachable
     const float circles = (_7p_quadruple_circle ? 1.5 :
                            _7p_triple_circle    ? 1.0 :
                            _7p_double_circle    ? 0.5 : 0),
                 r = (1 + circles * 0.1) * delta_calibration_radius;
-    for (uint8_t axis = 1; axis < 13; ++axis) {
+    for (uint8_t axis = 1; axis <= 12; ++axis) {
       const float a = RADIANS(180 + 30 * axis);
       if (!position_is_reachable_xy(cos(a) * r, sin(a) * r)) {
         SERIAL_PROTOCOLLNPGM("?(M665 B)ed radius is implausible.");
@@ -209,7 +424,6 @@ void GcodeSuite::G33() {
       }
     }
   }
-  SERIAL_PROTOCOLLNPGM("G33 Auto Calibrate");
 
   stepper.synchronize();
   #if HAS_LEVELING
@@ -232,7 +446,17 @@ void GcodeSuite::G33() {
     endstops.not_homing();
   }
 
-  // print settings
+  if (auto_tune) {
+    #if ENABLED(PROBE_MANUALLY)
+      SERIAL_PROTOCOLLNPGM("A probe is needed for auto-tune");
+    #else
+      G33_auto_tune();
+    #endif
+    G33_CLEANUP();
+    return;
+  }
+
+  // Report settings
 
   const char *checkingac = PSTR("Checking... AC"); // TODO: Make translatable string
   serialprintPGM(checkingac);
@@ -240,78 +464,19 @@ void GcodeSuite::G33() {
   SERIAL_EOL();
   lcd_setstatusPGM(checkingac);
 
-  print_G33_settings(!_1p_calibration, _7p_calibration && towers_set);
+  print_G33_settings(_endstop_results, _angle_results);
 
   do {
 
     float z_at_pt[13] = { 0.0 };
 
-    test_precision = zero_std_dev_old != 999.0 ? (zero_std_dev + zero_std_dev_old) / 2 : zero_std_dev;
+    test_precision = zero_std_dev;
 
     iterations++;
 
     // Probe the points
 
-    if (!_0p_calibration){
-      if (!_7p_half_circle && !_7p_triple_circle) { // probe the center
-        #if ENABLED(PROBE_MANUALLY)
-          z_at_pt[0] += lcd_probe_pt(0, 0);
-        #else
-          z_at_pt[0] += probe_pt(dx, dy, stow_after_each, 1, false);
-          if (isnan(z_at_pt[0])) return G33_CLEANUP();
-        #endif
-      }
-      if (_7p_calibration) { // probe extra center points
-        for (int8_t axis = _7p_multi_circle ? 11 : 9; axis > 0; axis -= _7p_multi_circle ? 2 : 4) {
-          const float a = RADIANS(180 + 30 * axis), r = delta_calibration_radius * 0.1;
-          #if ENABLED(PROBE_MANUALLY)
-            z_at_pt[0] += lcd_probe_pt(cos(a) * r, sin(a) * r);
-          #else
-            z_at_pt[0] += probe_pt(cos(a) * r + dx, sin(a) * r + dy, stow_after_each, 1);
-            if (isnan(z_at_pt[0])) return G33_CLEANUP();
-          #endif
-        }
-        z_at_pt[0] /= float(_7p_double_circle ? 7 : probe_points);
-      }
-      if (!_1p_calibration) {  // probe the radius
-        bool zig_zag = true;
-        const uint8_t start = _4p_opposite_points ? 3 : 1,
-                       step = _4p_calibration ? 4 : _7p_half_circle ? 2 : 1;
-        for (uint8_t axis = start; axis < 13; axis += step) {
-          const float zigadd = (zig_zag ? 0.5 : 0.0),
-                      offset_circles = _7p_quadruple_circle ? zigadd + 1.0 :
-                                       _7p_triple_circle    ? zigadd + 0.5 :
-                                       _7p_double_circle    ? zigadd : 0;
-          for (float circles = -offset_circles ; circles <= offset_circles; circles++) {
-            const float a = RADIANS(180 + 30 * axis),
-                        r = delta_calibration_radius * (1 + circles * (zig_zag ? 0.1 : -0.1));
-            #if ENABLED(PROBE_MANUALLY)
-              z_at_pt[axis] += lcd_probe_pt(cos(a) * r, sin(a) * r);
-            #else
-              z_at_pt[axis] += probe_pt(cos(a) * r + dx, sin(a) * r + dy, stow_after_each, 1);
-              if (isnan(z_at_pt[axis])) return G33_CLEANUP();
-            #endif
-          }
-          zig_zag = !zig_zag;
-          z_at_pt[axis] /= (2 * offset_circles + 1);
-        }
-      }
-      if (_7p_intermed_points) // average intermediates to tower and opposites
-        for (uint8_t axis = 1; axis < 13; axis += 2)
-          z_at_pt[axis] = (z_at_pt[axis] + (z_at_pt[axis + 1] + z_at_pt[(axis + 10) % 12 + 1]) / 2.0) / 2.0;
-    }
-
-    float S1 = z_at_pt[0],
-          S2 = sq(z_at_pt[0]);
-    int16_t N = 1;
-    if (!_1p_calibration) // std dev from zero plane
-      for (uint8_t axis = (_4p_opposite_points ? 3 : 1); axis < 13; axis += (_4p_calibration ? 4 : 2)) {
-        S1 += z_at_pt[axis];
-        S2 += sq(z_at_pt[axis]);
-        N++;
-      }
-    zero_std_dev_old = zero_std_dev;
-    zero_std_dev = round(SQRT(S2 / N) * 1000.0) / 1000.0 + 0.00001;
+    zero_std_dev = probe_G33_points(z_at_pt, probe_points, towers_set, stow_after_each);
 
     // Solve matrices
 
@@ -325,9 +490,24 @@ void GcodeSuite::G33() {
 
       float e_delta[ABC] = { 0.0 }, r_delta = 0.0, t_delta[ABC] = { 0.0 };
       const float r_diff = delta_radius - delta_calibration_radius,
-                  h_factor = (1.00 + r_diff * 0.001) / 6.0,                                       // 1.02 for r_diff = 20mm
-                  r_factor = (-(1.75 + 0.005 * r_diff + 0.001 * sq(r_diff))) / 6.0,               // 2.25 for r_diff = 20mm
-                  a_factor = (66.66 / delta_calibration_radius) / (iterations == 1 ? 16.0 : 2.0); // 0.83 for cal_rd = 80mm
+                  h_factor = 1 / 6.0 *
+                    #ifdef H_FACTOR
+                      (H_FACTOR),                                       // Set in Configuration.h
+                    #else
+                      (1.00 + r_diff * 0.001),                          // 1.02 for r_diff = 20mm
+                    #endif
+                  r_factor = 1 / 6.0 *
+                    #ifdef R_FACTOR
+                      -(R_FACTOR),                                      // Set in Configuration.h
+                    #else
+                      -(1.75 + 0.005 * r_diff + 0.001 * sq(r_diff)),    // 2.25 for r_diff = 20mm
+                    #endif
+                  a_factor = 1 / 6.0 *
+                    #ifdef A_FACTOR
+                      (A_FACTOR);                                       // Set in Configuration.h
+                    #else
+                      (66.66 / delta_calibration_radius);               // 0.83 for cal_rd = 80mm
+                    #endif
 
       #define ZP(N,I) ((N) * z_at_pt[I])
       #define Z6(I) ZP(6, I)
@@ -341,15 +521,11 @@ void GcodeSuite::G33() {
 
       switch (probe_points) {
         case 0:
-          #if DISABLED(PROBE_MANUALLY)
-            test_precision = 0.00; // forced end
-          #endif
+          test_precision = 0.00; // forced end
           break;
 
         case 1:
-          #if DISABLED(PROBE_MANUALLY)
-            test_precision = 0.00; // forced end
-          #endif
+          test_precision = 0.00; // forced end
           LOOP_XYZ(axis) e_delta[axis] = Z1(0);
           break;
 
@@ -375,9 +551,9 @@ void GcodeSuite::G33() {
           r_delta         = (Z6(0) - Z1(1) - Z1(5) - Z1(9) - Z1(7) - Z1(11) - Z1(3)) * r_factor;
 
           if (towers_set) {
-            t_delta[A_AXIS] = (       - Z2(5) + Z2(9)         - Z2(11) + Z2(3)) * a_factor;
-            t_delta[B_AXIS] = ( Z2(1)         - Z2(9) + Z2(7)          - Z2(3)) * a_factor;
-            t_delta[C_AXIS] = (-Z2(1) + Z2(5)         - Z2(7) + Z2(11)        ) * a_factor;
+            t_delta[A_AXIS] = (       - Z4(5) + Z4(9)         - Z4(11) + Z4(3)) * a_factor;
+            t_delta[B_AXIS] = ( Z4(1)         - Z4(9) + Z4(7)          - Z4(3)) * a_factor;
+            t_delta[C_AXIS] = (-Z4(1) + Z4(5)         - Z4(7) + Z4(11)        ) * a_factor;
             e_delta[A_AXIS] += (t_delta[B_AXIS] - t_delta[C_AXIS]) / 4.5;
             e_delta[B_AXIS] += (t_delta[C_AXIS] - t_delta[A_AXIS]) / 4.5;
             e_delta[C_AXIS] += (t_delta[A_AXIS] - t_delta[B_AXIS]) / 4.5;
@@ -395,11 +571,14 @@ void GcodeSuite::G33() {
       home_offset[Z_AXIS] = zh_old;
       COPY(delta_tower_angle_trim, ta_old);
     }
+
     if (verbose_level != 0) {                                    // !dry run
       // normalise angles to least squares
-      float a_sum = 0.0;
-      LOOP_XYZ(axis) a_sum += delta_tower_angle_trim[axis];
-      LOOP_XYZ(axis) delta_tower_angle_trim[axis] -= a_sum / 3.0;
+      if (_angle_results) {
+        float a_sum = 0.0;
+        LOOP_XYZ(axis) a_sum += delta_tower_angle_trim[axis];
+        LOOP_XYZ(axis) delta_tower_angle_trim[axis] -= a_sum / 3.0;
+      }
 
       // adjust delta_height and endstops by the max amount
       const float z_temp = MAX3(delta_endstop_adj[A_AXIS], delta_endstop_adj[B_AXIS], delta_endstop_adj[C_AXIS]);
@@ -411,30 +590,13 @@ void GcodeSuite::G33() {
 
     // print report
 
-    if (verbose_level != 1) {
-      SERIAL_PROTOCOLPGM(".    ");
-      print_signed_float(PSTR("c"), z_at_pt[0]);
-      if (_4p_towers_points || _7p_calibration) {
-        print_signed_float(PSTR("   x"), z_at_pt[1]);
-        print_signed_float(PSTR(" y"), z_at_pt[5]);
-        print_signed_float(PSTR(" z"), z_at_pt[9]);
-      }
-      if (!_4p_opposite_points) SERIAL_EOL();
-      if ((_4p_opposite_points) || _7p_calibration) {
-        if (_7p_calibration) {
-          SERIAL_CHAR('.');
-          SERIAL_PROTOCOL_SP(13);
-        }
-        print_signed_float(PSTR("  yz"), z_at_pt[7]);
-        print_signed_float(PSTR("zx"), z_at_pt[11]);
-        print_signed_float(PSTR("xy"), z_at_pt[3]);
-        SERIAL_EOL();
-      }
-    }
+    if (verbose_level != 1)
+      print_G33_results(z_at_pt, _tower_results, _opposite_results);
+
     if (verbose_level != 0) {                                    // !dry run
       if ((zero_std_dev >= test_precision && iterations > force_iterations) || zero_std_dev <= calibration_precision) {  // end iterations
         SERIAL_PROTOCOLPGM("Calibration OK");
-        SERIAL_PROTOCOL_SP(36);
+        SERIAL_PROTOCOL_SP(32);
         #if DISABLED(PROBE_MANUALLY)
           if (zero_std_dev >= test_precision && !_1p_calibration)
             SERIAL_PROTOCOLPGM("rolling back.");
@@ -452,7 +614,7 @@ void GcodeSuite::G33() {
         else
           sprintf_P(&mess[15], PSTR("%03i.x"), (int)round(zero_std_dev_min));
         lcd_setstatus(mess);
-        print_G33_settings(!_1p_calibration, _7p_calibration && towers_set);
+        print_G33_settings(_endstop_results, _angle_results);
         serialprintPGM(save_message);
         SERIAL_EOL();
       }
@@ -463,18 +625,18 @@ void GcodeSuite::G33() {
         else
           sprintf_P(mess, PSTR("No convergence"));
         SERIAL_PROTOCOL(mess);
-        SERIAL_PROTOCOL_SP(36);
+        SERIAL_PROTOCOL_SP(32);
         SERIAL_PROTOCOLPGM("std dev:");
         SERIAL_PROTOCOL_F(zero_std_dev, 3);
         SERIAL_EOL();
         lcd_setstatus(mess);
-        print_G33_settings(!_1p_calibration, _7p_calibration && towers_set);
+        print_G33_settings(_endstop_results, _angle_results);
       }
     }
     else {                                                       // dry run
       const char *enddryrun = PSTR("End DRY-RUN");
       serialprintPGM(enddryrun);
-      SERIAL_PROTOCOL_SP(39);
+      SERIAL_PROTOCOL_SP(35);
       SERIAL_PROTOCOLPGM("std dev:");
       SERIAL_PROTOCOL_F(zero_std_dev, 3);
       SERIAL_EOL();
@@ -490,7 +652,8 @@ void GcodeSuite::G33() {
     }
 
     endstops.enable(true);
-    home_delta();
+    if (!home_delta())
+      return;
     endstops.not_homing();
 
   }