{ "cells": [ { "cell_type": "markdown", "id": "31e30c68-b19e-4b0b-a1a9-9174eef5b6d3", "metadata": { "nbsphinx": "hidden", "tags": [] }, "source": [ "This notebook is part of the *orix* documentation https://orix.readthedocs.io. Links to the documentation won’t work from the notebook." ] }, { "cell_type": "markdown", "id": "b82f51da-5413-4e00-bfed-c6f53f7502ee", "metadata": {}, "source": [ "# Inverse pole figures\n", "\n", "In this tutorial we will express in which crystal direction $\\mathbf{t} = \\left$ any sample direction $\\mathbf{v} = (x, y, z)$, or any [Vector3d](../reference/generated/orix.vector.Vector3d.rst), is parallel to.\n", "Formally, this so-called inverse pole figure (IPF) shows the directional distribution of a set of reference vectors in terms of a fixed crystallographic frame Nolze (2015).\n", "By assigning a unique colour to each crystal direction in a symmetry's fundamental sector (the IPF), we can colour orientations Nolze and Hielscher (2016).\n", "\n", "While orix has descriptions for all 32 crystallographic point groups $S$, so far only the fundamental sectors of the eleven Laue group symmetries have been defined.\n", "Thus, the Laue class (in the left and center column in the table below) of the point groups (center and right column) will be used when plotting the IPF and colouring orientations.\n", "\n", "| Schoenflies | Laue | Non-centrosymmetric |\n", "| ----------- | --------------- | -------------------------- |\n", "| $Ci$ | $\\bar{1}$ | $1$ |\n", "| $C2h$ | $2/m$ | $2$, $m$ |\n", "| $D2h$ | $mmm$ | $222$, $2mm$ |\n", "| $S6$ | $\\bar{3}$ | $3$ |\n", "| $D3d$ | $\\bar{3}m$ | $32$, $3m$ |\n", "| $C4h$ | $4/m$ | $4$, $\\bar{4}$ |\n", "| $D4h$ | $4/mmm$ | $422$, $\\bar{4}2m$, $4mm$ |\n", "| $C6h$ | $6/m$ | $6$, $\\bar{6}$ |\n", "| $D6h$ | $6/mmm$ | $622$, $\\bar{6}2m$, $6mmm$ |\n", "| $Th$ | $m\\bar{3}$ | $23$ |\n", "| $Oh$ | $m\\bar{3}m$ | $432$, $\\bar{4}3m$ |\n", "\n", "Only the colour keys used in EDAX TSL OIM Analysis are supported, via the [IPFColorKeyTSL](../reference/generated/orix.plot.IPFColorKeyTSL.rst).\n", "\n", "The IPF functionality in the excellent software [MTEX](https://mtex-toolbox.github.io/Documentation.html) and the two above referenced works have been the basis for the IPF functionality in orix." ] }, { "cell_type": "code", "execution_count": null, "id": "35bbfb73-8d49-4d24-bd88-cf7e81592f9e", "metadata": {}, "outputs": [], "source": [ "%matplotlib inline\n", "\n", "import matplotlib.pyplot as plt\n", "import numpy as np\n", "\n", "from orix import plot, sampling\n", "from orix.crystal_map import Phase\n", "from orix.quaternion import Orientation, symmetry\n", "from orix.vector import Vector3d\n", "\n", "\n", "# We'll want our plots to look a bit larger than the default size\n", "new_params = {\n", " \"figure.facecolor\": \"w\",\n", " \"figure.figsize\": (20, 7),\n", " \"lines.markersize\": 10,\n", " \"font.size\": 15,\n", " \"axes.grid\": True,\n", "}\n", "plt.rcParams.update(new_params)" ] }, { "cell_type": "markdown", "id": "790fa928-6db3-4515-b844-bb735c98c579", "metadata": {}, "source": [ "## Minimal example of plotting the IPF colour key" ] }, { "cell_type": "code", "execution_count": null, "id": "63ae1026-b81d-4c60-836a-9f2b7647ddda", "metadata": { "nbsphinx-thumbnail": { "tooltip": "Sample directions expressed in terms of the crystallographic frame" }, "tags": [ "nbsphinx-thumbnail" ] }, "outputs": [], "source": [ "plot.IPFColorKeyTSL(symmetry.D3d).plot()" ] }, { "cell_type": "markdown", "id": "e60892ca-5425-473b-a1c2-e84fce94810d", "metadata": {}, "source": [ "## Plot colour keys\n", "\n", "Here we plot the colour keys for the eleven Laue groups.\n", "An example demonstrating how an IPF color key can be added to a custom plot is shown in the [Crystal map tutorial](crystal_map.ipynb#Plotting)." ] }, { "cell_type": "code", "execution_count": null, "id": "357988b9-c58d-4f76-a893-ec24f85088b8", "metadata": {}, "outputs": [], "source": [ "pg_laue = [\n", " symmetry.Ci,\n", " symmetry.C2h,\n", " symmetry.D2h,\n", " symmetry.S6,\n", " symmetry.D3d,\n", " symmetry.C4h,\n", " symmetry.D4h,\n", " symmetry.C6h,\n", " symmetry.D6h,\n", " symmetry.Th,\n", " symmetry.Oh,\n", "]" ] }, { "cell_type": "code", "execution_count": null, "id": "a9069847-9e5e-4a89-a359-a077fc43874e", "metadata": {}, "outputs": [], "source": [ "for pg in pg_laue:\n", " ipfkey = plot.IPFColorKeyTSL(pg)\n", " ipfkey.plot()" ] }, { "cell_type": "markdown", "id": "2d39a57c-7363-4a0e-a8cb-dccb64ed94c7", "metadata": {}, "source": [ "## Scatter plot\n", "\n", "Here we show which crystal direction $\\left$ the sample reference frame directions $(x, y, Z)$ point in.\n", "We can turn off the grid by updating the Matplotlib parameters with `plt.rcParams[\"axes.grid\"] = False`, or passing `grid=False` to [Orientation.scatter()](../reference/generated/orix.quaternion.Orientation.scatter.rst), which we use to plot the IPFs." ] }, { "cell_type": "code", "execution_count": null, "id": "f88fb7f1-51a8-4238-977b-d39d87eba01a", "metadata": {}, "outputs": [], "source": [ "v = Vector3d([[1, 0, 0], [0, 1, 0], [0, 0, 1]])\n", "kwargs = {\"projection\": \"ipf\", \"direction\": v}\n", "\n", "for S in pg_laue:\n", " O = Orientation.from_euler([325, 48, 163], S, degrees=True)\n", " O.scatter(**kwargs)" ] }, { "cell_type": "markdown", "id": "9e502318-9bb3-44da-a707-a739b8858ca9", "metadata": {}, "source": [ "## Coloring orientations\n", "\n", "Here we show how to obtain an RGB colour for each crystal orientation given a sample direction, using [IPFColorKeyTSL.orientation2color](../reference/generated/orix.plot.IPFColorKeyTSL.orientation2color.rst).\n", "We then plot them in the IPF. Note that we could instead plot them in a map, say when the orientations come from a crystallographic mapping experiment.\n", "See the [crystal map](crystal_map.rst) tutorial for an example." ] }, { "cell_type": "code", "execution_count": null, "id": "64e454e2-34a8-4e45-b704-925e79e36e97", "metadata": {}, "outputs": [], "source": [ "plt.rcParams[\"axes.grid\"] = False\n", "\n", "O2 = Orientation.random(1000)\n", "for S in pg_laue:\n", " ipfkey = plot.IPFColorKeyTSL(S)\n", " O2.symmetry = ipfkey.symmetry\n", " rgb_z = ipfkey.orientation2color(O2)\n", " O2.scatter(\"ipf\", c=rgb_z, direction=ipfkey.direction)" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "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.11.6" } }, "nbformat": 4, "nbformat_minor": 5 }