ImageD11 package¶
Subpackages¶
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')¶
-
class
ImageD11.ImageD11options.
FileType
(mode='r')¶ Bases:
object
-
class
ImageD11.ImageD11options.
GvectorFileType
(mode='r')¶
-
class
ImageD11.ImageD11options.
HdfFileType
(mode='r')¶
-
class
ImageD11.ImageD11options.
ImageFileType
(mode='r')¶
-
class
ImageD11.ImageD11options.
ParameterFileType
(mode='r')¶
-
class
ImageD11.ImageD11options.
SplineFileType
(mode='r')¶
-
class
ImageD11.ImageD11options.
UbiFileType
(mode='r')¶
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
- 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 :
-
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)
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 qhttp://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
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 wavelengthwavelength = 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 wavelengthwavelength = 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 centreomega 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 weightUpdated 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 centreomega 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.
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