10 STEM DPC Analysis of Magnetic Sample.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Analysing magnetic materials using STEM-DPC\n",
    "\n",
    "This notebook shows how to use the `pyXem` library to analyse 4-D scanning transmission electron microscopy (STEM) data, specifical magnetic materials using differential phase contrast (DPC). For more information about this imaging method, see the Wikipedia article on Scanning Transmission Electron Microscopy, which has a subsection on DPC: https://en.wikipedia.org/wiki/Scanning_transmission_electron_microscopy#Differential_phase_contrast\n",
    "\n",
    "The data we'll be looking at there is from the paper **Strain Anisotropy and Magnetic Domains in Embedded Nanomagnets**, and is STEM data recorded on a Merlin fast pixelated electron detector system, where the objective lens has been turned off.\n",
    "This allows for magnetic information to be extracted, by carefully mapping the beam shifts.\n",
    "\n",
    "More documentation about pyXem is found at https://pyxem.github.io/pyxem-website/\n",
    "\n",
    "Journal article:\n",
    "* **Strain Anisotropy and Magnetic Domains in Embedded Nanomagnets**\n",
    "* Nord, M., Semisalova, A., Kákay, A., Hlawacek, G., MacLaren, I., Liersch, V., Volkov, O. M., Makarov, D., Paterson, G. W., Potzger, K., Lindner, J., Fassbender, J., McGrouther, D., Bali, R.\n",
    "* Small 2019, 15, 1904738. https://doi.org/10.1002/smll.201904738\n",
    "\n",
    "The full dataset and scripts used in analysing this data is found at Zenodo: https://zenodo.org/record/3466591"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Importing libraries"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The first step is setting the plotting toolkit"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "nbval-skip"
    ]
   },
   "outputs": [],
   "source": [
    "%matplotlib qt5"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Then import the library itself"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import hyperspy.api as hs\n",
    "import pyxem as pxm"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Downloading data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that is a fairly large file of about 5 GB, so downloading it might take a while."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "nbval-skip"
    ]
   },
   "outputs": [],
   "source": [
    "import os.path\n",
    "try:\n",
    "    os.mkdir('datasets')\n",
    "except:\n",
    "    pass\n",
    "if not os.path.exists('datasets/fe60al40_stripe_pattern.hspy'):\n",
    "    import urllib.request\n",
    "    urllib.request.urlretrieve('https://zenodo.org/record/4312960/files/fe60al40_stripe_pattern.hspy', 'datasets/fe60al40_stripe_pattern.hspy')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Loading data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As the file is pretty big, we load it lazily."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s = hs.load(\"datasets/fe60al40_stripe_pattern.hspy\", lazy=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s = s.inav[0:200, :]  # Run this cell to reduce the size of the signal, which will make most of the processing much faster"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Plotting the data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Lets first have a look at the data, using the `slider` navigator."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s.plot(navigator='slider')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The direct beam is at (x, y) = (128, 125). Lets create a navigation signal from a single pixel in the direct beam."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s_nav = s.isig[128, 125].T\n",
    "s_nav.compute()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This was pretty quick, thanks to how the chunking is for this dataset: (16, 16, 16, 16). So only about 0.4% of the full dataset had to be loaded into memory.\n",
    "\n",
    "Plotting it, we can see the nanocrystalline structure of the material."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s_nav.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And then we can use it to navigate over the full dataset."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s.plot(navigator=s_nav)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "What we're interested in here is the movement of the centre beam. To make the plotting more responsive, we can grab just the centre of the diffraction patterns.\n",
    "\n",
    "Move the navigator up and down in the centre of the sample, and look for small changes in the beam position."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s_crop = s.isig[128-16:128+16, 128-16:128+16]\n",
    "s_crop.plot(navigator=s_nav)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "It can be really tricky noticing any systematic beam shifts. A better way is transposing the navigation dimensions."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Navigation in the detector dimensions"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A powerful feature in pyXem, which is inherited from HyperSpy, is the ability to \"flip\" the navigation dimensions.\n",
    "So instead of navigating across the probe position, you can navigate over the detector positions.\n",
    "\n",
    "Firstly, lets make a navigation signal for that, by averaging several of the diffraction patterns.\n",
    "Here we use the `s_crop` signal from earlier."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s_diff_nav = s_crop.inav[100:110, 100:110].mean(axis=(0, 1))\n",
    "s_diff_nav.compute()\n",
    "s_diff_nav.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Then we transpose `s_crop`, and use `s_diff_nav` for navigation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s_transpose = s_crop.transpose()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now, move the navigator across the direct beam, and look for the magnetic stripe in the signal plot."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s_transpose.plot(navigator=s_diff_nav)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This visualizes the shift of the beam, which is caused by the beam passing through the ferromagnetic domains in the material.\n",
    "\n",
    "However, it is not very quantitative. So lets try to extract the beam shifts using center of mass."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Extracting the beam shift"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "There are many ways of extracting the beam shift. One way is using center of mass, which is a fairly simple method.\n",
    "\n",
    "Again, we use `s_crop`, since we're only interested in the direct beam."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s_crop_com = s_crop.center_of_mass()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This returns a `DPCSignal2D` class, which will be explored more later. What we need to know is that is it basically a HyperSpy `Signal2D` class, where the x-beam shifts are in the first navigation index (`s.inav[0]`), while the y-shifts are in the second navigation index (`s.inav[1]`).\n",
    "\n",
    "Plotting it, we can see that the ferromagnetic domains are much more visible!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s_crop_com.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "However, we're also getting contrast from the other effects, such as structural effects. Since the sample is nanocrystalline, some of the grains will be close to some zone axis, giving heavy Bragg scattering. While the Bragg spots themselves won't be visible at such low scattering angles as we have in `s_crop`, it will still change the intensity distribution _within_ the direct beam. Essentially, the direct beam becomes non-uniform, which will have an effect similarly to beam shift.\n",
    "\n",
    "One way of reducing this is by using thresholding and masking. However, we first need to find reasonable values for these.\n",
    "\n",
    "For this, we use `threshold_and_mask` on a subset of the dataset."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s_subset = s_crop.inav[:20, :20]\n",
    "s_subset.compute()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s_subset.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`threshold_and_mask` is a useful way to see what exactly what the pre-processing of the center of mass algorithm does, as we want to extract the beam shifts themselves, not intensity distributions _within_ the disk, or diffraction spots _outside_ the disk.\n",
    "\n",
    "It works by getting the mean of the intensity inside the masked area, times the threshold. Then, any pixel lower or equal than that value is set to zero, while any value above that value is set to one. Ideally, this should remove the influence of other diffraction effects, and non-uniform direct beam."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s_threshold_mask = s_subset.threshold_and_mask(threshold=3., mask=(16, 16, 16))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s_threshold_mask.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Having found some optimized parameters, we can use this to run the center of mass processing again."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s_crop_com_threshold = s_crop.center_of_mass(threshold=3., mask=(16, 16, 16), lazy_result=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s_crop_com_threshold.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This seems to have improved it a little bit.\n",
    "\n",
    "For doing more advanced processing of this type of data, see https://fpdpy.gitlab.io/fpd/, which has better edge detection algorithms in the `fpd.fpd_processing.phase_correlation` function. For a complete example on how to use this, see https://zenodo.org/record/3466591/files/d001_get_dpc_raw_signal.py which processes the same type of dataset as used in this example."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Correcting d-scan"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "With the beam shift extracted, we will remove the effects of impure beam shift (d-scan).\n",
    "This is due to various instrument misalignments, and leads to a change in beam position in the probe plane becoming a shift of the beam in the detector plane.\n",
    "Luckily, in most situations, the d-scan is linear across the dataset, meaning it can be removed using a simple plane subtraction.\n",
    "\n",
    "However, for the full version of this dataset (if you didn't crop it earlier), a simple plane subtraction is not enough to remove all the effects of this d-scan.\n",
    "See `fpd.ransac_tools.ransac_im_fit` for more advanced ways of removing it, for example in https://zenodo.org/record/3466591/files/d001_get_dpc_raw_signal.py"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s_dpc = s_crop_com_threshold.correct_ramp(corner_size=0.05)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s_dpc.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Visualization"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we can visualize the signal as a magnitude and direction maps: `get_color_signal`, `get_magnitude_signal` and `get_color_image_with_indicator`.\n",
    "\n",
    "The two former returns a HyperSpy signal, while the latter interfaces directly with the matplotlib backend making it more customizable."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s_color = s_dpc.get_color_signal()\n",
    "s_color.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`get_magnitude_signal` gives the magnitude of the beam shift vector. Which can be useful for visualizing the domain walls."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s_magnitude = s_dpc.get_magnitude_signal()\n",
    "s_magnitude.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The `get_color_image_with_indicator` method has a large degree of customizability, which is useful when making images for presentations, posters or articles.\n",
    "\n",
    "By default it returns a matplotlib figure object, which can be saved directly"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig = s_dpc.get_color_image_with_indicator()\n",
    "fig.savefig(\"dpc_image.png\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "One way quick-and-dirty of suppressing the contrast from the nanocrystalline grains is by blurring.\n",
    "This suppresses the contrast from the rapidly varying grain, leading to an enhancing of the contrast from the large features, like the ferromagnetic domains."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s_dpc_blur = s_dpc.gaussian_blur(1.2)\n",
    "s_dpc_blur.get_color_signal().plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Bivariate histogram\n",
    "\n",
    "Lastly, we can visualize the 2-dimensional histogram of the beam shifts.\n",
    "\n",
    "Tip: press the H key on your keyboard, to tune the intensity in the image."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s_hist = s_dpc_blur.get_bivariate_histogram()\n",
    "s_hist.plot(cmap='viridis')"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.9.1+"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}