battery: R-free against the deposition without refinement only, one protocol for rotation and stills
The batteries now score the deposited model against a merge one way: model_check.py, REFMAC's first cycle before any shift (bulk solvent and anisotropic scale by default), the depositor's free set, d_min_used = max(ours, deposited), compared against the deposited data through the same protocol. - model_check.py: the restrained-refinement protocol is removed (refine, refined_shared, rename_misfit_ligands, REFINE_CYCLES, the refined parameter, --no-refine). Both data sets are now scored on their common reflections only, so a low-resolution or anisotropic cut of either applies to both (measured on three still sets: <= 0.001 change). - battery.py / paper_plots.py: refmac_refined_* and refmac_refined_ratio are gone; refmac_rfree_ratio is the headline; paper_plots draws the first-cycle pair only (--rfree removed). - still_battery.py: dimple is replaced by model_check.check on the same p.mtz. R_free rugnux is rfree_depflags, the target is rfree_depdata (no longer min(dimple, published)); the published R-free stays as an information column. The dimple cache, free-flag writer, deposited-MTZ writer and the reference verb are removed. This also retires dimple's freerflag COMPLETE renumbering, which scored R-free on depositor work reflections whenever our merge had hkl the deposition lacks (model_check carries the depositor's flags per hkl and marks ours-only hkl as work, and with the common-hkl restriction there are none). `report RUN --out DIR` rescores an older run's merges without rerunning rugnux. Refinement cost about 6.4 h per full rotation battery; the 0-cycle check takes seconds per set. Verified: 9V3D still merge reproduces the 0-cycle common-hkl numbers (0.2461 ours, 0.2367 deposited). Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_015CeEGgpKdfkCTJb3YNMJcu
This commit is contained in:
1 parent
ccf81183cc
commit
4d8c2cc031
6 files changed
+93
-355
No files matched your search
+22
-57
@@ -529,19 +529,22 @@ with the deposition, that needs CCP4: `tools/battery/model_check.py` scores each
|
||||
(`p.mtz`) against the deposited model, unmodified (ligands, waters and hydrogens kept; only atoms of
|
||||
unknown element, UNX, are left out because REFMAC stops on them): REFMAC's rigid-body mode,
|
||||
R-factors taken before any shift, the data reindexed into the model's setting and reduced to the
|
||||
deposited space group. The depositor's own structure factors go through the same setting choice
|
||||
(they can sit in the other branch of a merohedral ambiguity) and are scored by the same protocol.
|
||||
deposited space group. The model is never refined. The depositor's own structure factors go through
|
||||
the same setting choice (they can sit in the other branch of a merohedral ambiguity) and are scored
|
||||
by the same protocol.
|
||||
|
||||
All REFMAC numbers use the same high-resolution limit, `d_min_used` = max(our d_min, the deposited
|
||||
d_min), i.e. the lower of the two resolutions, so neither data set is scored past what the other
|
||||
reaches. The four R-free fields, all of the same deposited model:
|
||||
Both data sets are scored on the same reflections: the hkl present in both, at `d_min_used` =
|
||||
max(our d_min, the deposited d_min), i.e. the lower of the two resolutions, so neither data set is
|
||||
scored past what the other reaches and a cut either one made (low-resolution limit, anisotropic
|
||||
cut) applies to both. Bulk solvent and the overall anisotropic scale are REFMAC's defaults, fitted
|
||||
to each data set. The four R-free fields, all of the same deposited model:
|
||||
|
||||
| field | data | free set | reads as |
|
||||
|---|---|---|---|
|
||||
| `refmac_rfree` | ours | ours | close to a work R: most of our free reflections were work reflections to the depositor's refinement |
|
||||
| `refmac_rfree_depflags` | ours | the depositor's | **the fair number for our data** |
|
||||
| `refmac_rfree_depdata` | the depositor's | the depositor's | the baseline: same model, same protocol, the data the model was refined against |
|
||||
| `refmac_rfree_ratio` | | | `refmac_rfree_depflags / refmac_rfree_depdata`; below 1 = our data fit the model better; the report plots it |
|
||||
| `refmac_rfree_ratio` | | | `refmac_rfree_depflags / refmac_rfree_depdata`: **the headline**; below 1 = our data fit the model better; the report plots it |
|
||||
|
||||
`refmac_rwork` is R on our own work reflections. The depositor-data numbers carry a home
|
||||
advantage (the model was refined against those data), and `refmac_rfree_depdata` is not the
|
||||
@@ -560,39 +563,6 @@ set (`refmac_reason` says why). The runner records `refmac_rfree`, `refmac_rwork
|
||||
A `--model-check` run records the REFMAC version in its manifest (`refmac`), as does every recheck,
|
||||
and the report lists it.
|
||||
|
||||
**Refined, shared free set (second protocol).** The first-cycle numbers above score the model as
|
||||
deposited at the common resolution limit. The second protocol refines the same deposited model
|
||||
against each data set at that data set's OWN resolution and compares R-free on a free set both
|
||||
share:
|
||||
|
||||
- same starting model, same REFMAC protocol for both data sets: restrained refinement, 10 cycles,
|
||||
automatic weight, isotropic B for every atom (a deposited ANISOU is often TLS-derived, and
|
||||
refining it atom by atom diverged), twin refinement where the entry declares twinning; a ligand
|
||||
whose code the monomer library describes with other atoms is renamed to a code the library
|
||||
lacks, so REFMAC makes its restraints from the model's own coordinates (`ligands_renamed` in
|
||||
the check's output);
|
||||
- each data set over its whole resolution range; its work set is every reflection it has except
|
||||
the shared free set;
|
||||
- shared free set = the depositor's free reflections present in BOTH data sets (after the change
|
||||
of basis), which takes the lower of the two limits in every direction, anisotropic cuts and
|
||||
rejected reflections included; depositor free reflections outside it are removed from both;
|
||||
- R-free is REFMAC's own: each data set's free set IS the shared set. It is not recomputed from
|
||||
the output MTZ, whose F carries REFMAC's anisotropic scale - on strongly anisotropic data that
|
||||
scale multiplies the noise of the weak direction and the recomputed R comes out far too high.
|
||||
|
||||
| field | meaning |
|
||||
|---|---|
|
||||
| `refmac_refined_rfree_shared` | our data, refined, R-free on the shared set - **the refined number for our data** |
|
||||
| `refmac_refined_rfree_shared_depdata` | the depositor's data, same protocol, same shared set |
|
||||
| `refmac_refined_ratio` | the ratio of the two; below 1 = our data give the better refined model |
|
||||
| `refmac_refined_rwork`, `refmac_refined_rwork_depdata` | REFMAC's R-work after refinement, each data set (their own work sets and resolutions: not comparable with each other) |
|
||||
| `refmac_refined_rfree_refmac`, `refmac_refined_rfree_refmac_depdata` | REFMAC's own R-free, for reference |
|
||||
| `refmac_shared_free_n` | size of the shared free set |
|
||||
| `refmac_refined_reason` | why the refined protocol gave nothing, when it did not (the first-cycle numbers stand) |
|
||||
|
||||
The refined numbers need the depositor's structure factors with a free set, like the ratio; they
|
||||
take minutes to an hour per set. `model_check.py --no-refine` skips them.
|
||||
|
||||
### Rechecking an existing run
|
||||
|
||||
`battery.py recheck RUN [--only a,b] [--jobs N]` reruns the model check and the depositor-data
|
||||
@@ -758,7 +728,7 @@ was run `--unforced`), so its forced XDS-arm rows are left out the same way.
|
||||
| `rfree`, `rwork`, `cc_model`, `model_fit`, `rfree_deposited`, `rfree_ratio` | any row with a model: R_FREE, R_WORK, CC_MODEL_OVERALL and MODEL_FIT as rugnux reports them (placement-only, own free set: trend fields), the published R-free, and rfree / rfree_deposited |
|
||||
| `anom_sigma` | any row with a model: ANOMALOUS_SCATTERER_MEAN_SIGMA, the anomalous difference map (F+ - F- on the model's phases) read at every model atom from phosphorus up (S of Met/Cys, metals, Cl, I, ...) and averaged, in map sigma - how much anomalous signal the merge carries; null where the merge kept no Bijvoet split or the model has no such atom |
|
||||
| `rmodel`, `rmodel_shell_scaled`, `radial_misfit` | any row with a model: R_MODEL, R_MODEL_SHELL_SCALED and MODEL_RADIAL_MISFIT - R over all the reflections, with and without one free scale per resolution shell, and how big that rescale was. `rmodel_shell_scaled` is the R that compares between runs (below) |
|
||||
| `refmac_rfree`, `refmac_rwork`, `refmac_rfree_depflags`, `refmac_rfree_depdata`, `refmac_rfree_ratio`, `refmac_twin`, `refmac_untwinned`, `refmac_status`, `refmac_reason`, `refmac_refined_*`, `refmac_shared_free_n` | the REFMAC check (`--model-check`, open arm; see its section for which R-free is which) |
|
||||
| `refmac_rfree`, `refmac_rwork`, `refmac_rfree_depflags`, `refmac_rfree_depdata`, `refmac_rfree_ratio`, `refmac_twin`, `refmac_untwinned`, `refmac_status`, `refmac_reason` | the REFMAC check (`--model-check`, open arm; see its section for which R-free is which) |
|
||||
| `dep_status`, `dep_reason`, `dep_kind`, `dep_d_min`, `dep_n_common`, `dep_cc_delta_all`, `dep_cc_delta_outer`, `dep_beyond_cc`, `dep_beyond_d` | open arm with a model: our merge against the depositor's data per shell (`depdata_check.py`, above) |
|
||||
| `shx_*` | small molecule with a `cod` structure: the SHELXL check (`shelx_check.py`, above) |
|
||||
| `wall_s`, `rugnux_wall_s`, `elapsed_s`, `exit_code`, `gpu_others`, `rugnux_verdict` | timing (`wall_s` is rugnux's own WALL_TIME, not the time spent waiting for the GPU queue), the exit code, other GPU processes seen before the set, and rugnux's own verdict line |
|
||||
@@ -911,11 +881,9 @@ on the server (`20261011-1236_549939_1.0.0-rc.175`, and `20261014-0141_6e9b89_1.
|
||||
|
||||
## Figures for publication
|
||||
|
||||
`paper_plots.py RUN [RUN-REDO ...] --out DIR [--time-runs RUN ...] [--uncontended RUN] [--rfree refined|first-cycle]`
|
||||
draws three open-arm figures from a run's `results.json`: R-free of the deposited model after the REFMAC
|
||||
check (`--model-check`) on the depositor's data against Rugnux's data (by default the refined protocol:
|
||||
each data set at its own resolution, R-free on the shared free set; `--rfree first-cycle` for the
|
||||
first-cycle numbers at the common resolution), the deposited high-resolution limit against Rugnux's cut, and a
|
||||
`paper_plots.py RUN [RUN-REDO ...] --out DIR [--time-runs RUN ...] [--uncontended RUN]`
|
||||
draws three open-arm figures from a run's `results.json`: R-free of the deposited model by the REFMAC
|
||||
check (`--model-check`) on the depositor's data against Rugnux's data, the deposited high-resolution limit against Rugnux's cut, and a
|
||||
histogram of the processing time per dataset with the median marked. Later runs fill in sets an
|
||||
earlier one did not run (a `-redo`). `--time-runs` takes the time histogram from other runs - one
|
||||
whose sets had the GPU to themselves - and `summary.txt` lists the provenance of both. Beside the PDF/PNG figures it writes the
|
||||
@@ -949,11 +917,10 @@ once and its row says what the merge is worth and what it cost:
|
||||
| facility, beamline, detector, dataset, PDB | from the manifest |
|
||||
| images | frames in the input files (rugnux's `IMAGES_PROCESSED`); images (deposition) the deposition's own count, where it is known |
|
||||
| indexed %, d_min, R_split | rugnux's report (`INDEXING_RATE`, `INCLUDE_RESOLUTION_RANGE`, `R_SPLIT`: whole frames go to one half, as in CrystFEL) |
|
||||
| R_free rugnux (dimple) | `dimple` of the deposited model against our merge, reindexed into the deposition's space group, with the deposition's own free set, at the deposition's resolution |
|
||||
| R_free deposition (dimple) | the same protocol on the deposited data (the entry's `-sf.cif`: `intensity_meas`, else the mean of the anomalous pair, else `F_meas_au` squared) |
|
||||
| R_free published | the entry's own R-free (dark and ground-state entries only), from its own refinement program |
|
||||
| target | which of the two is lower (`dimple` or `published`): the R-free to beat |
|
||||
| difference, quality | ours minus the target; win where ours is lower by more than 0.01, draw within 0.01, loss otherwise |
|
||||
| R_free rugnux | the deposited model against our merge, not refined: the rotation battery's REFMAC check (`model_check.py`, see above), i.e. `refmac_rfree_depflags` - REFMAC's first cycle, the deposition's own free set, the common reflections at the lower of the two resolutions |
|
||||
| R_free deposition | the same protocol on the deposited data (`refmac_rfree_depdata`): the R-free to beat |
|
||||
| R_free published | the entry's own R-free (dark and ground-state entries only), from its own refinement program; for information |
|
||||
| difference, quality | ours minus the deposition's; win where ours is lower by more than 0.01, draw within 0.01, loss otherwise |
|
||||
| time, ms/image | wall clock of the rugnux command, which ends with the model validation |
|
||||
| I/O floor | the input's bytes over `--io-gbps` (default 5 GB/s): the time reading alone would take from an SSD |
|
||||
| page cache | share of the input already in memory when the run started (`fincore`): cold, warm or partial |
|
||||
@@ -969,17 +936,15 @@ still_battery.py run --rugnux build/rugnux/rugnux [--only pal90,sacla236 | --tie
|
||||
still_battery.py compare RUN_A RUN_B
|
||||
still_battery.py report RUN [--out DIR]
|
||||
still_battery.py list
|
||||
still_battery.py reference [--only pal90,sacla236]
|
||||
still_battery.py keep RUN --as workstation|server
|
||||
```
|
||||
|
||||
The deposited-data reference is computed once per entry and resolution and kept in the PDB cache
|
||||
(`<pdb_cache>/<pdb>-dimple-<d_cut>/reference.json`, with the dimple run beside it); `run` and
|
||||
`report` compute a missing one on the way (so `report` of an older run gains the column), and
|
||||
`reference` computes or lists them all up front.
|
||||
Both R-free values come from one `model_check.py` call per set (seconds), in `work/<set>/model_check/`.
|
||||
`report RUN --out DIR` of a run scored by an older protocol rescores its merges into `DIR` (rugnux is
|
||||
not run again).
|
||||
|
||||
A run writes `<runs_root>/<YYYYMMDD-HHMM>_<rugnux git>_still-<label>/` with `manifest.json`,
|
||||
`results.json`, `report.md` and `work/<set>/` (the rugnux outputs and the dimple run), made
|
||||
`results.json`, `report.md` and `work/<set>/` (the rugnux outputs and the REFMAC check), made
|
||||
read-only when done; `keep` copies the three files into `results/<role>/<run>/`, like a rotation run.
|
||||
`--tier smoke` runs the sets marked so in the manifest on their first few thousand frames. Data
|
||||
paths come from the `"still"` section of the site config:
|
||||
@@ -990,7 +955,7 @@ paths come from the `"still"` section of the site config:
|
||||
|
||||
`private_manifest` names a local-only manifest of user data (`--private`); such a run gets its own
|
||||
run directory and `report_PRIVATE.md`, and `keep` refuses it. Needs `gemmi` (Python), CCP4's
|
||||
`dimple` and `reindex`, and `fincore` (util-linux) for the page-cache column.
|
||||
`refmac5` and `ctruncate`, and `fincore` (util-linux) for the page-cache column.
|
||||
|
||||
Data preparation for the `_int16` sets is `still/float32_to_int16.py` (see the stills data page).
|
||||
`still/rugnux_stills_ab.py` is the older A/B harness: it compares two argument sets on a fixed
|
||||
|
||||
@@ -496,12 +496,9 @@ def run_one(run_dir, binary, e, opts):
|
||||
|
||||
|
||||
# what the REFMAC model check (model_check.py, --model-check) adds to an open-arm row
|
||||
REFMAC_REFINED = ("refined_rfree_shared", "refined_rfree_shared_depdata", "refined_rwork",
|
||||
"refined_rwork_depdata", "refined_rfree_refmac", "refined_rfree_refmac_depdata",
|
||||
"shared_free_n", "refined_reason")
|
||||
REFMAC_KEYS = ("refmac_rfree", "refmac_rwork", "refmac_rfree_depflags", "refmac_rfree_depdata",
|
||||
"refmac_rfree_ratio", "refmac_twin", "refmac_untwinned", "refmac_status",
|
||||
"refmac_reason", "refmac_refined_ratio") + tuple("refmac_" + k for k in REFMAC_REFINED)
|
||||
"refmac_reason")
|
||||
|
||||
|
||||
def check_model(e, wd):
|
||||
@@ -518,19 +515,15 @@ def check_model(e, wd):
|
||||
except Exception as ex: # a failing check must not lose the set's processing result
|
||||
return dict(out, refmac_status="error", refmac_reason=f"{type(ex).__name__}: {ex}"[:300])
|
||||
out.update({"refmac_" + k: res.get(k) for k in ("rfree", "rwork", "rfree_depflags", "rfree_depdata",
|
||||
"untwinned") + REFMAC_REFINED},
|
||||
"untwinned")},
|
||||
refmac_twin=" ".join(res.get("twin_declared") or []) or None,
|
||||
refmac_status=res.get("status"),
|
||||
refmac_reason=res.get("reason") or res.get("baseline_note"))
|
||||
# The fair ratio: our merge and the depositor's structure factors, both scored on the
|
||||
# The headline: our merge and the depositor's structure factors, both scored on the
|
||||
# depositor's free set by the same protocol. Our own free set was mostly work reflections
|
||||
# for the depositor, so R_free on it reads low.
|
||||
if out["refmac_rfree_depflags"] and out["refmac_rfree_depdata"]:
|
||||
out["refmac_rfree_ratio"] = round(out["refmac_rfree_depflags"] / out["refmac_rfree_depdata"], 4)
|
||||
# the same after refinement of each data set at its own resolution, on the shared free set
|
||||
if out["refmac_refined_rfree_shared"] and out["refmac_refined_rfree_shared_depdata"]:
|
||||
out["refmac_refined_ratio"] = round(out["refmac_refined_rfree_shared"]
|
||||
/ out["refmac_refined_rfree_shared_depdata"], 4)
|
||||
return out
|
||||
|
||||
|
||||
|
||||
+13
-103
@@ -27,13 +27,13 @@ processing (Rugnux is never given the model). For one set it:
|
||||
data go through the same setting choice as Rugnux's (they can be indexed in the other branch
|
||||
of a merohedral ambiguity); where only intensities were deposited, amplitudes come from
|
||||
ctruncate's French-Wilson;
|
||||
* scores both data sets on their common reflections only: the hkl present in both, at
|
||||
d_min_used = max(Rugnux's d_min, the deposited d_min), so a cut either data set made (a
|
||||
low-resolution limit, an anisotropic cut) applies to both;
|
||||
* scores an entry that declares twinning (_pdbx_reflns_twin) with REFMAC's twin refinement,
|
||||
both data sets alike, and keeps the untwinned numbers beside it;
|
||||
* and, as a second protocol (refined_*), refines the same deposited model against each data set
|
||||
at its OWN resolution (restrained, same protocol for both) with R-free on a SHARED free set:
|
||||
the depositor's free reflections present in both data sets (see refined_shared).
|
||||
both data sets alike, and keeps the untwinned numbers beside it.
|
||||
|
||||
The first-cycle numbers are taken at d_min_used = max(Rugnux's d_min, the deposited d_min).
|
||||
The model is never refined: every number is REFMAC's first cycle, before any shift.
|
||||
|
||||
On free sets: the depositor's model was refined against every reflection except THEIR free set,
|
||||
so most of Rugnux's free reflections were work reflections to it. `rfree` (Rugnux's own flags)
|
||||
@@ -64,7 +64,6 @@ PDB_ID = re.compile(r"^[0-9][a-z0-9]{3}$")
|
||||
MAX_ATOMS = 200000 # REFMAC time grows with atoms x reflections
|
||||
CELL_TOL = 0.05 # Niggli cells compared loosely on purpose; the R score decides
|
||||
ANGLE_TOL = 3.0
|
||||
REFINE_CYCLES = 10 # restrained refinement of the shared-free-set protocol
|
||||
SYNTHETIC_FREE = 0.05 # free fraction given to deposited structure factors that carry no free set
|
||||
|
||||
|
||||
@@ -368,86 +367,6 @@ def refmac(workdir, tag, mtz, xyz, d_min, free_label="FreeR_flag", twin=False):
|
||||
return rw[0], rf[0]
|
||||
|
||||
|
||||
def refine(workdir, tag, mtz, xyz, twin=False):
|
||||
"""R-work, R-free after REFMAC's restrained refinement of the model, over the data's whole
|
||||
resolution range: REFINE_CYCLES cycles, automatic weight, isotropic B for every atom (a
|
||||
deposited ANISOU is often TLS-derived, and refining it atom by atom diverged: -LL rising every
|
||||
cycle at 1.4 A on a model deposited at 1.7 A); a ligand without a dictionary entry gets REFMAC's
|
||||
own restraints instead of stopping the run."""
|
||||
rw, rf = run_refmac(workdir, tag, mtz, xyz,
|
||||
f"labin FP=F SIGFP=SIGF FREE=FreeR_flag\nncyc {REFINE_CYCLES}\nweight auto\n"
|
||||
"refi bref iso\nmake newligand continue\nmake check none\n"
|
||||
+ ("twin\n" if twin else "") + "end\n")
|
||||
return rw[-1], rf[-1]
|
||||
|
||||
|
||||
def rename_misfit_ligands(workdir, tag, xyz):
|
||||
"""A model whose ligand carries a CCD code that REFMAC's monomer library describes with other
|
||||
atoms ("atom ... is absent in the library": a code reused for another compound, or a model
|
||||
built with its own) stops a restrained refinement. Such a residue is given a code the library
|
||||
does not have, so REFMAC makes its restraints from the model's own coordinates, as for any new
|
||||
ligand. Returns the renamed model and a note, or None when that is not the failure."""
|
||||
log = open(os.path.join(workdir, tag + ".log"), errors="replace").read()
|
||||
names = sorted(set(re.findall(r"ERROR : atom :\S+\s+(\S+)\s.*is absent in the library", log)))
|
||||
if not names:
|
||||
return None
|
||||
lib = ccp4_env().get("CLIBD_MON", "")
|
||||
free = (c for c in ("".join(t) for t in itertools.product("DQX", "0123456789", "0123456789"))
|
||||
if not os.path.exists(os.path.join(lib, c[0].lower(), c + ".cif")))
|
||||
new = {n: next(free) for n in names}
|
||||
st = gemmi.read_structure(xyz)
|
||||
for ch in st[0]:
|
||||
for res in ch:
|
||||
if res.name in new:
|
||||
res.name = new[res.name]
|
||||
out = os.path.join(workdir, "model_ligands_renamed.pdb")
|
||||
st.write_pdb(out)
|
||||
return out, ", ".join(f"{k}->{v}" for k, v in new.items())
|
||||
|
||||
|
||||
def asu_key(hkl, sg):
|
||||
"""hkl_key of the reflections mapped into gemmi's asymmetric unit (REFMAC writes its output
|
||||
in its own convention, which is not always the same one)."""
|
||||
asu, ops = gemmi.ReciprocalAsu(sg), sg.operations()
|
||||
return hkl_key(np.array([asu.to_asu([int(x) for x in h], ops)[0] for h in hkl]).reshape(-1, 3))
|
||||
|
||||
|
||||
def refined_shared(workdir, data, dep, cell_m, sg_m, xyz, twin):
|
||||
"""Both data sets refined from the deposited model by the same protocol, each over its own
|
||||
resolution range, and R-free on a SHARED free set: the depositor's free reflections present in
|
||||
both data sets. Depositor free reflections outside that intersection are removed from both, so
|
||||
neither refinement uses a free reflection as a work reflection; every other reflection of a
|
||||
data set is a work reflection of that data set."""
|
||||
hkl, f, s = data
|
||||
hkl_d, f_d, s_d, flags_d = dep
|
||||
k, k_d = asu_key(hkl, sg_m), asu_key(hkl_d, sg_m)
|
||||
shared = np.intersect1d(k_d[flags_d == 0], k)
|
||||
keep_d = (flags_d != 0) | np.isin(k_d, shared)
|
||||
free = np.isin(k, shared)
|
||||
out = {"shared_free_n": int(len(shared))}
|
||||
tags = (("", hkl, f, s, free, "data_refine.mtz"),
|
||||
("_depdata", hkl_d[keep_d], f_d[keep_d], s_d[keep_d], flags_d[keep_d] == 0,
|
||||
"deposited_refine.mtz"))
|
||||
for suffix, h, fo, so, fr, name in tags:
|
||||
mtz = os.path.join(workdir, name)
|
||||
make_mtz(h, fo, so, np.where(fr, 0, 1), cell_m, sg_m).write_to_file(mtz)
|
||||
tag = "refine" + (suffix or "_rugnux")
|
||||
try:
|
||||
r = refine(workdir, tag, mtz, xyz, twin)
|
||||
except RuntimeError:
|
||||
renamed = rename_misfit_ligands(workdir, tag, xyz)
|
||||
if renamed is None:
|
||||
raise
|
||||
xyz, out["ligands_renamed"] = renamed
|
||||
r = refine(workdir, tag, mtz, xyz, twin)
|
||||
out["refined_rwork" + suffix], out["refined_rfree_refmac" + suffix] = r
|
||||
# REFMAC's own R-free is over exactly the shared set: it is both data sets' free set. Not
|
||||
# recomputed from the output MTZ - its F carries REFMAC's anisotropic scale, which on
|
||||
# strongly anisotropic data blows up the weak direction's noise and inflates R
|
||||
out["refined_rfree_shared" + suffix] = out["refined_rfree_refmac" + suffix]
|
||||
return {k_: (round(v, 4) if isinstance(v, float) else v) for k_, v in out.items()}
|
||||
|
||||
|
||||
def read_model(xyz):
|
||||
"""The deposited model with its strict-NCS copies generated: an entry can hold part of the
|
||||
asymmetric unit and the operators that make the rest (_struct_ncs_oper 'generate', MTRIX without
|
||||
@@ -611,7 +530,7 @@ def choose_setting(mtz, fc_hkl, fc, cell_m, sg_m, d_lo, d_hi):
|
||||
|
||||
# ------------------------------------------------------------------ the check
|
||||
|
||||
def check(mtz_path, pdb_id, workdir, cache_dir=DEFAULT_CACHE, baseline=True, refined=True):
|
||||
def check(mtz_path, pdb_id, workdir, cache_dir=DEFAULT_CACHE, baseline=True):
|
||||
t0 = time.time()
|
||||
res = {"status": "failed", "reason": "", "rwork": None, "rfree": None,
|
||||
"rfree_deposited": None, "rwork_deposited": None, "d_min_used": None, "sg_used": None,
|
||||
@@ -619,9 +538,6 @@ def check(mtz_path, pdb_id, workdir, cache_dir=DEFAULT_CACHE, baseline=True, ref
|
||||
"rfree_depflags": None,
|
||||
"rwork_depdata": None, "rfree_depdata": None, "baseline_note": None,
|
||||
"depdata_kind": None, "reindex_op_depdata": None, "twin_declared": None,
|
||||
"refined_rwork": None, "refined_rfree_shared": None, "refined_rfree_refmac": None,
|
||||
"refined_rwork_depdata": None, "refined_rfree_shared_depdata": None,
|
||||
"refined_rfree_refmac_depdata": None, "shared_free_n": None, "refined_reason": None,
|
||||
"sg_data": None, "d_min_data": None, "d_min_deposited": None, "cell_diff_pct": None,
|
||||
"setting_r": None, "setting_r_next": None, "fraction_indexed": None,
|
||||
"model_modification": "none", "protocol": "REFMAC R-factors of the model as deposited, bulk solvent, no refinement"}
|
||||
@@ -692,8 +608,12 @@ def check(mtz_path, pdb_id, workdir, cache_dir=DEFAULT_CACHE, baseline=True, ref
|
||||
if dep_scored:
|
||||
res["reindex_op_depdata"] = op_string(dep_scored[0][2])
|
||||
dep = dep_scored[0][3]
|
||||
# both data sets scored on their common reflections only
|
||||
common = np.isin(hkl_key(hkl), hkl_key(dep[0]))
|
||||
hkl, f, s, fl = hkl[common], f[common], s[common], fl[common]
|
||||
dep = tuple(x[np.isin(hkl_key(dep[0]), hkl_key(hkl))] for x in dep)
|
||||
_, i, j = np.intersect1d(hkl_key(hkl), hkl_key(dep[0]), return_indices=True)
|
||||
dep_flags = np.full(len(hkl), 1) # absent from the deposited set -> work
|
||||
dep_flags = np.full(len(hkl), 1)
|
||||
dep_flags[i] = dep[3][j]
|
||||
else:
|
||||
dep, res["baseline_note"] = None, "deposited structure factors too few to place"
|
||||
@@ -723,14 +643,6 @@ def check(mtz_path, pdb_id, workdir, cache_dir=DEFAULT_CACHE, baseline=True, ref
|
||||
if res["twin_declared"]:
|
||||
res["protocol"] += ", REFMAC twin refinement (the entry declares twinning)"
|
||||
res["untwinned"] = score(False, "_untwinned")
|
||||
if dep is not None and refined:
|
||||
try:
|
||||
res.update(refined_shared(workdir, (hkl, f, s), dep, cell_m, sg_m, xyz,
|
||||
bool(res["twin_declared"])))
|
||||
except Exception as e: # the first-cycle numbers stand without it
|
||||
if os.environ.get("MODEL_CHECK_RAISE"):
|
||||
raise
|
||||
res["refined_reason"] = f"{type(e).__name__}: {e}"
|
||||
return done("ok")
|
||||
except Exception as e: # one bad entry must not stop a battery
|
||||
if os.environ.get("MODEL_CHECK_RAISE"):
|
||||
@@ -746,12 +658,10 @@ def main():
|
||||
ap.add_argument("--cache", default=DEFAULT_CACHE)
|
||||
ap.add_argument("--no-baseline", action="store_true",
|
||||
help="skip the depositor's structure factors")
|
||||
ap.add_argument("--no-refine", action="store_true",
|
||||
help="first-cycle R factors only (skip the refined shared-free-set protocol)")
|
||||
a = ap.parse_args()
|
||||
wd = a.workdir or os.path.join(os.path.dirname(os.path.abspath(a.mtz)), "model_check_" + a.pdb_id)
|
||||
print(json.dumps(check(a.mtz, a.pdb_id, os.path.abspath(wd), a.cache, baseline=not a.no_baseline,
|
||||
refined=not a.no_refine), indent=1))
|
||||
print(json.dumps(check(a.mtz, a.pdb_id, os.path.abspath(wd), a.cache, baseline=not a.no_baseline),
|
||||
indent=1))
|
||||
|
||||
|
||||
if __name__ == "__main__":
|
||||
|
||||
@@ -19,9 +19,8 @@ with embedded fonts and the same figure as a 600 d.p.i. PNG. gnuplot and pdftopp
|
||||
on PATH; gemmi reads the facility and beamline from each entry's cached mmCIF (_diffrn_source).
|
||||
|
||||
Left figure: R-free of the deposited model after the battery's REFMAC check on the depositor's data
|
||||
(x) and on Rugnux's data (y). By default the refined protocol (model_check.py: the same model
|
||||
refined against each data set at its own resolution, R-free over the free reflections both data sets
|
||||
share); --rfree first-cycle draws the first-cycle rigid-body numbers at the common resolution instead.
|
||||
(x) and on Rugnux's data (y): the model as deposited, not refined, REFMAC's first cycle on the
|
||||
depositor's free set and both data sets' common reflections (model_check.py).
|
||||
A run without REFMAC values (the check is opt-in, --model-check) gets no R-free figure; the rest are drawn.
|
||||
Right figure: the deposited high-resolution limit (x) against Rugnux's own cut
|
||||
(y), log axes so the low-resolution entries show. Small molecules (no deposited model) and sets
|
||||
@@ -201,8 +200,6 @@ def main():
|
||||
ap = argparse.ArgumentParser(description=__doc__.split("\n\n")[0])
|
||||
ap.add_argument("runs", nargs="+", help="battery run directories or names (later ones fill in)")
|
||||
ap.add_argument("--out", required=True, help="output directory")
|
||||
ap.add_argument("--rfree", choices=["refined", "first-cycle"], default="refined",
|
||||
help="which REFMAC R-free pair to draw (default: refined, on the shared free set)")
|
||||
ap.add_argument("--time-runs", nargs="+", help="runs to draw the time histogram from "
|
||||
"(default: the runs above)")
|
||||
ap.add_argument("--time-arms", nargs="+", default=["open"],
|
||||
@@ -234,14 +231,11 @@ def main():
|
||||
excluded.append((s, "both", "verdict fail: " + str(r.get("reason"))))
|
||||
continue
|
||||
where[s] = source(r["model"])
|
||||
if a.rfree == "refined":
|
||||
dep, rug = r.get("refmac_refined_rfree_shared_depdata"), r.get("refmac_refined_rfree_shared")
|
||||
else:
|
||||
dep, rug = r.get("refmac_rfree_depdata"), r.get("refmac_rfree_depflags")
|
||||
dep, rug = r.get("refmac_rfree_depdata"), r.get("refmac_rfree_depflags")
|
||||
if dep and rug:
|
||||
rfree.append((dep, rug, s))
|
||||
else:
|
||||
excluded.append((s, "rfree", "no REFMAC pair: " + str(r.get("refmac_refined_reason") or r.get("refmac_reason")
|
||||
excluded.append((s, "rfree", "no REFMAC pair: " + str(r.get("refmac_reason")
|
||||
or r.get("dep_reason") or "depositor data not scored")))
|
||||
if r.get("d_min_ref") and r.get("d_min"):
|
||||
dmin.append((r["d_min_ref"], r["d_min"], s))
|
||||
|
||||
+50
-174
@@ -1,21 +1,22 @@
|
||||
#!/usr/bin/env python3
|
||||
"""Rugnux still battery: run serial (XFEL stills) datasets, refine the deposited model against each
|
||||
"""Rugnux still battery: run serial (XFEL stills) datasets, score the deposited model against each
|
||||
merge, report, compare.
|
||||
|
||||
The rotation battery (battery.py) scores a lattice and a space group. For stills both are given
|
||||
(-C/-S, the way a serial user runs), so this battery reports what the merge is worth and what it
|
||||
costs instead: indexing rate, R_split, the R-free of the deposited model refined against our merge
|
||||
(dimple, with the deposition's own free set) beside the R-free of the same protocol on the deposited
|
||||
data (the reference the verdict is taken against; the published R-free is shown for information), and
|
||||
the wall time of the rugnux command with the input's I/O floor next to it.
|
||||
costs instead: indexing rate, R_split, the R-free of the deposited model against our merge, not
|
||||
refined (model_check.py: REFMAC's first cycle, the deposition's own free set) beside the R-free of
|
||||
the same protocol on the deposited data (the target the verdict is taken against; the published
|
||||
R-free is shown for information), and the wall time of the rugnux command with the input's I/O floor
|
||||
next to it.
|
||||
|
||||
Usage (README.md in this directory, section "The still battery", is the guide):
|
||||
still_battery.py run --rugnux BIN [--label L] [--only a,b | --tier smoke] [--threads N]
|
||||
[--cold] [--no-refine] [--io-gbps G]
|
||||
[--cold] [--io-gbps G]
|
||||
still_battery.py compare RUN_A RUN_B
|
||||
still_battery.py report RUN [--out DIR] re-render a run's report.md
|
||||
still_battery.py report RUN [--out DIR] re-render a run's report.md (a run scored by an
|
||||
older protocol is rescored: give --out)
|
||||
still_battery.py list still runs under the runs root
|
||||
still_battery.py reference [--only a,b] compute (or show) the deposited-data dimple references
|
||||
still_battery.py keep RUN --as workstation|server copy report.md, results.json, manifest.json
|
||||
into results/<role>/<run>/ (committed); refuses a
|
||||
private run
|
||||
@@ -42,10 +43,7 @@ HERE = os.path.dirname(os.path.abspath(__file__))
|
||||
sys.path.insert(0, HERE)
|
||||
|
||||
import battery # noqa: E402 binary_info, hardware_info, run_logged, write_json, make_read_only, ...
|
||||
import model_check # noqa: E402 fetch, ccp4_env, RCSB
|
||||
|
||||
import gemmi # noqa: E402
|
||||
import numpy as np # noqa: E402
|
||||
import model_check # noqa: E402 check, fetch, RCSB
|
||||
|
||||
MANIFEST = os.path.join(HERE, "still_manifest.json")
|
||||
KEPT_FILES = ("report.md", "results.json", "manifest.json")
|
||||
@@ -55,9 +53,8 @@ DRAW = 0.01 # R-free difference to the target still called a draw
|
||||
COLUMNS = [("facility", "facility"), ("beamline", "beamline"), ("detector", "detector"),
|
||||
("dataset", "set"), ("PDB", "pdb"), ("images", "images"), ("images (deposition)", "images_deposited"),
|
||||
("indexed %", "indexed_pct"), ("d_min", "d_min"), ("R_split", "r_split"),
|
||||
("R_free rugnux (dimple)", "r_free"), ("R_free deposition (dimple)", "r_free_dimple_deposited"),
|
||||
("R_free published", "r_free_deposited"), ("target", "r_free_target_from"),
|
||||
("difference", "r_free_diff"), ("quality", "quality"),
|
||||
("R_free rugnux", "r_free"), ("R_free deposition", "r_free_depdata"),
|
||||
("R_free published", "r_free_deposited"), ("difference", "r_free_diff"), ("quality", "quality"),
|
||||
("time s", "wall_s"), ("ms/image", "ms_per_image"), ("I/O floor s", "io_floor_s"),
|
||||
("page cache", "cache")]
|
||||
|
||||
@@ -166,136 +163,25 @@ def num(s):
|
||||
return None
|
||||
|
||||
|
||||
def free_flags(sf_cif, out_mtz, free):
|
||||
"""An MTZ of the deposition's free set (FreeR_flag = 1 free, 0 working) from its -sf.cif. The free
|
||||
set is `status` 'f' where the entry has it, else the pdbx_r_free_flag value held by 1-15% of the
|
||||
reflections; the manifest row may name it instead ({"column": ..., "value": ...})."""
|
||||
rb = gemmi.as_refln_blocks(gemmi.cif.read(sf_cif))[0]
|
||||
h, k, l = (rb.make_int_array("index_" + c, 0) for c in "hkl")
|
||||
col = val = None
|
||||
if free:
|
||||
col, val = free["column"], str(free["value"])
|
||||
else:
|
||||
for c in ("status", "pdbx_r_free_flag"):
|
||||
tag = rb.block.find_values("_refln." + c)
|
||||
if not len(tag):
|
||||
continue
|
||||
vals = [gemmi.cif.as_string(v) for v in tag]
|
||||
counts = {v: vals.count(v) for v in set(vals)}
|
||||
cands = [v for v, n in counts.items() if 0.01 <= n / len(vals) <= 0.15]
|
||||
if c == "status" and "f" in cands:
|
||||
col, val = c, "f"
|
||||
break
|
||||
if c == "pdbx_r_free_flag" and cands:
|
||||
col, val = c, min(cands, key=lambda v: counts[v])
|
||||
break
|
||||
if col is None:
|
||||
return None
|
||||
flags = [1.0 if gemmi.cif.as_string(v) == val else 0.0 for v in rb.block.find_values("_refln." + col)]
|
||||
mtz = gemmi.Mtz(with_base=True)
|
||||
mtz.spacegroup = rb.spacegroup
|
||||
mtz.cell = rb.cell
|
||||
mtz.add_dataset("free")
|
||||
mtz.add_column("FreeR_flag", "I")
|
||||
mtz.set_data(np.column_stack([h, k, l, flags]).astype(np.float32))
|
||||
mtz.write_to_file(out_mtz)
|
||||
return f"{col}={val}"
|
||||
|
||||
|
||||
def deposited_model(pdb, cache):
|
||||
"""The deposited coordinates: PDB format where the RCSB serves it, else mmCIF (entries too large or
|
||||
too new for the PDB format exist only as mmCIF; rugnux and dimple read both)."""
|
||||
too new for the PDB format exist only as mmCIF; rugnux reads both)."""
|
||||
return (model_check.fetch(model_check.RCSB + pdb + ".pdb", os.path.join(cache, pdb + ".pdb"))
|
||||
or model_check.fetch(model_check.RCSB + pdb + ".cif", os.path.join(cache, pdb + ".cif")))
|
||||
|
||||
|
||||
def refine(wd, e, cache):
|
||||
"""dimple of the deposited model against our merge, reindexed into the deposition's group, with the
|
||||
deposition's free set, cut at the deposition's resolution. (R_work, R_free, note)."""
|
||||
pdb = e["pdb"].lower()
|
||||
xyz = deposited_model(pdb, cache)
|
||||
sf = model_check.fetch(model_check.RCSB + pdb + "-sf.cif", os.path.join(cache, pdb + "-sf.cif"))
|
||||
if not (xyz and sf) or not os.path.exists(os.path.join(wd, "p.mtz")):
|
||||
return None, None, "no model, structure factors or merge"
|
||||
free_mtz = os.path.join(wd, "deposited_free.mtz")
|
||||
free = free_flags(sf, free_mtz, e.get("free"))
|
||||
if free is None:
|
||||
return None, None, "no free set in the deposition"
|
||||
env = model_check.ccp4_env()
|
||||
if xyz.endswith(".cif"): # dimple takes PDB-format coordinates only
|
||||
subprocess.run(["gemmi", "convert", xyz, "deposited.pdb"], cwd=wd, env=env, capture_output=True)
|
||||
xyz = os.path.join(wd, "deposited.pdb")
|
||||
subprocess.run(["reindex", "hklin", "p.mtz", "hklout", "p_sg.mtz"], cwd=wd, env=env,
|
||||
input=f"SYMM '{e['sg']}'\nEND\n", capture_output=True, text=True)
|
||||
with open(os.path.join(wd, "dimple.out"), "w") as log:
|
||||
subprocess.run(["dimple", "--free-r-flags", free_mtz, "--freecolumn", "FreeR_flag=1", "-I", "IMEAN",
|
||||
"--reso", str(e["d_cut"]), "p_sg.mtz", xyz, "dimple"],
|
||||
cwd=wd, env=env, stdout=log, stderr=subprocess.STDOUT)
|
||||
log = os.path.join(wd, "dimple", "dimple.log")
|
||||
m = re.findall(r"R/Rfree\s+\S+\s+->\s+([\d.]+)/([\d.]+)", open(log).read()) if os.path.exists(log) else []
|
||||
if not m:
|
||||
return None, None, "dimple failed"
|
||||
return float(m[-1][0]), float(m[-1][1]), f"free set {free}"
|
||||
|
||||
|
||||
def deposited_mtz(sf_cif, out_mtz):
|
||||
"""IMEAN/SIGIMEAN of the deposition from its -sf.cif: intensity_meas where deposited, else the mean
|
||||
of the anomalous pair pdbx_I_plus/minus, else F_meas_au squared."""
|
||||
rb = gemmi.as_refln_blocks(gemmi.cif.read(sf_cif))[0]
|
||||
cols = rb.column_labels()
|
||||
if "intensity_meas" in cols:
|
||||
i, s = rb.make_float_array("intensity_meas"), rb.make_float_array("intensity_sigma")
|
||||
elif "pdbx_I_plus" in cols:
|
||||
ip, sp, im, sm = (rb.make_float_array(c) for c in ("pdbx_I_plus", "pdbx_I_plus_sigma",
|
||||
"pdbx_I_minus", "pdbx_I_minus_sigma"))
|
||||
n = np.isfinite(ip).astype(int) + np.isfinite(im)
|
||||
with np.errstate(invalid="ignore"):
|
||||
i = (np.nan_to_num(ip) + np.nan_to_num(im)) / n
|
||||
s = np.sqrt(np.nan_to_num(sp) ** 2 + np.nan_to_num(sm) ** 2) / n
|
||||
else:
|
||||
f, sf = rb.make_float_array("F_meas_au"), rb.make_float_array("F_meas_sigma_au")
|
||||
i, s = f ** 2, 2 * f * sf
|
||||
ok = np.isfinite(i)
|
||||
mtz = gemmi.Mtz(with_base=True)
|
||||
mtz.spacegroup = rb.spacegroup
|
||||
mtz.cell = rb.cell
|
||||
mtz.add_dataset("deposited")
|
||||
mtz.add_column("IMEAN", "J")
|
||||
mtz.add_column("SIGIMEAN", "Q")
|
||||
mtz.set_data(np.column_stack([rb.make_miller_array()[ok], i[ok], s[ok]]).astype(np.float32))
|
||||
mtz.write_to_file(out_mtz)
|
||||
|
||||
|
||||
def deposited_reference(e, cache):
|
||||
"""R-free of the deposited data through the same dimple protocol as our merge (refine()): what the
|
||||
deposition reaches under 8 REFMAC cycles, so the like-for-like reference rather than the published
|
||||
R-free. Computed once per entry and resolution and kept in the PDB cache. {r_work, r_free, note}."""
|
||||
pdb = e["pdb"].lower()
|
||||
wd = os.path.join(cache, f"{pdb}-dimple-{e['d_cut']}")
|
||||
done = os.path.join(wd, "reference.json")
|
||||
if os.path.exists(done):
|
||||
return json.load(open(done))
|
||||
sf = model_check.fetch(model_check.RCSB + pdb + "-sf.cif", os.path.join(cache, pdb + "-sf.cif"))
|
||||
shutil.rmtree(wd, ignore_errors=True)
|
||||
os.makedirs(wd)
|
||||
deposited_mtz(sf, os.path.join(wd, "p.mtz"))
|
||||
r_work, r_free, note = refine(wd, e, cache)
|
||||
ref = {"r_work": r_work, "r_free": r_free, "note": note}
|
||||
if r_free is not None:
|
||||
battery.write_json(done, ref)
|
||||
return ref
|
||||
|
||||
|
||||
def target(r):
|
||||
"""The R-free to beat: the lower of the deposition's under our dimple protocol and its published
|
||||
one. (value, "dimple" | "published"), or (None, None) without either."""
|
||||
cands = [(r.get("r_free_dimple_deposited"), "dimple"), (r.get("r_free_deposited"), "published")]
|
||||
cands = [c for c in cands if c[0] is not None]
|
||||
return min(cands) if cands else (None, None)
|
||||
def score(wd, mtz, e, cache):
|
||||
"""R-free of the deposited model, not refined, against our merge and against the deposited data
|
||||
through the same protocol (model_check.py: REFMAC's first cycle, the deposition's free set, both
|
||||
data sets' common reflections at the lower of the two resolutions)."""
|
||||
res = model_check.check(mtz, e["pdb"], os.path.join(wd, "model_check"), cache)
|
||||
return {"r_work": res["rwork"], "r_free": res["rfree_depflags"], "r_free_depdata": res["rfree_depdata"],
|
||||
"d_min_used": res["d_min_used"],
|
||||
"score_note": res["reason"] or res["baseline_note"] or res["model_modification"]}
|
||||
|
||||
|
||||
def quality(r_free, goal):
|
||||
"""win / draw / loss of our dimple R-free against the target."""
|
||||
"""win / draw / loss of our R-free against the deposited data's."""
|
||||
if r_free is None or goal is None:
|
||||
return None
|
||||
d = r_free - goal
|
||||
@@ -337,9 +223,8 @@ def run_one(run_dir, binary, e, st, opts):
|
||||
cc_half_split=num(keys.get("CC_HALF_SPLIT")), r_free_model_check=num(keys.get("R_FREE")))
|
||||
if rc != 0 or not keys:
|
||||
return dict(row, verdict="crash" if rc else "fail", reason=battery.failure_note(os.path.join(wd, "run.log")))
|
||||
if opts["refine"] and e.get("pdb"):
|
||||
row["r_work"], row["r_free"], row["refine_note"] = refine(wd, e, st["pdb_cache"])
|
||||
row["r_free_dimple_deposited"] = deposited_reference(e, st["pdb_cache"])["r_free"]
|
||||
if e.get("pdb"):
|
||||
row.update(score(wd, os.path.join(wd, "p.mtz"), e, st["pdb_cache"]))
|
||||
return dict(row, verdict="done", reason="")
|
||||
|
||||
|
||||
@@ -350,7 +235,7 @@ def fmt(v, key):
|
||||
return "-"
|
||||
if key == "r_split":
|
||||
return f"{100 * v:.1f}%"
|
||||
if key in ("r_free", "r_free_deposited", "r_free_dimple_deposited"):
|
||||
if key in ("r_free", "r_free_deposited", "r_free_depdata"):
|
||||
return f"{v:.3f}"
|
||||
if key == "r_free_diff":
|
||||
return f"{v:+.3f}"
|
||||
@@ -361,19 +246,23 @@ def fmt(v, key):
|
||||
return str(v)
|
||||
|
||||
|
||||
def render(run_dir, cache):
|
||||
"""report.md of a run. A run from before the deposited-data reference existed gets it from the
|
||||
PDB cache (computed there if missing)."""
|
||||
def render(run_dir, cache, src=None):
|
||||
"""report.md of a run. A row scored by an older protocol (no r_free_depdata) is rescored from the
|
||||
merge in src's work/<set>/ (src: the run itself unless the report is written elsewhere), into
|
||||
run_dir's work/<set>/model_check/, and results.json is rewritten."""
|
||||
man = json.load(open(os.path.join(run_dir, "manifest.json")))
|
||||
res = json.load(open(os.path.join(run_dir, "results.json")))
|
||||
sets = {e["id"]: e for e in man["sets"]}
|
||||
stale = [r for r in res if r.get("verdict") == "done" and r.get("pdb") and "r_free_depdata" not in r]
|
||||
for r in stale:
|
||||
r.update(score(os.path.join(run_dir, "work", r["set"]),
|
||||
os.path.join(src or run_dir, "work", r["set"], "p.mtz"), sets[r["set"]], cache))
|
||||
if stale:
|
||||
battery.write_json(os.path.join(run_dir, "results.json"), res)
|
||||
for r in res:
|
||||
if r.get("r_free") is not None and r.get("r_free_dimple_deposited") is None:
|
||||
r["r_free_dimple_deposited"] = deposited_reference(sets[r["set"]], cache)["r_free"]
|
||||
goal, r["r_free_target_from"] = target(r)
|
||||
if r.get("r_free") is not None and goal is not None:
|
||||
r["r_free_diff"] = r["r_free"] - goal
|
||||
r["quality"] = quality(r.get("r_free"), goal)
|
||||
if r.get("r_free") is not None and r.get("r_free_depdata") is not None:
|
||||
r["r_free_diff"] = r["r_free"] - r["r_free_depdata"]
|
||||
r["quality"] = quality(r.get("r_free"), r.get("r_free_depdata"))
|
||||
b = man["binary"]
|
||||
hw = man.get("hardware") or {}
|
||||
lines = [f"# Rugnux still battery - {os.path.basename(os.path.normpath(run_dir))}", "",
|
||||
@@ -383,10 +272,11 @@ def render(run_dir, cache):
|
||||
f"{', '.join(g['name'] for g in hw.get('gpus', [])) or 'none'}",
|
||||
f"- start {man['start']}, end {man.get('end') or 'unfinished'}; subset: {man.get('subset') or 'all sets'}"
|
||||
+ ("; ABORTED" if man.get("aborted") else ""),
|
||||
f"- R_free rugnux (dimple): dimple of the deposited model against our merge, with the deposition's "
|
||||
f"free set, at its resolution; R_free deposition (dimple): the same protocol on the deposited data; "
|
||||
f"R_free published: the entry's own value (phenix etc.); target: the lower of the two, the one "
|
||||
f"named; difference: ours minus the target; quality: win where ours is lower by more than {DRAW}, "
|
||||
f"- R_free rugnux: the deposited model against our merge, not refined (REFMAC's first cycle, bulk "
|
||||
f"solvent and anisotropic scale, the deposition's free set, the common reflections at the lower of "
|
||||
f"the two resolutions; model_check.py); R_free deposition: the same protocol on the deposited data, "
|
||||
f"the target; R_free published: the entry's own value, for information; difference: ours minus "
|
||||
f"the deposition's; quality: win where ours is lower by more than {DRAW}, "
|
||||
f"draw within {DRAW}, loss otherwise. Time: wall clock of the rugnux "
|
||||
f"command (to model validation); I/O floor: input bytes / {man['options']['io_gbps']} GB/s; page "
|
||||
f"cache: share of the input already in memory when the run started.", "",
|
||||
@@ -417,7 +307,7 @@ def cmd_run(a, st):
|
||||
os.makedirs(os.path.dirname(binary))
|
||||
shutil.copy2(a.rugnux, binary)
|
||||
info.update(copied_from=os.path.realpath(a.rugnux), sha256=battery.sha256(binary))
|
||||
opts = {"threads": a.threads, "timeout": a.timeout, "tier": a.tier, "refine": not a.no_refine,
|
||||
opts = {"threads": a.threads, "timeout": a.timeout, "tier": a.tier,
|
||||
"io_gbps": a.io_gbps, "cold": a.cold}
|
||||
man = {"kind": "still", "label": label, "private": a.private, "subset": a.only or (f"tier {a.tier}" if a.tier else None),
|
||||
"binary": info, "runner_git": battery.git(HERE, "rev-parse", "--short", "HEAD") or "unknown",
|
||||
@@ -472,22 +362,13 @@ def cmd_list(a, st):
|
||||
|
||||
def cmd_report(a, st):
|
||||
src = run_path(st, a.run)
|
||||
dest = src
|
||||
if a.out:
|
||||
tmp = os.path.join(a.out, os.path.basename(os.path.normpath(src)))
|
||||
os.makedirs(tmp, exist_ok=True)
|
||||
dest = os.path.join(a.out, os.path.basename(os.path.normpath(src)))
|
||||
os.makedirs(dest, exist_ok=True)
|
||||
for f in ("manifest.json", "results.json"):
|
||||
shutil.copyfile(os.path.join(src, f), os.path.join(tmp, f))
|
||||
src = tmp
|
||||
print(render(src, st["pdb_cache"]))
|
||||
|
||||
|
||||
def cmd_reference(a, st):
|
||||
print("| set | PDB | d_cut | R_work | R_free deposition (dimple) | R_free published |")
|
||||
print("|---|---|---|---|---|---|")
|
||||
for e in select(st, a.only, None, a.private):
|
||||
ref = deposited_reference(e, st["pdb_cache"])
|
||||
print(f"| {e['id']} | {e['pdb']} | {e['d_cut']} | {ref['r_work']} | {ref['r_free']} | "
|
||||
f"{e.get('r_free_deposited')} |", flush=True)
|
||||
shutil.copyfile(os.path.join(src, f), os.path.join(dest, f))
|
||||
print(render(dest, st["pdb_cache"], src))
|
||||
|
||||
|
||||
def cmd_keep(a, st):
|
||||
@@ -516,7 +397,6 @@ def main():
|
||||
r.add_argument("--cold", action="store_true", help="evict each set's input from the page cache before its run")
|
||||
r.add_argument("--timeout", type=int, default=6 * 3600)
|
||||
r.add_argument("--build-dir")
|
||||
r.add_argument("--no-refine", action="store_true", help="skip the dimple refinement")
|
||||
r.add_argument("--io-gbps", type=float, default=IO_GBPS)
|
||||
c = sub.add_parser("compare")
|
||||
c.add_argument("run_a")
|
||||
@@ -525,16 +405,12 @@ def main():
|
||||
p.add_argument("run")
|
||||
p.add_argument("--out")
|
||||
sub.add_parser("list")
|
||||
f = sub.add_parser("reference")
|
||||
f.add_argument("--only")
|
||||
f.add_argument("--private", action="store_true")
|
||||
k = sub.add_parser("keep")
|
||||
k.add_argument("run")
|
||||
k.add_argument("--as", dest="role", required=True, choices=["workstation", "server"])
|
||||
a = ap.parse_args()
|
||||
st = load_site(a.site)
|
||||
{"run": cmd_run, "compare": cmd_compare, "list": cmd_list, "report": cmd_report, "reference": cmd_reference,
|
||||
"keep": cmd_keep}[a.cmd](a, st)
|
||||
{"run": cmd_run, "compare": cmd_compare, "list": cmd_list, "report": cmd_report, "keep": cmd_keep}[a.cmd](a, st)
|
||||
|
||||
|
||||
if __name__ == "__main__":
|
||||
|
||||
@@ -1,4 +1,4 @@
|
||||
{"_comment": "The still battery (still_battery.py). inputs: globs relative to the site's still data_root, each glob's matches sorted by name; {data} in args is the data root. Every set runs with its deposition's cell and space group (-C/-S); r_free_deposited is the entry's published R-free (dark/ground-state entries only); d_cut is the deposition's resolution, used for the dimple refinement. Where to get the data: docs/EXTERNAL_TEST_DATA_STILLS.md.",
|
||||
{"_comment": "The still battery (still_battery.py). inputs: globs relative to the site's still data_root, each glob's matches sorted by name; {data} in args is the data root. Every set runs with its deposition's cell and space group (-C/-S); r_free_deposited is the entry's published R-free (dark/ground-state entries only); d_cut is the deposition's resolution. Where to get the data: docs/EXTERNAL_TEST_DATA_STILLS.md.",
|
||||
"sets": [
|
||||
{"id": "sacla243", "cxidb": 243, "facility": "SACLA", "beamline": "BL2", "detector": "MPCCD (8 panels, 50 um)", "format": "Cheetah-SACLA HDF5", "pdb": "31JT", "sg": "P 2 21 21", "d_cut": 1.63, "r_free_deposited": 0.1943, "images_deposited": 52151, "inputs": ["sacla_243/raw_all/run*.h5"], "args": ["-C", "39.36,79.16,145.40,90,90,90", "-S", "P 2 21 21"], "tiers": {}},
|
||||
{"id": "sacla237", "cxidb": 237, "facility": "SACLA", "beamline": "BL2", "detector": "MPCCD (8 panels, 50 um)", "format": "Cheetah-SACLA HDF5", "pdb": "9S06", "sg": "P 21 21 21", "d_cut": 2.05, "r_free_deposited": 0.219, "images_deposited": null, "inputs": ["sacla_237/raw/run13063[7-9]-*.h5", "sacla_237/raw/run13064[0-2]-*.h5"], "args": ["-C", "64.80,71.09,207.26,90,90,90", "-S", "P212121"], "tiers": {}},
|
||||
|
||||
Reference in new issue
Block a user