still_process fails spotfinding on image(s) in middle of range

[ndevenish]

Using https://zenodo.org/record/55058 as a test data set:

$ dials.stills_process test_1/test_1_0400.cbf.bz2
....
Error spotfinding test_1_0400.cbf 3

However, normal dials import works:

$ dials.import test_1/test_1_0400.cbf.bz2 && dials.find_spots imported.expt
...
Found 19 strong pixels on image 400
...
Saved 3 reflections to strong.refl

[phyy-nx]

Debugging now.

dials.stills_process test_1/test_1_0400.cbf.bz2 squash_errors=False

Error spotfinding test_1_0400.cbf 3
Traceback (most recent call last):
  File "/net/dials/raid1/aaron/dialsbuild/build/../modules/dials/command_line/stills_process.py", line 1482, in <module>
    script.run()
  File "/net/dials/raid1/aaron/dialsbuild/build/../modules/dials/command_line/stills_process.py", line 647, in run
    do_work(0, iterable)
  File "/net/dials/raid1/aaron/dialsbuild/build/../modules/dials/command_line/stills_process.py", line 490, in do_work
    processor.process_experiments(item[0], item[1])
  File "/net/dials/raid1/aaron/dialsbuild/build/../modules/dials/command_line/stills_process.py", line 820, in process_experiments
    observed = self.find_spots(experiments)
  File "/net/dials/raid1/aaron/dialsbuild/build/../modules/dials/command_line/stills_process.py", line 902, in find_spots
    observed = flex.reflection_table.from_observations(experiments, self.params)
  File "/net/dials/raid1/aaron/dialsbuild/modules/dials/array_family/flex_ext.py", line 173, in from_observations
    return find_spots(experiments)
  File "/net/dials/raid1/aaron/dialsbuild/modules/dials/algorithms/spot_finding/finder.py", line 776, in __call__
    assert missed.count(True) == 0, missed.count(True)
AssertionError: Please report this error to dials-support@lists.sourceforge.net: 3

[phyy-nx]

Thanks, bugfixed. Here are the results for the whole dataset:

dials.stills_process ../data/* mp.nproc=64 known_symmetry.space_group=P41212 known_symmetry.unit_cell=57.782,57.782,150.004,90,90,90
cctbx.xfel.plot_run_stats_from_experiments .

[graeme-winter]

This has broken handling of image sequences - maybe not adequately tested though as requires big-ish data.

230 experiments in an image sequence - with this change id == 0 for all -

Silver-Surfer tst-sp-fix :) $ dials.show strong.refl 
DIALS (2018) Acta Cryst. D74, 85-97. https://doi.org/10.1107/S2059798317017235
The following parameters have been modified:

input {
  reflections = strong.refl
}


Reflection list contains 2578 reflections
+------------------------+------------------------+--------------------------+--------------------------+
| Column                 | min                    | max                      | mean                     |
|------------------------+------------------------+--------------------------+--------------------------|
| flags                  | 32                     | 32                       | 32                       |
| id                     | 0                      | 0                        | 0                        |
| intensity.sum.value    | 7.0                    | 11676.0                  | 333.6                    |
| intensity.sum.variance | 7.0                    | 11676.0                  | 333.6                    |
| n_signal               | 3.0                    | 393.0                    | 11.122963537626067       |
| panel                  | 0                      | 0                        | 0                        |
| shoebox                |                        |                          |                          |
| summed I               | 7.0                    | 11676.0                  | 333.6                    |
| N pix                  | 3.0                    | 1643.0                   | 28.6                     |
| N valid foreground pix | 3.0                    | 393.0                    | 11.1                     |
| xyzobs.px.value        | 78.56, 42.50, 0.50     | 2448.67, 2442.50, 229.50 | 1231.04, 1211.16, 111.53 |
| xyzobs.px.variance     | 0.0833, 0.0833, 0.0833 | 0.1483, 0.1483, 0.0833   | 0.0905, 0.0929, 0.0833   |
+------------------------+------------------------+--------------------------+--------------------------+

should be

Silver-Surfer tst-sp-fix :) $ dials.show strong.refl 
DIALS (2018) Acta Cryst. D74, 85-97. https://doi.org/10.1107/S2059798317017235
The following parameters have been modified:

input {
  reflections = strong.refl
}


Reflection list contains 2578 reflections
+------------------------+------------------------+--------------------------+--------------------------+
| Column                 | min                    | max                      | mean                     |
|------------------------+------------------------+--------------------------+--------------------------|
| flags                  | 32                     | 32                       | 32                       |
| id                     | 0                      | 229                      | 111                      |
| intensity.sum.value    | 7.0                    | 11676.0                  | 333.6                    |
| intensity.sum.variance | 7.0                    | 11676.0                  | 333.6                    |
| n_signal               | 3.0                    | 393.0                    | 11.122963537626067       |
| panel                  | 0                      | 0                        | 0                        |
| shoebox                |                        |                          |                          |
| summed I               | 7.0                    | 11676.0                  | 333.6                    |
| N pix                  | 3.0                    | 1643.0                   | 28.6                     |
| N valid foreground pix | 3.0                    | 393.0                    | 11.1                     |
| xyzobs.px.value        | 78.56, 42.50, 0.50     | 2448.67, 2442.50, 229.50 | 1231.04, 1211.16, 111.53 |
| xyzobs.px.variance     | 0.0833, 0.0833, 0.0833 | 0.1483, 0.1483, 0.0833   | 0.0905, 0.0929, 0.0833   |
+------------------------+------------------------+--------------------------+--------------------------+

i.e. adding not is_still() negated the purpose of the loop.

[graeme-winter]

For the moment I have reverted the commit - will add a proper test

[graeme-winter]

@phyy-nx to explain thinking here - the spot finding works old-style on the image set then we iterate through this and assign spots to the correct experiments afterwards - this allows plenty of speed and parallelism as well as fast start up time as their is only one image set, indexed into by the scan.

[graeme-winter]

N.B. test fails with the change above:

Silver-Surfer command_line :) [master] $ pytest --regression test_find_spots.py::test_find_spots_from_imported_as_grid
Test session starts (platform: darwin, Python 3.6.7, pytest 4.6.4, pytest-sugar 0.9.2)
rootdir: /Users/graeme/svn/cctbx/modules/dials, inifile: pytest.ini
plugins: dials-data-2.0.69, sugar-0.9.2, mock-2.0.0, forked-1.1.2, xdist-1.31.0
collecting ... 

―――――――――――――――――――― test_find_spots_from_imported_as_grid ―――――――――――――――――――――

dials_data = <DataFetcher: /Users/graeme/svn/cctbx/build/dials_data>
tmpdir = local('/private/var/folders/hs/gvssfcd910s0jczjhjczj6h80000gn/T/pytest-of-graeme/pytest-12/test_find_spots_from_imported_0')

    def test_find_spots_from_imported_as_grid(dials_data, tmpdir):
        """First run import to generate an imported.expt and use this."""
        _ = procrunner.run(
            ["dials.import", "oscillation=0,0"]
            + [
                f.strpath for f in dials_data("centroid_test_data").listdir("centroid*.cbf")
            ],
            working_directory=tmpdir.strpath,
        )
    
        result = procrunner.run(
            [
                "dials.find_spots",
                tmpdir.join("imported.expt").strpath,
                "output.reflections=spotfinder.refl",
                "output.shoeboxes=True",
                "algorithm=dispersion",
            ],
            working_directory=tmpdir.strpath,
        )
        assert not result.returncode and not result.stderr
        assert tmpdir.join("spotfinder.refl").check(file=1)
    
        reflections = flex.reflection_table.from_file(tmpdir / "spotfinder.refl")
    
>       assert len(set(reflections["id"])) == 9, len(set(reflections["id"]))
E       AssertionError: 1
E       assert 1 == 9
E        +  where 1 = len({0})
E        +    where {0} = set(<scitbx_array_fami...ect at 0x1118ee870>)

test_find_spots.py:167: AssertionError
----------------------------- Captured stdout call -----------------------------
Warning: Ambiguous parameter definition: oscillation = 0,0
Best matches:
  input.tolerance.scan.oscillation
  geometry.scan.oscillation
Assuming geometry.scan.oscillation was intended.
DIALS (2018) Acta Cryst. D74, 85-97. https://doi.org/10.1107/S2059798317017235
DIALS 2.dev.1431-gb88160f3c
Warning: Ambiguous parameter definition: oscillation = 0,0
Best matches:
  input.tolerance.scan.oscillation
  geometry.scan.oscillation
Assuming geometry.scan.oscillation was intended.
The following parameters have been modified:

input {
  experiments = <image files>
}
geometry {
  scan {
    oscillation = 0 0
  }
}


Applying input geometry in the following order:
  1. Manual geometry

--------------------------------------------------------------------------------
  format: <class 'dxtbx.format.FormatCBFMiniPilatus.FormatCBFMiniPilatus'>
  num images: 9
  sequences:
    still:    1
    sweep:    0
  num stills: 0
--------------------------------------------------------------------------------
Writing experiments to imported.expt
DIALS (2018) Acta Cryst. D74, 85-97. https://doi.org/10.1107/S2059798317017235
DIALS 2.dev.1431-gb88160f3c
The following parameters have been modified:

output {
  reflections = "spotfinder.refl"
}
spotfinder {
  threshold {
    algorithm = *dispersion dispersion_extended
  }
}
input {
  experiments = /private/var/folders/hs/gvssfcd910s0jczjhjczj6h80000gn/T/pytest-of-graeme/pytest-12/test_find_spots_from_imported_0/imported.expt
}

Setting spotfinder.filter.min_spot_size=3
Configuring spot finder from input parameters
--------------------------------------------------------------------------------
Finding strong spots in imageset 0
--------------------------------------------------------------------------------


Finding spots in image 1 to 9...
Setting chunksize=9
Extracting strong pixels from images
 Using multiprocessing with 1 parallel job(s)

Found 1276 strong pixels on image 1
Found 1372 strong pixels on image 2
Found 1351 strong pixels on image 3
Found 1346 strong pixels on image 4
Found 1336 strong pixels on image 5
Found 1354 strong pixels on image 6
Found 1329 strong pixels on image 7
Found 1348 strong pixels on image 8
Found 1310 strong pixels on image 9

Extracted 3516 spots
Removed 1697 spots with size < 3 pixels
Removed 0 spots with size > 1000 pixels
Calculated 1819 spot centroids
Calculated 1819 spot intensities
Filtered 1812 of 1819 spots by peak-centroid distance

Histogram of per-image spot count for imageset 8:
1812 spots found on 8 images (max 301 / bin)
*      *
*      *
*      *
********
********
********
********
********
********
********
1      8

--------------------------------------------------------------------------------
Saved 1812 reflections to spotfinder.refl

 test/command_line/test_find_spots.py ⨯                          100% ██████████

Results (5.59s):
       1 failed
         - test/command_line/test_find_spots.py:142 test_find_spots_from_imported_as_grid

[graeme-winter]

Trivial reproducer of error:

dials.stills_process squash_errors=0 convert_sequences_to_stills=1 ins_1d_grid_004?.cbf.gz

exam.zip

=>

Grey-Area work :) $ dials.stills_process squash_errors=0 convert_sequences_to_stills=1 ins_1d_grid_004?.cbf.gz
Have 10 files
DIALS (2018) Acta Cryst. D74, 85-97. https://doi.org/10.1107/S2059798317017235
The following parameters have been modified:

dispatch {
  squash_errors = 0
}
geometry {
  convert_sequences_to_stills = 1
}

Loading files...
Loading ins_1d_grid_0040.cbf.gz
********************************************************************************
Finding Strong Spots
********************************************************************************
Setting spotfinder.filter.min_spot_size=3
Configuring spot finder from input parameters
--------------------------------------------------------------------------------
Finding strong spots in imageset 0
--------------------------------------------------------------------------------


Finding spots in image 1 to 1...
Setting chunksize=1
Extracting strong pixels from images
 Using multiprocessing with 1 parallel job(s)

Found 803 strong pixels on image 1

Extracted 46 spots
Removed 7 spots with size < 3 pixels
Removed 0 spots with size > 1000 pixels
Calculated 39 spot centroids
Calculated 39 spot intensities
Filtered 39 of 39 spots by peak-centroid distance
Error spotfinding ins_1d_grid_0040.cbf 39
Traceback (most recent call last):
  File "/Users/graeme/svn/cctbx/build/../modules/dials/command_line/stills_process.py", line 1482, in <module>
    script.run()
  File "/Users/graeme/svn/cctbx/build/../modules/dials/command_line/stills_process.py", line 647, in run
    do_work(0, iterable)
  File "/Users/graeme/svn/cctbx/build/../modules/dials/command_line/stills_process.py", line 563, in do_work
    processor.process_experiments(tag, experiments)
  File "/Users/graeme/svn/cctbx/build/../modules/dials/command_line/stills_process.py", line 820, in process_experiments
    observed = self.find_spots(experiments)
  File "/Users/graeme/svn/cctbx/build/../modules/dials/command_line/stills_process.py", line 902, in find_spots
    observed = flex.reflection_table.from_observations(experiments, self.params)
  File "/Users/graeme/svn/cctbx/modules/dials/array_family/flex_ext.py", line 173, in from_observations
    return find_spots(experiments)
  File "/Users/graeme/svn/cctbx/modules/dials/algorithms/spot_finding/finder.py", line 776, in __call__
    assert missed.count(True) == 0, missed.count(True)
AssertionError: Please report this error to dials-support@lists.sourceforge.net: 39

So we can at least look at this very easily.

[ndevenish]

I've just had a user encounter this error again, do we have any ideas on fixing at the moment?

[phyy-nx]

My original fix broke grid scans and so was reverted. We also need a reproducer of that error.

[graeme-winter]

My original fix broke grid scans and so was reverted. We also need a reproducer of that error.

@phyy-nx reproducer included above at #1180 (comment)

[phyy-nx]

Just wanted to note, working with @mewall I added a commit to the lunus_stills_process branch that fixes this issue in kind of a kludgy way. Comments invited on commit c570d44. Specifically, @graeme-winter, can you verify that grid scan processing is not borked on the lunus_stills_process branch?

[phyy-nx]

The test @graeme-winter added to prevent regression of grid scan spotfinding
passes on the lunus_stills_process branch using commit c570d44:

libtbx.pytest ./test/command_line/test_find_spots.py::test_find_spots_from_imported_as_grid

Further, the example dataset above (#1180 (comment)) also works using commit c570d44. Command I used which indexed a couple images (didn't try and optimize it):

dials.stills_process squash_errors=1 convert_sequences_to_stills=1 ins_1d_grid_004?.cbf.gz known_symmetry.unit_cell=69,69,69,109,109,109

I'd like to cherry pick that commit into a new PR, but I'd like some review here first. I imagine the flag I added might be better served as a phil parameter set by default on dials.stills_process. @graeme-winter?

Thanks!

[graeme-winter]

Just wanted to note, working with @mewall I added a commit to the lunus_stills_process branch that fixes this issue in kind of a kludgy way. Comments invited on commit c570d44. Specifically, @graeme-winter, can you verify that grid scan processing is not borked on the lunus_stills_process branch?

Kinda kludgy is true, peppering with is_stills, but if it fixes stills_process & we can somehow mark the flag as a temporary-until-we-actually-fix dxtbx am not completely hostile. Would like to think on it some though...