Add cryo-EM density rendering, showcase, and release-check

CHANGELOG.md

@@ -2,6 +2,14 @@
 
 ## Unreleased
 
+- Added a `density` subcommand to `pymol_agent.py`: render a model in cryo-EM
+  density (real `--map` or `--simulate`d gaussian density), with the mesh carved
+  around the model or a `--residue` selection and the contour level taken from
+  the map's sigma via `map_info.py`.
+- Added `SHOWCASE.md`: copy-pasteable, reproducible figures and analyses on
+  canonical public structures (every command is a real, tested code path).
+- Added a `make release-check` target: pre-publish hygiene sweep for tracked
+  structures/maps/media/secret files and large files.
 - Added `interface_report.py`: zero-dependency protein-protein interface residue
   analysis between chains (the chain-chain analog of `pocket_report.py`).
 - Added a `pocket` subcommand to `pymol_agent.py` for one-command annotated

Makefile

@@ -1,4 +1,18 @@
-.PHONY: test
+.PHONY: test release-check
 
 test:
 	python3 -m unittest -v
+
+# Pre-publish hygiene sweep: no tracked structures/maps/media/secrets, no large
+# files, no obvious secret strings. Run before pushing public changes.
+release-check:
+	@echo "Checking for tracked data/media/secret files..."
+	@! git ls-files | grep -iE '\.(pdb|cif|mmcif|bcif|ent|map|mrc|mrcs|ccp4|mp4|mov|gif|webm|fasta|fa|fastq|pt|pth|npy|npz|pem|key)$$' | grep -vE '^tests/fixtures/' \
+		|| { echo "ERROR: forbidden files tracked (see above)"; exit 1; }
+	@echo "Checking for large files (>25MB)..."
+	@! find . -path ./.git -prune -o -type f -size +25M -print | grep . \
+		|| { echo "ERROR: large files present (see above)"; exit 1; }
+	@echo "Scanning for secret-like strings..."
+	@! git grep -nIE '(AKIA[0-9A-Z]{16}|ghp_[A-Za-z0-9]{36}|sk-[A-Za-z0-9]{20,}|xox[baprs]-|BEGIN (RSA|OPENSSH|EC|PGP) PRIVATE KEY)' -- . ':(exclude)Makefile' \
+		|| { echo "ERROR: possible secrets (see above)"; exit 1; }
+	@echo "release-check passed."

README.md

@@ -51,6 +51,7 @@ Rosetta-oriented protein design guidance without building a custom plugin.
 - Zero-dependency protein-protein interface residue analysis between chains
 - HELIX-record injection for CA-only backbones (RFdiffusion / Genie designs) so cartoons render correctly
 - MRC/CCP4 map inspection with sigma-based contour-level suggestions
+- Cryo-EM density-fit figures (model in carved density mesh; works with a real map or simulated density)
 - AlphaFold DB fetch with confidence interpretation and pLDDT coloring
 - RCSB PDB fetch for experimental coordinates, metadata, and biological assembly mmCIF
 - UniProt lookup for resolving gene/protein names before AlphaFold fetches
@@ -173,8 +174,13 @@ python3 scripts/chimerax_agent.py align reference.pdb mobile.pdb   # structure a
 python3 scripts/chimerax_rest.py render structure.pdb out.png      # GPU render via managed REST
 python3 scripts/add_helix_records.py model.pdb --json              # fix CA-only backbone cartoons
 python3 scripts/map_info.py map.mrc --json                         # cryo-EM contour levels
+python3 scripts/pymol_agent.py density model.pdb fit.png --simulate # model in (simulated) density
 ```
 
+See [`SHOWCASE.md`](SHOWCASE.md) for end-to-end, copy-pasteable examples on
+canonical public structures (HIV protease, ACE2–RBD, hemoglobin, AlphaFold p53,
+and more).
+
 ## Layout
 
 ```text
@@ -223,7 +229,8 @@ Helper scripts emit machine-readable JSON, with human-readable text as a fallbac
 ## Test
 
 ```bash
-make test
+make test            # run the unit tests
+make release-check   # pre-publish hygiene sweep (no tracked structures/maps/media/secrets, no large files)
 ```
 
 ## Safety And Privacy

SHOWCASE.md

@@ -0,0 +1,111 @@
+# Proteus Showcase
+
+Reproducible figures and analyses you can generate with Proteus, using only
+public structures. Each entry lists the exact commands and what you should see.
+
+Outputs (`.pdb`, `.png`, `.mp4`) are git-ignored, so running these from the repo
+root never dirties the repository. Every command here is a real code path in the
+helper scripts, so this doubles as a set of living examples.
+
+> Requires Python 3.10+. Rendering entries need PyMOL (or ChimeraX for the GPU
+> path); analysis entries are stdlib-only. Run `python3 scripts/proteus_doctor.py`
+> to see what's available on your machine. Turntable movies additionally need
+> `ffmpeg`.
+
+## 1. Binding pocket — HIV-1 protease + indinavir (1HSG)
+
+```bash
+python3 scripts/fetch_pdb.py 1HSG --format pdb
+python3 scripts/pymol_agent.py pocket 1HSG.pdb hiv_protease_pocket.png --label
+```
+
+**Shows:** the indinavir inhibitor nested in the protease active site, with the
+surrounding pocket residues as labeled sticks and polar contacts dashed. The
+JSON reports the ligand atom count and number of pocket residues.
+
+## 2. Protein–protein interface — SARS-CoV-2 RBD + ACE2 (6M0J)
+
+```bash
+python3 scripts/fetch_pdb.py 6M0J --format pdb
+python3 scripts/interface_report.py 6M0J.pdb --chains A,E --cutoff 4.5 --json
+python3 scripts/pymol_agent.py render 6M0J.pdb ace2_rbd.png --color chain
+```
+
+**Shows:** the interface residues on each side of the ACE2 (A) / spike RBD (E)
+contact — the hotspot of SARS-CoV-2 recognition — plus a chain-colored render of
+the complex.
+
+## 3. Quaternary structure — hemoglobin tetramer (4HHB)
+
+```bash
+python3 scripts/fetch_pdb.py 4HHB --format pdb
+python3 scripts/pymol_agent.py render 4HHB.pdb hemoglobin.png --color chain
+python3 scripts/pymol_agent.py spin 4HHB.pdb hemoglobin_spin.mp4 --frames 60 --color chain
+```
+
+**Shows:** the four globin subunits in distinct colors, and a 360° turntable
+movie of the assembly.
+
+## 4. Protein–DNA assembly — nucleosome core particle (1AOI)
+
+```bash
+python3 scripts/fetch_pdb.py 1AOI --format pdb
+python3 scripts/pymol_agent.py render 1AOI.pdb nucleosome.png --color chain
+```
+
+**Shows:** the histone octamer wrapped by duplex DNA, each chain colored
+separately.
+
+## 5. AlphaFold confidence — human p53 (P04637)
+
+```bash
+python3 scripts/uniprot_lookup.py TP53 --gene-exact --json
+python3 scripts/fetch_alphafold.py P04637 --pae --json
+python3 scripts/pae_report.py AF-P04637-F1_pae.json --json
+python3 scripts/pymol_agent.py render AF-P04637-F1.pdb p53_confidence.png --color plddt
+```
+
+**Shows:** the prediction colored by the official AlphaFold pLDDT scheme (blue
+high-confidence core, orange low-confidence disordered regions), with the PAE
+report flagging the inter-domain uncertainty.
+
+## 6. Predicted vs experimental / conformational change — adenylate kinase
+
+```bash
+python3 scripts/fetch_pdb.py 4AKE --format pdb   # open conformation
+python3 scripts/fetch_pdb.py 1AKE --format pdb   # closed conformation
+python3 scripts/compare_structures.py 4AKE.pdb 1AKE.pdb --per-residue --json
+```
+
+**Shows:** the CE-alignment RMSD plus per-residue deviations that pinpoint the
+LID and NMP domain motions between the open and closed states.
+
+## 7. Cryo-EM density fit — sidechain in density
+
+```bash
+# Simulated density (no map download): does this sidechain sit in its density?
+python3 scripts/pymol_agent.py density 1HSG.pdb asp25_density.png \
+    --simulate --residue "chain A and resi 25"
+
+# With a real map, the contour level is chosen from the map's sigma:
+python3 scripts/map_info.py your_map.mrc --json
+python3 scripts/pymol_agent.py density model.pdb fit.png --map your_map.mrc --residue "chain A and resi 25"
+```
+
+**Shows:** the selected residue as sticks carved inside its density mesh, with
+the rest of the model as transparent cartoon for context.
+
+## 8. Publication presets and illustration style
+
+```bash
+python3 scripts/pymol_agent.py render 4HHB.pdb hb_soft.png --preset soft --color chain
+python3 scripts/pymol_agent.py render 1HSG.pdb hiv_illustration.png --preset illustration
+```
+
+**Shows:** the same structures in the neutral-background `soft` preset and the
+outlined `illustration` (molecular-illustration) style.
+
+---
+
+For the full command surface, see [`SKILL.md`](SKILL.md). For tool-specific
+recipes and gotchas, see [`references/`](references/).

SKILL.md

@@ -106,7 +106,7 @@ debugging, patching, or the help text is insufficient for the task.
 | `scripts/pocket_report.py` | Zero-dep ligand pocket contacts from PDB/PDB ID | `python3 scripts/pocket_report.py 1HSG --json` |
 | `scripts/interface_report.py` | Zero-dep protein-protein interface residues between chains | `python3 scripts/interface_report.py 1BRS --chains A,D --json` |
 | `scripts/compare_structures.py` | PyMOL CE alignment + optional per-residue deviations | `python3 scripts/compare_structures.py ref.pdb mobile.pdb --json` |
-| `scripts/pymol_agent.py` | Headless PyMOL driver (info, render, **pocket figure, spin movie**) | `python3 scripts/pymol_agent.py render structure.pdb out.png --color plddt` |
+| `scripts/pymol_agent.py` | Headless PyMOL driver (info, render, **pocket figure, density fit, spin movie**) | `python3 scripts/pymol_agent.py render structure.pdb out.png --color plddt` |
 | `scripts/chimerax_agent.py` | Headless ChimeraX driver (analysis, `--nogui`) | `python3 scripts/chimerax_agent.py run "open 1ubq; info chains #1"` |
 | `scripts/chimerax_rest.py` | Managed ChimeraX REST GUI render (GPU) + turntable | `python3 scripts/chimerax_rest.py render structure.pdb out.png --color plddt` |
 | `scripts/add_helix_records.py` | Add HELIX records to CA-only backbones so cartoons render | `python3 scripts/add_helix_records.py model.pdb --json` |
@@ -349,10 +349,17 @@ python3 scripts/add_helix_records.py model.pdb -o model_ss.pdb --json
 # Then render model_ss.pdb; in PyMOL also: set cartoon_trace_atoms, 1
 ```
 
-### Cryo-EM Contour Level
+### Cryo-EM Density Fit
 ```bash
-# Sigma-based level for `volume`/`isomesh` (absolute level differs per map).
+# Sigma-based contour level (absolute level differs per map).
 python3 scripts/map_info.py map.mrc --json   # -> suggested_level at 1/1.5/2/3 sigma
+
+# Render a model in density — the mesh is carved around the model, which avoids
+# the whole-map contour stall in headless PyMOL (gotcha 20).
+python3 scripts/pymol_agent.py density model.pdb fit.png --map map.mrc --residue "chain A and resi 25"
+
+# No map yet? Simulate gaussian density from the model (predicted/designed structures).
+python3 scripts/pymol_agent.py density model.pdb fit.png --simulate
 ```
 
 ## Good Demo Proteins
@@ -414,6 +421,7 @@ For multi-step workflows, write a summary JSON report at the end with:
 | Render via ChimeraX GPU (REST) | `python3 scripts/chimerax_rest.py render file.pdb out.png` |
 | Fix a CA-only backbone for cartoons | `python3 scripts/add_helix_records.py model.pdb` |
 | Pick a cryo-EM contour level | `python3 scripts/map_info.py map.mrc --json` |
+| Render a model in cryo-EM density | `python3 scripts/pymol_agent.py density model.pdb out.png --map map.mrc` |
 | Fetch an AlphaFold prediction | `python3 scripts/fetch_alphafold.py UNIPROT_ID --pae --json` |
 | Align two structures (ChimeraX) | `python3 scripts/chimerax_agent.py align ref.pdb mobile.pdb` |
 | Measure SASA | `python3 scripts/chimerax_agent.py sasa file.pdb` |

references/pymol.md

@@ -316,6 +316,14 @@ and can hang in headless PyMOL (the headless build lacks the VTKm accelerator) 
 do whole-map surfaces in ChimeraX instead. Use `scripts/map_info.py` to choose a
 sigma-based contour level (absolute levels differ per map).
 
+`pymol_agent.py density` packages this: `density model.pdb out.png --map map.mrc
+--residue "chain A and resi 25"` carves the mesh around the residue and zooms on
+it; `--simulate` generates the gaussian density when there is no experimental map.
+
+> Object naming: `model` is a reserved word in PyMOL's selection language —
+> loading an object as `model` makes later selections fail with "Invalid
+> selection". Use a name like `struct` or `mol` instead.
+
 ## Publication Settings
 
 Standard block for high-quality renders:

scripts/pymol_agent.py

@@ -446,6 +446,94 @@ def render_pocket(pdb_path: str, output: str, ligand: str = "organic", radius: f
     return _verify_png(result, abs_out)
 
 
+def render_density(pdb_path: str, output: str, map_path: str = None, simulate: bool = False,
+                   level: float = None, n_sigma: float = 2.0, carve: float = 2.5,
+                   residue: str = None, color: str = "chain", width: int = 1200,
+                   height: int = 900, preset: str = "publication") -> dict:
+    """Render a model fitted in cryo-EM density (real map or simulated).
+
+    The density mesh is carved around the model (or around `residue`), which
+    avoids the whole-map contour stall in headless PyMOL. With `residue`, shows
+    that selection as sticks and zooms on it — the "is this sidechain supported
+    by density?" figure. Without a map, `simulate` generates a gaussian density
+    from the model so the same path works for predicted/designed structures.
+
+    The contour level defaults to mean + n_sigma * sigma from the map (via
+    map_info), or 1.0 for simulated density.
+    """
+    # Validate inputs before checking for the tool, so argument errors are
+    # deterministic whether or not PyMOL is installed.
+    if not os.path.isfile(pdb_path):
+        return {"status": "error", "error": f"Structure file not found: {pdb_path}"}
+    if not simulate and not map_path:
+        return {"status": "error", "error": "Provide --map, or use --simulate to generate density from the model."}
+    if map_path and not os.path.isfile(map_path):
+        return {"status": "error", "error": f"Map file not found: {map_path}"}
+    if not PYMOL:
+        return {"status": "error", "error": "PyMOL not found. Install it or set PYMOL_BIN."}
+    try:
+        color = _validate_pymol_color(color)
+    except ValueError as exc:
+        return {"status": "error", "error": str(exc)}
+
+    # Contour level: explicit > sigma-from-map > simulated default.
+    auto_level = None
+    if level is None:
+        if simulate:
+            level = 1.0
+        else:
+            try:
+                import map_info
+                _, _, mean, sigma = map_info.read_map_stats(map_path)
+                level = round(mean + n_sigma * sigma, 4)
+                auto_level = {"mean": round(mean, 4), "sigma": round(sigma, 4), "n_sigma": n_sigma}
+            except Exception:
+                level = 1.0
+
+    abs_pdb = os.path.abspath(pdb_path)
+    abs_out = os.path.abspath(output)
+    sel = "(%s)" % residue if residue else None
+    carve_target = sel if sel else "struct"
+
+    lines = [f'cmd.load({_py_literal(abs_pdb)}, "struct")']
+    if simulate:
+        lines.append('cmd.map_new("emap", "gaussian", 1.0, "struct", 5)')
+    else:
+        lines.append(f'cmd.load({_py_literal(os.path.abspath(map_path))}, "emap")')
+    lines += [
+        f'cmd.isomesh("dens", "emap", {level}, {_py_literal(carve_target)}, carve={carve})',
+        'cmd.color("gray70", "dens")',
+        'cmd.hide("everything", "struct")',
+        'cmd.show("cartoon", "struct")',
+        _color_script(color, "struct"),
+    ]
+    if sel:
+        lines += [
+            'cmd.set("cartoon_transparency", 0.6, "struct")',
+            f'cmd.show("sticks", {_py_literal(sel)})',
+            f'util.cnc({_py_literal(sel)})',
+        ]
+    lines.append(_preset_script(preset))
+    if sel:
+        lines += [f'cmd.orient({_py_literal(sel)})', f'cmd.zoom({_py_literal(sel)}, 6)']
+    else:
+        lines.append('cmd.orient("struct")')
+    lines += [
+        f'cmd.ray({width}, {height})',
+        f'cmd.png({_py_literal(abs_out)})',
+        f'_output["data"]["rendered"] = {_py_literal(abs_out)}',
+    ]
+    result = run_pymol_script("\n".join(lines), timeout=300)
+    result = _verify_png(result, abs_out)
+    if result.get("status") == "ok":
+        result["data"]["level"] = level
+        result["data"]["carve"] = carve
+        result["data"]["simulated"] = simulate
+        if auto_level:
+            result["data"]["auto_level"] = auto_level
+    return result
+
+
 def main():
     parser = argparse.ArgumentParser(
         description="PyMOL headless agent helper — run commands, inspect structures, render images.",
@@ -488,6 +576,21 @@ def main():
     p_pocket.add_argument("--label", action="store_true", help="Label pocket residues (resn+resi)")
     p_pocket.add_argument("--preset", default="publication", choices=["publication", "illustration", "soft"])
 
+    # density
+    p_density = sub.add_parser("density", help="Render a model in cryo-EM density (real map or --simulate)")
+    p_density.add_argument("pdb", help="Path to model (PDB/CIF)")
+    p_density.add_argument("output", nargs="?", default="/tmp/pymol_density.png", help="Output PNG path")
+    p_density.add_argument("--map", dest="map_path", help="MRC/CCP4 density map (omit with --simulate)")
+    p_density.add_argument("--simulate", action="store_true", help="Simulate gaussian density from the model")
+    p_density.add_argument("--level", type=float, help="Contour level (default: map sigma, or 1.0 simulated)")
+    p_density.add_argument("--n-sigma", type=float, default=2.0, help="Sigma multiple for auto level (default: 2)")
+    p_density.add_argument("--carve", type=float, default=2.5, help="Carve radius around the model in Angstroms")
+    p_density.add_argument("--residue", help="Residue selection to show as sticks and zoom on")
+    p_density.add_argument("--color", default="chain", help="Model color: spectrum, chain, bfactor, plddt, or a color")
+    p_density.add_argument("--width", type=int, default=1200)
+    p_density.add_argument("--height", type=int, default=900)
+    p_density.add_argument("--preset", default="publication", choices=["publication", "illustration", "soft"])
+
     # spin
     p_spin = sub.add_parser("spin", help="Render a 360-degree turntable movie (frames -> ffmpeg)")
     p_spin.add_argument("pdb", help="Path to structure file")
@@ -512,6 +615,10 @@ def main():
     elif args.command == "pocket":
         result = render_pocket(args.pdb, args.output, args.ligand, args.radius,
                                args.width, args.height, args.label, args.preset)
+    elif args.command == "density":
+        result = render_density(args.pdb, args.output, args.map_path, args.simulate,
+                                args.level, args.n_sigma, args.carve, args.residue,
+                                args.color, args.width, args.height, args.preset)
     elif args.command == "spin":
         result = render_spin(args.pdb, args.output, args.frames, args.width, args.height,
                              args.style, args.color, args.preset, args.fps)

tests/test_scripts.py

@@ -171,6 +171,20 @@ def test_pymol_pocket_missing_file(self):
         result = pymol_agent.render_pocket("tests/fixtures/does_not_exist.pdb", "out.png")
         self.assertEqual(result["status"], "error")
 
+    def test_pymol_density_requires_map_or_simulate(self):
+        # Missing model file
+        self.assertEqual(
+            pymol_agent.render_density("tests/fixtures/missing.pdb", "out.png", simulate=True)["status"],
+            "error")
+        # Neither --map nor --simulate
+        result = pymol_agent.render_density("tests/fixtures/tiny.pdb", "out.png")
+        self.assertEqual(result["status"], "error")
+        self.assertIn("simulate", result["error"])
+        # Map path that does not exist
+        self.assertEqual(
+            pymol_agent.render_density("tests/fixtures/tiny.pdb", "out.png", map_path="nope.mrc")["status"],
+            "error")
+
     def test_map_info_sigma_from_synthetic_mrc(self):
         nx = ny = nz = 4
         vals = [float(i % 7) for i in range(nx * ny * nz)]