Add unpush command with primitive canonicalization

bu_parallel_upgrade.txt

@@ -1,324 +0,0 @@
-bu_parallel recursive scheduling upgrade
-========================================
-
-Status
-------
-
-This document proposes a replacement for the recursive scheduling and thread
-ID bookkeeping in src/libbu/parallel.c.  It records the constraints behind the
-current implementation, the recommended C++17 design, and an incremental plan
-for implementing and validating it.
-
-
-Confirmed source of the Valgrind warnings
------------------------------------------
-
-The following declaration is the source of both observed Valgrind messages:
-
-    static struct parallel_info mapping[MAX_PSW*MAX_PSW];
-
-MAX_PSW is currently 4096 and struct parallel_info is 32 bytes on this LP64
-build.  The array therefore contributes exactly 512 MiB to libbu's BSS.  The
-complete libbu BSS is 563,028,096 bytes in the observed build.
-
-Preloading only libbu.so.20 into /bin/true under Valgrind reproduces both
-"set address range perms: large range" warnings while reporting zero bytes in
-use at exit.  The dynamic loader emits the warnings while establishing memory
-permissions on libbu's unusually large RW ELF segment.  They are not reports
-of heap leaks.
-
-
-Why the array is squared
-------------------------
-
-Commit 145c80ce922fae38105c01aacae1da690e0955c1 deliberately changed the
-mapping from MAX_PSW entries to MAX_PSW*MAX_PSW entries.  It was an explicitly
-temporary attempt to accommodate recursive bu_parallel() calls.
-
-The present algorithm uses a parallel_info record both as a live thread record
-and as the concurrency scope for a nested invocation.  If an outer invocation
-starts N callbacks, each callback can call bu_parallel(..., 0) and inherit an
-independent limit of N.  The implementation can consequently create about N*N
-simultaneous parent and child records rather than sharing one limit of N.
-
-There are two additional historical constraints:
-
-* A parent ID cannot be released when a nested bu_parallel() returns if the
-  parent callback is still alive.  It may retain per-ID state and may invoke
-  bu_parallel() again.  Commit 988887e815 documents this case.
-
-* bu_parallel_id() values are used throughout BRL-CAD to index MAX_PSW-sized
-  lock-free resource arrays.  The value may be arbitrary, but it must remain in
-  the range [0, MAX_PSW).  Commit af04ca7bb4 documents this requirement.
-
-The squared mapping postpones exhaustion for one level of recursion, but it
-does not solve arbitrary nesting and it permits IDs that violate the public
-MAX_PSW array contract.  Increasing MAX_PSW from 1024 to 4096 in commit
-9aa03acf836 made the 48/64-core recursive test fit, but also increased the
-squared mapping from about 32 MiB to 512 MiB.
-
-Reducing the mapping to MAX_PSW without changing the scheduling model is not a
-valid fix.  It restores the failure that the 2016 change was trying to avoid.
-
-
-Required invariants
--------------------
-
-The replacement must maintain the following invariants:
-
-1. Inside a callback, func_cpu_id and bu_parallel_id() identify the same
-   logical resource slot.
-
-2. Every live logical callback context has an ID in [0, MAX_PSW), and no two
-   simultaneously live contexts share an ID.
-
-3. A callback retains its logical ID across a nested bu_parallel() call and
-   until the callback itself returns.  Nested callbacks receive other IDs.
-
-4. bu_parallel() remains synchronous: it does not return until all callbacks
-   requested by that invocation have completed.
-
-5. ncpu == 0 inherits the enclosing executor's concurrency width instead of
-   multiplying it at every recursive level.
-
-6. Recursive fork/join cannot deadlock merely because all physical workers are
-   currently executing parent callbacks.
-
-7. Explicit requests are bounded by MAX_PSW.  Exhaustion must produce a clear,
-   deterministic failure before any out-of-range ID is exposed.
-
-8. Bookkeeping memory is proportional to actual workers, active calls, and
-   nesting depth.  It must not be proportional to MAX_PSW squared.
-
-9. Existing C callers and the exported C ABI remain compatible.
-
-
-Recommended design: a C++17 fork/join executor
-----------------------------------------------
-
-C++17 makes it practical to replace the platform-specific recursive thread
-creation and manual counters with a small fork/join executor using std::thread,
-std::mutex, std::condition_variable, atomics, and RAII.
-
-The preferred organization is to put the scheduler in a C++ implementation
-file and retain C-callable entry points.  It should not be necessary to convert
-all of parallel.c immediately.  CPU discovery and platform helpers can remain
-in C while the current parallel_cpp11thread.cpp is replaced or superseded by a
-new parallel_executor.cpp.  If moving bu_parallel() cleanly requires renaming
-parallel.c to parallel.cpp, its exported functions must retain C linkage.
-
-Do not enable the existing CPP11THREAD implementation as-is.  It does not
-implement BRL-CAD logical IDs, affinity behavior, recursive limits, nested
-fork/join progress, or the MAX_PSW resource contract.
-
-The executor should contain these conceptual components:
-
-1. Physical executor
-
-   An outermost bu_parallel() creates an executor with a bounded number of
-   physical workers.  An explicit ncpu selects the requested width, clamped to
-   MAX_PSW.  A zero value uses bu_avail_cpus().  A nested invocation reuses the
-   executor already recorded in thread-local state and never creates N more OS
-   threads for each parent callback.
-
-   Starting with an executor owned by each outermost call avoids introducing a
-   process-lifetime thread pool, complicated fork behavior, and new shutdown
-   ordering.  A persistent process-wide executor can be considered later if
-   measurements justify it.
-
-2. Logical ID pool
-
-   Maintain a process-wide pool of exactly MAX_PSW logical resource IDs.  IDs
-   are acquired immediately before a callback begins and released immediately
-   after that callback returns.  A free list or bitmap avoids scanning a large
-   mapping table.
-
-   A nested callback running on the same physical std::thread as its suspended
-   parent must acquire a different logical ID.  An RAII guard saves the
-   physical thread's previous thread-local ID, installs the child's ID, and
-   restores the parent ID when the child returns.  The parent ID remains
-   reserved for the lifetime of the suspended parent stack.
-
-   If nesting consumes all MAX_PSW logical IDs, bu_parallel() must fail with a
-   diagnostic rather than wait forever or index beyond a caller's resource
-   array.  The API currently returns void, so the existing bu_bomb-style hard
-   failure is the compatible behavior unless a future API adds error returns.
-
-3. Parallel call state
-
-   Each invocation has dynamically allocated call state containing at least:
-
-       callback and data
-       requested callback count
-       next unclaimed callback
-       number of incomplete callbacks
-       parent/root executor association
-       completion mutex/condition
-       failure state
-
-   Callback work should be claimed lazily.  Do not allocate or enqueue N*N job
-   objects for a recursive N-by-N call.  A bounded number of runnable tokens
-   can refer to one call state while workers atomically claim its remaining
-   callback invocations.  Memory then remains proportional to executor width
-   and the number of active calls.
-
-4. Helping wait
-
-   A non-worker caller may sleep on the call's completion condition.  A worker
-   making a nested bu_parallel() call must not simply block, since all workers
-   could otherwise block waiting for children that have no thread on which to
-   run.
-
-   Instead, a waiting worker repeatedly executes available work until its
-   nested call completes.  This is the standard fork/join "help while waiting"
-   rule.  It keeps the fixed worker set making progress at arbitrary nesting
-   depth and replaces the current one-second polling loop.
-
-   Every re-entered callback gets its own logical ID guard, even when it runs
-   on the same physical std::thread as a suspended parent.
-
-5. Affinity
-
-   Affinity belongs to physical workers, not logical callback IDs.  Set it once
-   when a worker starts, using a physical worker ordinal.  Logical IDs may be
-   larger than the available CPU count and may change as one worker helps
-   nested calls.
-
-6. Callback ID semantics
-
-   The callback's func_cpu_id remains a logical resource ID and must match
-   bu_parallel_id().  It is not a stable callback ordinal and is not guaranteed
-   to start at zero or be contiguous.  Documentation that promises increasing
-   callback numbers starting at zero should be corrected.  Code needing a work
-   ordinal should claim work from its shared data, as existing self-dispatching
-   callbacks generally do.
-
-
-Why this fixes the recursive case
----------------------------------
-
-For an N-way outer call in which every callback makes an N-way nested call:
-
-* The executor still has N physical workers, not N*N workers.
-* All N*N nested callback invocations are performed, but workers claim them
-  incrementally while waiting for their own nested joins.
-* Each suspended outer callback retains its logical ID.
-* Each currently executing nested callback obtains another logical ID.
-* The usual two-level case needs on the order of 2*N live logical IDs rather
-  than N*N mapping entries.
-* Deeper nesting consumes logical IDs according to live stack depth.  The hard
-  MAX_PSW bound remains explicit and safe instead of silently corrupting a
-  per-ID array.
-
-
-Optional short-term containment
--------------------------------
-
-If the loader warnings must be removed before the executor is ready, replace
-the monolithic squared array with stable, lazily allocated blocks while
-preserving its current logical capacity.  For example, a static top-level
-array of block pointers can address blocks of 1024 parallel_info records,
-allocating a block only when an ID in that block is first needed.  Individual
-blocks must never move because thread_context currently retains pointers to
-parallel_info records.
-
-A free list or high-water hint should accompany the blocks so ID allocation
-does not scan the entire MAX_PSW*MAX_PSW address space.  This would reduce
-libbu's initial BSS and eliminate the two Valgrind loader warnings without
-reintroducing the immediate recursive-capacity regression.
-
-This is containment only.  It retains the possibility of returning an ID
-greater than or equal to MAX_PSW and does not bound recursive thread creation.
-It must not be presented as the correctness fix or allowed to delay the ID and
-scheduler work above.
-
-
-Implementation plan
--------------------
-
-Phase 1: characterize and strengthen tests
-
-* Add instrumentation available in test builds to report the maximum logical
-  ID, number of simultaneously live logical contexts, physical worker count,
-  and nesting depth.
-* Extend test_parallel with two- and three-level recursion.
-* Test repeated nested calls from one parent, preserving the case described by
-  commit 988887e815.
-* Exercise ncpu values 1, 2, the available CPU count, 48, and 64 where the host
-  permits them without making every test host create thousands of threads.
-* Add concurrent independent outermost bu_parallel() calls from ordinary
-  application threads.
-* Assert that every observed callback ID is less than MAX_PSW and that
-  func_cpu_id equals bu_parallel_id().
-
-Phase 2: implement the executor behind an internal switch
-
-* Add the C++17 executor and logical ID pool without deleting the legacy path.
-* Keep both implementations selectable at build time or by a test-only switch
-  so identical workloads can be compared.
-* Use RAII for logical IDs, TLS restoration, mutex ownership, worker joining,
-  and cleanup after thread creation failures.
-* Catch std::system_error from worker creation, stop accepting work, join all
-  workers already created, and issue the normal libbu fatal diagnostic.
-
-Phase 3: implement recursive helping and lazy work claims
-
-* Make nested calls reuse their current executor.
-* Implement help-while-waiting and completion notification.
-* Generate callback work lazily so recursive fan-out does not allocate a job
-  object per callback.
-* Confirm that nested explicit ncpu requests control the number of callback
-  invocations but do not create a second unbounded set of OS threads.
-
-Phase 4: switch the default and remove legacy mapping
-
-* Run the complete libbu, librt, MGED, and ray-tracing regression suites using
-  the new executor.
-* Stress under ThreadSanitizer and run the recursive tests repeatedly to expose
-  bookkeeping races.
-* Compare representative ray-tracing and BREP-preparation performance against
-  the legacy implementation.
-* Make the executor the default, then remove parallel_mapping(),
-  parallel_wait_for_slot(), the squared mapping, and obsolete platform thread
-  creation branches.
-* Remove the Valgrind warning FIXME once libbu's BSS no longer contains the
-  squared array, and verify by preloading libbu into /bin/true under Valgrind.
-
-Phase 5: documentation and cleanup
-
-* Update bu/parallel.h to describe arbitrary bounded logical IDs, executor
-  reuse during recursion, and explicit failure at MAX_PSW exhaustion.
-* Update the TODO item concerning recursive bu_parallel limits and known races.
-* Decide whether the old parallel_cpp11thread.cpp should be replaced, renamed,
-  or removed to avoid two competing implementations.
-* Consider a persistent executor only as a separate performance project after
-  fork behavior, library teardown, and application concurrency are specified.
-
-
-Validation checklist
---------------------
-
-The upgrade is ready to replace the legacy path only when all of the following
-are true:
-
-* No callback observes an ID outside [0, MAX_PSW).
-* No simultaneously live callback contexts share an ID.
-* A parent retains its ID before and after every nested call.
-* Two- and three-level recursive tests complete without deadlock.
-* Recursive ncpu == 0 never increases the physical worker count.
-* Multiple independent outer calls complete without ID collisions.
-* Thread creation and logical ID exhaustion fail cleanly.
-* ThreadSanitizer reports no scheduler bookkeeping races.
-* Existing libbu and application regression tests pass.
-* Representative performance is not materially worse than the legacy path.
-* libbu no longer contains the 512 MiB mapping BSS.
-* Preloading libbu under Valgrind produces neither large-range warning.
-
-
-Non-goals
----------
-
-This upgrade should not initially attempt to remove MAX_PSW-sized resource
-arrays throughout BRL-CAD, change the public callback signature, introduce a
-permanent global thread pool, or guarantee a callback ordinal through
-func_cpu_id.  Those are separate API and architecture projects.

doc/asciidoc/system/mann/CMakeLists.txt

@@ -371,6 +371,7 @@ set(
   translate.adoc
   tree.adoc
   unhide.adoc
+  unpush.adoc
   units.adoc
   v2m_point.adoc
   vars.adoc

doc/asciidoc/system/mann/unpush.adoc

@@ -0,0 +1,140 @@
+= UNPUSH(nged)
+:doctype: manpage
+:mansource: BRL-CAD
+:manmanual: BRL-CAD MGED Commands
+
+== NAME
+unpush - infer and restore shared geometry from pushed object trees
+
+== SYNOPSIS
+
+[source]
+----
+unpush [-D] [-L] [-v] [-m rigid|similarity|affine] object [object ...]
+----
+
+
+[[description]]
+== DESCRIPTION
+
+The *unpush* command examines the trees rooted at the specified _objects_ and
+looks for primitive shapes and combinations that differ only by a recoverable
+placement transform.  It can replace those copies with shared objects and move
+their differing placements into combination leaf matrices.  This is a
+geometric inference operation: it does not require provenance or other
+metadata from an earlier *push*, *xpush*, or *npush* operation.
+
+Candidate primitives are first canonicalized into a stable local frame.
+Candidates are accepted as a group only after their attributes and geometry
+match within the database tolerance and the proposed replacement has been
+reconstructed under every leaf matrix that will be changed.  Combination and
+region candidates are analyzed from the leaves upward.  Their boolean tree
+structure, operand order, child identities and relative matrices, region
+properties, material and rendering fields, and custom attributes must agree.
+Boolean operands are not reordered.
+
+The command currently canonicalizes halfspaces, ellipsoids and spheres, tori,
+elliptical tori, right parabolic and hyperbolic cylinders, elliptical
+paraboloids and hyperboloids, particles, TGC and REC shapes, ARB4 through ARB8
+shapes, BOT meshes, and superellipsoids.  Other primitive types are reported
+as unsupported and are left unchanged.
+Invalid supported shapes are reported as failures and are also left unchanged.
+
+
+[[modes]]
+== TRANSFORM MODES
+
+The *-m* option limits the class of transforms that may be recovered:
+
+* *rigid* permits translation and rotation.
+* *similarity* additionally permits uniform scale.
+* *affine* permits primitive canonicalizers to consider nonuniform scale and
+  shear.
+
+Dry-run analysis defaults to *affine*.  Database writes constrain an affine
+request to *similarity*, because BRL-CAD combination matrices do not accept
+general shear.  The report explicitly notes this constraint.
+
+
+[[options]]
+== OPTIONS
+
+*-D*, *--dry-run*::
+Analyze and report possible sharing without modifying the database.  Running a
+dry run first is strongly recommended.
+
+*-L*, *--local*::
+Restrict changed references to the selected trees.  Objects still referenced
+outside those trees are retained, and descendants of a selected combination
+whose definition is itself exposed outside the selection are not changed.
+Without *-L*, *unpush* may update matching references elsewhere in the database
+so a redundant object can be removed globally.
+
+*-m* _mode_, *--mode* _mode_::
+Select *rigid*, *similarity*, or *affine* transform inference as described
+above.
+
+*-v*, *--verbosity*::
+Increase reporting detail.  Repeat the option to include replacement matrices
+for primitive references.
+
+*-h*, *--help*::
+Print command help.
+
+
+[[writes]]
+== DATABASE CHANGES AND RECOVERY
+
+Before writing, *unpush* accounts for database references and prepares every
+affected object in memory.  It writes new canonical primitives first, then
+updates parent combinations and top-level wrappers, reads the results back, and
+validates the exact rewritten leaves.  Superseded objects are deleted only
+after the updated reference counts show that they are unused.  If an object
+write or validation fails, the command restores the saved object definitions
+and removes newly written canonical primitives when it is safe to do so.
+
+Selected top-level primitive names are preserved by replacing those primitives
+with one-leaf wrapper combinations.  A selected top-level non-region
+combination may likewise become a one-leaf combination referencing the shared
+representative.  Top-level regions are never nested inside wrapper regions;
+such consolidation is deferred, although safe sharing below those regions may
+still be restored.
+
+The recovery process does not copy the complete database, which keeps it
+practical for multi-gigabyte models.  It protects against errors detected by
+the running command, but it is not a crash-atomic database transaction.  Keep
+an external backup when modifying valuable data, particularly when storage or
+system reliability is uncertain.
+
+
+[[examples]]
+== EXAMPLES
+
+.Analyze an assembly without changing it.
+[example]
+====
+mged>**unpush -D -v -m similarity vehicle**
+====
+
+.Restore sharing only within two selected assemblies.
+[example]
+====
+mged>**unpush -L left_wheel right_wheel**
+====
+
+
+[[see_also]]
+== SEE ALSO
+
+*npush*(nged), *push*(nged), *xpush*(nged), *facetize*(nged)
+
+[[author]]
+== AUTHOR
+
+BRL-CAD Team
+
+
+[[bug_reports]]
+== BUG REPORTS
+
+Reports of bugs or problems should be submitted via electronic mail to devs@brlcad.org

doc/asciidoc/toplevel/mann_index.adoc

@@ -444,6 +444,7 @@ Interactive MGED and GED commands. See also link:man1_index.html[Command-Line To
 == U
 
 * link:mann/unhide.html[unhide]
+* link:mann/unpush.html[unpush]
 * link:mann/units.html[units]
 
 == V

include/bg/pca.h

@@ -21,7 +21,7 @@
 /*----------------------------------------------------------------------*/
 /** @addtogroup bg_pca
  *
- * Principle Component Analysis
+ * Principal Component Analysis
  *
  * Calculates an XYZ coordinate system such that it aligns with the largest
  * variations in the supplied data.  Intuitively, it "aligns" with the
@@ -61,6 +61,9 @@
 #define BG_PCA_H
 
 #include "common.h"
+
+#include <stddef.h>
+
 #include "vmath.h"
 #include "bg/defines.h"
 
@@ -69,7 +72,7 @@ __BEGIN_DECLS
 
 /**
  * @brief
- * Perform a Principle Component Analysis on a set of points.
+ * Perform a Principal Component Analysis on a set of points.
  *
  * Outputs are a center point and XYZ vectors for the coordinate system.
  *
@@ -86,9 +89,19 @@ __BEGIN_DECLS
  */
 BG_EXPORT extern int bg_pca(point_t *c, vect_t *xaxis, vect_t *yaxis, vect_t *zaxis, size_t npnts, const point_t *pnts);
 
+/**
+ * Perform PCA and also report the singular values of the centered point
+ * matrix, ordered to correspond to xaxis, yaxis, and zaxis.  A zero singular
+ * value identifies a dimension in which the input point set has no extent.
+ * The values scale with the square root of the number of input points.
+ */
+BG_EXPORT extern int bg_pca_svd(point_t *c, vect_t *xaxis, vect_t *yaxis,
+				vect_t *zaxis, vect_t *singular_values,
+				size_t npnts, const point_t *pnts);
+
 __END_DECLS
 
-#endif  /* BG_PLANE_H */
+#endif  /* BG_PCA_H */
 /** @} */
 /*
  * Local Variables:

include/rt/func.h

@@ -178,6 +178,20 @@ RT_EXPORT extern int rt_obj_make(const char *label, const point_t origin, double
  */
 RT_EXPORT extern int rt_obj_xform(struct rt_db_internal *op, const mat_t mat, struct rt_db_internal *ip, int release, struct db_i *dbip);
 
+/**
+ * Factor an object's placement into a canonical primitive and a matrix mapping
+ * that primitive back to the input geometry.
+ *
+ * The caller must pass an initialized, empty @p canonical and owns its contents
+ * on success.  Applying @p canonical_to_input to @p canonical with ft_mat must
+ * reproduce geometry equivalent to @p input within @p tol.
+ */
+RT_EXPORT extern int rt_obj_canonicalize(struct rt_db_internal *canonical,
+					 mat_t canonical_to_input,
+					 const struct rt_db_internal *input,
+					 const struct bn_tol *tol,
+					 enum rt_canonicalize_mode mode);
+
 /**
  * obtain parameters for an object in libpc form
  */

include/rt/functab.h

@@ -62,6 +62,27 @@ __BEGIN_DECLS
 
 struct rt_piecestate; /* forward declaration for ft_piece_shot / ft_piece_hitsegs */
 
+/**
+ * Transformation classes a primitive canonicalizer may factor out of an
+ * object's stored parameters.  Each successive mode permits all operations
+ * allowed by the preceding mode.
+ */
+enum rt_canonicalize_mode {
+    RT_CANONICALIZE_RIGID = 0,
+    RT_CANONICALIZE_SIMILARITY,
+    RT_CANONICALIZE_AFFINE
+};
+
+/**
+ * Return values for primitive canonicalization.  Unsupported is distinct from
+ * an error so callers may safely skip primitives without an implementation.
+ */
+enum rt_canonicalize_status {
+    RT_CANONICALIZE_OK = 0,
+    RT_CANONICALIZE_ERROR,
+    RT_CANONICALIZE_UNSUPPORTED
+};
+
 /**
  * This needs to be at the end of the raytrace.h header file, so that
  * all the structure names are known.  The "union record" and "struct
@@ -332,6 +353,22 @@ struct rt_functab {
     int (*ft_validate)(struct bu_vls *error_msg, const struct rt_db_internal *ip, const struct bn_tol *tol);
 #define RTFUNCTAB_FUNC_VALIDATE_CAST(_func) ((int (*)(struct bu_vls *, const struct rt_db_internal *, const struct bn_tol *))((void (*)(void))_func))
 
+    /**
+     * Factor primitive geometry into a deterministic local representation and
+     * a matrix mapping that representation back to the input geometry.
+     *
+     * On success, applying @p canonical_to_input to @p canonical with ft_mat
+     * must reproduce geometry equivalent to @p input within @p tol.  The
+     * callback initializes and owns the contents of @p canonical.  It must not
+     * modify @p input or access database state.
+     */
+    int (*ft_canonicalize)(struct rt_db_internal *canonical,
+			   mat_t canonical_to_input,
+			   const struct rt_db_internal *input,
+			   const struct bn_tol *tol,
+			   enum rt_canonicalize_mode mode);
+#define RTFUNCTAB_FUNC_CANONICALIZE_CAST(_func) ((int (*)(struct rt_db_internal *, mat_t, const struct rt_db_internal *, const struct bn_tol *, enum rt_canonicalize_mode))((void (*)(void))_func))
+
 };
 
 /**

regress/ged/CMakeLists.txt

@@ -19,6 +19,13 @@ if(TARGET libged)
   distclean("${CMAKE_CURRENT_BINARY_DIR}/regress-ged_push.log")
   distclean("${CMAKE_CURRENT_BINARY_DIR}/ged_push.g")
 
+  brlcad_addexec(regress_unpush unpush.cpp "libged;libwdb;librt;libbg;libbu" TEST_USESDATA)
+  if(TARGET regress_unpush)
+    add_dependencies(regress_unpush ged_plugins)
+  endif(TARGET regress_unpush)
+  brlcad_add_test(NAME regress-ged_unpush COMMAND regress_unpush)
+  brlcad_regression_test(regress-ged_unpush "regress_unpush" TEST_DEFINED)
+
   # Phase D regression: coplanar window-frame subtraction pattern.
   # Exercises the instance-aware variant planning pass in the Manifold
   # facetize pipeline (plan.cpp + _booltree_leaf_tess variant resolution).
@@ -44,6 +51,7 @@ cmakefiles(
   mater.c
   regdef.c
   push.cpp
+  unpush.cpp
   regress_facetize.cpp
   regress_lint_raytrace.cpp
   d2push_tests.g

src/libbg/pca.cpp

@@ -52,8 +52,9 @@
 #include "bg/pca.h"
 
 // Use SVD algorithm from Soderkvist to fit a plane to vertex points
-extern "C" int
-bg_pca(point_t *c, vect_t *xa, vect_t *ya, vect_t *za, size_t npnts, const point_t *pnts)
+static int
+bg_pca_impl(point_t *c, vect_t *xa, vect_t *ya, vect_t *za,
+	    vect_t *singular_values, size_t npnts, const point_t *pnts)
 {
     if (!c || !xa || !ya || !za || npnts == 0 || !pnts)
 	return BRLCAD_ERROR;
@@ -74,7 +75,9 @@ bg_pca(point_t *c, vect_t *xa, vect_t *ya, vect_t *za, size_t npnts, const point
     }
 
     // 3.  Perform SVD
-    Eigen::JacobiSVD<Eigen::MatrixXd> svd(A, Eigen::ComputeThinU);
+    // Full U supplies a complete 3D frame for one- and two-point inputs as
+    // well as ordinary 3D point sets.
+    Eigen::JacobiSVD<Eigen::MatrixXd> svd(A, Eigen::ComputeFullU);
 
     // 4.  Extract the vectors from the U matrix
     vect_t xaxis, yaxis, zaxis;
@@ -93,10 +96,33 @@ bg_pca(point_t *c, vect_t *xa, vect_t *ya, vect_t *za, size_t npnts, const point
     VMOVE(*xa, xaxis);
     VMOVE(*ya, yaxis);
     VMOVE(*za, zaxis);
+    if (singular_values) {
+	VSETALL(*singular_values, 0.0);
+	for (Eigen::Index i = 0; i < svd.singularValues().size() && i < 3; i++)
+	    (*singular_values)[i] = svd.singularValues()(i);
+    }
 
     return BRLCAD_OK;
 }
 
+
+extern "C" int
+bg_pca_svd(point_t *c, vect_t *xa, vect_t *ya, vect_t *za,
+	   vect_t *singular_values, size_t npnts, const point_t *pnts)
+{
+    if (!singular_values)
+	return BRLCAD_ERROR;
+
+    return bg_pca_impl(c, xa, ya, za, singular_values, npnts, pnts);
+}
+
+
+extern "C" int
+bg_pca(point_t *c, vect_t *xa, vect_t *ya, vect_t *za, size_t npnts, const point_t *pnts)
+{
+    return bg_pca_impl(c, xa, ya, za, NULL, npnts, pnts);
+}
+
 /** @} */
 // Local Variables:
 // tab-width: 8

src/libbg/tests/CMakeLists.txt

@@ -27,6 +27,9 @@ brlcad_addexec(bg_test "${CMAKE_CURRENT_BINARY_DIR}/bg_test.c" "${BG_TEST_LIBS}"
 # TODO - convert this into the scheme used by other tests below
 brlcad_addexec(bgtester bgtester.c "${BG_TEST_LIBS}" TEST)
 
+brlcad_addexec(bg_pca pca.c "${BG_TEST_LIBS}" TEST)
+brlcad_add_test(NAME bg_pca COMMAND bg_pca)
+
 #  ************ aabb_ray.c tests *************
 
 brlcad_addexec(bg_aabb_ray aabb_ray.c "${BG_TEST_LIBS}" TEST)
@@ -492,6 +495,7 @@ brlcad_add_test(NAME bg_geogram_optimizer COMMAND bg_geogram_optimizer)
 cmakefiles(
   bg_test.c.in
   geogram_optimizer.cpp
+  pca.c
   plane_dist.c
   plane_isect.c
   plane_pt.c

src/libbg/tests/pca.c

@@ -0,0 +1,145 @@
+/*                         P C A . C
+ * BRL-CAD
+ *
+ * Copyright (c) 2026 United States Government as represented by
+ * the U.S. Army Research Laboratory.
+ *
+ * This library is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * version 2.1 as published by the Free Software Foundation.
+ *
+ * This library is distributed in the hope that it will be useful, but
+ * WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with this file; see the file named COPYING for more
+ * information.
+ */
+
+#include "common.h"
+
+#include <math.h>
+
+#include "bu/app.h"
+#include "bu/log.h"
+#include "bg/pca.h"
+
+
+#define PCA_TEST_EPSILON 1.0e-9
+
+
+static int
+near_value(fastf_t a, fastf_t b)
+{
+    fastf_t scale = fmax(1.0, fmax(fabs(a), fabs(b)));
+
+    return fabs(a - b) <= PCA_TEST_EPSILON * scale;
+}
+
+
+static int
+test_box(void)
+{
+    point_t points[8];
+    point_t center;
+    vect_t xaxis;
+    vect_t yaxis;
+    vect_t zaxis;
+    vect_t singular_values;
+    size_t index = 0;
+
+    for (int xsign = -1; xsign <= 1; xsign += 2) {
+	for (int ysign = -1; ysign <= 1; ysign += 2) {
+	    for (int zsign = -1; zsign <= 1; zsign += 2) {
+		VSET(points[index], 4.0 + 3.0 * xsign,
+		     -2.0 + 2.0 * ysign, 7.0 + zsign);
+		index++;
+	    }
+	}
+    }
+
+    if (bg_pca_svd(&center, &xaxis, &yaxis, &zaxis, &singular_values,
+	8, (const point_t *)points) != BRLCAD_OK) {
+	bu_log("box PCA failed\n");
+	return 1;
+    }
+
+    if (!near_value(center[X], 4.0) || !near_value(center[Y], -2.0) ||
+	!near_value(center[Z], 7.0) ||
+	!near_value(singular_values[X], 3.0 * sqrt(8.0)) ||
+	!near_value(singular_values[Y], 2.0 * sqrt(8.0)) ||
+	!near_value(singular_values[Z], sqrt(8.0)) ||
+	!near_value(MAGNITUDE(xaxis), 1.0) ||
+	!near_value(MAGNITUDE(yaxis), 1.0) ||
+	!near_value(MAGNITUDE(zaxis), 1.0) ||
+	!near_value(VDOT(xaxis, yaxis), 0.0) ||
+	!near_value(VDOT(xaxis, zaxis), 0.0) ||
+	!near_value(VDOT(yaxis, zaxis), 0.0)) {
+	bu_log("box PCA returned incorrect center, frame, or singular values\n");
+	return 1;
+    }
+
+    return 0;
+}
+
+
+static int
+test_single_point(void)
+{
+    const point_t point = {3.0, -5.0, 8.0};
+    point_t center;
+    vect_t xaxis;
+    vect_t yaxis;
+    vect_t zaxis;
+    vect_t singular_values;
+
+    if (bg_pca_svd(&center, &xaxis, &yaxis, &zaxis, &singular_values,
+	1, (const point_t *)&point) != BRLCAD_OK) {
+	bu_log("single-point PCA failed\n");
+	return 1;
+    }
+
+    if (!near_value(center[X], point[X]) || !near_value(center[Y], point[Y]) ||
+	!near_value(center[Z], point[Z]) || !VNEAR_ZERO(singular_values, PCA_TEST_EPSILON) ||
+	!near_value(MAGNITUDE(xaxis), 1.0) ||
+	!near_value(MAGNITUDE(yaxis), 1.0) ||
+	!near_value(MAGNITUDE(zaxis), 1.0)) {
+	bu_log("single-point PCA returned invalid output\n");
+	return 1;
+    }
+
+    return 0;
+}
+
+
+int
+main(int UNUSED(argc), const char *argv[])
+{
+    point_t center;
+    vect_t xaxis;
+    vect_t yaxis;
+    vect_t zaxis;
+
+    bu_setprogname(argv[0]);
+
+    if (test_box() || test_single_point())
+	return 1;
+    if (bg_pca_svd(&center, &xaxis, &yaxis, &zaxis, NULL, 0, NULL) != BRLCAD_ERROR) {
+	bu_log("invalid PCA input was not rejected\n");
+	return 1;
+    }
+
+    return 0;
+}
+
+/*
+ * Local Variables:
+ * mode: C
+ * tab-width: 8
+ * indent-tabs-mode: t
+ * c-file-style: "stroustrup"
+ * End:
+ * ex: shiftwidth=4 tabstop=8
+ */

src/libgcv/plugins/iges/tests/test_iges_orientation.cpp

@@ -188,8 +188,7 @@ test_uncertain()
     require(result.shells.empty() || (result.shells.front().orientation == ShellOrientation::Indeterminate &&
 	!result.shells.front().corrected), "invalid parameter-loop senses caused an incorrect face repair");
     auto flat = box();
-    ON_Xform flatten;
-    flatten.DiagonalTransformation(1.0);
+    ON_Xform flatten = ON_Xform::DiagonalTransformation(1.0);
     flatten[2][2] = 1.0e-13;
     require(flat->Transform(flatten), "could not construct thin shell");
     result = repair_brep_orientation(*flat);

src/libged/CMakeLists.txt

@@ -328,6 +328,7 @@ set(
   LIBGED_SOURCES
   ${GED_CMD_SRC}
   ${GED_STATIC_REGISTRATION_CPP}
+  canonicalize.cpp
   edit_buf.cpp
   exec.cpp
   ged_init.cpp
@@ -355,6 +356,8 @@ set(
   wireframe_eval.c
 )
 
+cmakefiles(canonicalize_private.h)
+
 set_property(
   SOURCE ged_init.cpp
   APPEND

src/libged/canonicalize_private.h

@@ -0,0 +1,109 @@
+/*                 C A N O N I C A L I Z E _ P R I V A T E . H
+ * BRL-CAD
+ *
+ * Copyright (c) 2026 United States Government as represented by
+ * the U.S. Army Research Laboratory.
+ *
+ * This library is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * version 2.1 as published by the Free Software Foundation.
+ */
+
+#ifndef LIBGED_CANONICALIZE_PRIVATE_H
+#define LIBGED_CANONICALIZE_PRIVATE_H
+
+#include "common.h"
+
+#include "bn/tol.h"
+#include "ged/defines.h"
+#include "raytrace.h"
+
+__BEGIN_DECLS
+
+/* These helpers are private libged interfaces shared by commands that need to
+ * recognize or reconstruct canonical primitive geometry. */
+GED_EXPORT fastf_t _ged_matrix_roundoff_tolerance(fastf_t a, fastf_t b,
+	fastf_t database_limit);
+
+GED_EXPORT int _ged_matrices_numerically_equal(const mat_t a, const mat_t b,
+	const struct bn_tol *tol);
+
+GED_EXPORT int _ged_canonical_geometry_equal(
+	const struct rt_db_internal *a, const struct rt_db_internal *b,
+	const struct bn_tol *tol);
+
+GED_EXPORT int _ged_canonical_geometry_metric(
+	const struct rt_db_internal *object, const struct bn_tol *tol,
+	fastf_t *metric, fastf_t *metric_tolerance);
+
+GED_EXPORT int _ged_primitive_geometry_equal(
+	const struct rt_db_internal *a, const struct rt_db_internal *b,
+	const struct bn_tol *tol);
+
+/* Allocate and populate output with input transformed by matrix.  Output must
+ * be empty and remains empty on error. */
+GED_EXPORT int _ged_transform_primitive(struct rt_db_internal *output,
+	const mat_t matrix, const struct rt_db_internal *input);
+
+GED_EXPORT int _ged_transformed_geometry_equal(
+	const struct rt_db_internal *a, const mat_t a_matrix,
+	const struct rt_db_internal *b, const mat_t b_matrix,
+	const struct bn_tol *tol);
+
+__END_DECLS
+
+#ifdef __cplusplus
+
+#include <map>
+#include <string>
+#include <vector>
+
+struct GedCanonicalChildInfo {
+    size_t identity = 0;
+    mat_t placement = MAT_INIT_IDN;
+};
+
+struct GedCanonicalCombinationLeaf {
+    size_t child_identity = 0;
+    mat_t effective_matrix = MAT_INIT_IDN;
+    mat_t relative_matrix = MAT_INIT_IDN;
+};
+
+struct GedCanonicalCombinationRecord {
+    struct directory *dp = nullptr;
+    mat_t placement = MAT_INIT_IDN;
+    mat_t identity_placement = MAT_INIT_IDN;
+    mat_t representative_to_input = MAT_INIT_IDN;
+    std::string bucket;
+    std::vector<GedCanonicalCombinationLeaf> leaves;
+    fastf_t metric = 0.0;
+    fastf_t metric_tolerance = 0.0;
+    size_t identity = 0;
+    bool region = false;
+    bool valid = false;
+};
+
+struct GedCanonicalCombinationGroup {
+    std::vector<size_t> records;
+};
+
+struct GedCanonicalCombinationAnalysis {
+    std::vector<GedCanonicalCombinationRecord> records;
+    std::vector<GedCanonicalCombinationGroup> groups;
+    size_t failures = 0;
+};
+
+/* Identify equivalent combination DAG nodes bottom-up.  Primitive children
+ * supply identities and canonical-to-input placements established by the
+ * caller's selected transform mode.  Database metadata participates in
+ * identity when include_metadata is true; geometry-only consumers such as
+ * facetize deliberately ignore it. */
+GED_EXPORT int _ged_canonical_combination_analysis(
+    GedCanonicalCombinationAnalysis *analysis, struct db_i *dbip,
+    const std::vector<struct directory *> &combinations,
+    const std::map<struct directory *, GedCanonicalChildInfo> &primitive_children,
+    const struct bn_tol *tol, bool include_metadata);
+
+#endif
+
+#endif /* LIBGED_CANONICALIZE_PRIVATE_H */

src/libged/npush/CMakeLists.txt

@@ -1,12 +1,14 @@
-set(PLUGIN_LIBS libged libbu)
+set(PLUGIN_LIBS libged librt libbg libbu)
 
 ged_plugin_library(ged-npush
   npush.cpp
+  unpush.cpp
   LIBS ${PLUGIN_LIBS}
 )
 
 cmakefiles(
   CMakeLists.txt
+  unpush.cpp
 )
 
 # Local Variables:

src/libged/npush/npush.cpp

@@ -1708,11 +1708,14 @@ ged_npush_core(struct ged *gedp, int argc, const char *argv[])
 
 #include "../include/plugin.h"
 
+extern "C" int ged_unpush_core(struct ged *gedp, int argc, const char *argv[]);
+
 #define GED_NPUSH_COMMANDS(X, XID) \
     X(npush, ged_npush_core, GED_CMD_DEFAULT) \
+    X(unpush, ged_unpush_core, GED_CMD_DEFAULT) \
 
 GED_DECLARE_COMMAND_SET(GED_NPUSH_COMMANDS)
-GED_DECLARE_PLUGIN_MANIFEST("libged_npush", 1, GED_NPUSH_COMMANDS)
+GED_DECLARE_PLUGIN_MANIFEST("libged_npush", 2, GED_NPUSH_COMMANDS)
 
 // Local Variables:
 // tab-width: 8

src/libged/tests/CMakeLists.txt

@@ -4,6 +4,9 @@ if(TARGET ged_test_tops)
 endif(TARGET ged_test_tops)
 brlcad_add_test(NAME ged_test_tops_moss COMMAND ged_test_tops "${CMAKE_CURRENT_SOURCE_DIR}/ged_tests.g")
 
+brlcad_addexec(ged_test_canonicalize canonicalize.cpp libged TEST)
+brlcad_add_test(NAME ged_canonicalize COMMAND ged_test_canonicalize)
+
 brlcad_addexec(ged_test_rm ged_test_rm.c "libged;libwdb;librt;libbu" TEST)
 if(TARGET ged_test_rm)
   add_dependencies(ged_test_rm ged_plugins)
@@ -122,6 +125,7 @@ add_subdirectory(edit)
 # in a separate file, since its bad input examples stand an excellent chance of
 # breaking commands.
 cmakefiles(
+  canonicalize.cpp
   lint_tests.g
   rt_search_tests.g
   search_tests.g

src/libged/tests/canonicalize.cpp

@@ -0,0 +1,127 @@
+/*                       C A N O N I C A L I Z E . C P P
+ * BRL-CAD
+ *
+ * Copyright (c) 2026 United States Government as represented by
+ * the U.S. Army Research Laboratory.
+ *
+ * This library is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * version 2.1 as published by the Free Software Foundation.
+ */
+
+#include "common.h"
+
+#include <limits>
+
+#include "bu/app.h"
+#include "rt/func.h"
+
+#include "../canonicalize_private.h"
+
+
+namespace {
+
+void
+make_ell(struct rt_db_internal *intern, const point_t center,
+	 const vect_t a, const vect_t b, const vect_t c)
+{
+    RT_DB_INTERNAL_INIT(intern);
+    intern->idb_major_type = DB5_MAJORTYPE_BRLCAD;
+    intern->idb_minor_type = ID_ELL;
+    intern->idb_meth = &OBJ[ID_ELL];
+    BU_ALLOC(intern->idb_ptr, struct rt_ell_internal);
+    auto *ell = static_cast<struct rt_ell_internal *>(intern->idb_ptr);
+    ell->magic = RT_ELL_INTERNAL_MAGIC;
+    VMOVE(ell->v, center);
+    VMOVE(ell->a, a);
+    VMOVE(ell->b, b);
+    VMOVE(ell->c, c);
+}
+
+
+int
+test_matrix_comparison(const struct bn_tol *tol)
+{
+    mat_t a = MAT_INIT_IDN;
+    mat_t b = MAT_INIT_IDN;
+    const fastf_t epsilon = std::numeric_limits<fastf_t>::epsilon();
+
+    b[0] += 16.0 * epsilon;
+    if (!_ged_matrices_numerically_equal(a, b, tol))
+	return 1;
+
+    b[0] = a[0] + 0.5 * tol->perp;
+    if (_ged_matrices_numerically_equal(a, b, tol))
+	return 1;
+
+    b[0] = a[0];
+    b[3] = 0.5 * tol->dist;
+    return _ged_matrices_numerically_equal(a, b, tol) ? 1 : 0;
+}
+
+
+int
+test_geometry_comparison(const struct bn_tol *tol)
+{
+    const point_t center = VINIT_ZERO;
+    const vect_t x_axis = {2.0, 0.0, 0.0};
+    const vect_t y_axis = {0.0, 3.0, 0.0};
+    const vect_t z_axis = {0.0, 0.0, 4.0};
+    struct rt_db_internal a;
+    struct rt_db_internal equivalent;
+
+    make_ell(&a, center, x_axis, y_axis, z_axis);
+    make_ell(&equivalent, center, y_axis, x_axis, z_axis);
+    int failed = _ged_canonical_geometry_equal(&a, &equivalent, tol) ||
+	!_ged_primitive_geometry_equal(&a, &equivalent, tol);
+    rt_db_free_internal(&equivalent);
+    rt_db_free_internal(&a);
+    return failed;
+}
+
+
+int
+test_transform(const struct bn_tol *tol)
+{
+    const point_t center = VINIT_ZERO;
+    const vect_t x_axis = {2.0, 0.0, 0.0};
+    const vect_t y_axis = {0.0, 3.0, 0.0};
+    const vect_t z_axis = {0.0, 0.0, 4.0};
+    struct rt_db_internal input;
+    struct rt_db_internal transformed;
+    mat_t matrix = MAT_INIT_IDN;
+
+    make_ell(&input, center, x_axis, y_axis, z_axis);
+    matrix[0] = 0.0;
+    matrix[1] = -1.0;
+    matrix[4] = 1.0;
+    matrix[5] = 0.0;
+    MAT_DELTAS(matrix, 11.0, -7.0, 5.0);
+    RT_DB_INTERNAL_INIT(&transformed);
+    if (_ged_transform_primitive(&transformed, matrix, &input) != BRLCAD_OK) {
+	rt_db_free_internal(&input);
+	return 1;
+    }
+
+    int failed = !_ged_transformed_geometry_equal(&input, matrix,
+	&transformed, bn_mat_identity, tol);
+    rt_db_free_internal(&transformed);
+    rt_db_free_internal(&input);
+    return failed;
+}
+
+} /* namespace */
+
+
+int
+main(int UNUSED(argc), const char *argv[])
+{
+    const struct bn_tol tol = BN_TOL_INIT_TOL;
+
+    bu_setprogname(argv[0]);
+    int failures = test_matrix_comparison(&tol) +
+	test_geometry_comparison(&tol) + test_transform(&tol);
+    if (failures)
+	bu_log("libged canonical geometry tests failed: %d\n", failures);
+    return failures ? 1 : 0;
+}

src/librt/CMakeLists.txt

@@ -161,6 +161,7 @@ set(
   primitives/ars/ars_mirror.c
   primitives/bot/bot.c
   primitives/bot/bot_brep.cpp
+  primitives/bot/bot_canonicalize.cpp
   primitives/bot/bot_edge.c
   primitives/bot/bot_mirror.c
   primitives/bot/bot_oriented_bbox.cpp
@@ -176,6 +177,7 @@ set(
   primitives/bspline/bspline.cpp
   primitives/bspline/bspline_brep.cpp
   primitives/bspline/bspline_mirror.c
+  primitives/canonicalize.c
   primitives/crofton.cpp
   primitives/cline/cline.c
   primitives/cline/cline_brep.cpp
@@ -228,6 +230,7 @@ set(
   primitives/nmg/nmg_mirror.c
   primitives/nmg/nmg_tri_mc.c
   primitives/obj_adjust.c
+  primitives/obj_canonicalize.c
   primitives/obj_curve.c
   primitives/obj_describe.c
   primitives/obj_export.c
@@ -325,6 +328,7 @@ set(
   material/materialX.h
   cut_hlbvh.h
   prcomb.c
+  primitives/canonicalize_private.h
   primitives/bot/bot_edge.h
   primitives/bot/decimate_private.h
   primitives/bot/bot_wireframe.cpp
@@ -482,6 +486,7 @@ brlcad_addlib(librt "${LIBRT_SOURCES}" "${RT_INCLUDE_DIRS}" "${RT_LOCAL_INCLUDE_
     primitives/bot/edbot.c
     primitives/bot/bot.c
     primitives/bot/bot_brep.cpp
+    primitives/bot/bot_canonicalize.cpp
     primitives/bot/bot_edge.c
     primitives/bot/bot_mirror.c
     primitives/bot/bot_oriented_bbox.cpp

src/librt/primitives/arb8/arb8.c

@@ -1970,6 +1970,233 @@ rt_arb_mat(struct rt_db_internal *rop, const mat_t mat, const struct rt_db_inter
 }
 
 
+static int
+arb_affine_frame(mat_t canonical_to_input, mat_t input_to_canonical,
+		 const struct rt_arb_internal *arb, const struct bn_tol *tol)
+{
+    vect_t first, second, third;
+    vect_t cross;
+    fastf_t first_mag = 0.0;
+    fastf_t cross_mag = 0.0;
+    int first_index = -1;
+    int second_index = -1;
+    int third_index = -1;
+
+    for (int i = 1; i < 8; i++) {
+	VSUB2(first, arb->pt[i], arb->pt[0]);
+	first_mag = MAGNITUDE(first);
+	if (first_mag > tol->dist) {
+	    first_index = i;
+	    break;
+	}
+    }
+    if (first_index < 0)
+	return 0;
+
+    for (int i = first_index + 1; i < 8; i++) {
+	fastf_t second_mag;
+	VSUB2(second, arb->pt[i], arb->pt[0]);
+	second_mag = MAGNITUDE(second);
+	if (second_mag <= tol->dist)
+	    continue;
+	VCROSS(cross, first, second);
+	cross_mag = MAGNITUDE(cross);
+	if (cross_mag > first_mag * second_mag * tol->perp) {
+	    second_index = i;
+	    break;
+	}
+    }
+    if (second_index < 0)
+	return 0;
+
+    for (int i = second_index + 1; i < 8; i++) {
+	fastf_t third_mag;
+	VSUB2(third, arb->pt[i], arb->pt[0]);
+	third_mag = MAGNITUDE(third);
+	if (third_mag <= tol->dist)
+	    continue;
+	if (fabs(VDOT(third, cross)) > third_mag * cross_mag * tol->perp) {
+	    third_index = i;
+	    break;
+	}
+    }
+    if (third_index < 0)
+	return 0;
+
+    VSUB2(first, arb->pt[first_index], arb->pt[0]);
+    VSUB2(second, arb->pt[second_index], arb->pt[0]);
+    VSUB2(third, arb->pt[third_index], arb->pt[0]);
+    MAT_IDN(canonical_to_input);
+    canonical_to_input[0] = first[X];
+    canonical_to_input[4] = first[Y];
+    canonical_to_input[8] = first[Z];
+    canonical_to_input[1] = second[X];
+    canonical_to_input[5] = second[Y];
+    canonical_to_input[9] = second[Z];
+    canonical_to_input[2] = third[X];
+    canonical_to_input[6] = third[Y];
+    canonical_to_input[10] = third[Z];
+    MAT_DELTAS_VEC(canonical_to_input, arb->pt[0]);
+
+    return bn_mat_inverse(input_to_canonical, canonical_to_input);
+}
+
+
+static int
+arb_rigid_frame(mat_t canonical_to_input, mat_t input_to_canonical,
+		const struct rt_arb_internal *arb, const struct bn_tol *tol,
+		enum rt_canonicalize_mode mode)
+{
+    vect_t xaxis, yaxis, zaxis;
+    fastf_t longest_sq = 0.0;
+    fastf_t widest_sq = 0.0;
+    fastf_t highest = 0.0;
+    fastf_t uniform_scale = 1.0;
+    int x_index = -1;
+    int y_index = -1;
+    int z_index = -1;
+
+    for (int i = 1; i < 8; i++) {
+	vect_t delta;
+	fastf_t distance_sq;
+	VSUB2(delta, arb->pt[i], arb->pt[0]);
+	distance_sq = MAGSQ(delta);
+	if (distance_sq > longest_sq + tol->dist_sq) {
+	    longest_sq = distance_sq;
+	    x_index = i;
+	}
+    }
+    if (x_index < 0 || longest_sq <= tol->dist_sq)
+	return 0;
+    VSUB2(xaxis, arb->pt[x_index], arb->pt[0]);
+    VUNITIZE(xaxis);
+
+    for (int i = 1; i < 8; i++) {
+	vect_t delta, perpendicular;
+	fastf_t width_sq;
+	VSUB2(delta, arb->pt[i], arb->pt[0]);
+	VJOIN1(perpendicular, delta, -VDOT(delta, xaxis), xaxis);
+	width_sq = MAGSQ(perpendicular);
+	if (width_sq > widest_sq + tol->dist_sq) {
+	    widest_sq = width_sq;
+	    VMOVE(yaxis, perpendicular);
+	    y_index = i;
+	}
+    }
+    if (y_index < 0 || widest_sq <= tol->dist_sq)
+	return 0;
+    VUNITIZE(yaxis);
+    VCROSS(zaxis, xaxis, yaxis);
+
+    for (int i = 1; i < 8; i++) {
+	vect_t delta;
+	fastf_t height;
+	VSUB2(delta, arb->pt[i], arb->pt[0]);
+	height = fabs(VDOT(delta, zaxis));
+	if (height > highest + tol->dist) {
+	    highest = height;
+	    z_index = i;
+	}
+    }
+    if (z_index < 0 || highest <= tol->dist)
+	return 0;
+    {
+	vect_t delta;
+	VSUB2(delta, arb->pt[z_index], arb->pt[0]);
+	if (VDOT(delta, zaxis) < 0.0) {
+	    VREVERSE(yaxis, yaxis);
+	    VREVERSE(zaxis, zaxis);
+	}
+    }
+
+    if (mode == RT_CANONICALIZE_SIMILARITY) {
+	fastf_t diameter_sq = 0.0;
+	for (int i = 0; i < 7; i++) {
+	    for (int j = i + 1; j < 8; j++) {
+		vect_t delta;
+		VSUB2(delta, arb->pt[j], arb->pt[i]);
+		diameter_sq = fmax(diameter_sq, MAGSQ(delta));
+	    }
+	}
+	uniform_scale = sqrt(diameter_sq);
+	if (uniform_scale <= tol->dist)
+	    return 0;
+    }
+
+    MAT_IDN(canonical_to_input);
+    canonical_to_input[0] = xaxis[X];
+    canonical_to_input[4] = xaxis[Y];
+    canonical_to_input[8] = xaxis[Z];
+    canonical_to_input[1] = yaxis[X];
+    canonical_to_input[5] = yaxis[Y];
+    canonical_to_input[9] = yaxis[Z];
+    canonical_to_input[2] = zaxis[X];
+    canonical_to_input[6] = zaxis[Y];
+    canonical_to_input[10] = zaxis[Z];
+    MAT_DELTAS_VEC(canonical_to_input, arb->pt[0]);
+    canonical_to_input[15] = 1.0 / uniform_scale;
+    canonical_to_input[3] /= uniform_scale;
+    canonical_to_input[7] /= uniform_scale;
+    canonical_to_input[11] /= uniform_scale;
+
+    return bn_mat_inverse(input_to_canonical, canonical_to_input);
+}
+
+
+C_DECL int
+rt_arb_canonicalize(struct rt_db_internal *canonical,
+		    mat_t canonical_to_input,
+		    const struct rt_db_internal *input,
+		    const struct bn_tol *tol,
+		    enum rt_canonicalize_mode mode)
+{
+    const struct rt_arb_internal *aip;
+    struct rt_arb_internal *caip;
+    mat_t input_to_canonical;
+    int frame_ok;
+    int validation_issues = 0;
+
+    if (!canonical || !canonical_to_input || !input || !tol)
+	return RT_CANONICALIZE_ERROR;
+    if (mode < RT_CANONICALIZE_RIGID || mode > RT_CANONICALIZE_AFFINE)
+	return RT_CANONICALIZE_ERROR;
+    if (input->idb_type != ID_ARB8)
+	return RT_CANONICALIZE_ERROR;
+
+    aip = (const struct rt_arb_internal *)input->idb_ptr;
+    RT_ARB_CK_MAGIC(aip);
+    for (int i = 0; i < 8; i++) {
+	if (!isfinite(aip->pt[i][X]) || !isfinite(aip->pt[i][Y]) ||
+	    !isfinite(aip->pt[i][Z]))
+	    return RT_CANONICALIZE_ERROR;
+    }
+    (void)rt_arb_validate(NULL, aip, tol, &validation_issues);
+    /* Historic databases contain usable non-standard vertex encodings.  The
+     * index-stable factorization preserves those encodings exactly. */
+    if (validation_issues & ~RT_ARB_VALIDATE_NONSTANDARD)
+	return RT_CANONICALIZE_ERROR;
+
+    if (mode == RT_CANONICALIZE_AFFINE)
+	frame_ok = arb_affine_frame(canonical_to_input, input_to_canonical, aip, tol);
+    else
+	frame_ok = arb_rigid_frame(canonical_to_input, input_to_canonical,
+	    aip, tol, mode);
+    if (!frame_ok)
+	return RT_CANONICALIZE_ERROR;
+
+    canonical->idb_major_type = DB5_MAJORTYPE_BRLCAD;
+    canonical->idb_minor_type = ID_ARB8;
+    canonical->idb_meth = &OBJ[ID_ARB8];
+    BU_ALLOC(canonical->idb_ptr, struct rt_arb_internal);
+    caip = (struct rt_arb_internal *)canonical->idb_ptr;
+    caip->magic = RT_ARB_INTERNAL_MAGIC;
+    for (int i = 0; i < 8; i++)
+	MAT4X3PNT(caip->pt[i], input_to_canonical, aip->pt[i]);
+
+    return RT_CANONICALIZE_OK;
+}
+
+
 /**
  * Import an arb from the db5 format and convert to the internal
  * structure.  Code duplicated from rt_arb_import4() with db5 help from

src/librt/primitives/bot/bot_canonicalize.cpp

@@ -0,0 +1,283 @@
+/*           B O T _ C A N O N I C A L I Z E . C P P
+ * BRL-CAD
+ *
+ * Copyright (c) 2026 United States Government as represented by
+ * the U.S. Army Research Laboratory.
+ *
+ * This library is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * version 2.1 as published by the Free Software Foundation.
+ *
+ * This library is distributed in the hope that it will be useful, but
+ * WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with this file; see the file named COPYING for more
+ * information.
+ */
+
+#include "common.h"
+
+#include <cmath>
+
+#include "bg/pca.h"
+#include "raytrace.h"
+#include "rt/primitives/bot.h"
+
+
+namespace {
+
+constexpr size_t FRAME_AXIS_COUNT = 3;
+constexpr size_t FRAME_PERMUTATION_COUNT = 6;
+
+const size_t frame_permutations[FRAME_PERMUTATION_COUNT][FRAME_AXIS_COUNT] = {
+    {0, 1, 2},
+    {0, 2, 1},
+    {1, 0, 2},
+    {1, 2, 0},
+    {2, 0, 1},
+    {2, 1, 0}
+};
+
+const int frame_permutation_parity[FRAME_PERMUTATION_COUNT] = {
+    1, -1, -1, 1, 1, -1
+};
+
+
+bool
+bot_input_valid(const struct rt_bot_internal *bot)
+{
+    if (!bot || bot->num_vertices < 3 || !bot->vertices ||
+	(bot->num_faces && !bot->faces))
+	return false;
+
+    for (size_t i = 0; i < bot->num_vertices; i++) {
+	if (VINVALID(&bot->vertices[3 * i]))
+	    return false;
+    }
+    for (size_t i = 0; i < 3 * bot->num_faces; i++) {
+	if (bot->faces[i] < 0 || static_cast<size_t>(bot->faces[i]) >= bot->num_vertices)
+	    return false;
+    }
+
+    return true;
+}
+
+
+void
+proper_pca_frame(vect_t axes[FRAME_AXIS_COUNT])
+{
+    vect_t proper_z;
+
+    VCROSS(proper_z, axes[X], axes[Y]);
+    VUNITIZE(proper_z);
+    VMOVE(axes[Z], proper_z);
+}
+
+
+bool
+pca_frame_is_unique(const vect_t singular_values, size_t point_count,
+		    const struct bn_tol *tol)
+{
+    fastf_t singular_tolerance = tol->dist * std::sqrt(static_cast<fastf_t>(point_count));
+
+    return !NEAR_EQUAL(singular_values[X], singular_values[Y], singular_tolerance) &&
+	!NEAR_EQUAL(singular_values[Y], singular_values[Z], singular_tolerance);
+}
+
+
+bool
+anchor_frame(vect_t axes[FRAME_AXIS_COUNT], const point_t center,
+	     const struct rt_bot_internal *bot, const struct bn_tol *tol)
+{
+    vect_t offset;
+    vect_t residual;
+    fastf_t farthest_sq = -1.0;
+    fastf_t widest_sq = -1.0;
+    size_t farthest = 0;
+    size_t widest = 0;
+
+    for (size_t i = 0; i < bot->num_vertices; i++) {
+	VSUB2(offset, &bot->vertices[3 * i], center);
+	fastf_t distance_sq = MAGSQ(offset);
+	if (distance_sq > farthest_sq + tol->dist_sq) {
+	    farthest_sq = distance_sq;
+	    farthest = i;
+	}
+    }
+    if (farthest_sq <= tol->dist_sq)
+	return false;
+
+    VSUB2(axes[X], &bot->vertices[3 * farthest], center);
+    VUNITIZE(axes[X]);
+
+    for (size_t i = 0; i < bot->num_vertices; i++) {
+	VSUB2(offset, &bot->vertices[3 * i], center);
+	VJOIN1(residual, offset, -VDOT(offset, axes[X]), axes[X]);
+	fastf_t residual_sq = MAGSQ(residual);
+	if (residual_sq > widest_sq + tol->dist_sq) {
+	    widest_sq = residual_sq;
+	    widest = i;
+	}
+    }
+    if (widest_sq <= tol->dist_sq)
+	return false;
+
+    VSUB2(offset, &bot->vertices[3 * widest], center);
+    VJOIN1(axes[Y], offset, -VDOT(offset, axes[X]), axes[X]);
+    VUNITIZE(axes[Y]);
+    VCROSS(axes[Z], axes[X], axes[Y]);
+    VUNITIZE(axes[Z]);
+    return true;
+}
+
+
+void
+frame_matrix(mat_t matrix, const vect_t axes[FRAME_AXIS_COUNT],
+	     const size_t permutation[FRAME_AXIS_COUNT], const int signs[FRAME_AXIS_COUNT],
+	     const point_t center)
+{
+    mat_t rotate;
+    mat_t translate;
+
+    MAT_IDN(rotate);
+    for (size_t row = 0; row < FRAME_AXIS_COUNT; row++)
+	VSCALE(&rotate[4 * row], axes[permutation[row]], signs[row]);
+    MAT_IDN(translate);
+    MAT_DELTAS_VEC_NEG(translate, center);
+    bn_mat_mul(matrix, rotate, translate);
+}
+
+
+fastf_t
+comparison_coordinate(fastf_t coordinate, const struct bn_tol *tol)
+{
+    return NEAR_ZERO(coordinate, tol->dist) ? 0.0 : coordinate;
+}
+
+
+bool
+candidate_is_less(const mat_t candidate, const mat_t current,
+		  const struct rt_bot_internal *bot, const struct bn_tol *tol)
+{
+    point_t candidate_point;
+    point_t current_point;
+
+    for (size_t i = 0; i < bot->num_vertices; i++) {
+	MAT4X3PNT(candidate_point, candidate, &bot->vertices[3 * i]);
+	MAT4X3PNT(current_point, current, &bot->vertices[3 * i]);
+	for (size_t coordinate = 0; coordinate < ELEMENTS_PER_POINT; coordinate++) {
+	    fastf_t candidate_value = comparison_coordinate(candidate_point[coordinate], tol);
+	    fastf_t current_value = comparison_coordinate(current_point[coordinate], tol);
+	    if (candidate_value < current_value - tol->dist)
+		return true;
+	    if (candidate_value > current_value + tol->dist)
+		return false;
+	}
+    }
+
+    return false;
+}
+
+}
+
+
+extern "C" int
+rt_bot_canonicalize(struct rt_db_internal *canonical,
+		    mat_t canonical_to_input,
+		    const struct rt_db_internal *input,
+		    const struct bn_tol *tol,
+		    enum rt_canonicalize_mode mode)
+{
+    const struct rt_bot_internal *bot;
+    struct rt_bot_internal *canonical_bot;
+    struct rt_db_internal transformed;
+    point_t center;
+    vect_t axes[FRAME_AXIS_COUNT];
+    vect_t singular_values;
+    mat_t candidate;
+    mat_t input_to_canonical;
+    bool have_candidate = false;
+
+    if (!canonical || !canonical_to_input || !input || !tol)
+	return RT_CANONICALIZE_ERROR;
+    if (mode < RT_CANONICALIZE_RIGID || mode > RT_CANONICALIZE_AFFINE)
+	return RT_CANONICALIZE_ERROR;
+    if (input->idb_type != ID_BOT)
+	return RT_CANONICALIZE_ERROR;
+
+    bot = static_cast<const struct rt_bot_internal *>(input->idb_ptr);
+    RT_BOT_CK_MAGIC(bot);
+    if (!bot_input_valid(bot))
+	return RT_CANONICALIZE_ERROR;
+
+    /* rt_bot_mat currently reconstructs only the common per-vertex normal
+     * layout.  Decline other layouts until that transform contract is fixed. */
+    if ((bot->bot_flags & RT_BOT_HAS_SURFACE_NORMALS) &&
+	bot->num_normals != bot->num_vertices)
+	return RT_CANONICALIZE_UNSUPPORTED;
+
+    if (bg_pca_svd(&center, &axes[X], &axes[Y], &axes[Z], &singular_values,
+	bot->num_vertices, reinterpret_cast<const point_t *>(bot->vertices)) != BRLCAD_OK)
+	return RT_CANONICALIZE_ERROR;
+
+    if (pca_frame_is_unique(singular_values, bot->num_vertices, tol)) {
+	proper_pca_frame(axes);
+    } else if (!anchor_frame(axes, center, bot, tol)) {
+	return RT_CANONICALIZE_ERROR;
+    }
+
+    /* PCA axis signs are arbitrary, and equivalent parameterizations may
+     * permute axes.  Search all proper signed permutations and retain the
+     * lexicographically smallest vertex sequence.  This preserves BOT index
+     * topology while making pushed copies choose the same local frame. */
+    for (size_t permutation = 0; permutation < FRAME_PERMUTATION_COUNT; permutation++) {
+	for (int xsign = -1; xsign <= 1; xsign += 2) {
+	    for (int ysign = -1; ysign <= 1; ysign += 2) {
+		int signs[FRAME_AXIS_COUNT] = {
+		    xsign,
+		    ysign,
+		    frame_permutation_parity[permutation] * xsign * ysign
+		};
+		frame_matrix(candidate, axes, frame_permutations[permutation], signs, center);
+		if (!have_candidate || candidate_is_less(candidate, input_to_canonical, bot, tol)) {
+		    MAT_COPY(input_to_canonical, candidate);
+		    have_candidate = true;
+		}
+	    }
+	}
+    }
+
+    canonical_bot = rt_bot_dup(bot);
+    if (!canonical_bot)
+	return RT_CANONICALIZE_ERROR;
+
+    RT_DB_INTERNAL_INIT(&transformed);
+    transformed.idb_major_type = DB5_MAJORTYPE_BRLCAD;
+    transformed.idb_minor_type = ID_BOT;
+    transformed.idb_meth = &OBJ[ID_BOT];
+    transformed.idb_ptr = canonical_bot;
+    if (OBJ[ID_BOT].ft_mat(&transformed, input_to_canonical, input) != BRLCAD_OK) {
+	rt_bot_internal_free(canonical_bot);
+	return RT_CANONICALIZE_ERROR;
+    }
+
+    bn_mat_inv(canonical_to_input, input_to_canonical);
+    canonical->idb_major_type = DB5_MAJORTYPE_BRLCAD;
+    canonical->idb_minor_type = ID_BOT;
+    canonical->idb_meth = &OBJ[ID_BOT];
+    canonical->idb_ptr = canonical_bot;
+    return RT_CANONICALIZE_OK;
+}
+
+/*
+ * Local Variables:
+ * mode: C++
+ * tab-width: 8
+ * indent-tabs-mode: t
+ * c-file-style: "stroustrup"
+ * End:
+ * ex: shiftwidth=4 tabstop=8
+ */

src/librt/primitives/canonicalize.c

@@ -0,0 +1,93 @@
+/*                       C A N O N I C A L I Z E . C
+ * BRL-CAD
+ *
+ * Copyright (c) 2026 United States Government as represented by
+ * the U.S. Army Research Laboratory.
+ *
+ * This library is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * version 2.1 as published by the Free Software Foundation.
+ */
+
+#include "common.h"
+
+#include <math.h>
+
+#include "bn/mat.h"
+#include "vmath.h"
+
+#include "canonicalize_private.h"
+
+
+int
+rt_make_canonical_similarity_frame(mat_t canonical_to_input,
+	fastf_t *canonical_lengths, const point_t origin,
+	const vect_t input_x, const vect_t input_z, const fastf_t *lengths,
+	size_t length_count, const struct bn_tol *tol,
+	enum rt_canonicalize_mode mode)
+{
+    vect_t x_axis;
+    vect_t y_axis;
+    vect_t z_axis;
+    fastf_t x_length;
+    fastf_t z_length;
+    fastf_t projection;
+    fastf_t uniform_scale = mode == RT_CANONICALIZE_RIGID ? 1.0 : 0.0;
+    mat_t orientation;
+    mat_t scale;
+    mat_t oriented_scale;
+    mat_t translate;
+
+    if (!canonical_to_input || !canonical_lengths || !origin || !input_x ||
+	!input_z || !lengths || !length_count || !tol ||
+	mode < RT_CANONICALIZE_RIGID || mode > RT_CANONICALIZE_AFFINE)
+	return 0;
+    if (!isfinite(origin[X]) || !isfinite(origin[Y]) ||
+	!isfinite(origin[Z]))
+	return 0;
+
+    x_length = MAGNITUDE(input_x);
+    z_length = MAGNITUDE(input_z);
+    if (!isfinite(x_length) || !isfinite(z_length) ||
+	x_length <= tol->dist || z_length <= tol->dist)
+	return 0;
+    VSCALE(x_axis, input_x, 1.0 / x_length);
+    VSCALE(z_axis, input_z, 1.0 / z_length);
+    projection = VDOT(x_axis, z_axis);
+    if (!NEAR_ZERO(projection, tol->perp))
+	return 0;
+
+    /* Remove accepted input roundoff so write-mode matrices are proper
+     * rotations rather than barely-sheared approximations. */
+    VJOIN1(x_axis, x_axis, -projection, z_axis);
+    VUNITIZE(x_axis);
+    VCROSS(y_axis, z_axis, x_axis);
+    VUNITIZE(y_axis);
+
+    for (size_t i = 0; i < length_count; i++) {
+	if (!isfinite(lengths[i]) || lengths[i] <= tol->dist)
+	    return 0;
+	if (mode != RT_CANONICALIZE_RIGID)
+	    uniform_scale = fmax(uniform_scale, lengths[i]);
+    }
+    for (size_t i = 0; i < length_count; i++)
+	canonical_lengths[i] = lengths[i] / uniform_scale;
+
+    MAT_IDN(orientation);
+    orientation[0] = x_axis[X];
+    orientation[4] = x_axis[Y];
+    orientation[8] = x_axis[Z];
+    orientation[1] = y_axis[X];
+    orientation[5] = y_axis[Y];
+    orientation[9] = y_axis[Z];
+    orientation[2] = z_axis[X];
+    orientation[6] = z_axis[Y];
+    orientation[10] = z_axis[Z];
+    MAT_IDN(scale);
+    scale[15] = 1.0 / uniform_scale;
+    bn_mat_mul(oriented_scale, orientation, scale);
+    MAT_IDN(translate);
+    MAT_DELTAS_VEC(translate, origin);
+    bn_mat_mul(canonical_to_input, translate, oriented_scale);
+    return 1;
+}

src/librt/primitives/canonicalize_private.h

@@ -0,0 +1,33 @@
+/*                 C A N O N I C A L I Z E _ P R I V A T E . H
+ * BRL-CAD
+ *
+ * Copyright (c) 2026 United States Government as represented by
+ * the U.S. Army Research Laboratory.
+ *
+ * This library is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * version 2.1 as published by the Free Software Foundation.
+ */
+
+#ifndef LIBRT_PRIMITIVES_CANONICALIZE_PRIVATE_H
+#define LIBRT_PRIMITIVES_CANONICALIZE_PRIVATE_H
+
+#include "common.h"
+
+#include <stddef.h>
+
+#include "rt/functab.h"
+
+
+/* Construct a right-handed frame whose local X and Z axes follow two
+ * perpendicular primitive directions.  Scalar-radius primitives cannot
+ * represent a nonuniform affine transform in ft_mat, so affine requests use
+ * the strongest exactly reconstructable subset: a similarity transform. */
+int
+rt_make_canonical_similarity_frame(mat_t canonical_to_input,
+	fastf_t *canonical_lengths, const point_t origin,
+	const vect_t input_x, const vect_t input_z, const fastf_t *lengths,
+	size_t length_count, const struct bn_tol *tol,
+	enum rt_canonicalize_mode mode);
+
+#endif /* LIBRT_PRIMITIVES_CANONICALIZE_PRIVATE_H */

src/librt/primitives/ehy/ehy.c

@@ -162,6 +162,7 @@
 #include "raytrace.h"
 
 #include "../../librt_private.h"
+#include "../canonicalize_private.h"
 
 
 struct ehy_specific {
@@ -2118,6 +2119,48 @@ rt_ehy_mat(struct rt_db_internal *rop, const mat_t mat, const struct rt_db_inter
     return BRLCAD_OK;
 }
 
+
+C_DECL int
+rt_ehy_canonicalize(struct rt_db_internal *canonical,
+		    mat_t canonical_to_input,
+		    const struct rt_db_internal *input,
+		    const struct bn_tol *tol,
+		    enum rt_canonicalize_mode mode)
+{
+    const struct rt_ehy_internal *eip;
+    struct rt_ehy_internal *ceip;
+    fastf_t lengths[4];
+    fastf_t canonical_lengths[4];
+
+    if (!canonical || !canonical_to_input || !input || !tol ||
+	input->idb_type != ID_EHY)
+	return RT_CANONICALIZE_ERROR;
+    eip = (const struct rt_ehy_internal *)input->idb_ptr;
+    RT_EHY_CK_MAGIC(eip);
+    lengths[0] = MAGNITUDE(eip->ehy_H);
+    lengths[1] = eip->ehy_r1;
+    lengths[2] = eip->ehy_r2;
+    lengths[3] = eip->ehy_c;
+    if (!rt_make_canonical_similarity_frame(canonical_to_input,
+	    canonical_lengths, eip->ehy_V, eip->ehy_Au, eip->ehy_H,
+	    lengths, 4, tol, mode))
+	return RT_CANONICALIZE_ERROR;
+
+    canonical->idb_major_type = DB5_MAJORTYPE_BRLCAD;
+    canonical->idb_minor_type = ID_EHY;
+    canonical->idb_meth = &OBJ[ID_EHY];
+    BU_ALLOC(canonical->idb_ptr, struct rt_ehy_internal);
+    ceip = (struct rt_ehy_internal *)canonical->idb_ptr;
+    ceip->ehy_magic = RT_EHY_INTERNAL_MAGIC;
+    VSETALL(ceip->ehy_V, 0.0);
+    VSET(ceip->ehy_H, 0.0, 0.0, canonical_lengths[0]);
+    VSET(ceip->ehy_Au, 1.0, 0.0, 0.0);
+    ceip->ehy_r1 = canonical_lengths[1];
+    ceip->ehy_r2 = canonical_lengths[2];
+    ceip->ehy_c = canonical_lengths[3];
+    return RT_CANONICALIZE_OK;
+}
+
 int
 rt_ehy_functab_validate(struct bu_vls *error_msg, const struct rt_db_internal *ip, const struct bn_tol *tol)
 {

src/librt/primitives/ell/ell.c

@@ -1423,6 +1423,156 @@ rt_ell_mat(struct rt_db_internal *rop, const mat_t mat, const struct rt_db_inter
     return BRLCAD_OK;
 }
 
+
+struct ell_canonical_axis {
+    vect_t dir;
+    fastf_t length;
+    size_t input_index;
+};
+
+
+static void
+ell_canonical_axis_swap(struct ell_canonical_axis *a, struct ell_canonical_axis *b)
+{
+    struct ell_canonical_axis tmp;
+
+    VMOVE(tmp.dir, a->dir);
+    tmp.length = a->length;
+    tmp.input_index = a->input_index;
+
+    VMOVE(a->dir, b->dir);
+    a->length = b->length;
+    a->input_index = b->input_index;
+
+    VMOVE(b->dir, tmp.dir);
+    b->length = tmp.length;
+    b->input_index = tmp.input_index;
+}
+
+
+static void
+ell_canonical_axes_sort(struct ell_canonical_axis axes[3], const struct bn_tol *tol)
+{
+    for (size_t i = 0; i < 2; i++) {
+	for (size_t j = i + 1; j < 3; j++) {
+	    if (axes[j].length > axes[i].length + tol->dist ||
+		(NEAR_EQUAL(axes[j].length, axes[i].length, tol->dist) &&
+		 axes[j].input_index < axes[i].input_index))
+		ell_canonical_axis_swap(&axes[i], &axes[j]);
+	}
+    }
+}
+
+
+static void
+ell_canonical_orientation(mat_t orientation, const struct ell_canonical_axis axes[3])
+{
+    MAT_IDN(orientation);
+
+    orientation[0] = axes[0].dir[X];
+    orientation[4] = axes[0].dir[Y];
+    orientation[8] = axes[0].dir[Z];
+
+    orientation[1] = axes[1].dir[X];
+    orientation[5] = axes[1].dir[Y];
+    orientation[9] = axes[1].dir[Z];
+
+    orientation[2] = axes[2].dir[X];
+    orientation[6] = axes[2].dir[Y];
+    orientation[10] = axes[2].dir[Z];
+}
+
+
+C_DECL int
+rt_ell_canonicalize(struct rt_db_internal *canonical,
+		    mat_t canonical_to_input,
+		    const struct rt_db_internal *input,
+		    const struct bn_tol *tol,
+		    enum rt_canonicalize_mode mode)
+{
+    const struct rt_ell_internal *eip;
+    struct rt_ell_internal *ceip;
+    struct ell_canonical_axis axes[3];
+    const fastf_t *input_axes[3];
+    vect_t cross;
+    mat_t orientation;
+    mat_t scale_mat;
+    mat_t oriented_scale;
+    mat_t translate;
+    fastf_t canonical_lengths[3];
+    fastf_t uniform_scale = 1.0;
+
+    if (!canonical || !canonical_to_input || !input || !tol)
+	return RT_CANONICALIZE_ERROR;
+    if (mode < RT_CANONICALIZE_RIGID || mode > RT_CANONICALIZE_AFFINE)
+	return RT_CANONICALIZE_ERROR;
+    if (input->idb_type != ID_ELL && input->idb_type != ID_SPH)
+	return RT_CANONICALIZE_ERROR;
+
+    eip = (const struct rt_ell_internal *)input->idb_ptr;
+    RT_ELL_CK_MAGIC(eip);
+
+    input_axes[0] = eip->a;
+    input_axes[1] = eip->b;
+    input_axes[2] = eip->c;
+    for (size_t i = 0; i < 3; i++) {
+	axes[i].length = MAGNITUDE(input_axes[i]);
+	if (axes[i].length <= tol->dist)
+	    return RT_CANONICALIZE_ERROR;
+	VSCALE(axes[i].dir, input_axes[i], 1.0 / axes[i].length);
+	axes[i].input_index = i;
+    }
+
+    if (!NEAR_ZERO(VDOT(axes[0].dir, axes[1].dir), tol->perp) ||
+	!NEAR_ZERO(VDOT(axes[0].dir, axes[2].dir), tol->perp) ||
+	!NEAR_ZERO(VDOT(axes[1].dir, axes[2].dir), tol->perp))
+	return RT_CANONICALIZE_ERROR;
+
+    ell_canonical_axes_sort(axes, tol);
+
+    /* Axis signs do not alter ellipsoid geometry.  Select the sign of the
+     * third axis so the placement matrix is a proper rotation. */
+    VCROSS(cross, axes[0].dir, axes[1].dir);
+    if (VDOT(cross, axes[2].dir) < 0.0)
+	VREVERSE(axes[2].dir, axes[2].dir);
+
+    if (mode == RT_CANONICALIZE_AFFINE) {
+	VSETALL(canonical_lengths, 1.0);
+	for (size_t i = 0; i < 3; i++)
+	    VSCALE(axes[i].dir, axes[i].dir, axes[i].length);
+	ell_canonical_orientation(orientation, axes);
+	MAT_IDN(translate);
+	MAT_DELTAS_VEC(translate, eip->v);
+	bn_mat_mul(canonical_to_input, translate, orientation);
+    } else {
+	if (mode == RT_CANONICALIZE_SIMILARITY)
+	    uniform_scale = axes[0].length;
+	for (size_t i = 0; i < 3; i++)
+	    canonical_lengths[i] = axes[i].length / uniform_scale;
+
+	ell_canonical_orientation(orientation, axes);
+	MAT_IDN(scale_mat);
+	scale_mat[15] = 1.0 / uniform_scale;
+	bn_mat_mul(oriented_scale, orientation, scale_mat);
+	MAT_IDN(translate);
+	MAT_DELTAS_VEC(translate, eip->v);
+	bn_mat_mul(canonical_to_input, translate, oriented_scale);
+    }
+
+    canonical->idb_major_type = DB5_MAJORTYPE_BRLCAD;
+    canonical->idb_minor_type = input->idb_minor_type;
+    canonical->idb_meth = &OBJ[input->idb_minor_type];
+    BU_ALLOC(canonical->idb_ptr, struct rt_ell_internal);
+    ceip = (struct rt_ell_internal *)canonical->idb_ptr;
+    ceip->magic = RT_ELL_INTERNAL_MAGIC;
+    VSETALL(ceip->v, 0.0);
+    VSET(ceip->a, canonical_lengths[0], 0.0, 0.0);
+    VSET(ceip->b, 0.0, canonical_lengths[1], 0.0);
+    VSET(ceip->c, 0.0, 0.0, canonical_lengths[2]);
+
+    return RT_CANONICALIZE_OK;
+}
+
 /**
  * Import an ellipsoid/sphere from the database format to the internal
  * structure.  Apply modeling transformations as well.

src/librt/primitives/epa/epa.c

@@ -159,6 +159,7 @@
 #include "raytrace.h"
 
 #include "../../librt_private.h"
+#include "../canonicalize_private.h"
 
 static int epa_is_valid(struct rt_epa_internal *epa);
 
@@ -2044,6 +2045,46 @@ rt_epa_mat(struct rt_db_internal *rop, const mat_t mat, const struct rt_db_inter
     return BRLCAD_OK;
 }
 
+
+C_DECL int
+rt_epa_canonicalize(struct rt_db_internal *canonical,
+		    mat_t canonical_to_input,
+		    const struct rt_db_internal *input,
+		    const struct bn_tol *tol,
+		    enum rt_canonicalize_mode mode)
+{
+    const struct rt_epa_internal *eip;
+    struct rt_epa_internal *ceip;
+    fastf_t lengths[3];
+    fastf_t canonical_lengths[3];
+
+    if (!canonical || !canonical_to_input || !input || !tol ||
+	input->idb_type != ID_EPA)
+	return RT_CANONICALIZE_ERROR;
+    eip = (const struct rt_epa_internal *)input->idb_ptr;
+    RT_EPA_CK_MAGIC(eip);
+    lengths[0] = MAGNITUDE(eip->epa_H);
+    lengths[1] = eip->epa_r1;
+    lengths[2] = eip->epa_r2;
+    if (!rt_make_canonical_similarity_frame(canonical_to_input,
+	    canonical_lengths, eip->epa_V, eip->epa_Au, eip->epa_H,
+	    lengths, 3, tol, mode))
+	return RT_CANONICALIZE_ERROR;
+
+    canonical->idb_major_type = DB5_MAJORTYPE_BRLCAD;
+    canonical->idb_minor_type = ID_EPA;
+    canonical->idb_meth = &OBJ[ID_EPA];
+    BU_ALLOC(canonical->idb_ptr, struct rt_epa_internal);
+    ceip = (struct rt_epa_internal *)canonical->idb_ptr;
+    ceip->epa_magic = RT_EPA_INTERNAL_MAGIC;
+    VSETALL(ceip->epa_V, 0.0);
+    VSET(ceip->epa_H, 0.0, 0.0, canonical_lengths[0]);
+    VSET(ceip->epa_Au, 1.0, 0.0, 0.0);
+    ceip->epa_r1 = canonical_lengths[1];
+    ceip->epa_r2 = canonical_lengths[2];
+    return RT_CANONICALIZE_OK;
+}
+
 /**
  * Import an EPA from the database format to the internal format.
  * Apply modeling transformations as well.

src/librt/primitives/eto/eto.c

@@ -1931,6 +1931,117 @@ rt_eto_mat(struct rt_db_internal *rop, const mat_t mat, const struct rt_db_inter
 }
 
 
+C_DECL int
+rt_eto_canonicalize(struct rt_db_internal *canonical,
+		    mat_t canonical_to_input,
+		    const struct rt_db_internal *input,
+		    const struct bn_tol *tol,
+		    enum rt_canonicalize_mode mode)
+{
+    const struct rt_eto_internal *eip;
+    struct rt_eto_internal *ceip;
+    vect_t normal, horizontal, xaxis, yaxis;
+    mat_t orientation;
+    mat_t scale_mat;
+    mat_t oriented_scale;
+    mat_t translate;
+    fastf_t normal_mag;
+    fastf_t major_mag;
+    fastf_t vertical_mag;
+    fastf_t horizontal_mag;
+    fastf_t horizontal_sq;
+    fastf_t uniform_scale = 1.0;
+
+    if (!canonical || !canonical_to_input || !input || !tol)
+	return RT_CANONICALIZE_ERROR;
+    if (mode < RT_CANONICALIZE_RIGID || mode > RT_CANONICALIZE_AFFINE)
+	return RT_CANONICALIZE_ERROR;
+    if (input->idb_type != ID_ETO)
+	return RT_CANONICALIZE_ERROR;
+
+    eip = (const struct rt_eto_internal *)input->idb_ptr;
+    RT_ETO_CK_MAGIC(eip);
+    normal_mag = MAGNITUDE(eip->eto_N);
+    major_mag = MAGNITUDE(eip->eto_C);
+    if (!isfinite(eip->eto_V[X]) || !isfinite(eip->eto_V[Y]) ||
+	!isfinite(eip->eto_V[Z]) || !isfinite(normal_mag) ||
+	!isfinite(major_mag) || !isfinite(eip->eto_r) ||
+	!isfinite(eip->eto_rd) || normal_mag <= tol->dist ||
+	major_mag <= tol->dist || eip->eto_r <= tol->dist ||
+	eip->eto_rd <= tol->dist || eip->eto_rd > major_mag)
+	return RT_CANONICALIZE_ERROR;
+
+    VSCALE(normal, eip->eto_N, 1.0 / normal_mag);
+    vertical_mag = VDOT(eip->eto_C, normal);
+    if (vertical_mag < 0.0) {
+	VREVERSE(normal, normal);
+	vertical_mag = -vertical_mag;
+    }
+    VJOIN1(horizontal, eip->eto_C, -vertical_mag, normal);
+    horizontal_sq = MAGSQ(horizontal);
+    horizontal_mag = sqrt(fmax(0.0, horizontal_sq));
+
+    if (horizontal_mag > SMALL_FASTF) {
+	VSCALE(xaxis, horizontal, 1.0 / horizontal_mag);
+    } else {
+	/* A vertical C leaves rotation around N unconstrained.  Project the
+	 * least-aligned world axis so canonical input still produces an identity
+	 * placement while other orientations remain numerically stable. */
+	vect_t reference = VINIT_ZERO;
+	if (fabs(normal[X]) <= fabs(normal[Y]) && fabs(normal[X]) <= fabs(normal[Z]))
+	    reference[X] = 1.0;
+	else if (fabs(normal[Y]) <= fabs(normal[Z]))
+	    reference[Y] = 1.0;
+	else
+	    reference[Z] = 1.0;
+	VJOIN1(xaxis, reference, -VDOT(reference, normal), normal);
+	VUNITIZE(xaxis);
+	VSETALL(horizontal, 0.0);
+	horizontal_mag = 0.0;
+    }
+    VCROSS(yaxis, normal, xaxis);
+    VUNITIZE(yaxis);
+
+    /* rt_eto_mat can apply only rigid and uniform-scale placements to r and
+	 * rd.  AFFINE therefore intentionally uses the strongest transform that
+	 * this primitive can reconstruct exactly: a similarity. */
+    if (mode != RT_CANONICALIZE_RIGID)
+	uniform_scale = eip->eto_r;
+
+    MAT_IDN(orientation);
+    orientation[0] = xaxis[X];
+    orientation[4] = xaxis[Y];
+    orientation[8] = xaxis[Z];
+    orientation[1] = yaxis[X];
+    orientation[5] = yaxis[Y];
+    orientation[9] = yaxis[Z];
+    orientation[2] = normal[X];
+    orientation[6] = normal[Y];
+    orientation[10] = normal[Z];
+    MAT_IDN(scale_mat);
+    scale_mat[15] = 1.0 / uniform_scale;
+    bn_mat_mul(oriented_scale, orientation, scale_mat);
+    MAT_IDN(translate);
+    MAT_DELTAS_VEC(translate, eip->eto_V);
+    bn_mat_mul(canonical_to_input, translate, oriented_scale);
+
+    canonical->idb_major_type = DB5_MAJORTYPE_BRLCAD;
+    canonical->idb_minor_type = ID_ETO;
+    canonical->idb_meth = &OBJ[ID_ETO];
+    BU_ALLOC(canonical->idb_ptr, struct rt_eto_internal);
+    ceip = (struct rt_eto_internal *)canonical->idb_ptr;
+    ceip->eto_magic = RT_ETO_INTERNAL_MAGIC;
+    VSETALL(ceip->eto_V, 0.0);
+    VSET(ceip->eto_N, 0.0, 0.0, 1.0);
+    VSET(ceip->eto_C, horizontal_mag / uniform_scale, 0.0,
+	vertical_mag / uniform_scale);
+    ceip->eto_r = eip->eto_r / uniform_scale;
+    ceip->eto_rd = eip->eto_rd / uniform_scale;
+
+    return RT_CANONICALIZE_OK;
+}
+
+
 /**
  * Import a eto from the database format to the internal format.
  * Apply modeling transformations at the same time.

src/librt/primitives/half/half.c

@@ -653,6 +653,61 @@ rt_hlf_mat(struct rt_db_internal *rop, const mat_t mat, const struct rt_db_inter
     return 0;			/* OK */
 }
 
+
+C_DECL int
+rt_hlf_canonicalize(struct rt_db_internal *canonical,
+		    mat_t canonical_to_input,
+		    const struct rt_db_internal *input,
+		    const struct bn_tol *tol,
+		    enum rt_canonicalize_mode mode)
+{
+    struct rt_half_internal *hip;
+    struct rt_half_internal *chip;
+    vect_t normal;
+    vect_t zaxis = VINIT_ZERO;
+    point_t plane_point;
+    mat_t rotate;
+    mat_t translate;
+    fastf_t normal_mag;
+    fastf_t plane_dist;
+
+    if (!canonical || !canonical_to_input || !input || !tol)
+	return RT_CANONICALIZE_ERROR;
+    if (mode < RT_CANONICALIZE_RIGID || mode > RT_CANONICALIZE_AFFINE)
+	return RT_CANONICALIZE_ERROR;
+    if (input->idb_type != ID_HALF)
+	return RT_CANONICALIZE_ERROR;
+
+    hip = (struct rt_half_internal *)input->idb_ptr;
+    RT_HALF_CK_MAGIC(hip);
+
+    VMOVE(normal, hip->eqn);
+    normal_mag = MAGNITUDE(normal);
+    if (normal_mag <= tol->dist)
+	return RT_CANONICALIZE_ERROR;
+
+    VSCALE(normal, normal, 1.0 / normal_mag);
+    plane_dist = hip->eqn[W] / normal_mag;
+    VSCALE(plane_point, normal, plane_dist);
+
+    VSET(zaxis, 0.0, 0.0, 1.0);
+    bn_mat_fromto(rotate, zaxis, normal, tol);
+    MAT_IDN(translate);
+    MAT_DELTAS_VEC(translate, plane_point);
+    bn_mat_mul(canonical_to_input, translate, rotate);
+
+    canonical->idb_major_type = DB5_MAJORTYPE_BRLCAD;
+    canonical->idb_minor_type = ID_HALF;
+    canonical->idb_meth = &OBJ[ID_HALF];
+    BU_ALLOC(canonical->idb_ptr, struct rt_half_internal);
+    chip = (struct rt_half_internal *)canonical->idb_ptr;
+    chip->magic = RT_HALF_INTERNAL_MAGIC;
+    VSET(chip->eqn, 0.0, 0.0, 1.0);
+    chip->eqn[W] = 0.0;
+
+    return RT_CANONICALIZE_OK;
+}
+
 C_DECL int
 rt_hlf_import5(struct rt_db_internal *ip, const struct bu_external *ep, register const fastf_t *mat, const struct db_i *dbip)
 {

src/librt/primitives/obj_canonicalize.c

@@ -0,0 +1,68 @@
+/*                 O B J _ C A N O N I C A L I Z E . C
+ * BRL-CAD
+ *
+ * Copyright (c) 2026 United States Government as represented by
+ * the U.S. Army Research Laboratory.
+ *
+ * This library is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * version 2.1 as published by the Free Software Foundation.
+ *
+ * This library is distributed in the hope that it will be useful, but
+ * WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with this file; see the file named COPYING for more
+ * information.
+ */
+
+#include "common.h"
+
+#include "raytrace.h"
+#include "rt/func.h"
+
+
+int
+rt_obj_canonicalize(struct rt_db_internal *canonical,
+		    mat_t canonical_to_input,
+		    const struct rt_db_internal *input,
+		    const struct bn_tol *tol,
+		    enum rt_canonicalize_mode mode)
+{
+    const struct rt_functab *ft;
+
+    if (!canonical || !canonical_to_input || !input || !tol)
+	return RT_CANONICALIZE_ERROR;
+
+    RT_CK_DB_INTERNAL(canonical);
+    RT_CK_DB_INTERNAL(input);
+    BN_CK_TOL(tol);
+
+    if (canonical->idb_ptr || canonical->idb_meth ||
+	canonical->idb_major_type != -1 || canonical->idb_minor_type != -1)
+	return RT_CANONICALIZE_ERROR;
+
+    if (mode < RT_CANONICALIZE_RIGID || mode > RT_CANONICALIZE_AFFINE)
+	return RT_CANONICALIZE_ERROR;
+
+    if (input->idb_minor_type < 0 || input->idb_minor_type > ID_MAX_SOLID)
+	return RT_CANONICALIZE_ERROR;
+
+    ft = &OBJ[input->idb_minor_type];
+    if (!ft->ft_canonicalize)
+	return RT_CANONICALIZE_UNSUPPORTED;
+
+    return ft->ft_canonicalize(canonical, canonical_to_input, input, tol, mode);
+}
+
+/*
+ * Local Variables:
+ * mode: C
+ * tab-width: 8
+ * indent-tabs-mode: t
+ * c-file-style: "stroustrup"
+ * End:
+ * ex: shiftwidth=4 tabstop=8
+ */

src/librt/primitives/part/part.c

@@ -191,6 +191,7 @@
 #include "raytrace.h"
 #include "nmg.h"
 #include "../../librt_private.h"
+#include "../canonicalize_private.h"
 
 
 struct part_specific {
@@ -1978,6 +1979,110 @@ rt_part_mat(struct rt_db_internal *rop, const mat_t mat, const struct rt_db_inte
 
 }
 
+
+C_DECL int
+rt_part_canonicalize(struct rt_db_internal *canonical,
+		     mat_t canonical_to_input,
+		     const struct rt_db_internal *input,
+		     const struct bn_tol *tol,
+		     enum rt_canonicalize_mode mode)
+{
+    const struct rt_part_internal *pip;
+    struct rt_part_internal *cpip;
+    point_t origin;
+    vect_t axis;
+    vect_t x_axis;
+    fastf_t vertex_radius;
+    fastf_t height_radius;
+    fastf_t height_length;
+    fastf_t maximum_radius;
+    fastf_t uniform_scale;
+
+    if (!canonical || !canonical_to_input || !input || !tol ||
+	input->idb_type != ID_PARTICLE || mode < RT_CANONICALIZE_RIGID ||
+	mode > RT_CANONICALIZE_AFFINE)
+	return RT_CANONICALIZE_ERROR;
+    pip = (const struct rt_part_internal *)input->idb_ptr;
+    RT_PART_CK_MAGIC(pip);
+    if (!isfinite(pip->part_V[X]) || !isfinite(pip->part_V[Y]) ||
+	!isfinite(pip->part_V[Z]) || !isfinite(pip->part_H[X]) ||
+	!isfinite(pip->part_H[Y]) || !isfinite(pip->part_H[Z]) ||
+	!isfinite(pip->part_vrad) || !isfinite(pip->part_hrad) ||
+	pip->part_vrad < 0.0 || pip->part_hrad < 0.0)
+	return RT_CANONICALIZE_ERROR;
+
+    height_length = MAGNITUDE(pip->part_H);
+    maximum_radius = fmax(pip->part_vrad, pip->part_hrad);
+    if (!isfinite(height_length) || maximum_radius <= tol->dist)
+	return RT_CANONICALIZE_ERROR;
+    uniform_scale = mode == RT_CANONICALIZE_RIGID ? 1.0 :
+	fmax(height_length, maximum_radius);
+
+    canonical->idb_major_type = DB5_MAJORTYPE_BRLCAD;
+    canonical->idb_minor_type = ID_PARTICLE;
+    canonical->idb_meth = &OBJ[ID_PARTICLE];
+    BU_ALLOC(canonical->idb_ptr, struct rt_part_internal);
+    cpip = (struct rt_part_internal *)canonical->idb_ptr;
+    cpip->part_magic = RT_PART_INTERNAL_MAGIC;
+    cpip->part_type = pip->part_type;
+    VSETALL(cpip->part_V, 0.0);
+
+    /* Imported particle spheres have H cleared by rt_part_mat.  Do not turn a
+     * merely short cone into a different solid just because its height falls
+     * below the caller's modeling tolerance; the frame helper will reject
+     * nonzero axes that are too short to orient reliably. */
+    if (height_length <= SMALL_FASTF) {
+	mat_t scale;
+	mat_t translate;
+	MAT_IDN(scale);
+	scale[15] = 1.0 / uniform_scale;
+	MAT_IDN(translate);
+	MAT_DELTAS_VEC(translate, pip->part_V);
+	bn_mat_mul(canonical_to_input, translate, scale);
+	VSETALL(cpip->part_H, 0.0);
+	cpip->part_vrad = maximum_radius / uniform_scale;
+	cpip->part_hrad = cpip->part_vrad;
+	cpip->part_type = RT_PARTICLE_TYPE_SPHERE;
+	return RT_CANONICALIZE_OK;
+    }
+
+    if (pip->part_vrad <= pip->part_hrad) {
+	VMOVE(origin, pip->part_V);
+	VMOVE(axis, pip->part_H);
+	vertex_radius = pip->part_vrad;
+	height_radius = pip->part_hrad;
+    } else {
+	VADD2(origin, pip->part_V, pip->part_H);
+	VREVERSE(axis, pip->part_H);
+	vertex_radius = pip->part_hrad;
+	height_radius = pip->part_vrad;
+    }
+    /* Canonical input must retain an identity placement.  Other axial
+     * directions may use any deterministic perpendicular because PART is
+     * rotationally symmetric about H. */
+    if (NEAR_ZERO(axis[X], tol->perp) &&
+	NEAR_ZERO(axis[Y], tol->perp) && axis[Z] > 0.0)
+	VSET(x_axis, 1.0, 0.0, 0.0);
+    else
+	bn_vec_ortho(x_axis, axis);
+    {
+	const fastf_t frame_lengths[2] = {height_length, maximum_radius};
+	fastf_t canonical_frame_lengths[2];
+	if (!rt_make_canonical_similarity_frame(canonical_to_input,
+		canonical_frame_lengths, origin, x_axis, axis, frame_lengths, 2,
+		tol, mode)) {
+	    rt_db_free_internal(canonical);
+	    return RT_CANONICALIZE_ERROR;
+	}
+    }
+
+    VSET(cpip->part_H, 0.0, 0.0, height_length / uniform_scale);
+    cpip->part_vrad = vertex_radius / uniform_scale;
+    cpip->part_hrad = height_radius / uniform_scale;
+    return RT_CANONICALIZE_OK;
+}
+
+
 C_DECL int
 rt_part_import5(struct rt_db_internal *ip, const struct bu_external *ep, register const fastf_t *mat, const struct db_i *dbip)
 {

src/librt/primitives/rhc/rhc.c

@@ -174,6 +174,7 @@
 #include "raytrace.h"
 
 #include "../../librt_private.h"
+#include "../canonicalize_private.h"
 
 static int
 rhc_is_valid(struct rt_rhc_internal *rhc);
@@ -1932,6 +1933,47 @@ rt_rhc_mat(struct rt_db_internal *rop, const mat_t mat, const struct rt_db_inter
     return BRLCAD_OK;
 }
 
+
+C_DECL int
+rt_rhc_canonicalize(struct rt_db_internal *canonical,
+		    mat_t canonical_to_input,
+		    const struct rt_db_internal *input,
+		    const struct bn_tol *tol,
+		    enum rt_canonicalize_mode mode)
+{
+    const struct rt_rhc_internal *rip;
+    struct rt_rhc_internal *crip;
+    fastf_t lengths[4];
+    fastf_t canonical_lengths[4];
+
+    if (!canonical || !canonical_to_input || !input || !tol ||
+	input->idb_type != ID_RHC)
+	return RT_CANONICALIZE_ERROR;
+    rip = (const struct rt_rhc_internal *)input->idb_ptr;
+    RT_RHC_CK_MAGIC(rip);
+    lengths[0] = MAGNITUDE(rip->rhc_B);
+    lengths[1] = MAGNITUDE(rip->rhc_H);
+    lengths[2] = rip->rhc_r;
+    lengths[3] = rip->rhc_c;
+    if (!rt_make_canonical_similarity_frame(canonical_to_input,
+	    canonical_lengths, rip->rhc_V, rip->rhc_B, rip->rhc_H, lengths,
+	    4, tol, mode))
+	return RT_CANONICALIZE_ERROR;
+
+    canonical->idb_major_type = DB5_MAJORTYPE_BRLCAD;
+    canonical->idb_minor_type = ID_RHC;
+    canonical->idb_meth = &OBJ[ID_RHC];
+    BU_ALLOC(canonical->idb_ptr, struct rt_rhc_internal);
+    crip = (struct rt_rhc_internal *)canonical->idb_ptr;
+    crip->rhc_magic = RT_RHC_INTERNAL_MAGIC;
+    VSETALL(crip->rhc_V, 0.0);
+    VSET(crip->rhc_B, canonical_lengths[0], 0.0, 0.0);
+    VSET(crip->rhc_H, 0.0, 0.0, canonical_lengths[1]);
+    crip->rhc_r = canonical_lengths[2];
+    crip->rhc_c = canonical_lengths[3];
+    return RT_CANONICALIZE_OK;
+}
+
 /**
  * Import an RHC from the database format to the internal format.
  * Apply modeling transformations as well.

src/librt/primitives/rpc/rpc.c

@@ -172,6 +172,7 @@
 #include "raytrace.h"
 
 #include "../../librt_private.h"
+#include "../canonicalize_private.h"
 
 #if defined(HAVE_ASINH) && !defined(HAVE_DECL_ASINH)
 extern double asinh(double x);
@@ -1719,6 +1720,45 @@ rt_rpc_mat(struct rt_db_internal *rop, const mat_t mat, const struct rt_db_inter
     return BRLCAD_OK;
 }
 
+
+C_DECL int
+rt_rpc_canonicalize(struct rt_db_internal *canonical,
+		    mat_t canonical_to_input,
+		    const struct rt_db_internal *input,
+		    const struct bn_tol *tol,
+		    enum rt_canonicalize_mode mode)
+{
+    const struct rt_rpc_internal *rip;
+    struct rt_rpc_internal *crip;
+    fastf_t lengths[3];
+    fastf_t canonical_lengths[3];
+
+    if (!canonical || !canonical_to_input || !input || !tol ||
+	input->idb_type != ID_RPC)
+	return RT_CANONICALIZE_ERROR;
+    rip = (const struct rt_rpc_internal *)input->idb_ptr;
+    RT_RPC_CK_MAGIC(rip);
+    lengths[0] = MAGNITUDE(rip->rpc_B);
+    lengths[1] = MAGNITUDE(rip->rpc_H);
+    lengths[2] = rip->rpc_r;
+    if (!rt_make_canonical_similarity_frame(canonical_to_input,
+	    canonical_lengths, rip->rpc_V, rip->rpc_B, rip->rpc_H, lengths,
+	    3, tol, mode))
+	return RT_CANONICALIZE_ERROR;
+
+    canonical->idb_major_type = DB5_MAJORTYPE_BRLCAD;
+    canonical->idb_minor_type = ID_RPC;
+    canonical->idb_meth = &OBJ[ID_RPC];
+    BU_ALLOC(canonical->idb_ptr, struct rt_rpc_internal);
+    crip = (struct rt_rpc_internal *)canonical->idb_ptr;
+    crip->rpc_magic = RT_RPC_INTERNAL_MAGIC;
+    VSETALL(crip->rpc_V, 0.0);
+    VSET(crip->rpc_B, canonical_lengths[0], 0.0, 0.0);
+    VSET(crip->rpc_H, 0.0, 0.0, canonical_lengths[1]);
+    crip->rpc_r = canonical_lengths[2];
+    return RT_CANONICALIZE_OK;
+}
+
 /**
  * Import an RPC from the database format to the internal format.
  * Apply modeling transformations as well.

src/librt/primitives/superell/superell.c

@@ -998,11 +998,125 @@ rt_superell_mat(struct rt_db_internal *rop, const mat_t mat, const struct rt_db_
     MAT4X3VEC(top->a, mat, sa);
     MAT4X3VEC(top->b, mat, sb);
     MAT4X3VEC(top->c, mat, sc);
+    top->n = tip->n;
+    top->e = tip->e;
 
     return BRLCAD_OK;
 }
 
 
+static void
+superell_canonical_orientation(mat_t orientation, const vect_t a,
+			       const vect_t b, const vect_t c)
+{
+    MAT_IDN(orientation);
+    orientation[0] = a[X];
+    orientation[4] = a[Y];
+    orientation[8] = a[Z];
+    orientation[1] = b[X];
+    orientation[5] = b[Y];
+    orientation[9] = b[Z];
+    orientation[2] = c[X];
+    orientation[6] = c[Y];
+    orientation[10] = c[Z];
+}
+
+
+C_DECL int
+rt_superell_canonicalize(struct rt_db_internal *canonical,
+			 mat_t canonical_to_input,
+			 const struct rt_db_internal *input,
+			 const struct bn_tol *tol,
+			 enum rt_canonicalize_mode mode)
+{
+    static const fastf_t maximum_curvature = 10000.0;
+    const struct rt_superell_internal *sip;
+    struct rt_superell_internal *csip;
+    const fastf_t *input_axes[3];
+    vect_t axes[3];
+    vect_t cross;
+    fastf_t lengths[3];
+    fastf_t canonical_lengths[3];
+    fastf_t uniform_scale = 1.0;
+    mat_t orientation;
+    mat_t scale;
+    mat_t oriented_scale;
+    mat_t translate;
+
+    if (!canonical || !canonical_to_input || !input || !tol)
+	return RT_CANONICALIZE_ERROR;
+    if (mode < RT_CANONICALIZE_RIGID || mode > RT_CANONICALIZE_AFFINE)
+	return RT_CANONICALIZE_ERROR;
+    if (input->idb_type != ID_SUPERELL)
+	return RT_CANONICALIZE_ERROR;
+
+    sip = (const struct rt_superell_internal *)input->idb_ptr;
+    RT_SUPERELL_CK_MAGIC(sip);
+    if (!isfinite(sip->v[X]) || !isfinite(sip->v[Y]) ||
+	!isfinite(sip->v[Z]) || !isfinite(sip->n) || !isfinite(sip->e) ||
+	sip->n <= tol->dist || sip->e <= tol->dist ||
+	sip->n > maximum_curvature || sip->e > maximum_curvature)
+	return RT_CANONICALIZE_ERROR;
+
+    input_axes[0] = sip->a;
+    input_axes[1] = sip->b;
+    input_axes[2] = sip->c;
+    for (size_t i = 0; i < 3; i++) {
+	lengths[i] = MAGNITUDE(input_axes[i]);
+	if (!isfinite(lengths[i]) || lengths[i] <= tol->dist)
+	    return RT_CANONICALIZE_ERROR;
+	VSCALE(axes[i], input_axes[i], 1.0 / lengths[i]);
+    }
+    if (!NEAR_ZERO(VDOT(axes[0], axes[1]), tol->perp) ||
+	!NEAR_ZERO(VDOT(axes[0], axes[2]), tol->perp) ||
+	!NEAR_ZERO(VDOT(axes[1], axes[2]), tol->perp))
+	return RT_CANONICALIZE_ERROR;
+
+    /* Axis signs do not change a superellipsoid.  Keep the placement a
+     * proper transform by selecting the third sign from the input frame. */
+    VCROSS(cross, axes[0], axes[1]);
+    if (VDOT(cross, axes[2]) < 0.0)
+	VREVERSE(axes[2], axes[2]);
+
+    if (mode == RT_CANONICALIZE_AFFINE) {
+	VSETALL(canonical_lengths, 1.0);
+	for (size_t i = 0; i < 3; i++)
+	    VSCALE(axes[i], axes[i], lengths[i]);
+	superell_canonical_orientation(orientation, axes[0], axes[1], axes[2]);
+	MAT_IDN(translate);
+	MAT_DELTAS_VEC(translate, sip->v);
+	bn_mat_mul(canonical_to_input, translate, orientation);
+    } else {
+	if (mode == RT_CANONICALIZE_SIMILARITY)
+	    uniform_scale = fmax(lengths[0], fmax(lengths[1], lengths[2]));
+	for (size_t i = 0; i < 3; i++)
+	    canonical_lengths[i] = lengths[i] / uniform_scale;
+	superell_canonical_orientation(orientation, axes[0], axes[1], axes[2]);
+	MAT_IDN(scale);
+	scale[15] = 1.0 / uniform_scale;
+	bn_mat_mul(oriented_scale, orientation, scale);
+	MAT_IDN(translate);
+	MAT_DELTAS_VEC(translate, sip->v);
+	bn_mat_mul(canonical_to_input, translate, oriented_scale);
+    }
+
+    canonical->idb_major_type = DB5_MAJORTYPE_BRLCAD;
+    canonical->idb_minor_type = ID_SUPERELL;
+    canonical->idb_meth = &OBJ[ID_SUPERELL];
+    BU_ALLOC(canonical->idb_ptr, struct rt_superell_internal);
+    csip = (struct rt_superell_internal *)canonical->idb_ptr;
+    csip->magic = RT_SUPERELL_INTERNAL_MAGIC;
+    VSETALL(csip->v, 0.0);
+    VSET(csip->a, canonical_lengths[0], 0.0, 0.0);
+    VSET(csip->b, 0.0, canonical_lengths[1], 0.0);
+    VSET(csip->c, 0.0, 0.0, canonical_lengths[2]);
+    csip->n = sip->n;
+    csip->e = sip->e;
+
+    return RT_CANONICALIZE_OK;
+}
+
+
 /**
  * Import an superellipsoid/sphere from the database format to the
  * internal structure.  Apply modeling transformations as wsuperell.

src/librt/primitives/table.cpp

@@ -83,6 +83,7 @@ extern "C" {
     extern int rt_##name##_labels(struct rt_point_labels *pl, int pl_max, const mat_t xform, const struct rt_db_internal *ip, const struct bn_tol *utol); \
     extern const char *rt_##name##_keypoint(point_t *pt, const char *keystr, const mat_t xform, const struct rt_db_internal *ip, const struct bn_tol *tol); \
     extern int rt_##name##_mat(struct rt_db_internal *op, const mat_t mat, const struct rt_db_internal *ip); \
+    extern int rt_##name##_canonicalize(struct rt_db_internal *canonical, mat_t canonical_to_input, const struct rt_db_internal *input, const struct bn_tol *tol, enum rt_canonicalize_mode mode); \
     extern int rt_##name##_perturb(struct rt_db_internal **oip, const struct rt_db_internal *ip, int planar_only, fastf_t factor); \
     extern int rt_##name##_scene_obj(struct bv_scene_obj *vhead, struct directory *dp, struct db_i *dbip, const struct bg_tess_tol *ttol, const struct bn_tol *tol, const struct bview *info) \
 
@@ -246,7 +247,8 @@ const struct rt_functab OBJ[] = {
 	NULL,
 	NULL,
 	NULL,
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -298,7 +300,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_tor_mat),
 	RTFUNCTAB_FUNC_PERTURB_CAST(rt_tor_perturb), /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	RTFUNCTAB_FUNC_CANONICALIZE_CAST(rt_tor_canonicalize) /* canonicalize */
     },
 
     {
@@ -350,7 +353,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_tgc_mat),
 	RTFUNCTAB_FUNC_PERTURB_CAST(rt_tgc_perturb), /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	RTFUNCTAB_FUNC_VALIDATE_CAST(rt_tgc_functab_validate) /* validate */
+	RTFUNCTAB_FUNC_VALIDATE_CAST(rt_tgc_functab_validate), /* validate */
+	RTFUNCTAB_FUNC_CANONICALIZE_CAST(rt_tgc_canonicalize) /* canonicalize */
     },
 
     {
@@ -402,7 +406,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_ell_mat),
 	RTFUNCTAB_FUNC_PERTURB_CAST(rt_ell_perturb), /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	RTFUNCTAB_FUNC_VALIDATE_CAST(rt_ell_functab_validate) /* validate */
+	RTFUNCTAB_FUNC_VALIDATE_CAST(rt_ell_functab_validate), /* validate */
+	RTFUNCTAB_FUNC_CANONICALIZE_CAST(rt_ell_canonicalize)
     },
 
     {
@@ -454,7 +459,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_arb_mat),
 	RTFUNCTAB_FUNC_PERTURB_CAST(rt_arb_perturb), /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	RTFUNCTAB_FUNC_VALIDATE_CAST(rt_arb_functab_validate) /* validate */
+	RTFUNCTAB_FUNC_VALIDATE_CAST(rt_arb_functab_validate), /* validate */
+	RTFUNCTAB_FUNC_CANONICALIZE_CAST(rt_arb_canonicalize) /* canonicalize */
     },
 
     {
@@ -506,7 +512,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_ars_mat),
 	NULL, /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -558,7 +565,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_hlf_mat),
 	RTFUNCTAB_FUNC_PERTURB_CAST(rt_hlf_perturb), /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	RTFUNCTAB_FUNC_CANONICALIZE_CAST(rt_hlf_canonicalize)
     },
 
     {
@@ -610,7 +618,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_tgc_mat),
 	RTFUNCTAB_FUNC_PERTURB_CAST(rt_tgc_perturb), /* perturb - REC shares rt_tgc_internal */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	RTFUNCTAB_FUNC_CANONICALIZE_CAST(rt_tgc_canonicalize) /* canonicalize */
     },
 
     {
@@ -662,7 +671,8 @@ const struct rt_functab OBJ[] = {
 	NULL, /* mat */
 	NULL, /* perturb */
 	NULL  /* scene_obj */,
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -714,7 +724,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_nurb_mat),
 	NULL, /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -766,7 +777,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_ell_mat),
 	RTFUNCTAB_FUNC_PERTURB_CAST(rt_ell_perturb), /* perturb - SPH shares rt_ell_internal */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	RTFUNCTAB_FUNC_CANONICALIZE_CAST(rt_ell_canonicalize)
     },
 
     {
@@ -818,7 +830,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_nmg_mat),
 	NULL, /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -870,7 +883,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_ebm_mat),
 	NULL, /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -922,7 +936,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_vol_mat),
 	NULL, /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -974,7 +989,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_arbn_mat),
 	RTFUNCTAB_FUNC_PERTURB_CAST(rt_arbn_perturb), /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -1026,7 +1042,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_pipe_mat),
 	NULL, /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -1078,7 +1095,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_part_mat),
 	RTFUNCTAB_FUNC_PERTURB_CAST(rt_part_perturb), /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	RTFUNCTAB_FUNC_CANONICALIZE_CAST(rt_part_canonicalize) /* canonicalize */
     },
 
     {
@@ -1130,7 +1148,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_rpc_mat),
 	RTFUNCTAB_FUNC_PERTURB_CAST(rt_rpc_perturb), /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	RTFUNCTAB_FUNC_VALIDATE_CAST(rt_rpc_functab_validate) /* validate */
+	RTFUNCTAB_FUNC_VALIDATE_CAST(rt_rpc_functab_validate), /* validate */
+	RTFUNCTAB_FUNC_CANONICALIZE_CAST(rt_rpc_canonicalize) /* canonicalize */
     },
 
     {
@@ -1182,7 +1201,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_rhc_mat),
 	RTFUNCTAB_FUNC_PERTURB_CAST(rt_rhc_perturb), /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	RTFUNCTAB_FUNC_VALIDATE_CAST(rt_rhc_functab_validate) /* validate */
+	RTFUNCTAB_FUNC_VALIDATE_CAST(rt_rhc_functab_validate), /* validate */
+	RTFUNCTAB_FUNC_CANONICALIZE_CAST(rt_rhc_canonicalize) /* canonicalize */
     },
 
     {
@@ -1234,7 +1254,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_epa_mat),
 	RTFUNCTAB_FUNC_PERTURB_CAST(rt_epa_perturb), /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	RTFUNCTAB_FUNC_VALIDATE_CAST(rt_epa_functab_validate) /* validate */
+	RTFUNCTAB_FUNC_VALIDATE_CAST(rt_epa_functab_validate), /* validate */
+	RTFUNCTAB_FUNC_CANONICALIZE_CAST(rt_epa_canonicalize) /* canonicalize */
     },
 
     {
@@ -1286,7 +1307,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_ehy_mat),
 	RTFUNCTAB_FUNC_PERTURB_CAST(rt_ehy_perturb), /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	RTFUNCTAB_FUNC_VALIDATE_CAST(rt_ehy_functab_validate) /* validate */
+	RTFUNCTAB_FUNC_VALIDATE_CAST(rt_ehy_functab_validate), /* validate */
+	RTFUNCTAB_FUNC_CANONICALIZE_CAST(rt_ehy_canonicalize) /* canonicalize */
     },
 
     {
@@ -1338,7 +1360,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_eto_mat),
 	RTFUNCTAB_FUNC_PERTURB_CAST(rt_eto_perturb), /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	RTFUNCTAB_FUNC_VALIDATE_CAST(rt_eto_functab_validate) /* validate */
+	RTFUNCTAB_FUNC_VALIDATE_CAST(rt_eto_functab_validate), /* validate */
+	RTFUNCTAB_FUNC_CANONICALIZE_CAST(rt_eto_canonicalize) /* canonicalize */
     },
 
     {
@@ -1390,7 +1413,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_grp_mat),
 	NULL, /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -1442,7 +1466,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_joint_mat),
 	NULL, /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -1494,7 +1519,8 @@ const struct rt_functab OBJ[] = {
 	NULL, /* mat */
 	NULL, /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -1546,7 +1572,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_dsp_mat),
 	NULL, /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -1598,7 +1625,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_sketch_mat),
 	NULL, /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -1650,7 +1678,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_extrude_mat),
 	NULL, /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -1702,7 +1731,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_submodel_mat),
 	NULL, /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -1754,7 +1784,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_cline_mat),
 	NULL, /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -1806,7 +1837,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_bot_mat),
 	NULL, /* perturb */
 	NULL   /* scene_obj */,
-	RTFUNCTAB_FUNC_VALIDATE_CAST(rt_bot_functab_validate) /* validate */
+	RTFUNCTAB_FUNC_VALIDATE_CAST(rt_bot_functab_validate), /* validate */
+	RTFUNCTAB_FUNC_CANONICALIZE_CAST(rt_bot_canonicalize) /* canonicalize */
     },
 
     {
@@ -1858,7 +1890,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_comb_mat),
 	NULL, /* perturb */
 	NULL   /* scene_obj */,
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -1912,7 +1945,8 @@ const struct rt_functab OBJ[] = {
 	NULL, /* mat */
 	NULL, /* perturb */
 	NULL  /* scene_obj */,
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -1964,7 +1998,8 @@ const struct rt_functab OBJ[] = {
 	NULL, /* mat */
 	NULL, /* perturb */
 	NULL  /* scene_obj */,
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -2018,7 +2053,8 @@ const struct rt_functab OBJ[] = {
 	NULL, /* mat */
 	NULL, /* perturb */
 	NULL  /* scene_obj */,
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -2070,7 +2106,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_superell_mat),
 	RTFUNCTAB_FUNC_PERTURB_CAST(rt_superell_perturb), /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	RTFUNCTAB_FUNC_CANONICALIZE_CAST(rt_superell_canonicalize) /* canonicalize */
     },
 
     {
@@ -2122,7 +2159,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_metaball_mat),
 	NULL, /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -2174,7 +2212,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_brep_mat),
 	NULL, /* perturb */
 	NULL   /* scene_obj */,
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -2226,7 +2265,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_hyp_mat),
 	RTFUNCTAB_FUNC_PERTURB_CAST(rt_hyp_perturb), /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -2278,7 +2318,8 @@ const struct rt_functab OBJ[] = {
 	NULL, /* mat */
 	NULL, /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -2330,7 +2371,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_revolve_mat),
 	NULL, /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -2382,7 +2424,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_pnts_mat),
 	NULL, /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -2434,7 +2477,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_annot_mat),
 	NULL, /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -2486,7 +2530,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_hrt_mat),
 	RTFUNCTAB_FUNC_PERTURB_CAST(rt_hrt_perturb), /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
 
@@ -2539,7 +2584,8 @@ const struct rt_functab OBJ[] = {
 	RTFUNCTAB_FUNC_MAT_CAST(rt_datum_mat),
 	NULL, /* perturb */
 	RTFUNCTAB_FUNC_SCENE_OBJ_CAST(rt_generic_scene_obj),
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
 
@@ -2592,7 +2638,8 @@ const struct rt_functab OBJ[] = {
 	NULL, /* mat */
 	NULL, /* perturb */
 	NULL  /* scene_obj */,
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -2644,7 +2691,8 @@ const struct rt_functab OBJ[] = {
 	NULL, /* mat */
 	NULL, /* perturb */
 	NULL  /* scene_obj */,
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     },
 
     {
@@ -2696,7 +2744,8 @@ const struct rt_functab OBJ[] = {
 	NULL, /* mat */
 	NULL, /* perturb */
 	NULL  /* scene_obj */,
-	NULL /* validate */
+	NULL, /* validate */
+	NULL /* canonicalize */
     }
 };
 

src/librt/primitives/tgc/tgc.c

@@ -1757,6 +1757,248 @@ rt_tgc_mat(struct rt_db_internal *rop, const mat_t mat, const struct rt_db_inter
 }
 
 
+static void
+tgc_swap_vectors(vect_t a, vect_t b)
+{
+    vect_t tmp;
+
+    VMOVE(tmp, a);
+    VMOVE(a, b);
+    VMOVE(b, tmp);
+}
+
+
+static int
+tgc_ratio_pair_greater(fastf_t a0, fastf_t a1, fastf_t b0, fastf_t b1,
+		       fastf_t tolerance)
+{
+    fastf_t ahi = fmax(a0, a1);
+    fastf_t alo = fmin(a0, a1);
+    fastf_t bhi = fmax(b0, b1);
+    fastf_t blo = fmin(b0, b1);
+
+    if (ahi > bhi + tolerance)
+	return 1;
+    if (bhi > ahi + tolerance)
+	return 0;
+    return alo > blo + tolerance;
+}
+
+
+C_DECL int
+rt_tgc_canonicalize(struct rt_db_internal *canonical,
+		    mat_t canonical_to_input,
+		    const struct rt_db_internal *input,
+		    const struct bn_tol *tol,
+		    enum rt_canonicalize_mode mode)
+{
+    const struct rt_tgc_internal *tip;
+    struct rt_tgc_internal *ctip;
+    point_t v;
+    vect_t h, a, b, c, d;
+    vect_t u, w, normal;
+    fastf_t mag_h, mag_a, mag_b, mag_c, mag_d;
+    fastf_t ratio_c, ratio_d;
+    fastf_t h_a, h_b, h_w;
+    fastf_t uniform_scale = 1.0;
+    fastf_t ratio_tolerance;
+    int a_sign = 1;
+    int b_sign = 1;
+    int a_projection_zero;
+    int b_projection_zero;
+    int reverse_ends = 0;
+
+    if (!canonical || !canonical_to_input || !input || !tol)
+	return RT_CANONICALIZE_ERROR;
+    if (mode < RT_CANONICALIZE_RIGID || mode > RT_CANONICALIZE_AFFINE)
+	return RT_CANONICALIZE_ERROR;
+    if (input->idb_type != ID_TGC && input->idb_type != ID_REC)
+	return RT_CANONICALIZE_ERROR;
+
+    tip = (const struct rt_tgc_internal *)input->idb_ptr;
+    RT_TGC_CK_MAGIC(tip);
+    VMOVE(v, tip->v);
+    VMOVE(h, tip->h);
+    VMOVE(a, tip->a);
+    VMOVE(b, tip->b);
+    VMOVE(c, tip->c);
+    VMOVE(d, tip->d);
+
+    mag_h = MAGNITUDE(h);
+    mag_a = MAGNITUDE(a);
+    mag_b = MAGNITUDE(b);
+    mag_c = MAGNITUDE(c);
+    mag_d = MAGNITUDE(d);
+    if (!isfinite(v[X]) || !isfinite(v[Y]) || !isfinite(v[Z]) ||
+	!isfinite(mag_h) || !isfinite(mag_a) || !isfinite(mag_b) ||
+	!isfinite(mag_c) || !isfinite(mag_d) || mag_h <= tol->dist)
+	return RT_CANONICALIZE_ERROR;
+    if ((mag_a <= tol->dist && mag_c <= tol->dist) ||
+	(mag_b <= tol->dist && mag_d <= tol->dist))
+	return RT_CANONICALIZE_ERROR;
+    if ((mag_a <= tol->dist || mag_b <= tol->dist) &&
+	(mag_c <= tol->dist || mag_d <= tol->dist))
+	return RT_CANONICALIZE_ERROR;
+    if (input->idb_type == ID_REC) {
+	vect_t difference;
+	VSUB2(difference, a, c);
+	if (MAGNITUDE(difference) > tol->dist)
+	    return RT_CANONICALIZE_ERROR;
+	VSUB2(difference, b, d);
+	if (MAGNITUDE(difference) > tol->dist ||
+	    !NEAR_ZERO(VDOT(h, a) / (mag_h * mag_a), tol->perp) ||
+	    !NEAR_ZERO(VDOT(h, b) / (mag_h * mag_b), tol->perp))
+	    return RT_CANONICALIZE_ERROR;
+    }
+
+    /* A complete end supplies the affine basis.  When both ends are
+     * complete, select between the two reciprocal taper descriptions so
+     * reversing an otherwise identical TGC does not change its canonical
+     * form. */
+    ratio_tolerance = fmax(tol->perp, SMALL_FASTF);
+    if (mag_a <= tol->dist || mag_b <= tol->dist) {
+	reverse_ends = 1;
+    } else if (mag_c > tol->dist && mag_d > tol->dist) {
+	fastf_t taper_c = mag_c / mag_a;
+	fastf_t taper_d = mag_d / mag_b;
+	int forward_greater = tgc_ratio_pair_greater(taper_c, taper_d,
+	    1.0 / taper_c, 1.0 / taper_d, ratio_tolerance);
+	int reverse_greater = tgc_ratio_pair_greater(1.0 / taper_c,
+	    1.0 / taper_d, taper_c, taper_d, ratio_tolerance);
+
+	reverse_ends = forward_greater;
+	if (!forward_greater && !reverse_greater &&
+	    mode != RT_CANONICALIZE_AFFINE)
+	    reverse_ends = tgc_ratio_pair_greater(mag_c, mag_d,
+		mag_a, mag_b, tol->dist);
+    }
+    if (reverse_ends) {
+	VADD2(v, v, h);
+	VREVERSE(h, h);
+	tgc_swap_vectors(a, c);
+	tgc_swap_vectors(b, d);
+	mag_a = MAGNITUDE(a);
+	mag_b = MAGNITUDE(b);
+	mag_c = MAGNITUDE(c);
+	mag_d = MAGNITUDE(d);
+    }
+
+    ratio_c = mag_c / mag_a;
+    ratio_d = mag_d / mag_b;
+
+    /* Axis interchange is a parameterization change, not a geometry
+     * change.  Taper ratio is affine-invariant, so use it as the primary
+     * ordering key. */
+    if (ratio_d > ratio_c + ratio_tolerance ||
+	(mode != RT_CANONICALIZE_AFFINE &&
+	 NEAR_EQUAL(ratio_c, ratio_d, ratio_tolerance) &&
+	 mag_b > mag_a + tol->dist)) {
+	tgc_swap_vectors(a, b);
+	tgc_swap_vectors(c, d);
+	{
+	    fastf_t tmp = mag_a;
+	    mag_a = mag_b;
+	    mag_b = tmp;
+	    tmp = mag_c;
+	    mag_c = mag_d;
+	    mag_d = tmp;
+	    tmp = ratio_c;
+	    ratio_c = ratio_d;
+	    ratio_d = tmp;
+	}
+    }
+
+    VSCALE(u, a, 1.0 / mag_a);
+    VSCALE(w, b, 1.0 / mag_b);
+    if (!NEAR_ZERO(VDOT(u, w), tol->perp))
+	return RT_CANONICALIZE_ERROR;
+    if (mag_c > tol->dist && !NEAR_EQUAL(VDOT(u, c) / mag_c, 1.0, tol->perp))
+	return RT_CANONICALIZE_ERROR;
+    if (mag_d > tol->dist && !NEAR_EQUAL(VDOT(w, d) / mag_d, 1.0, tol->perp))
+	return RT_CANONICALIZE_ERROR;
+
+    h_a = VDOT(h, u);
+    h_b = VDOT(h, w);
+    VCROSS(normal, u, w);
+    h_w = VDOT(h, normal);
+    if (fabs(h_w) <= mag_h * tol->perp)
+	return RT_CANONICALIZE_ERROR;
+
+    /* A/C and B/D may each be negated without changing the TGC.  Use the
+     * height projections to resolve those signs without referring to world
+     * coordinates, which keeps the result invariant under placement. */
+    a_projection_zero = fabs(h_a) <= tol->dist;
+    b_projection_zero = fabs(h_b) <= tol->dist;
+    if (!a_projection_zero && h_a < 0.0)
+	a_sign = -1;
+    if (!b_projection_zero && h_b < 0.0)
+	b_sign = -1;
+    if (a_projection_zero && a_sign * b_sign * h_w < 0.0)
+	a_sign = -1;
+    else if (b_projection_zero && a_sign * b_sign * h_w < 0.0)
+	b_sign = -1;
+    if (a_sign < 0)
+	VREVERSE(u, u);
+    if (b_sign < 0)
+	VREVERSE(w, w);
+    VCROSS(normal, u, w);
+
+    MAT_IDN(canonical_to_input);
+    MAT_DELTAS_VEC(canonical_to_input, v);
+    if (mode == RT_CANONICALIZE_AFFINE) {
+	canonical_to_input[0] = u[X] * mag_a;
+	canonical_to_input[4] = u[Y] * mag_a;
+	canonical_to_input[8] = u[Z] * mag_a;
+	canonical_to_input[1] = w[X] * mag_b;
+	canonical_to_input[5] = w[Y] * mag_b;
+	canonical_to_input[9] = w[Z] * mag_b;
+	canonical_to_input[2] = h[X];
+	canonical_to_input[6] = h[Y];
+	canonical_to_input[10] = h[Z];
+    } else {
+	if (mode == RT_CANONICALIZE_SIMILARITY)
+	    uniform_scale = fmax(mag_h, fmax(mag_a, fmax(mag_b, fmax(mag_c, mag_d))));
+	canonical_to_input[0] = u[X];
+	canonical_to_input[4] = u[Y];
+	canonical_to_input[8] = u[Z];
+	canonical_to_input[1] = w[X];
+	canonical_to_input[5] = w[Y];
+	canonical_to_input[9] = w[Z];
+	canonical_to_input[2] = normal[X];
+	canonical_to_input[6] = normal[Y];
+	canonical_to_input[10] = normal[Z];
+	canonical_to_input[15] = 1.0 / uniform_scale;
+	canonical_to_input[3] /= uniform_scale;
+	canonical_to_input[7] /= uniform_scale;
+	canonical_to_input[11] /= uniform_scale;
+    }
+
+    canonical->idb_major_type = DB5_MAJORTYPE_BRLCAD;
+    canonical->idb_minor_type = input->idb_minor_type;
+    canonical->idb_meth = &OBJ[input->idb_minor_type];
+    BU_ALLOC(canonical->idb_ptr, struct rt_tgc_internal);
+    ctip = (struct rt_tgc_internal *)canonical->idb_ptr;
+    ctip->magic = RT_TGC_INTERNAL_MAGIC;
+    VSETALL(ctip->v, 0.0);
+    if (mode == RT_CANONICALIZE_AFFINE) {
+	VSET(ctip->h, 0.0, 0.0, 1.0);
+	VSET(ctip->a, 1.0, 0.0, 0.0);
+	VSET(ctip->b, 0.0, 1.0, 0.0);
+	VSET(ctip->c, ratio_c, 0.0, 0.0);
+	VSET(ctip->d, 0.0, ratio_d, 0.0);
+    } else {
+	VSET(ctip->h, VDOT(h, u) / uniform_scale,
+	    VDOT(h, w) / uniform_scale, VDOT(h, normal) / uniform_scale);
+	VSET(ctip->a, mag_a / uniform_scale, 0.0, 0.0);
+	VSET(ctip->b, 0.0, mag_b / uniform_scale, 0.0);
+	VSET(ctip->c, mag_c / uniform_scale, 0.0, 0.0);
+	VSET(ctip->d, 0.0, mag_d / uniform_scale, 0.0);
+    }
+
+    return RT_CANONICALIZE_OK;
+}
+
+
 /**
  * Import a TGC from the database format to the internal format.
  * Apply modeling transformations as well.

src/librt/primitives/tor/tor.c

@@ -1840,6 +1840,68 @@ rt_tor_mat(struct rt_db_internal *rop, const mat_t mat, const struct rt_db_inter
     return BRLCAD_OK;
 }
 
+
+C_DECL int
+rt_tor_canonicalize(struct rt_db_internal *canonical,
+		    mat_t canonical_to_input,
+		    const struct rt_db_internal *input,
+		    const struct bn_tol *tol,
+		    enum rt_canonicalize_mode mode)
+{
+    const struct rt_tor_internal *tip;
+    struct rt_tor_internal *ctip;
+    vect_t normal;
+    vect_t zaxis = VINIT_ZERO;
+    mat_t rotate;
+    mat_t scale_mat;
+    mat_t rotated_scale;
+    mat_t translate;
+    fastf_t normal_mag;
+    fastf_t uniform_scale = 1.0;
+
+    if (!canonical || !canonical_to_input || !input || !tol)
+	return RT_CANONICALIZE_ERROR;
+    if (mode < RT_CANONICALIZE_RIGID || mode > RT_CANONICALIZE_AFFINE)
+	return RT_CANONICALIZE_ERROR;
+    if (input->idb_type != ID_TOR)
+	return RT_CANONICALIZE_ERROR;
+
+    tip = (const struct rt_tor_internal *)input->idb_ptr;
+    RT_TOR_CK_MAGIC(tip);
+    normal_mag = MAGNITUDE(tip->h);
+    if (normal_mag <= tol->dist || tip->r_a <= tol->dist || tip->r_h <= tol->dist)
+	return RT_CANONICALIZE_ERROR;
+
+    VSCALE(normal, tip->h, 1.0 / normal_mag);
+    VSET(zaxis, 0.0, 0.0, 1.0);
+    bn_mat_fromto(rotate, zaxis, normal, tol);
+
+    if (mode != RT_CANONICALIZE_RIGID)
+	uniform_scale = tip->r_a;
+    MAT_IDN(scale_mat);
+    scale_mat[15] = 1.0 / uniform_scale;
+    bn_mat_mul(rotated_scale, rotate, scale_mat);
+    MAT_IDN(translate);
+    MAT_DELTAS_VEC(translate, tip->v);
+    bn_mat_mul(canonical_to_input, translate, rotated_scale);
+
+    canonical->idb_major_type = DB5_MAJORTYPE_BRLCAD;
+    canonical->idb_minor_type = ID_TOR;
+    canonical->idb_meth = &OBJ[ID_TOR];
+    BU_ALLOC(canonical->idb_ptr, struct rt_tor_internal);
+    ctip = (struct rt_tor_internal *)canonical->idb_ptr;
+    ctip->magic = RT_TOR_INTERNAL_MAGIC;
+    VSETALL(ctip->v, 0.0);
+    VSET(ctip->h, 0.0, 0.0, 1.0);
+    ctip->r_a = tip->r_a / uniform_scale;
+    ctip->r_h = tip->r_h / uniform_scale;
+    VSET(ctip->a, ctip->r_a, 0.0, 0.0);
+    VSET(ctip->b, 0.0, ctip->r_a, 0.0);
+    ctip->r_b = ctip->r_a;
+
+    return RT_CANONICALIZE_OK;
+}
+
 /**
  * Taken from the database record:
  * v vertex (point) of center of torus.

src/librt/tests/CMakeLists.txt

@@ -55,6 +55,9 @@ brlcad_addexec(rt_diff diff.c "${RT_TEST_LIBS}" TEST)
 brlcad_addexec(rt_default_object default_object.c "librt;libwdb;libbu;${M_LIBRARY}" TEST)
 brlcad_add_test(NAME rt_default_object COMMAND rt_default_object)
 
+brlcad_addexec(rt_canonicalize canonicalize.c "${RT_TEST_LIBS}" TEST)
+brlcad_add_test(NAME rt_canonicalize COMMAND rt_canonicalize)
+
 brlcad_addexec(rt_db_name db_name.c "librt;libbu" TEST)
 brlcad_add_test(NAME rt_db_name COMMAND rt_db_name)
 
@@ -216,6 +219,7 @@ set(
   binunif.c
   brep_io.cpp
   bound_internal.c
+  canonicalize.c
   crofton.c
   db_name.c
   default_object.c

src/libtclcad/commands.c

@@ -796,6 +796,7 @@ static struct to_cmdtab ged_cmds[] = {
     {"track",	(char *)0, TO_UNLIMITED, to_pass_through_func, ged_exec_track},
     {"tree",	(char *)0, TO_UNLIMITED, to_pass_through_func, ged_exec_tree},
     {"unhide",	(char *)0, TO_UNLIMITED, to_pass_through_func, ged_exec_unhide},
+    {"unpush",	(char *)0, TO_UNLIMITED, to_pass_through_func, ged_exec_unpush},
     {"units",	(char *)0, TO_UNLIMITED, to_pass_through_func, ged_exec_units},
     {"v2m_point",	"x y z", 5, to_view_func, ged_exec_v2m_point},
     {"vdraw",	(char *)0, TO_UNLIMITED, to_pass_through_and_refresh_func, ged_exec_vdraw},

src/mged/setup.c

@@ -457,6 +457,7 @@ static struct cmdtab mged_cmdtab[] = {
     {MGED_CMD_MAGIC, "translate", f_tr_obj, GED_FUNC_PTR_NULL, NULL},
     {MGED_CMD_MAGIC, "tree", cmd_ged_plain_wrapper, ged_exec_tree, NULL},
     {MGED_CMD_MAGIC, "unhide", cmd_ged_plain_wrapper, ged_exec_unhide, NULL},
+    {MGED_CMD_MAGIC, "unpush", cmd_ged_plain_wrapper, ged_exec_unpush, NULL},
     {MGED_CMD_MAGIC, "units", cmd_units, GED_FUNC_PTR_NULL, NULL},
     {MGED_CMD_MAGIC, "v2m_point", cmd_ged_plain_wrapper, ged_exec_v2m_point, NULL},
     {MGED_CMD_MAGIC, "vars", f_set, GED_FUNC_PTR_NULL, NULL},