From e381d2fd50973dc58a8019afb8f1620df3689eff Mon Sep 17 00:00:00 2001 From: Filip Leonarski Date: Tue, 8 Sep 2026 09:04:52 +0200 Subject: [PATCH] grid scan: review fixes - the ice channel sees its own spots, and a needle is not mirrored Two independent reviews of the merged grid-scan work. The findings that changed behaviour: The ice score's spot channel was fed a list the spot budget had already stripped. FilterSpotsByCount orders ice-band spots LAST when indexing is not to use them, so on a frame with more spots than the budget the ice spots are the first discarded - and the channel that exists for "ice arrives as discrete spots and leaves the radial profile flat" then read zero on exactly the frames it was written for. Probed at 3000 spots with 1200 on the hexagonal radii and a budget of 1000: 1.000 before the cap, 0.000 after. IceScore now takes d-spacings and is handed the list from before the cap. The viewer scaled the crystal box by the SIGNED grid step, where every other consumer takes the magnitude. On a negative step that mirrors the box - +30 deg drawn as -30 - and hands QRectF a negative width. rugnux --mode raster never put its settings on the experiment, so the indexing switch was read at its default while a deprecated per-run flag did the actual work; and RugnuxCommandLine emitted no --mode for Grid, so a raster job copied to a cluster ran the default mx - indexing, integrating and merging every cell of the raster. res_A is NaN where nothing in a blob measured a resolution, and nlohmann writes NaN as null, which the schema and the generated clients both reject. It is now left unset. The broker's configuration example named a key that does not exist (calibration, not calibration_settings); nlohmann ignores unknown keys, so a user copying it got a silently ignored block. The changelog had lost the rc.166 heading and 21 rc.167 entries to a bad edit of mine, and three entries had been filed under rc.166. Also: a warning where mode Grid meets a dataset with no grid scan, which was silent and indistinguishable from finding nothing; the viewer combo still named the retired ice_ring_score; and the claim that growth "cannot invent a crystal" was too strong - it cannot start a patch, but the cell count is read over the grown patch, so it does decide which patches pass. Co-Authored-By: Claude Opus 5 Claude-Session: https://claude.ai/code/session_01EFEJG6WBQv8th4UJFNe53N --- broker/OpenAPIConvert.cpp | 5 ++- common/JFJochMessages.h | 4 ++ docs/CHANGELOG.md | 28 +++++++++++- docs/CPU_DATA_ANALYSIS.md | 4 +- docs/JFJOCH_BROKER.md | 4 +- image_analysis/IceScore.cpp | 16 +++---- image_analysis/IceScore.h | 4 +- image_analysis/MXAnalysisWithoutFPGA.cpp | 2 +- image_analysis/spot_finding/SpotUtils.cpp | 7 +++ receiver/JFJochReceiver.cpp | 4 ++ receiver/JFJochReceiverFPGA.cpp | 2 +- rugnux/Rugnux.cpp | 2 +- rugnux/RugnuxCommandLine.cpp | 3 ++ rugnux/rugnux_cli.cpp | 4 +- tests/DetectionScoreTest.cpp | 49 +++++++++++++++++++-- viewer/JFJochViewerDatasetInfo.cpp | 2 +- viewer/image_viewer/JFJochGridScanImage.cpp | 7 ++- 17 files changed, 119 insertions(+), 28 deletions(-) diff --git a/broker/OpenAPIConvert.cpp b/broker/OpenAPIConvert.cpp index 2daabccf5..a5106e652 100644 --- a/broker/OpenAPIConvert.cpp +++ b/broker/OpenAPIConvert.cpp @@ -1323,7 +1323,10 @@ org::openapitools::server::model::Grid_scan_crystal Convert(const GridScanCrysta ret.setAngleDeg(input.angle_deg); ret.setScore(input.score); ret.setIceScore(input.ice_score); - ret.setResA(input.res_A); + // NaN means "nothing in the blob measured a resolution". nlohmann serialises NaN as null, + // which the schema and the generated clients reject, so the property is simply left unset. + if (std::isfinite(input.res_A)) + ret.setResA(input.res_A); ret.setNImages(input.n_images); return ret; } diff --git a/common/JFJochMessages.h b/common/JFJochMessages.h index 0bf2f0dbe..1b225a8a1 100644 --- a/common/JFJochMessages.h +++ b/common/JFJochMessages.h @@ -133,6 +133,10 @@ struct DataMessage { // unbounded above - the opposite convention to the two *_score members below, which is why the // name says ratio. Wire key "ice_ring_ratio", dataset /entry/MX/iceRingRatio; both readers also // accept the retired "ice_ring_score"/iceRingScore spellings. + // Spot d-spacings BEFORE the max-spot-count cap. FilterSpotsByCount orders ice-band spots LAST + // and so discards them first, which zeroes exactly the ice channel that exists for ice arriving + // as discrete spots. Transient: filled and consumed in the same process, never serialised. + std::vector spot_d_A_unfiltered; std::optional ice_ring_ratio; // How much of a single sweep's blind cone this orientation makes unrecoverable: 0 = one sweep // about the spindle reaches everything the point group can give, 1 = a short lattice row lies on diff --git a/docs/CHANGELOG.md b/docs/CHANGELOG.md index e77d0e85f..72a606008 100644 --- a/docs/CHANGELOG.md +++ b/docs/CHANGELOG.md @@ -7,8 +7,8 @@ * jfjoch_broker takes an explicit analysis mode - `MXStills` (the default), `Azint`, `Grid`, `PowderCalibration` or `None` - through `/config/analysis`, the `analysis` block of the configuration file and the Analysis panel of the web frontend. * The analysis mode takes precedence over `dataset_settings.spot_finding` and `spot_finding_settings.enable`, which are deprecated. * Rotation MX analysis is not available in jfjoch_broker, and is processed with rugnux. -* Each analysis method has its own settings endpoint, `/config/grid_scan_analysis` and `/config/calibration`; `/config/analysis` carries the mode alone. -* `PowderCalibration` computes the azimuthal profile on the CPU with a sector count of its own, at a few Hz. +* Each analysis method has its own settings endpoint, `/config/grid_scan_analysis` and `/config/calibration_settings`; `/config/analysis` carries the mode alone. +* `PowderCalibration` computes the azimuthal profile on the CPU with a sector count of its own. * A grid scan indexes each raster cell by default; turn it off with `indexing` in `/config/grid_scan_analysis`. * The crystals a grid scan found reach the CBOR end message as `grid_crystals`, the HDF5 master as `/entry/MX/crystals` and `GET /result/scan`, which also carries a powder-calibration result. * jfjoch_viewer draws a grid scan as one map - protein in orange, ice in cyan, mixed where both are present, over the shape of the loop - with a frame around each crystal found. @@ -17,6 +17,30 @@ * `rugnux --mode raster` reports the crystals in a stored grid scan - grid position, nearest image, measured size and orientation, score, ice score and resolution - as `_raster_report.txt` and `_raster.json`. * `rugnux --mode raster` takes `--beam-size`, `--raster-protein-threshold`, `--raster-grow-score`, `--raster-min-cells`, `--raster-decisive-score` and `--raster-max-crystals`. * The per-image ice ring measure is renamed `ice_ring_ratio` (`/entry/MX/iceRingRatio`, plot type `ice_ring_ratio`); the former `ice_ring_score` spelling is still read from older files and streams and still accepted by the API. +* `rugnux --model` reports CC(model, data) - the correlation of the merged intensities with the placed, scaled model - by resolution shell, on the same shells as CC1/2, with the reflection count and a significance for each. +* `rugnux --model` fits the model's scale, anisotropic B and bulk-solvent parameters on the working reflections only, so the R-free it reports is measured against a model no free reflection helped scale. +* The bulk-solvent parameters of `rugnux --model` are searched over their physically meaningful range instead of being fitted without bounds, so a model is never scaled with a solvent term that has silently switched itself off. +* The rigid-body placement of `rugnux --model` uses the same bounded bulk solvent as the reported fit, so a model is no longer placed against a target carrying a solvent term with no physical meaning. +* The rugnux results report opens with a summary - `VERDICT=` (`OK`, `WARNINGS`, `UNUSABLE`, `FAILED`), `VERDICT_TEXT=`, `PATHOLOGY_FLAGS=` with one closed-vocabulary code per condition that warned, and the `WARNING:` lines, which used to close the file - and the sections after it are renumbered 1-5 with no gaps. +* `rugnux --developer` writes the full results report - the pipeline-internal keys and the long explanations the default report now leaves out - and `--finalist-ledger` adds the evidence for every space group the search considered, not only the one it adopted. +* The results report warns when the merged data carry no usable signal and when too little of reciprocal space was measured inside the fitted resolution, and omits `FITTED_RESOLUTION` where the CC1/2 curve it is fitted on never falls off. +* rugnux detects translational pseudo-symmetry and reports it under the `PSEUDO_TRANSLATION` flag as `TNCS_DETECTED=` and the `TNCS_*` keys - a translation the merged data are exactly invariant under is reported as `UNDECLARED_LATTICE_TRANSLATION=` under `LATTICE_TRANSLATION` instead - and a detected pseudo-translation can no longer buy a false screw axis in the space-group search or hide a twin from the L-test (`L_TEST_VS_TNCS=`). +* The space-group search determines glide planes from zonal systematic absences, so a non-Sohncke space group such as P 2_1/c or Pbca is named where the run previously stopped at its Sohncke subgroup; `SOHNCKE_SPACE_GROUP=` carries the best Sohncke group beside it on every run that searched, and a centre of symmetry is never claimed. +* Where the cell metric carries more rotational symmetry than the Bravais class the indexer named, the extra rotations are put to the intensities and the space-group search is asked again on the metric's own cell - adopted only where the intensities confirm the higher symmetry - so a lattice that is nearly but not exactly hexagonal, or whose reduction landed in a sub-cell, still reaches its true point group. +* Systematic-absence calls rest on the evidence rather than on counts: a screw axis whose absent class the data show extinct is no longer refused because a handful of reflections in it read as present, and `SPACE_GROUP_ALTERNATIVES=` no longer drops a candidate that differs only on a zone the sweep never measured. +* A reference correlation measured on too few reflections is refused instead of scored zero, so a run given a reference MTZ is no longer reindexed on an operator that mapped almost everything outside the reference's coverage. +* A frame counts as indexed from 6 spots on its lattice rather than 9, so a weakly diffracting crystal whose frames cannot carry 9 is no longer refused the lattice it fits; `--min-indexed-spots` overrides it. +* `-C` accepts a known cell in any equivalent description - conventional or primitive, centred or not - instead of only the reduced primitive form, so a centred cell given the way it is published no longer makes the run report that it found no lattice. +* Each reflection is corrected for the sensor's quantum efficiency at the angle it meets the detector (attenuation lengths from the NIST tables, which also fixes the spot-width parallax term on CdTe) and for the attenuation of the flight path between the sample and its pixel; `--flight-path air|helium|vacuum` declares the medium - default air, since no file states it - and the report says what was assumed and what it was worth. The unmerged MTZ records the factors in new `QE` and `FLIGHT` columns beside `LP`, so raw counts are `I / LP * QE * FLIGHT`, and `_process.h5` in new optional `qe` and `flight` datasets. +* Rotation geometry post-refinement fits the crystal and the detector at once, against the observed spot positions and the observed rocking angles together, so the refined distance depends far less on how wrong the file's distance was. +* A coarsely sliced sweep integrates correctly: partials are joined into one rocking event by angle rather than by frame count, so two crossings of the Ewald sphere are no longer summed into one full, and at 0.5 degrees per image or coarser the per-frame geometry refinement accepts a spot whose miss the exposure's own rotation accounts for. +* `rugnux --mode scale` reports the detector tilt and direct beam of the geometry it re-scaled at, instead of zeros that read as a flat detector, and no longer warns that no image was indexed on a run whose lattice came from its input file. +* Every rotation run that determined a space group and merged reports what the mounting cost: `SPINDLE_LOST_UNIQUE_FRACTION=` is the fraction (0-1) of unique reflections the mounting made unmeasurable under the measured point group, also written to the master as `/entry/MX/spindleLostUniqueFraction` and what the mounting warning fires on; `SPINDLE_SYMMETRY_AXIS_ANGLE_DEG=` / `SPINDLE_SYMMETRY_AXIS_ORDER=` describe the mounting in the `--developer` report. +* Stills and grid scans carry a per-image `spindle_blind_fraction` - how much of a rotation sweep's blind cone this orientation would make unrecoverable, 0.5 and above calling for a second orientation - through the CBOR stream, HDF5 (`/entry/MX/spindleBlindFraction`), the plot and scan-result APIs, and the viewer and frontend plots; an absent value means the frame could not be assessed and is not a 0. +* The results report's `REPORT_VERSION` is 7. + +### 1.0.0-rc.166 + * `rugnux --model` treats the model as a hypothesis: it decides the enantiomorph and the indexing only where its R-work beats that of the same model in random orientations, and a model the data reject is still scored, placed and mapped, but leaves the reflection files byte for byte what a run with no model writes. * `rugnux --model` places the model against the data as a rigid body before scoring it, writes sigma_A-weighted 2mFo-DFc and mFo-DFc maps in place of the unweighted 2Fo-Fc and Fo-Fc, and writes the model as it was placed - `_model.cif`, and `_model.pdb` where the PDB format can express the cell - in the cell and space group of the reflection files beside it. * `rugnux` and `jfjoch_viewer` read PILATUS miniCBF sweeps natively, and open masters written at other facilities, including Eiger 1.x and third-party NXmx. diff --git a/docs/CPU_DATA_ANALYSIS.md b/docs/CPU_DATA_ANALYSIS.md index 3c9bca136..06e1e3cf3 100644 --- a/docs/CPU_DATA_ANALYSIS.md +++ b/docs/CPU_DATA_ANALYSIS.md @@ -53,8 +53,8 @@ The identifiers, the stored names and the wire names all follow it. The ice ring | Surface | Now | Retired spelling | |---|---|---| -| HDF5 | `/entry/MX/iceRingRatio`, `iceRingRatioMean` | `iceRingScore`, `iceRingScoreMean` — **still read**, so a file written before the rename still opens | -| CBOR | `ice_ring_ratio`, `ice_ring_ratio_mean` | `ice_ring_score`, `ice_ring_score_mean` — **still decoded**, so an older producer's stream still reads | +| HDF5 | `/entry/MX/iceRingRatio`, `iceRingRatioMean` | `iceRingScore`, `iceRingScoreMean` — **still read (the per-image array only - no `*Mean` scalar is read back by anything)**, so a file written before the rename still opens | +| CBOR | `ice_ring_ratio`, `ice_ring_ratio_mean` | `ice_ring_score`, `ice_ring_score_mean` — **still decoded**, so an older producer's stream still read (the per-image array only - no `*Mean` scalar is read back by anything)s | | REST `plot_type` | `ice_ring_ratio` | `ice_ring_score` — **still accepted**, and not scheduled for removal; both name the same plot | | CLI | `--ice-min-score` | *unchanged*, deliberately. It is a threshold on the ratio, and it is the one surface a user types | diff --git a/docs/JFJOCH_BROKER.md b/docs/JFJOCH_BROKER.md index 3f0c7460c..3e246903e 100644 --- a/docs/JFJOCH_BROKER.md +++ b/docs/JFJOCH_BROKER.md @@ -39,7 +39,7 @@ The mode takes precedence over `dataset_settings.spot_finding` and `spot_finding of which are deprecated. Each method carries its own settings, on its own endpoint: `/config/grid_scan_analysis` and -`/config/calibration`. `/config/analysis` carries the mode alone. +`/config/calibration_settings`. `/config/analysis` carries the mode alone. `PowderCalibration` computes the azimuthal profile on the CPU with a sector count of its own, so a calibration exposure runs at a few Hz rather than at the detector's frame rate. @@ -171,7 +171,7 @@ Example with all fields: "decisive_single_cell_score": 0.6, "indexing": true }, - "calibration": { + "calibration_settings": { "calibrant": "lab6", "method": "Rings" }, diff --git a/image_analysis/IceScore.cpp b/image_analysis/IceScore.cpp index 6257e4568..f35554ca2 100644 --- a/image_analysis/IceScore.cpp +++ b/image_analysis/IceScore.cpp @@ -258,21 +258,21 @@ float IceScoreRadial(const std::vector &profile, const std::vector return s * s / (s * s + S_HALF * S_HALF); } -float IceScoreSpots(const std::vector &spots, const std::vector &profile_count, +float IceScoreSpots(const std::vector &spot_d_A, const std::vector &profile_count, int32_t q_bins, const AzimuthalIntegrationSettings &settings, float half_width_q) { const int nq = std::max(q_bins, 0); const float low_q = settings.GetLowQ_recipA(); const float dq = settings.GetQSpacing_recipA(); - if (nq < 2 || !(dq > 0.0f) || !(half_width_q > 0.0f) || spots.empty() || profile_count.empty()) + if (nq < 2 || !(dq > 0.0f) || !(half_width_q > 0.0f) || spot_d_A.empty() || profile_count.empty()) return 0.0f; const std::vector count_q = FoldCountToQ(profile_count, nq); std::vector spot_q; - spot_q.reserve(spots.size()); - for (const auto &s: spots) - if (s.d_A > 0.0f) - spot_q.push_back(TWO_PI / s.d_A); + spot_q.reserve(spot_d_A.size()); + for (const float d: spot_d_A) + if (d > 0.0f) + spot_q.push_back(TWO_PI / d); std::ranges::sort(spot_q); std::vector band_q; @@ -326,7 +326,7 @@ float IceScoreSpots(const std::vector &spots, const std::vector &profile, const std::vector &profile_std, const std::vector &profile_count, int32_t q_bins, - const AzimuthalIntegrationSettings &settings, const std::vector &spots, + const AzimuthalIntegrationSettings &settings, const std::vector &spot_d_A, float half_width_q) { // The larger of the two, deliberately, and this is a decided question rather than an open one. // @@ -347,5 +347,5 @@ float IceScore(const std::vector &profile, const std::vector &prof // before scoring it (measured on the prototype: clean protein 20.5 -> 12.2%, ice 88.8 -> ~83%), not // to remove a channel. return std::max(IceScoreRadial(profile, profile_std, profile_count, q_bins, settings), - IceScoreSpots(spots, profile_count, q_bins, settings, half_width_q)); + IceScoreSpots(spot_d_A, profile_count, q_bins, settings, half_width_q)); } diff --git a/image_analysis/IceScore.h b/image_analysis/IceScore.h index 036508edd..903482691 100644 --- a/image_analysis/IceScore.h +++ b/image_analysis/IceScore.h @@ -58,12 +58,12 @@ float IceScoreRadial(const std::vector &profile, const std::vector // azimuthal bins; it is what makes the offsets comparable where the detector edge cuts a radius short. // half_width_q is the band half-width, the same one the spot finder marks ice rings with // (SpotFindingSettings::ice_ring_width_Q_recipA) - one width, measured over the corpus, not two. -float IceScoreSpots(const std::vector &spots, const std::vector &profile_count, +float IceScoreSpots(const std::vector &spot_d_A, const std::vector &profile_count, int32_t q_bins, const AzimuthalIntegrationSettings &settings, float half_width_q); // The score itself: whichever channel sees more. See the .cpp for why the maximum and not something // gentler, and for the lever to reach for if the rate ever has to come down. float IceScore(const std::vector &profile, const std::vector &profile_std, const std::vector &profile_count, int32_t q_bins, - const AzimuthalIntegrationSettings &settings, const std::vector &spots, + const AzimuthalIntegrationSettings &settings, const std::vector &spot_d_A, float half_width_q); diff --git a/image_analysis/MXAnalysisWithoutFPGA.cpp b/image_analysis/MXAnalysisWithoutFPGA.cpp index a0fefefd9..4a8df0472 100644 --- a/image_analysis/MXAnalysisWithoutFPGA.cpp +++ b/image_analysis/MXAnalysisWithoutFPGA.cpp @@ -298,7 +298,7 @@ void MXAnalysisWithoutFPGA::Analyze(DataMessage &output, // carry. Its second channel reads the spots, which SpotAnalyze has already put in the message. output.ice_score = IceScore(output.az_int_profile, output.az_int_profile_std, output.az_int_profile_count, integration.GetQBinCount(), - integration.Settings(), output.spots, + integration.Settings(), output.spot_d_A_unfiltered, spot_finding_settings.ice_ring_width_Q_recipA); } diff --git a/image_analysis/spot_finding/SpotUtils.cpp b/image_analysis/spot_finding/SpotUtils.cpp index 23b78ee77..68849a438 100644 --- a/image_analysis/spot_finding/SpotUtils.cpp +++ b/image_analysis/spot_finding/SpotUtils.cpp @@ -284,6 +284,13 @@ void SpotAnalyze(const DiffractionExperiment &experiment, output.resolution_estimate = GetResolution(spots_out); output.protein_score = ProteinScore(spots_out); + // The ice score's spot channel has to see the ice spots, and the cap below throws them away + // first, so its input is taken here rather than from output.spots. + output.spot_d_A_unfiltered.clear(); + output.spot_d_A_unfiltered.reserve(spots_out.size()); + for (const auto &s: spots_out) + output.spot_d_A_unfiltered.push_back(s.d_A); + // One decision drives both: if indexing is to use the ice-band spots, the spot budget must not // throw them away before it gets the chance. FilterSpotsByCount(spots_out, experiment.GetMaxSpotCount(), diff --git a/receiver/JFJochReceiver.cpp b/receiver/JFJochReceiver.cpp index e3be6bf6e..aabf32887 100644 --- a/receiver/JFJochReceiver.cpp +++ b/receiver/JFJochReceiver.cpp @@ -195,6 +195,10 @@ void JFJochReceiver::SendEndMessage() { // is often the shoulder of a better one two rows further down, and a scan exists precisely to see // the whole loop before choosing. Done here so one list reaches both the written file (through // this end message) and the API (through GetFinalStatistics). + if (experiment.GetAnalysisMode() == AnalysisMode::Grid && !experiment.GetGridScan().has_value()) + logger.Warning("Analysis mode is Grid but the dataset has no grid scan: images were scored, " + "but there is no map to find crystals in and the crystal list is empty"); + if (const auto grid = experiment.GetGridScan(); grid.has_value() && experiment.GetAnalysisMode() == AnalysisMode::Grid) { grid_scan_result = AnalyzeGridScan(scan_result.GetResult(), grid.value(), diff --git a/receiver/JFJochReceiverFPGA.cpp b/receiver/JFJochReceiverFPGA.cpp index 72fc23be8..3ff265d86 100644 --- a/receiver/JFJochReceiverFPGA.cpp +++ b/receiver/JFJochReceiverFPGA.cpp @@ -443,7 +443,7 @@ void JFJochReceiverFPGA::FrameTransformationThread(uint32_t threadid) { message.ice_score = IceScore(message.az_int_profile, message.az_int_profile_std, message.az_int_profile_count, experiment.GetAzimuthalIntegrationSettings().GetQBinCount(), - experiment.GetAzimuthalIntegrationSettings(), message.spots, + experiment.GetAzimuthalIntegrationSettings(), message.spot_d_A_unfiltered, spot_finding_settings.ice_ring_width_Q_recipA); scan_result.Add(message); diff --git a/rugnux/Rugnux.cpp b/rugnux/Rugnux.cpp index be95ee5cd..e21720e41 100644 --- a/rugnux/Rugnux.cpp +++ b/rugnux/Rugnux.cpp @@ -3456,7 +3456,7 @@ ProcessResult Rugnux::RunPipeline(RugnuxObserver *observer, bool write_output, b // ice score has. msg.ice_score = IceScore(msg.az_int_profile, msg.az_int_profile_std, msg.az_int_profile_count, mapping.GetQBinCount(), - mapping.Settings(), msg.spots, + mapping.Settings(), msg.spot_d_A_unfiltered, config_.spot_finding.ice_ring_width_Q_recipA); msg.run_number = experiment_.GetRunNumber(); msg.run_name = experiment_.GetRunName(); diff --git a/rugnux/RugnuxCommandLine.cpp b/rugnux/RugnuxCommandLine.cpp index d16c8c14b..fc2568ee1 100644 --- a/rugnux/RugnuxCommandLine.cpp +++ b/rugnux/RugnuxCommandLine.cpp @@ -63,6 +63,9 @@ std::string RugnuxCommandLine(const ProcessConfig &config, } else if (calibration) { args.emplace_back("--mode"); args.emplace_back("calibration"); + } else if (config.mode == AnalysisMode::Grid) { + args.emplace_back("--mode"); + args.emplace_back("raster"); } auto add = [&](const std::string &flag, const std::string &val) { diff --git a/rugnux/rugnux_cli.cpp b/rugnux/rugnux_cli.cpp index 3eecdc6f9..0a8ca1598 100644 --- a/rugnux/rugnux_cli.cpp +++ b/rugnux/rugnux_cli.cpp @@ -2199,8 +2199,10 @@ static int RunRugnux(int argc, char **argv) { // scores, and is worth having. config.write_process_h5 = write_process_h5_flag; + experiment.ImportGridScanAnalysisSettings(raster_analysis); + config.spot_finding.enable = true; - config.spot_finding.indexing = false; + config.spot_finding.indexing = raster_analysis.IsIndexing(); config.spot_finding.signal_to_noise_threshold = sigma_spot_finding; config.spot_finding.photon_count_threshold = photon_count_threshold_spot_finding; config.spot_finding.min_pix_per_spot = min_pix_per_spot.value_or(2); diff --git a/tests/DetectionScoreTest.cpp b/tests/DetectionScoreTest.cpp index 64830cf7a..7f89b3107 100644 --- a/tests/DetectionScoreTest.cpp +++ b/tests/DetectionScoreTest.cpp @@ -147,6 +147,18 @@ TEST_CASE("IceScoreRadial_AzimuthalProfileFoldsToTheSameAnswer") { == Catch::Approx(IceScoreRadial(flat, sigma, count, q_bins, settings))); } +namespace { + // IceScoreSpots reads d-spacings, not spots: it is fed the list from BEFORE the max-spot-count + // cap, which orders ice-band spots last and would otherwise discard the very spots it counts. + std::vector d_of(const std::vector &spots) { + std::vector d; + d.reserve(spots.size()); + for (const auto &s: spots) + d.push_back(s.d_A); + return d; + } +} + TEST_CASE("IceScoreSpots_ExcessOnTheIceRadii") { const auto settings = ice_settings(); const int q_bins = settings.GetQBinCount(); @@ -160,21 +172,21 @@ TEST_CASE("IceScoreSpots_ExcessOnTheIceRadii") { const float q = q_lo + (q_hi - q_lo) * static_cast(i) / 399.0f; even.push_back(spot(TWO_PI / q, 1000.0f)); } - CHECK(IceScoreSpots(even, count, q_bins, settings, ICE_W) < 0.5f); + CHECK(IceScoreSpots(d_of(even), count, q_bins, settings, ICE_W) < 0.5f); // The same frame with 10 extra spots planted on each hexagonal radius. std::vector with_ice = even; for (const float d: ICE_RING_RES_A) for (int i = 0; i < 10; i++) with_ice.push_back(spot(d, 1000.0f)); - CHECK(IceScoreSpots(with_ice, count, q_bins, settings, ICE_W) > 0.5f); + CHECK(IceScoreSpots(d_of(with_ice), count, q_bins, settings, ICE_W) > 0.5f); // Two spots that both happen to sit on a ring are not ice: the ratio term refuses them even // though the Poisson tail alone would not. std::vector two; two.push_back(spot(ICE_RING_RES_A[0], 1000.0f)); two.push_back(spot(ICE_RING_RES_A[1], 1000.0f)); - CHECK(IceScoreSpots(two, count, q_bins, settings, ICE_W) < 0.5f); + CHECK(IceScoreSpots(d_of(two), count, q_bins, settings, ICE_W) < 0.5f); CHECK(IceScoreSpots({}, count, q_bins, settings, ICE_W) == 0.0f); } @@ -203,5 +215,34 @@ TEST_CASE("IceScore_TakesTheStrongerChannel") { for (int i = 0; i < 10; i++) textured.push_back(spot(d, 1000.0f)); CHECK(IceScoreRadial(flat, sigma, count, q_bins, settings) == 0.0f); - CHECK(IceScore(flat, sigma, count, q_bins, settings, textured, ICE_W) > 0.5f); + CHECK(IceScore(flat, sigma, count, q_bins, settings, d_of(textured), ICE_W) > 0.5f); +} + +// The spot channel must be fed the list from BEFORE the max-spot-count cap. FilterSpotsByCount +// orders ice-band spots LAST when indexing is not to use them, so on a frame with more spots than +// the budget it discards the ice first - and the channel that exists for ice arriving as discrete +// spots then reads zero on exactly the frames it was written for. +TEST_CASE("IceScoreSpots_ReadsThePreCapList") { + const auto settings = ice_settings(); + const int q_bins = settings.GetQBinCount(); + const std::vector count(q_bins, 10000); + + std::vector spots; + const float q_lo = 1.3f, q_hi = 4.2f; + for (int i = 0; i < 1800; i++) { + const float q = q_lo + (q_hi - q_lo) * static_cast(i) / 1799.0f; + spots.push_back(spot(TWO_PI / q, 1000.0f)); + } + for (const float d: ICE_RING_RES_A) + for (int i = 0; i < 60; i++) + spots.push_back(spot(d, 1000.0f)); + MarkIceRings(spots, ICE_W); + + const std::vector before = d_of(spots); + FilterSpotsByCount(spots, 1000, true); + const std::vector after = d_of(spots); + + REQUIRE(before.size() > after.size()); + CHECK(IceScoreSpots(before, count, q_bins, settings, ICE_W) > 0.5f); + CHECK(IceScoreSpots(after, count, q_bins, settings, ICE_W) < 0.5f); } diff --git a/viewer/JFJochViewerDatasetInfo.cpp b/viewer/JFJochViewerDatasetInfo.cpp index dc3adfcea..584cb1be7 100644 --- a/viewer/JFJochViewerDatasetInfo.cpp +++ b/viewer/JFJochViewerDatasetInfo.cpp @@ -145,7 +145,7 @@ void JFJochViewerDatasetInfo::UpdateLabels() { if (dataset->experiment.GetGridScan()) combo_box->addItem("Protein / ice map", kCompositeMetric); combo_box->addItem("Background estimate", 0); - combo_box->addItem("Ice ring score", 14); + combo_box->addItem("Ice ring ratio", 14); combo_box->addItem("Spindle blind fraction", 15); combo_box->addItem("Resolution estimate", 7); combo_box->addItem("Spot count", 1); diff --git a/viewer/image_viewer/JFJochGridScanImage.cpp b/viewer/image_viewer/JFJochGridScanImage.cpp index 23d5f3587..12d24edc3 100644 --- a/viewer/image_viewer/JFJochGridScanImage.cpp +++ b/viewer/image_viewer/JFJochGridScanImage.cpp @@ -246,8 +246,11 @@ void JFJochGridScanImage::drawCrystalFrames() { QPen line(Qt::white, 2); line.setCosmetic(true); - const float step_x = settings->GetGridStepX_um(); - const float step_y = settings->GetGridStepY_um(); + // Magnitudes: the step sign says which way the stage moved, and the display grid has already + // applied it. Using the signed value here mirrors the box on a negative step, drawing +30 deg + // as -30 and handing QRectF a negative width. + const float step_x = std::fabs(settings->GetGridStepX_um()); + const float step_y = std::fabs(settings->GetGridStepY_um()); for (const auto &c : crystals) { if (grid_scan::IsRoundBlob(c.major_um, c.minor_um)) {