Issues with reading in multiple passes

[DHekstra]

With reference to /n/hekstra_lab/projects/laue-dials-tests/hemoglobin_TR/playground/tutorial-hemoglobin.ipynb, the following,

laue.initial_solution imported.expt \
                        indexer.refinement.parameterisation.auto_reduction.action=fix \
                        spotfinder.threshold.dispersion.gain=2 \
                        spotfinder.filter.d_min=1.6 \
                        indexer.indexing.known_symmetry.space_group=5 \
                        indexer.indexing.known_symmetry.unit_cell=93.25,43.98,83.5,90.0,122.03,90.0 

fails with a DialsIndexError("No suitable lattice could be found."). The cause of this appears to be related to the angles of the four passes. Specifically, I did not specify geometry.scan.oscillation during dials.import, According to dials.show imported.expt in that case, the code advances the angle of each image by one degree in each pass (instead of the true step size of 6 degrees), and subsequent laue.initial_solution fails as described above. Presumably the cause of the error is therefore that the indexing routine is being fed garbage frame geometry.

When we do specify geometry.scan.oscillation as 0,6 during dials.import, dials.show imported.expt produces four scans with

Scan:
    number of images:   31
    image range:   {1,31}
    oscillation:   {0,6}
    exposure time: 0

and laue.initial_solution seems to succeed (very confusing given my comment below) but leads to trouble when invoking sequence_to_stills.

The difficulty is that the first ON/OFF pass uses angles (0:6:180) and the second pass uses (183:6:363).

Part of the problem seems to be that there are no angles/no angles are read from the MarCCD header, so we need to specify them externally. Naively it seems that there are three possible solutions:

• somehow putting the angles in the headers of the MarCCD files (currently our MarCCD package doesn't support that but it could be extended)
• some other way of getting the DIALS import done, e.g. by doing two separate imports and combining them. I describe this in a comment below.
• use some DIALS function I do not know off to manipulate the MarCCD images.

[DHekstra]

When we limit ourselves to an ON and an OFF pass at the same angles, and specify geometry.scan.oscillation=0,6 during dials.import, laue.initial_solution indexes about 16% of reflections for both datasets, with >3,000 reflections each, but then fails with

Traceback (most recent call last):
  File "/n/hekstra_lab/people/dhekstra/conda_envs/laue-dials/bin/laue.initial_solution", line 8, in <module>
    sys.exit(run())
             ^^^^^
  File "/n/hekstra_lab/people/dhekstra/conda_envs/laue-dials/lib/python3.11/contextlib.py", line 81, in inner
    return func(*args, **kwds)
           ^^^^^^^^^^^^^^^^^^^
  File "/n/home12/dhekstra/ipython_notebooks/laue-dials/src/laue_dials/command_line/initial_solution.py", line 292, in run
    refined_expts, refined_refls = scan_varying_refine(
                                   ^^^^^^^^^^^^^^^^^^^^
  File "/n/home12/dhekstra/ipython_notebooks/laue-dials/src/laue_dials/algorithms/monochromatic.py", line 66, in scan_varying_refine
    expts_refined, refls_refined, _, _ = run_dials_refine(expts, refls, params)
                                         ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/n/hekstra_lab/people/dhekstra/conda_envs/laue-dials/lib/python3.11/site-packages/dials/command_line/refine.py", line 433, in run_dials_refine
    params.input.reflections = None
    ^^^^^^^^^^^^
AttributeError: Please report this error to [email protected]: 'scope_extract' object has no attribute 'input'

followed by a CalledProcessError

This happens when laue.initial_solution purports to be "Duplicating crystal model for scan-varying refinement of experiment 1".

[DHekstra]

The following work-around by suggested by @kmdalton in the DIALS Slack channel seems to circumvent any problems with laue.initial_solution but still fails when it comes to laue.sequence_to_stills.

The workaround is:

dials.import \
    geometry.beam.wavelength=1.04 \
    geometry.scan.oscillation=0.,6. \
    geometry.goniometer.axis=-1,0,0 \
    geometry.beam.polarization_fraction=0.99 \
    geometry.detector.panel.pixel=0.079,0.079 \
    input.template="/n/hekstra_lab/projects/laue-dials-tests/hemoglobin_TR/data/new_imgs/HbI5a_off_###.mccd" \
    input.template="/n/hekstra_lab/projects/laue-dials-tests/hemoglobin_TR/data/new_imgs/HbI5a_100ps_###.mccd" \
    output.experiments=scan_a.expt 

dials.import \
    geometry.beam.wavelength=1.04 \
    geometry.scan.oscillation=183.,6. \
    geometry.goniometer.axis=-1,0,0 \
    geometry.beam.polarization_fraction=0.99 \
    geometry.detector.panel.pixel=0.079,0.079 \
    input.template="/n/hekstra_lab/projects/laue-dials-tests/hemoglobin_TR/data/new_imgs/HbI5b_off_###.mccd" \
    input.template="/n/hekstra_lab/projects/laue-dials-tests/hemoglobin_TR/data/new_imgs/HbI5b_100ps_###.mccd" \
    output.experiments=scan_b.expt 

dials.find_spots scan_a.expt output.reflections=strong_scan_a.refl
dials.find_spots scan_b.expt output.reflections=strong_scan_b.refl

dials.combine_experiments \
    scan_a.expt \
    scan_b.expt \
    strong_scan_a.refl \
    strong_scan_b.refl \
    reference_from_experiment.beam=0 \
    reference_from_experiment.detector=0 \
    reference_from_experiment.goniometer=0  

but upon laue.sequence_to_stills monochromatic.*, the error is:

Traceback (most recent call last):
  File "/n/hekstra_lab/people/dhekstra/conda_envs/laue-dials/bin/laue.sequence_to_stills", line 8, in <module>
    sys.exit(run())
             ^^^^^
  File "/n/hekstra_lab/people/dhekstra/conda_envs/laue-dials/lib/python3.11/contextlib.py", line 81, in inner
    return func(*args, **kwds)
           ^^^^^^^^^^^^^^^^^^^
  File "/n/home12/dhekstra/ipython_notebooks/laue-dials/src/laue_dials/command_line/sequence_to_stills.py", line 229, in run
    total_reflections.extend(new_reflections[i])
                             ~~~~~~~~~~~~~~~^^^
IndexError: Please report this error to [email protected]: list index out of range

followed by a CalledProcessError.

[hkwang]

on laue-dials version 0.3.dev48+gacda0c2:

I run the following script, living in /n/holyscratch01/hekstra_lab/hwang/laue/e35_PDZ/ld_gain_sweep/multipass.sh. I attempt to use dials.combine_experiments to combine multiple passes of e35 PDZ2. The script fails at laue.sequence_to_stills.

################README###########################
# An example script for processing one timepoint and all passes of e35 PDZ2. This can be called as:
# sh multipass.sh <timepoint, e.g. 'off'> <gain> 
# the image files live at /n/holyscratch01/hekstra_lab/hwang/laue/e35_PDZ/images/e35${pass}_${TIME}_###.mccd"
# currently, the script breaks at `laue.sequence_to_stills`. Having this fixed would greatly help 
# processing of multipass datasets. 
#################################################

TIME=${1}
CELL='"65.2,39.45,38.9,90.000,117.45,90.000"'
GAIN=${2}
GAIN_STR=$(echo "$GAIN" | tr . ,)
N=4 # Max multiprocessing

# Make needed directories

mkdir gain_${GAIN_STR}
cd gain_${GAIN_STR}
mkdir dials_files_${TIME}
cd dials_files_${TIME}
eval "$(conda shell.bash hook)"
conda activate laue-dials
rm indexed.expt indexed.refl laue.initial_solution.log

#these start-angles and sweep angles are to be manually input by the user using information from their particular experimental design. 
declare -A START_ANGLES=( ["c"]=0 ["d"]=92 ["e"]=181 ["f"]=361.5)
declare -A OSCS=( ["c"]=2 ["d"]=2 ["e"]=2 ["f"]=1)

for pass in c d e f
do
    
    FILE_INPUT_TEMPLATE="/n/holyscratch01/hekstra_lab/hwang/laue/e35_PDZ/images/e35${pass}_${TIME}_###.mccd"
    # Import data into DIALS files
    dials.import geometry.scan.oscillation=${START_ANGLES[$pass]},${OSCS[$pass]}\
        geometry.goniometer.invert_rotation_axis=True \
        geometry.goniometer.axes=0,1,0 \
        geometry.beam.wavelength=1.04 \
        geometry.detector.panel.pixel_size=0.08854,0.08854 \
        input.template=$FILE_INPUT_TEMPLATE \
        output.experiments=imported_${pass}.expt 
    dials.find_spots imported_${pass}.expt output.reflections=strong_${pass}.refl 
done

dials.combine_experiments imported_*.expt strong_*.refl \
reference_from_experiment.beam=0 reference_from_experiment.detector=0 

# Get a monochromatic geometry model

laue.initial_solution combined.expt \
    indexer.refinement.parameterisation.auto_reduction.action=fix \
    spotfinder.spotfinder.threshold.dispersion.gain=$GAIN \
    indexer.indexing.known_symmetry.space_group=5 \
    indexer.indexing.known_symmetry.unit_cell=$CELL

# Split sequence into stills
laue.sequence_to_stills monochromatic.* 

The error message:

The following parameters have been modified:

input {
  experiments = monochromatic.expt
  reflections = monochromatic.refl
}

Building experiment lists.
Traceback (most recent call last):
  File "/n/home01/hwang6/mambaforge/envs/laue-dials_old/bin/laue.sequence_to_stills", line 8, in <module>
    sys.exit(run())
             ^^^^^
  File "/n/home01/hwang6/mambaforge/envs/laue-dials_old/lib/python3.11/contextlib.py", line 81, in inner
    return func(*args, **kwds)
           ^^^^^^^^^^^^^^^^^^^
  File "/n/home01/hwang6/mambaforge/envs/laue-dials_old/lib/python3.11/site-packages/laue_dials/command_line/sequence_to_stills.py", line 219, in run
    (new_experiments, new_reflections) = sequence_to_stills(
                                         ^^^^^^^^^^^^^^^^^^^
  File "/n/home01/hwang6/mambaforge/envs/laue-dials_old/lib/python3.11/site-packages/laue_dials/command_line/sequence_to_stills.py", line 134, in sequence_to_stills
    crystal=crystals[i_scan_point],
            ~~~~~~~~^^^^^^^^^^^^^^
IndexError: Please report this error to [email protected]: list index out of range

This appears to be the same error as above.

[PrinceWalnut]

The laue.sequence_to_stills errors should be fixed as of this commit. Can I get verification from @DHekstra and @hkwang that your datasets are analyzed properly using the most recent laue-dials version?

[DHekstra]

Could you provide a brief use example?

[DHekstra]

Attached my script for the hemoglobin example. This works in principle. The first ~12 frames of run a have poor geometry, but that is probably unrelated to what is described in this issue. Let's get feedback from @hkwang and then close this issue.

process.sh.txt

[hkwang]

Yes, this works in principle: on either the hemoglobin or e35_PDZ example, the sequence_to_stills.py fix finishes splitting the sequence into stills, and completes downstream processing, instead of throwing the previous bug. However, I still haven't yet scaled and merged a combined dataset due to

1. errors with indexing when the crystal slips in the e35_PDZ case
2. errors with the image headers in the hemoglobin case
   so I don't know that the sequence_to_stills.py fix generates time-resolved signal.

[hkwang]

After using the sequence_to_stills.py fix on the e35_PDZ dataset, where each pass is processed individually, the time-resolved signal looks different than previously, but not noticeably worse.