Thursday, 23 February 2017

Making masks

Yesterday I was asked how to mask a backstop region in DIALS processing, did not have a ready answer so (with the help of "dials.${program} -c -e2" and later pointers from colleagues) figured this out. For my own future reference here are my notes:

First load up the image, twiddle with the "brightness" and "color scheme" settings to allow you to easily see the region you want to mask:




Then press shift and left mouse button to draw out the region you want masking (only rectangles through this mechanism though the dials.generate_mask (command-line) tool can also handle circles. Should look like this:



Next select "save mask" which makes mask.pickle and exit. Reload the image viewer as before but this time also with mask=mask.pickle, select "show mask" and you should see something kin to:



Now, to use this start the usual dials.import with mask=mask.pickle added then everything else (find spots, index, refine, integrate etc.) as usual - when you look at the integration region you will see that there are no spots used on this area of the detector:



This is particularly useful when the rotation axis is not coincident with the backstop...

Friday, 14 October 2016

Scaling data processed with DIALS in SADABS

The chemical crystallography community like "sadabs" for scaling, so we have added the ability to export data in a format which sadabs can consume (strictly: this is an ersatz EvalCCD format.) Essentially following integration with DIALS (e.g. with xia2 or otherwise) run:

dials.export 1_integrated_experiments.json 1_integrated.pickle \
    sadabs.hklout=int1.sad run=1 format=sadabs
...
dials.export 2_integrated_experiments.json 2_integrated.pickle \
    sadabs.hklout=int2.sad run=2 format=sadabs

where you have runs 1...N in your data set (typical for small molecule crystallography work). You can then run sadabs as you usually would, and should see something like:

 Reading file int1.sad
 Reading file int2.sad
 Reading file int3.sad
 Reading file int4.sad

 Mean and maximum errors in direction cosine check function =   0.000   0.000
 The mean error should not exceed 0.005, and is usually caused by matrix
 changes during data processing.

 Approximate wavelength, cell and maximum 2-theta (from cosines etc.):
 0.69010    5.432    8.147   12.053   89.994   90.026   90.002   72.46

in the output, showing in this case that the answers were about as expected. 

Your mileage may vary, this has just been implemented, you will also want the output from:

dials.two_theta_refine p4p=sad.p4p 1_integrated_experiments.json 1_integrated.pickle ...

with the pickle files and json files above to get the unit cell and error for subsequent analysis. Took this data through sadabs, xprep, shelxt, shelxl and the answers seemed sensible...

Thursday, 25 February 2016

Processing Eiger data

Starting with a standard Eiger example data set from DECTRIS FTP site:

insu_with_bs_master.h5
insu_with_bs_data_000001.h5
...
insu_with_bs_data_000037.h5

Broadly we can ignore the 000 files and just run from the master file - everything else will come out by magic (or by HDF5 references). First import the data into DIALS - this just reads the master file and tries to make sense of it:

dials.import insu_with_bs_master.h5

Finding spots etc. works exactly as normal however sometimes there are some unmasked "test" pixels which should be taken into account via the hot pixel mask:

dials.find_spots datablock.json nproc=20 write_hot_mask=true spotfinder.filter.min_spot_size=3

These are the same as the usual commands:

dials.index datablock.json strong.pickle space_group=I23
dials.refine indexed.pickle experiments.json use_all_reflections=true \
scan_varying=true

But for integration pass in the mask:

dials.integrate refined_experiments.json refined.pickle nproc=20 \
mask=hot_mask_0.pickle

Now everything else is as per normal:

dials.export integrated.pickle integrated_experiments.json \
mtz.hklout=integrated.mtz

pointless hklin integrated.mtz hklout sorted.mtz 

aimless hklin sorted.mtz hklout scaled.mtz << eof
anomalous on
resolution 1.4
eof


Starting from scaled.mtz can phase by S-SAD with SHELX C/D/E which suggests this is about right, feedback welcome - dials-support@lists.sourceforge.net 

Friday, 13 November 2015

Drawing predictions

Sometimes it would be useful to see what DIALS thinks are the reflections on an image - this is easily achieved by indexing the data then asking to have the reflection profiles generated, then predicting, then plotting the predictions as follows - from indexing get an experiments file and indexed reflections, then:

Graemes-MacBook-Pro-2:index graeme$ dials.create_profile_model 15_experiments.json 15_indexed.pickle 
The following parameters have been modified:

input {
  experiments = 15_experiments.json
  reflections = 15_indexed.pickle
}

Removed invalid coordinates, 2513 remaining................................0.00s
Sigma B: 0.030353
Sigma M: 0.048908
Wrote experiments to experiments_with_profile_model.json...................0.01s

followed by:

Graemes-MacBook-Pro-2:index graeme$ dials.predict experiments_with_profile_model.json 
The following parameters have been modified:

input {
  experiments = experiments_with_profile_model.json
}

Saved 7392 reflections to predicted.pickle.................................0.01s

and finally:

dials.image_viewer 14_datablock.json predicted.pickle

Job done - you can now see the reflections! This will allow you to look to see that the reflection model matches the observed spots.

Thursday, 8 October 2015

Using Refined Geometry

From time to time you may be faced with a situation where the geometric description of the experiment was not quite right, and you know you have data where you can refine it but want to apply this refinement to another data set. With DIALS, you can. Here's how.

First following the usual recipes get to the point of having the "reference" data processed i.e. some refined_experiments.json file available. You can then fold this information back in with dials.import:


dials.import datablock=../demo/datablock.json \
reference_geometry=../demo/experiments_3.json 

that way when you use the geometry in index it's using the refined geometry. It also shows the neat import datablock from datablock feature - you could just as well have passed a path to the images here. If you have got this refined geometry from a previous look at this data then you can just recycle the strong spot pickle file and rerun indexing, as I did in previous post...

By way of an aside, if you are off to a beamline to record some hard to process data it is well worth taking some nice crystals along i.e. insulin, thaumatin, and recording a calibration set at the same distance and wavelength you will use for the hard experiments. BTW the xia2 equivalent is:

reference_geometry=../demo/experiments_3.json

So you don't even need to do the processing by hand...

Checking Indexing Symmetry

From time to time it may happen you record data with the wrong beam centre in the header (say) and then, if you have a large unit cell, you can have a lot of fun trying to get the data to process right... why? Because it may index perfectly well but with the wrong origin in reciprocal space.

Clearly you can resolve this by looking at the data - however it's much more fun to use DIALS!

OK, so we have run the usual dials.find_spots, dials.index shuffle & ended up with something which looks ... sensible:

Final refined crystal models:
model 1 (93307 reflections):
Crystal:
   Unit cell: (57.593, 57.688, 148.089, 90.579, 90.455, 89.758)
   Space group: P 1
   U matrix:  {{ 0.3287,  0.1200, -0.9368},
               { 0.0693, -0.9923, -0.1027},
               {-0.9419, -0.0312, -0.3344}}
   B matrix:  {{ 0.0174,  0.0000,  0.0000},
               {-0.0001,  0.0173,  0.0000},
               { 0.0001,  0.0002,  0.0068}}
   A = UB:    {{ 0.0056,  0.0019, -0.0063},
               { 0.0013, -0.0172, -0.0007},
               {-0.0164, -0.0006, -0.0023}}

Saving refined experiments to experiments.json

Saving refined reflections to indexed.pickle

However we know here that the beam centre is wrong - how to work out the right centre? Look for the inversion symmetry in the intensities of the strong spots as follows:

dials.check_indexing_symmetry experiments.json indexed.pickle \
grid_search_scope=3 symop_threshold=0.7

This will check the symmetry operations implied by the input experiment (in this case P1; nothing to see here) and since we have set the grid search scope it will also compute the CC on the -h,-k,-l operator taking in to account shifts in h,k,l of -3 to +3 => 343 permutations; it will then print (in this case) those which give a CC > 0.7 and also the 0,0,0 for reference - here we find that the data are indeed misindexed:

symop_threshold = 0.7
grid_search_scope = 3
input {
 experiments = experiments.json
 reflections = indexed.pickle
}

Checking HKL origin:

dH dK dL   Nref    CC
0  0  0  26890 0.393
0  0  2  27849 0.965

So you can then apply this shift to the indexed reflections and rerun the refinement - this should sort out your geometry nicely. You do this with:

dials.reindex indexed.pickle hkl_offset=0,0,2 output.reflections=reindexed.pickle

If this happened to me, I would take the refined geometry and rerun the indexing (that's for another post). In the case I was working on here this improved the #reflections indexed from 93307 to 269577 - worthwhile!

Happy twiddling!

What is this blog?

This is an unofficial blog showing some of the fun things which can be done with the DIALS software (from http://dials.diamond.ac.uk or http://dials.lbl.gov). This does not form part of the official documentation, and readers should be warned that sometimes things will go "out of date" as the software evolves. However it's also fun to see what can be done, right? I'll post stuff here when I do interesting things which I think others may find useful.