Address code-review findings on the 2.x PRs

CHANGELOG.md

@@ -2,6 +2,19 @@
 
 All notable changes to DBSTEP. Versions follow [semantic versioning](https://semver.org/).
 
+## Unreleased
+
+### Fixed
+- SDF energy fields holding plain integers (a relative energy of `0`, rounded values) are accepted for `--boltzmann`; the integer guard now applies only to xyz comment lines.
+- `--nometals` removes metals through the same path as `--noH`: spec atoms are renumbered (a metal chosen as atom1 becomes a zero-radius ghost, keeping the Sterimol alignment) and per-atom metadata stays aligned, so `--decompose --nometals` works.
+- `--decompose` with a scan starting at R = 0 paired contributions with the wrong radius rows.
+- `--tensor --save` on multi-frame files or `--residue all` wrote every result to the same `.npy`; the frame index and residue are now part of the filename.
+- `all_frames(frames=...)` no longer changes the caller's options object; `frames=0` from Python selects frame 0 instead of all frames.
+- `--residue all --atom1 N` uses the named atom to pick residues (previously only `--atom` was consulted); a default atom2 can no longer coincide with atom1 (e.g. `--atom CB`), and an explicit coincidence is an error.
+
+### Added
+- `path` column (the input path as given) in `results`, the CSV and the contributions CSV, so same-named files from different folders stay distinguishable.
+
 ## 2.1.0 — 2026-09-25
 
 ### Added

dbstep/Dbstep.py

@@ -123,6 +123,18 @@ def __init__(self, *args, **kwargs):
 						print("   Cutoff {:.2f} Ang around atom1: keeping {} of {} atoms".format(self.cutoff, self.n_atoms_kept, self.n_atoms_total))
 			elif not options.quiet:
 				print("   Note: --cutoff is ignored for density cube input (the grid comes from the cube file)")
+		# --nometals: drop metal atoms through the same ghost/renumber path as --noH, so spec atoms
+		# and per-atom metadata stay consistent (a metal chosen as atom1 becomes a zero-radius ghost)
+		if options.no_metals and options.surface == "vdw":
+			is_metal = np.isin(mol.ATOMTYPES, list(metals))
+			if is_metal.any():
+				spec = [options.spec_atom_1] + list(options.spec_atom_2) + ([int(options.atom3)] if options.atom3 else [])
+				new_spec = mol.exclude_mask(is_metal, spec)
+				options.spec_atom_1, options.spec_atom_2 = new_spec[0], new_spec[1:1 + len(options.spec_atom_2)]
+				if options.atom3:
+					options.atom3 = new_spec[-1]
+				if options.verbose:
+					print("   Excluded {} metal atom(s) (--nometals)".format(int(is_metal.sum())))
 		# the structure actually measured (after selection and crop), before translation/rotation
 		self.atoms, self.coords = np.array(mol.ATOMTYPES), np.array(mol.CARTESIANS)
 		self.spec_atoms = [options.spec_atom_1] + list(options.spec_atom_2)
@@ -162,6 +174,11 @@ def __init__(self, *args, **kwargs):
 			}
 			if options.save:
 				save_base = name if isinstance(name, str) else "tensor"
+				structure_idx = getattr(options, "structure", None)
+				if structure_idx is not None:
+					save_base += "_frame{}".format(structure_idx)
+				if self.residue_label:
+					save_base += "_" + self.residue_label.replace(":", "").replace(" ", "_")
 				np.save(save_base + "_tensor.npy", self.tensor)
 				if not options.quiet:
 					print("   Tensor saved to {}_tensor.npy (shape: {})".format(save_base, self.tensor.shape))
@@ -217,14 +234,6 @@ def _assign_surface(self, mol, file, options, origin):
 			if options.sterimol or options.volume:
 				mol.CARTESIANS = calculator.translate_mol(mol, options, origin)
 
-			# Remove metals when --nometals is specified (iterate in reverse to avoid index shifting)
-			if options.no_metals:
-				for i in range(len(mol.ATOMTYPES) - 1, -1, -1):
-					if mol.ATOMTYPES[i] in metals:
-						mol.ATOMTYPES = np.delete(mol.ATOMTYPES, i)
-						mol.CARTESIANS = np.delete(mol.CARTESIANS, i, axis=0)
-						mol.RADII = np.delete(mol.RADII, i)
-
 			# Determine grid bounds from molecule extent
 			[x_min, x_max, y_min, y_max, z_min, z_max, xyz_max] = sterics.max_dim(mol.CARTESIANS, mol.RADII, options)
 			if options.gridsize and options.verbose:
@@ -399,6 +408,8 @@ def _compute(self, mol, file, options, origin, occ_grid, occ_vol, point_tree, gr
 			if options.volume:
 				if rad == 0:
 					bur_vol, bur_shell = 0.0, 0.0
+					if options.decompose:
+						contributions_list.append({})
 				else:
 					if options.vshell:
 						strip_width = options.vshell
@@ -423,6 +434,7 @@ def _compute(self, mol, file, options, origin, occ_grid, occ_vol, point_tree, gr
 			# Record the result row (also written by --csv)
 			self.results.append({
 				"file": plain_name,
+				"path": file if isinstance(file, str) else "",
 				"frame": structure_idx if structure_idx is not None else "",
 				"structure": structure_label,
 				"residue": self.residue_label or "",
@@ -627,7 +639,7 @@ def contribution_rows(run):
 	for result, contributions in zip(run.results, per_radius):
 		for label, percent in contributions.items():
 			rows.append({"file": result["file"], "frame": result["frame"], "structure": result["structure"], "residue": result["residue"],
-						 "radius": result["radius"], "contributor": label, "percent_vbur": percent})
+						 "radius": result["radius"], "contributor": label, "percent_vbur": percent, "path": result.get("path", "")})
 	return rows
 
 
@@ -736,9 +748,13 @@ def all_frames(file, frames=None, **kwargs):
 		list of dbstep objects, one per frame (times residues for residue="all")
 	"""
 	options = kwargs["options"] if "options" in kwargs else set_options(kwargs)
+	previous = getattr(options, "frames", False)
 	if frames is not None:
 		options.frames = frames
-	return run_file(file, options)
+	try:
+		return run_file(file, options)
+	finally:
+		options.frames = previous
 
 
 def build_parser():

dbstep/ensemble.py

@@ -44,6 +44,14 @@ def number_in_text(text):
 	return float(match.group(0)) if match else None
 
 
+def field_number(text):
+	"""Number in an SDF data field: the whole value as a float (so "0" and "12" work), else the first float in the text."""
+	try:
+		return float(str(text).strip())
+	except ValueError:
+		return number_in_text(text)
+
+
 def _label(run):
 	return run.structure_name or run.results[0]["file"] if run.results else str(run.file)
 
@@ -61,13 +69,13 @@ def energy_of(run, tag=None):
 		key = str(tag).strip().lower()
 		if key not in lowered:
 			sys.exit(f"   No <{tag}> data field for structure '{_label(run)}'. Fields present: {', '.join(properties) or 'none'}")
-		value = number_in_text(lowered[key])
+		value = field_number(lowered[key])
 		if value is None:
 			sys.exit(f"   Data field <{tag}> of structure '{_label(run)}' is not a number: {lowered[key]!r}")
 		return value
 	for candidate in ENERGY_TAGS:
 		if candidate in lowered:
-			value = number_in_text(lowered[candidate])
+			value = field_number(lowered[candidate])
 			if value is not None:
 				return value
 	if "comment" in lowered:
@@ -115,6 +123,7 @@ def boltzmann_average(runs, tag=None, temperature=298.15, units="kcal"):
 		first = runs[0].results[i]
 		row = {
 			"file": first["file"],
+			"path": first.get("path", ""),
 			"frame": "",
 			"structure": "boltzmann",
 			"residue": first["residue"],

dbstep/selection.py

@@ -156,11 +156,12 @@ def find_atom_by_name(mol, mask, name):
 	return int(hits[0]) + 1
 
 
-def resolve_spec_atom(mol, mask, value, default_names, what):
+def resolve_spec_atom(mol, mask, value, default_names, what, avoid=None):
 	"""Turn an atom given by name ('CA') or by 1-indexed file position ('12'/12) into a 1-indexed position.
 
 	Names are looked up within the selected residue(s). With `value` unset, the first of
-	`default_names` present in the residue is used. Returns (index, name_used).
+	`default_names` present in the residue (and different from the 1-indexed atom `avoid`) is used.
+	Returns (index, name_used).
 	"""
 	if value not in (False, None, ""):
 		try:
@@ -176,7 +177,7 @@ def resolve_spec_atom(mol, mask, value, default_names, what):
 		return index, mol.METADATA["name"][index - 1]
 	for name in default_names:
 		index = find_atom_by_name(mol, mask, name)
-		if index is not None:
+		if index is not None and index != avoid:
 			return index, name
 	present = ", ".join(dict.fromkeys(mol.METADATA["name"][mask]))
 	sys.exit(f"   Could not pick a default atom for {what} (tried {', '.join(default_names)}). Atoms present: {present}")
@@ -240,7 +241,9 @@ def apply_residue_selection(mol, options, verbose=False):
 		atom2_values = [atom2_values]
 	atom2, atom2_names = [], []
 	for value in atom2_values:
-		index, name = resolve_spec_atom(mol, first_mask, value, _DEFAULT_ATOM2, "atom2")
+		index, name = resolve_spec_atom(mol, first_mask, value, _DEFAULT_ATOM2, "atom2", avoid=atom1)
+		if index == atom1:
+			sys.exit(f"   atom2 ({name}) is the same atom as atom1 ({atom1_name}); the Sterimol axis needs two different atoms")
 		atom2.append(index)
 		atom2_names.append(name)
 	if options.sterimol and atom2_names == ["N"] and verbose:
@@ -302,13 +305,16 @@ def list_residues(mol, options):
 	peptide backbone) that contain the atom1 name, in file order; waters and other hetero groups
 	are skipped and --chain restricts the chains."""
 	_require_metadata(mol, "--residue all")
-	atom_name = options.atom if options.atom not in (False, None, "") else "CA"
+	# atom1 may be given as --atom or as --atom1; either must be a name in all-residue mode
+	atom_name = options.atom if options.atom not in (False, None, "") else options.spec_atom_1
+	if atom_name in (False, None, ""):
+		atom_name = "CA"
 	try:
 		int(atom_name)
 	except (TypeError, ValueError):
 		pass
 	else:
-		sys.exit("   With --residue all, --atom must be an atom name such as CA, not a file index")
+		sys.exit("   With --residue all, atom1 must be an atom name such as CA (--atom CA), not a file index")
 	meta = mol.METADATA
 	names = np.char.upper(meta["name"].astype(str))
 	polymer = (meta["record"] == "ATOM") | polymer_hetatm_mask(mol)

dbstep/trajectory.py

@@ -43,7 +43,7 @@ def parse_frames(spec, n_frames):
 	Returns:
 		list of 0-based frame indices in the order given
 	"""
-	if spec in (False, None, ""):
+	if spec is False or spec is None or (isinstance(spec, str) and spec.strip() == ""):
 		return list(range(n_frames))
 	tokens = [str(t).strip() for t in spec] if isinstance(spec, (list, tuple)) else [t.strip() for t in str(spec).split(",")]
 	indices = []
@@ -75,6 +75,7 @@ def frame_indices(file, options):
 	"""Frames to run for a file: a list of 0-based indices, or [None] for a single-structure file without --frames."""
 	n_frames = count_frames(file)
 	frames = getattr(options, "frames", False)
-	if n_frames == 1 and not frames:
+	unset = frames is False or frames is None or (isinstance(frames, str) and frames.strip() == "")
+	if n_frames == 1 and unset:
 		return [None]
 	return parse_frames(frames, n_frames)

dbstep/writer.py

@@ -277,7 +277,7 @@ def xyz_export(file, mol):
 	log.close()
 
 
-CSV_COLUMNS = ["file", "frame", "structure", "residue", "atom1", "atom2", "radius", "mol_vol", "percent_vbur", "percent_sbur", "bmin", "bmax", "L", "population"]
+CSV_COLUMNS = ["file", "frame", "structure", "residue", "atom1", "atom2", "radius", "mol_vol", "percent_vbur", "percent_sbur", "bmin", "bmax", "L", "population", "path"]
 
 
 def csv_export(path, rows):
@@ -289,7 +289,7 @@ def csv_export(path, rows):
 			out.writerow({key: row.get(key, "") for key in CSV_COLUMNS})
 
 
-CONTRIBUTION_COLUMNS = ["file", "frame", "structure", "residue", "radius", "contributor", "percent_vbur"]
+CONTRIBUTION_COLUMNS = ["file", "frame", "structure", "residue", "radius", "contributor", "percent_vbur", "path"]
 
 
 def contributions_csv_export(path, rows):

docs/plans/void-featurization.md

@@ -39,10 +39,12 @@ describe the shape of the void as a whole.
 
 ### Boundary
 
-- **Sphere-bounded (consistent with %V_bur):** empty voxels within radius R that are connected
-  to the centre by a flood fill through empty space (`scipy.ndimage.label`). Disconnected gaps
-  behind residues are excluded. Free volume is the complement of what `--exclude-self --vbur`
-  reports today.
+- **Sphere-bounded:** empty voxels within radius R that are connected to the centre by a flood
+  fill through empty space (`scipy.ndimage.label`). Disconnected gaps behind residues are
+  excluded, so this is the *connected-pocket* volume. It is smaller than or equal to the plain
+  free volume, the complement of what `--exclude-self --vbur` reports today, which counts every
+  empty voxel in the sphere including gaps the centre cannot reach. Report both: their difference
+  is itself informative (buried voids near the site).
 - **Morphologically closed (fpocket-like):** dilate the occupancy with a large probe (about
   8 Å), erode back, which seals the mouth; then flood-fill from the centre inside the sealed
   envelope. Pocket volume becomes independent of R.

tests/test_decompose.py

@@ -92,7 +92,7 @@ def test_cli_prints_and_writes_contributions(tmp_path):
 	assert f"Residue contributions written to {contributions}" in text
 	with open(contributions, newline="") as f:
 		rows = list(csv.DictReader(f))
-	assert list(rows[0].keys()) == ["file", "frame", "structure", "residue", "radius", "contributor", "percent_vbur"]
+	assert list(rows[0].keys()) == ["file", "frame", "structure", "residue", "radius", "contributor", "percent_vbur", "path"]
 	assert {row["contributor"] for row in rows} >= {"A:3 ALA", "A:2 ALA", "A:4 ALA", "A:101 HOH"}
 	table = float(re.search(r"ala5\.pdb A:3 ALA\s+CA\s+3\.50\s+[\d.]+\s+([\d.]+)", text).group(1))
 	assert sum(float(row["percent_vbur"]) for row in rows) == pytest.approx(table, abs=0.006)

tests/test_residue_all.py

@@ -110,7 +110,7 @@ def test_csv_all_residues_matches_table(tmp_path):
 	text = run_cli(ALA5, "--residue", "all", "-b", "-s", "--csv", str(out))
 	rows = read_csv(out)
 	assert len(rows) == 5
-	assert list(rows[0].keys()) == ["file", "frame", "structure", "residue", "atom1", "atom2", "radius", "mol_vol", "percent_vbur", "percent_sbur", "bmin", "bmax", "L", "population"]
+	assert list(rows[0].keys()) == ["file", "frame", "structure", "residue", "atom1", "atom2", "radius", "mol_vol", "percent_vbur", "percent_sbur", "bmin", "bmax", "L", "population", "path"]
 	assert f"Results written to {out}" in text
 	for row in rows:
 		printed = re.search(rf"ala5\.pdb {re.escape(row['residue'])}\s+CA\s+CB\s+3\.50\s+([\d.]+)\s+([\d.]+)\s+[\d.]+\s+([\d.]+)\s+([\d.]+)\s+([\d.]+)", text)