Merge 8854090 into fb07613

pyxem/utils/expt_utils.py

@@ -102,6 +102,16 @@ def _polar2cart(r, theta):
     x = r * np.sin(theta)
     return x, y
 
+def _capped_gaussian(x,prefactor,sigma,center):
+    """
+    Worker functions for : subpixel_beam_finder()
+    Creates a 1D Gaussian model, and then truncates all values
+    greate than 1 to unity. There is no discretisation here.
+    """
+    model_value = prefactor*np.exp(-(x-center)**2/(2*sigma**2))
+    model_value[np.greater(model_value,np.ones_like(model_value))] = 1 #finite dynamic range
+    return np.ravel(model_value) #convert this to 1D
+
 def radial_average(z, center,cython=True):
     """Calculate the radial profile by azimuthal averaging about a specified
     center.
@@ -459,6 +469,34 @@ def refine_beam_position(z, start, radius):
 
     return c
 
+def subpixel_beam_finder(single_pattern,half_shape):
+    """
+    This routine is designed to find, to sub-pixel accuracy, the 
+    center of a pattern that has saturated a detector. The noise
+    model is Gaussian. The input dp should be approximately (within about 3 
+    pixels) centered. Use one of the other centering methods to achieve this.
+    
+    single_pattern : numpy array : A dp (for .map purposes)
+    half_shape     : int         : An int for the approx center of the pattern
+    """
+
+    # Set up 
+    hs = half_shape #readability
+    size = np.int(half_shape/7) #this prevents fitting to the far out noise
+    pattern = single_pattern[hs-size:hs+size,hs-size:hs+size]
+
+    if np.max(pattern) > 1:
+        raise ValueError('Patterns should be normalised to max intensity')
+
+    #fitting
+    from scipy.optimize import curve_fit as cf
+    i_array = np.array([np.arange(hs-size,hs+size),np.arange(hs-size,hs+size),np.arange(hs-size,hs+size)])
+    x_center = cf(_capped_gaussian,i_array,np.ravel(pattern[size:size+3,:].T),p0=[2,5,hs])[0][2]
+    y_center = cf(_capped_gaussian,i_array,np.ravel(pattern[:,size:size+3].T),p0=[2,5,hs])[0][2]
+
+    return np.asarray([x_center,y_center])
+
+
 def enhance_gauss_sauvola(z, sigma_blur, sigma_enhance, k, window_size, threshold, morph_opening=True):
     z = z.astype(np.float64)
     im1 = ndi.gaussian_filter(z, sigma=sigma_blur, mode='mirror')

pyxem_examples/marker_attachment_testbook.ipynb

@@ -2,7 +2,7 @@
  "cells": [
   {
    "cell_type": "code",
-   "execution_count": 43,
+   "execution_count": 2,
    "metadata": {},
    "outputs": [],
    "source": [
@@ -13,34 +13,25 @@
     "import pymatgen as pmg\n",
     "from matplotlib import pyplot as plt\n",
     "from pymatgen.transformations.standard_transformations import RotationTransformation\n",
-    "from pyxem.indexation_generator import IndexationGenerator\n",
+    "from pyxem.generators.indexation_generator import IndexationGenerator\n",
     "from pyxem.utils.sim_utils import peaks_from_best_template\n",
     "from pyxem.utils.plot import generate_marker_inputs_from_peaks\n",
     "from scipy.constants import pi"
    ]
   },
   {
    "cell_type": "code",
-   "execution_count": 44,
+   "execution_count": 3,
    "metadata": {},
-   "outputs": [
-    {
-     "name": "stdout",
-     "output_type": "stream",
-     "text": [
-      "The autoreload extension is already loaded. To reload it, use:\n",
-      "  %reload_ext autoreload\n"
-     ]
-    }
-   ],
+   "outputs": [],
    "source": [
     "%load_ext autoreload\n",
     "%autoreload 2"
    ]
   },
   {
    "cell_type": "code",
-   "execution_count": 45,
+   "execution_count": 4,
    "metadata": {
     "collapsed": true
    },
@@ -57,15 +48,13 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 46,
-   "metadata": {
-    "collapsed": true
-   },
+   "execution_count": 6,
+   "metadata": {},
    "outputs": [],
    "source": [
     "size = 256\n",
     "radius=1.5\n",
-    "ediff = pxm.ElectronDiffractionCalculator(300., 0.025)\n",
+    "ediff = pxm.DiffractionGenerator(300., 0.025)\n",
     "\n",
     "rotaxis = [0, 0, 1]\n",
     "thetas = np.arange(0, 46, 45)\n",
@@ -87,7 +76,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 47,
+   "execution_count": 7,
    "metadata": {
     "collapsed": true
    },
@@ -102,38 +91,39 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 48,
+   "execution_count": 24,
    "metadata": {},
    "outputs": [
     {
      "name": "stderr",
      "output_type": "stream",
      "text": [
-      "                                                 \r"
+      "                                                  \r"
      ]
     }
    ],
    "source": [
-    "edc = pxm.ElectronDiffractionCalculator(300, 0.025)\n",
-    "diff_gen = pxm.DiffractionLibraryGenerator(edc)\n",
+    "diff_gen = pxm.DiffractionLibraryGenerator(ediff)\n",
     "struc_lib = dict()\n",
     "struc_lib['si'] = (silicon, rot_list)\n",
     "struc_lib['ga'] = (gall,rot_list)\n",
     "library = diff_gen.get_diffraction_library(struc_lib,\n",
     "                                            calibration=1.2/128,\n",
-    "                                            reciprocal_radius=1.5,\n",
-    "                                            representation='euler')"
+    "                                            reciprocal_radius=1.1,\n",
+    "                                            representation='euler',\n",
+    "                                            half_shape=(256/2,256/2),\n",
+    "                                            with_direct_beam=False)"
    ]
   },
   {
    "cell_type": "code",
-   "execution_count": 49,
+   "execution_count": 25,
    "metadata": {},
    "outputs": [
     {
      "data": {
       "application/vnd.jupyter.widget-view+json": {
-       "model_id": "675d0008b8d348db87981eb4292e9051",
+       "model_id": "80070396142b40ae945fe0344323e56e",
        "version_major": 2,
        "version_minor": 0
       },
@@ -143,13 +133,6 @@
      },
      "metadata": {},
      "output_type": "display_data"
-    },
-    {
-     "name": "stdout",
-     "output_type": "stream",
-     "text": [
-      "\n"
-     ]
     }
    ],
    "source": [
@@ -160,13 +143,13 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 50,
+   "execution_count": 26,
    "metadata": {},
    "outputs": [
     {
      "data": {
       "application/vnd.jupyter.widget-view+json": {
-       "model_id": "69d040b131d74bbd8c103ce0370d1f12",
+       "model_id": "b6d864e1bad44b5b9b500d113e72c927",
        "version_major": 2,
        "version_minor": 0
       },
@@ -176,13 +159,6 @@
      },
      "metadata": {},
      "output_type": "display_data"
-    },
-    {
-     "name": "stdout",
-     "output_type": "stream",
-     "text": [
-      "\n"
-     ]
     }
    ],
    "source": [
@@ -191,7 +167,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 51,
+   "execution_count": 27,
    "metadata": {
     "collapsed": true
    },
@@ -222,7 +198,7 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython3",
-   "version": "3.6.2"
+   "version": "3.6.3"
   }
  },
  "nbformat": 4,

tests/test_diffraction_simulation.py

@@ -84,15 +84,14 @@ def check_pattern_equivilance(p1,p2,coords_only=False):
 fake_pattern   = get_pattern(Mscope,fake_cubic_I)
 larger_pattern = get_pattern(Mscope,larger_cubic_I)
 
-
 def test_casual_formal():
     # Checks that Pymatgen understands that these are the same structure
     assert formal_cubic_I == casual_cubic_I
 
 def test_casual_formal_in_simulation():
     ## Checks that are simulations also realise that
     check_pattern_equivilance(formal_pattern,casual_pattern)
-    
+
 def test_systematic_absence():
     ## Cubic I thus each peak must have indicies that sum to an even number
     assert np.all(np.sum(formal_pattern.indices,axis=1) % 2 == 0)