assembly#

class cryocat.analysis.assembly.BlockDefinition(sites, pairing=(('site', 'site'),), flip_site=1, symmetry=None)#

Bases: object

Geometry and contact rules for one class of building block.

Parameters:
sitestuple[ContactSite, …]

Ordered, counter-clockwise sequence of contact sites. At least one site is required.

pairingtuple[tuple[str, str], …], default=((“site”, “site”),)

Allowed contact pairs as (type_a, type_b) tuples. Both type strings must appear in the site types of sites.

flip_siteint, default=1

1-based index of the site used to define the polarity-flip axis. The flip rotation is a 180° rotation around the in-plane part of this site’s vector.

foldint or None, default=None

Cyclic symmetry order. When set, must equal n_sites. Used by PleomorphicSurface.get_sites_as_motl() to delegate to split_in_asymmetric_subunits() so that the returned positions also carry per-symmetry-copy orientations. Set automatically by cyclic(); None for non-cyclic definitions such as those created with microtubule().

Raises:
ValueError

If sites is empty, any site vector has the wrong length or zero in-plane part, flip_site is out of range, any pairing type is missing from the site types, the sites are not listed in counter-clockwise azimuth order, or fold is set but does not equal n_sites.

classmethod cyclic(symmetry, site_shift, site_type='site')#

Create a C_n-symmetric block with equally-spaced sites.

Site k (1-based) is site_shift rotated about +z by 360·(k-1)/n degrees — the same convention as split_in_asymmetric_subunits().

Parameters:
symmetrySymmetry

Cyclic symmetry specifier: an integer n, or a string "Cn" (e.g. "C3"). Non-cyclic groups (D, T, O, I) raise ValueError.

site_shiftTripletLike

Shift from the block origin to the first contact site, in voxels, expressed in the block’s local frame. Typically about half the centre-to-centre distance, pointing along one leg.

site_typestr, default=”site”

Type label for all sites.

Returns:
BlockDefinition

A frozen BlockDefinition with a single site (the shift) and symmetry set to "C{n}".

Examples

>>> BlockDefinition.cyclic(3, [-5, 0, 0])      # C3, arm along −x
>>> BlockDefinition.cyclic("C6", [-3, 0, -1])  # C6, tilted arm
property effective_n_sites#

Number of sites after symmetry expansion.

For cyclic definitions (created via cyclic()), returns the cyclic order n. For non-cyclic definitions, returns n_sites.

flip_site = 1#
classmethod microtubule(lateral_length, axial_length, lateral_elevation=0.0)#

Create a 4-site block definition for a microtubule protofilament.

The block frame convention is: x points toward the right neighbouring protofilament, y points toward the plus end of the protofilament, and z points radially outward from the tube axis. The four sites are placed at:

  • site 1 — lateral_right: (lc, 0, −ls)

  • site 2 — plus: (0, axial_length, 0)

  • site 3 — lateral_left: (−lc, 0, −ls)

  • site 4 — minus: (0, −axial_length, 0)

where lc = lateral_length × cos(lateral_elevation) and ls = lateral_length × sin(lateral_elevation).

Parameters:
lateral_lengthfloat

Distance from the block origin to each lateral site in voxels.

axial_lengthfloat

Distance from the block origin to each axial (plus/minus) site.

lateral_elevationfloat, default=0.0

Elevation (degrees) of the lateral sites below the xy-plane in the block frame. Positive values tilt the sites toward the tube axis (negative z).

Returns:
BlockDefinition

Four-site definition with pairing=(("plus", "minus"), ("lateral_right", "lateral_left")) and flip_site=2 (plus end).

property n_sites#

Number of contact sites stored in this definition.

pairing = (('site', 'site'),)#
property site_types#

Ordered tuple of site-type labels.

site_vectors()#

Return site vectors as an (n_sites, 3) float array.

sites#
symmetry = None#
class cryocat.analysis.assembly.BlockLayer(blocks, block_definition, *, block_type_column='geom1', block_id_column='geom3', affiliation_column='object_id', site_index_column='geom1', site_type_column='geom2', assembly_id_column='geom4', face_id_column='geom4', face_size_column='geom5', source_block_count_column='geom2', pixel_size=1.0, tomo_id_column='tomo_id', max_cycle_length=25, bonding_tolerance=None, tip_tolerance=None, angle_tolerance=None, merge_tolerance=None, clean_distance=None)#

Bases: object

Block layer bundling blocks motl, block definition, column names, and scale.

Validates inputs once; pass the result to PleomorphicSurface as block_layer=. Use PleomorphicSurface.from_blocks() instead of constructing BlockLayer directly unless you already have a BlockDefinition.

Parameters:
blocksMotlSource

Block particle list. Loaded with Motl.load; the input is never modified.

block_definitionBlockDefinition or dict[float, BlockDefinition]

Geometry for one block type, or a mapping from block-type float values to per-type definitions.

block_type_columnMotlColumn, default=”class”

Column in blocks that identifies the block class when block_definition is a dict.

block_id_columnMotlColumn, default=”geom3”

Column used to identify blocks in site-motl output.

affiliation_columnMotlColumn, default=”object_id”

Column whose value identifies which assembly each block belongs to. Blocks with different values in this column are never connected, even when their sites are within max_distance.

site_index_columnMotlColumn, default=”geom1”

Column written with the 1-based site index.

site_type_columnMotlColumn, default=”geom2”

Column written with the integer site-type code.

assembly_id_columnMotlColumn, default=”geom4”

Column written with the assembly (connected-component) id in the faces motl produced by get_faces_as_motl().

face_id_columnMotlColumn, default=”geom5”

Column written with the face id in both the faces motl (get_faces_as_motl()) and the missing-block motl (get_missing_block_motl()).

face_size_columnMotlColumn, default=”score”

Column written with the face size (vertex count) in the faces motl.

source_block_count_columnMotlColumn, default=”subtomo_mean”

Column written with the number of distinct source blocks in the gaps motl (get_gaps_as_motl()).

pixel_sizefloat, default=1.0

Ångström-per-voxel scale factor.

tomo_id_columnMotlColumn, default=”tomo_id”

Column identifying the tomogram each block belongs to.

bonding_tolerancefloat or None, default=None

Half-width (voxels) of the acceptance window around the nominal centre-to-centre bonding distance 2 × |site_shift|. None derives to 0.4 × |site_shift| at call time.

tip_tolerancefloat or None, default=None

Radius tolerance (voxels) for accepting a leg tip as near the candidate centre. None derives to 0.3 × |site_shift| at call time.

angle_tolerancefloat or None, default=None

Maximum allowed angular residual (degrees) for the Cn pairwise spacing check. None derives to 0.125 × (360 / n) at call time (giving 15° for C3, 7.5° for C6, 1.15° for C39). Raises ValueError if the resolved value is >= 0.5 × (360 / n).

merge_tolerancefloat or None, default=None

Half-width (voxels) used as the default merge distance in infer_missing_blocks() when merge_distance is None. None derives to 0.5 × |site_shift|.

clean_distancefloat or None, default=None

Default clean distance (voxels) for infer_missing_blocks(). None derives to 0.4 × |site_shift| at call time.

Raises:
ValueError

If block_definition is None, if blocks has non-unique subtomo_id values, or if block_definition is a dict and the block-type column contains values without a matching definition.

class cryocat.analysis.assembly.Chain(traced_motl, pixel_size=1.0, column_name='tomo_id', chain_id_col='object_id', order_id_col='geom2', step_dist_col='geom4')#

Bases: object

Generic linear-chain analysis on a traced motl.

A traced motl is a motl whose object_id (or other store_idx1) identifies the chain that each particle belongs to, and whose geom2 (or other store_idx2) gives the position within the chain. geom4 (or store_dist) holds the distance to the next particle in the chain. These columns are produced by nnana.trace_chains().

Use the class constructor when you already have a traced motl; use from_motls() / from_motl() when you have raw entry/exit motls and want tracing performed in one step.

Parameters:
traced_motlstr or Motl
pixel_sizefloat, default=1.0
column_namestr, default=’tomo_id’
chain_id_colstr, default=’object_id’
order_id_colstr, default=’geom2’
step_dist_colstr, default=’geom4’
add_traced_info(input_motl, output_motl_path=None, sort_by_subtomo=True, occupancy_id='geom1')#

Copy chain columns from the traced motl onto input_motl.

The columns occupancy_id, order_id_col, step_dist_col, and chain_id_col are transferred by matching subtomo_id values. If occupancy has not yet been computed it is computed first.

Parameters:
input_motlMotlSource

Target motl that will receive the chain annotations.

output_motl_pathPathOrStr, optional

Path to save the annotated motl.

sort_by_subtomobool, default=True

Sort both motls by subtomo_id before copying to ensure correct row alignment.

occupancy_idMotlColumn, default=’geom1’

Column that holds (or will hold) the chain-length value.

Returns:
Motl

A new Motl with chain columns populated.

Raises:
ValueError

When input_motl contains different subtomogram IDs than the traced motl.

classmethod from_motl(motl, max_distance, min_distance=0, column_name='tomo_id', pixel_size=1.0, output_motl=None, chain_id_col='object_id', order_id_col='geom2', step_dist_col='geom4')#

Build a Chain by tracing a single motl (single-site mode).

Useful for structures where each particle has only one binding site, such as nucleosomes in a chromatin chain. Passes the same motl as both entry and exit to nnana.trace_chains().

Parameters:
motlMotlSource

Particle list to trace.

max_distancefloat

Maximum allowed step distance (in voxels).

min_distancefloat, default=0

Minimum allowed step distance.

column_nameMotlColumn, default=’tomo_id’

Column used to group particles before tracing.

pixel_sizefloat, default=1.0

Pixel size in Å.

output_motlPathOrStr, optional

Path to save the traced motl.

chain_id_colMotlColumn, default=’object_id’

Column that receives the chain identifier.

order_id_colMotlColumn, default=’geom2’

Column that receives the within-chain position index.

step_dist_colMotlColumn, default=’geom4’

Column that receives the step distance.

Returns:
Chain
classmethod from_motls(motl_entry, motl_exit, max_distance, min_distance=0, column_name='tomo_id', pixel_size=1.0, output_motl=None, chain_id_col='object_id', order_id_col='geom2', step_dist_col='geom4')#

Build a Chain by tracing an entry/exit motl pair.

Calls nnana.trace_chains() on motl_entry and motl_exit and wraps the resulting traced motl in a Chain instance.

Parameters:
motl_entryMotlSource

Entry-site particle list.

motl_exitMotlSource

Exit-site particle list.

max_distancefloat

Maximum allowed step distance (in voxels) between successive entry/exit pairs.

min_distancefloat, default=0

Minimum allowed step distance.

column_nameMotlColumn, default=’tomo_id’

Column used to group particles before tracing.

pixel_sizefloat, default=1.0

Pixel size in Å; stored on the instance for later distance scaling.

output_motlPathOrStr, optional

Path to save the traced motl.

chain_id_colMotlColumn, default=’object_id’

Column that receives the chain identifier.

order_id_colMotlColumn, default=’geom2’

Column that receives the within-chain position index.

step_dist_colMotlColumn, default=’geom4’

Column that receives the step distance.

Returns:
Chain
get_chain_stats(min_chain_size=2)#

Per-step statistics across all chains.

Parameters:
min_chain_sizeint, default=2

Skip chains shorter than this.

Returns:
pandas.DataFrame

Columns: chain_size, distance, coord_x/y/z, coord_rx/ry/rz, angular_distance, rot_x/y/z, phi/theta/psi, type.

get_class_chain_occupancies(mode='mp', occupancy_id='geom1', class_col='class')#

Return per-class chain-occupancy counts broken down by chain type.

Parameters:
mode{‘mp’, ‘mdp’}, default=’mp’

Breakdown resolution:

'mp'

Two categories — monomers (chain length 1) vs. polysomes (chain length > 1).

'mdp'

Three categories — monomers, disomes (length 2), and polysomes (length > 2).

occupancy_idMotlColumn, default=’geom1’

Column that holds chain-length values. Computed automatically if not yet present.

class_colMotlColumn, default=’class’

Column used to group particles by class.

Returns:
pandas.DataFrame

For mode='mp': columns class, particle_number, chain_type, percentage. For mode='mdp': columns class, particle_number, chain_type.

Raises:
ValueError

When mode is not 'mp' or 'mdp'.

get_occupancy(occupancy_id='geom1', output_motl=None)#

Write the chain length (occupancy) per particle into occupancy_id.

Each particle receives the length of its chain, i.e. the maximum within-chain position index in order_id_col. The result is stored in self.traced_motl in place.

Parameters:
occupancy_idMotlColumn, default=’geom1’

Column name that receives the chain-length value.

output_motlPathOrStr, optional

Path to save the updated motl.

Returns:
Motl

The updated self.traced_motl.

get_step_stats()#

Per-step angular and normal-vector distances along each chain.

Within each chain particles are sorted by order_id_col and orientation distances are computed for every consecutive pair (particle i → particle i + 1). Calls nnana.angular_distances() with rotation_type="all" on the two paired ZXZ-angle arrays.

Grouping and sorting use self.column_name, self.chain_id_col, and self.order_id_col; these come from the Chain instance and are never hardcoded. A chain of one particle contributes no rows.

Returns:
pandas.DataFrame

One row per step; n particles → n − 1 rows.

<column_name>

—

value of self.column_name (e.g. tomo_id)

chain_id

—

value of self.chain_id_col

step

—

order_id_col value of the upstream particle

step_dist

motl’s own units

step_dist_col of the upstream particle; NaN when the column is absent

angular_distance

degrees

SO(3) geodesic distance between orientations

cone_distance

degrees

angle between z-axes (normal-vector distance)

in_plane_distance

degrees

in-plane rotation component

Notes

Sorting is always by order_id_col within each group — never by DataFrame row order. Duplicate order_id_col values within one chain leave the order among tied particles undefined (stable sort, arbitrary tie-break). Gaps in order_id_col are harmless.

class cryocat.analysis.assembly.CnComplex(motl, symmetry, *, affiliation_column='object_id', order_column='geom1', tomo_id_column='tomo_id', center_method='circle_fit')#

Bases: SymmetricComplex

Cyclic Cn-symmetric ring assembly.

Extends SymmetricComplex with methods specific to cyclic symmetry: subunit ordering, affiliation clustering, occupancy analysis, and diameter computation.

Parameters:
motlMotlSource

Subunit particle list.

symmetrySymmetry

Cyclic fold, e.g. "C8" or 8. Dihedral or Platonic symmetries raise ValueError.

affiliation_columnMotlColumn, default=’object_id’

Column that identifies which object each particle belongs to.

order_columnMotlColumn, default=’geom1’

Column that assign_subunit_order() writes cyclic indices into.

tomo_id_columnMotlColumn, default=’tomo_id’

Column that identifies the tomogram.

center_method{‘circle_fit’, ‘barycentric’}, default=’circle_fit’

Algorithm used by get_centers_as_motl() and related helpers. 'circle_fit' falls back to barycentric when the fit fails.

Raises:
ValueError

When symmetry is not a cyclic Cn group.

add_component(motl=None, name='basket', *, assign_affiliation=False, assign_su_order=False, copy_su_order=False, reference_ring=0, affiliation_max_distance=None, order_max_distance=None)#

Register a non-ring component and optionally affiliate and order it.

The component is stored in self._components[name] and never added to self._ring_motls, so it is invisible to all ring-fitting, merge, and statistics methods.

Behaviour is driven only by explicit flags — no inference from column presence or absence.

Parameters:
motlMotlSource

Component particle list. Loaded with cryomotl.Motl.load(). None raises ValueError.

namestr, default 'basket'

Key under which the component is stored in self._components.

assign_affiliationbool, default False

When True, run NN-based affiliation against the reference ring within each tomo_id, overwriting any existing values. When False, the affiliation column is left untouched.

assign_su_orderbool, default False

When True, run chain ordering, overwriting any existing values. Mutually exclusive with copy_su_order.

copy_su_orderbool, default False

When True, copy order_column from the nearest reference-ring particle within the same (tomo_id, affiliation_column) pore. Requires affiliation to be set (either pre-existing or via assign_affiliation=True). Mutually exclusive with assign_su_order.

reference_ringint, default 0

0-based ring index used as the affiliation/copy reference. For NPC: 0=CR, 1=IR, 2=NR, 3=gp210. A warning is printed (once per call) when the default (0) is used and affiliation NN or copy is performed.

affiliation_max_distancefloat or None

Maximum voxel distance for affiliation NN assignment. Required when assign_affiliation=True.

order_max_distancefloat or None

Absolute chain-link cutoff (voxels) for ordering or copying. Required when assign_su_order=True or copy_su_order=True.

assign_subunit_order(ref_direction=None)#

Assign 1-based cyclic subunit indices into self.order_column.

For every ring group (as defined by _ring_group_columns): fits the angular grid to all particles and assigns each one the nearest grid position (1 … n), choosing the chirality that minimises duplicates then residual. Writes results into self.motl.df[self.order_column].

Parameters:
ref_directionndarray, optional

Global 3-vector that defines azimuth zero after projection onto each ring plane. When supplied, all rings use the same external reference so cross-ring index 1 points in a consistent direction. When None (default) the first particle in each ring group sets the reference.

Returns:
dict

Mapping of group key tuples (tomo_id, affiliation_id) to the per-object maximum grid residual in degrees.

Notes

Object centres are determined by _compute_object_center(), which respects self.center_method. Modifies self.motl.df in place. The @gui_exposed(returns="none") decorator means the GUI ignores the return value; Python callers can use it to filter distorted rings.

property central_angle#

Angle between adjacent subunits (360 / n degrees).

central_angles(gaps='holey', inward_axis=None, include_flagged=False, flag_column='geom5')#

Compute central angles between neighbouring subunit pairs.

For each ring group two methods are evaluated for every consecutive pair (query particle -> nearest in ring order):

  • Positional: |atan2(n * (u x v), u * v)| where u = S1 - C and v = S2 - C, with C the fitted ring centre and n the ring normal from SVD of displacement vectors.

  • Orientational: same formula applied to the inward-axis vectors obtained by rotating inward_axis by each subunit’s stored orientation. Does not depend on the accuracy of the circle fit for the centre position.

All reported angle columns are non-negative (absolute values). The underlying signed values are retained in angle_pos_signed and angle_ori_signed for callers that need directionality (e.g. to detect CW/CCW consistency across a ring).

Call assign_subunit_order() before this method to populate self.order_column.

Parameters:
gaps{‘holey’, ‘full’}, default=’holey’

Pairing mode. 'holey' pairs every present subunit with the next present one (by sorted position); 'full' skips pairs whose target index (order + 1 mod n) is absent from the ring.

inward_axisarray-like of shape (3,) or None, default=None

Local axis pointing toward the ring centre in the particle frame. None uses [-1, 0, 0], the convention in _compute_object_center() and assign_subunit_order().

include_flaggedbool, default=False

When False, pores with flag_column == 1.0 are excluded from the output.

flag_columnstr, default 'geom5'

Column checked for the per-pore ambiguity flag.

Returns:
pandas.DataFrame

One row per neighbouring pair with columns:

tomo_id, object_id, qp_subtomo_id, nn_subtomo_id, qp_idx, nn_idx, idx_diff, angle_pos, angle_ori, angle_diff, angle_pos_per_pos, angle_ori_per_pos, dev_pos, dev_ori, angle_pos_signed, angle_ori_signed

All angle_* and dev_* columns are in degrees. angle_pos and angle_ori are non-negative (magnitudes). angle_diff = |angle_pos - angle_ori| (unsigned disagreement). dev_pos and dev_ori are signed deviations from the ideal idx_diff x central_angle (positive = wider than ideal, negative = narrower). angle_pos_signed / angle_ori_signed retain the sign from atan2; their sign depends on the SVD normal orientation and is arbitrary between rings – use only for within-ring consistency checks, never for cross-ring comparison.

Raises:
ValueError

If affiliation_column is absent from self.motl.df.

ValueError

If order_column is absent from self.motl.df.

ValueError

If a ring group has duplicate order_column values.

check_for_duplicates()#

Return ring groups whose order_column has a within-group duplicate.

Groups by the ring-group columns (tomo_id, affiliation_column, and _ring_column when present) and collects every group where any order_column value appears more than once.

Returns:
pandas.DataFrame

One row per offending group. Columns: tomo_id, object_id, and ring (only when a ring column is part of the grouping). Empty when no duplicates are found.

circumference(*, pixel_size=None)#

Per-object circumference derived from the circumradius.

Computes 2 π × circumradius × pixel_size for each object. The circumradius is estimated via the Pratt circle fit (≥ 4 particles) or falls back to the mean particle–to–barycenter distance.

Parameters:
pixel_sizefloat or None, default None

Ångström-per-voxel scale factor. When None, uses self.pixel_size (defaulting to 1.0 when unset).

Returns:
pandas.DataFrame

Columns tomo_id, object_id, circumference.

Raises:
ValueError

If affiliation_column is absent from self.motl.df.

diameter(*, pixel_size=1.0, store_column='geom4')#

Compute the mean diameter for each object.

For even n with order_column present, opposite-subunit pairs (i, i + n//2) are matched (1-based, matching assign_subunit_order()’s convention). For odd n or when order_column is absent, the diameter is derived from the circumradius (2 × circumradius × pixel_size).

Parameters:
pixel_sizefloat, default=1.0

Ångström-per-voxel scale factor applied to all distances.

store_columnMotlColumn, default=’geom4’

Column in the returned motl_out that carries each object’s mean diameter on every row; NaN for objects with no result.

Returns:
summary_dfpandas.DataFrame

One row per (tomo_id, object_id) with columns tomo_id, object_id, mean_diameter, n_pairs. n_pairs is 0 for the circumradius fallback.

motl_outMotl

Copy of self.motl with store_column populated.

Raises:
ValueError

If affiliation_column is absent from self.motl.df.

Warns:
UserWarning

When the circumradius fallback is used (odd n, missing order_column, or even n but no pairs could be matched).

get_centers_as_motl()#

Return a Motl with one centre particle per object per tomogram.

Overrides the barycentric base implementation: uses the Pratt circle fit (for center_method='circle_fit') with automatic fallback to barycentric when the fit fails.

Returns:
Motl

One row per (tomogram, object) pair. tomo_id holds the tomogram identifier and object_id holds the affiliation value. All other columns are zero-filled.

get_component_stats(name='basket', *, pixel_size=None, degrees=True, gaps='full')#

Per-particle geometry of a registered component relative to its reference ring.

All geometry is measured against the reference ring’s centre and pore axis — never a component-fitted centre. Positions are scaled by self.pixel_size.

Parameters:
namestr, default 'basket'

Key used in add_component().

degreesbool, default True

Return angles in degrees when True, radians otherwise.

gaps'full' or 'holey', default 'full'

Passed to nnana.NearestNeighbors.ordered_pairs(). Ignored when the component has no order_column.

Returns:
pandas.DataFrame

One row per component particle. Columns: tomo_id, affiliation_column, subtomo_id, radius, axial_offset, distance_to_centre, tilt_angle. order_column, central_angle_prev, central_angle_next, and interior_angle are added when the component carries an order_column.

get_object_stats(*, pixel_size=None)#

Comprehensive per-object statistics table.

Composes occupancy(), circumference(), diameter(), and centre/radius computation into one row per (tomo_id, object_id) group.

Parameters:
pixel_sizefloat or None, default None

Ångström-per-voxel scale factor for distance columns. When None, uses self.pixel_size (defaulting to 1.0 when unset).

Returns:
pandas.DataFrame

One row per (tomo_id, object_id). Columns:

tomo_id, object_id

Group identifiers.

n_present, occupancy, missing

From occupancy().

x, y, z

Object centre coordinates (voxels).

radius

Circumradius (voxels).

circumference

From circumference().

mean_diameter, n_pairs

From diameter().

mean_angle_pos, std_angle_pos, var_angle_pos

Mean, standard deviation, and variance of the per-position positional central angle (degrees). Present and non-NaN only when assign_subunit_order() has been called and the ring has ≥ 2 subunits.

mean_angle_ori, std_angle_ori, var_angle_ori

Same statistics for the orientational central angle.

n_angle_pairs

Number of pairs contributing to the angle statistics.

Raises:
ValueError

If affiliation_column is absent from self.motl.df.

get_pore_frame_coords(*, rings=None)#

Per-subunit coordinates in each pore’s reference frame, phase-aligned on order 1.

For each ring within a pore the centre and axis are computed identically to NPC.get_subunit_stats() (via _compute_object_center() and euler_angles_to_normals()). Each subunit position is projected onto the plane perpendicular to the ring axis, giving in-plane Cartesian x, y and polar r, az coordinates.

Phase alignment is done per pore: the mean azimuth of order-1 subunits across all rings in the pore is subtracted from every azimuth, placing order 1 at ≈ 0°. x and y are recomputed from r and the aligned az.

Parameters:
ringslist of float, optional

Restrict output to these ring ids (applies only after merge, when a ring column is present). None (default) returns all rings.

Returns:
pd.DataFrame

One row per subunit with columns: tomo_id_column, affiliation_column, order_column, ring, r (pixels × pixel_size), az (degrees), x (pixels × pixel_size), y (pixels × pixel_size).

property interior_angle#

Interior angle of the regular Cn polygon ((n-2)*180 / n degrees).

renumber_subunits(object_id, mapping, *, tomo_id=None, flag_column='geom5')#

Manually correct the subunit order for one object.

Applies mapping to order_column values across all rings of the object, validates the result, then clears the ambiguity flag so the object re-enters the per-subunit analyses.

Parameters:
object_idint or float

Affiliation ID of the object to correct.

mappingdict or tuple of (old_values, new_values)

Maps old order_column value → new value. Accept a dict {old: new, …} or a pair of equal-length sequences ([old, …], [new, …]). Values not in the mapping are unchanged.

tomo_idint, float, or None, optional

Required when object_id appears in more than one tomogram.

flag_columnstr, default 'geom5'

Column cleared (set to 0.0) on success.

Raises:
ValueError

If mapping is not injective (duplicate target values).

ValueError

If object_id is ambiguous and tomo_id is not given.

ValueError

If the remapped values would create within-ring duplicates.

ring_summary(*, pixel_size=None)#

One-row-per-ring summary of radius, occupancy, pore count, and subunit count.

Radius is the mean circumradius across pores (via _circumradius_for_group()) multiplied by pixel_size. Occupancy is the mean per-pore fraction from occupancy().

Parameters:
pixel_sizefloat or None, default None

Ångström-per-voxel scale factor. None uses self.pixel_size (defaulting to 1.0).

Returns:
pandas.DataFrame

Columns: ring, radius, occupancy, n_pores, n_subunits. One row per ring; single-ring complexes return a one-row frame with ring=1.

class cryocat.analysis.assembly.ContactSite(vector, site_type='site')#

Bases: object

A single contact site in the block frame.

Parameters:
vectortuple[float, float, float]

Displacement from the block origin to the site, in voxels, expressed in the block’s local coordinate frame.

site_typestr, default=”site”

Logical type label used when defining which site pairs are allowed to contact each other in a BlockDefinition.

Raises:
ValueError

If vector does not have length 3 or if its in-plane part (x, y) is zero (which would leave the site without a defined azimuth).

site_type = 'site'#
vector#
class cryocat.analysis.assembly.DnComplex(motl, symmetry, *, affiliation_column='object_id', order_column='geom1', tomo_id_column='tomo_id', center_method='circle_fit')#

Bases: CnComplex

Dihedral Dn-symmetric structure modelled as two stacked Cn rings.

Extends CnComplex with ring-splitting, ring-aware subunit ordering (1 … n for the top ring, n+1 … 2n for the bottom ring), and inter-ring metrics (axial spacing, rotational twist).

Parameters:
motlMotlSource

Subunit particle list.

symmetrySymmetry

Dihedral fold, e.g. "D6" or 6 (integer folds are accepted and treated as Dn). Non-dihedral symmetries raise ValueError.

affiliation_columnMotlColumn, default=’object_id’

Column that identifies which object each particle belongs to.

order_columnMotlColumn, default=’geom1’

Column that assign_subunit_order() writes subunit indices into.

tomo_id_columnMotlColumn, default=’tomo_id’

Column that identifies the tomogram.

center_method{‘circle_fit’, ‘barycentric’}, default=’circle_fit’

Algorithm used by centre-computation helpers.

Raises:
ValueError

When symmetry is not a dihedral Dn group.

Notes

n_subunits equals 2 * fold (full dihedral group size). n (inherited from CnComplex via _cyclic_setup()) equals fold — the per-ring subunit count.

Ring 0 is the ring whose subunits have a higher mean axial coordinate along _split_axis (the “top” ring). Ring 1 is the “bottom” ring. After assign_subunit_order(), ring 0 subunits receive indices 1 … n and ring 1 subunits receive indices n+1 … 2n.

assign_subunit_order()#

Assign 1-based subunit indices across both rings.

Calls split_rings() when the ring column is absent from self.motl.df. Then delegates per-ring cyclic ordering to CnComplex.assign_subunit_order() (indices 1 … n within each ring). Finally offsets ring 1 indices by self.n so that the full range is 1 … 2n (ring 0 first, ring 1 second).

Modifies self.motl.df in place.

get_object_stats(*, pixel_size=None)#

Comprehensive per-object statistics for dihedral structures.

Composes per-object occupancy, ring spacing, inter-ring twist, and per-ring diameter/circumference (averaged to one row per object).

Parameters:
pixel_sizefloat or None, default None

Ångström-per-voxel scale factor for distance columns. When None, uses self.pixel_size (defaulting to 1.0 when unset).

Returns:
pandas.DataFrame

One row per (tomo_id, object_id). Columns:

tomo_id, object_id

Group identifiers.

n_present, occupancy, missing

From occupancy().

x, y, z

Object barycentre (voxels).

radius

Mean circumradius across both rings (voxels).

ring_spacing

Axial distance between rings (scaled by pixel_size).

inter_ring_twist

Rotational twist between rings (degrees).

circumference

Mean per-ring circumference (scaled by pixel_size).

mean_diameter, n_pairs

Mean per-ring diameter and total pair count.

Raises:
ValueError

If affiliation_column is absent from self.motl.df.

inter_ring_twist(*, degrees=True)#

Rotational twist between the two rings for each object.

Thin wrapper over _ring_centre_twist() that passes the dihedral axis, preserving the n-fold circular-mean phase computation. Auto-splits the rings on first call if split_rings() has not been called explicitly.

Returns the signed n-fold phase difference in [−180/n, 180/n) degrees (0 = rotational alignment). An eclipsed arrangement gives 0; a perfectly staggered arrangement (offset by exactly 180°/n) lies at the fold boundary and returns −180/n.

Parameters:
degreesbool, default=True

Return twist in degrees when True, radians when False.

Returns:
pandas.DataFrame

Columns tomo_id, object_id, inter_ring_twist.

ring_spacing(*, pixel_size=None)#

Axial distance between the two rings for each object.

Thin wrapper over _ring_centre_spacing() that passes the dihedral axis, preserving the axial-projection behaviour. Auto-splits the rings on first call if split_rings() has not been called explicitly.

Parameters:
pixel_sizefloat or None, default None

Ångström-per-voxel scale factor. When None, uses self.pixel_size (defaulting to 1.0 when unset).

Returns:
pandas.DataFrame

Columns tomo_id, object_id, ring_spacing.

split_rings(*, ring_column='geom5', axis=(0.0, 0.0, 1.0))#

Partition subunits into two axial rings and label them 0 / 1.

For each object, projects every subunit’s position relative to the object barycentre along axis. Subunits with a non-negative projection (higher axial coordinate) are labelled ring 0; those with a negative projection are labelled ring 1.

The result is written into motl.df[ring_column] and self._ring_group_columns is updated to [tomo_id_column, affiliation_column, ring_column].

Parameters:
ring_columnMotlColumn, default=’geom5’

Column to write ring labels (0 or 1) into.

axisarray-like of shape (3,), default=(0, 0, 1)

Splitting axis. Need not be normalised.

Returns:
cryomotl.Motl

self.motl with ring_column populated in place.

class cryocat.analysis.assembly.IcosahedralComplex(motl, *, affiliation_column='object_id', order_column='geom1', tomo_id_column='tomo_id')#

Bases: PolyhedralComplex

Icosahedral (I, order 60) motl-analysis complex.

class cryocat.analysis.assembly.NPC(motl, symmetry=None, *, radius=None, affiliation_column='object_id', order_column='geom1', tomo_id_column='tomo_id', center_method='circle_fit', ring_column='geom3', pixel_size=1.0)#

Bases: CnComplex

NPC-specific extensions of CnComplex.

Inherits all single-ring methods (centre computation, subunit ordering, and object merging) from CnComplex. The methods below are NPC-specific: orientation unification, multi-ring assembly, and opposite-subunit diameter analysis.

Typical workflow for a multi-ring NPC:

  1. cluster_subunits_to_rings() — trace subunits into rings and merge nearby rings (call once per ring motl).

  2. unify_nn_orientations() — flip ambiguous orientations (call with ring_index to target one ring in the multi-motl list).

  3. merge() — unify object_id across rings, stamp ring column, and reconcile IR orientations. One-way irreversible.

assign_subunit_order(ref_direction=None, flag_column='geom5', ambiguity_tol=None, *, radius=None, entry_shift=None, exit_shift=None, max_distance_factor=1.65)#

Assign 1-based subunit indices per ring-object via chain tracing.

For each ring in self._ring_motls, groups particles by (tomo_id_column, affiliation_column) and calls _chain_order_for_ring_object(). The cross-ring alignment step is deferred to _finalize_subunit_order() (called by merge()).

Parameters:
ref_directionnumpy.ndarray or None

Forwarded to the pooled-azimuth fallback inside _chain_order_for_ring_object().

flag_columnstr, default 'geom5'

Column written with 1.0 for objects whose ordering used the fallback or involved multiple chains. Created in each ring motl if absent.

ambiguity_tolfloat or None

Forwarded to the fallback. Defaults to central_angle / 4.

radiusfloat or None

Circumradius of the NPC ring in voxels, used only to compute the default tangential entry/exit shift R × sin(π/n). Needed only when entry_shift/exit_shift are not both given; falls back to self.radius, and a ValueError is raised if both are None.

entry_shift3-tuple of float or None

Local-frame shift vector for the entry throwaway motl. When given together with exit_shift, overrides the auto-computed value for every object and no radius is required.

exit_shift3-tuple of float or None

Local-frame shift vector for the exit throwaway motl. When given together with entry_shift, overrides the auto-computed value for every object and no radius is required.

max_distance_factorfloat, default 1.65

Per-ring chain-link cutoff = max_distance_factor × median nearest-neighbour spacing of the ring’s subunits (voxels). Default 1.65 spans dilated-adjacent pairs up to ~1.65× the ring’s own spacing while staying below the distance to a non-adjacent subunit.

Returns:
cryomotl.Motl

Combined motl of all rings with order_column, _ring_column, and flag_column populated.

central_angles(*args, **kwargs)#

Compute central angles between neighbouring subunit pairs.

For each ring group two methods are evaluated for every consecutive pair (query particle -> nearest in ring order):

  • Positional: |atan2(n * (u x v), u * v)| where u = S1 - C and v = S2 - C, with C the fitted ring centre and n the ring normal from SVD of displacement vectors.

  • Orientational: same formula applied to the inward-axis vectors obtained by rotating inward_axis by each subunit’s stored orientation. Does not depend on the accuracy of the circle fit for the centre position.

All reported angle columns are non-negative (absolute values). The underlying signed values are retained in angle_pos_signed and angle_ori_signed for callers that need directionality (e.g. to detect CW/CCW consistency across a ring).

Call assign_subunit_order() before this method to populate self.order_column.

Parameters:
gaps{‘holey’, ‘full’}, default=’holey’

Pairing mode. 'holey' pairs every present subunit with the next present one (by sorted position); 'full' skips pairs whose target index (order + 1 mod n) is absent from the ring.

inward_axisarray-like of shape (3,) or None, default=None

Local axis pointing toward the ring centre in the particle frame. None uses [-1, 0, 0], the convention in _compute_object_center() and assign_subunit_order().

include_flaggedbool, default=False

When False, pores with flag_column == 1.0 are excluded from the output.

flag_columnstr, default 'geom5'

Column checked for the per-pore ambiguity flag.

Returns:
pandas.DataFrame

One row per neighbouring pair with columns:

tomo_id, object_id, qp_subtomo_id, nn_subtomo_id, qp_idx, nn_idx, idx_diff, angle_pos, angle_ori, angle_diff, angle_pos_per_pos, angle_ori_per_pos, dev_pos, dev_ori, angle_pos_signed, angle_ori_signed

All angle_* and dev_* columns are in degrees. angle_pos and angle_ori are non-negative (magnitudes). angle_diff = |angle_pos - angle_ori| (unsigned disagreement). dev_pos and dev_ori are signed deviations from the ideal idx_diff x central_angle (positive = wider than ideal, negative = narrower). angle_pos_signed / angle_ori_signed retain the sign from atan2; their sign depends on the SVD normal orientation and is arbitrary between rings – use only for within-ring consistency checks, never for cross-ring comparison.

Raises:
ValueError

If affiliation_column is absent from self.motl.df.

ValueError

If order_column is absent from self.motl.df.

ValueError

If a ring group has duplicate order_column values.

circumference(*args, **kwargs)#

Per-object circumference derived from the circumradius.

Computes 2 π × circumradius × pixel_size for each object. The circumradius is estimated via the Pratt circle fit (≥ 4 particles) or falls back to the mean particle–to–barycenter distance.

Parameters:
pixel_sizefloat or None, default None

Ångström-per-voxel scale factor. When None, uses self.pixel_size (defaulting to 1.0 when unset).

Returns:
pandas.DataFrame

Columns tomo_id, object_id, circumference.

Raises:
ValueError

If affiliation_column is absent from self.motl.df.

classify_pores(*, gaps='full', range_min=12.0, vacancy_ratio=1.7)#

Classify each (pore, ring) by its central-angle pattern.

Reuses get_subunit_stats() (with gaps) for all per-subunit geometry. Classification is based on three thresholds only (the three keyword args):

  • "vacancy" — ring is incomplete (n_present < full_n) and max_gap ≥ vacancy_ratio × nominal where nominal = 360 / full_n. A vacancy gap can be a milling-truncation edge, not necessarily a missing subunit.

  • "shear" — within-ring range = max_gap − min_gap ≥ range_min.

  • "regular" — otherwise.

Cross-ring coherence columns are added per pore; they are NaN for pores with only one ring present:

  • coherence_corr — mean pairwise Spearman of the per-order central_angle_next vectors across rings.

  • argmax_agreement — fraction of rings whose widest-gap order equals the most common widest-gap order across the pore’s rings.

Parameters:
gapsstr, default "full"

Passed to get_subunit_stats().

range_minfloat, default 12.0

Minimum within-ring angular range (degrees) for "shear" label.

vacancy_ratiofloat, default 1.7

Widest-gap threshold as a multiple of the nominal spacing.

Returns:
pd.DataFrame

One row per (pore, ring) with columns: tomo_id_column, affiliation_column, ring, n_present, max_gap, min_gap, range, complete, type, n_rings, coherence_corr, argmax_agreement.

clean_per_object(*args, **kwargs)#

Keep the n best rows per object and drop the rest.

For each (tomo_id_column, affiliation_column) group: sort by column, keep the top-n rows according to keep, and discard the remainder. Objects that already have at most n rows are returned unchanged.

Parameters:
columnMotlColumn

Column to sort and filter by.

keep{‘high’, ‘low’}, default=’high’

'high' retains the n rows with the largest values (e.g. scores); 'low' retains the n rows with the smallest values (e.g. cone distances).

nint, optional

Number of rows to keep per object. Defaults to self.n_subunits.

Returns:
cryomotl.Motl

A copy of self.motl with over-occupied objects trimmed.

Raises:
ValueError

If affiliation_column is absent from self.motl.df.

static cluster_subunits_to_rings(input_motl, npc_radius=None, max_trace_distance=None, min_trace_distance=0, *, mask_size=None, entry_mask_coord=None, exit_mask_coord=None, entry_mask=None, exit_mask=None, tomo_id_column='tomo_id', affiliation_column='object_id')#

Cluster NPC subunit particles into rings.

Workflow:

  1. Build (or accept) spherical entry/exit masks.

  2. Re-centre the input motl to the entry and exit sub-particle positions.

  3. Trace entry/exit pairs into chains with Chain.from_motls().

  4. Copy chain annotations onto the original motl and merge nearby subunits with merge_subunits().

Mask handling is fully in-memory: when a coord + mask_size is given, cryocat.core.cryomask.spherical_mask() is called without output_path and the returned ndarray is forwarded to cryocat.core.cryomotl.Motl.recenter_to_subparticle() (whose input_map parameter accepts ndarrays via cryocat.core.cryomap.read()). No temporary files are written.

For each of the entry and exit sides, either the mask itself or the (coord + mask_size) pair must be supplied.

Parameters:
input_motlMotlSource

Subunit particle list. A Motl, a DataFrame, or a path to a motl file – cryocat.core.cryomotl.Motl.load() normalises all three.

npc_radiusfloat

Approximate NPC ring radius (voxels) used by merge_subunits().

max_trace_distancefloat

Maximum allowed step distance during chain tracing (voxels).

min_trace_distancefloat, default=0

Minimum allowed step distance during chain tracing.

mask_sizeTripletLike, optional

Box size for the in-memory entry / exit masks. Required when entry_mask / exit_mask are not supplied; ignored otherwise.

entry_mask_coordTripletLike, optional

Centre of the entry spherical mask (voxels). Required when entry_mask is not supplied.

exit_mask_coordTripletLike, optional

Centre of the exit spherical mask (voxels). Required when exit_mask is not supplied.

entry_maskMapSource, optional

User-provided entry mask. Path or ndarray. When supplied, entry_mask_coord / mask_size are ignored on the entry side.

exit_maskMapSource, optional

User-provided exit mask. Path or ndarray. When supplied, exit_mask_coord / mask_size are ignored on the exit side.

tomo_id_columnstr, default=’tomo_id’

Column holding the tomogram identifier.

affiliation_columnstr, default=’object_id’

Column holding the ring/object identifier used during sorting and merging.

Returns:
Motl

Motl with affiliation_column identifying each ring, geom1 holding ring occupancy, and geom2 the within-ring subunit index.

Raises:
ValueError

If, for either side, neither the mask nor the (coord + mask_size) pair was supplied.

static compute_diameter(input_motl, *, pixel_size=1.0, store_column='geom4', symmetry=8, tomo_id_column='tomo_id', affiliation_column='object_id')#

Compute the mean NPC diameter per ring using opposite-subunit pairs.

Matches subunit pairs (i, i + symmetry//2) using the 1-based index stored in geom2. Objects with no matching opposite pair are omitted from the summary and receive NaN in store_column. Unlike CnComplex.diameter(), no circumradius fallback is applied.

Parameters:
input_motlMotlSource

Particle list with NPC subunits. Requires affiliation_column for ring affiliation and geom2 for the 1-based subunit order within each ring.

pixel_sizefloat, default=1.0

Ångström-per-voxel scale factor applied to all distances.

store_columnMotlColumn, default=’geom4’

Column in the returned motl that carries each ring’s mean diameter. NaN for rings with no opposite-pair matches.

symmetryint, default=8

Rotational symmetry order. Determines the pair offset symmetry // 2.

tomo_id_columnstr, default=’tomo_id’

Column holding the tomogram identifier.

affiliation_columnstr, default=’object_id’

Column holding the ring/object identifier.

Returns:
summary_dfpandas.DataFrame

One row per (tomo_id_column, affiliation_column) that produced at least one opposite-subunit pair. Columns: tomo_id_column, affiliation_column, mean_diameter, n_pairs. Empty when no ring has matching pairs.

motl_outMotl

Copy of input_motl with store_column populated; NaN for rings without pairs.

create_affiliation(*args, **kwargs)#

Cluster subunit particles into objects and write affiliation labels.

Operates on a copy of self.motl and returns it with affiliation_column populated. For every object that ends up with more than self.n_subunits particles, outlier-normal particles are removed (see normals_threshold). If an object remains over-populated after cleaning its geom5 column is set to 1.0 and a one-line notice is printed.

Subunit ordering is a separate stage-2 call to assign_subunit_order(). Occupancy and normal deviations are computed on demand via occupancy() and normal_deviations().

Parameters:
method{‘radius’, ‘tracing’}, default=’radius’

Clustering strategy.

'radius'

Optionally shift particles along their local −x axis by shift, then run a self nearest-neighbour search within radius. Connected components of the NN graph become objects. Isolated particles (no NN within radius) are kept as singleton objects with unique affiliation_column values.

'tracing'

Optionally shift particles along their local −x axis by shift, then trace chains via nnana.trace_chains() with radius as the maximum link distance. Each chain becomes one object.

shiftfloat, optional

Magnitude of the local-frame shift along −x applied before clustering (voxels). When None no recentring is performed. Typical value: approximate ring radius.

radiusfloat

For method='radius': NN search radius (voxels). For method='tracing': maximum chain-link distance (voxels). Required.

normals_thresholdfloat or None, default=None

Cone-distance cutoff (degrees) applied only to over-populated objects (size > n_subunits). None uses the adaptive rule median(cone) + 1 × std(cone) per object.

occupancy_columnMotlColumn or None, default=None

Accepted but ignored. Use occupancy() on demand.

cone_distance_columnMotlColumn or None, default=None

Accepted but ignored. Use normal_deviations() on demand.

min_occupancyint or None, default=None

When set, objects whose size (after over-population cleaning) is below this value are removed. None keeps all objects.

Returns:
cryomotl.Motl

A new motl with affiliation_column and geom5 populated. geom5 is 0.0 for all particles; 1.0 marks objects that are still over-populated after cleaning.

Raises:
ValueError

If radius is None or method is unrecognised.

flagged_pores(flag_column='geom5')#

Return a summary of pores whose flag_column is set to 1.0.

Each row of the returned DataFrame corresponds to one (tomo_id, object_id, ring) group. The order_sequence column contains the sorted list of order_column values assigned to subunits in that ring, making it easy to spot missing or duplicated slots.

Parameters:
flag_columnstr, default 'geom5'

Column written by assign_subunit_order() / _finalize_subunit_order() to mark ambiguous pores.

Returns:
pandas.DataFrame

Columns: tomo_id, object_id, ring, order_sequence. Empty when no pores are flagged.

static get_centers_as_motl(tomo_motl, *, tomo_id=None, radius=None, tomo_id_column='tomo_id', affiliation_column='object_id')#

Return one centre particle per ring using the radius-shift estimator.

Shifts each subunit by (-radius, 0, 0) along its local X axis and averages the resulting positions to estimate the NPC ring centre. Unlike the inherited CnComplex.get_centers_as_motl(), this method works correctly for any ring occupancy including a single subunit.

Parameters:
tomo_motlMotlSource

Particle list for one tomogram (or all tomograms).

tomo_idfloat, optional

Tomogram identifier stored in the output motl. Defaults to the tomo_id_column value found on each ring’s particles.

radiusfloat, default=55.0

Approximate NPC ring radius in voxels.

tomo_id_columnstr, default=’tomo_id’

Column holding the tomogram identifier.

affiliation_columnstr, default=’object_id’

Column holding the ring/object identifier.

Returns:
Motl

One row per unique affiliation_column with the estimated ring centre in x, y, z.

get_object_stats(*, pixel_size=None)#

Comprehensive per-object statistics for NPC structures.

Pre-merge direct calls dispatch per-ring via per_ring(); post-merge calls run on the full merged motl, grouping by (tomo_id, object_id, ring_column) so each ring gets its own centre, radius and circumference row.

Uses geom.barycenter and _circumradius_for_group. Circumference = 2π × circumradius × pixel_size. Spacing and twist (per adjacent ring pair, not per ring) are merged on (tomo_id, object_id).

get_subunit_orientation(*, degrees=True, reference_ring=None, include_flagged=False, flag_column='geom5')#

Per-subunit orientation difference between ring pairs.

For each subunit order k matched across rings, reports the cone (out-of-plane normal deviation) and in-plane (roll) differences between each ring pair, via geom.cone_inplane_distance().

Requires merge() to have been called.

Parameters:
degreesbool, default=True

Return angles in degrees when True, radians otherwise.

reference_ringint or None, default=None

Ring whose subunit order defines the iteration. None uses the first ring present in each pore.

include_flaggedbool, default=False

When False, pores with flag_column == 1.0 are excluded.

flag_columnstr, default 'geom5'

Column checked for the per-pore ambiguity flag.

Returns:
pandas.DataFrame

One row per (pore, subunit). Columns:

  • tomo_id, affiliation_column, order_column

  • cone_{r0}_{r1} — cone (normal) angular difference for ring pair; NaN when the subunit is absent in either ring

  • inplane_{r0}_{r1} — in-plane (roll) angular difference; NaN when the subunit is absent in either ring

get_subunit_spacing(*, pixel_size=1.0, degrees=True, reference_ring=None, include_flagged=False, flag_column='geom5')#

Per-subunit spacing and twist between adjacent rings.

For each subunit order k, reports the Euclidean distance and azimuthal twist between the corresponding subunit in each pair of adjacent rings. The per-pore summary numbers from ring_spacing() and inter_ring_twist() remain available as aggregate columns.

Ring centres are computed with _compute_object_center() (respects center_method).

Requires merge() to have been called.

Parameters:
pixel_sizefloat, default=1.0

Voxel size in physical units; distances are multiplied by this.

degreesbool, default=True

Return twist angles in degrees when True, radians otherwise.

reference_ringint or None, default=None

Ring whose subunit order defines the iteration. None uses the first ring present in each pore.

include_flaggedbool, default=False

When False, pores with flag_column == 1.0 are excluded from the output.

flag_columnstr, default 'geom5'

Column checked for the per-pore ambiguity flag.

Returns:
pandas.DataFrame

One row per (pore, subunit). Columns:

  • tomo_id, object_id, order_column (subunit index)

  • spacing_{r0}_{r1} — Euclidean distance (voxels × pixel_size) between subunit k in ring r0 and ring r1 for each ring pair, including non-adjacent pairs

  • twist_{r0}_{r1} — signed azimuthal angle (degrees) of subunit k in ring r1 relative to its position in ring r0, measured from the ring-pair axis; one column per ring pair

get_subunit_stats(*, pixel_size=None, degrees=True, include_flagged=False, flag_column='geom5', gaps='full')#

Per-subunit geometry statistics.

Computes local geometry for every particle in the NPC motl. Prev/next neighbours within each ring are found with nnana.NearestNeighbors.ordered_pairs() (circular topology), which respects the assigned subunit order rather than distance. Ring centres use _compute_object_center() (respects center_method). Z-normals use geom.euler_angles_to_normals().

Parameters:
pixel_sizefloat, default=1.0

Voxel size; distances are multiplied by this.

degreesbool, default=True

Return angles in degrees when True, radians otherwise.

include_flaggedbool, default=False

When False, pores with flag_column == 1.0 are excluded from the output.

flag_columnstr, default 'geom5'

Column checked for the per-pore ambiguity flag.

gaps'full' or 'holey', default 'full'

Passed directly to nnana.NearestNeighbors.ordered_pairs().

'full' — a subunit pairs only with the subunit whose order value is exactly order ± 1 (circular). If that slot is absent (milling-truncated edge, internal gap), there is no neighbour on that side and the corresponding central_angle_* and interior_angle are NaN.

'holey' — a subunit pairs with the previous/next present subunit by order, skipping over absent slots. Every subunit in a ring with ≥ 2 present members gets both neighbours, so no gap-NaN arises; the angle across a skipped gap reflects the actual spacing, not the ideal 45°.

Returns:
pandas.DataFrame

One row per particle. Columns:

  • tomo_id, object_id, ring_column (post-merge), order_column, subtomo_id

  • distance_to_centre — distance from particle to ring centre (voxels × pixel_size)

  • tilt_angle — angle (°) between particle z-normal and ring mean normal

  • central_angle_prev, central_angle_next — unsigned angle (°) at ring centre between this particle and each order-adjacent neighbour; NaN when the required order-adjacent slot is absent and gaps is 'full'

  • interior_angle — unsigned vertex angle (°) at this particle: the angle between the two edge vectors meeting at this vertex, in [0°, 180°], direction-independent. NaN when either neighbour is absent. Cannot distinguish a genuine reflex concave vertex (it folds to ≤180°); real pores are convex so this is intentional.

  • shear — difference between the two flanking central angles (central_angle_next − central_angle_prev); ≈0 when the subunit sits equidistant from both neighbours, large when it is displaced toward one side. NaN when either flanking angle is NaN.

inter_ring_twist(*, degrees=True)#

Rotational twist between adjacent rings about the pore axis.

For each pore and every unordered ring pair, the pore axis is estimated as the unit vector from ring r0’s centre to ring r1’s centre. Produces n*(n-1)/2 columns for n rings. The rotational twist is the signed n-fold circular-mean phase difference folded into [−π/n, π/n) ([−180/n, 180/n) degrees); 0 means the rings are rotationally aligned.

Parameters:
degreesbool, default=True

Return twist in degrees when True, radians when False.

Returns:
pandas.DataFrame

One row per (tomo_id, object_id). Columns:

tomo_id, object_id

Pore identifiers.

twist_{r0}_{r1}

Rotational twist between ring r0 and ring r1; one column per ring pair, including non-adjacent pairs. NaN when one or both rings are absent or when ring centres coincide.

Raises:
ValueError

If called before merge().

merge(*, npc_radius=None, ring_order=None, distance_threshold=40, store_sid_column=None)#

Unify affiliation across rings, stamp ring column, and reconcile IR orientations.

This is a one-way irreversible operation. After the call, self.motl holds all rings as a single merged motl, self._rings_merged is True, and self._ring_group_columns is extended to include _ring_column (geom3).

Steps:

  1. Re-order rings according to ring_order (default: input order).

  2. Stamp ring index into _ring_column for every particle.

  3. Renumber affiliation_column sequentially so every ring starts fresh, then use nearest-neighbour pore matching (nnana.find_nn_indices) to unify affiliation_column across rings that belong to the same pore.

  4. Reconcile IR orientations via _align_to_reference(): flip any IR particle (ring index 1) whose z-normal is more than 90° away from the per-pore mean CR z-normal.

  5. Concatenate and store.

  6. If store_sid_column is given, save the original subtomo_id values into that column and call renumber_particles() so the merged motl has globally unique particle IDs.

Parameters:
npc_radiusfloat

Ring radius in voxels. Forwarded to get_centers_as_motl().

ring_orderlist of int or None

Permutation of range(len(rings)). Index 0 = CR, index 1 = IR, index 2 = NR. Defaults to input order.

distance_thresholdfloat, default=40

Maximum pore-centre distance (voxels) for two rings to be considered the same pore.

store_sid_columnstr or None, default=None

When set, the original subtomo_id values (per-ring particle IDs) are stored in this column before renumbering. Useful for tracing each merged particle back to its source.

static merge_rings(input_motls, npc_radius=None, distance_threshold=40, *, tomo_id_column='tomo_id', affiliation_column='object_id')#

Merge corresponding rings across multiple ring-motls.

Assigns sequential affiliation_column values across all motls, then for every pair of motls finds rings (by their estimated centres) that are closer than distance_threshold and unifies their affiliation_column entries so matched rings share the same identifier.

Parameters:
input_motlslist of MotlSource

At least two ring-motls to merge.

npc_radiusfloat

Ring radius in voxels, forwarded to get_centers_as_motl().

distance_thresholdfloat, default=40

Maximum centre-to-centre distance (voxels) for two rings from different motls to be considered the same NPC.

tomo_id_columnMotlColumn, default=”tomo_id”

Column used to group particles by tomogram.

affiliation_columnMotlColumn, default=”object_id”

Column used to identify which NPC pore a particle belongs to.

Returns:
list of Motl

The input motls with updated affiliation_column values so that matched rings share the same identifier.

Raises:
UserWarning

When input_motls is not a list or contains fewer than two items.

merge_subunits(*args, **kwargs)#

Merge near-duplicate objects whose centres are within radius.

For each tomogram:

  1. Compute per-object barycentric centres.

  2. Find object-centre pairs within radius using nnana.get_nn_within_distance().

  3. Re-assign affiliation_column of near objects to the first encountered partner.

  4. Recount occupancy into geom1 and recompute subunit order for all objects via assign_subunit_order().

Parameters:
radiusfloat, default=55

Distance threshold in voxels.

Notes

Modifies self.motl.df in place.

occupancy()#

Per-object (per-ring when merged) subunit occupancy.

Overrides SymmetricComplex.occupancy() to group by _ring_group_columns so that, after merge(), each (tomo_id, object_id, ring) group is reported separately.

Pre-merge direct calls dispatch via per_ring(); calls from inside per_ring() (_per_ring_active=True) run this body directly on the current self.motl.

per_ring(method_name, *args, **kwargs)#

Call method_name on each held ring motl and return aggregated results.

When rings are already merged, delegates directly to the instance method. Otherwise, walks the MRO (skipping NPC itself) to find the base-class implementation, temporarily swaps self.motl to each ring, calls the base implementation, and tags the result with a ring index in _ring_column (geom3). _per_ring_active is set to True during the loop so that routing guards on NPC overrides (e.g. occupancy()) skip re-dispatch when called from inside the base method. DataFrame results are concatenated; Motl results are returned as a list.

ring_ellipticity(*, pixel_size=None, flag_column='geom5', include_flagged=False)#

Ellipticity of each ring from a least-squares 2-fold radial fit.

For each ring, particles are projected onto the ring plane (defined by _compute_object_center() and the mean z-normal) and the in-plane radii fitted to r = c0 + c2·cos(2θ) + s2·sin(2θ) by least squares. The 2-fold amplitude A = hypot(c2, s2) characterises the ellipse.

Rings with fewer than 5 particles or azimuthal coverage below 180° are flagged (partial_arc == 1) and metric columns set to NaN, as partial arcs produce biased 2-fold estimates.

Parameters:
pixel_sizefloat, default=1.0

Voxel size; major and minor are multiplied by this.

flag_columnstr, default 'geom5'

Column checked for the per-pore ambiguity flag.

include_flaggedbool, default=False

When False, pores with flag_column == 1.0 are excluded.

Returns:
pandas.DataFrame

One row per ring group (tomo_id, object_id, ring_column when post-merge). Columns:

  • n_used — number of particles in the ring

  • ellipticity — fractional 2-fold amplitude A / c0 (0 = circle)

  • axis_ratio — minor/major ratio (c0 − A) / (c0 + A)

  • major, minor — semi-axes (voxels × pixel_size)

  • orientation — major-axis angle (°) in the ring plane, 0.5 · atan2(s2, c2)

  • partial_arc — 1 if coverage < 180° or n_used < 5, else 0

ring_spacing(*, pixel_size=None)#

Inter-ring spacing per NPC pore, computed from the merged motl.

For every unordered pair of rings, computes the 3-D Euclidean distance between the mean coordinates of the two ring centres for each pore and scales by pixel_size. Produces n*(n-1)/2 columns for n rings.

Parameters:
pixel_sizefloat, default=1.0

Ångström-per-voxel scale factor applied to the raw voxel distances.

Returns:
pandas.DataFrame

One row per (tomo_id, object_id). Columns:

tomo_id, object_id

Pore identifiers.

spacing_{r0}_{r1}

Distance (Å when pixel_size given) between the centre of ring r0 and ring r1; one column per ring pair, including non-adjacent pairs. NaN when one or both ring centres are absent for a pore.

Raises:
ValueError

If called before merge().

split_by_ring(result)#

Split result into per-ring sub-DataFrames using _ring_column.

Returns a list of length len(_ring_motls), one element per ring. Each element contains only rows whose _ring_column value equals that ring’s index. If result has no ring column (e.g. single-ring NPC or a pre-dispatch call), the whole DataFrame is returned as a one-element list.

This is the library-side of the send-to-motl splitting contract: the caller maps each returned DataFrame to the corresponding source motl pool entry. No new matching logic is needed because the ring column records the origin of every particle.

step_statistics(*args, **kwargs)#
unify_nn_orientations(dist_threshold=10000, ring_index=0, *, reference=None)#

Flip orientations so that neighbouring subunits point consistently.

When reference is supplied, delegates to _align_to_reference(): any particle in self.motl whose z-normal is more than 90° away from the per-object mean z-normal of reference is flipped 180° about z. dist_threshold and ring_index are ignored in this branch.

Without reference, traces particles into chains via nnana.trace_chains(), then walks each chain and applies a 180° rotation (srot.from_euler("zxz", [0,180,0])) whenever the cone angle (geom.cone_distance) between successive subunits exceeds 90°. Updates self.motl in place, and when multiple ring motls are held, also updates self._ring_motls[ring_index].

Parameters:
dist_thresholdfloat, default=10000

Maximum nearest-neighbour distance for tracing (voxels). Ignored when reference is given.

ring_indexint, default=0

Which ring motl to operate on when the NPC holds multiple unmerged rings. Must be 0 after merge() is called (raises otherwise). Ignored when reference is given.

referencecryomotl.Motl or None, default=None

Reference motl that defines the expected orientation direction. When provided, each particle is compared against the per-object mean z-normal of reference rather than against its chain neighbours.

Raises:
ValueError

If called post-merge with ring_index != 0 (chain-tracing path only; the reference path never raises this). After merge() there is a single combined motl; a non-zero ring_index has no target and would silently operate on the wrong data.

update_ring_motls(flag_column='geom5')#

Propagate corrected subunit order from the merged motl back to _ring_motls.

After manual corrections to the merged motl (e.g. fixing flagged pores with renumber_subunits()), call this to sync the updated order_column values back into each source ring motl so they can be written out. Uses the original subtomo_id values stored in _store_sid_column during merge() to match each merged row to its source ring motl row.

Parameters:
flag_columnstr, default 'geom5'

When present in the merged motl, its values are also propagated back so the ring motls record which particles were flagged.

Raises:
ValueError

If merge() has not been called, or if store_sid_column was None during merge() (making back-tracking impossible).

write_ring_motls(path_prefix, suffix='.em')#

Write each ring motl to {path_prefix}{ring_index}{suffix}.

Parameters:
path_prefixstr

Directory + condition string, e.g. "/data/control_". Produces /data/control_0.em, /data/control_1.em, …

suffixstr, default '.em'

File extension.

class cryocat.analysis.assembly.OctahedralComplex(motl, *, affiliation_column='object_id', order_column='geom1', tomo_id_column='tomo_id')#

Bases: PolyhedralComplex

Octahedral (O, order 24) motl-analysis complex.

class cryocat.analysis.assembly.ParametricSurface(quadrics, column_name='object_id')#

Bases: object

Wrapper around QuadricsM for the ellipsoid particle-assignment workflow.

Mirrors the static-method interface of the old PleomorphicSurface ellipsoid methods but as a proper class with instance state.

Parameters:
quadricsQuadricsM

Already-constructed container of analytic surfaces.

column_nameMotlColumn, default=’object_id’

Column name used as the surface-object identifier (one fitted surface per unique (tomo_id, column_name) group). Per the column-naming convention this is cryocat._types.MotlColumn.

assign_affiliation_distance_based(input_motl, output_path=None, unassigned_value=None)#

Assign each particle to the nearest surface centre.

Parameters:
input_motlMotlSource

Particles to assign.

output_pathPathOrStr, optional

Path to save the result.

unassigned_valuefloat, optional

When provided, only particles whose current self.column_name value equals this value are re-assigned; the rest are kept unchanged.

Returns:
Motl

Motl with self.column_name updated.

assign_affiliation_intersection_based(input_motl, output_path=None, keep_unassigned=True)#

Assign each particle to the surface it points toward (ray casting).

A ray is cast along the negated particle normal. The particle is labelled with the identifier of the surface whose intersection is closest along that ray. Particles that lie inside a surface or have no valid intersection receive -1.

Parameters:
input_motlMotlSource

Particles to assign. Euler angles are used to derive normals.

output_pathPathOrStr, optional

Path to save the result.

keep_unassignedbool, default=True

When False, particles whose self.column_name is -1 are removed.

Returns:
Motl

Motl with self.column_name updated.

static assign_affiliation_mask_based(input_motl, object_motl, tomo_dim, shell_size, column_name='object_id', output_path=None, radius_offset=0.0, motl_radius_id='geom5')#

Assign each particle to a surface object using a spherical-shell mask.

Parameters:
input_motlMotlSource

Particles to assign.

object_motlMotlSource

Surface-object positions (one row per object).

tomo_dimTomoDimensions

Tomogram dimensions; normalized via ioutils.dimensions_load().

shell_sizeint

Thickness of the spherical shell mask in voxels.

column_nameMotlColumn, default=’object_id’

Identifier column on input_motl that receives the assigned object id.

output_pathPathOrStr, optional

Path to save the result.

radius_offsetfloat, default=0.0

Constant added to object_motl[motl_radius_id] to pad the shell.

motl_radius_idMotlColumn, default=’geom5’

Column in object_motl holding each object’s radius.

Returns:
Motl
clean_by_normals(input_motl, compute_normals=True, normals_id='geom4', threshold=None, output_path=None)#

Remove particles whose orientation deviates too far from the surface normal.

Parameters:
input_motlMotlSource

Source particles.

compute_normalsbool, default=True

Recompute the angle-to-normal column before filtering.

normals_idMotlColumn, default=’geom4’

Column holding the angle-to-normal values.

thresholdfloat, optional

Maximum allowed angle (degrees). Defaults to one standard deviation.

output_pathPathOrStr, optional

Path to save the result.

Returns:
Motl
clean_by_radius(input_motl, threshold=None, output_path=None)#

Remove particles that lie too far from the mean ellipsoid radius.

Parameters:
input_motlMotlSource

Source particles.

thresholdfloat, optional

Half-width of the allowed distance band. Defaults to one standard deviation.

output_pathPathOrStr, optional

Path to save the result.

Returns:
Motl
compute_intersection(input_motl)#

Compute ray-ellipsoid intersection distances for each particle.

For every particle a ray is cast along -euler_angles_to_normals and the two intersection distances with the assigned ellipsoid are returned.

Parameters:
input_motlMotlSource

Particles grouped by self.column_name.

Returns:
pandas.DataFrame

Columns: subtomo_id, <self.column_name>, d1, d2.

compute_normals_angle(input_motl, store_column_name='geom4', output_path=None)#

Compute the angle between each particle’s orientation and the ellipsoid radial normal.

The radial normal is the vector from the fitted ellipsoid centre to the particle. The stored value is the angle (degrees) between that vector and the particle’s orientation normal.

Parameters:
input_motlMotlSource

Particles grouped by self.column_name.

store_column_nameMotlColumn, default=’geom4’

Column that receives the angle values.

output_pathPathOrStr, optional

Path to save the result.

Returns:
Motl

Input motl with store_column_name populated.

compute_point_surface_distance(input_motl, output_path=None, store_column_name='geom4')#

Compute the shortest distance from each particle to its assigned surface.

Parameters:
input_motlMotlSource

Particles with self.column_name already assigned.

output_pathPathOrStr, optional

Path to save the result.

store_column_nameMotlColumn, default=’geom4’

Column that receives the distance values.

Returns:
Motl

Input motl with store_column_name populated.

static create_spherical_oversampling(input_motl, motl_radius_id, sampling_distance, sampling_angle=360.0, output_path=None)#

Generate oversampled particles on a sphere around each input particle.

Parameters:
input_motlMotlSource

Source particles (one sphere per row).

motl_radius_idMotlColumn

Column holding the sphere radius for each particle.

sampling_distancefloat

Angular sampling step forwarded to geom.sample_cone().

sampling_anglefloat, default=360.0

Half-opening angle of the sampling cone. 360 samples the full sphere.

output_pathPathOrStr, optional

Path to save the result.

Returns:
Motl
classmethod from_csv(path, surface_type='ellipsoid', column_name='object_id')#

Load analytic surface parameters from a CSV file.

Parameters:
pathPathOrStr

Path to a CSV produced by write_out().

surface_type{‘ellipsoid’}, default=’ellipsoid’

Quadric type to materialise.

column_nameMotlColumn, default=’object_id’

Surface-object identifier column the file is keyed on.

Returns:
ParametricSurface
classmethod from_motl(input_motl, surface_type='ellipsoid', column_name='object_id')#

Fit analytic surfaces to particle groups and return a ParametricSurface.

Parameters:
input_motlMotlSource

Input particle list. One surface is fitted per unique (tomo_id, column_name) group.

surface_type{‘ellipsoid’}, default=’ellipsoid’

Quadric type; currently only ellipsoid is supported.

column_nameMotlColumn, default=’object_id’

Column used to group particles.

Returns:
ParametricSurface
write_out(output_path)#

Write the surface parameter table to output_path as CSV.

Parameters:
output_pathPathOrStr

Destination CSV file path.

class cryocat.analysis.assembly.PleomorphicSurface(surface=None, *, block_layer=None)#

Bases: object

Pleomorphic lattice assembly: an envelope layer and/or a block layer.

The two layers are independent.

  • Envelope layer (_surface): a Mesh or OrientedPointCloud representing the outer surface. Access it via surface; check its presence with has_envelope.

  • Block layer (blocks): a Motl together with a BlockDefinition (or dict thereof) describing the building-block geometry. Access it via blocks; check its presence with has_blocks.

Either layer can be absent, but at least one is required.

Blocks-only example (no envelope, using the plain-parameter builder):

ps = PleomorphicSurface.from_blocks(
    motl, symmetry="C3", site_length=9.0, site_elevation=-10.0
)
ps.connect(max_distance=4.0)
mesh = ps.envelope_from_faces()
ps.surface = mesh          # attach the envelope for polarity queries

For the GUI entry point see from_blocks(). For direct construction see __init__().

adopt_recovered(motl, diagnostics=None, *, min_neighbours=None, max_angular_rms=None, max_normal_dev=None)#

Adopt recovered particles into the original lattice.

Merges the accepted subset of motl (from recover_orientations()) into blocks, tags them with origin 1 in _origin_column, and rebuilds the contact graph so that scope='original' reflects the corrected baseline. Existing blocks retain origin 1.

Any inferred blocks (_added, origin 2) are cleared because they were computed from the pre-adoption baseline. "Inference state cleared." is printed if they existed.

Parameters:
motlcryomotl.Motl

Output of recover_orientations().

diagnosticspd.DataFrame or None, default None

Diagnostics DataFrame returned by recover_orientations(). When provided, min_neighbours, max_angular_rms, and max_normal_dev are applied against its n_neighbours, angular_rms, and normal_spread columns respectively. When None, those filters are not applied.

min_neighboursint or None, default None

Keep rows whose n_neighbours (from diagnostics) is at least this value. Requires diagnostics.

max_angular_rmsfloat or None, default None

Keep rows whose angular_rms (from diagnostics) is at most this value. Requires diagnostics.

max_normal_devfloat or None, default None

Keep rows whose normal_spread (from diagnostics) is at most this value. Requires diagnostics.

Returns:
cryomotl.Motl

The adopted subset, tagged with origin 1 in _origin_column.

Raises:
ValueError

If connect() has not been called, or if any filter parameter is given without diagnostics.

annotate_with_envelope(tomo_id=None)#

Annotate blocks with their distance to the envelope surface.

Requires both an envelope surface (self.surface) and a block layer (self.blocks). The envelope is never built internally; the caller is responsible for attaching one, e.g.:

ps_clean = PleomorphicSurface(blocks=motl, block_definition=bd, ...)
ps_clean.connect(max_distance=4.0)
mesh = ps_clean.envelope_from_faces()
ps = PleomorphicSurface(mesh, blocks=motl, block_definition=bd, ...)
annot = ps.annotate_with_envelope()

Do not pass the envelope back to the same instance that produced it; envelope_from_faces uses block centres as vertices, so the resulting mesh is self-referentially tied to those same blocks and cannot give a meaningful distance query.

Parameters:
tomo_idfloat or None, default=None

Restrict blocks to this tomogram. The surface itself is never filtered; the full envelope is queried for the selected blocks.

Returns:
pandas.DataFrame

Columns: tomo_id, block_id, envelope_distance, closest_x, closest_y, closest_z, primitive_id, normal_angle, mean_curvature, gaussian_curvature.

Raises:
ValueError

If no block layer or no envelope surface is present.

apply_vertex_mask(mask, inplace=False)#

Return a surface containing only vertices where mask is True.

Pass one of the boolean masks returned by separate_surfaces():

  • For surface_type='closed': pass inner_mask or outer_mask.

  • For surface_type='planar': pass surface1_mask or surface2_mask.

Example:

inner_mask, outer_mask = ps.separate_surfaces(surface_type='closed')
inner = ps.apply_vertex_mask(inner_mask)

s1_mask, s2_mask = ps.separate_surfaces(surface_type='planar')
half1 = ps.apply_vertex_mask(s1_mask)
Parameters:
masknp.ndarray

Boolean array of shape (N,) aligned with self.surface.vertices.

inplacebool, default=False

If True, modify this instance in place and return None. If False, return a new PleomorphicSurface.

Returns:
PleomorphicSurface or None

New instance if inplace=False, else None.

property block_pixel_size#

Ångström-per-voxel scale factor for block-layer computations.

Block-layer methods (contact distances, face centroids, site annotations) multiply voxel-space coordinates by this value. It never overwrites the mesh’s own scale; surface measurements use self.surface.pixel_size directly.

browse_objects(*, plot_faces=True, scope='original')#

Browse objects interactively in a Jupyter notebook.

Shows one object at a time using plot_object(), with Previous / Next buttons and a dropdown to jump directly to an object. Requires ipywidgets; if it is not installed a clear ImportError is raised with an install hint.

Parameters:
plot_facesbool, default=True

Passed to plot_object().

scope{'original', 'inferred', 'all'}, default=``’original’``

Passed to plot_object().

check_object_grouping(scope='original')#

Diagnose whether each affiliation holds one physical assembly.

An affiliation whose contact graph is disconnected contains multiple independent pieces (n_components > 1). An affiliation whose blocks have contacts reaching into another affiliation is physically joined to it (n_cross_affiliation_contacts > 0). Use regroup_objects() to act on these findings.

Parameters:
scope{'original', 'inferred', 'all'}, default=``’original’``

'original' uses the real-block contact graph. 'all' uses the combined lattice + inferred block graph. 'inferred' checks grouping for inferred faces only.

Returns:
pandas.DataFrame

One row per (tomo_id, affiliation). Columns:

  • tomo_id — tomogram identifier

  • object_id (or the configured affiliation column) — affiliation value

  • n_blocks — total blocks in this affiliation

  • n_contacts — undirected contacts (matched half-edges / 2) within the affiliation

  • n_components — connected components of the within-affiliation contact graph (> 1 means the affiliation contains multiple pieces)

  • euler_characteristic — V − E + F for this affiliation

  • component_sizes — blocks per component, sorted descending

  • n_cross_affiliation_contacts — half-edges leaving this affiliation to a block in a different affiliation of the same tomogram; each undirected cross contact contributes 1 to both affiliations involved

Raises:
ValueError

If connect() has not been called for scope.

clean_by_angle(max_angle_deg, reference_normal, signed=False)#

Remove points whose normal deviates more than max_angle_deg from a given axis.

Companion to clean_by_normals(); both delegate to the same engine (apply_normals_mask()) but this variant cleans against a caller-supplied direction rather than the mean normal.

Parameters:
max_angle_degfloat

Maximum allowed angle (degrees) between a point’s normal and reference_normal. Points exceeding this threshold are removed.

reference_normalnp.ndarray (3,)

Reference direction/axis to measure each point’s normal against.

signedbool, default=False

If False (default), antiparallel normals count as aligned. If True, a normal pointing opposite reference_normal is treated as a 180° deviation (directional cleaning).

Returns:
PleomorphicSurface

self (modified in-place).

clean_by_normals(max_angle_deg=90.0)#

Remove points whose normal deviates more than max_angle_deg from the mean direction.

Cleans against the mean normal of the surface. To clean against a fixed axis/direction instead, use clean_by_angle(). To split by orientation relative to a reference point (e.g. centroid), use separate_surfaces() with surface_type='closed'.

Parameters:
max_angle_degfloat, default=90.0

Maximum allowed angle (degrees) between a point’s normal and the mean normal of the surface. Points exceeding this threshold are removed.

Returns:
PleomorphicSurface

self (modified in-place).

compute_normals(**kwargs)#

Delegates to Mesh.compute_normals() or OrientedPointCloud.compute_normals().

Parameters:
**kwargs

Mesh: no keywords are used; unknown keys are silently consumed.

OrientedPointCloud:

  • knn : int, default=30 — neighbors for normal estimation.

  • orient_normals : bool, default=True — orient normals consistently.

  • tangent_plane_knn : int, default=50 — neighbors for normal orientation.

  • inplace : bool, default=True — update in place; if False, wrap and return a new PleomorphicSurface.

Returns:
PleomorphicSurface

self when the delegate updates in place; a new wrapper when the delegate returns a copy (point cloud with inplace=False).

connect(max_distance, *, scope='original')#

Build and cache the contact graph from the block layer.

Matches sites across blocks within max_distance voxels, assigns each block to a connected-component assembly, and computes faces as the minimum cycle basis of the block contact graph (weighted by centre-to-centre distance).

After calling this method the following getters become available: get_contact_stats(), get_face_stats(), get_block_stats(), get_assembly_stats().

Parameters:
max_distancefloat

Upper bound on site tip-to-tip distance (voxels). Two contact points bond when they are within this distance of each other.

scopestr, default=``’original’``

'original' connects using blocks only and stores the result in the instance’s primary contact-graph slots (_site_table, _faces). 'all' is accepted for backward compatibility but has no lasting effect; the combined real + added graph is recomputed on demand whenever a scope='all' getter is called.

Raises:
ValueError

If no block layer is present, or scope is not 'original' or 'all'.

Notes

Matching algorithm:

  1. Two contact points from different blocks are candidates when their tip-to-tip distance is ≤ max_distance.

  2. Candidates are sorted by increasing tip distance and committed greedily: a pair is accepted only when both legs are still free. Each leg therefore pairs with at most one partner.

convex_hull()#

Return the convex hull as a PleomorphicSurface wrapping a Mesh.

Works for both Mesh and OrientedPointCloud backing surfaces. The per-hull statistics (volume, surface area, etc.) returned by Open3D are discarded; call Mesh.convex_hull() directly if you need them.

Returns:
PleomorphicSurface

New wrapper whose backing surface is a Mesh of the convex hull.

crop(bbox, inplace=False)#

Delegate to Mesh.crop() / OrientedPointCloud.crop().

Parameters:
bboxopen3d.geometry.AxisAlignedBoundingBox or dict

Bounding box for cropping. When a dict, must have 'min_bound' and 'max_bound' keys.

inplacebool, default=False

If True, modify in place and return None. If False, return a new wrapper.

Returns:
PleomorphicSurface or None

Wrapped surface when inplace=False; None when inplace=True.

distance_to_pointcloud(target, method='nn_unoriented', max_distance=None, ray_length=None, reverse_normals=False, bidirectional=False, one_hit_per_target=False, knn_radius=10.0, return_stats=True)#

Compute distance from this surface to another point cloud surface. - If source is Mesh: always uses raycasting - If source is OrientedPointCloud search nearest neighbours (unoriented or along normals)

Parameters:
targetPleomorphicSurface or OrientedPointCloud

Target surface. Wrapped targets are unwrapped internally; the concrete target must be an OrientedPointCloud.

method{‘nn_unoriented’, ‘nn’, ‘nearest’, ‘nn_oriented’, ‘ray’}, default=’nn_unoriented’

Distance computation method (only used if source is OrientedPointCloud). 'nn_unoriented' (aliases 'nn', 'nearest') — KDTree nearest neighbour. 'nn_oriented' (alias 'ray') — ray cast along normals.

max_distancefloat, optional

Maximum distance threshold. Points beyond this distance are excluded.

ray_lengthfloat, optional

For raycasting: maximum ray length. If None, uses infinite rays. If max_distance is set and ray_length is None, ray_length = max_distance.

reverse_normalsbool, default=False

For normal: if True, cast rays opposite to normal direction

bidirectionalbool, default=False

For Mesh sources, cast along both normal directions and keep the closer hit.

one_hit_per_targetbool, default=False

For mesh sources: keep only the closest mesh vertex per target particle (deduplicate distance_to_pointcloud hits).

knn_radiusfloat, default=10.0

For point cloud search along normals: search radius for finding points along ray trajectory

return_statsbool, default=True

If True, return stats dictionary. If False, return only distances array.

Returns:
dict or np.ndarray
If return_details=True, returns dictionary with:

‘distances’: np.ndarray (N,) - distance for each source point ‘closest_points’: np.ndarray (N, 3) - coordinates of closest/hit points ‘closest_indices’: np.ndarray (N,) - indices in target point cloud (-1 if no hit) ‘hit_mask’: np.ndarray (N,) - boolean mask of successful matches ‘closest_normals’: np.ndarray (N, 3) - normals at closest points (if available) ‘stats’: dict with min, max, mean, median, std distance statistics

If return_stats=False, returns only distances array (N,)

distance_to_points(target, compute_occupancy=True, compute_signed=False, return_closest_points=False)#

Compute distance from a point to a Mesh or an OrientedPointCloud surface.

Parameters:
targetnp.ndarray

Query points as (N, 3) array

compute_occupancybool, default=True

Compute occupancy (inside/outside). Only for Mesh.

compute_signedbool, default=False

Compute signed distance instead of unsigned. Only for Mesh. If True, compute_occupancy is automatically enabled.

return_closest_pointsbool, default=False

Return closest surface points and triangle/point IDs

Returns:
dict

Dictionary containing: - ‘distances’: unsigned or signed distances for each point - ‘distance_type’: ‘signed’ or ‘unsigned’ - ‘n_total’: total number of query points

If compute_occupancy=True and Mesh: - ‘occupancy’: binary array (1=inside, 0=outside) - ‘inside_mask’: boolean mask for inside points - ‘outside_mask’: boolean mask for outside points - ‘n_inside’: number of points inside - ‘n_outside’: number of points outside

If return_closest_points=True: - ‘closest_points’: closest points on surface (N, 3) - ‘primitive_ids’: triangle IDs (Mesh) or point IDs (PointCloud) (N,) - ‘closest_distances’: distances to closest points (same as ‘distances’ for unsigned)

Raises:
TypeError

If trying to compute occupancy/signed distance for non-Mesh surface

envelope_from_blocks(radius=None, assembly_id=None, tomo_id=None, min_blocks=4, scope='original')#

Build one triangulated mesh envelope per object using block positions.

A ball-pivoting reconstruction is applied to the block-centre point cloud of each object. Block orientations are used directly as normals; no normal re-estimation is performed. When radius is None, the pivot radius is estimated as the mean nearest-neighbour distance among the blocks of each object (via from_ball_pivoting()).

Objects with fewer than min_blocks blocks, objects whose meshing raises, and meshes with fewer than 4 faces are skipped with a warning. A ValueError is raised only when no object could be meshed.

The result is also stored as _envelopes.

Parameters:
radiusfloat or None, default=None

Ball-pivoting radius. None auto-estimates per object.

assembly_idint or None, default=None

Restrict to this assembly. None processes every assembly.

tomo_idfloat or None, default=None

Restrict to this tomogram. None processes every tomogram.

min_blocksint, default=4

Minimum number of blocks to attempt meshing an object.

scope{'original', 'inferred', 'all'}, default=``’original’``

Block pool to use (see _blocks_for_scope()).

Returns:
dict of {(tomo_id, assembly_id): Mesh}

One mesh per object.

Raises:
ImportError

If open3d is not available (required by from_ball_pivoting()).

ValueError

If no blocks remain after filtering, or no object could be meshed.

envelope_from_faces(assembly_id=None, tomo_id=None, scope='original')#

Build one triangulated mesh envelope per (tomo_id, assembly_id) object.

Each object’s mesh is built from its face polygon network: block positions become vertices, face centroids become hub vertices, and each face is fan-triangulated from its centroid. Normals are computed on the returned meshes.

The result is also stored as _envelopes so that get_shape_stats() can use it without rebuilding.

Parameters:
assembly_idint or None, default=None

Restrict to this assembly. None processes every assembly.

tomo_idfloat or None, default=None

Restrict to this tomogram. None processes every tomogram.

scope{'original', 'inferred', 'all'}, default=``’original’``

Block pool to use (see _blocks_for_scope()).

Returns:
dict of {(tomo_id, assembly_id): Mesh}

One mesh per object.

Raises:
ValueError

If connect() has not been called, or no faces remain after filtering.

property envelopes#

Per-object mesh envelopes keyed by (tomo_id, assembly_id), or None.

extract_region(indices, element='triangles', preserve_curvatures=True)#

Extract an indexed subregion from the wrapped surface.

For meshes, element='triangles' extracts a triangle submesh and preserves per-vertex curvature fields by default. For point clouds, element='points' extracts selected points; element='mask' treats indices as a boolean mask.

Parameters:
indicesnp.ndarray

Integer indices of the elements to keep, or a boolean mask when element='mask'.

element{‘triangles’, ‘points’, ‘mask’}, default=’triangles’

Which surface primitive indices refers to:

  • 'triangles' (Mesh only): select by triangle index.

  • 'points' (OrientedPointCloud only): select by point index.

  • 'mask' (OrientedPointCloud only): boolean selection mask.

preserve_curvaturesbool, default=True

Mesh only — when True, per-vertex curvature fields are copied to the extracted submesh. Ignored for point clouds.

Returns:
PleomorphicSurface

New wrapper containing only the extracted elements.

Raises:
ValueError

If element is not valid for the wrapped surface’s representation (e.g. 'triangles' on an OrientedPointCloud).

TypeError

If the wrapped surface is neither Mesh nor OrientedPointCloud.

property faces#

Triangle connectivity for mesh-backed surfaces.

flip_normals(inplace=True, **kwargs)#

Delegate normal-direction flipping to the wrapped surface.

Parameters:
inplacebool, default=True

If True, modify the wrapped surface in place and return self. If False, return a new PleomorphicSurface wrapping a flipped copy.

**kwargs

Mesh only:

  • flip_faces : bool, default=True — also reverse triangle winding so that normals recomputed from faces keep the flipped orientation.

OrientedPointCloud: no extra keywords are accepted; passing any raises TypeError.

Returns:
PleomorphicSurface or None

self when inplace=True; a new wrapper when inplace=False.

classmethod from_blocks(blocks, symmetry='C3', *, site_shift, symmetry_column=None, affiliation_column='object_id', pixel_size=1.0, tomo_id_column='tomo_id')#

Build a PleomorphicSurface from a block motl without constructing a BlockDefinition directly.

This is the GUI entry point and the preferred programmatic path. Builds a BlockLayer internally and passes it to __init__(). Example:

PleomorphicSurface.from_blocks("trimers.em", symmetry="C3", site_shift=[-5, 0, 0])

For a single symmetry (symmetry_column=None), one BlockDefinition.cyclic() definition is created with symmetry and site_shift and applied to all blocks.

For mixed folds (symmetry_column set), the fold order for each block is read from that column of blocks. Every unique value must be a positive integer (floats like 3.0 are accepted; non-integer floats like 2.5 raise ValueError). One BlockDefinition.cyclic() definition is created per unique fold, using the shared site_shift.

The result has no envelope layer; use envelope_from_faces() after connect() if one is needed.

Parameters:
blocksMotlSource

Block particle list (path or Motl).

symmetrySymmetry, default=”C3”

Cyclic symmetry specifier for the single-symmetry path: an integer n or a string like "C3". Ignored when symmetry_column is set.

site_shiftTripletLike

Shift from each block origin to its first contact site, in voxels, expressed in the block’s local frame. Typically about half the centre-to-centre distance, pointing along one leg.

symmetry_columnMotlColumn or None, default=None

If set, read the fold order from this column of blocks instead of using symmetry. All unique values must be positive integers.

affiliation_columnMotlColumn, default=”object_id”

Column whose value identifies which assembly each block belongs to. Blocks with different values are never connected by connect().

pixel_sizefloat, default=1.0

Pixel size passed to PleomorphicSurface.

tomo_id_columnMotlColumn, default=”tomo_id”

Column identifying the tomogram each block belongs to.

Returns:
PleomorphicSurface

A blocks-only surface (no envelope).

Raises:
ValueError

If symmetry_column contains non-integer or less-than-1 values.

Notes

blocks is loaded with load(), which copies a Motl instance; the input is never modified. unify_polarity(), store_block_stats() and store_block_stat() change only this object’s copy (angles, or the requested columns; row order and index are kept). Use get_blocks_as_motl() to obtain the current state, e.g. unified orientations for averaging or for building a new assembly.

get_assembly_stats(scope='original')#

Return one row per (tomo_id, assembly_id).

All measurements are derived from the contact graph that connect() builds; no assumed ideal lattice is required.

Parameters:
scope{'original', 'inferred', 'all'}, default=``’original’``

'original' uses the real-block contact graph. 'all' uses the combined lattice + inferred block graph. 'inferred' counts only assemblies that have at least one inferred face.

Returns:
pandas.DataFrame

Fixed columns: tomo_id, assembly_id, n_blocks, n_contacts, n_faces, n_boundary_half_edges, n_boundary_blocks (blocks with at least one unmatched site), closed, euler_characteristic, angle_deficit_sum, majority_face_size (mode of observed face sizes, NaN on tie or no faces), majority_face_defect (sum (mode − size) / mode over all faces; positive when the assembly contains more faces smaller than the background — e.g. pentagons in a hexagonal lattice close a sphere and give +2; negative when it contains faces larger than the background). Dynamic columns: one n_faces_<m> per distinct observed face size (size) and one n_degree_<k> per distinct vertex degree found in the data.

Raises:
ValueError

If connect() has not been called for scope.

Notes

euler_characteristic (V − E + F) is the closure measurement. A connected closed shell gives 2; an open sheet gives a different value. majority_face_defect measures face-topology deviation from the observed background face size. For a sphere built from hexagons with twelve pentagons the background is six-sided, the deviation is +2. A geodesic dome of all triangles has no face-size defect (0.0); its closure is captured by euler_characteristic = 2 and twelve degree-5 vertices in the n_degree_5 column. When the face-size distribution has no unique mode (equal split between two or more sizes), both majority_face_size and majority_face_defect are NaN.

get_block_stats(scope='original')#

Return one row per block.

Parameters:
scope{'original', 'inferred', 'all'}, default=``’original’``

'original' uses the real-block contact graph. 'all' uses the combined lattice + inferred block graph and adds an origin column (1 = lattice, 2 = inferred). 'inferred' returns only blocks that participate in inferred faces, also with origin.

Returns:
pandas.DataFrame

Fixed columns: tomo_id, block_id, block_type, n_sites, degree, complete, assembly_id, face_signature, angle_sum, angle_deficit, face_sizes, n_faces. Dynamic columns (one per distinct face size present in the lattice): n_faces_<m> — how many faces of size m border this block. The set of m values is observed from the data; no size is hardcoded. When scope is 'all' or 'inferred': additionally origin.

Raises:
ValueError

If connect() has not been called for scope.

get_blocks_as_motl(scope='original')#

Return the block layer as an independent Motl.

Parameters:
scope{“original”, “inferred”, “all”}, default=”original”

"original" returns blocks only. "inferred" returns the inferred blocks from _added. "all" returns both concatenated: original IDs are preserved; inferred IDs continue from max(subtomo_id) + 1.

Returns:
cryomotl.Motl

Independent copy of the requested block set.

Raises:
ValueError

If no block layer is present, or scope="inferred"/"all" is requested but no inferred blocks exist.

get_bond_elevations(scope='original')#

Return one directed row per bonded site pair with the elevation angle.

For each bond A → B, the elevation is the signed angle between the bond vector (from A’s block centre to B’s block centre) and the plane perpendicular to A’s z-axis (A’s normal direction in world coordinates). Positive means B is on the same side as A’s normal; negative means the opposite side. The angle is not folded to the positive half-space.

Because elevation is a directed quantity (A looking at B differs from B looking at A when their orientations differ), each undirected bond produces two rows — one for each ordered pair. Rows where unit_id == partner_id are never produced.

connect() must be called before this method.

Parameters:
scope{'original', 'inferred', 'all'}, default=``’original’``
Returns:
pandas.DataFrame

Columns:

  • tomo_id — tomogram identifier

  • assembly_id

  • unit_id — subtomo_id of the source block

  • partner_id — subtomo_id of the target block

  • elevation_deg — signed elevation angle in degrees

  • face_sizes — tuple of unique face sizes bordering the unit

  • n_faces — number of unique faces bordering the unit

  • n_faces_<m> (one per observed face size) — face count by size

Raises:
ValueError

If connect() has not been called.

get_connected_triangles(triangle_id, max_hops=1)#

Return edge-connected neighboring triangle IDs for a mesh-backed surface.

Parameters:
triangle_idint

Seed triangle index.

max_hopsint, default=1

Number of edge-traversal steps from the seed. max_hops=1 returns only immediate face-neighbors; higher values expand the region.

get_contact_stats(scope='original')#

Return one row per matched half-edge (each contact appears twice).

Parameters:
scope{'original', 'inferred', 'all'}, default=``’original’``

'original' uses the real-block contact graph. 'all' uses the combined lattice + inferred block graph. 'inferred' returns only contacts that participate in inferred faces.

Returns:
pandas.DataFrame

Columns: tomo_id, block_id, site, site_type, partner_block_id, partner_site, assembly_id, face_id, n_candidates, site_distance, block_distance, bend_angle, label_offset, mismatch_x, mismatch_y, mismatch_z.

Raises:
ValueError

If connect() has not been called for scope, or scope is 'inferred' and infer_missing_blocks() has not been called.

get_curvature_directions()#

Return per-vertex principal curvature direction vectors for a mesh-backed surface.

Returns:
np.ndarray, shape (N, 3, 2)

Direction vectors at each vertex: [:, :, 0] is the first principal direction (k1), [:, :, 1] is the second (k2).

get_curvedness()#

Return per-vertex curvedness for a mesh-backed surface.

Returns:
np.ndarray, shape (N,)

Curvedness C = sqrt((k1^2 + k2^2) / 2) at each vertex, in [0, inf).

get_face_filling_summary(other)#

Return per-size face counts before and after block filling.

Parameters:
otherPleomorphicSurface

A surface produced by with_inferred_blocks() from this surface.

Returns:
pandas.DataFrame

Columns: size, n_before, n_after, n_with_inferred. One row per face size that appears in either surface.

get_face_inference()#

Per-face inference table: merged-face detection and missing block count.

Runs _face_hole_analysis() on every closed face and returns the diagnostic result. All decisions are derived from interior-angle geometry; no assumed ideal lattice is required. The face normal used for interior-angle projection is the normalised mean of the block z-vectors (each block’s orientation applied to [0, 0, 1]), matching get_face_stats().

get_face_stats() reports the measured face topology; this method reports the inferred interpretation. get_missing_block_motl() uses the same detector to produce positions and orientations for the inferred missing blocks.

Returns:
pandas.DataFrame

Columns per face:

  • tomo_id, object_id — identity, from the face record.

  • face_id — face identifier.

  • size — observed boundary length (number of half-edges).

  • is_merged — True when the median interior angle deviates from the ideal n-gon value by more than 20°.

  • n_reflex — count of reflex vertices (turn angle < −30° in the minority winding direction); the raw signal driving n_missing.

  • n_missing — inferred missing trivalent blocks (n_reflex // 3); 0 when is_merged is False or the reflex count is not a multiple of 3.

  • n_faces_recovered — faces that would be restored: 0 when n_missing == 0; 1 + 2 * n_missing otherwise (Euler: each re-inserted unit adds ΔV=+1, ΔE=+3 → ΔF=+2, Δχ=0).

  • recovered_face_size — majority observed face size for this assembly (the expected size of each recovered face); NaN when n_missing == 0 or the assembly size distribution has no unique mode.

Raises:
ValueError

If connect() has not been called.

get_face_stats(scope='original')#

Return one row per closed face.

Faces reflect connectivity only — the undirected block contact graph. A block whose orientation is inverted relative to its neighbours keeps its bonds, so it does not change the face count or topology (χ is unaffected). Use get_inverted_blocks() to identify such blocks before or after connecting.

All columns are derived directly from the contact graph; no inference is performed here. For merged-face detection and missing block counts use get_face_inference().

Parameters:
scope{'original', 'inferred', 'all'}, default=``’original’``

'original' uses the real-block contact graph. 'all' uses the combined lattice + inferred block graph and adds the per-face column n_added (count of inferred blocks per face). 'inferred' returns only faces that exist in 'all' but not in 'original', also with n_added.

Returns:
pandas.DataFrame

Columns: tomo_id, face_id, assembly_id, size, n_unique_blocks, block_ids, centroid_x, centroid_y, centroid_z, normal_x, normal_y, normal_z, n_vertices (= size), n_boundary_blocks (face blocks with at least one unmatched site), n_inferred_blocks (blocks whose origin column equals 2; 0 when no inferred blocks are present), planarity (RMS out-of-plane deviation of cycle blocks divided by the cycle radius; 0 for a perfectly flat face, larger for a cycle that wanders off the shell), centroid_depth (distance from the face centroid to the nearest block in the assembly, divided by the cycle radius; small for a face lying on the shell, large for a shortcut cycle crossing the interior). When scope is 'all' or 'inferred': additionally n_added (count of origin-2 blocks in the face).

Radius is the mean Euclidean distance of the cycle’s blocks from their centroid. The face normal is the normalised mean of the block z-vectors (each block’s orientation applied to [0, 0, 1]).

Raises:
ValueError

If connect() has not been called, or scope is 'inferred' and infer_missing_blocks() has not been called.

get_face_vertices_as_motl(size=None, *, vertex_class=None, scope='original')#

Return a Motl with one row per (block, face) pair.

Intended for face-opening subtomogram averaging: each output particle sits at a block position with its normal unchanged, but the in-plane angle rotated so the direction from the block to its face’s centroid falls on the particle X-axis. Averaging all output particles then superposes the openings at a consistent orientation.

A block bounding two qualifying faces appears twice with different in-plane angles — that is intentional, not a duplicate.

Column assignments use the BlockLayer constructor parameters (defaults in parentheses):

  • tomo_id_column (tomo_id) — tomogram identifier

  • affiliation_column (object_id) — assembly affiliation

  • block_id_column (geom3) — subtomo_id of the source block in blocks

  • face_id_column (geom4) — 1-based face id within the assembly

  • face_size_column (geom5) — face vertex count

Parameters:
sizeint or None, default=None

Restrict to faces with exactly size vertices (e.g. size=6 for hexagons). None includes every face.

vertex_classdict of {int: int} or None, default=None

Exact vertex-class filter: {face_size: count, ...}. A block is included only when its set of face memberships matches every entry exactly — {5: 1, 6: 2} keeps only blocks that border exactly one pentagon and exactly two hexagons, with no other faces. A partial specification such as {5: 1} means exactly one pentagon and nothing else (not “one pentagon and any remainder”). None keeps all blocks.

scope{'original', 'inferred', 'all'}, default=``’original’``

'original' uses the real-block contact graph. 'all' includes inferred blocks. 'inferred' returns only blocks from inferred faces.

Returns:
cryomotl.Motl

One row per (block, face) pair, renumbered from 1.

Raises:
ValueError

If connect() has not been called, or if any two of block_id_column, face_id_column, face_size_column share the same column name.

get_faces_as_motl(size=None, *, scope='original')#

Return a Motl with one row per closed face.

Column assignments are controlled by BlockLayer constructor parameters (defaults in parentheses):

  • tomo_id_column (tomo_id) — tomogram identifier

  • affiliation_column (object_id) — assembly affiliation

  • face_id_column (geom4) — 1-based face id within the assembly

  • face_size_column (geom5) — face size (number of vertices)

  • x, y, z — face centroid

  • phi, theta, psi — Euler angles from face normal

Parameters:
sizeint or None, default=None

When set, return only faces with exactly size vertices (e.g. size=5 for pentagons). None returns all faces.

scope{'original', 'inferred', 'all'}, default=``’original’``

'original' uses the real-block contact graph. 'all' uses the combined lattice + inferred block graph. 'inferred' returns only faces that exist in 'all' but not in 'original'. Use the n_added column returned by get_face_stats() to count inferred blocks per face.

get_gaps_as_motl(cluster_radius, min_blocks=None, *, scope='original')#

Predict missing-block positions from unmatched contact sites.

Parameters:
cluster_radiusfloat

Maximum distance (voxels) between two gap points to be in one cluster.

min_blocksint or None, default=None

Minimum number of distinct source blocks in a cluster to keep it. Defaults to the highest derivable ideal degree from the block definitions (3, 4, or 6 based on symmetry), else 3.

scope{'original', 'inferred', 'all'}, default=``’original’``

'original' uses the real-block contact graph. 'all' uses the combined lattice + inferred block graph. 'inferred' uses only contacts in inferred faces.

Returns:
motlcryomotl.Motl

One row per kept cluster. Column assignments are controlled by BlockLayer constructor parameters (defaults in parentheses):

  • tomo_id_column (tomo_id) — tomogram identifier

  • affiliation_column (object_id) — assembly affiliation

  • x, y, z — cluster centroid

  • phi, theta, psi — Euler from mean z of source blocks

diagnosticspd.DataFrame

One row per kept cluster, indexed by gap_id (matching subtomo_id in motl). Columns:

  • gap_id — 1-based cluster id

  • n_source_blocks — distinct source blocks contributing to this cluster

get_gaussian_curvature()#

Return per-vertex Gaussian curvature for a mesh-backed surface.

Returns:
np.ndarray, shape (N,)

Gaussian curvature K = k1 * k2 at each vertex.

get_inverted_blocks(max_block_distance, *, scope='original')#

Return blocks whose orientation disagrees with a majority of their neighbours.

Uses the same proximity graph and z-normal comparison as unify_polarity() (reference='neighbours'). Does not modify the motl.

A block is considered inverted when n_disagree > n_agree, i.e. more than half of its proximity neighbours have a negative dot product with its z-normal (z_block · z_neighbour < 0).

Parameters:
max_block_distancefloat

Maximum centre-to-centre distance (voxels) for two blocks to be considered neighbours. Identical definition to unify_polarity().

scope{'original', 'inferred', 'all'}, default=``’original’``

'original' uses blocks only. 'all' uses lattice + inferred blocks combined. 'inferred' uses only inferred blocks.

Returns:
pandas.DataFrame

One row per inverted block, empty when none are found.

tomo_id

tomogram identifier

affiliation_column

assembly affiliation (from affiliation_column)

subtomo_id

particle identifier

n_neighbours

total proximity neighbours

n_agree

neighbours with z·z > 0

n_disagree

neighbours with z·z < 0

angle_to_mean_normal_deg

angle (°) between this block’s z-normal and the mean z-normal of its neighbours

Raises:
ValueError

If no block layer is present, or scope is invalid.

get_inverted_blocks_as_motl(max_block_distance, *, scope='original')#

Return inverted blocks as a Motl.

Calls get_inverted_blocks() and filters the selected block pool to those rows. Returns an empty Motl when no inverted blocks are found.

Parameters:
max_block_distancefloat

Passed directly to get_inverted_blocks().

scope{'original', 'inferred', 'all'}, default=``’original’``

Passed directly to get_inverted_blocks().

Returns:
cryomotl.Motl

Subset of the selected block pool containing only the inverted blocks, reset to a clean 1-based subtomo_id index.

Raises:
ValueError

If no block layer is present, or scope is invalid.

get_mean_curvature()#

Return per-vertex mean curvature for a mesh-backed surface.

Returns:
np.ndarray, shape (N,)

Mean curvature H = (k1 + k2) / 2 at each vertex.

get_missing_block_motl(*, scope='original')#

Infer positions of missing blocks from oversized face cycles.

A face cycle of length L in an assembly whose majority face size is m encloses k = L/m − 1 missing blocks when L is an exact multiple of m and L > m. Faces whose length is not a multiple of m are skipped (rim cycles, irregular holes).

Position: centroid of the three cycle vertices that bond to the missing block. For the j-th missing block (j = 0 … k−1) the bonding triplet is taken at equally-spaced positions around the cycle: round(j·L/k), round(j·L/k + L/3), round(j·L/k + 2·L/3) (all modulo L).

Orientation: mean SO(3) rotation of the three bonding blocks, computed with scipy.spatial.transform.Rotation.mean().

Parameters:
scope{'original', 'inferred', 'all'}, default=``’original’``

'original' uses the real-block contact graph. 'all' uses the combined lattice + inferred block graph. 'inferred' uses only inferred faces.

Returns:
cryomotl.Motl

One row per inferred missing block. Column assignments follow BlockLayer constructor parameters (defaults in parentheses):

  • tomo_id_column (tomo_id) — from the owning face

  • affiliation_column (object_id) — assembly id

  • face_id_column (geom2) — face_id of the owning face

  • x, y, z — inferred 3-D position

  • phi, theta, psi — ZXZ Euler of averaged neighbours

All other Motl columns are 0. An empty Motl is returned when no eligible faces are found.

Raises:
ValueError

If connect() has not been called for the requested scope.

get_neighboring_triangles(triangle_id, method='edge-connected', **kwargs)#

Get neighboring triangles (Mesh only).

Parameters:
triangle_idint

ID of the seed triangle.

method{‘edge-connected’, ‘radius’}, default=’edge-connected’

Traversal strategy. 'edge-connected' — topological edge walk. 'radius' — distance-based search (requires radius kwarg).

**kwargs

Additional parameters: - For ‘q’: max_hops (int, default=1) - For ‘radius’: radius (float, required), use_kdtree (bool, default=True)

Returns:
set or dict

For ‘edge-connected’: set of triangle IDs For ‘radius’: dict with ‘neighbor_ids’, ‘distances’, ‘seed_centroid’, ‘n_neighbors’

Raises:
TypeError

If surface is not a Mesh

ValueError

If invalid method or missing required parameters

get_point_neighborhoods(seed_point_ids, radii)#

Expand seed points on an oriented point cloud using surface radii.

Delegates to OrientedPointCloud.get_point_neighborhoods(). The returned dict always includes "hit points"; for each radius r, adds "r <= {r} nm" and annulus bands "{r_inner} < r <= {r_outer} nm" between consecutive radii.

Parameters:
seed_point_idsnp.ndarray

Integer indices of the seed points to expand from.

radiiArrayLike

Expansion radii in the same units as the point-cloud coordinates.

Returns:
dict[str, np.ndarray]

Keys: "hit points", "r <= {r} nm" for each radius, and "{r_inner} < r <= {r_outer} nm" for each consecutive pair. Values are sorted integer arrays of point indices.

get_points_within_distance(target, threshold)#

Find points within a distance threshold from a Mesh or an OrientedPointCloud surface.

Parameters:
targetnp.ndarray

Query points as (N, 3) array

thresholdfloat

Distance threshold

Returns:
dict

Dictionary containing: - ‘mask’: boolean array indicating points within threshold - ‘distances’: unsigned distances for all points - ‘indices’: indices of points within threshold - ‘within_points’: coordinates of points within threshold - ‘n_within’: number of points within threshold - ‘n_total’: total number of query points

get_principal_curvatures()#

Return per-vertex principal curvatures for a mesh-backed surface.

Returns:
np.ndarray, shape (N, 2)

Columns are the two principal curvature values k1 and k2 at each vertex.

get_shape_index()#

Return per-vertex shape index for a mesh-backed surface.

Returns:
np.ndarray, shape (N,)

Shape index S = (2/pi) * arctan2(k1 + k2, k1 - k2) at each vertex, in [-1, 1].

get_shape_stats()#

Return shape statistics for each object envelope in _envelopes.

Call envelope_from_faces() or envelope_from_blocks() first to populate _envelopes. Curvatures are computed internally on each mesh.

Returns:
pandas.DataFrame

One row per object. Columns:

  • tomo_id

  • assembly_id

  • scope — block pool used when the envelopes were built

  • radius — 1 / median(|H|) where H is the mean curvature at each vertex; equivalent to the characteristic curvature radius

  • sphere_fraction — Σ(K_i · area_i) / (4π) where K_i is the Gaussian curvature and area_i the Voronoi area at vertex i; equals 1 for a topological sphere, < 1 for partial coverage

  • area — total mesh surface area

  • b_a — ratio of the second-largest to largest semi-axis of the fitted ellipsoid (range 0–1; 1 = oblate/prolate degenerate, near 1 = near-sphere)

  • c_a — ratio of the smallest to largest semi-axis (range 0–1)

  • long_axis_x/y/z — unit vector along the longest ellipsoid axis

  • short_axis_x/y/z — unit vector along the shortest ellipsoid axis

  • fit_residual — RMS spread of per-vertex curvature radii around the median radius

Raises:
ValueError

If _envelopes is None (no envelope has been built).

get_sites_as_motl(scope='original')#

Return a particle list with one row per contact site of every block.

For each block in the pool selected by scope, each contact site defined in block_definition is placed at c + R.apply(v_k), where c is the block centre, R its rotation, and v_k the k-th site vector.

When a BlockDefinition has symmetry set (created via cyclic()), this method delegates to split_in_asymmetric_subunits(), which also rotates the output angles by the corresponding symmetry operation so that each site row carries the block rotation composed with the k-th cyclic rotation R_k = R_z(360*(k-1)/n).

For definitions without symmetry (non-symmetric / typed sites), expand_motl() is used and angles stay as the block rotation (orientation="keep").

The returned motl is independent of the internal block layer: later calls to connect() and all statistics methods read rotations from blocks directly, not from the returned motl.

Parameters:
scope{'original', 'inferred', 'all'}, default=``’original’``

'original' uses blocks only. 'all' uses lattice + inferred blocks combined. 'inferred' uses only inferred blocks.

Returns:
cryomotl.Motl

One row per (block, site) combination. The following columns carry values written by this method:

  • block_id_column (default geom3) — source block’s subtomo_id

  • site_index_column (default geom1) — 1-based site index (CCW order)

  • site_type_column (default geom2) — integer site-type code from site_type_codes

  • affiliation_column (default object_id) — assembly affiliation copied from source block (propagated when affiliation_column != block_id_column)

  • x, y, z — site position in voxels

All other columns are inherited from the parent block rows unchanged.

Raises:
ValueError

If no block layer is present, or scope is invalid.

get_surface_area()#

Return total surface area of a mesh-backed surface.

get_surface_type(as_labels=False)#

Return per-vertex categorical surface type for a mesh-backed surface.

Parameters:
as_labelsbool, default=False

If True, return string labels (e.g. "cap"); otherwise integer category codes (-1 flat, 0 cup .. 8 cap).

Returns:
np.ndarray, shape (N,)

Surface type per vertex.

get_triangle_neighborhoods(seed_triangle_ids, radii, use_kdtree=True)#

Expand seed triangles on the mesh using centroid-distance radii.

Always includes the "hit triangles" key. For each radius r in radii, adds a cumulative "r <= {r} nm" key and an annulus band "{r_inner} < r <= {r_outer} nm" between consecutive radii.

Parameters:
seed_triangle_idsnp.ndarray

Integer indices of the seed triangles to expand from.

radiiArrayLike

Expansion radii in the same units as the mesh coordinates. Each entry produces a cumulative shell and (for consecutive pairs) an annulus band.

use_kdtreebool, default=True

If True, use a KD-tree for centroid lookups; otherwise use brute-force search.

Returns:
dict[str, np.ndarray]

Keys: "hit triangles", "r <= {r} nm" for each radius, and "{r_inner} < r <= {r_outer} nm" for each consecutive pair. Values are sorted integer arrays of triangle indices.

get_triangles_within_radius(triangle_id, radius, use_kdtree=True)#

Return triangle-neighborhood query result for a mesh-backed surface.

Parameters:
triangle_idint

Seed triangle index.

radiusfloat

Maximum centroid-to-centroid distance for a triangle to be included.

use_kdtreebool, default=True

Use a KDTree for fast radius queries (recommended for large meshes).

property has_blocks#

True when a block layer is present.

property has_envelope#

True when an envelope surface (Mesh or OrientedPointCloud) is present.

property has_envelopes#

True when per-object envelopes have been built.

infer_missing_blocks(min_support=2, tip_tolerance=None, angle_tolerance=None, max_rounds=5, site_shift=None, clean_distance=None, min_real_support=0)#

Infer missing-block positions from unpaired legs via Cn symmetry fitting.

A missing block leaves its neighbours’ legs unpaired. For each unpaired leg the candidate centre is placed two site-shift lengths along the leg from the owning block’s centre. A candidate is kept only when at least min_support unpaired legs agree and their tips form a Cn-consistent arrangement.

The symmetry axis (z_c) is taken from the mean z-vector of the supporting real blocks, which carries the correct sign from the particles’ own orientations. The in-plane angle is then fitted in that frame. The tip-plane angle (angle between the SVD-derived tip-plane normal and the supporting blocks’ mean z-vector) is returned in the diagnostics DataFrame; large values flag geometrically suspect candidates.

Duplicate candidates (closer than clean_distance) are suppressed; the one with lower angular RMS wins; support breaks ties. Iteration continues until no new candidates are found or max_rounds is reached; later rounds carry higher iteration numbers in the output.

When block_definition is a dict each candidate position is tried against every definition; the type with the lowest angular RMS wins. The winning type is recorded in block_type_column and the runner-up RMS is returned in the diagnostics DataFrame.

After all rounds a final cleaning pass is run over the complete accepted set (all groups combined). Real blocks always beat any candidate; among candidates the priority is: earlier round, then higher support, then lower angular RMS. The number removed by this final pass is printed to stdout.

When the returned motl is non-empty the accepted candidates are also stored in _added (tagged with origin 2 in _origin_column). A subsequent call replaces the inferred portion of _added. Call connect(scope='all') to build a contact graph that includes both real and inferred blocks.

Parameters:
min_supportint, default=2

Minimum number of unpaired legs that must support a candidate. 1 is ambiguous (may be a genuine lattice edge); 2 is the practical minimum.

tip_tolerancefloat or None, default=None

Tolerance (voxels) for accepting a leg tip as near the candidate centre. None falls back to self.tip_tolerance; if that is also None, derives to 0.3 × |site_shift|.

angle_tolerancefloat or None, default=None

Maximum angular residual (degrees) for the Cn pairwise spacing check. Only applied to Cn-symmetric block definitions. None falls back to self.angle_tolerance; if that is also None, derives to 0.125 × (360 / n). Raises ValueError when the resolved value is >= 0.5 × (360 / n).

max_roundsint, default=5

Cap on iteration rounds.

site_shiftTripletLike or None, default=None

Override the shift vector used for candidate placement and support collection. None uses the block definition’s own first-site vector. Supply the corrected or mirrored vector when the contact is asymmetric — e.g. if the definition uses [45, -3, 0] but the partner’s tip sits at [45, +3, 0]. The magnitude determines the site length; the direction governs the in-plane orientation fit and virtual-leg placement.

clean_distancefloat or None, default=None

Deduplication distance (voxels) applied at two moments: (1) after each round, candidates closer than this to a better-ranked candidate or real block are discarded; (2) after all rounds, a final global pass removes any remaining duplicates. Higher support wins; lower angular RMS breaks ties. None falls back to self.clean_distance; if that is also None, derives to 0.4 × |site_shift|.

min_real_supportint, default=0

Minimum number of real-block legs a candidate must have to appear in the returned motl. 0 returns everything; 1 excludes pure-virtual candidates; 2 requires at least two confirmations from actual data.

Returns:
motlcryomotl.Motl

One row per inferred block. Carries position, orientation, affiliation, block type, and origin. Columns: tomo_id_column (tomo_id), affiliation_column (object_id), x, y, z, phi, theta, psi, block_type_column (default class), _origin_column (default geom5, value 2). All other columns are zero.

diagnosticspd.DataFrame

One row per inferred block, keyed by subtomo_id (matches the returned motl). Columns:

  • subtomo_id

  • n_support — total unpaired legs supporting the candidate (real + virtual); this is the value used by min_support

  • n_real_support — legs from real blocks; filtered by min_real_support

  • angular_rms — angular RMS of the winning type’s orientation fit (degrees)

  • iteration — round that accepted this candidate (1 = first)

  • block_type — type code of the winning definition

  • runner_up_rms — angular RMS of the runner-up definition (degrees); inf when only one definition was tried; a small margin (winner ≈ runner-up) means the type assignment is ambiguous

  • tip_plane_angle — angle (degrees) between the SVD-derived tip-plane normal and the supporting blocks’ mean z-vector; 0 for non-Cn blocks; large values flag geometrically suspect candidates

Raises:
ValueError

If connect() has not been called or no block definition is set.

intersection_data(result, query_type=None, min_distance_source_target=None, max_distance_source_target=None, source_id_name='source_id', target_id_name='target_id', include_curvatures=True, surface_radii=None, surface_element=None, surface_seeds='auto', use_kdtree=True)#

Turn raw intersection or distance-query output into analysis-ready tables.

Works with results from ray_intersections() (mesh or point cloud target) and distance_to_pointcloud(). Optional surface_radii expansion grows regions around hit sites on self (triangles/vertices on meshes, points on oriented point clouds).

Parameters:
resultdict

Output of ray_intersections() or distance_to_pointcloud().

query_type{‘ray’, ‘distance_to_pointcloud’}, optional

Inferred from result when omitted. The body also accepts a small set of aliases ('rays', 'distance', 'pointcloud', 'distance_to_point_cloud') for backward compatibility.

min_distance_source_target, max_distance_source_targetfloat, optional

Keep hits whose source-target distance lies in this interval.

source_id_name, target_id_namestr

Column names for source and target indices in the hit table.

include_curvaturesbool, default=True

Attach curvature columns for mesh-backed self.

surface_radiisequence of float, optional

Radii for cumulative regions around hit seeds on self.

surface_element{‘triangles’, ‘vertices’, ‘points’}, optional

Defaults: ray+mesh → triangles; distance+mesh → vertices; point cloud → points.

surface_seeds{‘auto’, ‘hit_sources’, ‘hit_targets’}, default=’auto’

Which hit IDs seed expansion.

use_kdtreebool, default=True

Passed to mesh triangle expansion.

Returns:
dict
  • hits: hit table

  • regions: region name → index arrays (if surface_radii set)

  • region_summary: per-region summary table

  • triangle_curvatures: per-triangle arrays (mesh-backed self)

invalidate_caches()#

Invalidate cached geometry on the wrapped surface (mesh ray scene, neighbor trees, etc.).

property is_mesh#

True when the backing geometry has triangle connectivity (Mesh).

property is_point_cloud#

True when the backing geometry is discrete samples (OrientedPointCloud).

property normals#

DiscreteSurface normals.

oversample(**kwargs)#

Delegate to oversample on surface; mesh and point-cloud semantics differ.

Parameters:
**kwargs

Mesh (Mesh.oversample()):

  • oversample_factor : float, optional — desired factor increase in vertices. Defaults to 1.0 (no change) when both this and point_spacing are None.

  • point_spacing : float, optional — desired spacing between sampled points (same units as mesh coordinates). Uses two-pass Poisson-disk calibration.

  • poisson_init_factor : int, default=5 — initial candidate factor for Poisson-disk sampling (larger → more uniform distribution).

OrientedPointCloud (OrientedPointCloud.oversample()):

  • oversample_factor : float, optional — desired factor increase in points. Defaults to 1.0 when both this and point_spacing are None.

  • point_spacing : float, optional — desired spacing; uses greedy Poisson-disk sampling to enforce spacing directly.

  • random_seed : int, optional — seed for reproducible sampling.

Returns:
PleomorphicSurface

New wrapper around the resampled surface.

property pixel_size#

Ångström-per-voxel scale for block-layer computations.

Alias for block_pixel_size. The mesh’s own scale is stored on the Mesh object and is never overwritten by this value.

plot_object(object_id, *, plot_faces=True, scope='original')#

Return a Plotly figure showing all blocks of one object, coloured by provenance.

One trace is drawn per origin class present in the data: "Real" (origin 1) and "Inferred" (origin 2). When plot_faces is True, face centroids are overlaid as a second set of traces — one per distinct face edge-count — using diamond markers larger than the block markers.

The figure is returned without being displayed. Pass it to fig.show() or use browse_objects() for an interactive notebook browser.

Parameters:
object_idfloat

Value of affiliation_column that identifies the object.

plot_facesbool, default=True

When True, overlay face centroids coloured by edge count.

scope{'original', 'inferred', 'all'}, default=``’original’``

'original' draws lattice blocks only. 'all' includes inferred blocks as well. 'inferred' draws only inferred faces.

Returns:
plotly.graph_objects.Figure
ray_intersections(rays, one_hit_per_target=False, knn_radius=10.0, return_orientations=False, target_orientation='normal')#

Compute ray intersections with this surface.

For meshes, uses Open3D’s exact raycasting. For oriented point clouds, uses KDTree-based nearest neighbor search along ray trajectories.

Parameters:
raysnp.ndarray, shape (N, 6)

Ray array where each row is [origin_x, origin_y, origin_z, dir_x, dir_y, dir_z]

one_hit_per_targetbool, default=False

If True and multiple rays are supplied, add shortest_distance and shortest_indices for the global shortest hit across all rays.

knn_radiusfloat, default=10.0

For OrientedPointCloud only: maximum search radius for finding nearest points along ray trajectory

return_orientationsbool, default=False

If True, compute relative orientations between ray directions and surface orientations at hit points. Returns additional fields: - ‘ray_directions’: normalized ray direction vectors - ‘surface_orientations’: orientation vectors at hit points - ‘angles_deg’: angles in degrees between ray and surface orientation - ‘dot_products’: dot products (cosine of angle)

target_orientation{‘normal’, ‘principal_1’, ‘principal_2’} or callable, default=’normal’

Which orientation to use for comparison. Options:

For OrientedPointCloud:
  • ‘normal’: Use the normals stored in the point cloud (for filaments, these are the axis/tangent directions)

For Mesh:
  • ‘normal’: Use surface normals at hit points (default)

  • ‘principal_1’: Use first principal curvature direction

  • ‘principal_2’: Use second principal curvature direction

Custom function:
  • A callable that takes (surface, primitive_ids, hit_points) and returns an (N, 3) array of orientation vectors, where N is the number of hits. Example: lambda surf, ids, pts: surf.get_curvature_directions()[ids, :, 0]

Returns:
dict

Always contains:

  • t_hit (R,): ray travel distance to the intersection; inf for misses.

  • primitive_ids (R,): triangle index (Mesh) or point index (OrientedPointCloud) of the hit; -1 for misses.

  • hit_points (R, 3): 3-D coordinates of hit points; NaN for misses.

  • primitive_normals (R, 3): surface normals at hit points (Mesh only); NaN for misses or when the backing surface is an OrientedPointCloud.

  • geometry_ids (R,): geometry identifier (Mesh only).

When return_orientations=True, also adds:

  • ray_directions (R, 3): normalized ray direction vectors.

  • surface_orientations (R, 3): surface orientation vectors at hit points; NaN for misses.

  • angles_deg (R,): angle in degrees between the ray and surface orientation; NaN for misses.

  • dot_products (R,): cosine of that angle; NaN for misses.

classmethod read(input_path, method='mesh', **kwargs)#

Create a wrapped surface from common on-disk inputs.

Parameters:
input_pathstr or Path

Input file path.

method{‘mesh’, ‘mesh_curvatures’, ‘mesh_from_mrc’, ‘point_cloud’, ‘point_cloud_from_mrc’, ‘point_cloud_from_motl’}, default=’mesh’

Loader to use:

  • “mesh”: geometry-only triangle mesh via Mesh.read()

  • “mesh_curvatures”: VTP triangle mesh with curvature fields via Mesh.read_curvatures()

  • “mesh_from_mrc”: segmentation-to-mesh via Mesh.from_mrc()

  • “point_cloud”: oriented point cloud via OrientedPointCloud.read()

  • “point_cloud_from_mrc”: segmentation-to-point-cloud via OrientedPointCloud.from_mrc()

  • “point_cloud_from_motl”: oriented point cloud from a motl file via OrientedPointCloud.from_motl(). Pass group_by (e.g. 'object_id' or 'subtomo_id') to split the motl into one wrapped surface per unique value; returns a dict of PleomorphicSurface in that case, otherwise a single PleomorphicSurface.

**kwargs

Forwarded to the selected loader. Accepted keywords depend on method:

“mesh” and “mesh_curvatures”:

  • units : str, optional — coordinate units ('nm', 'angstrom', 'pixel', …).

“mesh_from_mrc”:

  • transpose : bool, default=True — transpose the segmentation array on load.

  • labels_dict : dict, optional — map label names to integer values; binary if None.

  • level : float, default=0.5 — marching-cubes iso-level.

  • pixel_size : float, default=1.0 — voxel size for coordinate scaling.

  • smooth_sigma : float, optional — Gaussian pre-smooth sigma.

  • step_size : int, default=1 — marching-cubes step size.

“point_cloud”:

  • recompute_normals : bool, default=False — recompute even if file has normals.

  • knn : int, default=30 — neighbors for normal estimation.

  • orient_normals : bool, default=True — orient normals consistently.

  • tangent_plane_knn : int, default=50 — neighbors for normal orientation.

“point_cloud_from_mrc”:

  • labels_dict : dict, optional — map label names to integer values; binary if None.

  • pixel_size : float or array-like, optional — voxel size.

  • compute_normals : bool, default=True — estimate normals after extraction.

  • knn : int, default=30 — neighbors for normal estimation.

  • orient_normals : bool, default=True — orient normals consistently.

  • tangent_plane_knn : int, default=50 — neighbors for normal orientation.

  • transpose : bool, default=True — transpose the segmentation array on load.

  • smooth_sigma : float, optional — Gaussian pre-smooth sigma.

“point_cloud_from_motl”:

  • group_by : MotlColumn, optional — motl column to split on (e.g. 'object_id', 'subtomo_id'). If given and >1 unique value exists, returns a dict of wrapped surfaces keyed by that column’s value.

  • recompute_normals : bool, default=False — recompute normals from geometry.

  • knn : int, default=30 — neighbors for normal estimation.

  • orient_normals : bool, default=True — orient normals consistently.

  • tangent_plane_knn : int, default=50 — neighbors for normal orientation.

Returns:
PleomorphicSurface or dict[Any, PleomorphicSurface]

Wrapped surface loaded from input_path. For method="point_cloud_from_motl" with group_by splitting into multiple groups, a dict of wrapped surfaces keyed by group value.

recover_orientations(motl, *, min_neighbours=2, bonding_tolerance=None, angle_tolerance=None)#

Recover orientations for particles removed from blocks.

Fits each candidate’s orientation from the surrounding reference blocks in blocks (the cleaned population). Candidates never support each other — only real reference blocks are used. Results are stored in _added tagged with origin 3 in _origin_column. Running again replaces the origin-3 portion.

The input motl may be either the discarded particles (all rows are treated as candidates) or the full population (rows whose subtomo_id is not in blocks become candidates).

When block_definition is a dict, each candidate’s definition is resolved from its block_type_column value. A candidate whose type has no definition is skipped and included in unfitted; a warning naming the missing types is printed.

Parameters:
motlcryomotl.Motl

Particles to recover.

min_neighboursint, default=2

Minimum number of distinct reference-block neighbours a candidate needs for a valid fit.

bonding_tolerancefloat or None, default=None

Half-width (voxels) of the acceptance window around the nominal centre-to-centre bonding distance 2 × |site_shift|. Overrides bonding_tolerance for this call only. None falls back to self.bonding_tolerance; if that is also None, derives to 0.4 × |site_shift|.

angle_tolerancefloat or None, default=None

Maximum allowed angular residual (degrees) for the Cn pairwise spacing check. Overrides angle_tolerance for this call only. None falls back to self.angle_tolerance; if that is also None, derives to 0.125 × (360 / n). Raises ValueError when the resolved value is >= 0.5 × (360 / n).

Returns:
motlcryomotl.Motl

One row per successfully fitted candidate. Carries position, orientation, affiliation, block type, and origin. Columns: tomo_id_column (tomo_id), affiliation_column (object_id), x, y, z, phi, theta, psi, block_type_column (default class), _origin_column (default geom5, value 3). All other columns are zero.

diagnosticspd.DataFrame

One row per candidate (fitted and unfitted), keyed by subtomo_id. Columns:

  • subtomo_id

  • status — "fitted" or "unfitted"

  • n_neighbours — number of accepted reference neighbours

  • angular_rms — angular RMS of the Cn check (degrees); NaN for candidates that failed before the Cn check

  • normal_spread — RMS spread of the accepted neighbours’ normals around their mean (degrees); NaN for unfitted

  • tip_plane_angle — angle (degrees) between the SVD-derived tip-plane normal and the neighbours’ mean z-vector; NaN for unfitted candidates

  • reason — empty string for fitted candidates; for unfitted: "no_reference_blocks", "undefined_type", "insufficient_neighbours", or "cn_check_failed"

Raises:
ValueError

If no block definition or block layer is present.

refine_normals(radius_hit=3.0, batch_size=2000, n_iter=1, mask=None, logger=None, inplace=True, **kwargs)#

Refine normals on the wrapped surface by neighborhood averaging.

Delegates to Mesh.refine_normals() or OrientedPointCloud.refine_normals() (both inherit DiscreteSurface.refine_normals()).

Parameters:
radius_hitfloat, default=3.0

Neighborhood radius for normal averaging, in mesh/point-cloud units.

batch_sizeint, default=2000

Batch size for spatial neighbor queries.

n_iterint, default=1

Number of refinement passes.

masknp.ndarray, optional

Boolean mask of vertices/samples to update. If None, all are refined.

loggerlogging.Logger, optional

Logger passed through to the delegate.

inplacebool, default=True

If True, update the wrapped surface in place. If False, return a new wrapper.

**kwargs

Additional keyword arguments forwarded to the delegate.

Returns:
PleomorphicSurface

self when inplace=True; a new wrapper when inplace=False.

regroup_objects(split=True, merge=True, original_column='original_object_id')#

Return a new motl with affiliations reassigned so each holds one connected piece.

Acts on blocks without modifying it. The original affiliation value is preserved in original_column so the change is traceable.

Parameters:
splitbool, default=True

When True, each connected component of the within-affiliation contact graph is promoted to its own affiliation. Multi-piece affiliations are broken apart.

mergebool, default=True

When True, affiliations joined by cross-affiliation contacts are folded into one. The block-level contact graph (including cross-affiliation edges) determines which affiliations are physically connected.

original_columnstr, default=”original_object_id”

Name of the column added to the returned motl to record the original affiliation value.

Returns:
cryomotl.Motl

Deep copy of blocks with affiliation_column rewritten. Affiliation values are assigned sequentially from 1 within each tomogram, ordered by the smallest block subtomo_id in the group. The returned motl is independent of self; blocks is unchanged.

Raises:
ValueError

If connect() has not been called.

Notes

split=False, merge=False returns an unchanged copy (plus the original_column).

remove_nonfinite_vertices(inplace=True, **kwargs)#

Remove NaN/Inf vertices or point samples from the wrapped surface.

For meshes, affected faces are also dropped and vertex connectivity is remapped.

Parameters:
inplacebool, default=True

If True, modify the wrapped surface in place and return self. If False, return a new PleomorphicSurface wrapping a repaired copy.

**kwargs
  • recompute_normals : bool — recompute normals after filtering. Default is True for Mesh, False for OrientedPointCloud.

Returns:
PleomorphicSurface

self when inplace=True; a new wrapper around the repaired surface when inplace=False.

save(output_path, format=None, **kwargs)#

Save the wrapped surface.

If format is None, the wrapped surface may infer it from output_path. Additional keyword arguments are forwarded to the concrete surface save method.

Parameters:
output_pathPathOrStr

Output file path.

format{‘ply’, ‘vtp’, ‘motl’, ‘em’}, optional

File format. 'ply' and 'vtp' are Mesh formats; 'motl'/'em' are OrientedPointCloud formats. If None, inferred from the file suffix.

**kwargs

Forwarded to the concrete save method. Accepted keywords depend on format:

Mesh (format='ply' or 'vtp'):

  • include_curvatures : bool, default=False — embed per-vertex curvature scalars and principal-direction vectors; requires curvatures to have been computed and format='vtp'.

OrientedPointCloud — PLY (format='ply'):

  • write_ascii : bool, default=False — write ASCII rather than binary PLY.

OrientedPointCloud — MOTL / EM (format='motl' or 'em'):

  • input_dict : dict, optional — extra motive-list columns to fill.

  • subtomo_ids : array-like, shape (N,), optional — per-point subtomogram IDs; sequential IDs are assigned when None.

  • tomo_id : int, float, or array-like, optional — tomogram ID; scalar applies to all points, array assigns per-point IDs.

Returns:
None
separate_surfaces(surface_type='closed', threshold_angle=90.0, reference_point=None)#

Separate the two halves of a Mesh or OrientedPointCloud surface.

Parameters:
surface_type{‘closed’, ‘planar’}, default=’closed’

Strategy for separation:

'closed'

For enclosed volumes (vesicles, organelles). Classifies each vertex by whether its normal points toward or away from a reference point (default: centroid). Returns (inner_mask, outer_mask). Calls separate_closed_surface().

'planar'

For surfaces with lower curvature or flatter geometry. Uses PCA on normals to find the axis of greatest normal spread and splits by projection sign. Returns (surface1_mask, surface2_mask) — no inherent inner/outer meaning; inspect spatially to assign labels. Calls separate_planar_surface().

threshold_anglefloat, default=90.0

Angle threshold in degrees. Only used when surface_type='closed'.

reference_pointnp.ndarray (3,), optional

Reference point for 'closed' separation. Defaults to centroid. Ignored for 'planar'.

Returns:
mask1, mask2(N,) bool ndarray each

For 'closed': (inner_mask, outer_mask). For 'planar': (surface1_mask, surface2_mask). Pass either mask directly to apply_vertex_mask().

store_block_stat(stat, column)#

Write one block-stats column into the blocks motl and return the updated blocks.

Parameters:
stat{“n_sites”, “degree”, “complete”, “assembly_id”, “angle_sum”, “angle_deficit”}

Name of the get_block_stats() column to store.

columnMotlColumn

Motl column to write the values into.

Returns:
cryomotl.Motl

Deep copy of blocks after writing.

Raises:
ValueError

If stat is not one of the accepted column names.

store_block_stats(columns)#

Write selected block-stats columns into the blocks motl.

Parameters:
columnsdict[str, MotlColumn]

Mapping of get_block_stats() column name → motl column to write into. Only numeric columns are accepted.

Raises:
ValueError

If a requested column is non-numeric (e.g. face_signature) or not a recognised get_block_stats() column.

property surface#

The wrapped Mesh or OrientedPointCloud.

Raises:
ValueError

If this instance has no envelope (was created with blocks only).

unify_polarity(max_block_distance=None, reference='neighbours', *, _annotate_fn=None)#

Align block polarities so neighbours point the same way.

Works per tomogram. For reference='neighbours' or 'centroid', uses a BFS walk over the block-proximity graph (edges between blocks whose centres are within max_block_distance): when a newly visited block has z_n · z_m < 0 relative to its BFS parent, it is flipped by the 180° rotation around the in-plane part of its flip_site vector.

For reference='centroid', each connected component is additionally aligned so that its blocks point outward from the component centroid (Σ z · (c - g) > 0).

For reference='envelope', an envelope surface must be attached (via ps.surface = ...). max_block_distance may be None.

Parameters:
max_block_distancefloat or None, default=None

Maximum centre-to-centre distance (voxels) for two blocks to be considered neighbours. Required for 'neighbours' and 'centroid'; ignored for 'envelope'.

reference{“neighbours”, “centroid”, “envelope”}, default=”neighbours”

Alignment reference.

_annotate_fnCallable or None, default=None

Private override for the envelope annotation callable. When None and reference='envelope', falls back to annotate_with_envelope().

Returns:
int

Number of blocks whose final orientation differs from the initial one.

Raises:
ValueError

If no block layer is present, reference is invalid, max_block_distance is None for a BFS-based reference, or reference='envelope' is requested but no envelope is attached.

property units#

Coordinate units stored on the wrapped surface.

property vertices#

DiscreteSurface vertices / points.

class cryocat.analysis.assembly.PolyhedralComplex(motl, *, affiliation_column='object_id', order_column='geom1', tomo_id_column='tomo_id')#

Bases: SymmetricComplex

Abstract base for Platonic-solid (T/O/I) motl-analysis complexes.

Do not instantiate this class directly. Use one of the concrete subclasses:

Subclasses declare two class attributes:

_solid

The cryocat.utils.geom solid class (e.g. geom.Icosahedron).

_symmetry

The symmetry letter "T", "O", or "I".

All logic lives here; subclasses only override those two attributes.

angular_dissimilarity(rotations_1, rotations_2)#

Symmetry-aware dissimilarity of paired particle orientations.

Compares the corners of the complex’s solid turned by each pair of rotations (see cryocat.utils.geom.Polyhedron.angular_dissimilarity()); 0 means the two particles look the same up to the complex’s symmetry, and higher values mean more different orientations.

Uses self.solid when geometry has been fitted via fit_geometry(), so the corners follow the symmetry axes of the reference map; otherwise the canonical solid is used (reference assumed canonically oriented), as in feature_vectors() and symmetry_group().

Parameters:
rotations_1RotationLike

First orientation of each pair: one rotation or a stack of N, e.g. motl_a.get_rotations().

rotations_2RotationLike

Second orientation of each pair, same number as rotations_1 (e.g. motl_b.get_rotations()). A single rotation on either side is compared with every rotation on the other side.

Returns:
np.ndarray

(N,) dissimilarities in radians.

Raises:
ValueError

If either input is empty, or the numbers of rotations differ and neither input holds exactly one.

assign_subunit_order()#

Assign 1-based subunit indices ordered by x→y→z (ascending).

center = None#
expand(*, mode='vertices', project_to_sphere=False, radius=1.0, shift_vecs=None, original_id_col='object_id', order_id_col='geom1', output_motl_type='emmotl', output_path=None, **output_kwargs)#

Expand each particle into subparticles at polyhedral feature positions.

Parameters:
mode{“vertices”, “edges”, “faces”}, default=”vertices”

Feature type used when shift_vecs is None.

project_to_spherebool, default=False

Project features to the circumscribed sphere (passed to feature_vectors()).

radiusfloat, default=1.0

Fallback solid radius when no geometry is fitted and shift_vecs is None.

shift_vecsnp.ndarray of shape (N, 3), optional

Explicit shift vectors. When None, vectors are derived from feature_vectors().

original_id_colMotlColumn, default=”object_id”

Column that stores the subtomo_id of the source particle.

order_id_colMotlColumn, default=”geom1”

Column that stores the 0-based subunit extraction index.

output_motl_typeMotlType, default=”emmotl”

Format of the returned/written motive list.

output_pathPathOrStr, optional

Write path. No file is written when None.

**output_kwargs

Forwarded to cryocat.core.cryomotl.motl_converter_kwargs().

Returns:
MotlSource

Expanded motive list in the requested format.

Raises:
ValueError

If shift_vecs has wrong shape or a required column is missing.

Notes

The number of subparticles per particle is the solid’s vertex, edge or face count (icosahedron: 12, 30 or 20), not the group order (60 for "I") used by split_in_asymmetric_units(). Corners, edge midpoints and face centres lie on symmetry axes, so the group’s rotations only produce order / n distinct places for them (table in the cryocat.utils.symmetry module Notes).

feature_vectors(mode='vertices', project_to_sphere=False, radius=1.0)#

Return feature centres for the fitted solid (or a canonical one at radius).

Uses self.solid when geometry has been fitted via fit_geometry(); otherwise falls back to the canonical solid at radius.

Parameters:
mode{“vertices”, “edges”, “faces”}, default=”vertices”

Which feature centres to return.

project_to_spherebool, default=False

When True, rescale each vector to the circumscribed sphere radius.

radiusfloat, default=1.0

Fallback radius when no geometry is fitted.

Returns:
np.ndarray

Array of shape (N, 3).

fit_geometry(markers, reference_map, center=None)#

Fit self.solid and self.center from two non-collinear markers.

Reads the first two marker positions from markers, resolves the box centre from reference_map when center is None, and sets self.solid = self._solid.from_vectors(v1 - c, v2 - c). Pixel size is read from the map and used to convert marker coordinates (Å) to voxels before building the solid.

classmethod recover_features(input_cmm_file, input_map, *, center=None, mode='vertices', project_to_sphere=False, output_cmm_file=None)#

Recover polyhedral feature coordinates from two marker positions.

Must be called on a concrete subclass (TetrahedralComplex.recover_features(…), etc.).

Parameters:
input_cmm_filePathOrStr

Path to a .cmm file with markers for two non-collinear vertices.

input_mapPathOrStr

Map used to prepare the marker file (supplies pixel size and box dimensions).

centerTripletLike, optional

Centre of the solid in map voxels. Defaults to the box centre.

mode{“vertices”, “edges”, “faces”}, default=”vertices”

Feature type to recover.

project_to_spherebool, default=False

Project features to the circumscribed sphere.

output_cmm_filePathOrStr, optional

If given, write recovered features to this .cmm file.

Returns:
feature_vecnp.ndarray of shape (N, 3)

Vectors from the solid centre to each feature.

features_coordsnp.ndarray of shape (N, 3)

Feature positions in the map box (in Å).

Raises:
TypeError

When called on PolyhedralComplex directly rather than on a concrete subclass.

ValueError

If mode is invalid.

solid = None#
split_in_asymmetric_units(xyz_shift, *, output_motl_type='emmotl', output_path=None, **output_kwargs)#

Split each particle into its asymmetric units, in the fitted frame.

Calls cryocat.core.cryomotl.Motl.split_in_asymmetric_subunits() with the orientation of symmetry_group(), so the subunits are placed around the symmetry axes of the reference map rather than the canonical ones. This works for references in any orientation and uses the same frame as expand(). Without fit_geometry(), the reference is assumed to be canonically oriented.

Parameters:
xyz_shiftArrayLike

Position of the reference subunit in the reference map, in voxels relative to the map centre (the frame of feature_vectors()). It should lie off every symmetry axis; otherwise copies overlap.

output_motl_typeMotlType, default=”emmotl”

Format of the returned/written motive list.

output_pathPathOrStr, optional

Write path. No file is written when None.

**output_kwargs

Forwarded to cryocat.core.cryomotl.motl_converter_kwargs().

Returns:
MotlSource

Expanded motive list with 12, 24 or 60 subunits per particle.

Notes

One subunit per group rotation (60 for "I"), which differs from the 12/30/20 subparticles expand() places on the icosahedron’s vertices/edges/faces; see the cryocat.utils.symmetry module Notes.

symmetry_group()#

Return the complex’s symmetry group in its fitted orientation.

After fit_geometry(), the group is derived from self.solid via cryocat.utils.symmetry.SymmGroup.from_polyhedron(), so its rotations match the symmetry axes of the reference map (checked for consistency with the fitted solid). Before fitting, the canonical group is returned, matching the canonical solid used by feature_vectors().

Returns:
cryocat.utils.symmetry.SymmGroup

12, 24 or 60 rotations for T, O or I; .rotation holds the fitted orientation (identity if not fitted).

vertex_subunit_vectors(subunit_markers, *, symmetry=None, reference=None, output_cmm_file=None)#

Place a ring of subunits around every vertex of the fitted solid.

One marker on one subunit near a vertex is turned about that vertex’s axis into a ring of n evenly spaced copies (local Cn), and the ring is copied to every vertex with cryocat.utils.geom.Polyhedron.transport_points(). Requires fit_geometry() first, which supplies the solid, its centre and the pixel size.

Parameters:
subunit_markersPathOrStr

Marker file (e.g. ChimeraX .cmm, Å) placed on the reference map used in fit_geometry(); only the first marker is used.

symmetrySymmetry, optional

Cyclic symmetry of the ring around each vertex, e.g. "C5" or 5. Defaults to the fold of the vertex axis (T: 3, O: 4, I: 5).

referenceint, optional

Index of the vertex the marked subunit belongs to. Defaults to the vertex closest in angle to the marker.

output_cmm_filePathOrStr, optional

If given, write the subunit positions (Å) to this .cmm file.

Returns:
vecsnp.ndarray of shape (V * n, 3)

Vectors in voxels from the solid’s centre, vertex by vertex (rows i * n to i * n + n - 1 belong to vertex i); ready for expand(shift_vecs=vecs).

coordsnp.ndarray of shape (V * n, 3)

The same positions as map coordinates in Å.

Raises:
ValueError

If no geometry has been fitted, the marker file is empty, or reference is not a valid vertex index.

NotImplementedError

If symmetry is not cyclic.

Warns:
UserWarning

If n is not a multiple of the vertex fold (the copies then depend on the neighbour ordering and do not share the solid’s symmetry), or if the marker lies on the vertex axis (all ring copies overlap).

Notes

  • With the default n (the vertex fold), the result is the same set of positions as the symmetry group’s orbit of the marker (self.symmetry_group().orbit(shift)), i.e. the positions used by split_in_asymmetric_units(); the orientations assigned later by expand() are radial instead.

  • expand() sorts the shift vectors, so its order column does not follow the vertex-by-vertex order of vecs.

Examples

>>> ico = IcosahedralComplex("complexes.em")
>>> ico.fit_geometry("two_vertices.cmm", "reference.mrc")
>>> vecs, coords = ico.vertex_subunit_vectors("one_subunit.cmm")
>>> vecs.shape                      # 12 vertices x 5 subunits
(60, 3)
>>> sub_motl = ico.expand(shift_vecs=vecs)
write_features_cmm(output_path, mode='vertices', project_to_sphere=False)#

Write feature vectors as ChimeraX marker coordinates to a .cmm file.

Requires that fit_geometry() has been called first.

Parameters:
output_pathPathOrStr

Destination .cmm file path.

mode{“vertices”, “edges”, “faces”}, default=”vertices”

Feature type to write.

project_to_spherebool, default=False

Project features to the circumscribed sphere before writing.

Raises:
ValueError

When no geometry has been fitted yet.

class cryocat.analysis.assembly.SymmetricComplex(motl, symmetry, *, affiliation_column='object_id', order_column='geom1', tomo_id_column='tomo_id')#

Bases: object

Base class for point-group symmetric multi-subunit complexes.

Encapsulates a motl of subunit particles that form a symmetric complex (cyclic, dihedral, or Platonic) and provides per-object centre computation, orientation unification, and geometric statistics that are independent of the specific symmetry type.

Parameters:
motlMotlSource

Subunit particle list.

symmetrySymmetry

Symmetry specifier, e.g. "C8", "D6", "T", "O", "I", or a bare integer (interpreted as cyclic).

affiliation_columnMotlColumn, default=’object_id’

Column that identifies which object each particle belongs to.

order_columnMotlColumn, default=’geom1’

Column that subunit-ordering methods write indices into.

tomo_id_columnMotlColumn, default=’tomo_id’

Column that identifies the tomogram.

assign_subunit_order()#

Assign subunit indices into self.order_column.

Subclasses must override this method with a symmetry-specific implementation.

Raises:
NotImplementedError

Always; subclasses define subunit ordering.

clean_per_object(column, keep='high', *, n=None)#

Keep the n best rows per object and drop the rest.

For each (tomo_id_column, affiliation_column) group: sort by column, keep the top-n rows according to keep, and discard the remainder. Objects that already have at most n rows are returned unchanged.

Parameters:
columnMotlColumn

Column to sort and filter by.

keep{‘high’, ‘low’}, default=’high’

'high' retains the n rows with the largest values (e.g. scores); 'low' retains the n rows with the smallest values (e.g. cone distances).

nint, optional

Number of rows to keep per object. Defaults to self.n_subunits.

Returns:
cryomotl.Motl

A copy of self.motl with over-occupied objects trimmed.

Raises:
ValueError

If affiliation_column is absent from self.motl.df.

create_affiliation(method='radius', *, shift=None, radius=None, normals_threshold=None, occupancy_column=None, cone_distance_column=None, min_occupancy=None)#

Cluster subunit particles into objects and write affiliation labels.

Operates on a copy of self.motl and returns it with affiliation_column populated. For every object that ends up with more than self.n_subunits particles, outlier-normal particles are removed (see normals_threshold). If an object remains over-populated after cleaning its geom5 column is set to 1.0 and a one-line notice is printed.

Subunit ordering is a separate stage-2 call to assign_subunit_order(). Occupancy and normal deviations are computed on demand via occupancy() and normal_deviations().

Parameters:
method{‘radius’, ‘tracing’}, default=’radius’

Clustering strategy.

'radius'

Optionally shift particles along their local −x axis by shift, then run a self nearest-neighbour search within radius. Connected components of the NN graph become objects. Isolated particles (no NN within radius) are kept as singleton objects with unique affiliation_column values.

'tracing'

Optionally shift particles along their local −x axis by shift, then trace chains via nnana.trace_chains() with radius as the maximum link distance. Each chain becomes one object.

shiftfloat, optional

Magnitude of the local-frame shift along −x applied before clustering (voxels). When None no recentring is performed. Typical value: approximate ring radius.

radiusfloat

For method='radius': NN search radius (voxels). For method='tracing': maximum chain-link distance (voxels). Required.

normals_thresholdfloat or None, default=None

Cone-distance cutoff (degrees) applied only to over-populated objects (size > n_subunits). None uses the adaptive rule median(cone) + 1 × std(cone) per object.

occupancy_columnMotlColumn or None, default=None

Accepted but ignored. Use occupancy() on demand.

cone_distance_columnMotlColumn or None, default=None

Accepted but ignored. Use normal_deviations() on demand.

min_occupancyint or None, default=None

When set, objects whose size (after over-population cleaning) is below this value are removed. None keeps all objects.

Returns:
cryomotl.Motl

A new motl with affiliation_column and geom5 populated. geom5 is 0.0 for all particles; 1.0 marks objects that are still over-populated after cleaning.

Raises:
ValueError

If radius is None or method is unrecognised.

get_centers_as_motl()#

Return a Motl with one barycentric centre particle per object per tomogram.

Iterates over all tomograms in self.motl, groups by affiliation_column within each tomogram, and returns the barycentric centre of each group.

Returns:
Motl

One row per (tomogram, object) pair. tomo_id holds the tomogram identifier and object_id holds the affiliation value. All other columns are zero-filled.

merge_subunits(radius=55)#

Merge near-duplicate objects whose centres are within radius.

For each tomogram:

  1. Compute per-object barycentric centres.

  2. Find object-centre pairs within radius using nnana.get_nn_within_distance().

  3. Re-assign affiliation_column of near objects to the first encountered partner.

  4. Recount occupancy into geom1 and recompute subunit order for all objects via assign_subunit_order().

Parameters:
radiusfloat, default=55

Distance threshold in voxels.

Notes

Modifies self.motl.df in place.

normal_deviations()#

Per-particle angle to per-object consensus z-axis with object stats.

For each particle, reports the angle (degrees) between its z-normal (from phi/theta/psi via geom.euler_angles_to_normals()) and the mean z-normal of its object. Also reports per-object median and std of those angles.

Returns:
pandas.DataFrame

Columns: tomo_id, object_id, subtomo_id, deviation, object_median, object_std.

occupancy()#

Per-object subunit occupancy.

For every (tomo_id, object_id) group, counts present subunits and computes the fraction of the expected n_subunits. When order_column is populated, the missing indices (1 … n_subunits) are also reported.

Returns:
pandas.DataFrame

Columns:

tomo_id, object_id

Group identifiers.

n_present

Number of particles in the group.

occupancy

n_present / self.n_subunits.

missing

Sorted list of 1-based subunit indices absent from order_column (empty list when fully occupied); None when order_column is not in the motl.

Raises:
ValueError

If affiliation_column is absent from self.motl.df.

class cryocat.analysis.assembly.TetrahedralComplex(motl, *, affiliation_column='object_id', order_column='geom1', tomo_id_column='tomo_id')#

Bases: PolyhedralComplex

Tetrahedral (T, order 12) motl-analysis complex.

cryocat.analysis.assembly.complex_centers(motl, *, affiliation_column='object_id', tomo_id_column='tomo_id', weights=None)#

Return one barycentric centre particle per (tomogram, object) group.

Parameters:
motlMotlSource

Particle list.

affiliation_columnMotlColumn, default=’object_id’

Column that identifies which object each particle belongs to.

tomo_id_columnMotlColumn, default=’tomo_id’

Column that identifies the tomogram.

weightsArrayLike, optional

Per-particle weights forwarded to geom.barycenter().

Returns:
cryomotl.Motl

One row per (tomo_id, affiliation) pair. tomo_id and object_id carry the group identifiers; all other columns are zero-filled.

cryocat.analysis.assembly.expand_motl(motl, shift_vecs, *, original_id_col='object_id', order_id_col='geom1', sort_vectors=True, orientation='radial', start_index=0)#

Expand a particle list by placing one copy at each shift vector.

For every shift vector v_k and every source particle with rotation R and centre c, a new particle is placed at c + R.apply(v_k). When orientation is "radial", the new particle’s z-axis is aligned with v_k and a random in-plane rotation (phi) is drawn. When orientation is "keep", the original orientation is left unchanged.

Parameters:
motlMotlSource

Input particle list.

shift_vecsnumpy.ndarray, shape (M, 3)

Site displacement vectors in the particle frame (voxels).

original_id_colMotlColumn, default=”object_id”

Column in which to store the source particle’s subtomo_id.

order_id_colMotlColumn, default=”geom1”

Column in which to store the site index (start_index + position).

sort_vectorsbool, default=True

When True, vectors are sorted by (x, y, z) before expansion so that order_id_col values follow a deterministic order.

orientation{“radial”, “keep”}, default=”radial”

"radial" rotates each copy so that its z-axis points along the shift vector and assigns a random phi. "keep" copies the source orientation without modification.

start_indexint, default=0

Value assigned to order_id_col for the first shift vector. Subsequent vectors receive start_index + 1, start_index + 2, …

Returns:
cryomotl.Motl

Expanded particle list with M × N rows (M vectors × N source particles), sorted by (original_id_col, order_id_col) and renumbered.

Raises:
ValueError

If orientation is not a recognised string.

Notes

When orientation="radial" and the shift vector is collinear with the z-axis, the cross-product used to build the rotation is zero. The guard is: ‖v × ẑ‖ < 1e-12 → identity rotation for +z, 180° around x for −z. This avoids NaN Euler angles for shifts along ±z.

cryocat.analysis.assembly.rays_from_motl(motl, pixel_size, reverse_direction=False, ray_length=None)#

Build a ray array from a Motl’s particle positions and orientations.

Coordinates are scaled by pixel_size to convert from voxels to physical units (e.g. nm) without modifying the input Motl. Normal vectors are derived from the particles’ Euler angles via euler_angles_to_normals() (zxz convention, z-axis as the particle’s forward axis).

Parameters:
motlcryomotl.Motl

Particle list. Coordinates must be in voxels.

pixel_sizefloat

Physical units per voxel (e.g. nm/voxel).

reverse_directionbool, default=False

When True, negate the normals before building ray directions so rays fly away from the surface rather than toward it.

ray_lengthfloat, optional

Forwarded to construct_rays(). None gives effectively infinite rays (magnitude 1e10).

Returns:
numpy.ndarray

Shape (N, 6); columns [ox, oy, oz, dx, dy, dz].

cryocat.analysis.assembly.trace_faces(partner, block, site, n_sites)#

Assign face IDs to half-edges by walking the contact graph.

A face is a closed cycle obtained by repeatedly crossing a contact and advancing to the next site counter-clockwise on the partner block. Open paths (hitting an unmatched site or a previously visited half-edge) are labelled −1 (boundary).

Parameters:
partnernumpy.ndarray, shape (H,)

Index of the matching half-edge for each half-edge, or −1 if unmatched. Indices are 0-based positions in the calling array.

blocknumpy.ndarray, shape (H,)

Block ID for each half-edge.

sitenumpy.ndarray, shape (H,)

1-based site index for each half-edge.

n_sitesnumpy.ndarray, shape (H,)

Number of sites on the block of each half-edge.

Returns:
numpy.ndarray, shape (H,), dtype int

face[h] is −1 for boundary half-edges, or a positive integer identifying the closed face that contains h. Face IDs are assigned in traversal order (first closed face found gets id 1, second gets 2, …).

Notes

The traversal rule is: after crossing to half-edge p = partner[h], advance to the next CCW site on p’s block: lookup[(block[p], site[p] % n_sites[p] + 1)].

The “next” map is injective on matched half-edges, so every half-edge lies on exactly one cycle (closed face) or one path (boundary).