Analysis (mpes.analysis)

Data analysis pipeline including background removal, segmentation and fitting

@author: R. Patrick Xian

class mpes.analysis.BoundedArea(image=None, shape=None, subimage=None)

Bounded area object from a parametric equation.

property mask

Subimage attribute as mask

setBoundary(pmz='linear', boundtype='>', **kwds)

Add bound to grid to redefine subgrid.

Parameters

pmz: str | ‘linear’

Parametrization (pmz) of the decision boundary (‘linear’ or ‘circular’).

boundtype: str | ‘>’

Bound region specification (‘>’ or ‘<’).

**kwds: keyword arguments

property subgrid

Substituent pixel coordinates of the image.

toMask(inbound=1, exbound=0)

Generate a scaled mask from existing shape.

Parameters

inbound: float | 1

Value for the pixels within the boundary.

exbound: float | 0

Value for the pixels outside the boundary.

Return

modmask: 2d array

Modified mask as a 2d array.

view(origin='lower', cmap='terrain_r', axes=True, **kwds)

Display the current mask.

Parameters

origin: str | ‘lower’

Location of the image origin.

cmap: str | ‘terrain_r’

Color map

axes: bool | True

Axes visibility option in plot.

**kwds: keyword arguments

Additional arguments for matplotlib.pyplot.imshow().

class mpes.analysis.EnergyCalibrator(biases=None, files=[], folder=[], file_sorting=True, traces=None, tof=None)

Electron binding energy calibration workflow.

addFeatures(ranges, refid=0, traces=None, infer_others=False, mode='append', **kwds)

Select or extract the equivalent landmarks (e.g. peaks) among all traces.

Parameters

ranges: list/tuple

Collection of feature detection ranges, within which an algorithm (i.e. 1D peak detector) with look for the feature.

refid: int | 0

Index of the reference trace (EDC).

traces: 2D array | None

Collection of energy dispersion curves (EDCs).

infer_others: bool | True

Option to infer the feature detection range in other traces (EDCs) from a given one.

mode: str | ‘append’

Specification on how to change the feature ranges (‘append’ or ‘replace’).

**kwds: keyword arguments

Dictionarized keyword arguments for trace alignment (See self.findCorrespondence())

calibrate(refid=0, ret=['coeffs'], **kwds)

Calculate the functional mapping between time-of-flight and the energy scale using optimization methods.

Parameters

refid: int | 0

The reference trace index (an integer).

ret: list | [‘coeffs’]

Options for return values (see mpes.analysis.calibrateE()).

**kwds: keyword arguments

See available keywords for mpes.analysis.calibrateE().

property dup

The duplication number.

featureExtract(ranges=None, traces=None, **kwds)

Select or extract the equivalent landmarks (e.g. peaks) among all traces.

Parameters

ranges: list/tuple | None

Range in each trace to look for the peak feature, [start, end].

traces: 2D array | None

Collection of 1D spectra to use for calibration.

**kwds: keyword arguments

See available keywords in mpes.analysis.peaksearch().

static findCorrespondence(sig_still, sig_mov, method='dtw', **kwds)

Determine the correspondence between two 1D traces by alignment.

Parameters

sig_still, sig_mov: 1D array, 1D array

Input 1D signals.

method: str | ‘dtw’

Method for 1D signal correspondence detection (‘dtw’ or ‘ptw’).

**kwds: keyword arguments

See available keywords for the following functions, (1) fastdtw.fastdtw() (when method=='dtw') (2) ptw.ptw.timeWarp() (when method=='ptw')

Return

pathcorr: list

Pixel-wise path correspondences between two input 1D arrays (sig_still, sig_mov).

property nfiles

The number of loaded files.

normalize(**kwds)

Normalize the spectra along an axis.

Parameters

**kwds: keyword arguments

See the keywords for mpes.utils.normspec().

property nranges

The number of specified feature ranges.

property ntraces

The number of loaded/calculated traces.

read(form='h5', tracename='', tofname='ToF')

Read traces (e.g. energy dispersion curves) from files.

Parameters

form: str | ‘h5’

Format of the files (‘h5’ or ‘mat’).

tracename: str | ‘’

Name of the group/attribute corresponding to the trace.

tofname: str | ‘ToF’

Name of the group/attribute corresponding to the time-of-flight.

saveParameters(form='h5', save_addr='./energy')

Save all the attributes of the workflow instance for later use (e.g. energy scale conversion).

Parameters

form: str | ‘h5’

The file format to save the attributes in (‘h5’/’hdf5’ or ‘mat’).

save_addr: str | ‘./energy’

The filename to save the files with.

view(traces, segs=None, peaks=None, show_legend=True, ret=False, display=True, backend='matplotlib', linekwds={}, linesegkwds={}, scatterkwds={}, legkwds={}, **kwds)

Display a plot showing line traces with annotation.

Parameters

traces: 2d array

Matrix of traces to visualize.

segs: list/tuple

Segments to be highlighted in the visualization.

peaks: 2d array

Peak positions for labelling the traces.

ret: bool

Return specification.

backend: str | ‘matplotlib’

Backend specification, choose between ‘matplotlib’ (static) or ‘bokeh’ (interactive).

linekwds: dict | {}

Keyword arguments for line plotting (see matplotlib.pyplot.plot()).

scatterkwds: dict | {}

Keyword arguments for scatter plot (see matplotlib.pyplot.scatter()).

legkwds: dict | {}

Keyword arguments for legend (see matplotlib.pyplot.legend()).

**kwds: keyword arguments

keyword

data type

meaning

maincolor

str

labels

list

Labels for each curve

xaxis

1d array

x (horizontal) axis values

title

str

Title of the plot

legend_location

str

Location of the plot legend

class mpes.analysis.Model(func, xvar, name=None)

Class of fitting curve models.

_costfunc(inits, xv, form='original')

Define the cost function of the optimization process.

fit(data, inits, method='leastsq', **fitkwds)

Run the optimization.

model_eval(params)

Evaluate the fitting model with given parameters.

static normalize(data)

Normalize n-dimensional data.

partial_eval(params, part=0)

Evaluate parts of a composite fitting model.

class mpes.analysis.MomentumCorrector(image, rotsym=6)

Momentum distortion correction and momentum calibration workflow.

_featureUpdate(center_det='centroidnn', **kwds)

Update selected features.

_imageUpdate()

Update distortion-corrected images.

applyDeformation(image, ret=True, **kwds)

Apply the deformation field to a specified image slice.

Parameters

image: 2D array

Image slice to apply the deformation.

ret: bool | True

Option to return the image after deformation.

**kwds: keyword arguments
rdeform, cdeform

2D array, 2D array | self.rdeform_field, self.cdeform_field Row- and column-ordered deformation fields.

interp_order

int | 1 Interpolation order.

others

See scipy.ndimage.map_coordinates().

calcGeometricDistances()

Calculate geometric distances involving the center and the vertices. Distances calculated include center-vertex and nearest-neighbor vertex-vertex distances.

calcSymmetryScores(symtype='rotation')

Calculate the symmetry scores from geometric quantities.

Paramters

symtype: str | ‘rotation’

Type of symmetry.

calibrate(image, point_from, point_to, dist, ret='coeffs', **kwds)

Calibration of the momentum axes. Obtain all calibration-related values, return only the ones requested.

Parameters

image: 2d array

Image slice to construct the calibration function.

point_from, point_to: list/tuple, list/tuple

Pixel coordinates of the two special points in (row, col) ordering.

dist: float

Distance between the two selected points in inverse Angstrom.

ret: str | ‘coeffs’

Specification of return values (‘axes’, ‘extent’, ‘coeffs’, ‘grid’, ‘func’, ‘all’).

**kwds: keyword arguments

See arguments in mpes.analysis.calibrateE().

Return

Specified calibration parameters in a dictionary.

coordinateTransform(type, keep=False, ret=False, interp_order=1, mapkwds={}, **kwds)

Apply a pixel-wise coordinate transform to an image.

Parameters

type: str

Type of deformation to apply to image slice.

keep: bool | False

Option to keep the specified coordinate transform.

ret: bool | False

Option to return transformed image slice.

interp_order: int | 1

Interpolation order for filling in missed pixels.

mapkwds: dict | {}

Additional arguments passed to scipy.ndimage.map_coordinates().

**kwds: keyword arguments

Additional arguments in specific deformation field. See symmetrize.sym module.

correct(axis, use_composite_transform=False, update=False, use_deform_field=False, updatekwds={}, **kwds)

Apply a 2D transform to a stack of 2D images (3D) along a specific axis.

Parameters

axis: int

Axis for slice selection.

use_composite_transform: bool | False

Option to use the composite transform involving the rotation.

update: bool | False

Option to update the existing figure attributes.

use_deform_field: bool | False

Option to use deformation field for distortion correction.

**kwds: keyword arguments

keyword

data type

meaning

image

2d array

3D image for correction

dfield

list/tuple

row and column deformation field

warping

2d array

2D transform correction matrix

featureExtract(image, direction='ccw', type='points', center_det='centroidnn', symscores=True, **kwds)
Extract features from the selected 2D slice.

Currently only point feature detection is implemented.

Parameters

image: 2d array

The image slice to extract features from.

direction: str | ‘ccw’

The circular direction to reorder the features in (‘cw’ or ‘ccw’).

type: str | ‘points’

The type of features to extract.

center_det: str | ‘centroidnn’

Specification of center detection method (‘centroidnn’, ‘centroid’, None).

**kwds: keyword arguments

Extra keyword arguments for symmetrize.pointops.peakdetect2d().

property features

Dictionary of detected features for the symmetrization process. self.features is a derived attribute from existing ones.

static getWarpFunction(**kwds)

Construct warping function to apply to other datasets. # TODO: turn this into a fully operational method.

importBinningParameters(parp)

Import parameters of binning used for correction image from parallelHDF5Processor Class instance

Parameters

parp: instance of the ParallelHDF5Processor class

Import parameters used for creation of the distortion-corrected image.

intensityTransform(type='rot_sym', **kwds)

Apply pixel-wise intensity transform.

Parameters

type: str | ‘rot_sym’

Type of intensity transform.

**kwds: keyword arguments

linWarpEstimate(weights=1, 1, 1, optfunc='minimize', optmethod='Nelder-Mead', ret=True, warpkwds={}, **kwds)

Estimate the homography-based deformation field using landmark correspondences.

Parameters

weights: tuple/list/array | (1, 1, 1)

Weights added to the terms in the optimizer. The terms are assigned to the cost functions of (1) centeredness, (2) center-vertex symmetry, (3) vertex-vertex symmetry, respectively.

optfunc, optmethod: str/func, str | ‘minimize’, ‘Nelder-Mead’

Name of the optimizer function and the optimization method. See description in mpes.analysis.sym.target_set_optimize().

ret: bool | True

Specify if returning the corrected image slice.

warpkwds: dictionary | {}

Additional arguments passed to symmetrize.sym.imgWarping().

**kwds: keyword arguments

keyword

data type

meaning

niter

int

Maximum number of iterations

landmarks

list/array

Symmetry landmarks selected for registration

fitinit

tuple/list

Initial conditions for fitting

Return

Corrected 2D image slice (when ret=True is specified in the arguments).

resetDeformation(**kwds)

Reset the deformation field.

rotate(angle='auto', ret=False, **kwds)

Rotate 2D image in the homogeneous coordinate.

Parameters

angle: float/str

Angle of rotation (specify ‘auto’ to use automated estimation).

ret: bool | False

Return specification (True/False)

**kwds: keyword arguments

keyword

data type

meaning

image

2d array

2D image for correction

center

tuple/list

pixel coordinates of the image center

scale

float

scaling factor in rotation

See symmetrize.sym.sym_pose_estimate() for other keywords.

saveImage(form='tiff', save_addr='./', dtyp='float32', **kwds)

Save the distortion-corrected dataset (image only, without axes).

Parameters

form: str | ‘tiff’

File format for saving the corrected image (‘tiff’ or ‘mat’).

save_addr: str | ‘./’

The address to save the file at.

dtyp: str | ‘float32’

Data type (in case conversion if needed).

**kwds: keyword arguments

See keywords from tifffile.imsave().

saveParameters(form='h5', save_addr='./momentum')

Save all the attributes of the workflow instance for later use (e.g. momentum scale conversion, reconstructing the warping map function).

Parameters

form: str | ‘h5’

File format to for saving the parameters (‘h5’/’hdf5’, ‘mat’)

save_addr: str | ‘./momentum’

The address for the to be saved file.

selectSlice2D(selector, axis=2)

Select (hyper)slice from a (hyper)volume.

Parameters

selector: slice object/list/int

Selector along the specified axis to extract the slice (image). Use the construct slice(start, stop, step) to select a range of images and sum them. Use an integer to specify only a particular slice.

axis: int | 2

Axis along which to select the image.

splineWarpEstimate(image, include_center=True, fixed_center=True, iterative=False, interp_order=1, update=False, ret=False, **kwds)

Estimate the spline deformation field using thin plate spline registration.

Parameters

image: 2D array

Image slice to be corrected.

include_center: bool | True

Option to include the image center/centroid in the registration process.

fixed_center: bool | True

Option to have a fixed center during registration-based symmetrization.

iterative: bool | False

Option to use the iterative approach (may not work in all cases).

interp_order: int | 1

Order of interpolation (see scipy.ndimage.map_coordinates()).

update: bool | False

Option to keep the spline-deformed image as corrected one.

ret: bool | False

Option to return corrected image slice.

**kwds: keyword arguments
landmarks

list/array | self.pouter_ord Landmark positions (row, column) used for registration.

new_centers

dict | {} User-specified center positions for the reference and target sets. {‘lmkcenter’: (row, col), ‘targcenter’: (row, col)}

property symscores

Dictionary of symmetry-related scores.

static transform(points, transmat)

Coordinate transform of a point set in the (row, column) formulation.

Parameters

points: list/array

Cartesian pixel coordinates of the points to be transformed.

transmat: 2D array

The transform matrix.

Return

Transformed point coordinates.

update(content, **kwds)

Update specific attributes of the class.

Parameters

content: str | ‘all’

‘feature’ = update only feature attributes

‘image’ = update only image-related attributes

‘all’ = update both feature and image-related attributes

**kwds: keyword arguments

Extra keyword arguments passed into self._featureUpdate().

updateDeformation(rdeform, cdeform, reset=False, **kwds)

Update the deformation field.

Parameters

rdeform, cdeform: 2D array, 2D array

Row- and column-ordered deformation fields.

reset: bool | False

Option to reset the deformation field.

**kwds: keyword arguments

See mpes.analysis.MomentumCorrector.resetDeformation().

view(origin='lower', cmap='terrain_r', figsize=4, 4, points={}, annotated=False, display=True, backend='matplotlib', ret=False, imkwds={}, scatterkwds={}, crosshair=False, radii=[50, 100, 150], crosshair_thickness=1, **kwds)

Display image slice with specified annotations.

Parameters

origin: str | ‘lower’

Figure origin specification (‘lower’ or ‘upper’).

cmap: str | ‘terrain_r’

Colormap specification.

figsize: tuple/list | (4, 4)

Figure size.

points: dict | {}

Points for annotation.

annotated: bool | False

Option for annotation.

display: bool | True

Display option when using bokeh to render interactively.

backend: str | ‘matplotlib’

Visualization backend specification. :’matplotlib’: use static display rendered by matplotlib. :’bokeh’: use interactive display rendered by bokeh.

ret: bool | False

Option to return figure and axis objects.

imkwd: dict | {}

Keyword arguments for matplotlib.pyplot.imshow().

crosshair: bool | False

Display option to plot circles around center self.pcent. Works only in bokeh backend.

radii: list | [50,100,150]

Radii of circles to plot when crosshair optin is activated.

crosshair_thickness: int | 1

Thickness of crosshair circles.

**kwds: keyword arguments

General extra arguments for the plotting procedure.

mpes.analysis._datacheck_peakdetect(x_axis, y_axis)

Input format checking

mpes.analysis._rotate2d(image, center, angle, scale=1)

2D matrix scaled rotation carried out in the homogenous coordinate.

Parameters

image: 2d array

Image matrix.

center: tuple/list

Center of the image (row pixel, column pixel).

angle: numeric

Angle of image rotation.

scale: numeric | 1

Scale of image rotation.

Returns

image_rot: 2d array

Rotated image matrix.

rotmat: 2d array

Rotation matrix in the homogeneous coordinate system.

mpes.analysis._signedmask(imr, imc, maskr, maskc, sign)

Generate a binary mask using the masked coordinates.

Parameters

imr, imc: int

Row and column size of the image.

maskr, maskc: 1D array

Row and column coordinates of the masked pixels.

sign: int/str

Value of the masked region, (0, 1, ‘nan’, or ‘xnan’).

Return

mask: 2D array

Mask matrix.

mpes.analysis.applyWarping(imgstack, axis, hgmat)

Apply warping transform for a stack of images along an axis.

Parameters

imgstack: 3D array

Image stack before warping correction.

axis: int

Axis to iterate over to apply the transform.

hgmat: 2D array

Homography matrix.

Return

imstack_transformed: 3D array

Stack of images after correction for warping.

mpes.analysis.apply_mask_along(arr, mask, axes=None)

Apply a mask in a low dimensional slice throughout a high-dimensional array.

Parameters

arr: nD array

Multidimensional array for masking.

mask: nD array

Mask to apply.

axes: list/tuple of int | None

The axes to apply the mask to.

Return maskedarr: nD array

Masked multidimensional array.

mpes.analysis.bandpath_map(bsvol, pathr=None, pathc=None, path_coords=None, eaxis=2, method='analog')

Extract band diagram map from 2D/3D data.

Parameters

bsvol: 2D/3D array

Volumetric band structure data.

pathr, pathc: 1D array | None, None

Row and column pixel coordinates along the band path (ignored if path_coords is given).

path_coords: 2D array | None

Combined row and column pixel coordinates of the band path.

eaxis: int | 2

Energy axis index.

method: str | ‘analog’

Method for generating band path map (‘analog’ or ‘digital’). :’analog’: Using an interpolation scheme to calculate the exact pixel values. :’digital’: Using only the approximating pixel values (Bresenham’s algorithm).

Return

bpm: 2D array

Band path map (BPM) sampled from the volumetric data.

mpes.analysis.blocknorm(data, mavg_axis=0, blockwidth=1)

Block-thresholding 2D data.

Parameters

data: ndarray

Data to normalize.

mavg_axis: int | 0

Axis to move the block along.

blockwidth: int | 1

Width of the moving block.

Return

datanorm: ndarray

Block-normalized data.

mpes.analysis.bootstrapfit(data, axval, model, params, axis=0, dfcontainer=None, pbar=False, pbenv='classic', ret='all', **kwds)

Line-by-line fitting via bootstrapping fitted parameters from one line to the next.

Parameters

data: ndarray

Data used in fitting.

axval: list/numeric array

Value for the axis.

model: lmfit Model object

The fitting model.

params: lmfit Parameters object

Initial guesses for fitting parameters.

axis: int | 0

The axis of the data to fit.

dfcontainer: pandas DataFrame | None

Dataframe container for the fitting parameters.

pbar: bool | False

Progress bar condition.

pbenv: str | ‘classic’

Progress bar environment (‘classic’ for generic version, ‘notebook’ for notebook compatible version).

**kwds: keyword arguments

keyword

data type

meaning

maxiter

int

maximum iteration per fit (default = 20)

concat

bool

concatenate the fit parameters to DataFrame input

False (default) = no concatenation, use an empty DataFrame to start True = with concatenation to input DataFrame

bgremove

bool

toggle for background removal (default = True)

flipped

bool

toggle for fitting start position

(if flipped, fitting start from the last line)

limpropagate

bool

verbose

bool

toggle for output message (default = False)

Returns

df_fit: pandas DataFrame

Dataframe container populated with obtained fitting parameters.

data_nobg: ndarray

Background-removed (Shirley-type) traces.

mpes.analysis.build_dynamic_matrix(fitparams, display_range=slice(None, None, None), pre_t0_range=slice(None, 1, None))

Construct the dynamic matrix from the fitting results. For each fitting parameter, construct time-dependent value, time-dependent absolute and relative changes.

Parameters

fitparams: 3D ndarray

fitting output

display_range: slice object | slice(None, None, None)

display time range of the fitting parameters (default = full range)

pre_t0_range: slice object | slice(None, 1, None)

time range regarded as before time-zero

Returns

dyn_matrix: 4D ndarray

calculated dynamic matrix

mpes.analysis.calibrateE(pos, vals, order=3, refid=0, ret='func', E0=None, Eref=None, t=None, aug=1, method='lstsq', **kwds)

Energy calibration by nonlinear least squares fitting of spectral landmarks on a set of (energy dispersion curves (EDCs). This amounts to solving for the coefficient vector, a, in the system of equations T.a = b. Here T is the differential drift time matrix and b the differential bias vector, and assuming that the energy-drift-time relationship can be written in the form, E = sum_n (a_n * t**n) + E0

Parameters

pos: list/array

Positions of the spectral landmarks (e.g. peaks) in the EDCs.

vals: list/array

Bias voltage value associated with each EDC.

order: int | 3

Polynomial order of the fitting function.

refid: int | 0

Reference dataset index, varies from 0 to vals.size - 1.

ret: str | ‘func’

Return type, including ‘func’, ‘coeffs’, ‘full’, and ‘axis’ (see below).

E0: float | None

Constant energy offset.

t: numeric array | None

Drift time.

aug: int | 1

Fitting dimension augmentation (1=no change, 2=double, etc).

Returns

pfunc: partial function

Calibrating function with determined polynomial coefficients (except the constant offset).

ecalibdict: dict

A dictionary of fitting parameters including the following, :coeffs: Fitted polynomial coefficients (the a’s). :offset: Minimum time-of-flight corresponding to a peak. :Tmat: the T matrix (differential time-of-flight) in the equation Ta=b. :bvec: the b vector (differential bias) in the fitting Ta=b. :axis: Fitted energy axis.

mpes.analysis.calibrateK(img, pxla, pxlb, k_ab=None, kcoorda=None, kcoordb=[0.0, 0.0], equiscale=False, ret=['axes'])

Momentum axes calibration using the pixel positions of two symmetry points (a and b) and the absolute coordinate of a single point (b), defaulted to [0., 0.]. All coordinates should be specified in the (row_index, column_index) format. See the equiscale option for details on the specifications of point coordinates.

Parameters

img: 2D array

An energy cut of the band structure.

pxla, pxlb: list/tuple/1D array

Pixel coordinates of the two symmetry points (a and b). Point b has the default coordinates [0., 0.] (see below).

k_ab: float | None

The known momentum space distance between the two symmetry points.

kcoorda: list/tuple/1D array | None

Momentum coordinates of the symmetry point a.

kcoordb: list/tuple/1D array | [0., 0.]

Momentum coordinates of the symmetry point b (krow, kcol), default to k-space center.

equiscale: bool | False

Option to adopt equal scale along both the row and column directions. :True: Use a uniform scale for both x and y directions in the image coordinate system. This applies to the situation where the points a and b are (close to) parallel with one of the two image axes. :False: Calculate the momentum scale for both x and y directions separately. This applies to the situation where the points a and b are sufficiently different in both x and y directions in the image coordinate system.

ret: list | [‘axes’]

Return type specification, options include ‘axes’, ‘extent’, ‘coeffs’, ‘grid’, ‘func’, ‘all’.

Returns

k_row, k_col: 1D array

Momentum coordinates of the row and column.

axis_extent: list

Extent of the two momentum axis (can be used directly in imshow).

k_rowgrid, k_colgrid: 2D array

Row and column mesh grid generated from the coordinates (can be used directly in pcolormesh).

mpes.analysis.circmask(img, rcent, ccent, rad, sign=1, ret='mask', **kwds)

Use a circular binary mask to cover an image.

Parameters

img: 2D array

Input image to be masked.

rcent: float

Row center position.

ccent: float

Column center position.

rad: float

Radius of circle.

sign: int/str | 1

Value of the masked region (0, 1, ‘nan’ or ‘xnan’). ‘xnan’ means the masked region is 1 and the other region nan.

ret: str | ‘mask’

Return type (‘mask’, ‘masked_image’)

kwds: keyword arguments

Return

cmask or cmask*img: 2D array

Mask only or masked image

mpes.analysis.curvature2d(image, cx=1, cy=1)

Implementation of 2D curvature calculation. The formula follows Zhang et al. Rev. Sci. Instrum. 82, 043712 (2011).

Parameters

image: 2D array

2D image obtained from measurement.

cx, cy: numeric, numeric | 1, 1

Scaling parameters in x and y directions.

mpes.analysis.func_add(*funcs)

Addition of an arbitray number of functions

Parameters

*funcs: list/tuple

functions to combine

Returns

funcsum: function

functional sum

mpes.analysis.func_update(func, suffix='')

Attach a suffix to parameter names and their instances in the expression of a function

Parameters

func: function

input function

suffix: str | ‘’

suffix to attach to parameter names

Returns

params: list of str

updated function parameters

expr: str

updated function expression

mpes.analysis.gaussian(feval=False, vardict=None)

1D Gaussian lineshape model. Returns numerical values if feval=True.

Parameters

feval: bool | False

Option to evaluate function.

vardict: dict | None

Dictionary containing values for the variables named as follows (as dictionary keys).

amp function amplitude or scaling factor.

xvar x values (energy values in a lineshape).

ctr center position.

sig standard deviation.

mpes.analysis.gradn(array, axes, **kwds)

Calculate nth-order gradients of the array along different directions.

Parameters

array: numpy array

N-dimensional matrix for calculating the gradient.

axes: int/list/tuple/1D array

A sequence of axes (from first to last) to calculate the gradient. When input a single integer, the gradient is calculated along that particular axis. For example, the 4th-order mixed gradient d4f/(dxdydxdy) requires the sequence (1, 0, 1, 0).

**kwds: keyword arguments

See numpy.gradient().

mpes.analysis.image_interpolator(image, iptype='RGI')

Construction of an image interpolator.

Parameters

image: 2D array

2D image for interpolation.

iptype: str | ‘RGI’

Type of the interpolator.

Return

interp: interpolator instance

Instance of an interpolator.

mpes.analysis.interp_slice(data, pathr=None, pathc=None, path_coords=None, iptype='RGI')

Slicing 2D/3D data through interpolation.

Parameters

data: 2D/3D array

Data array for slicing.

pathr, pathc: list/tuple/array, list/tuple/array

Row and column coordinates of the interpolated path.

path_coords: array

Cartesian coordinates of the interpolated path.

iptype: str | ‘RGI’

Type of interpolator.

mpes.analysis.line_generator(A, B, npoints, endpoint=True, ret='separated')

Generate intermediate points in a line segment AB (A to B) given endpoints.

Parameters

A, B: tuple/list, tuple/list

Pixel coordinates of the endpoints of the line segment.

npoints: numeric

Number of points in the line segment.

endpoint: bool | True

Option to include the endpoint (B) in the line coordinates.

ret: str | ‘separated’

Option to return line coordinates (‘separated’ or ‘joined’).

mpes.analysis.peakdetect1d(y_axis, x_axis=None, lookahead=200, delta=0)

Function for detecting local maxima and minima in a signal. Discovers peaks by searching for values which are surrounded by lower or larger values for maxima and minima respectively

Converted from/based on a MATLAB script at: http://billauer.co.il/peakdet.html

Parameters

y_axis: list

A list containing the signal over which to find peaks

x_axis: list | None

A x-axis whose values correspond to the y_axis list and is used in the return to specify the position of the peaks. If omitted an index of the y_axis is used.

lookahead: int | 200

distance to look ahead from a peak candidate to determine if it is the actual peak ‘(samples / period) / f’ where ‘4 >= f >= 1.25’ might be a good value

delta: numeric | 0

this specifies a minimum difference between a peak and the following points, before a peak may be considered a peak. Useful to hinder the function from picking up false peaks towards to end of the signal. To work well delta should be set to delta >= RMSnoise * 5.

Returns

max_peaks: list

positions of the positive peaks

min_peaks: list

positions of the negative peaks

mpes.analysis.peakdetect2d(img, method='daofind', **kwds)

Peak-like feature detection in a 2D image.

Parameters

img: 2D array

Image matrix.

method: str | ‘daofind’

Detection method (‘daofind’ or ‘maxlist’).

**kwds: keyword arguments

Arguments passed to the specific methods chosen.

daofind

See astropy.stats.sigma_clipped_stats() and photutils.detection.DAOStarFinder().

sigma: float | 5.0

Standard deviation of the clipping Gaussian.

fwhm: float | 3.0

FWHM of the convoluting Gaussian kernel.

threshfactor: float | 8

Intensity threshold for background-foreground separation (foreground is above threshold).

maxlist

See skimage.feature.peak_local_max().

mindist: float | 10

Minimal distance between two local maxima.

numpeaks: int | 7

Maximum number of detected peaks.

Return

pks: 2D array

Pixel coordinates of detected peaks, in (column, row) ordering.

mpes.analysis.peaksearch(traces, tof, ranges=None, method='range-limited', pkwindow=3, plot=False)

Detect a list of peaks in the corresponding regions of multiple EDCs

Parameters

traces: 2D array

Collection of EDCs.

tof: 1D array

Time-of-flight values.

ranges: list of tuples/lists | None

List of ranges for peak detection in the format [(LowerBound1, UpperBound1), (LowerBound2, UpperBound2), ….].

method: str | ‘range-limited’

Method for peak-finding (‘range-limited’ and ‘alignment’).

pkwindow: int | 3

Window width of a peak(amounts to lookahead in mpes.analysis.peakdetect1d).

plot: bool | False

Specify whether to display a custom plot of the peak search results.

Returns

pkmaxs: 1D array

Collection of peak positions.

mpes.analysis.perspectiveWarping(img, landmarks, targs, ret='image')

Perform image warping based on a generic affine transform (homography).

Parameters

img: 2D array

Input image (distorted).

landmarks: list/array

List of pixel positions of the reference points.

targs: list/array

List of pixel positions of the target points.

Returns

imgaw: 2D array

Image after affine warping.

maw: 2D array

Homography matrix for the tranform.

mpes.analysis.points2path(pointsr, pointsc, method='analog', npoints=None, ret='separated')

Calculate ordered pixel cooridnates along a path defined by specific intermediate points. The approach constructs the path using a set of line segments bridging the specified points, therefore it is also able to trace the sequence indices of these special points.

Parameters

pointsr, pointsc: list/tuple/array

The row and column pixel coordinates of the special points along the sampling path.

method: str | ‘analog’

Method of sampling.

npoints: list/tuple | None

Number of points along each segment.

ret: str | ‘separated’

Specify if return combined (‘combined’) or separated (‘separated’) row and column coordinates.

Returns

polyr, polyc: 1D array

Pixel coordinates along the path traced out sequentially.

pid: 1D array

Pointwise indices of the special lpoints.

mpes.analysis.rangeConvert(x, xrng, pathcorr)

Convert value range using a pairwise path correspondence (e.g. obtained from time warping techniques).

Parameters

x: 1D array

Values of the x axis (e.g. time-of-flight values).

xrng: list/tuple

Boundary value range on the x axis.

pathcorr: list/tuple

Path correspondence between two 1D arrays in the following form, [(id_1_trace_1, id_1_trace_2), (id_2_trace_1, id_2_trace_2), …]

Return

xrange_trans: tuple

Transformed range according to the path correspondence.

mpes.analysis.rectmask(img, rcent, ccent, shift, direction='row', sign=1, ret='mask', **kwds)

Use a rectangular binary mask to cover an image.

Parameters

img: 2D array

Input image to be masked

rcent: float

Row center position

ccent: float

Column center position

shift: int/list of int

Pixel shifts

direction: str | ‘row’

Direction to apply the shift to, ‘row’ or ‘column’ indicates row-wise or column-wise shift for generating the rectangular mask

sign: int/str | 1

Value of the masked region (0, 1, ‘nan’ or ‘xnan’). ‘xnan’ means the masked region is 1 and the other region nan.

ret: str | ‘mask’

Return type (‘mask’, ‘masked_image’)

**kwds: keyword arguments

Return

cmask or cmask*img: 2D array

Mask only or masked image

mpes.analysis.regionExpand(mask, **kwds)

Expand the region of a binarized image around a line position

Parameters

mask: numeric binarized 2D array

the mask to be expanded

**kwds: keyword arguments

keyword

data type

meaning

method

str

method of choice (‘offset’, ‘growth’)

value

numeric

value to be assigned to the masked

linecoords

2D array

contains x and y positions of the line

axoffsets

tuple/list

[downshift upshift] pixel number

clipbounds

tuple/list

bounds in the clipping direction

selem

ndarray

structuring element

Return

mask: numeric 2D array

modified mask (returns the original mask if insufficient arguments are provided for the chosen method for region expansion)

mpes.analysis.ridgeDetect(mask, method='mask_mean_y', **kwds)

Detect the band ridges using selected methods.

Parameters

mask: numeric 2D array

the 2D integer-valued mask with labeled bands

method: str

the method used for ridge detection ‘mask_mean_y’: mean mask position along y direction (default) ‘mask_mean_x’: mean mask position along x direction

**kwds: keyword arguments

keyword

data type

meaning

x

int/float

x axis coordinates

y

int/float

y axis coordinates

Return

ridges: list of dataframes

the ridge coordinates

mpes.analysis.segment2d(img, nbands=1, **kwds)

Electronic band segmentation using local thresholding and connected component labeling

Parameters

img: 2D numeric array

the 2D matrix to segment

nbands: int

number of electronic bands

**kwds: keyword arguments

Return

imglabeled: 2D numeric array

labeled mask

mpes.analysis.shirley(x, y, tol=1e-05, maxiter=20, explicit=False, warning=False)

Calculate the 1D best Shirley-Proctor-Sherwood background S for a dataset (x, y). A. Proctor, P. M. A. Sherwood, Anal. Chem. 54 13 (1982). The function is adapted from Kane O’Donnell’s routine 1. Finds the biggest peak 2. Use the minimum value on either side of this peak as the terminal points of the Shirley background. 3. Iterate over the process within maximum allowed iteration (maxiter) to reach the tolerance level (tol).

Parameters

x: 1D numeric array

The photoelectron energy axis.

y: 1D numeric array

The photoemission intensity axis.

tol: float | 1e-5

The fitting tolerance.

maxiter: int | 20

The maximal iteration.

explicit: bool | False

Option for explicit display of iteration number.

warning: bool | False

Option to display of warnings during calculation.

Return

sbg: 1D numeric array

Calculated Shirley background.

mpes.analysis.shirley2d(x, y, tol=1e-05, maxiter=20, explicit=False, warning=False)

2D Shirley background removal

Parameters

x: 1D numeric array

Photoemission energy axis.

y: 2D numeric array

Photoemission intensity matrix.

tol: float | 1e-5

The fitting tolerance.

maxiter: int | 20

The maximal iteration.

explicit: bool | False

Option for explicit display of iteration number.

warning: bool | False

Option to display of warnings during calculation.

mpes.analysis.vertexGenerator(center, fixedvertex=None, cvd=None, arot=None, nside=None, direction=- 1, scale=1, diagdir=None, ret='all', rettype='float32')

Generation of the vertices of symmetric polygons.

Parameters

center: (int, int)

Pixel positions of the symmetry center (row pixel, column pixel).

fixedvertex: (int, int) | None

Pixel position of the fixed vertex (row pixel, column pixel).

cvd: numeric | None

Center-vertex distance.

arot: float | None

Spacing in angle of rotation.

nside: int | None

The total number of sides for the polygon.

direction: int | -1

Direction of angular rotation (1 = counterclockwise, -1 = clockwise)

scale: float | 1

Radial scaling factor.

diagdir: str | None

Diagonal direction of the polygon (‘x’ or ‘y’).

ret: str | ‘all’

Return type. Specify ‘all’ returns all vertices, specify ‘generated’ returns only the generated ones (without the fixedvertex in the argument).

Return

vertices: 2D array

Collection of generated vertices.

mpes.analysis.voigt(feval=False, vardict=None)

1D Voigt lineshape model. Returns numerical values if feval=True.

Parameters

feval: bool | False

Option to evaluate function.

vardict: dict | None

Dictionary containing values for the variables named as follows (as dictionary keys).

amp function amplitude or scaling factor.

xvar x values (energy values in a lineshape).

ctr center position.

sig standard deviation of the Gaussian component.

gam linewidth of the Lorentzian component.