Skip to content

Fix documentation that contradicts the code - #62

Merged
HazenBabcock merged 13 commits into
ZhuangLab:masterfrom
HazenBabcock:fix-docs
Sep 8, 2026
Merged

Fix documentation that contradicts the code#62
HazenBabcock merged 13 commits into
ZhuangLab:masterfrom
HazenBabcock:fix-docs

Conversation

@HazenBabcock

@HazenBabcock HazenBabcock commented Sep 8, 2026

Copy link
Copy Markdown
Member

Twelve places where the documentation contradicts the code. One commit per
issue.

Nothing here changes behaviour. Every commit that touches a C file carries a
byte-identical object file as evidence, so "comment only" is demonstrated
rather than asserted.

Stale formats and names

  • a518b134 The Spliner README listed find_peaks_fista.py (renamed to
    find_peaks_decon.py, and no longer FISTA-only) and offset_to_Z.py (a
    case error), and never listed spline_info.py. The README and the
    directory now agree exactly, 21 entries each way.
  • 57cda576, 67fdf934 measure_psf.py and psf_localizations.py still
    described their input as a .bin molecule list. Both read HDF5 through
    saH5Py, and measure_psf.py's usage line contradicted its own --bin
    help text three lines below it.
  • b38b1f4c mFitAnscombe() documented a var parameter it has never
    taken.
  • 54e25262 Three more comment blocks naming renamed parameters: cmpr for
    compression twice in mlem_sparse.c, and filter for flt in
    matched_filter.c, where filter is also the struct type in the same
    signature.

Wrong descriptions

  • b19b1b04 mFitCalcErrFWLS() was described as the data weighted
    version, copied from mFitCalcErrDWLS() above it. That distinction is why
    both exist: DWLS divides by the measurement, FWLS by the fit. Also records
    why only FWLS can fail on a negative fit.
  • 77eef16f The ZCENTER convergence annotation said 0.01 where the
    multiplier gives 0.0001. The other six comment/multiplier pairs map
    exactly, which settles which of the two is wrong.
  • 275812c8 ZScaler's comment had the divisibility backwards relative to
    its own assertion and failure message.
  • 95e55d47 Four lines across three spline accessors named the wrong
    derivative; two of the three 3D accessors claimed to return the y
    derivative.
  • 5fd11214 The maxima finder's assertions named the opposite axis to the
    sizes it sets four lines above.

Missing conditions

  • 47e67fb1 "This library is thread safe.. Pretty sure.." is now the actual
    condition: safe per splineData, unsafe shared, because computeDelta2D/3D
    write scratch results into the struct that the accessors then read back.
  • 5d0f3db5 The width-vs-z curve has seven terms and the Z fitting model
    uses five, dropping any C and D silently. Now stated in all three places a
    reader can land. The default --fit_order 2 produces exactly five, so this
    only bites a deliberately higher-order calibration.
  • 4deff621 DriftFromFile documented z as nanometres; it is added straight
    onto h5_data["z"], which is microns. The diagnostics write the same array
    into both drift.txt and z_offset.txt, and the latter is documented as
    microns.

Full test suite passes, 277 tests.

HazenBabcock and others added 13 commits September 8, 2026 13:29
@HazenBabcock @claude
Make the Spliner README match the directory
a518b13 
Cross-checking every entry against the files on disk turned up four
problems:

  find_peaks_fista.py   listed, renamed to find_peaks_decon.py by
                        2a3a8fa, "we are no longer limited to just
                        FISTA". The description went stale with the
                        name: the module now has both a FISTAPeakFinder
                        and an ADMMPeakFinder.

  offset_to_Z.py        listed, the file is offset_to_z.py.

  spline_info.py        on disk since 01/18, never listed. It prints a
                        spline's details, which is not the same thing as
                        the listed print_psf.py, a debugging tool that
                        dumps PSF values.

  measure_psf_beads.py  described as taking "a text for input rather
                        than a molecule list file (.bin file)". The
                        contrast is still right, it does read bead
                        positions from a text file, but measure_psf.py
                        has read HDF5 through saH5Py for years.

With these fixed the README and the directory agree exactly, 21 entries
each way, nothing listed that is missing and nothing present that is
unlisted.

No test run: this is a text file that nothing in the test suite or the
packaging reads.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
@HazenBabcock @claude
Stop measure_psf.py calling its input a .bin file
57cda57 
Two Insight3 era references left in a script that has read HDF5 for
years:

    parser = argparse.ArgumentParser(description = 'Measure PSF given a
        movie, a list.bin file and (optionally) a z_offset file')

    # Load dax file, z offset file and molecule list file.

measurePSF() opens its localizations with saH5Py.SAH5Py(), and the
--bin argument's own help text already says "The name of the
localizations HDF5 file". So the usage line a user sees first
contradicts the help text three lines below it, and names a file format
the project has left.

Renaming the --bin argument itself would be a user visible change to the
command line, so that is left alone. This is the wording only.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
@HazenBabcock @claude
Stop psf_localizations.py calling its files molecule lists
67fdf93 
The same Insight3 era wording as 57cda57, in the multiplane script that
does the equivalent job. The module doc string described its input as "a
molecule list" and its output as "molecule lists", and the --frame help
text referred to "the frame in .bin file".

psfLocalizations() reads with saH5Py.SAH5Py() and writes with
saH5Py.saveLocalizations(), producing files named <base>_c1_psf.hdf5 and
so on, one per channel. Nothing about it has involved a .bin file for
years.

Also fixes "Giving a mapping file" to "Given", and says what the outputs
actually are, which is one localization file per channel rather than an
unspecified number of lists.

The --bin argument name is left alone, as renaming it would change the
command line.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
@HazenBabcock @claude
Stop mFitAnscombe documenting a parameter it does not take
b38b1f4 
The comment block listed two parameters:

     * x - Signal / function value.
     * var - (Gaussian) variance.
     */
    double mFitAnscombe(double x)

There is no var. The function has taken a single argument since it was
introduced in d0ab0dd (2018-05-04), so this was never a leftover from a
signature change, it was wrong when it was written. A generalized
Anscombe transform does take a variance, which is presumably where it
came from, but this one is the plain 2*sqrt(x + 3/8) and the code
matches that.

Comment only, and provably so: multi_fit.o is byte identical before and
after, md5 4b8058df50f0d1b692237950820df4aa.

Checked by scanning every comment block in every .c file in the tree
against the parameter list of the function below it. Excluding the
"Return" lines, which are return value documentation rather than
parameters, mFitAnscombe was one of four mismatches and is now absent
from the scan.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
@HazenBabcock @claude
Name the right parameters in three more comment blocks
54e2526 
The scan that found the mFitAnscombe defect in b38b1f4 turned up three
others, all the same shape but from renames rather than a parameter that
never existed:

  decon_storm/mlem_sparse.c   backwardCompressed() and
                              backwardCompressedFixedBg() both document
                              "cmpr", but take "compression". cmpr is a
                              local inside each, 1.0/(1.0 + compression),
                              so the comment describes the derived value
                              rather than the argument.

  sa_library/matched_filter.c cleanup() documents "filter", but takes
                              "flt". Worse than a plain mismatch, since
                              filter is also the name of the struct type
                              in its own signature, void cleanup(filter
                              *flt).

Corrects the names only. Both mlem_sparse functions also leave their
first parameter, dec, undocumented, but writing documentation that was
never there is a different job from fixing documentation that is wrong.

Comment only. Verified by compiling each file at the parent commit and
at this one and comparing the objects, which are byte identical:
matched_filter 0242ef0177f1e28e78e3288fd3ef06c9, mlem_sparse
65c6e99e4335d80848df5e5b9a6830f9. Both have to be compiled from the same
filename for this to mean anything, since gcc embeds the source name.

Neither has a standing object file to compare against, which is worth
knowing. matched_filter.c is compiled straight into libmatched_filter.so
with no intermediate .o kept, and decon_storm is not in SConstruct at
all, so mlem_sparse.c is never built as part of the package. That
directory is largely MATLAB.

Note that matched_filter.c has CRLF line endings, unlike almost
everything else here, so it was edited in binary. Its CRLF count is
unchanged at 239.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
@HazenBabcock @claude
Describe mFitCalcErrFWLS as the fit weighted version
b19b1b0 
Its comment read "The data weighted least squares version of the error
function", copied verbatim from mFitCalcErrDWLS immediately above it.
The distinction those words name is the entire reason both functions
exist:

    DWLS   err_i[k] = di*di/fit_data->x_data[k];    the measurement
    FWLS   err_i[k] = di*di/fi;                     the fit

Also says what follows from that. FWLS divides by the fit, so it checks
for a negative fit and returns 1 when it finds one, incrementing
n_neg_fi. DWLS divides by the measurement and needs no such guard. A
reader given the DWLS description has no way to understand why only one
of the two can fail.

The other three variants, ALS, LS and DWLS, all describe themselves
correctly, so FWLS was the only one.

Comment only: multi_fit.o is byte identical before and after, md5
4b8058df50f0d1b692237950820df4aa.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
@HazenBabcock @claude
Fix the ZCENTER delta in the convergence comment
77eef16 
Every threshold in mFitDeltaConvergence() is annotated with the delta it
represents at the default tolerance of 1e-6. Six of the seven agree with
their multiplier, and they agree exactly:

    HEIGHT      10000 * tol   0.01     comment 0.01
    BACKGROUND  10000 * tol   0.01     comment 0.01
    XCENTER      1000 * tol   0.001    comment 0.001
    YCENTER      1000 * tol   0.001    comment 0.001
    ZCENTER       100 * tol   0.0001   comment 0.01     <-
    XWIDTH       1000 * tol   0.001    comment 0.001
    YWIDTH       1000 * tol   0.001    comment 0.001

So 100 means 0.0001 and the comment is the thing that is wrong, not the
multiplier. The units say the same: x and y converge at 0.001 pixels,
which is 0.1nm at a 100nm pixel, and z at 0.0001 microns is also 0.1nm.
Believing the comment would make the z criterion a hundred times looser
than the lateral ones for no stated reason.

Changing the multiplier instead would be a change to when fitting stops,
so it would need measuring rather than reasoning, and does not belong in
a documentation branch. Note that none of this runs today in any case:
DELTA_CONVERGENCE defaults to 0, so the whole function is unreachable.

Comment only: multi_fit.o is byte identical before and after, md5
4b8058df50f0d1b692237950820df4aa.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
@HazenBabcock @claude
Say which of z_range and z_step divides the other
275812c 
The comment above ZScaler's assertion had the relationship backwards:

    # Assert that the z_step size is a multiple of the z_range.
    assert ((int(z_range*1.0e+3) % int(z_step*1.0e+3)) == 0), \
           "The z range must be a multiple of the z step."

z_range % z_step tests whether the range is a multiple of the step,
which is what the assertion message says, and the comment says the
reverse. Exercised directly:

    range=0.6 step=0.2   accepted   600 % 200 = 0
    range=0.5 step=0.3   rejected   500 % 300 = 200
    range=0.6 step=0.4   rejected   600 % 400 = 200

so a step larger than the range cannot be what is meant, and in fact the
assertion two lines above already rejects that case outright.

The code and the failure message were both right, only the comment was
wrong, which is the awkward case: someone debugging a failure reads the
comment first and it tells them to change the wrong number.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
@HazenBabcock @claude
State the condition under which the spline library is thread safe
47e67fb 
The header said:

     *  3. This library is thread safe.. Pretty sure..

Which is the least useful form the claim could take. It is unqualified,
so it reads as a guarantee, and it is hedged, so it invites someone to
test it. The truth is conditional and the condition matters.

computeDelta2D() and computeDelta3D() write their results into the
delta_f, delta_dxf, delta_dyf and delta_dzf arrays, which live in the
splineData struct rather than on the stack. The accessors then read them
back. So the evaluation of a point is a write followed by a read through
shared state, and two threads evaluating the same splineData will
overwrite each other's deltas in between.

Per splineData it is safe, and initSpline2D() and initSpline3D() copy
the coefficients into their own allocation rather than aliasing the
caller's, so two splineData share nothing at all. Nothing is broken
today, storm-analysis scales with multiprocessing rather than threads.

Also fixes note 1 in the same block, "but this whether or not this works
correctly has not been tested", which lost a word somewhere.

Comment only: cubic_spline.o is byte identical before and after, md5
b0b37b59a7cd4270e7014183ed032699.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
@HazenBabcock @claude
Say that the Z model uses five of the seven width terms
5d0f3db 
getWidthParams() returns seven values, wo, c, d, A, B, C and D, and
z_calibration.doFit() will fit all four of the higher order terms
(assert(n_additional < 5)). The consumers then disagree:

    sa_utilities/fitz.c        7 terms, uses params[5] and params[6]
    sa_utilities/fitz_c.py     5 terms, calcSxSy() stops at params[4]
    sa_library/dao_fit.h       5 terms, double wx_z_params[5]

daoInitializeZ() copies with for(i=0;i<5;i++), so a calibration with C
and D terms is applied in full when z is calculated after fitting for
the "3d" model, and silently truncated by the Z fitting model and by
calcSxSy(). Nothing warns.

Either could be the right design, the Z model may well not want the
extra freedom, so this does not change any of it. What it changes is
that you no longer have to read three files to find out. The note goes
in all three places a reader can land:

  - parameters.py, on the shared doc string for all seven parameters,
    which is where they are actually documented.
  - dao_fit.h, on the [5] arrays, which is where the truncation happens.
  - fitz_c.py, on calcSxSy(), the other five term consumer.

Worth recording that the default is safe: z_calibration takes
--fit_order 2, which produces exactly five terms, so this only bites
someone who deliberately fits a higher order curve.

Comment only: dao_fit.o is byte identical before and after, md5
b2dc5b13076ae6d61574c86d6345c5ad.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
@HazenBabcock @claude
Say that DriftFromFile takes z in microns
4deff62 
The doc string said nanometers. drift() adds the third column straight
onto h5_data["z"], which the PSF classes document as microns, so the
file has to be microns too.

The callers settle it. diagnostics/spliner_configure_common.py and
diagnostics/measure_psf both build

    dz = numpy.arange(-spline_z_range, spline_z_range + 0.001, 0.01)

and write that same array into two files: column 2 of drift.txt, and
column 1 of z_offset.txt. The z offset file is documented, in
measure_psf.py's own help text, as "the z offset of the frame relative
to the focal plane in microns". One array cannot be in two different
units, so drift.txt is microns.

The values agree. spline_z_range is 0.75, which as microns is a plus or
minus 750nm stack in 10nm steps, and as nanometers would be plus or
minus 0.75nm in steps of ten picometers.

Following the doc string would put every z offset out by 1000x. Also
says what the file's three columns are, which was not written down
anywhere.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
@HazenBabcock @claude
Name the right axis in the maxima finder's assertions
5fd1121 
findMaxima() sets the sizes correctly and then complains about the wrong
axis four lines later:

    self.flm_data.xsize = int(images[0].shape[1])
    self.flm_data.ysize = int(images[0].shape[0])
    ...
    assert (image.shape[0] == ...), "Unexpected image x size!"
    assert (image.shape[1] == ...), "Unexpected image y size!"

The assignments are right, and match the C, which indexes as
images[zi][yi*xsize + xi], so x is the fast axis and xsize is shape[1].
The messages say the opposite.

Message text only. The assertions compare an image against taken, so
they are checking the right things and would fire correctly. But anyone
who trips one is told the opposite of which axis is wrong, in the module
that defines the x/y convention for the whole peak finder, and this
package has two conventions in play: the simulator works in [x, y] order
internally and transposes on write, see the comment in
simulator/camera.py.

Note that admm/admm_math.py also labels shape[0] as X, but it is
self-consistent about it, mx is shape[0] and my is shape[1] throughout,
and it describes matrix cells rather than images, so it is left alone.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
@HazenBabcock @claude
Name the right derivative in the spline accessor comments
95e55d4 
Each d?fAt?D() comment block opens with a summary and closes with a
"Return" line naming which derivative it gives. Four of those were
wrong, in three of the five functions:

    dxfAt3D   Return - The derivative in y.            should be x
    dyfAt3D   "the derivative of the in spline y"      word order
    dzfAt3D   "the derivative of the in spline z"      word order
    dzfAt3D   Return - The derivative in y.            should be z

so of the five accessors only dxfAt2D and dyfAt2D described themselves
correctly, and two of the three 3D ones claimed to return the y
derivative.

These are the accessors underneath the Spliner Jacobian, so someone
checking the derivative terms against the comments finds two of them
apparently returning the wrong quantity, which is exactly the kind of
doubt that costs an afternoon.

The parallel d?fSpline?D() family, which takes double coordinates rather
than integer ones, is correct throughout and is left alone.

Comment only: cubic_spline.o is byte identical before and after, md5
b0b37b59a7cd4270e7014183ed032699.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
@HazenBabcock
HazenBabcock merged commit 4d25fa0 into ZhuangLab:master Sep 8, 2026
5 checks passed
@HazenBabcock
HazenBabcock deleted the fix-docs branch September 8, 2026 22:10
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Reviewers

No reviews

Assignees

No one assigned

Labels

None yet

Projects

None yet

Milestone

No milestone

Development

Successfully merging this pull request may close these issues.

1 participant

@HazenBabcock