API documentation

Classes

Atom lattice

class atomap.atom_lattice.Atom_Lattice(image=None, name='', sublattice_list=None)
Parameters:
  • image (2D NumPy array, optional) –
  • sublattice_list (list of sublattice object, optional) –
image

2D NumPy array

x_position

list of floats – x positions for all sublattices.

y_position

list of floats – y positions for all sublattices.

get_sublattice(sublattice_id)

Get a sublattice object from either sublattice index or name.

get_sublattice_atom_list_on_image(image=None, add_numbers=False, markersize=20)
integrate_column_intensity(method='Voronoi', max_radius='Auto', data_to_integrate=None)

Integrate signal around the atoms in the atom lattice.

See atomap.tools.integrate for more information about the parameters.

Parameters:
  • method (string) – Voronoi or Watershed
  • max_radius (int, optional) –
  • data_to_integrate (NumPy array, HyperSpy signal or array-like) – Works with 2D, 3D and 4D arrays, so for example an EEL spectrum image can be used.
Returns:

Return type:

i_points, i_record, p_record

Examples

>>> import atomap.api as am
>>> al = am.dummy_data.get_simple_atom_lattice_two_sublattices()
>>> i_points, i_record, p_record = al.integrate_column_intensity()

See also

tools.integrate()

plot(image=None, add_numbers=False, markersize=20, **kwargs)

Plot all atom positions for all sub lattices on the image data.

The Atom_Lattice.image0 is used as the image. For the sublattices, sublattice._plot_color is used as marker color. This color is set when the sublattice is initialized, but it can also be changed.

Parameters:
  • image (2D NumPy array, optional) –
  • add_numbers (bool, default False) – Plot the number of the atom beside each atomic position in the plot. Useful for locating misfitted atoms.
  • markersize (number, default 20) – Size of the atom position markers
  • **kwargs – Addition keywords passed to HyperSpy’s signal plot function.

Examples

>>> import atomap.api as am
>>> import atomap.testing_tools as tt
>>> test_data = tt.MakeTestData(50, 50)
>>> import numpy as np
>>> test_data.add_atom_list(np.arange(5, 45, 5), np.arange(5, 45, 5))
>>> atom_lattice = test_data.atom_lattice
>>> atom_lattice.plot()

Change sublattice colour and color map

>>> atom_lattice.sublattice_list[0]._plot_color = 'green'
>>> atom_lattice.plot(cmap='viridis')

See also

get_sublattice_atom_list_on_image()
get HyperSpy signal with atom positions as markers. More customizability.
save(filename=None, overwrite=False)

Save the AtomLattice object as an HDF5 file. This will store the information about the individual atomic positions, like the pixel position, sigma and rotation. The input image(s) and modified version of these will also be saved.

Parameters:
  • filename (string, optional) –
  • overwrite (bool, default False) –

Examples

>>> from numpy.random import random
>>> import atomap.api as am
>>> sl = am.dummy_data.get_simple_cubic_sublattice()
>>> atom_lattice = am.Atom_Lattice(random((9, 9)), "test", [sl])
>>> atom_lattice.save("test.hdf5", overwrite=True)

Loading the atom lattice:

>>> import atomap.api as am
>>> atom_lattice1 = am.load_atom_lattice_from_hdf5("test.hdf5")
signal
x_position
y_position

Sublattice

class atomap.sublattice.Sublattice(atom_position_list, image, original_image=None, name='', color='red', pixel_size=1.0)
Parameters:
  • atom_position_list (NumPy array) – In the form [[x0, y0], [x1, y1], [x2, y2], ... ]
  • image (2D NumPy array) – A HyperSpy signal with 2 dimensions can also be used directly.
  • original_image (2D NumPy array, optional) –
  • name (string, default '') –
  • color (string, optional) – Plotting color, default red.
  • pixel_size (float, optional) – Scaling number, default 1.
image

2D NumPy array, or 2D array-like.

x_position

list of floats

y_position

list of floats

sigma_x

list of floats

sigma_y

list of floats

sigma_average

list of floats

rotation

list of floats – In radians. The rotation of the axes of each 2D-Gaussian relative to the horizontal axes. For the rotation of the ellipticity, see rotation_ellipticity.

ellipticity

list of floats

rotation_ellipticity

list of floats – In radians, the rotation between the “x-axis” and the major axis of the ellipse. Basically giving the direction of the ellipticity.

signal

HyperSpy 2D signal

name

string

Examples

>>> import numpy as np
>>> from atomap.sublattice import Sublattice
>>> atom_positions = [[2, 2], [2, 4], [4, 2], [4, 4]]
>>> image_data = np.random.random((7, 7))
>>> sublattice = Sublattice(atom_positions, image_data)
>>> s_sublattice = sublattice.get_atom_list_on_image()
>>> s_sublattice.plot()

More atom positions

>>> x, y = np.mgrid[0:100:10j, 0:100:10j]
>>> x, y = x.flatten(), y.flatten()
>>> atom_positions = np.dstack((x, y))[0]
>>> image_data = np.random.random((100, 100))
>>> sublattice = Sublattice(atom_positions, image_data, color='yellow',
...     name='the sublattice')
>>> sublattice.get_atom_list_on_image(markersize=50).plot()
atom_amplitude_gaussian2d
atom_amplitude_max_intensity
atom_positions
construct_zone_axes(atom_plane_tolerance=0.5, zone_axis_para_list=False)

Constructs the zone axes for an atomic resolution image.

The zone axes are found by finding the 15 nearest neighbors for each atom position in the sublattice, and finding major translation symmetries among the nearest neighbours. Only unique zone axes are kept, and “bad” ones are removed.

After finding the zone axes, atom planes are constructed.

Parameters:
  • atom_plane_tolerance (scalar, default 0.5) – When constructing the atomic planes, the method will try to locate the atoms by “jumping” one zone vector, and seeing if there is an atom with the pixel_separation times atom_plane_tolerance. So this value should be increased the atomic planes are non-continuous and “split”.
  • zone_axis_para_list (parameter list or bool, default False) – A zone axes parameter list is used to name and index the zone axes. See atomap.process_parameters for more info. Useful for automation.

Example

>>> import atomap.api as am
>>> sublattice = am.dummy_data.get_simple_cubic_sublattice()
>>> sublattice
<Sublattice,  (atoms:400,planes:0)>
>>> sublattice.construct_zone_axes()
>>> sublattice
<Sublattice,  (atoms:400,planes:4)>

See also

atom_finding_refining.construct_zone_axes_from_sublattice()

ellipticity
find_missing_atoms_from_zone_vector(zone_vector)

Returns a list of coordinates between atoms given by a zone vector.

These coordinates are given by the mid-point between adjacent atoms in the atom planes with the given zone_vector.

Parameters:zone_vector (tuple) – Zone vector for the atom planes where the new atoms are positioned between the atoms in the sublattice.
Returns:xy_coordinates
Return type:List of tuples

Example

>>> import atomap.api as am
>>> sublattice_A = am.dummy_data.get_simple_cubic_sublattice()
>>> sublattice_A.construct_zone_axes()
>>> zone_axis = sublattice_A.zones_axis_average_distances[0]
>>> B_pos = sublattice_A.find_missing_atoms_from_zone_vector(
...                       zone_axis)
find_nearest_neighbors(nearest_neighbors=9, leafsize=100)

Find the nearest neighbors for all the atom position objects. Needs to be run before doing any position refinements.

Parameters:
  • nearest_neighbors (int, default 9) –
  • leafsize (int, default 100) –
find_sublattice_intensity_from_masked_image(image_data=None, radius=2)

Find the image intensity of each atomic column in the sublattice.

The intensity of the atomic column is given by the average intensity of the pixels inside an area within a radius in pixels from each atom position.

Parameters:
  • image_data (2-D NumPy array, optional) – Image data for plotting. If none is given, will use the original_image.
  • radius (int, optional, default 2) – The radius in pixels away from the atom centre pixels, determining the area that shall not be masked. The default radius is 3 pixels.

Examples

>>> import atomap.api as am
>>> sublattice = am.dummy_data.get_simple_cubic_sublattice()
>>> sublattice.find_sublattice_intensity_from_masked_image(
...     sublattice.image, 3)
>>> intensity = sublattice.intensity_mask
get_all_atom_planes_by_zone_vector(zone_vector_list=None, image=None, add_numbers=True, color='red')

Get a overview of atomic planes for some or all zone vectors.

Parameters:
  • zone_vector_list (optional) – List of zone vectors for visualizing atomic planes. Default is visualizing all the zone vectors.
  • image (2D Array, optional) –
  • add_numbers (bool, optional, default True) – If True, will the number of the atom plane at the end of the atom plane line. Useful for finding the index of the atom plane.
  • color (string, optional, default red) – The color of the lines and text used to show the atom planes.
Returns:

  • HyperSpy signal2D object if given a single zone vector,
  • list of HyperSpy signal2D if given a list (or none) of zone vectors.

Examples

Getting a list signals showing the atomic planes for all the zone vectors

>>> import atomap.api as am
>>> sublattice = am.dummy_data.get_simple_cubic_sublattice()
>>> sublattice.construct_zone_axes()
>>> s = sublattice.get_all_atom_planes_by_zone_vector()
>>> s.plot()

Several zone vectors

>>> zone_vec_list = sublattice.zones_axis_average_distances[0:3]
>>> s = sublattice.get_all_atom_planes_by_zone_vector(zone_vec_list)
>>> s.plot()

Different image

>>> from numpy.random import random
>>> im = random((9, 9))
>>> s = sublattice.get_all_atom_planes_by_zone_vector(image=im)
>>> s.plot()
get_atom_angles_from_zone_vector(zone_vector0, zone_vector1, degrees=False)

Calculates for each atom in the sub lattice the angle between the atom, and the next atom in the atom planes in zone_vector0 and zone_vector1. Default will return the angles in radians.

Parameters:
  • zone_vector0 (tuple) – Vector for the first zone.
  • zone_vector1 (tuple) – Vector for the second zone.
  • degrees (bool, optional) – If True, will return the angles in degrees. Default False.
get_atom_column_amplitude_max_intensity(image=None, percent_to_nn=0.4)

Finds the maximal intensity for each atomic column.

Finds the maximal image intensity of each atomic column inside an area covering the atomic column.

Results are stored in each Atom_Position object as amplitude_max_intensity, which can most easily be accessed in through the sublattice object (see the examples below).

Parameters:
  • image (NumPy 2D array, default None) – Uses original_image by default.
  • percent_to_nn (float, default 0.4) – Determines the boundary of the area surrounding each atomic column, as fraction of the distance to the nearest neighbour.

Example

>>> import atomap.api as am
>>> sublattice = am.dummy_data.get_simple_cubic_sublattice()
>>> sublattice.find_nearest_neighbors()
>>> sublattice.get_atom_column_amplitude_max_intensity()
>>> intensity_list = sublattice.atom_amplitude_max_intensity
get_atom_distance_difference_from_zone_vector(zone_vector)

Get distance difference between atoms in atoms planes belonging to a zone axis.

Parameters:zone_vector (tuple) – Zone vector for the system
get_atom_distance_difference_line_profile(zone_vector, atom_plane, invert_line_profile=False, interpolate_value=50)
get_atom_distance_difference_map(zone_vector_list=None, atom_plane_list=None, prune_outer_values=False, invert_line_profile=False, add_zero_value_sublattice=None, upscale_map=2)
get_atom_distance_list_from_zone_vector(zone_vector)

Get distance between each atom and the next atom in an atom plane given by the zone_vector. Returns the x- and y-position, and the distance between the atom and the monolayer. The position is set between the atom and the monolayer.

For certain zone axes where there is several monolayers between the atoms in the atom plane, there will be some artifacts created. For example, in the perovskite (110) projection, the atoms in the <111> atom planes are separated by 3 monolayers.

To avoid this problem use the function get_monolayer_distance_list_from_zone_vector.

Parameters:zone_vector (tuple) – Zone vector for the system
get_atom_distance_map(zone_vector_list=None, atom_plane_list=None, prune_outer_values=False, invert_line_profile=False, add_zero_value_sublattice=None, upscale_map=2)
get_atom_list_atom_amplitude_gauss2d_range(amplitude_range)
get_atom_list_between_four_atom_planes(par_atom_plane1, par_atom_plane2, ort_atom_plane1, ort_atom_plane2)

Get a slice of atoms between two pairs of atom planes. Each pair must belong to the zone vector.

Parameters:
  • par_atom_plane2 (par_atom_plane1,) –
  • ort_atom_plane2 (ort_atom_plane1,) –
Returns:

atom_list

Return type:

list of Atom_Position objects

get_atom_list_on_image(atom_list=None, image=None, color=None, add_numbers=False, markersize=20)

Plot atom positions on the image data.

Parameters:
  • atom_list (list of Atom objects, optional) – Atom positions to plot. If no list is given, will use the atom_list.
  • image (2-D NumPy array, optional) – Image data for plotting. If none is given, will use the original_image.
  • color (string, optional) –
  • add_numbers (bool, default False) – Plot the number of the atom beside each atomic position in the plot. Useful for locating misfitted atoms.
  • markersize (number, default 20) – Size of the atom position markers
Returns:

The atom positions as permanent markers stored in the metadata.

Return type:

HyperSpy 2D-signal

Examples

>>> import atomap.api as am
>>> sublattice = am.dummy_data.get_simple_cubic_sublattice()
>>> sublattice.construct_zone_axes()
>>> s = sublattice.get_atom_list_on_image()
>>> s.plot()

Number all the atoms

>>> s = sublattice.get_atom_list_on_image(add_numbers=True)
>>> s.plot()

Plot a subset of the atom positions

>>> atoms = sublattice.atom_list[0:20]
>>> s = sublattice.get_atom_list_on_image(
...     atom_list=atoms, add_numbers=True)
>>> s.plot(cmap='viridis')

Saving the signal as HyperSpy HDF5 file, which saves the atom positions as permanent markers.

>>> s = sublattice.get_atom_list_on_image()
>>> s.save("sublattice_atom_positions.hdf5", overwrite=True)
get_atom_plane_slice_between_two_planes(atom_plane1, atom_plane2)

Get a list of atom planes between two atom planes. Both these atom planes must belong to the same zone vector.

The list will include the two atom planes passed to the function.

Parameters:atom_plane2 (atom_plane1,) –
Returns:atom_plane_slice
Return type:list of Atom Plane objects
get_atom_planes_on_image(atom_plane_list, image=None, add_numbers=True, color='red')

Get atom_planes signal as lines on the image.

Parameters:
  • atom_plane_list (list of atom_plane objects) – atom_planes to be plotted on the image contained in
  • image (2D Array, optional) –
  • add_numbers (bool, optional, default True) – If True, will the number of the atom plane at the end of the atom plane line. Useful for finding the index of the atom plane.
  • color (string, optional, default red) – The color of the lines and text used to show the atom planes.
Returns:

Return type:

HyperSpy signal2D object

Examples

>>> import atomap.api as am
>>> sublattice = am.dummy_data.get_simple_cubic_sublattice()
>>> sublattice.construct_zone_axes()
>>> zone_vector = sublattice.zones_axis_average_distances[0]
>>> atom_planes = sublattice.atom_planes_by_zone_vector[zone_vector]
>>> s = sublattice.get_atom_planes_on_image(atom_planes)
>>> s.plot()
get_ellipticity_line_profile(atom_plane, invert_line_profile=False, interpolate_value=50)

Returns a line profile of the atoms ellipticity.

This gives the ellipticity as a function of the distance from a given atom plane (interface).

Parameters:
  • atom_plane (Atomap AtomPlane object) – The plane which is defined as the 0-point in the spatial dimension.
  • invert_line_profile (bool, optional, default False) – Passed to _get_property_line_profile(). If True, will invert the x-axis values.
  • interpolate_value (int, default 50) – Passed to _get_property_line_profile(). The amount of data points between in monolayer, due to HyperSpy signals not supporting non-linear axes.
Returns:

signal

Return type:

HyperSpy Signal1D

Example

>>> import atomap.api as am
>>> sublattice = am.dummy_data.get_simple_cubic_sublattice()
>>> sublattice.construct_zone_axes()
>>> zone = sublattice.zones_axis_average_distances[1]
>>> plane = sublattice.atom_planes_by_zone_vector[zone][4]
>>> s_elli_line = sublattice.get_ellipticity_line_profile(plane)
get_ellipticity_map(upscale_map=2.0, atom_plane_list=None)

Get a HyperSpy signal showing the magnitude of the ellipticity for the sublattice.

Parameters:
  • upscale_map (number, default 2.0) – Amount of upscaling compared to the original image given to Atomap. Note, a high value here can greatly increase the memory use for large images.
  • atom_plane_list (List of Atomap AtomPlane object, optional) – If a list of AtomPlanes are given, the plane positions will be added to the signal as permanent markers. Which can be visualized using s.plot(plot_markers=True). Useful for showing the location of for example an interface.
Returns:

Return type:

HyperSpy 2D signal

Examples

>>> from numpy.random import random
>>> import atomap.api as am
>>> sublattice = am.dummy_data.get_simple_cubic_sublattice()
>>> for atom in sublattice.atom_list:
...     atom.sigma_x, atom.sigma_y = 0.5*random()+1, 0.5*random()+1
>>> s_elli = sublattice.get_ellipticity_map()
>>> s_elli.plot()

Include an atom plane, which is added to the signal as a marker

>>> sublattice.construct_zone_axes()
>>> atom_plane = [sublattice.atom_plane_list[10]]
>>> s_elli = sublattice.get_ellipticity_map(atom_plane_list=atom_plane)
>>> s_elli.plot()
get_ellipticity_vector(image=None, atom_plane_list=None, vector_scale=1.0, color='red')

Visualize the ellipticity and direction of the atomic columns using markers in a HyperSpy signal.

Parameters:
  • image (2-D NumPy array, optional) –
  • atom_plane_list (list of AtomPlane instances) –
  • vector_scale (scaling of the vector) –
  • color (string) –
Returns:

  • HyperSpy 2D-signal with the ellipticity vectors as
  • permanent markers

Examples

>>> import atomap.api as am
>>> sublattice = am.dummy_data.get_simple_cubic_sublattice()
>>> for atom in sublattice.atom_list:
...     atom.sigma_x, atom.sigma_y = 1., 1.2
>>> s = sublattice.get_ellipticity_vector(vector_scale=20)
>>> s.plot()
get_fingerprint_1d(pixel_radius=100)

Produce a distance fingerprint of the sublattice.

Example

>>> import atomap.api as am
>>> sublattice = am.dummy_data.get_simple_cubic_sublattice()
>>> fp = sublattice.get_fingerprint_1d()
>>> import matplotlib.pyplot as plt
>>> cax = plt.plot(fp, marker='o')
get_fingerprint_2d(pixel_radius=100)

Produce a distance and direction fingerprint of the sublattice.

Parameters:pixel_radius (int, default 100) –
Returns:cluster_array
Return type:NumPy array

Example

>>> import atomap.api as am
>>> sublattice = am.dummy_data.get_simple_cubic_sublattice()
>>> fp = sublattice.get_fingerprint_2d()
>>> import matplotlib.pyplot as plt
>>> cax = plt.scatter(fp[:,0], fp[:,1], marker='o')
get_model_image(image_shape=None, show_progressbar=True)

Generate an image of the atomic positions from the atom positions Gaussian parameter’s.

Parameters:
  • image_shape (tuple (x, y), optional) – Size of the image generated. Note that atoms might be outside the image if (x, y) is too small.
  • show_progressbar (default True) –
Returns:

Return type:

signal, HyperSpy 2D signal

get_monolayer_distance_line_profile(zone_vector, atom_plane, invert_line_profile=False, interpolate_value=50)

Finds the distance between each atom and the next monolayer, and the distance to atom_plane (the 0-point). The monolayers belong to the zone vector zone_vector. For more information on the definition of monolayer distance, check sublattice.get_monolayer_distance_list_from_zone_vector()

Parameters:
  • zone_vector (tuple) – Zone vector for the monolayers for which the separation will be found
  • atom_plane (Atomap AtomPlane object) – Passed to get_monolayer_distance_list_from_zone_vector(). 0-point in the line profile.
  • invert_line_profile (bool, optional, default False) – Passed to _get_property_line_profile(). If True, will invert the x-axis values.
  • interpolate_value (int, default 50) – Passed to _get_property_line_profile(). The amount of data points between in monolayer, due to HyperSpy signals not supporting non-linear axes.
Returns:

Return type:

HyperSpy signal1D

Example

>>> from numpy.random import random
>>> import atomap.api as am
>>> sublattice = am.dummy_data.get_simple_cubic_sublattice()
>>> for atom in sublattice.atom_list:
...     atom.sigma_x, atom.sigma_y = 0.5*random()+1, 0.5*random()+1
>>> sublattice.construct_zone_axes()
>>> zone = sublattice.zones_axis_average_distances[0]
>>> plane = sublattice.atom_planes_by_zone_vector[zone][0]
>>> s = sublattice.get_monolayer_distance_line_profile(zone,plane)
get_monolayer_distance_list_from_zone_vector(zone_vector)

Get distance between each atom and the next monolayer given by the zone_vector. Returns the x- and y-position, and the distance between the atom and the monolayer. The position is set between the atom and the monolayer.

The reason for finding the distance between monolayer, instead of directly between the atoms is due to some zone axis having having a rather large distance between the atoms in one atom plane. For example, in the perovskite (110) projection, the atoms in the <111> atom planes are separated by 3 monolayers. This can give very bad results.

To get the distance between atoms, use the function get_atom_distance_list_from_zone_vector.

Parameters:zone_vector (tuple) – Zone vector for the system
get_monolayer_distance_map(zone_vector_list=None, atom_plane_list=None, upscale_map=2)
get_nearest_neighbor_directions(pixel_scale=True, neighbors=None)

Get the vector to the nearest neighbors for the atoms in the sublattice. Giving information similar to a FFT of the image, but for real space.

Useful for seeing if the peakfinding and symmetry finder worked correctly. Potentially useful for doing structure fingerprinting.

Parameters:
  • pixel_scale (bool, optional. Default True) – If True, will return coordinates in pixel scale. If False, will return in data scale (pixel_size).
  • neighbors (int, optional) – The number of neighbors returned for each atoms. If no number is given, will return all the neighbors, which is typically 9 for each atom. As given when running the symmetry finder.
Returns:

Position – (x_position, y_position). Where both are NumPy arrays.

Return type:

tuple

Examples

>>> import atomap.api as am
>>> sublattice = am.dummy_data.get_simple_cubic_sublattice()
>>> sublattice.construct_zone_axes()
>>> x_pos, y_pos = sublattice.get_nearest_neighbor_directions()
>>> import matplotlib.pyplot as plt
>>> cax = plt.scatter(x_pos, y_pos)

With all the keywords

>>> x_pos, y_pos = sublattice.get_nearest_neighbor_directions(
...     pixel_scale=False, neighbors=3)
get_nearest_neighbor_directions_all()

Like get_nearest_neighbour_directions(), but considers all other atoms (instead of the typical 9) as neighbors from each atom.

This method also does not require atoms to have the atom.nearest_neighbor_list parameter populated with sublattice.find_nearest_neighbors().

Without the constraint of looking at only n nearest neighbours, blazing fast internal NumPy functions can be utilized to calculate directions. However, memory usage will grow quadratically with the number of atomic columns. E.g.: 1000 atomic columns will require ~8MB of memory. 10,000 atomic columns will require ~800MB of memory. 100,000 atomic columns will throw a MemoryError exception on most machines.

Returns:nn_position_array – In the form [x_pos, y_pos].
Return type:NumPy array

Examples

>>> import numpy as np
>>> import atomap.api as am
>>> sublattice = am.dummy_data.get_simple_cubic_sublattice()
>>> x_pos, y_pos = sublattice.get_nearest_neighbor_directions_all()
>>> mask = np.sqrt(x_pos**2 + y_pos**2) < 3
>>> import matplotlib.pyplot as plt
>>> cax = plt.scatter(x_pos[mask], y_pos[mask])
get_position_history(image=None, color='red', add_numbers=False, markersize=20, show_progressbar=True)

Plot position history of each atom positions on the image data.

Parameters:
  • image (2-D NumPy array, optional) – Image data for plotting. If none is given, will use the original_image.
  • color (string, default 'red') –
  • add_numbers (bool, default False) – Plot the number of the atom beside each atomic position in the plot. Useful for locating misfitted atoms.
  • markersize (number, default 20) – Size of the atom position markers
  • show_progressbar (default True) –
Returns:

The atom positions as permanent markers stored in the metadata.

Return type:

HyperSpy 2D-signal

get_property_map(x_list, y_list, z_list, atom_plane_list=None, add_zero_value_sublattice=None, upscale_map=2)

Returns an interpolated signal map of a property.

Maps the property in z_list. The map is interpolated, and the upscale factor of the interpolation can be set (default 2).

Parameters:
  • y_list (x_list,) – Lists of x and y positions
  • z_list (list of numbers) – List of properties for positions. z[0] is the property of the position at x[0] and y[0]. z_list must have the same length as x_list and y_list.
  • atom_plane_list (list of Atomap Atom_Plane objects, default None) – List of atom planes which will be added as markers in the signal.
  • add_zero_value_sublattice (Sublattice) – The sublattice for A-cations in a perovskite oxide can be included here, when maps of oxygen tilt patterns are made. The principle is that another sublattice can given as a parameter, in order to add zero-value points in the map. This means that he atom positions in this sublattice will be included in the map, and the property at these coordinate is set to 0. This helps in the visualization of tilt patterns. For example, the “tilt” property” outside oxygen octahedron in perovskite oxides is 0. The A-cations are outside the octahedra.
  • upscale_map (int, default 2) – Upscaling factor for the interpolated map.
Returns:

interpolated_signal

Return type:

HyperSpy Signal2D

Examples

Example with ellipticity as property.

>>> import atomap.api as am
>>> sublattice = am.dummy_data.get_simple_cubic_sublattice()
>>> s = sublattice.get_property_map(
...                    sublattice.x_position,
...                    sublattice.y_position,
...                    sublattice.ellipticity)
>>> s.plot()
get_zone_vector_index(zone_vector_id)

Find zone vector index from zone vector name

get_zone_vector_mean_angle(zone_vector)

Get the mean angle between the atoms planes with a specific zone vector and the horizontal axis. For each atom plane the angle between all the atoms, its neighbor and horizontal axis is calculated. The mean of these angles for all the atom planes is returned.

integrate_column_intensity(method='Voronoi', max_radius='Auto', data_to_integrate=None)

Integrate signal around the atoms in the sublattice

If the sublattice is a part of an Atom_Lattice object, this function will only take into account the atoms contained in this sublattice. To get the integrated intensity for all the sublattices, use the integrate_column_intensity in the Atom_Lattice object.

See atomap.tools.integrate for more information about the parameters.

Parameters:
  • method (string) – Voronoi or Watershed
  • max_radius (int, optional) –
  • data_to_integrate (NumPy array, HyperSpy signal or array-like) – Works with 2D, 3D and 4D arrays, so for example an EEL spectrum image can be used.
Returns:

Return type:

i_points, i_record, p_record

Examples

>>> import atomap.api as am
>>> sl = am.dummy_data.get_simple_cubic_sublattice(image_noise=True)
>>> i_points, i_record, p_record = sl.integrate_column_intensity()

See also

tools.integrate()

intensity_mask
mask_image_around_sublattice(image_data=None, radius=2)

Returns a HyperSpy signal containing a masked image. The mask covers the area around each atom position in the sublattice, from a given radius away from the center position of the atom. This radius is given in pixels.

Parameters:
  • image_data (2-D NumPy array, optional) – Image data for plotting. If none is given, will use the original_image.
  • radius (int, optional) – The radius in pixels away from the atom centre pixels, determining the area that shall not be masked. The default radius is 2 pixels.
Returns:

masked_signal

Return type:

HyperSpy 2D-signal

Examples

>>> import atomap.api as am
>>> sublattice = am.dummy_data.get_simple_cubic_sublattice()
>>> s = sublattice.mask_image_around_sublattice(radius=3)
>>> s.plot()
plot(color=None, add_numbers=False, markersize=20, **kwargs)

Plot all atom positions in the sublattice on the image data.

The sublattice.original_image attribute is used as the image.

Parameters:
  • color (string, optional) – Color of the atom positions. If none is specific the value set in sublattice._plot_color is used.
  • add_numbers (bool, default False) – Plot the number of the atom beside each atomic position in the plot. Useful for locating misfitted atoms.
  • markersize (number, default 20) – Size of the atom position markers
  • **kwargs – Addition keywords passed to HyperSpy’s signal plot function.

Examples

>>> import atomap.api as am
>>> sublattice = am.dummy_data.get_simple_cubic_sublattice()
>>> sublattice.plot()

Setting color and color map

>>> sublattice.plot(color='green', cmap='viridis')

See also

get_atom_list_on_image()
get HyperSpy signal with atom positions as markers. More customizability.
plot_ellipticity_map(**kwargs)

Plot the magnitude of the ellipticity.

Parameters:**kwargs – Addition keywords passed to HyperSpy’s signal plot function.

Examples

>>> import atomap.api as am
>>> sublattice = am.dummy_data.get_fantasite_sublattice()
>>> sublattice.construct_zone_axes()
>>> sublattice.refine_atom_positions_using_2d_gaussian()
>>> sublattice.plot_ellipticity_map()
plot_ellipticity_vectors(save=False)

Get a quiver plot of the rotation and ellipticity for the sublattice.

sublattice.refine_atom_positions_using_2d_gaussian has to be run at least once before using this function. If the sublattice hasn’t been refined with 2D-Gaussians, the value for ellipticity and rotation are default, 1 and 0 respectively. When sigma_x and sigma_y are equal (circle) the ellipticity is 1. To better visualize changes in ellipticity, the 0-point for ellipticity in the plot is set to circular atomic columns.

Parameters:save (bool) – If true, the figure is saved as ‘vector_field.png’.

Examples

>>> import atomap.api as am
>>> sublattice = am.dummy_data.get_fantasite_sublattice()
>>> sublattice.construct_zone_axes()
>>> sublattice.refine_atom_positions_using_2d_gaussian()
>>> sublattice.plot_ellipticity_vectors()
plot_planes(image=None, add_numbers=True, color='red', **kwargs)

Show the atomic planes for all zone vectors.

Parameters:
  • image (2D Array, optional) – If image None, the original image of the sublattice is used as background image.
  • add_numbers (bool, optional, default True) – If True, will the number of the atom plane at the end of the atom plane line. Useful for finding the index of the atom plane.
  • color (string, optional, default red) – The color of the lines and text used to show the atom planes.
  • **kwargs – Additional keywords passed to HyperSpy’s signal plot function.

Examples

>>> import atomap.api as am
>>> sublattice = am.dummy_data.get_simple_cubic_sublattice()
>>> sublattice.construct_zone_axes()
>>> sublattice.plot_planes()
refine_atom_positions_using_2d_gaussian(image_data=None, percent_to_nn=0.4, rotation_enabled=True, show_progressbar=True)

Use 2D-Gaussian fitting to refine the atom positions on the image data.

Parameters:
  • image_data (2D NumPy array, optional) – If this is not specified, original_image will be used.
  • percent_to_nn (float, default 0.4) – Percent to nearest neighbor. The function will find the closest nearest neighbor to the current atom position, and this value times percent_to_nn will be the radius of the mask centered on the atom position. Value should be somewhere between 0.01 (1%) and 1 (100%).
  • rotation_enabled (bool, default True) – If True, rotation will be enabled for the 2D Gaussians. This will most likely make the fitting better, but potentially at a cost of less robust fitting.
  • show_progressbar (default True) –

Example

>>> import atomap.api as am
>>> sublattice = am.dummy_data.get_simple_cubic_sublattice()
>>> sublattice.find_nearest_neighbors()
>>> sublattice.refine_atom_positions_using_2d_gaussian()

See also

refine_atom_positions_using_center_of_mass()

refine_atom_positions_using_center_of_mass(image_data=None, percent_to_nn=0.25, show_progressbar=True)

Use center of mass to refine the atom positions on the image data.

Parameters:
  • image_data (2D NumPy array, optional) – If this is not specified, original_image will be used.
  • percent_to_nn (float, default 0.25) – Percent to nearest neighbor. The function will find the closest nearest neighbor to the current atom position, and this value times percent_to_nn will be the radius of the mask centered on the atom position. Value should be somewhere between 0.01 (1%) and 1 (100%).
  • show_progressbar (default True) –

Example

>>> import atomap.api as am
>>> sublattice = am.dummy_data.get_simple_cubic_sublattice()
>>> sublattice.find_nearest_neighbors()
>>> sublattice.refine_atom_positions_using_center_of_mass()

See also

refine_atom_positions_using_2d_gaussian()

rotation
rotation_ellipticity
save_map_from_datalist(data_list, data_scale, atom_plane=None, dtype='float32', signal_name='datalist_map.hdf5')

data_list : numpy array, 4D

sigma_average
sigma_x
sigma_y
signal
x_position
y_position

Atom plane

class atomap.atom_plane.Atom_Plane(atom_list, zone_vector, atom_lattice)
Parameters:
  • atom_list (list of Atom_Position objects) –
  • zone_vector (tuple) –
  • atom_lattice (Atomap Atom_Lattice object) –
x_position

list of floats

y_position

list of floats

sigma_x

list of floats

sigma_y

list of floats

sigma_average

list of floats

rotation

list of floats

ellipticity

list of floats

amplitude_gaussian
amplitude_max_intensity
ellipticity
get_angle_to_horizontal_axis()

Get angle between atoms in the atom plane and horizontal axis.

get_atom_distance_list()
get_atom_index(check_atom)
get_closest_distance_and_angle_to_point(points_x, points_y, use_precalculated_line=False, plot_debug=False)

Return the smallest distances from each point in a list to the atom plane.

Parameters:
  • points_y (points_x,) – The x and y coordinates for the atoms.
  • use_precalculated_line (bool or list/array) – The coefficients [a, b] for the linear line y = ax + b to which the shortest distance should be found. By default False, in which the coefficients are found by fitting a straight line to self (atom plane).
  • plot_debug (bool, default False) – If true, a debug plot is saved. The plot is a 3D plot of x and y positions and distance to the plane.
Returns:

list of numbers

Return type:

list of the shortest distances.

get_closest_position_to_point(point_position, extend_line=False)
get_connecting_atom_planes(atom_plane, zone_vector)

Get the outer atom planes which connect self atom plane with another atom plane, through a specific atom plane direction.

The self atom plane, atom plane, and the two returned atom planes will span make a parallelogram.

Parameters:
  • atom_plane (Atomap atom_plane object) –
  • zone_vector (tuple) –
Returns:

Return type:

tuple, two atom plane objects

get_intersecting_atom_from_atom_plane(atom_plane)
get_net_distance_change_between_atoms()

Output [(x,y,z)]

get_slice_between_two_atoms(atom1, atom2)
get_x_position_list()
get_y_position_list()
position_distance_to_neighbor()

Get distance between all atoms and its next neighbor in the atom plane, and the position between these two atoms.

Returns:
Return type:Numpy array [x, y, distance]

Example

>>> from numpy.random import random
>>> from atomap.sublattice import Sublattice
>>> pos = [[x, y] for x in range(9) for y in range(9)]
>>> sublattice = Sublattice(pos, random((9, 9)))
>>> sublattice.construct_zone_axes()
>>> atom_plane = sublattice.atom_plane_list[10]
>>> pos_distance = atom_plane.position_distance_to_neighbor()
>>> x_pos = pos_distance[0]
>>> y_pos = pos_distance[1]
>>> distance = pos_distance[2]
rotation
rotation_ellipticity
sigma_average
sigma_x
sigma_y
sort_atoms_by_distance_to_point(point=(0, 0))
x_position
y_position

Atom position

class atomap.atom_position.Atom_Position(x, y, sigma_x=1.0, sigma_y=1.0, rotation=0.01, amplitude=1.0)

The Atom_Position class contain information about a single atom column.

Parameters:
  • x (float) –
  • y (float) –
  • sigma_x (float, optional) –
  • sigma_y (float, optional) –
  • rotation (float, optional) – In radians. The rotation of the axes of the 2D-Gaussian relative to the image axes. In other words: the rotation of the sigma_x relative to the horizontal axis (x-axis). This is different from the rotation_ellipticity, which is the rotation of the largest sigma in relation to the horizontal axis. For the rotation of the ellipticity, see rotation_ellipticity.
  • amplitude (float, optional) – Amplitude of Gaussian. Stored as amplitude_gaussian attribute.
ellipticity

float

rotation

float – The rotation of sigma_x axis, relative to the x-axis in radians. This value will always be between 0 and pi, as the elliptical 2D-Gaussian used here is always symmetric in the rotation direction, and the perpendicular rotation direction.

rotation_ellipticity

float – The rotation of the longest sigma, relative to the x-axis in radians.

Examples

>>> from atomap.atom_position import Atom_Position
>>> atom_position = Atom_Position(10, 5)

More parameters

>>> atom_pos = Atom_Position(10, 5, sigma_x=2, sigma_y=4, rotation=2)
as_gaussian()
calculate_max_intensity(image_data, percent_to_nn=0.4)

Find the maximum intensity of the atom.

The maximum is found within the the distance to the nearest neighbor times percent_to_nn.

can_atom_plane_be_reached_through_zone_vector(atom_plane, zone_vector)
ellipticity

Largest sigma divided by the shortest

find_atom_intensity_inside_mask(image_data, radius)

Find the average intensity inside a circle.

The circle is defined by the atom position, and the given radius (in pixels). The outside this area is covered by a mask. The average intensity is saved to self.intensity_mask.

get_angle_between_atoms(atom0, atom1=None)

Return the angle between atoms in radians.

Can either find the angle between self and two other atoms, or between another atom and the horizontal axis.

Parameters:
  • atom0 (Atom Position object) – The first atom.
  • atom1 (Atom Position object, optional) – If atom1 is not specified, the angle between itself, atom0 and the horizontal axis will be returned.
Returns:

Angle – Angle in radians

Return type:

float

Examples

>>> from atomap.atom_position import Atom_Position
>>> atom0 = Atom_Position(0, 0)
>>> atom1 = Atom_Position(1, 1)
>>> atom2 = Atom_Position(-1, 1)
>>> angle0 = atom0.get_angle_between_atoms(atom1, atom2)
>>> angle1 = atom0.get_angle_between_atoms(atom1)
get_angle_between_zone_vectors(zone_vector0, zone_vector1)

Return the angle between itself and the next atoms in the atom planes belonging to zone_vector0 and zone_vector1

get_atomic_plane_from_zone_vector(zone_vector)
get_center_position_com(image_data, percent_to_nn=0.4)
get_closest_neighbor()

Find the closest neighbor to an atom in the same sub lattice.

Returns:
Return type:Atomap atom_position object
get_ellipticity_rotation_vector()
get_ellipticity_vector()
get_index_in_atom_plane(atom_plane)
get_neighbor_atoms_in_atomic_plane_from_zone_vector(zone_vector)
get_next_atom_in_atom_plane(atom_plane)
get_next_atom_in_zone_vector(zone_vector)

Get the next atom in the atom plane belonging to zone vector.

get_pixel_difference(atom)

Vector between self and given atom

get_pixel_distance_from_another_atom(atom)
get_pixel_position()
get_position_convergence(distance_to_first_position=False)
get_previous_atom_in_atom_plane(atom_plane)
get_previous_atom_in_zone_vector(zone_vector)
get_rotation_vector()
is_in_atomic_plane(zone_direction)
pixel_distance_from_point(point=(0, 0))
refine_position_using_2d_gaussian(image_data, rotation_enabled=True, percent_to_nn=0.4, centre_free=True)

Use 2D Gaussian to refine the parameters of the atom position.

Parameters:
  • image_data (Numpy 2D array) –
  • rotation_enabled (bool, optional) – If True, the Gaussian will be able to rotate. Note, this can increase the chance of fitting failure. Default True.
  • percent_to_nn (float, optional) – The percent of the distance to the nearest neighbor atom in the same sub lattice. The distance times this percentage defines the mask around the atom where the Gaussian will be fitted. A smaller value can reduce the effect from neighboring atoms, but might also decrease the accuracy of the fitting due to less data to fit to. Default 0.4 (40%).
  • centre_free (bool, default True) – If True, the centre parameter will be free, meaning that the Gaussian can move.
refine_position_using_center_of_mass(image_data, percent_to_nn=0.4)
rotation

The rotation of the atom relative to the horizontal axis.

Given in radians. For the rotation of the ellipticity, see rotation_ellipticity.

rotation_ellipticity

Rotation between the “x-axis” and the major axis of the ellipse.

Rotation between the horizontal axis, and the longest part of the atom position, given by the longest sigma. Basically giving the direction of the ellipticity.

sigma_average
sigma_x
sigma_y

Process parameters

class atomap.process_parameters.ModelParametersBase
add_sublattice_config(sublattice_config_object)
get_sublattice_from_order(order_number)
number_of_sublattices
class atomap.process_parameters.GenericStructure

A Generic structure with one sublattice, the GenericSublattice.

Parameters:
  • peak_separation (None) – No peak separation is set, pixel_separation must be given separately in order to find initial atom positions.
  • name ('A structure') –
  • sublattice_list (list of Sublattices) – Contains one GenericSublattice
class atomap.process_parameters.PerovskiteOxide110

The Perovskite Oxide structure in the 110 projection.

Parameters:
  • name ('Perovskite 110') –
  • peak_separation (0.127 , distance in nm) – Approximately half the distance between the most intense atoms in the structure, used to get initial atom position by peak finding.
  • sublattice_list (list of Sublattices) – Contains 3 sublattices, for A, B and O: PerovskiteOxide110SublatticeACation, PerovskiteOxide110SublatticeBCation, PerovskiteOxide110SublatticeOxygen
class atomap.process_parameters.SrTiO3_110
class atomap.process_parameters.SublatticeParameterBase
class atomap.process_parameters.GenericSublattice

Process parameters for a generic sublattice

color

‘red’ , color of the markers indicating atom positions

image_type

0 – The image will not be inverted.

name

‘S0’ , name of the sublattice

sublattice_order

0 – The sublattice will get the order 0. Higher orders can be used for less intense sublattices in images with multiple sublattices.

zone_axis_list

list – Can have up to 11 zone axis with name from 0-10.

refinement_config

dict – Dict with configuration settings for the refinement of atom positions. 1st refinement : center-of-mass on image modified with background removal, PCA noise filtering and normalization. 2nd refinement : Center of mass on the original image. 3rd refinement : Fitting 2D Gaussians to the original image.

neighbor_distance

0.35 – Mask radius for fitting, set to 0.35 of the distance to nearest neighbor.

atom_subtract_config

list – Configuration for how to subtract higher order sublattices from the image. (Not really used in this structure, but included as an example.)

Examples

>>> import atomap.process_parameters as pp
>>> generic_sublattice = pp.GenericSublattice()
class atomap.process_parameters.PerovskiteOxide110SublatticeACation

Process parameters for the most intense atoms (typically the A-cations) in a Perovskite Oxide structure projected along the (110) direction.

color

‘blue’ , color of the markers indicating atom positions

image_type

0 – The image will not be inverted for finding atom positions.

name

‘A-cation’ , name of the sublattice

sublattice_order

0 – The most intense sublattice gets order 0.

zone_axis_list

list – A list of numbers and names for the zone axes in this projection.

refinement_config

dict – Dict with configuration settings for the refinement of atom positions. Two refinements by fitting 2D Gaussians to the original image.

neighbor_distance

0.35 – Mask radius for fitting set to 0.35 of nearest neighbor.

class atomap.process_parameters.PerovskiteOxide110SublatticeBCation

Process parameters for the second most intense atoms (typically the B-cations) in a Perovskite Oxide structure projected along the (110) direction.

color

‘green’ , color of the markers indicating atom positions

image_type

0 – The image will not be inverted.

name

‘B-cation’ , name of the sublattice

sublattice_order

1 – The sublattice will get the order 1. Higher orders can be used for less intense sublattices in images with multiple sublattices. Lower order is for more intense sublattices.

zone_axis_list

list – A list of numbers and names for the zone axes in this projection.

refinement_config

dict – Dict with configuration settings for the refinement of atom positions. 1st refinement : center-of-mass on the original image. 2nd refinement : Fitting 2D Gaussians to the original image.

neighbor_distance

0.25 – Mask radius for fitting set to 0.25 of nearest neighbor.

atom_subtract_config

list – Configuration for how to subtract higher order sublattices from the image. Subtracts a sublattice with name ‘A-cation’ from the original image.

class atomap.process_parameters.PerovskiteOxide110SublatticeOxygen

Process parameters for the least intense atoms (typically the Oxygen) in a Perovskite Oxide structure projected along the (110) direction.

color

‘red’ , color of the markers indicating atom positions

image_type

1 – The image will be inverted. Oxygen can be visible in ABF images. The image is inverted such that the atoms are bright dots on a dark background.

name

‘Oxygen’ , name of the sublattice

sublattice_order

2 – The sublattice will get the order 2, being the third most intense sublattice.

zone_axis_list

list – A list of numbers and names for the zone axes in this projection.

refinement_config

dict – Dict with configuration settings for the refinement of atom positions. 1st refinement : center-of-mass on the original image. 2nd refinement : Fitting 2D Gaussians to the original image.

neighbor_distance

0.25 – Mask radius for fitting set to 0.25 of nearest neighbor.

atom_subtract_config

list – Configuration for how to subtract higher order sublattices from the image. Subtracts first a sublattice with name ‘A-cation’ from the inverted image, then a sublattice with name ‘B-cation’.

Modules

Atom finding refining

atomap.atom_finding_refining.construct_zone_axes_from_sublattice(sublattice, atom_plane_tolerance=0.5, zone_axis_para_list=False)

Constructs zone axes for a sublattice.

The zone axes are constructed by finding the 15 nearest neighbors for each atom position in the sublattice, and finding major translation symmetries among the nearest neighbours. Only unique zone axes are kept, and “bad” ones are removed.

Parameters:
  • sublattice (Atomap Sublattice) –
  • atom_plane_tolerance (scalar, default 0.5) – When constructing the atomic planes, the method will try to locate the atoms by “jumping” one zone vector, and seeing if there is an atom with the pixel_separation times atom_plane_tolerance. So this value should be increased the atomic planes are non-continuous and “split”.
  • zone_axis_para_list (parameter list or bool, default False) – A zone axes parameter list is used to name and index the zone axes. See atomap.process_parameters for more info. Useful for automation.

Example

>>> import atomap.api as am
>>> sublattice = am.dummy_data.get_simple_cubic_sublattice()
>>> sublattice
<Sublattice,  (atoms:400,planes:0)>
>>> import atomap.atom_finding_refining as afr
>>> afr.construct_zone_axes_from_sublattice(sublattice)
>>> sublattice
<Sublattice,  (atoms:400,planes:4)>

See also

sublattice._make_translation_symmetry()
How unique zone axes are found
sublattice._remove_bad_zone_vectors()
How fragmented (“bad zone axis”) are identified and removed.
atomap.process_parameters()
more info on zone axes parameter list
atomap.atom_finding_refining.do_pca_on_signal(signal, pca_components=22)
atomap.atom_finding_refining.find_feature_density(image_data, separation_range=None, separation_step=1, plot_figure=False, plot_debug_figures=False)

Do peak finding with a varying amount of peak separation constrained. Gives a measure of feature density, and what peak separation should be used to find the initial sub-lattice.

Inspiration from the program Smart Align by Lewys Jones.

atomap.atom_finding_refining.find_features_by_separation(signal, separation_range, separation_step=1, threshold_rel=0.02, pca=False, subtract_background=False, normalize_intensity=False, show_progressbar=True)

Do peak finding with a varying amount of peak separation constrained.

Inspiration from the program Smart Align by Lewys Jones.

Parameters:
  • signal (HyperSpy 2D signal) –
  • separation_range (tuple) – Lower and upper end of minimum pixel distance between the features.
  • separation_step (int, optional) –
  • show_progressbar (bool, default True) –
Returns:

Return type:

tuple, (separation_list, peak_list)

atomap.atom_finding_refining.fit_atom_positions_gaussian(atom_list, image_data, rotation_enabled=True, percent_to_nn=0.4, mask_radius=None, centre_free=True)

Fit a list of Atom_Positions to an image using 2D Gaussians.

The results of the fitting will be saved in the Atom_Position objects themselves, and the old positions will be added to atom.old_pixel_x_list and atom.old_pixel_y_list.

This type of fitting is prone to errors, especially where the atoms are close together, and on noisy images. To reduce the odds of bad fitting the function will rerun the model with a slightly smaller mask_radius or percent_to_nn if: - Any of the fitted positions are outside the model region. - The sigma ratio (highest_sigma/lowest_sigma) is higher than 4.

This repeats 10 times, and if the fitting is still bad, center of mass will be used to set the center position for the atom.

Parameters:
  • atom_list (list of Atom_Position objects) – For example the members of a dumbbell.
  • image_data (NumPy 2D array) –
  • rotation_enabled (bool, optional) – If True (default), the 2D Gaussian will be able to rotate.
  • percent_to_nn (float, optional) –
  • mask_radius (float, optional) – Radius of the mask around each atom. If this is not set, the radius will be the distance to the nearest atom in the same sublattice times the percent_to_nn value. Note: if mask_radius is not specified, the Atom_Position objects must have a populated nearest_neighbor_list. This is normally done through the sublattice class, but can also be done manually. See below for an example how to do this.
  • centre_free (bool, optional) – If True (default), the Gaussian will be free to move. Setting this to False can be used to find better values for the other parameters like sigma and A, without the centre positions causing bad fitting.

See also

_make_model_from_atom_list(), _fit_atom_positions_with_gaussian_model(), atom_lattice.Dumbbell_Lattice(), positions()

Examples

>>> import numpy as np
>>> from atomap.atom_position import Atom_Position
>>> from atomap.atom_finding_refining import fit_atom_positions_gaussian

Fitting atomic columns one-by-one

>>> atom_list = [Atom_Position(2, 2), Atom_Position(4, 4)]
>>> image = np.zeros((9, 9))
>>> for atom_position in atom_list:
...     g_list = fit_atom_positions_gaussian(
...         atom_list=[atom_position], image_data=image, mask_radius=2)

Fitting two atoms together

>>> atom_list = [Atom_Position(2, 2), Atom_Position(4, 4)]
>>> image = np.zeros((9, 9))
>>> image[2, 2] = 1.
>>> image[4, 4] = 1.
>>> g_list = fit_atom_positions_gaussian(
...     atom_list=atom_list, image_data=image, mask_radius=2)

Not using mask_radius, populating the nearest_neighbor_list manually

>>> image = np.zeros((9, 9))
>>> image[2, 2] = 1.
>>> image[5, 5] = 1.
>>> atom0 = Atom_Position(2, 2, 0.5, 0.5)
>>> atom1 = Atom_Position(5, 5, 0.5, 0.5)
>>> atom0.nearest_neighbor_list = [atom1]
>>> atom1.nearest_neighbor_list = [atom0]
>>> g_list = fit_atom_positions_gaussian([atom0, atom1], image)
atomap.atom_finding_refining.get_atom_positions(signal, separation=5, threshold_rel=0.02, pca=False, subtract_background=False, normalize_intensity=False, remove_too_close_atoms=True)

Find the most intense features in a HyperSpy signal, where the features has to be separated by a minimum distance.

Parameters:
  • signal (HyperSpy 2D signal) –
  • separation (number) – Minimum separation between the features.
  • threshold_rel (float, default 0.02) –
  • pca (bool, default False) – Do PCA on the signal before doing the peak finding.
  • subtract_background (bool, default False) – Subtract the average background from the signal before doing the peak finding.
  • normalize_intensity (bool, default False) –
  • remove_too_close_atoms (bool, default True) – Will attempt to find and remove atoms which are too close to eachother, i.e. less than separation.
Returns:

Return type:

NumPy array, list of the most intense atom positions.

Example

>>> import numpy as np
>>> import atomap.api as am
>>> s = am.dummy_data.get_simple_cubic_signal()
>>> from atomap.atom_finding_refining import get_atom_positions
>>> atom_positions = get_atom_positions(s, 5)
>>> peak_x = atom_positions[:,0]
>>> peak_y = atom_positions[:,1]
atomap.atom_finding_refining.get_feature_separation(signal, separation_range=(5, 30), separation_step=1, pca=False, subtract_background=False, normalize_intensity=False, threshold_rel=0.02, show_progressbar=True)

Plot the peak positions on in a HyperSpy signal, as a function of peak separation.

Parameters:
  • signal (HyperSpy signal 2D) –
  • separation_range (tuple, optional, default (5, 30)) –
  • separation_step (int, optional, default 1) –
  • pca (bool, default False) –
  • subtract_background (bool, default False) –
  • normalize_intensity (bool, default False) –
  • threshold_rel (float, default 0.02) –
  • show_progressbar (bool, default True) –

Example

>>> import numpy as np
>>> import hyperspy.api as hs
>>> from atomap.atom_finding_refining import get_feature_separation
>>> s = hs.signals.Signal2D(np.random.random((500, 500)))
>>> s1 = get_feature_separation(s)
atomap.atom_finding_refining.normalize_signal(signal, invert_signal=False)
atomap.atom_finding_refining.refine_sublattice(sublattice, refinement_config_list, percent_to_nn)
atomap.atom_finding_refining.subtract_average_background(signal, gaussian_blur=30)

Input output (IO)

Module containing functions to save and load Atom_Lattice objects.

atomap.io.load_atom_lattice_from_hdf5(filename, construct_zone_axes=True)

Load an Atomap HDF5-file, restoring a saved Atom_Lattice.

Parameters:
  • filename (string) – Filename of the HDF5-file.
  • construct_zone_axes (bool) – If True, find relations between atomic positions by constructing atomic planes. Default True.
Returns:

Return type:

Atomap Atom_Lattice object

atomap.io.save_atom_lattice_to_hdf5(atom_lattice, filename, overwrite=False)

Store an Atom_Lattice object as an HDF5-file.

Parameters:
  • atom_lattice (Atomap Atom_Lattice object) –
  • filename (string) –
  • overwrite (bool, default False) –

Main

atomap.main.make_atom_lattice_from_image(s_image0, process_parameter=None, pixel_separation=None, s_image1=None, debug_plot=False)
atomap.main.run_image_filtering(signal, invert_signal=False)

Subtracts background, filters noise with PCA, and normalizes a signal.

Parameters:
  • signal (HyperSpy Signal2D) –
  • invert_signal (bool, default False) – Inverts the image to 1./signal.data
Returns:

filtered_signal

Return type:

HyperSpy Signal2D

Example

>>> import atomap.api as am
>>> s = am.dummy_data.get_simple_cubic_signal()
>>> from atomap.main import run_image_filtering
>>> s_new = run_image_filtering(s)

Plotting

atomap.plotting.get_rgb_array(angle, magnitude, rotation=0, angle_lim=None, magnitude_lim=None)
atomap.plotting.make_color_wheel(ax, rotation=0)
atomap.plotting.normalize_array(np_array, max_number=1.0)
atomap.plotting.plot_complex_image_map_line_profile_using_interface_plane(image_data, amplitude_data, phase_data, line_profile_amplitude_data, line_profile_phase_data, interface_plane, atom_plane_list=None, data_scale=1, atom_list=None, extra_marker_list=None, amplitude_image_lim=None, phase_image_lim=None, plot_title='', add_color_wheel=False, color_bar_markers=None, vector_to_plot=None, rotate_atom_plane_list_90_degrees=False, prune_outer_values=False, figname='map_data.jpg')
Parameters:
  • atom_list (list of Atom_Position instances) –
  • extra_marker_list (two arrays of x and y [[x_values], [y_values]]) –
atomap.plotting.plot_feature_density(separation_value_list, peakN_list, figname='feature_density.png')
atomap.plotting.plot_image_map_line_profile_using_interface_plane(image_data, heatmap_data_list, line_profile_data_list, interface_plane, atom_plane_list=None, data_scale=1, data_scale_z=1, atom_list=None, extra_marker_list=None, clim=None, plot_title='', vector_to_plot=None, rotate_atom_plane_list_90_degrees=False, prune_outer_values=False, figname='map_data.jpg')
Parameters:
  • atom_list (list of Atom_Position instances) –
  • extra_marker_list (two arrays of x and y [[x_values], [y_values]]) –
atomap.plotting.plot_line_profiles_from_parameter_input(parameter_list, parameter_name_list=None, invert_line_profiles=False, extra_line_marker_list=None, x_lim=False, figname='line_profile_list.jpg')
atomap.plotting.plot_vector_field(x_pos_list, y_pos_list, x_rot_list, y_rot_list, save=True)
atomap.plotting.plot_zone_vector_and_atom_distance_map(image_data, distance_data, atom_planes=None, distance_data_scale=1, atom_list=None, extra_marker_list=None, clim=None, atom_plane_marker=None, plot_title='', vector_to_plot=None, figsize=(10, 20), figname='map_data.jpg')
Parameters:
  • atom_list (list of Atom_Position instances) –
  • extra_marker_list (two arrays of x and y [[x_values], [y_values]]) –

Tools

class atomap.tools.Fingerprinter(cluster_algo=DBSCAN(algorithm='auto', eps=0.1, leaf_size=30, metric='euclidean', min_samples=10, n_jobs=1, p=None))

Produces a distance-fingerprint from an array of neighbor distance vectors.

To avoid introducing our own interface we’re going to use scikit-learn Estimator conventions to name the method, which produces our fingerprint, ‘fit’ and store our estimations as attributes with a trailing underscore in their names. http://scikit-learn.org/stable/developers/contributing.html#fitting http://scikit-learn.org/stable/developers/contributing.html#estimated-attributes)

fingerprint_

array, shape = (n_clusters,) – The main result. The contents of fingerprint_ can be described as the relative distances to neighbours in the generalized neighborhood.

cluster_centers_

array, shape = (n_clusters, n_dimensions) – The cluster center coordinates from which the fingerprint was produced.

cluster_algo_.labels_

array, shape = (n_points,) – Integer labels that denote which cluster each point belongs to.

fit(X, max_neighbors=150000)
Parameters:
  • X (array, shape = (n_points, n_dimensions)) – This array is typically a transpose of a subset of the returned value of sublattice.get_nearest_neighbor_directions_all()
  • max_neighbors (int, default 150000) – If the length of X is larger than max_neighbors, X will be reduced to max_neighbors. The selection is random. This is done to allow very large datasets to be processed, since having X too large causes the fitting to use too much memory.

Notes

More information about memory use: http://scikit-learn.org/stable/modules/clustering.html#dbscan

atomap.tools.array2signal1d(array, scale=1.0, offset=0.0)
atomap.tools.array2signal2d(numpy_array, scale=1.0, rotate_flip=False)
atomap.tools.calculate_angle(v1, v2)
atomap.tools.combine_clustered_positions_into_layers(data_list, layer_distance, combine_layers=True)

Combine clustered positions into groups.

Atoms with a similar distance for a line belong to the same plane parallel to this line. Atoms in data_list are grouped based on which plane they belong to. If there is only one atom in a layer, it will be disregarded as it gives a high uncertainty.

Parameters:
  • data_list (NumPy array) – An array where the distance from a point to a line is in input_data_list[:, 0], and the property of the point (atom) is in [:,1]
  • layer_distance (float) – The half width of a layer, used to determine which layer an atom belongs to.
  • combine_layers (bool, default True) – If True, the values for distance and property is averaged for each layer.
Returns:

  • A list, layer_list. If combine_layers is True, a list of the average
  • position and property of the points in the layer. If False, a nested
  • list where each element in layer_list contains a list the atoms (position
  • and property) in the layer.

atomap.tools.combine_clusters_using_average_distance(data_list, margin=0.5)
atomap.tools.dotproduct(v1, v2)
atomap.tools.find_atom_position_between_atom_planes(image, atom_plane0, atom_plane1, orthogonal_zone_vector, integration_width_percent=0.2, max_oxygen_sigma_percent=0.2)
atomap.tools.find_atom_positions_for_all_atom_planes(image, sublattice, parallel_zone_vector, orthogonal_zone_vector)
atomap.tools.find_atom_positions_for_an_atom_plane(image, atom_plane0, atom_plane1, orthogonal_zone_vector)
atomap.tools.find_average_distance_between_atoms(input_data_list, crop_start=3, crop_end=3, threshold=0.4)

Return the distance between monolayers.

Returns the maximal separation between two adjacent points in input_data_list[:,0], as a good approximation for monolayer separation.

Parameters:
  • input_data_list (NumPy array) – An array where the distance from a point to a line is in input_data_list[:, 0]
  • crop_end (crop_start,) – Before and after the index given by crop_start and crop_end, the data is ignored. By default 3, to ignore outliers.
  • threshold (float, default 0.4) – Atoms with a separation of more than the threshold times the largest atom separation will be counted as members of different planes.
Returns:

  • first_peak (float) – The monolayer separation.
  • monolayer_sep (array) – An array with monolayer separations
  • mean_separation (float) – The mean monolayer separation

atomap.tools.fliplr_points_and_signal(signal, x_array, y_array)

Horizontally flip a set of points and a HyperSpy signal.

For making sure both the image and points are flipped correctly.

Parameters:
  • signal (HyperSpy 2D signal) –
  • y_array (x_array,) –
Returns:

  • flipped_signal (HyperSpy signal)
  • flipped_x_array, flipped_y_array (NumPy arrays)

Examples

>>> import atomap.api as am
>>> import atomap.tools as to
>>> sublattice = am.dummy_data.get_distorted_cubic_sublattice()
>>> s = sublattice.get_atom_list_on_image()
>>> x, y = sublattice.x_position, sublattice.y_position
>>> s_flip, x_flip, y_flip = to.fliplr_points_and_signal(s, x, y)

Plotting the data

>>> import matplotlib.pyplot as plt
>>> fig, ax = plt.subplots()
>>> cax = ax.imshow(s_flip.data, origin='lower',
...           extent=s_flip.axes_manager.signal_extent)
>>> cpoints = ax.scatter(x_flip, y_flip)
atomap.tools.fliplr_points_around_signal_centre(signal, x_array, y_array)

Horizontally flip a set of points around the centre of a signal.

Parameters:
  • signal (HyperSpy 2D signal) –
  • y_array (x_array,) –
Returns:

flipped_x_array, flipped_y_array

Return type:

NumPy arrays

Examples

>>> import atomap.api as am
>>> import atomap.tools as to
>>> sublattice = am.dummy_data.get_distorted_cubic_sublattice()
>>> s = sublattice.get_atom_list_on_image()
>>> x, y = sublattice.x_position, sublattice.y_position
>>> x_rot, y_rot = to.fliplr_points_around_signal_centre(s, x, y)
atomap.tools.get_arbitrary_slice(image, start_point, end_point, width, debug_figname=None)
atomap.tools.get_atom_planes_square(sublattice, atom_plane1, atom_plane2, interface_atom_plane, zone_vector, debug_plot=False)
Parameters:
  • sublattice (Atomap sublattice object) –
  • atom_plane2 (atom_plane1,) –
atomap.tools.get_point_between_four_atoms(atom_list)
atomap.tools.get_point_between_two_atoms(atom_list)
atomap.tools.get_slice_between_four_atoms(image, start_atoms, end_atoms, width)
atomap.tools.get_slice_between_two_atoms(image, atom0, atom1, width)
atomap.tools.integrate(s, points_x, points_y, method='Voronoi', max_radius='Auto')

Given a spectrum image a set of points and a maximum outer radius, this function integrates around each point in an image, using either Voronoi cell or watershed segmentation methods.

Parameters:
  • s (HyperSpy signal or array-like object) – Assuming 2D, 3D or 4D dataset where the spatial dimensions are 2D and any remaining dimensions are spectral.
  • point_y (point_x,) – Detailed list of the x and y coordinates of each point of interest within the image.
  • method (string) – ‘Voronoi’ or ‘Watershed’
  • max_radius ({'Auto'} int) – A maximum outer radius for each Voronoi Cell. If a pixel exceeds this radius it will not be included in the cell. This allows analysis of a surface and particles. If ‘max_radius’ is left as ‘Auto’ then it will be set to the largest dimension in the image.
Returns:

  • integrated_intensity (NumPy array) – An array where dimension 0 is the same length as points, and subsequent subsequent dimension are energy dimensions.
  • intensity_record (HyperSpy signal, same size as s) – Each pixel/voxel in a particular segment or region has the value of the integration, value.
  • point_record (NumPy array, same size as image) – Image showing where each integration region is, pixels equating to point 0 (integrated_intensity[0]) all have value 0, all pixels equating to integrated_intensity[1] all have value 1 etc.

Examples

>>> import atomap.api as am
>>> from atomap.tools import integrate
>>> import hyperspy.api as hs
>>> sublattice = am.dummy_data.get_simple_cubic_sublattice(
...        image_noise=True)
>>> image = hs.signals.Signal2D(sublattice.image)
>>> i_points, i_record, p_record = integrate(
...        image,
...        points_x=sublattice.x_position,
...        points_y=sublattice.y_position, method='Voronoi')
>>> i_record.plot()

For a 3 dimensional dataset, with artificial EELS data

>>> s = am.dummy_data.get_eels_spectrum_survey_image()
>>> s_eels = am.dummy_data.get_eels_spectrum_map()
>>> peaks = am.get_atom_positions(s, separation=4)
>>> i_points, i_record, p_record = integrate(
...         s_eels, peaks[:, 0], peaks[:, 1], max_radius=3)

Note

Works in principle with 3D and 4D data sets but will quickly hit a memory error with large sizes.

atomap.tools.length(v)
atomap.tools.project_position_property_sum_planes(input_data_list, interface_plane, rebin_data=True)

Project 2D positions onto a 1D plane. The 2D positions are found as function of distance to the interface_plane. If rebin_data is True, the function will attempt to sum the positions belonging to the same plane. In this case, one will get the positions as a function of atomic plane from the interface_plane.

Note, positions will not be projected on to the interface_plane, but on a plane perpendicular to the interface_plane. The returned data will give the distance from the projected positions to the interface_plane.

Parameters:
  • input_data_list (Numpy array, [Nx3]) – Numpy array with positions and property value. Must be in the from [[x,y,z]], so that the x-positions are extracted using input_data_list[:,0]. y-positions using input_data_list[:,1]. Property value using input_data_list[:,2].
  • interface_plane (Atomap atom_plane object) –
  • rebin_data (bool, optional) – If True, will attempt to combine the data points which belong to the same atomic plane. The points which belong to the same plane will be averaged into a single value for each atomic plane. This will give the property value as a function distance to the interface_plane.
Returns:

Data list – Array in the form [[position, property]].

Return type:

NumPy Array

Example

>>> from numpy.random import random
>>> from atomap.sublattice import Sublattice
>>> pos = [[x, y] for x in range(9) for y in range(9)]
>>> sublattice = Sublattice(pos, random((9, 9)))
>>> sublattice.construct_zone_axes()
>>> x, y = sublattice.x_position, sublattice.y_position
>>> z = sublattice.ellipticity
>>> input_data_list = np.array([x, y, z]).swapaxes(0, 1)
>>> from atomap.tools import project_position_property_sum_planes
>>> plane = sublattice.atom_plane_list[10]
>>> data = project_position_property_sum_planes(input_data_list, plane)
>>> positions = data[:,0]
>>> property_values = data[:,1]
>>> cax = plt.plot(positions, property_values)
atomap.tools.remove_atoms_from_image_using_2d_gaussian(image, sublattice, percent_to_nn=0.4, show_progressbar=True)
Parameters:
  • image (NumPy 2D array) –
  • sublattice (Atomap sublattice object) –
  • percent_to_nn (float) –
  • show_progressbar (bool, default True) –
Returns:

subtracted_image

Return type:

NumPy 2D array

atomap.tools.rotate_points_and_signal(signal, x_array, y_array, rotation)

Rotate a set of points and a HyperSpy signal.

For making sure both the image and points are rotated correctly.

Parameters:
  • signal (HyperSpy 2D signal) –
  • y_array (x_array,) –
  • rotation (scalar) –
Returns:

  • rotated_signal (HyperSpy signal)
  • rotated_x_array, rotated_y_array (NumPy arrays)

Examples

>>> import atomap.api as am
>>> import atomap.tools as to
>>> sublattice = am.dummy_data.get_distorted_cubic_sublattice()
>>> s = sublattice.get_atom_list_on_image()
>>> x, y = sublattice.x_position, sublattice.y_position
>>> s_rot, x_rot, y_rot = to.rotate_points_and_signal(s, x, y, 30)

Plotting the data

>>> fig, ax = plt.subplots()
>>> cax = ax.imshow(s_rot.data, origin='lower',
...                 extent=s_rot.axes_manager.signal_extent)
>>> cpoints = ax.scatter(x_rot, y_rot)
atomap.tools.rotate_points_around_signal_centre(signal, x_array, y_array, rotation)

Rotate a set of points around the centre of a signal.

Parameters:
  • signal (HyperSpy 2D signal) –
  • y_array (x_array,) –
  • rotation (scalar) –
Returns:

rotated_x_array, rotated_y_array

Return type:

NumPy arrays

Examples

>>> import atomap.api as am
>>> import atomap.tools as to
>>> sublattice = am.dummy_data.get_distorted_cubic_sublattice()
>>> s = sublattice.get_atom_list_on_image()
>>> x, y = sublattice.x_position, sublattice.y_position
>>> x_rot, y_rot = to.rotate_points_around_signal_centre(s, x, y, 30)

Dummy data

atomap.dummy_data.get_distorted_cubic_signal(image_noise=False)

Generate a test image signal of a distorted cubic atomic structure.

Parameters:image_noise (default False) – If True, will add Gaussian noise to the image.
Returns:signal
Return type:HyperSpy 2D

Examples

>>> import atomap.api as am
>>> s = am.dummy_data.get_distorted_cubic_signal()
>>> s.plot()
atomap.dummy_data.get_distorted_cubic_sublattice(image_noise=False)

Generate a test sublattice of a distorted cubic atomic structure.

Parameters:image_noise (default False) – If True, will add Gaussian noise to the image.
Returns:sublattice
Return type:Atomap Sublattice

Examples

>>> import atomap.api as am
>>> sublattice = am.dummy_data.get_distorted_cubic_sublattice()
>>> sublattice.plot()
atomap.dummy_data.get_dumbbell_signal()
atomap.dummy_data.get_eels_spectrum_map(add_noise=True)

Get an artificial atomic resolution EELS map of LaMnO3

Containing the Mn-L23 and La-M54 edges.

Parameters:add_noise (bool) – If True, will add Gaussian noise to the spectra. Default True.
Returns:eels_map
Return type:HyperSpy EELSSpectrum

Example

>>> import atomap.api as am
>>> s_eels_map = am.dummy_data.get_eels_spectrum_map()
>>> s_eels_map.plot()

Not adding noise

>>> s_eels_map = am.dummy_data.get_eels_spectrum_map(add_noise=False)

See also

get_eels_spectrum_survey_image()
signal with same spatial dimensions
atomap.dummy_data.get_eels_spectrum_survey_image()

Get an artificial survey image of atomic resolution EELS map of LaMnO3

Returns:survey_image
Return type:HyperSpy Signal2D

Example

>>> import atomap.api as am
>>> s = am.dummy_data.get_eels_spectrum_survey_image()
>>> s.plot()

See also

get_eels_spectrum_map()
corresponding EELS map
atomap.dummy_data.get_fantasite()

Fantasite is a fantastic structure with several interesting structural variations.

It contains two sublattices, domains with elliptical atomic columns and tilt-patterns. This function returns a HyperSpy 2D signal.

Examples

>>> import atomap.api as am
>>> s = am.dummy_data.get_fantasite()
>>> s.plot()
Returns:fantasite_signal
Return type:HyperSpy 2D signal

See also

get_fantasite_sublattice()
get a sublattice object of the fantasite.
atomap.dummy_data.get_fantasite_atom_lattice()

Fantasite is a fantastic structure with several interesting structural variations.

It contains two sublattices, domains with elliptical atomic columns and tilt-patterns. This function returns an Atomap Atom_Lattice.

Examples

>>> import atomap.api as am
>>> atom_lattice = am.dummy_data.get_fantasite_atom_lattice()
>>> atom_lattice.plot()
atomap.dummy_data.get_fantasite_sublattice()

Fantasite is a fantastic structure with several interesting structural variations.

It contains two sublattices, domains with elliptical atomic columns and tilt-patterns. This function returns an Atomap sublattice.

Currently this function returns the two sublattices as one sublattice. To get these sublattices separately, see get_fantasite_atom_lattice

Examples

>>> import atomap.api as am
>>> sublattice = am.dummy_data.get_fantasite_sublattice()
>>> sublattice.plot()

See also

get_fantasite()
get a Signal2D object of the fantasite.
get_fantasite_atom_lattice()
get the atom lattice for fantasite, with both sublattices.
atomap.dummy_data.get_perovskite110_ABF_signal(image_noise=False)

Returns signals representing an ABF image of a perovskite <110>.

Parameters:image_noise (default False) – If True, will add Gaussian noise to the image.
Returns:signal
Return type:HyperSpy 2D,

Examples

>>> import atomap.api as am
>>> s_ABF = am.dummy_data.get_perovskite110_ABF_signal()
atomap.dummy_data.get_simple_atom_lattice_two_sublattices(image_noise=False)

Returns a simple atom_lattice with two sublattices.

Parameters:image_noise (bool, default False) –
Returns:simple_atom_lattice
Return type:Atom_Lattice object

Examples

>>> import atomap.api as am
>>> al = am.dummy_data.get_simple_atom_lattice_two_sublattices()
>>> al.plot()
atomap.dummy_data.get_simple_cubic_signal(image_noise=False)

Generate a test image signal of a simple cubic atomic structure.

Parameters:image_noise (default False) – If True, will add Gaussian noise to the image.
Returns:signal
Return type:HyperSpy 2D

Examples

>>> import atomap.api as am
>>> s = am.dummy_data.get_simple_cubic_signal()
>>> s.plot()
atomap.dummy_data.get_simple_cubic_sublattice(image_noise=False)

Generate a test sublattice of a simple cubic atomic structure.

Parameters:image_noise (default False) – If True, will add Gaussian noise to the image.
Returns:sublattice
Return type:Atomap Sublattice

Examples

>>> import atomap.api as am
>>> sublattice = am.dummy_data.get_simple_cubic_sublattice()
>>> sublattice.plot()
atomap.dummy_data.get_simple_heterostructure_signal(image_noise=True)
atomap.dummy_data.get_two_sublattice_signal()