Skip to content

Fetch, cache and display electron density and cryo-EM maps - #1134

Open
aalhossary wants to merge 8 commits into
biojava:masterfrom
aalhossary:aa/electron-density-maps
Open

Fetch, cache and display electron density and cryo-EM maps#1134
aalhossary wants to merge 8 commits into
biojava:masterfrom
aalhossary:aa/electron-density-maps

Conversation

@aalhossary

Copy link
Copy Markdown
Member

Closes #947.

This PR is about a new feature, therefore, I recommend releasing it as BioJava 7.3.0 rather than 6.2.7.
It adds org.biojava.nbio.structure.io.density, which downloads and caches density maps the way LocalPDBDirectory caches coordinate files, and a JmolPanel method that contours the result.

DensityMapCache cache = new DensityMapCache();
DensityMapResult map = cache.getDensityMap(new PdbId("1cbs"), DensityMapKind.TWO_FO_FC);
jmolPanel.loadDensityMap(map);

Why a chain of sources rather than one

edmaps.rcsb.org shut down in October 2024. What RCSB documents in its place —the map coefficients published with the wwPDB validation reports — are structure factors, not a sampled grid, and cannot be displayed without a Fourier transform.
Several other services do serve grids, and they differ enormously in size for the same entry, so rather than picking one, sources are tried in order until one answers:

order
X-ray RCSB density server → PDBe CCP4 → PDBe density server → wwPDB coefficients (off by default)
cryo-EM RCSB density server → PDBe density server → EMDB primary map

Measured sizes for the same data:

entry full archive density server
1cbs (X-ray) 2.1 MB (two CCP4 files) 210 kB at detail 0
EMD-0262 (6hu9) 116 MB 0.5–3.7 MB

wwPDB coefficients are supported for completeness but disabled by default, and DensityFileFormat.isJmolLoadable() lets a viewer refuse them rather than silently drawing nothing.

Design notes

  • A source that has nothing for an entry is skipped; any other transport failure aborts the chain. A network outage must never be reported as "this entry has no density". When every source is exhausted, NoDensityMapException carries the reason from each one, so the UI can explain rather than just fail.
  • A density server response contains both the 2Fo-Fc and Fo-Fc blocks, so the two kinds share one cache entry instead of downloading the identical file twice.
  • Cryo-EM entries are resolved through EMDB's search API, with RCSB as a fallback. That lookup also yields the depositors' recommended contour level, which is how an EM map should be contoured; it is attached to the result whichever source supplied the voxels. The experimental method is never inferred from resolution, which BioJava parses incorrectly for some cryo-EM entries (Resolution lost in cryo-EM #1000).
  • Ccp4Header checks for the MAP stamp at byte 208, so a server answering with an error page and HTTP 200 produces a clean cache miss rather than a corrupt entry.
  • A .meta sidecar fully describes each cached result, so LOCAL_ONLY is served without opening a connection.

Three Jmol behaviours established by experiment

Each of these changed the implementation, and none is documented:

  1. Option order in isosurface is load-bearing. With mesh nofill before the file name, Jmol accepts the command, reports no error, and draws nothing at all.
  2. A negative sigma does not contour at a negative level. Jmol reserves negative sigma internally, so sigma -3.0 silently contours at the default level — the intended red lobe of a difference map came out identical to the blue one. Difference maps are drawn as a single signed surface instead.
  3. Selecting the Fo-Fc block of a cached BinaryCIF file needs the marker in the file name. BCifDensityReader picks the block by testing whether the file name contains &diff=1, which normally arrives in the URL query string. A cached local file has no query string; appending the marker to the file URL fails, and the #diff=1 form the reader's own to-do comment suggests is ignored. Embedding it in the name works, so the cache exposes the difference map under a companion name, hard-linked to the same bytes.

Measured proof of (3): the plain file contours at cutoff 0.356 over a −1.31 to 3.78 range; the companion name contours at 0.374 over −0.69 to 0.85. Different data, so the marker genuinely selects the other block.

Testing

43 unit tests, entirely offline — the chain tests use stub providers, so "404 falls through" and "anything else aborts" are pinned exactly. 6 integration tests against the live services costing a couple of megabytes in total; the cryo-EM path is covered without downloading the 116 MB map by setting the size limit so the guard fires after resolution but before transfer.

Jmol version

No change: BioJava stays on net.sourceforge.jmol:jmol:14.31.10, which is the newest on Maven Central (October 2020). Every behaviour relied on here — MrcBinaryReader, BCifDensityReader, the &diff=1 selection, gzip sniffing — is identical in 14.31.10 and current Jmol, so downstream projects overriding the dependency with a newer build are unaffected.

@aalhossary

Copy link
Copy Markdown
Member Author

Worth discussing: how small is small enough?

This PR deliberately prefers the smallest adequate representation, and I would like a second opinion on where the line should sit.

The density servers return a downsampled grid, and the saving is not marginal. For 1cbs, detail 0 is 210 kB against 2.1 MB for the two full CCP4 files — and one request instead of two, since a single response carries both map kinds. For the cryo-EM map behind 6hu9 it is 0.5–3.7 MB against 116 MB. The current defaults are detail 3 for RCSB and detail 6 for PDBe, matching what each service's own viewer asks for.

Open questions:

  1. Is a downsampled grid adequate for what people actually use BioJava density for?
    For looking at a ligand pocket in a viewer I am confident it is. For anything quantitative — real-space correlation, occupancy refinement, automated density fit scoring — it may well not be, and someone doing that would want the full grid and would not necessarily notice they had been handed a coarser one. The result object does report its source, so the information is available, but "available if you check" is weaker than a sensible default.

  2. Should the default detail level be lower?
    Detail 0 for 1cbs is a tenth the size of detail 3 and still perfectly readable around a ligand. I kept the services' own defaults rather than choosing for them, but I am not attached to that.

  3. Should the size ceiling be the primary mechanism instead of the ordering?
    Right now the 256 MiB ceiling rarely fires: EMD-0262 is 116 MB and would sail under it. What actually avoids the 116 MB download is putting the density servers first. An alternative design would always prefer the full map and fall back only when it is genuinely huge, which is more predictable but much more expensive by default.

  4. Should a quantitative caller get a different default than a viewer?
    There is already allowNonRenderableFormats distinguishing the two audiences; a preferFullResolution flag would be a small addition if the distinction is real.

Happy to change any of this — the ordering is one list and a couple of setters.

@aalhossary

Copy link
Copy Markdown
Member Author

Reviewer note: This PR (when approved) should be merged after #1133; because the other PR is included, a prerequisite, and a dependency of this PR.

@josemduarte josemduarte left a comment

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thank you! I skimmed through it and it looks good. See only one comment about the PDB archives layout.

@josemduarte

Copy link
Copy Markdown
Contributor

And one thing to consider: since the tests are new, write them in junit5 (jupiter). Probably an easy task for AI.

Adds org.biojava.nbio.structure.io.density, which downloads and caches density
maps the same way LocalPDBDirectory caches coordinate files, and hands back a
File that a viewer can contour. Closes biojava#947.

Several services publish density for the PDB and they differ enormously in size
for the same entry, so rather than picking one, sources are tried in order until
one answers. The order is smallest-adequate-first, because the smallest form is
usually perfectly good to look at:

  X-ray:  RCSB density server -> PDBe CCP4 -> PDBe density server
          -> wwPDB map coefficients (disabled by default)
  cryo-EM: RCSB density server -> PDBe density server -> EMDB primary map

For 1cbs a density server slice is about a tenth the size of the equivalent pair
of CCP4 files. For the map behind 6hu9 it is 3.7 MB against a 106 MB primary
map. A size limit, 256 MiB by default, is checked against the size EMDB itself
reports before any of the body is transferred, and exceeding it is not an error:
the chain simply falls back to a smaller representation.

wwPDB map coefficients are supported for completeness, since they are the route
RCSB documents now that edmaps.rcsb.org has shut down, but they are structure
factors rather than a sampled grid and cannot be displayed without a Fourier
transform. They are therefore disabled by default, and DensityFileFormat carries
an isJmolLoadable() flag so that a viewer can refuse them rather than silently
drawing nothing.

Notes on the design:

* A source that has nothing for an entry is skipped and the next is tried, but
  any other transport failure aborts the chain. A network outage must never be
  reported as "this entry has no density". When every source is exhausted,
  NoDensityMapException carries the reason from each one, so a caller can say
  why rather than just that it failed.
* A density server response contains both the 2Fo-Fc and the Fo-Fc blocks, so
  the two kinds share one cache entry instead of downloading the identical file
  twice. Which block to read is a display-time decision.
* Cryo-EM entries are found through their EMDB identifier, looked up from EMDB's
  search API with RCSB as a fallback. That lookup also yields the contour level
  the depositors recommend, which is how an EM map should be contoured; it is
  attached to the result whichever source supplied the voxels. The experimental
  method is never inferred from resolution, which BioJava parses incorrectly for
  some cryo-EM entries (biojava#1000).
* Ccp4Header checks for the MAP stamp at byte 208, so a server that answers with
  an error page and HTTP 200 produces a clean cache miss rather than a corrupt
  cache entry.
* Cached results are fully described by a .meta sidecar, so LOCAL_ONLY requests
  are served without opening a connection.

DemoFetchElectronDensity exercises all three outcomes: an X-ray entry, a cryo-EM
entry resolved through EMDB, and an entry deposited without structure factors.
Adds loadDensityMap and clearDensityMaps to JmolPanel, the single point every
viewer in the module already goes through, so a map fetched by DensityMapCache
can be contoured with one call. Surfaces are given stable ids so they can be
addressed or removed individually, and the state that "Reset Display" restores
is re-saved afterwards; otherwise the button would silently discard the map the
user had just asked for.

Three things had to be established by experiment against Jmol 14.31.10 rather
than assumed, and each changed the implementation:

* Option order in the isosurface command is load-bearing. With "mesh nofill"
  placed before the file name, Jmol accepts the command, reports no error, and
  draws nothing whatsoever. It has to follow the file name.

* A negative sigma does not contour at a negative level. Jmol reserves negative
  sigma for its own internal signalling, so "sigma -3.0" silently contours at
  the default level instead: the intended red lobe of a difference map came out
  identical to the blue one. Difference maps are therefore drawn as a single
  signed surface, which also matches Jmol's own shortcut for them.

* A density server response carries both the 2FO-FC and FO-FC blocks, and Jmol
  chooses between them by testing whether the file NAME contains "&diff=1". That
  marker normally arrives in the URL query string, which a cached local file does
  not have. Appending it to the file URL does not work, and neither does the
  "#diff=1" form the reader's own comment suggests: both were measured and both
  returned the 2FO-FC block or nothing. Embedding the marker in the file name
  does work, so the cache exposes the difference map under a companion name,
  hard-linked to the same bytes where the filesystem allows it.

Verified by contouring the real cached files headlessly: the 2Fo-Fc map at 1
sigma gives cutoff 0.356 over a -1.31 to 3.78 range, and the difference map at 3
sigma gives 0.374 over -0.69 to 0.85 - a different data block, which is what
proves the marker works rather than merely being accepted.

A map that cannot be contoured without a Fourier transform is rejected with an
explanatory exception rather than producing an empty surface.
Puts the feature in reach from the alignment viewer's View menu. The fetch runs
on a SwingWorker: even the smallest source is a few hundred kilobytes and a
full-resolution map can be far larger, so fetching on the event dispatch thread
would freeze the window for the duration.

Requests are built with allowNonRenderableFormats(false), which keeps the map
coefficient source out of the chain automatically rather than relying on the
viewer to notice it cannot draw the result.

When nothing is available the dialog explains why rather than listing HTTP
codes: for an entry whose every source returns 404 the likely reason is that no
structure factors were deposited and there is no associated EMDB map, which is
worth saying plainly. The per-source detail is still shown underneath.

AbstractAlignmentJmol gains getFrame() and setStatus() so a listener in the
neighbouring package can own its dialogs and report progress; both were
previously reachable only as protected fields.

DemoShowElectronDensity displays 1CBS with both maps clipped around the bound
retinoic acid.
Forty-three unit tests run with no network at all. The chain tests use stub
providers, so the behaviour that matters can be pinned exactly: a 404 falls
through to the next source, a too-large map falls through as well, and anything
else aborts. That last one is the point of the design - reporting a dropped
connection as "this entry has no density" would be worse than failing.

The offline set also pins the two-character directory rule against both
spellings of an entry (1cbs and pdb_00001cbs must land in "cb", not "db"), that
every source and kind combination maps to a distinct file, and that a LOCAL_ONLY
request is served entirely from disk with every server pointed at a dead port.

Six integration tests exercise the real services for a couple of megabytes in
total. The cryo-EM path is covered without downloading the 116 MB map: the EMDB
entry is resolved, the author contour level checked against its known value, and
the size guard then declines the full map before any of its body is transferred.
The coefficient test corrupts a downloaded file afterwards to confirm the
ETag-derived MD5 actually catches it rather than merely being recorded.

Writing the header check turned up a real bug: isCcp4 rejected any file shorter
than a CCP4 header, but a gzipped map compresses to a small fraction of the
header it contains, so small EMDB maps would have been rejected as invalid. The
length shortcut is gone; reading decides it.
The <table summary="..."> form is obsolete: the summary attribute was removed
in HTML5, and javadoc has generated HTML5 since JDK 15, so doclint rejects it
whenever it is switched on. The build sets -Xdoclint:none so this never broke
CI, but it would surface in the release profile and the attribute does nothing
for accessibility any more.

A definition list suits a list of placeholders and their meanings better than a
two-column table in any case.
The project is migrating to jupiter, so tests added by this branch should
not arrive as JUnit 4.

Imports move to org.junit.jupiter.api, @before and @after become
@beforeeach and @AfterEach, and the seventeen assertions carrying a
message have it moved from the first argument to the last, which is where
JUnit 5 expects it.

The three assertEquals(expected, actual, delta) calls are left alone: the
third argument there is a floating point tolerance, not a message, and
that overload is unchanged between the two versions.

No pom changes: both modules already declare junit-jupiter-engine and
junit-jupiter-params.
The coefficients were fetched through the divided archive path, built from
the two-character hash and the identifier. That path stops being correct in
July 2027, when the PDB moves to extended identifiers and a per-entry layout
under data/entries/<hash>/<extended id>/, and every file name changes with
it. The documented download endpoint resolves an entry by file name instead,
so it survives the move untouched:

  https://files.wwpdb.org/validation/download/1cbs_validation_2fo-fc_map_coef.cif.gz

Verified against files.wwpdb.org, files.rcsb.org and files-beta.wwpdb.org:
all three serve it, in both the short and the extended spelling, and the
beta host already serves it from the new archive.

This was the only provider here that ever constructed a path. The density
servers, PDBe and the EMDB archive are addressed by identifier already, so
nothing else in the package is exposed to the transition.

The default host stays files.wwpdb.org rather than files-beta.wwpdb.org,
which is deliberate and worth recording. The wwPDB describes the beta host
as transitional: on 21 July 2027 the beta archive replaces the main one,
after which the beta URL is supported by redirection for three years.
Pointing at it would be the choice that has to be revisited, twice, while
the main host simply becomes the new archive. Nothing is given up by
preferring the durable name either: the two serve byte-identical files
today, checked across nine entries and agreeing even on which ones 404.
The beta host is kept as a constant precisely because it already holds the
post-2027 content, which makes it the way to test this endpoint against the
archive as it will be rather than as it is - and it resolves both spellings
there, so the endpoint's semantics survive the cutover.

The identifier spelling is left to PdbId. getId(true) yields the short form
where an entry has one and the extended form otherwise, which is what the
entries deposited after the four-character space is exhausted will need;
hard-coding either spelling would replace a rule that adapts with a constant
that does not.

Mirrors that publish directories rather than an endpoint stay reachable.
EBI is one - it offers no name-resolving endpoint at all - so the divided
templates remain as DIVIDED_*_TEMPLATE constants, and the layout that
arrives in 2027 is expressible as ENTRIES_*_TEMPLATE. Expressing the latter
needed the extended identifier inside a template, which nothing provided:
{pdbid_lc} yields whichever spelling PdbId chose, and the new tree needs the
extended one in two positions regardless. Hence {extid} in UrlTemplates.
A mirror of the new archive is now a configuration change rather than a
release.

The cache layout is deliberately unchanged. Its two-character directory is
a way of spreading files over directories, not a copy of the archive's own
layout, and getMiddleHash counts from the right hand end, so 1cbs and
pdb_00001cbs both land in cb - which is also the rule the wwPDB documents
for the new archive, confirmed against it: entries/cb/pdb_00001cbs/ resolves
and entries/bs/pdb_00001cbs/ does not. Nothing there needs attention in
2027. The javadoc now says why the cache is not laid out like the archive:
it cannot be a mirror, since the archive publishes structure factors and
coefficients but never grids, and writing density into a directory that is
meant to be an exact copy of upstream puts it at the mercy of the next
rsync --delete. Should that be revisited, every cached path is computed in
that one class.
@aalhossary
aalhossary force-pushed the aa/electron-density-maps branch from d2cdae2 to 510798c Compare August 30, 2026 15:17
@aalhossary

Copy link
Copy Markdown
Member Author

Comment drafted with the help of Claude

Thanks — changed, and you were right that it reaches past the one file.

Here are some decisions I took, for your review

The endpoint.

Coefficients now go through https://files.wwpdb.org/validation/download/1cbs_validation_2fo-fc_map_coef.cif.gz. All six host × spelling combinations (wwpdb / rcsb / files-beta, short / extended) return 200 and serve the same bytes, and it holds across entries — 3alb, 1m3q, 2atk, 3hbx, 8ug2 all 200, while 6hu9 (cryo-EM) and 1a03 (NMR) 404, which is the right answer and is what raises NoDensityMapException. This was the only provider that ever built a path; the density servers, PDBe and EMDB are addressed by identifier already, so nothing else here is exposed.

Kept files, not files-beta

Even though beta serves the same endpoint with newer content. The beta page says it replaces the main archive on 21 July 2027, "after [which] the beta archive URL will be supported with redirection for 3 years" — so beta is the name that would need revisiting twice, while files.wwpdb.org simply becomes the new archive. Nothing is lost by preferring it: the two are byte-identical today across the nine entries I compared. Beta is kept as a constant for the opposite reason — it is the post-2027 content, so it is how I checked the endpoint's semantics survive the cutover.

The spelling is left to PdbId.

As urlId() delegates to PdbId.getId(true), which returns the short form where an entry has one and the extended form otherwise, I'm relying on that deliberately.

{mid} stays; {extid} is new.

EBI publishes no name-resolving endpoint — three plausible flat forms all 404 — so mirrors still need directory paths, and the divided path survives as DIVIDED_*_TEMPLATE. That exposed a real gap: the 2027 layout needs the extended ID in two positions regardless of the entry's own spelling, and {pdbid_lc} cannot produce it. Hence {extid}, which makes the new archive a configuration line rather than a patch:

entries/{mid}/{extid}/validation_reports/{extid}_validation_2fo-fc_map_coef.cif.gz

kept as ENTRIES_*_TEMPLATE, tested against the beta tree's real paths.

Cache layout

Since it is the one thing we can't revise after release: I considered mirroring the per-entry tree and decided against it. A cache cannot become a mirror — the archive publishes structure factors and coefficients but never grids, so density is fetched whatever else is rsynced. And putting our files inside a tree whose contract is "an exact copy of upstream" misbehaves under the documented rsync --delete: excluded files are protected, so the usual symptom is rsync failing to remove the directory and reporting it every sync rather than silent loss, but it becomes real deletion under --delete-excluded. For a whole-archive sweep it is also ~195k directories holding one or two small files each, against ~1,300 today. If you would rather converge anyway, deciding later is cheap: every cached path is computed in DensityCacheLayout and nowhere else. That reasoning is now in its javadoc.

@aalhossary
aalhossary requested a review from josemduarte August 30, 2026 16:30
aalhossary added a commit to aalhossary/biojava that referenced this pull request Aug 30, 2026
Covers what is merged since 7.2.6 and what is open and expected to land:
the download and checksum work, the CATH and ECOD fixes, the contact and
ASA performance tweaks, electron density, and the JUnit 5 migration.

Five entries are for pull requests that are still open - biojava#1134, biojava#1147,
biojava#1148, biojava#1150 and biojava#1151 - and should be checked against what actually
merged before the release is tagged.
The flat /validation/download/ endpoint returns neither an ETag nor a
Content-Length, so no digest can be recorded from it. The divided archive
path does return the content MD5, but that path disappears at the July 2027
archive transition, and the beta archive returns neither header on any path
on files.wwpdb.org or files.rcsb.org.

So the digest is asserted when it was recorded and reported when it was not,
rather than being required. Size validation applies either way, since the
sidecar is written from the bytes actually read.
@aalhossary

Copy link
Copy Markdown
Member Author

Following up on the switch to /validation/download/: it works, but it costs the checksum, and the integration test caught it.

The divided archive path returns the content MD5 as its ETag; the flat endpoint returns neither an ETag nor a Content-Length:

files.wwpdb.org/pub/pdb/validation_reports/cb/1cbs/…   Content-Length: 219757   ETag: "f99fb9d9…5745cf"
files.wwpdb.org/validation/download/…                  (neither)
files.rcsb.org/validation/download/…                   (neither)

That is not an argument for reverting. Checking all four mirrors, both archives:

Content-Length ETag
current wwPDB / RCSB, static /pub yes content MD5
current wwPDB / RCSB, /download/ shortlinks no no
beta wwPDB / RCSB, any path no no
beta PDBe / PDBj, per-entry yes size-mtime

So the content-MD5 disappears at the July 2027 transition whatever URL form we choose — it is a property of the static file paths, not of the divided layout. Pinning the divided path would buy about ten months and then need redoing.

d5a8f1fb9 therefore asserts the digest when the server offered one and reports when it did not, rather than requiring it:

No MD5 recorded for https://files.wwpdb.org/validation/download/1cbs_validation_2fo-fc_map_coef.cif.gz
  - the server offered no usable ETag. Size validation still applies.

Size validation is unaffected — the .size sidecar is written from the bytes actually read, so it exists either way, and the corruption check at the end of that test still passes. What is lost is detecting a same-length corruption. All six integration tests pass.

If the headers come back, the assertion starts running again with no code change.

Reported upstream. The gap is wider than map coefficients: files.rcsb.org/download/1cbs.cif.gz — the coordinate shortlink, in production today — also returns no ETag, Content-Length, Last-Modified or Accept-Ranges, and neither does any beta path on those two hosts, while PDBe and PDBj beta mirrors return all four. That means no integrity check, no If-Modified-Since, and no resumable downloads, on precisely the endpoints the wwPDB news of 15 August recommends for automated workflows. I have written to wwPDB and RCSB about it, with @josemduarte copied.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

Fetch, cache, load, and view electron density maps

2 participants