ImageD11 package

Submodules

ImageD11.ImageD11_file_series module

ImageD11.ImageD11_file_series.get_options(parser)
ImageD11.ImageD11_file_series.get_series_from_hdf(hdf_file, dark=None, flood=None)
ImageD11.ImageD11_file_series.get_series_from_options(options, args)

Returns a file series thing - not a fabio one

This gives back a fabioimage object with dark and flood corrected data

ImageD11.ImageD11_file_series.get_series_from_stemnum(options, args, dark=None, flood=None)

Returns a file series thing - not a fabio one

ImageD11.ImageD11_file_series.getedfheader(filename)

Reads a header from an edf file in 1024 byte chunks. Assumes enclosing { } Returns string enclosed Adds a filename key at the top

ImageD11.ImageD11_file_series.motor_mne(hd)

expands the _mne and _pos header items of edf headers

ImageD11.ImageD11_file_series.series_from_fabioseries(fabioseries, dark, flood, options)

ImageD11.ImageD11_thread module

class ImageD11.ImageD11_thread.ImageD11_thread(myname='ImageD11_thread')

Bases: threading.Thread

Add a stopping mechanism for unhandled exceptions

ImageD11_stop_now()
run()

Method representing the thread’s activity.

You may override this method in a subclass. The standard run() method invokes the callable object passed to the object’s constructor as the target argument, if any, with sequential and keyword arguments taken from the args and kwargs arguments, respectively.

ImageD11.ImageD11options module

class ImageD11.ImageD11options.ColumnFileType(mode='r')

Bases: ImageD11.ImageD11options.FileType

class ImageD11.ImageD11options.FileType(mode='r')

Bases: object

class ImageD11.ImageD11options.GvectorFileType(mode='r')

Bases: ImageD11.ImageD11options.FileType

class ImageD11.ImageD11options.HdfFileType(mode='r')

Bases: ImageD11.ImageD11options.FileType

class ImageD11.ImageD11options.ImageFileType(mode='r')

Bases: ImageD11.ImageD11options.FileType

class ImageD11.ImageD11options.ParameterFileType(mode='r')

Bases: ImageD11.ImageD11options.FileType

class ImageD11.ImageD11options.SplineFileType(mode='r')

Bases: ImageD11.ImageD11options.FileType

class ImageD11.ImageD11options.UbiFileType(mode='r')

Bases: ImageD11.ImageD11options.FileType

ImageD11.blobcorrector module

class ImageD11.blobcorrector.correctorclass(argsplinefile, orientation='edf')

Bases: object

Applies a spatial distortion to a peak position using a fit2d splinefile

correct(xin, yin)

Transform x,y in raw image coordinates into x,y of an idealised image. Returns a tuple (x,y), expects a pair of floats as arguments

distort(xin, yin)

Distort a pair of points xnew, ynew to find where they would be in a raw image

Iterative algorithm…

make_pixel_lut(dims)

Generate an x and y image which maps the array indices into floating point array indices (to be corrected for pixel size later)

returns FIXME - check they are the right way around

add some sort of known splinefile testcase
make_pos_lut(dims)

Generate a look up table of pixel positions in microns # Cache the value in case of multiple calls returns …

readfit2dspline(name)

Reads a fit2d spline file into a scipy/fitpack tuple, tck A fairly long and dull routine…

test(xin, yin)

Checks that the correct and distort functions are indeed inversely related to each other

class ImageD11.blobcorrector.eiger_spatial(dxfile='/data/id11/nanoscope/Eiger/spatial_20210415_JW/e2dx.edf', dyfile='/data/id11/nanoscope/Eiger/spatial_20210415_JW/e2dy.edf')

Bases: object

pixel_lut()

returns (slow, fast) pixel postions of an image

class ImageD11.blobcorrector.perfect

Bases: ImageD11.blobcorrector.correctorclass

To use on previously corrected when there is no splinefile Allows pixel size etc to be set

correct(xin, yin)

Do nothing - just return the same values

make_pixel_lut(dims)

Generate an x and y image which maps the array indices into floating point array indices (to be corrected for pixel size later)

returns FIXME - check they are the right way around

add some sort of known splinefile testcase
splinefile = 'NO_CORRECTION_APPLIED'
xsize = 'UNKNOWN'
ysize = 'UNKNOWN'
ImageD11.blobcorrector.readfit2dfloats(filep, nfl)

Interprets a 5E14.7 formatted fortran line filep = file object (has readline method) nfl = the number of floats to read

ImageD11.cImageD11 module

ImageD11.cImageD11.fill_in_docstrings()
ImageD11.cImageD11.fix_doc(oldstring, to_be_added)

Adds a description of the function to the f2py string

ImageD11.cImageD11.put_incr(*a, **k)

redirects to put_incr64

ImageD11.cImageD11_docstrings module

Autogenerated from make_pyf.py Edit in _cImageD11.pyf please

ImageD11.columnfile module

ImageD11.columnfile.bench()

Compares the timing for reading with columfile versus np.loadtxt

ImageD11.columnfile.clean(str_lst)

trim whitespace from titles

ImageD11.columnfile.colfile2db(colfilename, dbname)

Read the columnfile into a database Ignores parameter metadata (not used yet)

ImageD11.columnfile.colfile_from_dict(c)

convert from a dictonary of numpy arrays

ImageD11.columnfile.colfile_from_hdf(hdffile, name=None, obj=None)

Read a columnfile from a hdf file FIXME TODO - add the parameters somewhere (attributes??)

ImageD11.columnfile.colfile_to_hdf(colfile, hdffile, name=None, compression=None, compression_opts=None)

Copy a columnfile into hdf file FIXME TODO - add the parameters somewhere (attributes??)

ImageD11.columnfile.colfileobj_to_hdf(cf, hdffile, name=None)

Save a columnfile into hdf file format FIXME TODO - add the parameters somewhere (attributes??)

class ImageD11.columnfile.columnfile(filename=None, new=False)

Bases: object

Class to represent an ascii file containing multiple named columns

addcolumn(col, name)

Add a new column col to the object with name “name” Overwrites in the case that the column already exists

bigarray
chkarray()

Ensure self.bigarray holds our attributes

copy()

Returns a (deep) copy of the columnfile

copyrows(rows)

Returns a copy of select rows of the columnfile

filter(mask)

mask is an nrows long array of true/false

get_bigarray()
getcolumn(name)

Gets data, if column exists

keys()
readfile(filename)

Reads in an ascii columned file

removerows(column_name, values, tol=0)

removes rows where self.column_name == values values is a list of values to remove column name should be in self.titles tol is for floating point (fuzzy) comparisons versus integer

reorder(indices)

Put array into the order given by indices … normally indices would come from np.argsort of something

set_attributes()

Set object vars to point into the big array

set_bigarray(ar)
setcolumn(col, name)

Add a new column col to the object with name “name” Overwrites in the case that the column already exists

setparameters(pars)

update the parameters

sortby(name)

Sort arrays according to column named “name”

updateGeometry(pars=None)
changing or not the parameters it (re)-computes:
xl,yl,zl = ImageD11.transform.compute_xyz_lab tth, eta = ImageD11.transform.compute_tth_eta gx,gy,gz = ImageD11.transform.compute_g_vectors
writefile(filename)

write an ascii columned file

ImageD11.columnfile.fillcols(lines, cols)
class ImageD11.columnfile.newcolumnfile(titles)

Bases: ImageD11.columnfile.columnfile

Just like a columnfile, but for creating new files

ImageD11.compute_fazit module

ImageD11.compute_fazit.get_options(parser)
ImageD11.compute_fazit.main()

A CLI interface

class ImageD11.compute_fazit.xydisp(splinefile=None, parfile=None)

Bases: object

compute_rad_arc()

This part needs more work - how to properly define the output after spd is patched For now we aim to make a triangle filling the same array size

# FIXME - this remains unclear to me, what are the output units supposed to be??

compute_tth_eta(dims)

Find the twotheta and azimuth images

required_pars = ['wavelength', 'distance']
write(stemname)

save the dx, dy images

ImageD11.correct module

ImageD11.correct.correct(data_object, dark=None, flood=None, do_median=False, monitorval=None, monitorcol=None, filterlist=[])

Does the dark and flood corrections Also PIL filters

ImageD11.eps_sig_solver module

ImageD11.eps_sig_solver.readubis(ubifile)

read ubifile and return a list of ubi arrays

class ImageD11.eps_sig_solver.solver(unitcell=None, ubis=None, crystal_symmetry=None, c11=None, c12=None, c13=None, c14=None, c15=None, c16=None, c22=None, c23=None, c24=None, c25=None, c26=None, c33=None, c34=None, c35=None, c36=None, c44=None, c45=None, c46=None, c55=None, c56=None, c66=None)

Bases: object

A class for getting strain and stress tensors

MVStiffness()
compute_write_eps_sig(outputfile)

Compute strain and stress in crystal and sample co-ordinates system

loadmap(filename)
loadpars(filename=None)
savepars(filename=None)
setunitcell(uc)

this is used to set all the unit cell elements in one shot

unitcell()
updateparameters()
ImageD11.eps_sig_solver.write_ubi_file(filename, ubilist)

save 3x3 matrices into file

ImageD11.fft_index_refac module

ImageD11.fft_index_refac.get_options(parser)
class ImageD11.fft_index_refac.grid(npx=128, mr=1.0, nsig=5)

Bases: object

fft()

Compute the Patterson

gv_to_grid_new(gv)

Put gvectors into our grid gv - ImageD11.indexing gvectors

peaksearch(peaksfile)

Peaksearch in the Patterson

props()

Print some properties of the Patterson

pv(v)

print vector

read_peaks(peaksfile)

Read in the peaks from a peaksearch

reduce(vecs)
slow_score()
ImageD11.fft_index_refac.refine_vector(v, gv, tol=0.25, ncycles=25, precision=1e-06)

Refine a (single) real space lattice vector Input

v : real space lattice vector gv : scattering vectors (no errors) tol : tolerance for initial assignment to peak ncycles : how much to refine precision : when to stop refining in v += s when |s|/|v| < precision
Returns
vref : refined real space lattice vector
Uses:
hr = v[0]*gv[i,0] + v[0]*gv[i,1] + v[1]*gv[i,2] 1 per peak hi = round(hr) = calc 1 per peak diff = calc - hr 1 per peak gradient = d(diff)/dv[i] = gv[j,i] 3 per peak matrix = gradient[i].gradient[j] 3x3 peaks^2 rhs = gradient . diff 3

Solve for shifts and iterates reducing tolerance as 1/(ncycles+1) Peaks that change index during refinement are removed Stops on :

no peaks reassigned ncycles runs out |shift| / |vec| < precision
ImageD11.fft_index_refac.test(options)

ImageD11.finite_strain module

class ImageD11.finite_strain.DeformationGradientTensor(ubi, ub0)

Bases: object

SVD

Returns the singular value decomposition of F

U

Returns the Busing and Levy U matrix relating ubi to ub0

VRS

Returns the Polar decomposition of F=V.R=R.S with R as a rotation matrix and V, S as symmetric

finite_strain_lab(m=0.5)

Returns the finite strain in the lab co-ordinate system if the polar decomposition gives : F = V.R = R.S = ui.bi.b0.u0 F.Ft removes u0 initial orientation effect Em = (V^m - I )/2m

finite_strain_ref(m=0.5)

Returns the finite strain in the reference co-ordinate system if the polar decomposition gives : F = V.R = R.S = ui.bi.b0.u0 Ft.F removes ui final orientation effect Em = (S^m - I )/2m

ImageD11.finite_strain.cell_to_B(cell)
ImageD11.finite_strain.e6_to_symm(e)

Follows the eps convention from xfab eps = [e11, e12, e13, e22, e23, e33]

ImageD11.finite_strain.symm_to_e6(m)

Follows the eps convention from xfab eps = [e11, e12, e13, e22, e23, e33]

ImageD11.grain module

ImageD11.grain.e6_to_symm(e)

Follows the eps convention from xfab eps = [e11, e12, e13, e22, e23, e33]

class ImageD11.grain.grain(ubi, translation=None, **kwds)

Bases: object

B
Rod

A Rodriguez vector. Length proportional to angle, direction is axis

U

The orientation matrix (U) from Busing and Levy

UB

The UB matrix from Busing and Levy columns are the reciprocal space lattice vectors

clear_cache()
eps_grain(dzero_cell, m=0.5)
dzero_cell can be another grain or cell parameters
[a,b,c,alpha,beta,gamma]

m is the exponent for the Seth-Hill finite strain tensors E = 1/2m (U^2m - I) (Negative m reverses the reference and lab systems).

Returns eps as a the 6 independent entries in the symmetric matrix
e11 e12 e13 e22 e23 e33

… in the grain reference system of dzero_cell

eps_grain_matrix(dzero_cell, m=0.5)
dzero_cell can be another grain or cell parameters
[a,b,c,alpha,beta,gamma]

m is the exponent for the Seth-Hill finite strain tensors E = 1/2m (U^2m - I) (Negative m reverses the reference and lab systems).

Returns eps as a symmetric matrix … in the grain reference system of dzero_cell

eps_sample(dzero_cell, m=0.5)
dzero_cell can be another grain or cell parameters:
[a,b,c,alpha,beta,gamma]

m is the exponent for the Seth-Hill finite strain tensors E = 1/2m (V^2m - I) (Negative m reverses the reference and lab systems).

Returns eps as a the 6 independent entries in the symmetric matrix
e11 e12 e13 e22 e23 e33

… in the sample system (z up, x along the beam at omega=0)

eps_sample_matrix(dzero_cell, m=0.5)
dzero_cell can be another grain or cell parameters:
[a,b,c,alpha,beta,gamma]

m is the exponent for the Seth-Hill finite strain tensors E = 1/2m (V^2m - I) (Negative m reverses the reference and lab systems).

Returns eps as a symmetric matrix … in the sample system (z up, x along the beam at omega=0)

mt

Metric tensor

rmt

Reciprocal metric tensor

set_ubi(ubi)

Update the orientation and clear cached values

u
ub
unitcell

a,b,c,alpha,beta,gamma

ImageD11.grain.read_grain_file(filename)

read ubifile and return a list of ubi arrays

ImageD11.grain.symm_to_e6(m)

Follows the eps convention from xfab eps = [e11, e12, e13, e22, e23, e33]

ImageD11.grain.write_grain_file(filename, list_of_grains)

ImageD11.grid_index_parallel module

ImageD11.grid_index_parallel.doindex(gve, x, y, z, w, gridpars)

Try to index some g-vectors

ImageD11.grid_index_parallel.domap(pars, colfile, grains, gridpars)

mapping function - does what makemap.py does, but in a function

ImageD11.grid_index_parallel.grid_index_parallel(fltfile, parfile, tmp, gridpars, translations)

fltfile containing peaks parfile containing instrument geometry and unit cell tmp - base name for scratch files and results gridpars : dictionary of control parameters (rings to use, etc) translations : list of translation positions to try

Runs a grid index algorithm using pythons multiprocessing module splits workload over processes (blocks of translations to each process) This thread should catch results via a queue

ImageD11.grid_index_parallel.initgrid(fltfile, parfile, tmp, gridpars)

Sets up a grid indexing by preparing the unitcell for indexing and checking the columns we want are in the colfile

ImageD11.grid_index_parallel.test_many_points(args)

Grid index - loop over points Places the results in a multiprocessing Queue

class ImageD11.grid_index_parallel.uniq_grain_list(symmetry, toldist, tolangle, grains=None)

Bases: object

Cope with finding the same grain over and over…

add(grains)
append_uniq(g)
ImageD11.grid_index_parallel.wrap_test_many_points(x)
ImageD11.grid_index_parallel.wrap_test_many_points_init(q)

This passes the q to the function Something from stackoverflow. It happens during the Pool initializer function

ImageD11.guicommand module

class ImageD11.guicommand.guicommand

Bases: object

Keeps a log of all commands issued - separates gui code from algorithmical code

execute(obj, command, *args, **kwds)

Pass in object as string [peakmerger|transformer|indexer] Pass in command as string, getattr(command) will be used Returns the return value of the function….

TODO : change this interface???
eg : works - returns True
you look for self.lastreturned
fails - returns False
you look for self.lasttraceback
getdata(obj, name)

Allows access to “live” data in the objects wrapped

By passing references back you can circumvent the cleanliness of the interface. Please dont.

Returns object.name

gethistory()

Returns the history of commands run by the gui commander

ImageD11.gv_general module

ImageD11.gv_general.angmod(a)
ImageD11.gv_general.axis_from_matrix(m)

Factory method to create a rotation_axis object from a direction cosine matrix http://en.wikipedia.org/wiki/Rotation_representation_%28mathematics%29

ImageD11.gv_general.chimat(c)
ImageD11.gv_general.chiwedge(chi=0.0, wedge=0.0)

in transform: g = omega . chi . wedge .k

ImageD11.gv_general.g_to_k(g, wavelength, axis=array([0., 0., 1.]), pre=None, post=None)
Get the k and angles given the g in
g = pre . rot(axis, angle) . post . k g = omega . chi . wedge . k

…where k will satisfy the Laue equations

The recipe is:
Find the components of g along and perpendicular to the rotation axis
co-ords are a0 = axis,
a1 = axis x g, a2 = a1 x axis = ( axis x g ) x axis

Rotated vector will be a0 + a1 sin(angle) + a2 cos(angle) Laue condition says that [incident beam].[k] = k.sin(theta)

= 2 sin^2(theta) / lambda = sin(t)/d = |g|.sin(t) = |g|*|g|*lambda / 2

Apply any post rotations to the incident beam Apply any pre-rotations to the g vector |g| = [a0 + a1 sin(angle) + a2 cos(angle) ] . [incident beam] => solve for angle

http://www.ping.be/~ping1339/gonio.htm a sin(u) + b cos(u) = c let: tan(u’) = -b/a and: A = a / cos(u’) Finally you’ll find:

A.sin(u-u’) = c A.cos(pi/2 - u - u’) = c
ImageD11.gv_general.k_to_g(k, angles, axis=None, pre=None, post=None)

Computes g = pre . rot(axis, angle) . post . k Typically in ImageD11 post = [wedge][chi] and pre = identity

since g = Omega.Chi.Wedge.k or k = Wedge-1.Chi-1.Omega-1.g that is to say omega is directly under the sample and wedge is at the bottom of the stack
class ImageD11.gv_general.rotation_axis(direction, angle=0.0)

Bases: object

Represents an axis - so far includes:
axis/angle matrix (direction cosines)

potentially can add quaternions and euler etc later?

rotate_vectors(vectors, angles=None)

Given a list of vectors, rotate them to new ones Angle is either self.angle, or 1 per vector, in degrees

http://www.ecse.rpi.edu/Homepages/wrf/Research/Short_Notes/rotation.html p’ = a . p a + (p - a . p a) cos q + a * p sin q

= p cos q + a . p a (1-cos q) + a * p sin q point p, axis a, angle q

http://mathworld.wolfram.com/RotationFormula.html r’ = r cos(t) + n(n.r)(1-cos(t)) + rxn sin(t)

rotate_vectors_inverse(vectors, angles=None)

Same as rotate vectors, but opposing sense

to_matrix()

Returns a 3x3 rotation matrix R = I3cos(t) + ww^T (1-cos(t)) - W sin(t) t = angle W = vector with hat operator = [ 0 -wz wy

wz 0 -wx
-wy wx 0 ]
ImageD11.gv_general.wedgechi(wedge=0.0, chi=0.0)

in transform: g = omega . chi . wedge .k

ImageD11.gv_general.wedgemat(w)

ImageD11.indexer module

ImageD11.indexer.deriv(xk, f, epsilon, *args)
ImageD11.indexer.get_tth(ds, wvln)

Convert 1/d to Bragg angle

class ImageD11.indexer.indexer(transformpars, colfile)

Bases: object

A class for searching for orientation matrices

NC = 0
assign(ubi, translation, tol)
assigntorings()

Assign the g-vectors to hkl rings that are in self.unitcell

choose(gnum, tol)
gof(p, *args)
loadcolfile(colfile)
pairs(hkl1, hkl2, cos_tol=0.02, hkl_tol=0.1)

We only look for reflection pairs matching a single hkl pairing

printringassign()
refine(ubi, translation=(0, 0, 0), inds=None, hkls=None, tol=0.05)

Fit ubi and translation ihkls = array of peak_index, h, k, l tol = hkl tolerance to use in fitting

reset()

when pars change or colfile changes etc

rings_2_use(rings=None, multimin=12)

Filter rings as having low multiplicity for indexing searches

Give rings = [list of rings to use] Or multimin = all rings with low multiplicity

search1d(gvec_id, hkl=None, angstart=-180, angend=180, nang=3600, tol=0.1)

gvec_id = an integer for a row of self.colfile hkl = hkl indices to assign to the peak (None

means guess from self.cf.ring)

This is inspired from Bernier’s fibre texture method (citation: https://github.com/HEXRD/hexrd/blob/master/hexrd/findorientations.py )

setcell()

Sets the unit cell parameters for indexing

tthcalc(hkls=None)

Computes the twotheta for a unit cell given a list of hkls or from the unit cell generating a list of hkls.

updatecolfile()
ImageD11.indexer.unit(a)

Normalise vector

ImageD11.indexing module

ImageD11.indexing.calc_drlv2(UBI, gv)

Get the difference from being integer hkls UBI: ubi matrix gv: list of g-vectors returns drlv2 = (h_calc - h_int)^2

class ImageD11.indexing.clogging

Bases: object

debug(*args)
error(*args)
info(*args)
log(*args)
warning(*args)
class ImageD11.indexing.indexer(unitcell=None, gv=None, cosine_tol=0.002, minpks=10, hkl_tol=0.01, ring_1=1, ring_2=2, ds_tol=0.005, wavelength=-1, uniqueness=0.5, eta_range=0.0, max_grains=100)

Bases: object

A class for searching for orientation matrices

assigntorings()

Assign the g-vectors to hkl rings

coverage()

Compute the expected coverage of reciprocal space use the min/max obs values of xp/yp/omega to work out what was measured in the scan? No lambda or

fight_over_peaks()

Get the best ubis from those proposed Use all peaks (ring assigned or not)

find()

Dig out the potential hits

friedelpairs(filename)

Attempt to identify Freidel pairs

Peaks must be assigned to the same powder ring Peaks will be the closest thing to being 180 degrees apart

getind(UBI, tol=None, drlv2tmp=None, labelstmp=None)

Returns the indices of peaks in self.gv indexed by matrix UBI

histogram_drlv_fit(UBI=None, bins=None)

Generate a histogram of |drlv| for a ubi matrix For use in validation of grains

loadpars(filename=None)
out_of_eta_range(eta)

decide if an eta is going to be kept

readgvfile(filename, quiet=False)
refine(UBI)

Refine an orientation matrix and rescore it.

From Paciorek et al Acta A55 543 (1999) UB = R H-1

where: R = sum_n r_n h_n^t H = sum_n h_n h_n^t r = g-vectors h = hkl indices
reset()

To get a really clean indexer just create a new one (e.g. via __init__) This was added for the gui to help it forget what happened before but keep parameters as they were set

saveindexing(filename, tol=None)

Save orientation matrices

FIXME : refactor this into something to do
grain by grain } peak by peak }
savepars(filename=None)
saveubis(filename)

Save the generated ubi matrices into a text file

score(UBI, tol=None)

Decide which are the best orientation matrices

score_all_pairs(n=None)

Generate all the potential pairs of rings and go score them too

scorethem(fitb4=False)

decide which trials listed in hits to keep

updateparameters()
ImageD11.indexing.indexer_from_colfile(colfile, **kwds)
ImageD11.indexing.mod_360(theta, target)

Find multiple of 360 to add to theta to be closest to target

ImageD11.indexing.myhistogram(data, bins)

The numpy histogram api was changed So here is an api that will not change It is based on that from the old Numeric manual

ImageD11.indexing.readubis(ubifile)

read ubifile and return a list of ubi arrays

ImageD11.indexing.refine(UBI, gv, tol, quiet=True)

Refine an orientation matrix and rescore it.

From Paciorek et al Acta A55 543 (1999)
UB = R H-1
where:
R = sum_n r_n h_n^t H = sum_n h_n h_n^t r = g-vectors h = hkl indices
ImageD11.indexing.ubi_fit_2pks(ubi, g1, g2)

Refine a ubi matrix so it matches the pair of g-vectors supplied (almost) accounts for cell parameters not being quite right

ImageD11.indexing.ubitoB(ubi)

give the B matrix from ubi

ImageD11.indexing.ubitoRod(ubi)

TODO Testcases!!!

ImageD11.indexing.ubitoU(ubi)

convert ubi to Grainspotter style U The convention is B as being triangular, hopefully as Busing and Levy TODO - make some testcases please!!

ImageD11.indexing.ubitocellpars(ubi)

convert ubi matrix to unit cell

ImageD11.indexing.write_ubi_file(filename, ubilist)

save 3x3 matrices into file

ImageD11.labelimage module

ImageD11.labelimage.flip1(x, y)

fast, slow to dety, detz

ImageD11.labelimage.flip2(x, y)

fast, slow to dety, detz

ImageD11.labelimage.flip3(x, y)

fast, slow to dety, detz

ImageD11.labelimage.flip4(x, y)

fast, slow to dety, detz

ImageD11.labelimage.flip5(x, y)

fast, slow to dety, detz

ImageD11.labelimage.flip6(x, y)

fast, slow to dety, detz

ImageD11.labelimage.flip7(x, y)

fast, slow to dety, detz

ImageD11.labelimage.flip8(x, y)

fast, slow to dety, detz

class ImageD11.labelimage.labelimage(shape, fileout=<_io.TextIOWrapper name='<stdout>' mode='w' encoding='UTF-8'>, spatial=<ImageD11.blobcorrector.perfect object>, flipper=<function flip2>, sptfile=<_io.TextIOWrapper name='<stdout>' mode='w' encoding='UTF-8'>)

Bases: object

For labelling spots in diffraction images

finalise()

Write out the last frame

format = ' %.4f %.4f %.4f %.0f %.4f %.4f %.4f %.4f %.4f %.4f %.4f %.4f %.4f %.4f %.4f %.4f %.0f %.0f %.4f %.0f %.0f %.0f %.0f %.4f %.4f %.4f %.4f %d %d %d\n'
labelpeaks(data, threshold)
measurepeaks(data, omega, blim=None)
mergelast()

Merge the last two images searches

output2dpeaks(file_obj)

Write something compatible with the old ImageD11 format which fabian is reading. This is called before mergelast, so we write self.npk/self.res

outputpeaks(peaks)

Peaks are in Numeric arrays nowadays

peaksearch(data, threshold, omega)

# Call the c extensions to do the peaksearch, on entry: # # data = 2D Numeric array (of your data) # threshold = float - pixels above this number are put into objects

titles = '# sc fc omega Number_of_pixels avg_intensity s_raw f_raw sigs sigf covsf sigo covso covfo sum_intensity sum_intensity^2 IMax_int IMax_s IMax_f IMax_o Min_s Max_s Min_f Max_f Min_o Max_o dety detz onfirst onlast spot3d_id\n'

ImageD11.lattice_reduction module

exception ImageD11.lattice_reduction.BadVectors

Bases: Exception

Raised by lattice class when vectors are coplanar or colinear

ImageD11.lattice_reduction.checkvol(v1, v2, v3, min_vec2=1e-18)

Check whether vectors are singular

ImageD11.lattice_reduction.cosangle_vec(ubi, v)

Angle between v in real and reciprocal space eg, is a* parallel to a or not?

ImageD11.lattice_reduction.find_lattice(vecs, min_vec2=1, n_try=None, test_vecs=None, tol=0.1, fraction_indexed=0.9, noisy=False)

vecs - vectors to use to generate the lattice min_vec2 - squared length of min_vec (units as vec) n_try - max number of vecs to test (clips vecs for generating) test_vecs - vectors to index with the unit cell (else use vecs) fraction_indexed - whether or not to return cells

ImageD11.lattice_reduction.fparl(x, y)

fraction of x parallel to y

ImageD11.lattice_reduction.get_options(parser)
ImageD11.lattice_reduction.iter3d(n)

Generate all possible unordered combinations of vectors i,j,k for i,j,k < n

This looping was rewritten thanks to: TAOCP V4 fascicle 3, section 7.2.1.3.

It gives much nicer ordering than previously, as is gradually expands down the list, instead of hammering the start

ImageD11.lattice_reduction.iter3d_old(n)

Generate all possible unordered combinations of vectors i,j,k for i,j,k < n

class ImageD11.lattice_reduction.lattice(v1, v2, v3, direction=None, min_vec2=1e-18)

Bases: object

Represents a 3D crystal lattice built from 3 vectors

flip(v)

See also __init__.__doc__

matrix(direction)
nearest(vecs)

Give back the nearest lattice point indices, in the same direction

remainders(vecs)

Difference between vecs and closest lattice points

score(vecs, tol=0.1, debug=False)

How many peaks have integer error less than tol?

withvec(x, direction='col')

Try to fit x into the lattice Make the remainder together with current vectors Index it as hkl indices whichever vector has the biggest projection is replaced remake the lattice with these 3 vectors

ImageD11.lattice_reduction.mod(x, y)

Returns the part of x with integer multiples of y removed assert that ||mod(x,y)|| <= ||x|| for all y

ImageD11.lattice_reduction.reduce(v1, v2, v3, min_vec2=1e-18)

Try to find a reduced lattice basis of v1, v2, v3

ImageD11.lattice_reduction.rsweep(vl)

One sweep subtracting each from other two This idea comes from Knuth TAOCP sec 3.3.4C

ImageD11.lattice_reduction.search_2folds(ubi)

Inspired by the Yvon Lepage’s method for finding lattice symmetry Check for 2 fold axes by measuring the directions between real and reciprocal vectors with the same indices. In the case of 2 fold axes they should be parallel

ImageD11.lattice_reduction.sortvec_len(vl)

Sorts according to length (d-s-u)

ImageD11.lattice_reduction.sortvec_xyz(vl)

For each choose between x, -x Sorts according to x then y then z components

ImageD11.license module

ImageD11.parameters module

ImageD11.peakmerge module

class ImageD11.peakmerge.peak(line, omega, threshold, num, tolerance)

Bases: object

Represents a peak

combine(other)

Combine this peak with another peak (eg merge!) # # If the thresholds are different, we see the same peaks # twice. We will just pick the one with the lower threshold # # Otherwise we try to average them #

class ImageD11.peakmerge.peakmerger(quiet='No')

Bases: object

The useful class - called by the gui to process peaksearch output into spots

filter()

Does nothing really, just makes a self.finalpeaks array

TODO: pass functions to filter/compress finalpeaks

getheaderinfo(key)

try to find “key” in the headers and return it

getheaderkeys()

things in the headers

harvestpeaks(numlim=None, omlim=None, thresholds=None)

Harvests the peaks from the images within a range of imagenumbers and/or omega eg: it actually reads the numbers now

mergepeaks()

Merges peaks together - if they overlap and are on adjacent frames

More complex than initially planned

readpeaks(filename, startom=0.0, omstep=1.0)

Read in the output from the peaksearching script Filename is the output file optionally startom and omstep fill in omega ONLY if missing in file

savepeaks(filename)

# Write out minimal information # list of xcorr,ycorr,omega, try to keep intensity now

setpixeltolerance(tolerance=2)

TODO ah… if only

class ImageD11.peakmerge.pkimage(name)

Bases: object

Contains the header information from an image Also pointers to positions of lines in peaksearch output file

addtoheader(headerline)

save header info in dict headerline is line starting with # in peaksearch output

otherheaderstuff(headerline)

Hard to parse header info - called by addtoheader

ImageD11.peakmerge.roundfloat(x, tol)

Return the float nearest to x stepsize tol

ImageD11.peaksearcher module

ImageD11.peaksearcher.get_help(usage=True)

return the help string for online help

ImageD11.peaksearcher.get_options(parser)

Add our options to a parser object

ImageD11.peaksearcher.peaksearch(filename, data_object, corrector, thresholds, labims)

filename : The name of the image file for progress info data_object : Fabio object containing data and header corrector : spatial and dark/flood , linearity etc

thresholds : [ float[1], float[2] etc]

labims : label image objects, one for each threshold

ImageD11.peaksearcher.peaksearch_driver(options, args)

To be called with options from command line

class ImageD11.peaksearcher.timer

Bases: object

msg(msg)
tick(msg='')
tock(msg='')

ImageD11.rc_array module

class ImageD11.rc_array.rc_array

Bases: numpy.ndarray

Row/Column array Represent a list of row or column vectors

check()

Ensure we have an rc_array which is well behaved Pattern assert(check(v)) should disappear in optimiser

flip(mat)

Flip row to col or col.row

inv()

Inverse matrix of self

nb_vector_axis()

The axis which has the n on it

norm2()

sum(v*v,axis=? for row or col

nvectors()
other_direction()

The one which is not self.direction

vector_axis()

The axis which has the 3 on it

ImageD11.refinegrains module

ImageD11.refinegrains.compute_lp_factor(colfile, **kwds)

We’ll assume for simplicity that wedge, chi, tilts are all zero Of course we should put that in sometime in the near future

Try using:

Table 6.2.1.1 in Volume C of international tables : (h) Rotation photograph of small crystal, volume V

1. Beam normal to axis rho =

rac{ N^2 e^4 lambda^3 V} { 4 pi m^2 c^4} .

rac{ 1 + cos^2{ heta} } { sin{2 heta } .

rac{ cos{ heta} { sqrt{ cos^2{psi} - sin^2{ heta} } }
p’ |F|^2

rho is Integrated reflection power ratio from a crystal element e, m Electronic charge and mass c Speed of light lambda Wavelength of radiation 2 heta Angle between incident and diffracted beams psi in (h), latitude of reciprocal-lattice point relative to axis of rotation V Volume of crystal, or of irradiated part of powder sample N Number of unit cells per unit volume F Structure factor of hkl reflection p’ Multiplicity factor for single-crystal methods

nb - what is this ? Hopefully would be 1 ?

From this we will take the (1+cos^2(2t))cos(2t)/sqrt(cos^2p-sin^2t)/sin2t

ImageD11.refinegrains.compute_total_intensity(colfile, indices, tth_range, ntrim=2, quiet=False)

Add up the intensity in colfile for peaks given in indices

ImageD11.refinegrains.cubic(cp)

a=b=c, alpha=beta=gamma=90

ImageD11.refinegrains.get_options(parser)
ImageD11.refinegrains.hexagonal(cp)

a=b,c, alpha=beta=90,gamma=120

ImageD11.refinegrains.lf(tth, eta)

INPUT: 2*theta and eta in degrees. Can be floats or numpy arrays.

OUTPUT: Lorentz scaling factor for intensity. Same data type as input.

EXPLANATION: Compute the Lorentz Factor defined for 3DXRD as

L( theta,eta ) = sin( 2*theta )*|sin( eta )|

This is verified by:

  • Kabsch, W. (1988). Evaluation of single-crystal X-ray di
raction data
from a position-sensitive detector
  • Poulsen, H. F. (2004). 3DXRD { a new probe for materials science. Roskilde: Riso National Laboratory
and can be derived with support of
  • Als-Nielsen, J. and McMorrow, D. (2017). Elements of Modern X-ray Physics

MODIFIED: 7 Feb 2019 by Axel Henningsson.

ImageD11.refinegrains.monoclinic_a(cp)
ImageD11.refinegrains.monoclinic_b(cp)
ImageD11.refinegrains.monoclinic_c(cp)
ImageD11.refinegrains.orthorhombic(cp)

a=b, c, 90,90,90

class ImageD11.refinegrains.refinegrains(tolerance=0.01, intensity_tth_range=(6.1, 6.3), latticesymmetry=<function triclinic>, OmFloat=True, OmSlop=0.25)

Bases: object

Class to refine the position and orientation of grains and also the detector parameters

applyargs(args)
assignlabels(quiet=False)

Fill out the appropriate labels for the spots

compute_gv(thisgrain, update_columns=False)

Makes self.gv refer be g-vectors computed for this grain in this scan

estimate_steps(gof, guess, steps)
fit(maxiters=100)

Fit the global parameters

generate_grains()
getgrains()
gof(args)

<drlv> for all of the grains in all of the scans

loadfiltered(filename)

Read in file containing filtered and merged peak positions

loadparameters(filename)
makeuniq(symmetry)

Flip orientations to a particular choice

Be aware that if the unit cell is distorted and you do this, you might have problems…

pars = {'cell__a': 4.1569162, 'cell__b': 4.1569162, 'cell__c': 4.1569162, 'cell_alpha': 90.0, 'cell_beta': 90.0, 'cell_gamma': 90.0, 'cell_lattice_[P,A,B,C,I,F,R]': 'P', 'chi': 0.0, 'distance': 7367.8452302, 'fit_tolerance': 0.5, 'o11': 1, 'o12': 0, 'o21': 0, 'o22': 1, 'omegasign': 1, 't_x': 33.0198146824, 't_y': 14.6384893741, 't_z': 0.0, 'tilt_x': 0.0, 'tilt_y': 0.0623952920101, 'tilt_z': 0.00995011461696, 'wavelength': 0.2646, 'wedge': 0.0, 'y_center': 732.950204632, 'y_size': 4.6, 'z_center': 517.007049626, 'z_size': 4.6}
printresult(arg)
readubis(filename)

Read ubi matrices from a text file

refine(ubi, quiet=True)

Fit the matrix without changing the peak assignments

refinepositions(quiet=True, maxiters=100)
refineubis(quiet=True, scoreonly=False)
reset_labels(scanname)
savegrains(filename, sort_npks=True)

Save the refined grains

saveparameters(filename)
set_translation(gr, sc)
stepsizes = {'chi': 0.0017453292519943296, 'distance': 200.0, 't_x': 0.2, 't_y': 0.2, 't_z': 0.2, 'tilt_x': 0.0017453292519943296, 'tilt_y': 0.0017453292519943296, 'tilt_z': 0.0017453292519943296, 'wavelength': 0.001, 'wedge': 0.0017453292519943296, 'y_center': 0.2, 'y_size': 0.2, 'z_center': 0.2, 'z_size': 0.2}
ImageD11.refinegrains.test_benoit()
ImageD11.refinegrains.test_nac()
ImageD11.refinegrains.tetragonal(cp)

a=b, c, 90,90,90

ImageD11.refinegrains.triclinic(cp)
ImageD11.refinegrains.trigonalH(cp)

a=b,c, alpha=beta=90,gamma=120

ImageD11.refinegrains.trigonalP(cp)

a=b=c, alpha=beta=gamma

ImageD11.rotdex module

ImageD11.rotdex.compute_Cgve(txyz, peaks_Cxyz, beam_Cxyz, wavelength)

Computes g-vectors from spots and beam in crystal frame txyz = (3,) origin position peaks_Cxyz = (3,n) spot position on rotated detector beam_Cxyz = (3,n) vectors of rotated incident beam

beam is already normalised to wavelength

wavelength = float radation, normalises length returns scattering vectors (g-vectors)

ImageD11.rotdex.compute_dgdt(txyz, peaks_Cxyz, beam_Cxyz, wavelength)

Compute the gvectors and the derivatives with respect to tx,ty,tz txyz = (3,) origin position peaks_Cxyz = (3,n) spot position on rotated detector beam_Cxyz = (3,n) vectors of rotated incident beam

beam is already normalised to wavelength

wavelength = float radation, normalises length returns scattering vectors (g-vectors) and d(g)/d(txyz)

ImageD11.rotdex.fit_ub_t(ub, translation, hkl, peaks_Cxyz, beam_Cxyz, wavelength)

Fits the ub and grain origin to a list of assigned peaks All unit cell and orientations parameters are free Runs 2 cycles (empirically this converges)

ub = (3,3) input UB matrix (so that h ~= UB.g) translation = (3,) grain origin hkl = (3,n) integer peak assignments peaks_Cxyz = (3,n) spot positions in crystal frame (getCxyz) beam_Cxyz = (3,n) beam directions in crystal frame (getCxyz) wavelength = float, radiation, normalises length returns fitted UB and translation

ImageD11.rotdex.fitagrain(gr, pars)
ImageD11.rotdex.getCxyz(colf, pars)

Find spot positions in crystal frame colf = columnfile (or other object) holding .sc, .fc, .omega pars = parameters for ImageD11 transform module

ImageD11.rotdex.main()
ImageD11.rotdex.main2()

ImageD11.rsv module

ImageD11.rsv.getbounds(vol, plane)

Returns the extent argument to use for pylab.imshow when plotting a plane

ImageD11.rsv.mem()

debug the memory usage

ImageD11.rsv.readvol(filename, savespace=False)

Read volume from a file returns an rsv object Take care to allocate and read to avoid temporaries

class ImageD11.rsv.rsv(dimensions, bounds, np, **kwds)

Bases: object

A reciprocal space volume

allocate_vol()

Allocates memory for a volume data

normalise(savespace=True)

Return the normalised but avoid divide by zero

plnames = {0: 0, 'h': 0, 'H': 0, 1: 1, 'k': 1, 'K': 1, 2: 2, 'l': 2, 'L': 2}
slice(plane, num)

return signal on plane index num

ImageD11.rsv.writenormedvol(vol, filename)

Write volume in vol to filename - save only the normalised to avoid using so much memory

Compress -1 is for the zeros, which there might be a lot of

ImageD11.rsv.writevol(vol, filename)

Write volume in vol to filename

Compress -1 is for the zeros, which there might be a lot of

ImageD11.rsv_mapper module

ImageD11.rsv_mapper.get_options(parser)

Command line interface for making a mapping Add our options to a parser object

ImageD11.rsv_mapper.main()

A user interface

class ImageD11.rsv_mapper.rsv_mapper(dims, pars, ubi, splinefile=None, np=16, border=10, omegarange=[0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99, 100, 101, 102, 103, 104, 105, 106, 107, 108, 109, 110, 111, 112, 113, 114, 115, 116, 117, 118, 119, 120, 121, 122, 123, 124, 125, 126, 127, 128, 129, 130, 131, 132, 133, 134, 135, 136, 137, 138, 139, 140, 141, 142, 143, 144, 145, 146, 147, 148, 149, 150, 151, 152, 153, 154, 155, 156, 157, 158, 159, 160, 161, 162, 163, 164, 165, 166, 167, 168, 169, 170, 171, 172, 173, 174, 175, 176, 177, 178, 179, 180, 181, 182, 183, 184, 185, 186, 187, 188, 189, 190, 191, 192, 193, 194, 195, 196, 197, 198, 199, 200, 201, 202, 203, 204, 205, 206, 207, 208, 209, 210, 211, 212, 213, 214, 215, 216, 217, 218, 219, 220, 221, 222, 223, 224, 225, 226, 227, 228, 229, 230, 231, 232, 233, 234, 235, 236, 237, 238, 239, 240, 241, 242, 243, 244, 245, 246, 247, 248, 249, 250, 251, 252, 253, 254, 255, 256, 257, 258, 259, 260, 261, 262, 263, 264, 265, 266, 267, 268, 269, 270, 271, 272, 273, 274, 275, 276, 277, 278, 279, 280, 281, 282, 283, 284, 285, 286, 287, 288, 289, 290, 291, 292, 293, 294, 295, 296, 297, 298, 299, 300, 301, 302, 303, 304, 305, 306, 307, 308, 309, 310, 311, 312, 313, 314, 315, 316, 317, 318, 319, 320, 321, 322, 323, 324, 325, 326, 327, 328, 329, 330, 331, 332, 333, 334, 335, 336, 337, 338, 339, 340, 341, 342, 343, 344, 345, 346, 347, 348, 349, 350, 351, 352, 353, 354, 355, 356, 357, 358, 359], maxpix=None, mask=None)

Bases: object

Maps images into a reciprocal space volume

Basic idea is for a fixed detector, scattering vectors in lab frame are fixed. When rotating sample, just rotate these into crystal frame and use as array indices into volume

add_image(om, data)

RSV = bounds of reciprocal space vol NR = dims of RSV k = scattering vector for image at om == 0 data = 2D image (dims == k.dims) SIG = signal MON = monitor

find_vol(border, omegarange)

find limiting volume The four image corners over 360 degrees

np is the number of pixels per hkl

returns (RSV, NR)
RSV min/max in reciprocal space (np*ubi).gv NR number of points == 1+RSV[i][1]-RSV[i][0]
make_k_vecs()

Generate the k vectors from the experiment parameters given in constructor

writevol(filename)

Save the volume in a hdf file

ImageD11.saintraw module

class ImageD11.saintraw.saintraw(filename=None)

Bases: object

condition_filter(name, func)

Remove the peaks according to condition

doc = '\nIHKL 3I4 HKL indices of the reflection \n\n #IMNP 3I4 PRESENT ONLY FOR MODULATED-STRUCTURE DATA\n MNP indices of the reflection. Unused\n indices (in 4- or 5-dimensional cases)\n are written as zero\n\nFI F8.x Integrated intensity in photons, background-\n subtracted, normalized to 1 min/deg, and\n Lp corrected. Number of places to right\n of the decimal point is adjusted to\n balance precision vs. overflow; exponential\n format is used for very large values.\n\nSI F8.x Estimated standard deviation of the\n intensity. Number of places to right\n of the decimal point is adjusted to\n balance precision vs. overflow; exponential\n format is used for very large values.\n\nIBATNO I4 Batch number (for scaling) assigned to\n the reflection\n\nCOSINES 6F8.5 Direction cosines of the incident and\n diffracted beams, relative to the \n unit cell axes. Order is XI,XD,YI,YD,\n ZI,ZD where XI is the X-component of\n the incident beam, XD is the X-component\n of the diffracted beam, etc. Used\n for absorption correction.\n\n\nMSTATUS I3 1X,Z2 Status mask reserved for flagging\n abnormal conditions. There are\n currently none defined.\n\nXO F7.2 Observed X-pixel coordinate of the\n intensity-weighted reflection centroid,\n in reduced pixels (scale of 512x512)\n\nYO F7.2 Observed Y-pixel coordinate of the\n intensity-weighted reflection centroid,\n in reduced pixels (scale of 512x512)\n\nZO F8.2 Observed frame number of the\n intensity-weighted reflection centroid\n\nXP F7.2 Predicted X-pixel of the reflection\n in reduced pixels (scale of 512x512)\n\nYP F7.2 Predicted Y-pixel of the reflection\n in reduced pixels (scale of 512x512)\n\nZP F8.2 Predicted frame number of the reflection\n\nCORLPAF F6.3 Multiplicative correction for Lorentz\n effect, polarization, air absorption,\n and detector faceplate absorption,\n already applied to integrated intensity,\n FI\n\nCORREL F5.2 The correlation coefficient between the\n 3-D profile observed for this reflection\n and the corresponding model 3-D profile\n\nACCUMTIME F7.2 Accumulated hours of exposure\n\nSWING F7.2 Detector swing angle in degrees\n\nROT F7.2 Scan axis (omega or phi) setting in\n degrees at which the reflection was\n observed\n\nIAXIS I2 Scan axis number (2=omega, 3=phi)\n\nISTL I5 SIN(THETA)/LAMBDA times 10,000 \n\nPKSUM I9 Total raw peak counts in photons, with no\n correction for background, normalization, \n or Lp\n\nBGAVG I7 Normally, average BG per pixel in \n photons * 1000. In data sets where the \n background was reported to be very large \n ( > 1024 photons after 1 min/deg \n normalization), the program issues a warning \n message during integration and the scale of\n 1000 is omitted.\n\nROTREL F7.2 Rotation of the scan axis (omega or phi)\n in degrees relative to the start of the\n run\n\nICRYST I4 Crystal number (for scaling)\n\nPKFRAC F6.3 Fraction of the profile volume nearest\n this HKL. 1-PKFRAC is the fraction of\n the intensity which had to be estimated\n from the model profiles because of\n overlap with neighboring spots\n\nIKEY I11 The unique sort key for the group of\n equivalents to which this HKL belongs.\n IKEY = 1048576 * (H+511) + 1024 *\n (K+511) + (L+511), where H,K and L are\n the HKL indices transformed to the base\n asymmetric unit.\n\nIMULT I3 The point-group multiplicity of this HKL\n\nCORLORENTZ F6.3 Lorentz correction\n\nXGEO F8.2 Spatially corrected X relative to beam\n center, in pixels\n\nYGEO F8.2 Spatially corrected Y relative to beam \n center, in pixels\n\nCHI F8.3 Chi setting angle, degrees\n\nOTHER_ANGLE F8.3 The remaining setting angle, degrees. This\n will be phi if scans were in omega, or omega\n if scans were in phi\n\nICOMPONENT I4 Twin component number in SHELXTL HKLF 5\n convention. In a multi-component overlap,\n ICOMPONENT is negated for all but the last\n record of the overlap\n\n'
formats = {}
helps = {}
parsedocs()

Parse the saint documentation for the Bruker format

read(filename)

Read an ascii formatted saint reflection file

sort(name)

Sort according to a column in self.data

take(order)

Put the peaks in the order given in order (indices)

titles = []
tocolumnfile()

Return a columnfile

write(filename)

Write an ascii formatted saint reflection file

ImageD11.scale module

class ImageD11.scale.scale(im1, threshold=None)

Bases: object

scale(im2)

Fill out RHS and solve returns the scale to apply to the image stored in the class

You probably want the scale to apply to the image you supply …use scale image for that

scaleimage(im2)

Return a copy of the image scaled to match the class

ImageD11.scale.scaleseries(target, stem, first, last, thresh=None, writeim=None)

Scale a series of [bruker] images to the target TODO - make it work with fabio file series

ImageD11.simplex module

class ImageD11.simplex.Simplex(testfunc, guess, increments, kR=-1, kE=2, kC=0.5)

Bases: object

accept_contracted_point()
accept_expanded_point()
accept_reflected_point()
calculate_errors_at_vertices()
contract_simplex()
expand_simplex()
minimize(epsilon=0.0001, maxiters=250, monitor=1)

Walks to the simplex down to a local minima. INPUTS —— epsilon convergence requirement maxiters maximum number of iterations monitor if non-zero, progress info is output to stdout

an array containing the final values lowest value of the error function number of iterations taken to get here

multiple_contract_simplex()
reflect_simplex()
ImageD11.simplex.main()
ImageD11.simplex.myfunc(args)

ImageD11.sparseframe module

class ImageD11.sparseframe.SparseScan(hname, scan, start=0, n=None, names=['row', 'col', 'intensity'], omeganames=['measurement/rot_center', 'measurement/rot', 'measurement/diffrz_center', 'measurement/diffrz'], dtynames=['measurement/dty_center', 'measurement/dty', 'measurement/diffty_center', 'measurement/diffty'])

Bases: object

cplabel(threshold=0, countall=True)

Label pixels using the connectedpixels assigment code Fills in:

self.nlabels = number of peaks per frame self.labels = peak labels (should be unique) self.total_labels = total number of peaks
if countall == True : labels all peaks from zero
== False : labels from 1 on each frame
getframe(i, SAFE=True)
lmlabel(threshold=0, countall=True, smooth=True)

Label pixels using the localmax assigment code Fills in:

self.nlabels = number of peaks per frame self.labels = peak labels (should be unique) self.total_labels = total number of peaks
if countall == True : labels all peaks from zero
== False : labels from 1 on each frame
moments()

Computes the center of mass in s/f/omega

ImageD11.sparseframe.from_data_cut(data, cut, header={}, detectormask=None)
ImageD11.sparseframe.from_data_mask(mask, data, header)

Create a sparse from a dense array

ImageD11.sparseframe.from_hdf_group(group)
ImageD11.sparseframe.overlaps(frame1, labels1, frame2, labels2)

figures out which label of self matches which label of other Assumes the zero label does not exist (background) Returns sparse array of: label in self (row) label in other (col) number of shared pixels (data)

class ImageD11.sparseframe.overlaps_linear(nnzmax=16384)

Bases: object

Memory caching object for the linear time algorithm to find peak overlaps

Given (row1, col1, label1) and (row2, col2, label2) it finds pixels where (row[i] == row2[i]) and (col1[i] == col2[i]) and returns (labels1[i], labels2[i], sum_pixels[i]) … so the number of overlapping pixels for that pair of labels

realloc()
class ImageD11.sparseframe.overlaps_matrix(npkmax=256)

Bases: object

Memory caching object for the quadratic time algorithm to find peak overlaps

Given (row1, col1, label1) and (row2, col2, label2) it finds pixels where (row[i] == row2[i]) and (col1[i] == col2[i]) and returns (labels1[i], labels2[i], sum_pixels[i]) … so the number of overlapping pixels for that pair of labels

This is easier to understand and faster for small number of peaks per frame

realloc()
ImageD11.sparseframe.sparse_connected_pixels(frame, label_name='connectedpixels', data_name='intensity', threshold=None)

frame = a sparse frame label_name = the array to save labels to in that frame data_name = an array in that frame threshold = float value or take data.threshold

class ImageD11.sparseframe.sparse_frame(row, col, shape, itype=<class 'numpy.uint16'>, pixels=None, SAFE=True)

Bases: object

Indices / shape mapping

This was developed for a single 2D frame
See SparseScan below for something aiming towards many frames
check(row, col, shape, itype, SAFE=True)

Ensure the index data makes sense and fits

is_sorted()

Tests whether the data are sorted into slow/fast order rows are slow direction columns are fast

mask(msk)

returns a subset of itself

reorder(order)

Put the pixels into a different order (in place)

set_pixels(name, values, meta=None)

Named arrays sharing these labels

sort()

Puts you into slow / fast looping order

sort_by(name)

Not sure when you would do this. For sorting by a peak labelling to get pixels per peak

threshold(threshold, name='intensity')

returns a new sparse frame with pixels > threshold

to_dense(data=None, out=None)

returns the full 2D image data = name in self.pixels or 1D array matching self.nnz Does not handle repeated indices e.g. obj.to_dense( obj.pixels[‘raw_intensity’] )

to_hdf_group(group)

Save a 2D sparse frame to a hdf group Makes 1 single frame per group

ImageD11.sparseframe.sparse_localmax(frame, label_name='localmax', data_name='intensity')
ImageD11.sparseframe.sparse_moments(frame, intensity_name, labels_name)

We rely on a labelling array carrying nlabel metadata (==labels.data.max())

ImageD11.sparseframe.sparse_smooth(frame, data_name='intensity')

ImageD11.sym_u module

ImageD11.sym_u.cubic()

P32

ImageD11.sym_u.find_uniq_hkls(hkls, grp, func=<function hklmax>)
ImageD11.sym_u.find_uniq_u(u, grp, debug=0, func=<function trace>)
ImageD11.sym_u.fmt(c)
ImageD11.sym_u.generate_group(*args)
ImageD11.sym_u.getgroup(s)

convert a user supplied string to a group … a little vague still

class ImageD11.sym_u.group(tol=1e-05)

Bases: object

An abstract mathematical finite(?) point rotation groups

additem(item)

add a new member

comp(x, y)

Compare two things for equality

isMember(item)

Decide if item is already in the group

makegroup()

ensure all items = op(x,y) are in group

op(x, y)

Normally multiplication ? Means of generating new thing from two others

ImageD11.sym_u.hexagonal()

P6/mmm 191

ImageD11.sym_u.hklmax(h, hmax=1000)
ImageD11.sym_u.m_from_string(s)

Creates a symmetry operator from a string

ImageD11.sym_u.m_to_string(m)

Creates a symmetry operator string from a matrix

ImageD11.sym_u.monoclinic_a()
ImageD11.sym_u.monoclinic_b()
ImageD11.sym_u.monoclinic_c()
ImageD11.sym_u.orthorhombic()

P222 16

ImageD11.sym_u.rhombohedralP()

R3 primitive

ImageD11.sym_u.test()
ImageD11.sym_u.tetragonal()

P4 75

class ImageD11.sym_u.trans_group(tol=1e-05)

Bases: ImageD11.sym_u.group

Translation group (eg crystal lattice)

FIXME - this is mostly done in lattice_reduction.py instead now

additem(x)

Do lattice reduction before adding as infinite group

isMember(x)

Decide if item is already in the group

mod(x, y)

Remove y from x to give smallest possible result Find component of x || to y and remove it

op(x, y)

Means of generating new thing from two others In this case add them and mod by group members

reduce(v)

Perform lattice reduction

ImageD11.sym_u.triclinic()
ImageD11.sym_u.trigonal()

P321 150

ImageD11.symops module

ImageD11.symops.checkop(h, k, l, op, axis)
ImageD11.symops.glide_plane(h, k, l, gtype, axis)
ImageD11.symops.glidetest(gtype, axis)
ImageD11.symops.lattice_centre(h, k, l, ctype)
ImageD11.symops.mirror_plane(h, k, l, axis)
ImageD11.symops.rotation_axis(h, k, l, rtype, axis)
ImageD11.symops.screw_axis(h, k, l, stype, axis)
ImageD11.symops.test_absence(h, k, l, sg)

ImageD11.threshold_image module

ImageD11.transform module

class ImageD11.transform.Ctransform(pars)

Bases: object

pnames = ('y_center', 'z_center', 'y_size', 'z_size', 'distance', 'wavelength', 'omegasign', 'tilt_x', 'tilt_y', 'tilt_z', 'o11', 'o12', 'o21', 'o22', 'wedge', 'chi')
reset()
sf2gv(sc, fc, omega, tx=0, ty=0, tz=0, out=None)
sf2xyz(sc, fc, tx=0, ty=0, tz=0, out=None)
xyz2gv(xyz, omega, tx=0, ty=0, tz=0, out=None)
class ImageD11.transform.PixelLUT(pars)

Bases: object

A look up table for a 2D image to store pixel-by-pixel values

pnames = ('y_center', 'z_center', 'y_size', 'z_size', 'distance', 'wavelength', 'omegasign', 'tilt_x', 'tilt_y', 'tilt_z', 'o11', 'o12', 'o21', 'o22', 'wedge', 'chi', 'dxfile', 'dyfile', 'spline', 'shape')
spatial(sraw, fraw)

applies a spatial distortion to sraw, fraw (for peak centroids)

ImageD11.transform.compute_g_from_k(k, omega, wedge=0, chi=0)

Compute g-vectors with cached k-vectors

ImageD11.transform.compute_g_vectors(tth, eta, omega, wvln, wedge=0.0, chi=0.0)
Generates spot positions in reciprocal space from
twotheta, wavelength, omega and eta

Assumes single axis vertical … unless a wedge angle is specified

ImageD11.transform.compute_grain_origins(omega, wedge=0.0, chi=0.0, t_x=0.0, t_y=0.0, t_z=0.0)

# print “Using translations t_x %f t_y %f t_z %f”%(t_x,t_y,t_z) # Compute positions of grains # expecting tx, ty, tz for each diffraction spot # # g = R . W . k # g - is g-vector w.r.t crystal # k is scattering vector in lab # so we want displacement in lab from displacement in sample # shift = W-1 R-1 crystal_translation # # R = ( cos(omega) , sin(omega), 0 ) # (-sin(omega) , cos(omega), 0 ) # ( 0 , 0 , 1 ) # # W = ( cos(wedge) , 0 , sin(wedge) ) # ( 0 , 1 , 0 ) # (-sin(wedge) , 0 , cos(wedge) ) # # C = ( 1 , 0 , 0 ) ??? Use eta0 instead # ( 0 , cos(chi) , sin(chi) ) ??? Use eta0 instead # ( 0 , -sin(chi) , cos(chi) ) ??? Use eta0 instead

ImageD11.transform.compute_k_vectors(tth, eta, wvln)

generate k vectors - scattering vectors in laboratory frame

ImageD11.transform.compute_lorentz_factors(tth, eta, omega, wavelength, wedge=0.0, chi=0.0)

From Kabsch 1988 J. Appl. Cryst. 21 619

Multiply the intensities by: Lorentz = | S.(u x So)| / |S|.|So| S = scattered vector So = incident vector u = unit vector along rotation axis

ImageD11.transform.compute_polarisation_factors(args)

From Kabsch 1988 J. Appl. Cryst. 21 619

DIVIDE the intensities by: <sin2 psi> = (1 - 2p) [ 1 - (n.S/|S|^2) ] + p { 1 + [S.S_0/(|S||S_0|)^2]^2}

p = degree of polarisation (sync = 1, tube = 0.5 , mono + tube in between)
or “probability of finding electric field vector in plane having
normal, n”

S = scattered vector S_0 = incident vector n = normal to polarisation plane, typically perpendicular to S_0

In ImageD11 we normally expect to find: x axis along the beam z axis being up, and parallel to the normal n mentioned above

ImageD11.transform.compute_sinsqth_from_xyz(xyz)

Computes sin(theta)**2 x,y,z = co-ordinates of the pixel in cartesian space

if you need grain translations then use func(peaks_xyz - compute_grain_origins(…) )

seems to be competitive with arctan2 (docs/sintheta_squared_geometry.ipynb)

returns sin(theta)**2

ImageD11.transform.compute_tth_eta(peaks, y_center=0.0, y_size=0.0, tilt_y=0.0, z_center=0.0, z_size=0.0, tilt_z=0.0, tilt_x=0.0, distance=0.0, o11=1.0, o12=0.0, o21=0.0, o22=-1.0, t_x=0.0, t_y=0.0, t_z=0.0, omega=None, wedge=0.0, chi=0.0, **kwds)

Finds x,y,z co-ordinates of peaks in the laboratory frame Computes tth/eta from these (in degrees)

kwds are not used (left for convenience if you have a parameter dict)

ImageD11.transform.compute_tth_eta_from_xyz(peaks_xyz, omega, t_x=0.0, t_y=0.0, t_z=0.0, wedge=0.0, chi=0.0, **kwds)

Peaks is a 3 d array of x,y,z peak co-ordinates crystal_translation is the position of the grain giving rise to a diffraction spot in x,y,z ImageD11 co-ordinates

x,y is with respect to the axis of rotation (usually also beam centre). z with respect to beam height, z centre

omega data are needed if crystal translations are used

computed via the arctan recipe.

returns tth/eta in degrees

ImageD11.transform.compute_tth_histo(tth, no_bins=100, weight=False, weights=None, **kwds)

Compute a histogram of tth values Uses numpy’s histogram rather that doing it by hand as above New feature: weight by something (peak intensity for instance), send true for weight and weights values

Returns a normalised histogram (should make this a probability and

For each datapoint, the corresponding histogram weight

Updated and modernized 2021-02-11 S. Merkel

ImageD11.transform.compute_xyz_from_tth_eta(tth, eta, omega, t_x=0.0, t_y=0.0, t_z=0.0, wedge=0.0, chi=0.0, **kwds)

Given the tth, eta and omega, compute the xyz on the detector

crystal_translation is the position of the grain giving rise to a diffraction spot

in x,y,z ImageD11 co-ordinates x,y with respect to axis of rotation and or beam centre ?? z with respect to beam height, z centre

omega data needed if crystal translations used

ImageD11.transform.compute_xyz_lab(peaks, y_center=0.0, y_size=0.0, tilt_y=0.0, z_center=0.0, z_size=0.0, tilt_z=0.0, tilt_x=0.0, distance=0.0, o11=1.0, o12=0.0, o21=0.0, o22=-1.0, **kwds)

Peaks is a 2 d array of x,y yc is the centre in y ys is the y pixel size ty is the tilt around y zc is the centre in z zs is the z pixel size tz is the tilt around z dist is the sample - detector distance detector_orientation is a matrix to apply to peaks arg to get ImageD11 convention

(( 0, 1),( 1, 0)) for ( y, x) ((-1, 0),( 0, 1)) for (-x, y) (( 0,-1),(-1, 0)) for (-y,-x)

etc…

kwds are not used (but lets you pass in a dict with other things in it)

ImageD11.transform.cross_product_2x2(a, b)

returns axb for two len(3) vectors a,b

ImageD11.transform.detector_rotation_matrix(tilt_x, tilt_y, tilt_z)

Return the tilt matrix to apply to peaks tilts are in radians typically applied to peaks rotating around beam center

ImageD11.transform.sinth2_sqrt_deriv(xyz)

sin(theta)**2 from xyz, and derivatives w.r.t x,y,z

ImageD11.transform.uncompute_g_vectors(g, wavelength, wedge=0.0, chi=0.0)

Given g-vectors compute tth,eta,omega assert uncompute_g_vectors(compute_g_vector(tth,eta,omega))==tth,eta,omega

ImageD11.transform.uncompute_one_g_vector(gv, wavelength, wedge=0.0)

Given g-vectors compute tth,eta,omega assert uncompute_g_vectors(compute_g_vector(tth,eta,omega))==tth,eta,omega

ImageD11.transformer module

class ImageD11.transformer.transformer(parfile=None, fltfile=None)

Bases: object

Handles the algorithmic, fitting and state information for fitting parameters to give experimental calibrations

addcellpeaks(limit=None)

Adds unit cell predicted peaks for fitting against Optional arg limit gives highest angle to use

Depends on parameters:
‘cell__a’,’cell__b’,’cell__c’, ‘wavelength’ # in angstrom ‘cell_alpha’,’cell_beta’,’cell_gamma’ # in degrees ‘cell_lattice_[P,A,B,C,I,F,R]’
addcolumn(col, name)

Return the data

applyargs(args)

for use with simplex/gof function, alter parameters

compute_histo(colname)

Compute the histogram over twotheta for peaks previous read in Filtering is moved to a separate function

colname is most-often “tth”

other parameters are set in the parameter object no_bins = number of bins weight_hist_intensities: True or False

False: histogram by number of measured peaks True: weight by peak intensities
compute_tth_eta()

Compute the twotheta and eta for peaks previous read in

compute_tth_histo()

Give hardwire access to tth

computegv()

Using tth, eta and omega angles, compute x,y,z of spot in reciprocal space

filter_min(col, minval)

and filter peaks out falling in bins with less than min_prob

fit(tthmin=0, tthmax=180)

Apply simplex to improve fit of obs/calc tth

get_variable_list()
getaxis()

Compute the rotation axis This handles omegasign in a more elegant way

getcols()
getcolumn(name)

Return the data

getvars()

decide what is refinable

gof(args)

Compute how good is the fit of obs/calc peak positions in tth

loadfileparameters(filename)

Read in beam center etc from file

loadfiltered(filename)

Read in 3D peaks from peaksearch

save_tth_his(filename, bins, hist)

Saves a 2 theta histogram into a file bins and histogram should have been calculated previously

savegv(filename)

Save g-vectors into a file Use crappy .ascii format from previous for now (testing)

saveparameters(filename)

Save beam center etc to file

setvars(varlist)

set the things to refine

setxyomcols(xname, yname, omeganame)
tth_entropy()

Compute the entropy of the two theta values … this may depend on the number of tth bins (perhaps?) … to be used for cell parameter free line straightening S = sum -p log p

updateparameters()
write_colfile(filename)

Save out the column file (all info stored we hope)

write_graindex_gv(filename)
write_pyFAI(filename, tthmin=0, tthmax=180)

Write file for Jerome Kieffer’s pyFAI fitting routine to use and run the refinment…

ImageD11.unitcell module

ImageD11.unitcell.A(h, k, l)
ImageD11.unitcell.B(h, k, l)
ImageD11.unitcell.BTmat(h1, h2, B, BI)

used for computing orientations

ImageD11.unitcell.C(h, k, l)
ImageD11.unitcell.F(h, k, l)
ImageD11.unitcell.I(h, k, l)
ImageD11.unitcell.P(h, k, l)
ImageD11.unitcell.R(h, k, l)
ImageD11.unitcell.cellfromstring(s)
ImageD11.unitcell.cosangles_many(ha, hb, gi)

Finds the cosines of angles between two lists of hkls in reciprocal metric gi

ImageD11.unitcell.cross(a, b)

a x b has length |a||b|sin(theta)

ImageD11.unitcell.degrees(x)
ImageD11.unitcell.filter_pairs(h1, h2, c2a, B, BI, tol=1e-05)

remove duplicate pairs for orientation searches h1 = reflections of ring1, N1 peaks h2 = reflections of ring2, N2 peaks c2a = cos angle between them, N1xN2 B = B matrix in reciprocal space BI = inverse in real space

ImageD11.unitcell.norm2(a)

Compute the unit 2 norm

ImageD11.unitcell.orient_BL(B, h1, h2, g1, g2)

Algorithm was intended to follow this one using 2 indexed reflections. W. R. Busing and H. A. Levy. 457. Acta Cryst. (1967). 22, 457

ImageD11.unitcell.radians(x)
ImageD11.unitcell.ubi_equiv(ubilist, ublist, tol=1e-08)

Two ubi are considered equivalent if they both index the peaks in the HKL0 array exactly

ImageD11.unitcell.unit(a)

Normalise vector a to unit length

class ImageD11.unitcell.unitcell(lattice_parameters, symmetry='P', verbose=0)

Bases: object

anglehkls(h1, h2)

Compute the angle between reciprocal lattice vectors h1, h2

ds(h)

computes 1/d for this hkl = hgh

getanglehkls(ring1, ring2)

Cache the previous pairs called for sorted by cos2angle

gethkls(dsmax)

Generate hkl list Argument dsmax is the d* limit (eg 1/d) Default of zero gives only the (000) reflection

assumes [h|k|l] < 200

gethkls_xfab(dsmax, spg)

Generate hkl list Argument dsmax is the d* limit (eg 1/d) Argument spg is the space group name, e.g. ‘R3-c’

makerings(limit, tol=0.001)

Makes a list of computed powder rings The tolerance is the difference in d* to decide if two peaks overlap

orient(ring1, g1, ring2, g2, verbose=0, crange=-1.0)

Compute an orientation matrix using cell parameters and the indexing of two reflections

Orientation matrix: A matrix such that h = A^1 g define 3 vectors t1,t2,t3 t1 is parallel to first peak (unit vector along g1) t2 is in the plane of both (unit vector along g1x(g1xg2)) t3 is perpendicular to both (unit vector along g1xg2)

tostring()

Write out a line containing unit cell information

ImageD11.unitcell.unitcell_from_parameters(pars)

ImageD11.write_graindex_gv module

ImageD11.write_graindex_gv.get_ds_string(g, ds_list)

Attempt to emulate the graindex (hkl) syntax. Obviously you would not want a (10,0,0) peak!

ImageD11.write_graindex_gv.make_ds_list(cell, limit=2.0)

Generates a list of d-spacings

ImageD11.write_graindex_gv.write_graindex_gv(outfilename, gv, tth, eta, omega, intensity, unitcell)

call with array of gvectors, tth, eta and omega vals

Using an indexing object to get the ring assignments

Module contents

ImageD11 - image analysis software for beamline ID11 at the ESRF

Copyright (C) 2005-2019 Jon Wright

This program is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation; either version 2 of the License, or (at your option) any later version.

This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.

You should have received a copy of the GNU General Public License along with this program; if not, write to the Free Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA