API documentation

Classes

PixelatedSTEM

class pixstem.pixelated_stem_class.PixelatedSTEM(*args, **kw)[source]
angular_mask(angle0, angle1, centre_x_array=None, centre_y_array=None)[source]

Get a bool array with True values between angle0 and angle1. Will use the (0, 0) point as given by the signal as the centre, giving an “angular” slice. Useful for analysing anisotropy in diffraction patterns.

Parameters:
  • angle1 (angle0,) –
  • centre_y_array (centre_x_array,) – Has to have the same shape as the navigation axis of the signal.
Returns:

mask_array – The True values will be the region between angle0 and angle1. The array will have the same dimensions as the signal.

Return type:

Numpy array

Examples

>>> import pixstem.dummy_data as dd
>>> s = dd.get_holz_simple_test_signal()
>>> s.axes_manager.signal_axes[0].offset = -25
>>> s.axes_manager.signal_axes[1].offset = -25
>>> mask_array = s.angular_mask(0.5*np.pi, np.pi)
angular_slice_radial_integration(angleN=20, centre_x=None, centre_y=None, slice_overlap=None, show_progressbar=True)[source]

Do radial integration of different angular slices. Useful for analysing anisotropy in round diffraction features, such as diffraction rings from polycrystalline materials or higher order laue zone rings.

Parameters:
  • angleN (int, default 20) – Number of angular slices. If angleN=4, each slice will be 90 degrees. The integration will start in the top left corner (0, 0) when plotting using s.plot(), and go clockwise.
  • centre_y (centre_x,) – If given as int, all the diffraction patterns will have the same centre position. Each diffraction pattern can also have different centre position, by passing a NumPy array with the same dimensions as the navigation axes. Note: in either case both x and y values must be given. If one is missing, both will be set from the signal (0., 0.) positions. If no values are given, the (0., 0.) positions in the signal will be used.
  • slice_overlap (float, optional) – Amount of overlap between the slices, given in fractions of angle slice (0 to 1). For angleN=4, each slice will be 90 degrees. If slice_overlap=0.5, each slice will overlap by 45 degrees on each side. The range of the slices will then be: (-45, 135), (45, 225), (135, 315) and (225, 45). Default off: meaning there is no overlap between the slices.
  • show_progressbar (bool) – Default True
Returns:

signal – With one more navigation dimensions (the angular slices) compared to the input signal.

Return type:

HyperSpy 1D signal

Examples

>>> import pixstem.dummy_data as dd
>>> s = dd.get_holz_simple_test_signal()
>>> s_com = s.center_of_mass(show_progressbar=False)
>>> s_ar = s.angular_slice_radial_integration(
...     angleN=10, centre_x=s_com.inav[0].data,
...     centre_y=s_com.inav[1].data, slice_overlap=0.2,
...     show_progressbar=False)
>>> s_ar.plot()
center_of_mass(threshold=None, mask=None, lazy_result=False, show_progressbar=True, chunk_calculations=None)[source]

Get the centre of the STEM diffraction pattern using center of mass. Threshold can be set to only use the most intense parts of the pattern. A mask can be used to exclude parts of the diffraction pattern.

Parameters:
  • threshold (number, optional) – The thresholding will be done at mean times this threshold value.
  • mask (tuple (x, y, r), optional) – Round mask centered on x and y, with radius r.
  • lazy_result (bool, optional) – If True, will not compute the data directly, but return a lazy signal. Default False
  • show_progressbar (bool, optional) – Default True
  • chunk_calculations (tuple, optional) – Chunking values when running the calculations.
Returns:

s_com – DPCSignal with beam shifts along the navigation dimension and spatial dimensions as the signal dimension(s).

Return type:

DPCSignal

Examples

With mask centered at x=105, y=120 and 30 pixel radius

>>> import pixstem.dummy_data as dd
>>> s = dd.get_disk_shift_simple_test_signal()
>>> mask = (25, 25, 10)
>>> s_com = s.center_of_mass(mask=mask, show_progressbar=False)
>>> s_color = s_com.get_color_signal()

Also threshold

>>> s_com = s.center_of_mass(threshold=1.5, show_progressbar=False)

Get a lazy signal, then calculate afterwards

>>> s_com = s.center_of_mass(lazy_result=True, show_progressbar=False)
>>> s_com.compute(progressbar=False)
correct_bad_pixels(bad_pixel_array, lazy_result=True, show_progressbar=True)[source]

Correct bad pixels by getting mean value of neighbors.

Note: this method is currently not very optimized with regards to memory use, so currently be careful when using it on large datasets.

Parameters:
  • bad_pixel_array (array-like) –
  • lazy_result (bool) – Default True.
  • show_progressbar (bool) – Default True
Returns:

signal_corrected

Return type:

PixelatedSTEM or LazyPixelatedSTEM

Examples

>>> s = ps.dummy_data.get_hot_pixel_signal()
>>> s_hot_pixels = s.find_hot_pixels(
...     show_progressbar=False, lazy_result=True)
>>> s_corr = s.correct_bad_pixels(s_hot_pixels)

Dead pixels

>>> s = ps.dummy_data.get_dead_pixel_signal()
>>> s_dead_pixels = s.find_dead_pixels(
...     show_progressbar=False, lazy_result=True)
>>> s_corr = s.correct_bad_pixels(s_dead_pixels)

Combine both dead pixels and hot pixels

>>> s_bad_pixels = s_hot_pixels + s_dead_pixels
>>> s_corr = s.correct_bad_pixels(s_dead_pixels)

See also

find_dead_pixels(), find_hot_pixels()

find_dead_pixels(dead_pixel_value=0, mask_array=None, lazy_result=False, show_progressbar=True)[source]

Find dead pixels in the diffraction images.

Parameters:
  • dead_pixel_value (scalar) – Default 0
  • mask_array (Boolean Numpy array) –
  • lazy_result (bool) – If True, return a lazy signal. If False, compute the result and return a non-lazy signal. Default False.
  • show_progressbar (bool) –
Returns:

s_dead_pixels – With dead pixels as True, rest as False.

Return type:

HyperSpy 2D signal

Examples

>>> s = ps.dummy_data.get_dead_pixel_signal()
>>> s_dead_pixels = s.find_dead_pixels(show_progressbar=False)

Using a mask array

>>> import numpy as np
>>> mask_array = np.zeros((128, 128), dtype=np.bool)
>>> mask_array[:, 100:] = True
>>> s = ps.dummy_data.get_dead_pixel_signal()
>>> s_dead_pixels = s.find_dead_pixels(
...     mask_array=mask_array, show_progressbar=False)

Getting a lazy signal as output

>>> s_dead_pixels = s.find_dead_pixels(
...     lazy_result=True, show_progressbar=False)

See also

find_hot_pixels(), correct_bad_pixels()

find_hot_pixels(threshold_multiplier=500, mask_array=None, lazy_result=True, show_progressbar=True)[source]

Find hot pixels in the diffraction images.

Note: this method will be default return a lazy signal, since the size of the returned signal is the same shape as the original signal. So for large datasets actually calculating computing the results can use a lot of memory.

In addition, this signal is currently not very optimized with regards to memory use, so be careful when using this method for large datasets.

Parameters:
  • threshold_multiplier (scalar) – Default 500
  • mask_array (Boolean Numpy array) –
  • lazy_result (bool) – If True, return a lazy signal. If False, compute the result and return a non-lazy signal. Default True.
  • show_progressbar (bool) –

Examples

>>> s = ps.dummy_data.get_hot_pixel_signal()
>>> s_hot_pixels = s.find_hot_pixels(show_progressbar=False)

Using a mask array

>>> import numpy as np
>>> mask_array = np.zeros((128, 128), dtype=np.bool)
>>> mask_array[:, 100:] = True
>>> s = ps.dummy_data.get_hot_pixel_signal()
>>> s_hot_pixels = s.find_hot_pixels(
...     mask_array=mask_array, show_progressbar=False)

Getting a non-lazy signal as output

>>> s_hot_pixels = s.find_hot_pixels(
...     lazy_result=False, show_progressbar=False)

See also

find_dead_pixels(), correct_bad_pixels()

flip_diffraction_x()[source]

Flip the dataset along the diffraction x-axis.

The function returns a new signal, but the data itself is a view of the original signal. So changing the returned signal will also change the original signal (and visa versa). To avoid changing the orignal signal, use the deepcopy method afterwards, but note that this requires double the amount of memory. See below for an example of this.

Returns:flipped_signal
Return type:PixelatedSTEM signal

Example

>>> s = ps.dummy_data.get_holz_simple_test_signal()
>>> s_flip = s.flip_diffraction_x()

To avoid changing the original object afterwards

>>> s_flip = s.flip_diffraction_x().deepcopy()
flip_diffraction_y()[source]

Flip the dataset along the diffraction y-axis.

The function returns a new signal, but the data itself is a view of the original signal. So changing the returned signal will also change the original signal (and visa versa). To avoid changing the original signal, use the deepcopy method afterwards, but note that this requires double the amount of memory. See below for an example of this.

Returns:flipped_signal
Return type:PixelatedSTEM signal

Example

>>> s = ps.dummy_data.get_holz_simple_test_signal()
>>> s_flip = s.flip_diffraction_y()

To avoid changing the original object afterwards

>>> s_flip = s.flip_diffraction_y().deepcopy()
plot(*args, **kwargs)[source]
radial_integration(centre_x=None, centre_y=None, mask_array=None, parallel=True, show_progressbar=True)[source]

Radially integrates a pixelated STEM diffraction signal.

Parameters:
  • centre_y (centre_x,) – If given as int, all the diffraction patterns will have the same centre position. Each diffraction pattern can also have different centre position, by passing a NumPy array with the same dimensions as the navigation axes. Note: in either case both x and y values must be given. If one is missing, both will be set from the signal (0., 0.) positions. If no values are given, the (0., 0.) positions in the signal will be used.
  • mask_array (Boolean NumPy array, optional) – Mask with the same shape as the signal.
  • parallel (bool, default True) – If True, run the processing on several cores. In most cases this should be True, but for debugging False can be useful.
  • show_progressbar (bool) – Default True
Returns:

Return type:

HyperSpy signal, one less signal dimension than the input signal.

Examples

>>> import pixstem.dummy_data as dd
>>> s = dd.get_holz_simple_test_signal()
>>> s_r = s.radial_integration(centre_x=25, centre_y=25,
...     show_progressbar=False)
>>> s_r.plot()

Using center_of_mass to find bright field disk position

>>> s = dd.get_disk_shift_simple_test_signal()
>>> s_com = s.center_of_mass(threshold=2, show_progressbar=False)
>>> s_r = s.radial_integration(
...     centre_x=s_com.inav[0].data, centre_y=s_com.inav[1].data,
...     show_progressbar=False)
>>> s_r.plot()
rotate_diffraction(angle, parallel=True, show_progressbar=True)[source]

Rotate the diffraction dimensions.

Parameters:
  • angle (scalar) – Clockwise rotation in degrees.
  • parallel (bool) – Default True
  • show_progressbar (bool) – Default True
Returns:

rotated_signal

Return type:

PixelatedSTEM class

Examples

>>> s = ps.dummy_data.get_holz_simple_test_signal()
>>> s_rot = s.rotate_diffraction(30, show_progressbar=False)
shift_diffraction(shift_x, shift_y, parallel=True, inplace=False, show_progressbar=True)[source]

Shift the diffraction patterns in a pixelated STEM signal.

The points outside the boundaries are set to zero.

Parameters:
  • shift_y (shift_x,) – If given as int, all the diffraction patterns will have the same shifts. Each diffraction pattern can also have different shifts, by passing a NumPy array with the same dimensions as the navigation axes.
  • parallel (bool) – If True, run the processing on several cores. In most cases this should be True, but for debugging False can be useful. Default True
  • inplace (bool) – If True (default), the data is replaced by the result. Useful when working with very large datasets, as this avoids doubling the amount of memory needed. If False, a new signal with the results is returned.
  • show_progressbar (bool) – Default True.
Returns:

shifted_signal

Return type:

PixelatedSTEM signal

Examples

>>> s = ps.dummy_data.get_disk_shift_simple_test_signal()
>>> s_c = s.center_of_mass(threshold=3., show_progressbar=False)
>>> s_c -= 25 # To shift the center disk to the middle (25, 25)
>>> s_shift = s.shift_diffraction(
...     s_c.inav[0].data, s_c.inav[1].data,
...     show_progressbar=False)
>>> s_shift.plot()
threshold_and_mask(threshold=None, mask=None, show_progressbar=True)[source]

Get a thresholded and masked of the signal.

Useful for figuring out optimal settings for the center_of_mass method.

Parameters:
  • threshold (number, optional) – The thresholding will be done at mean times this threshold value.
  • mask (tuple (x, y, r)) – Round mask centered on x and y, with radius r.
  • show_progressbar (bool) – Default True
Returns:

s_out

Return type:

PixelatedSTEM signal

Examples

>>> import pixstem.dummy_data as dd
>>> s = dd.get_disk_shift_simple_test_signal()
>>> mask = (25, 25, 10)
>>> s_out = s.threshold_and_mask(
...     mask=mask, threshold=2, show_progressbar=False)
>>> s_out.plot()
virtual_annular_dark_field(cx, cy, r_inner, r, lazy_result=False, show_progressbar=True)[source]

Get a virtual annular dark field signal.

Parameters:
  • cy (cx,) – x- and y-centre positions.
  • r_inner (float) – Inner radius.
  • r (float) – Outer radius.
  • lazy_result (bool, optional) – If True, will not compute the data directly, but return a lazy signal. Default False
  • show_progressbar (bool, default True) –
Returns:

virtual_adf_signal

Return type:

HyperSpy 2D signal

Examples

>>> s = ps.dummy_data.get_holz_heterostructure_test_signal()
>>> s_adf = s.virtual_annular_dark_field(
...     40, 40, 20, 40, show_progressbar=False)
>>> s_adf.plot()

Get a lazy signal, then compute

>>> s_adf = s.virtual_annular_dark_field(
...     40, 40, 20, 40, lazy_result=True, show_progressbar=False)
>>> s_adf.compute(progressbar=False)
>>> s_adf.plot()
virtual_bright_field(cx=None, cy=None, r=None, lazy_result=False, show_progressbar=True)[source]

Get a virtual bright field signal.

Can be sum the whole diffraction plane, or a circle subset. If any of the parameters are None, it will sum the whole diffraction plane.

Parameters:
  • cy (cx,) – x- and y-centre positions.
  • r (float, optional) – Outer radius.
  • lazy_result (bool, optional) – If True, will not compute the data directly, but return a lazy signal. Default False
  • show_progressbar (bool, optional) – Default True.
Returns:

virtual_bf_signal

Return type:

HyperSpy 2D signal

Examples

>>> s = ps.dummy_data.get_holz_heterostructure_test_signal()
>>> s_bf = s.virtual_bright_field(show_progressbar=False)
>>> s_bf.plot()

Sum a subset of the diffraction pattern

>>> s_bf = s.virtual_bright_field(40, 40, 10, show_progressbar=False)
>>> s_bf.plot()

Get a lazy signal, then compute

>>> s_bf = s.virtual_bright_field(
...     lazy_result=True, show_progressbar=False)
>>> s_bf.compute(progressbar=False)

DPCBaseSignal

class pixstem.pixelated_stem_class.DPCBaseSignal(*args, **kwargs)[source]

Signal for processing differential phase contrast (DPC) acquired using scanning transmission electron microscopy (STEM).

The signal assumes the data is 3 dimensions, where the two signal dimensions are the probe positions, and the navigation dimension is the x and y disk shifts.

The first navigation index (s.inav[0]) is assumed to the be x-shift and the second navigation is the y-shift (s.inav[1]).

DPCSignal1D

class pixstem.pixelated_stem_class.DPCSignal1D(*args, **kwargs)[source]

Signal for processing differential phase contrast (DPC) acquired using scanning transmission electron microscopy (STEM).

The signal assumes the data is 2 dimensions, where the signal dimension is the probe position, and the navigation dimension is the x and y disk shifts.

The first navigation index (s.inav[0]) is assumed to the be x-shift and the second navigation is the y-shift (s.inav[1]).

get_bivariate_histogram(histogram_range=None, masked=None, bins=200, spatial_std=3)[source]

Useful for finding the distribution of magnetic vectors(?).

Parameters:
  • histogram_range (tuple, optional) – Set the minimum and maximum of the histogram range. Default is setting it automatically.
  • masked (1-D numpy bool array, optional) – Mask parts of the data. The array must be the same size as the signal. The True values are masked. Default is not masking anything.
  • bins (integer, default 200) – Number of bins in the histogram
  • spatial_std (number, optional) – If histogram_range is not given, this value will be used to set the automatic histogram range. Default value is 3.
Returns:

s_hist

Return type:

Signal2D

DPCSignal2D

class pixstem.pixelated_stem_class.DPCSignal2D(*args, **kw)[source]

Signal for processing differential phase contrast (DPC) acquired using scanning transmission electron microscopy (STEM).

The signal assumes the data is 3 dimensions, where the two signal dimensions are the probe positions, and the navigation dimension is the x and y disk shifts.

The first navigation index (s.inav[0]) is assumed to the be x-shift and the second navigation is the y-shift (s.inav[1]).

correct_ramp(corner_size=0.05, only_offset=False, out=None)[source]

Subtracts a plane from the signal, useful for removing the effects of d-scan in a STEM beam shift dataset.

The plane is calculated by fitting a plane to the corner values of the signal. This will only work well when the property one wants to measure is zero in these corners.

Parameters:
  • corner_size (number, optional) – The size of the corners, as a percentage of the image’s axis. If corner_size is 0.05 (5%), and the image is 500 x 1000, the size of the corners will be (500*0.05) x (1000*0.05) = 25 x 50. Default 0.05
  • only_offset (bool, optional) – If True, will subtract a “flat” plane, i.e. it will subtract the mean value of the corners. Default False
  • out (optional, DPCSignal2D signal) –
Returns:

corrected_signal

Return type:

Signal2D

Examples

>>> s = ps.dummy_data.get_square_dpc_signal(add_ramp=True)
>>> s_corr = s.correct_ramp()
>>> s_corr.plot()

Only correct offset

>>> s_corr = s.correct_ramp(only_offset=True)
>>> s_corr.plot()
flip_axis_90_degrees(flips=1)[source]

Flip both the spatial and beam deflection axis

Will rotate both the image and the beam deflections by 90 degrees.

Parameters:flips (int, default 1) – Number of flips. The default (1) gives 90 degrees rotation. 2 gives 180, 3 gives 270, ...

Examples

>>> s = ps.dummy_data.get_stripe_pattern_dpc_signal()
>>> s
<DPCSignal2D, title: , dimensions: (2|50, 100)>
>>> s_rot = s.flip_axis_90_degrees()
>>> s_rot
<DPCSignal2D, title: , dimensions: (2|100, 50)>

Do several flips

>>> s_rot = s.flip_axis_90_degrees(2)
>>> s_rot
<DPCSignal2D, title: , dimensions: (2|50, 100)>
>>> s_rot = s.flip_axis_90_degrees(3)
>>> s_rot
<DPCSignal2D, title: , dimensions: (2|100, 50)>
gaussian_blur(sigma=2, output=None)[source]

Blur the x- and y-beam shifts.

Useful for reducing the effects of structural diffraction effects.

Parameters:
  • sigma (scalar, default 5) –
  • output (HyperSpy signal) –
Returns:

blurred_signal

Return type:

HyperSpy 2D Signal

Examples

>>> s = ps.dummy_data.get_square_dpc_signal(add_ramp=False)
>>> s_blur = s.gaussian_blur()

Different sigma

>>> s_blur = s.gaussian_blur(sigma=1.2)

Using the signal itself as output

>>> s.gaussian_blur(output=s)
>>> s.plot()
get_bivariate_histogram(histogram_range=None, masked=None, bins=200, spatial_std=3)[source]

Useful for finding the distribution of magnetic vectors(?).

Parameters:
  • histogram_range (tuple, optional) – Set the minimum and maximum of the histogram range. Default is setting it automatically.
  • masked (2-D numpy bool array, optional) – Mask parts of the data. The array must be the same size as the signal. The True values are masked. Default is not masking anything.
  • bins (integer, default 200) – Number of bins in the histogram
  • spatial_std (number, optional) – If histogram_range is not given, this value will be used to set the automatic histogram range. Default value is 3.
Returns:

s_hist

Return type:

HyperSpy Signal2D

Examples

>>> s = ps.dummy_data.get_stripe_pattern_dpc_signal()
>>> s_hist = s.get_bivariate_histogram()
>>> s_hist.plot()
get_color_image_with_indicator(phase_rotation=0, indicator_rotation=0, only_phase=False, autolim=True, autolim_sigma=4, scalebar_size=None, ax=None, ax_indicator=None)[source]

Make a matplotlib figure showing DPC contrast.

Parameters:
  • phase_rotation (float, default 0) – Changes the phase of the plotted data. Useful for correcting scan rotation.
  • indicator_rotation (float, default 0) – Changes the color wheel rotation.
  • only_phase (bool, default False) – If False, will plot both the magnitude and phase. If True, will only plot the phase.
  • autolim (bool, default True) –
  • autolim_sigma (float, default 4) –
  • scalebar_size (int, optional) –
  • ax (Matplotlib subplot, optional) –
  • ax_indicator (Matplotlib subplot, optional) – If None, generate a new subplot for the indicator. If False, do not include an indicator
  • s = ps.dummy_data.get_simple_dpc_signal() (>>>) –
  • fig = s.get_color_image_with_indicator() (>>>) –
  • fig.savefig("simple_dpc_test_signal.png") (>>>) –
  • plotting the phase (Only) –
  • fig = s.get_color_image_with_indicator(only_phase=True) (>>>) –
  • fig.savefig("simple_dpc_test_signal.png")
  • subplot as input (Matplotlib) –
  • import matplotlib.pyplot as plt (>>>) –
  • fig, ax = plt.subplots() (>>>) –
  • ax_indicator = fig.add_subplot(331) (>>>) –
  • fig_return = s.get_color_image_with_indicator( (>>>) –
  • scalebar_size=10, ax=ax, ax_indicator=ax_indicator) (..) –
get_color_signal(rotation=None, autolim=True, autolim_sigma=4)[source]

Get DPC image visualized using continuous color scale.

Converts the x and y beam shifts into an RGB array, showing the magnitude and direction of the beam shifts.

Useful for visualizing magnetic domain structures.

Parameters:
  • rotation (float, optional) – In degrees. Useful for correcting the mismatch between scan direction and diffraction pattern rotation.
  • autolim (bool, default True) –
  • autolim_sigma (float, default 4) –
Returns:

color_signal

Return type:

HyperSpy 2D RGB signal

Examples

>>> s = ps.dummy_data.get_simple_dpc_signal()
>>> s_color = s.get_color_signal()
>>> s_color.plot()

Rotate the beam shift by 30 degrees

>>> s_color = s.get_color_signal(rotation=30)

See also

get_color_signal()
Signal showing both phase and magnitude
get_phase_signal()
Signal showing the phase
get_magnitude_signal(autolim=True, autolim_sigma=4)[source]

Get DPC magnitude image visualized as greyscale.

Converts the x and y beam shifts into a magnitude map, showing the magnitude of the beam shifts.

Useful for visualizing magnetic domain structures.

Parameters:
  • autolim (bool, default True) –
  • autolim_sigma (float, default 4) –
Returns:

magnitude_signal

Return type:

HyperSpy 2D signal

Examples

>>> s = ps.dummy_data.get_simple_dpc_signal()
>>> s_magnitude = s.get_magnitude_signal()
>>> s_magnitude.plot()

See also

get_color_signal()
Signal showing both phase and magnitude
get_phase_signal()
Signal showing the phase
get_phase_signal(rotation=None)[source]

Get DPC phase image visualized using continuous color scale.

Converts the x and y beam shifts into an RGB array, showing the direction of the beam shifts.

Useful for visualizing magnetic domain structures.

Parameters:
  • rotation (float, optional) – In degrees. Useful for correcting the mismatch between scan direction and diffraction pattern rotation.
  • autolim (bool, default True) –
  • autolim_sigma (float, default 4) –
Returns:

phase_signal

Return type:

HyperSpy 2D RGB signal

Examples

>>> s = ps.dummy_data.get_simple_dpc_signal()
>>> s_color = s.get_phase_signal(rotation=20)
>>> s_color.plot()

See also

get_color_signal()
Signal showing both phase and magnitude
get_magnitude_signal()
Signal showing the magnitude
rotate_beam_shifts(angle)[source]

Rotate the beam shift vector.

Parameters:angle (float) – Clockwise rotation in degrees
Returns:shift_rotated_signal
Return type:DPCSignal2D

Example

Rotate beam shifts by 10 degrees clockwise

>>> s = ps.dummy_data.get_simple_dpc_signal()
>>> s_new = s.rotate_beam_shifts(10)
>>> s_new.plot()
rotate_data(angle, reshape=False)[source]

Rotate the scan dimensions by angle.

Parameters:angle (float) – Clockwise rotation in degrees
Returns:rotated_signal
Return type:DPCSignal2D

Example

Rotate data by 10 degrees clockwise

>>> s = ps.dummy_data.get_simple_dpc_signal()
>>> s_rot = s.rotate_data(10)
>>> s_rot.plot()

Modules

IO tools

pixstem.io_tools.load_binary_merlin_signal(filename, probe_x=None, probe_y=None, chunks=(32, 32, 32, 32), flyback_pixels=1, lazy_result=True)[source]

Temporary function for loading Merlin binary data.

This function will be replaced at some point, so do not rely on it for other functions!

pixstem.io_tools.load_dpc_signal(filename)[source]

Load a differential phase contrast style signal.

This function can both files saved directly using HyperSpy, and saved using this library. The only requirement is that the signal has one navigation dimension, with this one dimension having a size of two. The first navigation index is the x-shift, while the second is the y-shift. The signal dimension contains the spatial dimension(s), i.e. the probe positions.

The return signal depends on the dimensions of the input file: - If two signal dimensions: DPCSignal2D - If one signal dimension: DPCSignal1D - If zero signal dimension: DPCBaseSignal

Parameters:filename (string) –
Returns:dpc_signal – The type of return signal depends on the signal dimensions of the input file.
Return type:DPCBaseSignal, DPCSignal1D, DPCSignal2D

Examples

>>> import numpy as np
>>> s = ps.DPCSignal2D(np.random.random((2, 90, 50)))
>>> s.save("test_dpc_signal2d.hspy", overwrite=True)
>>> s_dpc = ps.load_dpc_signal("test_dpc_signal2d.hspy")
>>> s_dpc
<DPCSignal2D, title: , dimensions: (2|50, 90)>
>>> s_dpc.plot()

Saving a HyperSpy signal

>>> import hyperspy.api as hs
>>> s = hs.signals.Signal1D(np.random.random((2, 10)))
>>> s.save("test_dpc_signal1d.hspy", overwrite=True)
>>> s_dpc_1d = ps.load_dpc_signal("test_dpc_signal1d.hspy")
>>> s_dpc_1d
<DPCSignal1D, title: , dimensions: (2|10)>
pixstem.io_tools.load_ps_signal(filename, lazy=False, chunk_size=(16, 16), navigation_signal=None)[source]
Parameters:
  • filename (string) –
  • lazy (bool, default False) –
  • chunk_size (tuple, default (16, 16)) – Used if Lazy is True. Sets the chunk size of the signal in the navigation dimension. Higher number will potentially make the calculations be faster, but use more memory.
  • navigation_signal (Signal2D) –
pixstem.io_tools.signal_to_pixelated_stem(s)[source]

Make a PixelatedSTEM object from a HyperSpy signal.

This will retain both the axes information and the metadata. If the signal is lazy, the function will return LazyPixelatedSTEM.

Parameters:s (HyperSpy signal) – Should work for any HyperSpy signal.
Returns:pixelated_stem_signal
Return type:PixelatedSTEM or LazyPixelatedSTEM object

Examples

>>> import numpy as np
>>> import hyperspy.api as hs
>>> s = hs.signals.Signal2D(np.random.random((8, 11, 21, 13)))
>>> s.metadata.General.title = "test dataset"
>>> s
<Signal2D, title: test dataset, dimensions: (11, 8|13, 21)>
>>> from pixstem.io_tools import signal_to_pixelated_stem
>>> s_new = signal_to_pixelated_stem(s)
>>> s_new
<PixelatedSTEM, title: test dataset, dimensions: (11, 8|13, 21)>

Radial

pixstem.radial.fit_ellipses_to_signal(s, radial_signal_span_list, prepeak_range=None, angleN=20, show_progressbar=True)[source]
Parameters:
  • s (HyperSpy Signal2D) –
  • radial_signal_span_list (tuple) –
  • prepeak_range (tuple, optional) –
  • angleN (list or int, default 20) –
  • show_progressbar (bool, default True) –
Returns:

  • signal (HyperSpy 2D signal) – Fitted ellipse and fitting points are stored as HyperSpy markers in the metadata
  • xC, yC (floats) – Centre position
  • semi_len0, semi_len1 (floats) – Length of the two semi-axes.
  • rot (float) – Angle between semi_len0 and the positive x-axis. Since semi_len0, is not necessarily the longest semi-axis, the rotation will _not_ be between the major semi-axis and the positive x-axis. In radians, between 0 and pi. The rotation is clockwise, so at rot = 0.1 the ellipse will be pointing in the positive x-direction, and negative y-direction.
  • eccen (float) – Eccentricity of the ellipse

Examples

>>> import pixstem.make_diffraction_test_data as mdtd
>>> s = ps.PixelatedSTEM(np.zeros((200, 220)))
>>> s.axes_manager[0].offset, s.axes_manager[1].offset = -100, -110
>>> ellipse_ring = mdtd._get_elliptical_ring(s, 0, 0, 50, 70, 0.8)
>>> s.data += ellipse_ring
>>> from pixstem.radial import fit_ellipses_to_signal
>>> output = fit_ellipses_to_signal(
...     s, [(40, 80)], angleN=30, show_progressbar=False)
pixstem.radial.fit_single_ellipse_to_signal(s, radial_signal_span, prepeak_range=None, angleN=20, show_progressbar=True)[source]
Parameters:
  • s (HyperSpy Signal2D) –
  • radial_signal_span (tuple) –
  • prepeak_range (tuple, optional) –
  • angleN (int, default 20) –
  • show_progressbar (bool, default True) –
Returns:

  • signal (HyperSpy 2D signal) – Fitted ellipse and fitting points are stored as HyperSpy markers in the metadata
  • xC, yC (floats) – Centre position
  • semi_len0, semi_len1 (floats) – Length of minor and major semi-axis
  • rot (float) – Angle between semi_len0 and the positive x-axis. Since semi_len0, is not necessarily the longest semi-axis, the rotation will _not_ be between the major semi-axis and the positive x-axis. In radians, between 0 and pi. The rotation is clockwise, so at rot = 0.1 the ellipse will be pointing in the positive x-direction, and negative y-direction.
  • eccen (float) – Eccentricity of the ellipse
  • signal, xC, yC, semi0, semi1, rot, ecc

Examples

>>> import pixstem.make_diffraction_test_data as mdtd
>>> s = ps.PixelatedSTEM(np.zeros((200, 220)))
>>> s.axes_manager[0].offset, s.axes_manager[1].offset = -100, -110
>>> ellipse_ring = mdtd._get_elliptical_ring(s, 0, 0, 50, 70, 0.8)
>>> s.data += ellipse_ring
>>> from pixstem.radial import fit_single_ellipse_to_signal
>>> output = fit_single_ellipse_to_signal(
...     s, (40, 80), angleN=30, show_progressbar=False)
pixstem.radial.get_angle_image_comparison(s0, s1, angleN=12, mask_radius=None)[source]

Compare two images by overlaying one on the other in slices.

This function takes two images, extracts different angular slices and combines them into one image.

Useful for comparing two diffraction images, to see if the rings have the same radius.

Parameters:
  • s1 (s0,) – Both signals need to have the same shape, and no navigation dimensions.
  • angleN (int, default 12) – Number of angular slices.
  • mask_radius (int, optional) – Mask the centre of the image. The default is not to mask anything. Useful to mask the most intense parts of the diffraction pattern, so the less intense parts can be visualized.
Returns:

comparison_signal

Return type:

HyperSpy 2D

Examples

>>> from pixstem.make_diffraction_test_data import MakeTestData
>>> test_data0 = MakeTestData(300, 300)
>>> test_data0.add_ring(150, 150, 40)
>>> test_data1 = MakeTestData(300, 300)
>>> test_data1.add_ring(150, 150, 60)
>>> s0 = test_data0.signal
>>> s1 = test_data1.signal
>>> s0.axes_manager[0].offset, s0.axes_manager[1].offset = -150, -150
>>> s1.axes_manager[0].offset, s1.axes_manager[1].offset = -150, -150
>>> import pixstem.radial as ra
>>> s = ra.get_angle_image_comparison(s0, s1)
>>> s.plot()

Mask the inner parts

>>> s = ra.get_angle_image_comparison(s0, s1, mask_radius=10)
pixstem.radial.get_centre_position_list(s, steps, step_size)[source]

Returns a zip of x and y coordinates based on the offset of the center- point, and number of steps in each direction and a step-size given in pixels. Scale and offset taken from axes_manager.

pixstem.radial.get_coordinate_of_min(s)[source]

Returns the x and y values of the minimum in a signal.

pixstem.radial.get_optimal_centre_position(s, radial_signal_span, steps=3, step_size=1, angleN=8, show_progressbar=True)[source]

Find centre position of a ring by using angle sliced radial integration.

Takes signal s, radial span of feature used to determine the centre position. Radially integrates the feature in angleN segments, for each possible centre position. Models this integrated signal as a Gaussian and returns array of with the standard deviations of these Gaussians.

Note, the approximate centre position must be set using the offset parameter in the signal. The offset parameter is the negative of the position shown with s.plot(). So if the centre position in the image is x=52 and y=55: s.axes_manager.signal_axes[0].offset = -52 s.axes_manager.signal_axes[1].offset = -55 This can be checked by replotting the signal, to see if the centre position has the values x=0 and y=0.

Parameters:
  • s (HyperSpy 2D signal) – Only supports signals with no navigation dimensions.
  • radial_signal_span (tuple) – Range for finding the circular feature.
  • steps (int, default 3) – Number of steps in x/y direction to look for the optimal centre position. If the offset is (50, 55), and step_size 1: the positions x=(45, 46, 47, 48, ..., 55) and y=(50, 51, 52, ..., 60), will be checked.
  • step_size (int, default 1) –
  • angleN (int, default 8) – See angular_slice_radial_integration for information about angleN.
Returns:

Optimal centre position

Return type:

HyperSpy 2D signal

Examples

>>> import pixstem.dummy_data as dd
>>> s = dd.get_single_ring_diffraction_signal()
>>> s.axes_manager.signal_axes[0].offset = -105
>>> s.axes_manager.signal_axes[1].offset = -67
>>> import pixstem.radial as ra
>>> s_ocp = ra.get_optimal_centre_position(s, (35, 45), steps=2)
>>> centre_pos = ra.get_coordinate_of_min(s_ocp)
pixstem.radial.get_radius_vs_angle(signal, radial_signal_span, angleN=15, centre_x=None, centre_y=None, prepeak_range=None, show_progressbar=True)[source]

Get radius of a ring as a function of angle.

This is done by radially integrating angular slices of the image, followed by firstly fitting the background around the specified ring, then fitting a Gaussian to the remaining intensity in the ring (which gets reduced to a peak due to the radial integration).

Useful for finding if a ring is circular or elliptical. The centre position has be to set using the offset parameter in signal.axes_manager.

Parameters:
  • signal (HyperSpy 2D signal) – Diffraction image with calibrated centre point. This can be set manually, by looking at the diffraction pattern, or use a function like radial.get_optimal_centre_position
  • radial_signal_span (tuple) –
  • angleN (default 15) –
  • centre_y (centre_x,) –
  • prepeak_range (tuple, optional) –
  • show_progressbar (default True) –
Returns:

s_centre

Return type:

HyperSpy 1D signal

Examples

>>> import pixstem.dummy_data as dd
>>> s = dd.get_single_ring_diffraction_signal()
>>> s.axes_manager.signal_axes[0].offset = -105
>>> s.axes_manager.signal_axes[1].offset = -67
>>> import pixstem.radial as ra
>>> s_centre = ra.get_radius_vs_angle(s, (35, 45), show_progressbar=False)
pixstem.radial.refine_signal_centre_position(s, radial_signal_span, refine_step_size=None, **kwargs)[source]

Refine centre position of a diffraction signal by using round feature.

This function simply calls get_optimal_centre_position with a smaller and smaller step_size, giving a more accurate centre position.

The calculated centre position is stored in the offset value in the input signal.

Parameters:
  • s (HyperSpy 2D signal) – Approximate centre position must be set beforehand, see get_optimal_centre_position for more information.
  • radial_signal_span (tuple) –
  • refine_step_size (list, default [1, 0.5, 0.25]) –
  • **kwargs – Passed to get_optimal_centre_position

See also

get_optimal_centre_position()
finds the centre positions for a step_size

Example

>>> import pixstem.dummy_data as dd
>>> import pixstem.radial as ra
>>> s = dd.get_single_ring_diffraction_signal()
>>> s.axes_manager[0].offset, s.axes_manager[1].offset = -103, -69
>>> ra.refine_signal_centre_position(s, (32., 48.), angleN=4)
>>> x, y = s.axes_manager[0].offset, s.axes_manager[1].offset

Dummy data

pixstem.dummy_data.get_dead_pixel_signal(lazy=False)[source]

Get 2D PixelatedSTEM signal with a disk in the middle.

Has 4 pixels with value equal to 0, to simulate dead pixels.

Example

>>> s = ps.dummy_data.get_dead_pixel_signal()

Lazy signal

>>> s_lazy = ps.dummy_data.get_dead_pixel_signal(lazy=True)
pixstem.dummy_data.get_disk_shift_simple_test_signal(lazy=False)[source]

Get HyperSpy 2D signal with 2D navigation dimensions for DPC testing.

Probe size x/y (20, 20), and image size x/y (50, 50). Disk moves from 22-28 x/y.

Parameters:lazy (bool, default False) –
Returns:disk_shift_signal
Return type:PixelatedSTEM signal

Examples

>>> s = ps.dummy_data.get_disk_shift_simple_test_signal()
>>> s.plot()

Load as lazy

>>> s = ps.dummy_data.get_disk_shift_simple_test_signal(lazy=True)
pixstem.dummy_data.get_holz_heterostructure_test_signal(lazy=False)[source]

Get HyperSpy 2D signal with 2D navigation dimensions for HOLZ testing.

The centre, radius and intensity of the ring varies as a function of probe position. The disk centre position varies as a function of probe position.

Parameters:lazy (bool, default False) –
Returns:holz_signal
Return type:HyperSpy 2D signal

Example

>>> s = ps.dummy_data.get_holz_heterostructure_test_signal()
>>> s.plot()

Load as lazy

>>> s = ps.dummy_data.get_holz_heterostructure_test_signal(lazy=True)
pixstem.dummy_data.get_holz_simple_test_signal(lazy=False)[source]

Get HyperSpy 2D signal with 2D navigation dimensions for HOLZ testing.

Probe size x/y (20, 20), and image size x/y (50, 50). Contains a disk and a ring. The disk stays at x, y = 25, 25, with radius 2. The ring has a radius of 20, and moves from x, y = 24-26, 24-26.

Parameters:lazy (bool, default False) –
Returns:holz_signal
Return type:PixelatedSTEM signal

Examples

>>> s = ps.dummy_data.get_holz_simple_test_signal()
>>> s.plot()

Load as lazy

>>> s = ps.dummy_data.get_holz_simple_test_signal(lazy=True)
pixstem.dummy_data.get_hot_pixel_signal(lazy=False)[source]

Get 2D PixelatedSTEM signal with a disk in the middle.

Has 4 pixels with value equal to 50000, to simulate hot pixels.

Example

>>> s = ps.dummy_data.get_hot_pixel_signal()

Lazy signal

>>> s_lazy = ps.dummy_data.get_hot_pixel_signal(lazy=True)
pixstem.dummy_data.get_simple_dpc_signal()[source]

Get a simple DPCSignal2D with a zero point in the centre.

Example

>>> s = ps.dummy_data.get_simple_dpc_signal()
pixstem.dummy_data.get_single_ring_diffraction_signal()[source]

Get HyperSpy 2D signal with a single ring with centre position.

The ring has a centre at x=105 and y=67, and radius=40.

pixstem.dummy_data.get_square_dpc_signal(add_ramp=False)[source]

Get a 2D DPC signal resembling a Landau domain.

Parameters:add_ramp (bool, default False) – If True, will add a ramp in the beam shift across the image to emulate the effects of d-scan.
Returns:square_dpc_signal
Return type:DPCSignal2D

Examples

>>> s = ps.dummy_data.get_square_dpc_signal()
>>> s.plot()

Adding a ramp

>>> s = ps.dummy_data.get_square_dpc_signal(add_ramp=True)
>>> s.plot()
pixstem.dummy_data.get_stripe_pattern_dpc_signal()[source]

Get a 2D DPC signal with a stripe pattern.

The stripe pattern only has an x-component, with alternating left/right directions. There is a small a net moment in the positive x-direction (leftwards).

Returns:stripe_dpc_signal
Return type:DPCSignal2D

Example

>>> s = ps.dummy_data.get_stripe_pattern_dpc_signal()

Make diffraction test data

class pixstem.make_diffraction_test_data.Circle(xx, yy, x0, y0, r, intensity, scale, lw=None)[source]
centre_on_image(xx, yy)[source]
get_centre_pixel(xx, yy, scale)[source]

This function sets the indices for the pixels on which the centre point is. Because the centre point can sometimes be exactly on the boundary of two pixels, the pixels are held in a list. One list for x (self.centre_x_pixels) and one for y (self.centre_x_pixels). If the centre is outside the image, the lists will be empty.

mask_outside_r(scale)[source]
set_uniform_intensity()[source]
update_axis(xx, yy)[source]
class pixstem.make_diffraction_test_data.Disk(xx, yy, scale, x0, y0, r, intensity)[source]

Disk object, with outer edge of the ring at r

get_signal()[source]
set_centre_intensity()[source]

Sets the intensity of the centre pixels to I. Coordinates are self.z.circle[y, x], due to how numpy works.

update_axis(xx, yy)[source]
class pixstem.make_diffraction_test_data.MakeTestData(size_x=100, size_y=100, scale=1, default=False, blur=True, blur_sigma=1, downscale=True)[source]

MakeTestData is an object containing a generated test signal. The default signal is consisting of a Disk and concentric Ring, with the Ring being less intensive than the centre Disk. Unlimited number of Rings and Disks can be added separately.

Parameters:
  • size_y (size_x,) – The range of the x and y axis goes from 0 to size_x, size_y
  • scale (float, int) – The step size of the x and y axis
  • default (bool, default False) – If True, the default object should be generated. If false, Ring and Disk must be added separately by self.add_ring(), self.add_disk()
  • blur (bool, default True) – If True, do a Gaussian blur of the disk.
  • blur_sigma (int, default 1) – Sigma of the Gaussian blurring, if blur is True.
  • downscale (bool, default True) – Note: currently using downscale and adding a disk, will lead to the center of the disk being shifted. Ergo: the disk_x and disk_y will not be correct when downscale is True.
signal

hyperspy.signals.Signal2D – Test signal

z_list

list – List containing Ring and Disk objects added to the signal

downscale_factor

int – The data is upscaled before Circle is added, and similarly downscaled to return to given dimensions. This improves the quality of Circle

Examples

Default settings

>>> from pixstem.make_diffraction_test_data import MakeTestData
>>> test_data = MakeTestData()
>>> test_data.signal.plot()

More control

>>> test_data = MakeTestData(default=False)
>>> test_data.add_disk(x0=50, y0=50, r=10, intensity=30)
>>> test_data.add_ring(x0=45, y0=52, r=25, intensity=10)
>>> test_data.signal.plot()
add_disk(x0=50, y0=50, r=5, intensity=10)[source]
add_ring(x0=50, y0=50, r=20, intensity=10, lw_pix=0)[source]

Add a ring to the test data.

Parameters:
  • y0 (x0,) – Centre position of the ring
  • r (number, default 20) – Radius of the ring, defined as the distance from the centre to the middle of the line of the ring, which will be most intense after blurring.
  • intensity (number, default 10) – Pixel value of the ring. Note, this value will be lowered if blur or downscale is True
  • lw_pix (number, default 0) – Distance in pixels from radius to the outer and inner edge of the ring. Inner radius: r-lw, outer radius: r+lw. In total this gives a ring line width in pixels of 2*lw+1.
blur()[source]
downscale()[source]
generate_mesh()[source]
make_signal()[source]
set_downscale_factor(factor)[source]
set_signal_zero()[source]
to_signal()[source]
update_signal()[source]
class pixstem.make_diffraction_test_data.Ring(xx, yy, scale, x0, y0, r, intensity, lr)[source]

Ring object, with outer edge of the ring at r+lr, and inner r-lr. The radius of the ring is defined as in the middle of the line making up the ring.

get_signal()[source]
mask_inside_r(scale)[source]
update_axis(xx, yy)[source]
pixstem.make_diffraction_test_data.generate_4d_data(probe_size_x=10, probe_size_y=10, image_size_x=50, image_size_y=50, disk_x=25, disk_y=25, disk_r=5, disk_I=20, ring_x=25, ring_y=25, ring_r=20, ring_I=6, ring_lw=0, blur=True, blur_sigma=1, downscale=True, add_noise=False, noise_amplitude=1, lazy=False, lazy_chunks=None)[source]

Generate a test dataset containing a disk and diffraction ring. Useful for checking that radial integration algorithms are working properly.

The centre, intensity and radius position of the ring and disk can vary as a function of probe position, through the disk_x, disk_y, disk_r, disk_I, ring_x, ring_y, ring_r and ring_I arguments. In addition, the line width of the ring can be varied with ring_lw.

The ring can be deactivated by setting ring_x=None. The disk can be deactivated by setting disk_x=None.

Parameters:
  • probe_size_y (probe_size_x,) – Size of the navigation dimension.
  • image_size_y (image_size_x,) – Size of the signal dimension.
  • disk_y (disk_x,) – Centre position of the disk. Either integer or NumPy 2-D array. See examples on how to make them the correct size. To deactivate the ring, set disk_x=None.
  • disk_r (int or NumPy 2D-array, default 5) – Radius of the disk. Either integer or NumPy 2-D array. See examples on how to make it the correct size.
  • disk_I (int or NumPy 2D-array, default 20) – Intensity of the disk, for each of the pixels. So if I=30, the each pixel in the disk will have a value of 30. Note, this value will change if blur=True or downscale=True.
  • ring_y (ring_x,) – Centre position of the ring. Either integer or NumPy 2-D array. See examples on how to make them the correct size. To deactivate the ring, set ring_x=None.
  • ring_r (int or NumPy 2D-array, default 5) – Radius of the ring. Either integer or NumPy 2-D array. See examples on how to make it the correct size.
  • ring_I (int or NumPy 2D-array, default 6) – Intensity of the ring, for each of the pixels. So if I=5, each pixel in the ring will have a value of 5. Note, this value will change if blur=True or downscale=True.
  • ring_lw (int or NumPy 2D-array, default 0) – Line width of the ring. If ring_lw=1, the line will be 3 pixels wide. If ring_lw=2, the line will be 5 pixels wide.
  • blur (bool, default True) – If True, do a Gaussian blur of the disk.
  • blur_sigma (int, default 1) – Sigma of the Gaussian blurring, if blur is True.
  • downscale (bool, default True) – If True, use upscaling (then downscaling) to anti-alise the disk.
  • add_noise (bool, default False) – Add Gaussian random noise.
  • noise_amplitude (float, default 1) – The amplitude of the noise, if add_noise is True.
  • lazy (bool, default False) – If True, the signal will be lazy
  • lazy_chunks (tuple, optional) – Used if lazy is True, default (10, 10, 10, 10).
Returns:

signal – Signal with 2 navigation dimensions and 2 signal dimensions.

Return type:

HyperSpy Signal2D

Examples

>>> import pixstem.make_diffraction_test_data as mdtd
>>> s = mdtd.generate_4d_data()
>>> s.plot()

Using more arguments

>>> s = mdtd.generate_4d_data(probe_size_x=20, probe_size_y=30,
...         image_size_x=50, image_size_y=90,
...         disk_x=30, disk_y=70, disk_r=9, disk_I=30,
...         ring_x=35, ring_y=65, ring_r=20, ring_I=10,
...         blur=False, downscale=False)

Adding some Gaussian random noise

>>> s = mdtd.generate_4d_data(add_noise=True, noise_amplitude=3)

Different centre positions for each probe position. Note the size=(20, 10), and probe_x=10, probe_y=20: size=(y, x).

>>> import numpy as np
>>> disk_x = np.random.randint(5, 35, size=(20, 10))
>>> disk_y = np.random.randint(5, 45, size=(20, 10))
>>> disk_I = np.random.randint(50, 100, size=(20, 10))
>>> ring_x = np.random.randint(5, 35, size=(20, 10))
>>> ring_y = np.random.randint(5, 45, size=(20, 10))
>>> ring_r = np.random.randint(10, 15, size=(20, 10))
>>> ring_I = np.random.randint(1, 30, size=(20, 10))
>>> ring_lw = np.random.randint(1, 5, size=(20, 10))
>>> s = mdtd.generate_4d_data(probe_size_x=10, probe_size_y=20,
...         image_size_x=40, image_size_y=50, disk_x=disk_x, disk_y=disk_y,
...         disk_I=disk_I, ring_x=ring_x, ring_y=ring_y, ring_r=ring_r,
...         ring_I=ring_I, ring_lw=ring_lw)

Do not plot the disk

>>> s = mdtd.generate_4d_data(disk_x=None)

Do not plot the ring

>>> s = mdtd.generate_4d_data(ring_x=None)