Merge pull request #379 from jadball/geometry_fix

Geometry fixes

ImageD11/nbGui/S3DXRD/pbp_3_refinement.ipynb

@@ -201,14 +201,13 @@
     "# set up a refinement manager object\n",
     "\n",
     "y0 = 0.0\n",
-    "fpks = 0.95\n",
     "hkl_tol_origins = 0.05\n",
     "hkl_tol_refine = 0.1\n",
     "hkl_tol_refine_merged = 0.05\n",
     "ds_tol = 0.006\n",
     "ifrac = 6e-3\n",
     "\n",
-    "refine = PBPRefine(dset=ds, y0=y0, fpks=fpks, hkl_tol_origins=hkl_tol_origins, hkl_tol_refine=hkl_tol_refine, hkl_tol_refine_merged=hkl_tol_refine_merged, ds_tol=ds_tol, ifrac=ifrac, phase_name=phase_str)"
+    "refine = PBPRefine(dset=ds, y0=y0, hkl_tol_origins=hkl_tol_origins, hkl_tol_refine=hkl_tol_refine, hkl_tol_refine_merged=hkl_tol_refine_merged, ds_tol=ds_tol, ifrac=ifrac, phase_name=phase_str)"
    ]
   },
   {
@@ -394,7 +393,7 @@
     "        pmap = PBPMap(ds.pbpfile)\n",
     "        pmap.choose_best(min_unique)\n",
     "        \n",
-    "        refine = PBPRefine(dset=ds, y0=y0, fpks=fpks, hkl_tol_origins=hkl_tol_origins, hkl_tol_refine=hkl_tol_refine, hkl_tol_refine_merged=hkl_tol_refine_merged, ds_tol=ds_tol, ifrac=ifrac, phase_name=phase_str)\n",
+    "        refine = PBPRefine(dset=ds, y0=y0, ybeam=ybeam, hkl_tol_origins=hkl_tol_origins, hkl_tol_refine=hkl_tol_refine, hkl_tol_refine_merged=hkl_tol_refine_merged, ds_tol=ds_tol, ifrac=ifrac, phase_name=phase_str)\n",
     "        \n",
     "        refine.setmap(pmap)\n",
     "        refine.setpeaks(cf_2d)\n",

ImageD11/nbGui/S3DXRD/tomo_2_map.ipynb

@@ -63,6 +63,7 @@
     "from ImageD11.sinograms.sinogram import GrainSinogram, build_slice_arrays, write_slice_recon, read_h5, write_h5, get_2d_peaks_from_4d_peaks\n",
     "from ImageD11.sinograms.roi_iradon import run_iradon\n",
     "from ImageD11.sinograms.tensor_map import TensorMap\n",
+    "from ImageD11.sinograms.geometry import sino_shift_and_pad\n",
     "import ImageD11.sinograms.dataset\n",
     "import ImageD11.nbGui.nb_utils as utils\n",
     "\n",
@@ -133,18 +134,6 @@
     "is_half_scan = False"
    ]
   },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "metadata": {
-    "tags": []
-   },
-   "outputs": [],
-   "source": [
-    "if is_half_scan:\n",
-    "    ds.correct_bins_for_half_scan()"
-   ]
-  },
   {
    "cell_type": "code",
    "execution_count": null,
@@ -318,8 +307,8 @@
     "a[3].scatter(y, x, c=s)\n",
     "a[3].set(title='Number of 4D peaks', aspect='equal')\n",
     "\n",
-    "fig.supxlabel(\"<- Lab y (transverse)\")\n",
-    "fig.supylabel(\"Lab x (beam) ->\")\n",
+    "fig.supxlabel(\"<- Sample y (transverse)\")\n",
+    "fig.supylabel(\"Sample x (beam) ->\")\n",
     "\n",
     "for a in ax.ravel():\n",
     "    a.invert_xaxis()\n",
@@ -351,16 +340,19 @@
     "\n",
     "print('y0 is', y0)\n",
     "\n",
-    "# update the y0 for each grain with the median y0:\n",
+    "if is_half_scan:\n",
+    "    ds.correct_bins_for_half_scan(y0=y0)\n",
     "\n",
-    "for gs in grainsinos:\n",
-    "    gs.update_recon_parameters(y0=y0)\n",
+    "# try to automatically determine the sinogram shift and the padding from the y0 values\n",
+    "shift, pad = sino_shift_and_pad(y0, len(ds.ybincens), min(ds.ybincens), ds.ystep)\n",
     "\n",
-    "# the shift we have to apply to the reconstructions is equal to -y0/ystep (in integer coords)\n",
+    "print('shift is', shift)\n",
+    "print('pad is', pad)\n",
     "\n",
-    "shift = -y0/ds.ystep\n",
+    "# update the grainsinogram parameters accordingly:\n",
     "\n",
-    "print('shift is', shift)"
+    "for gs in grainsinos:\n",
+    "    gs.update_recon_parameters(y0=y0, shift=shift, pad=pad)"
    ]
   },
   {
@@ -392,7 +384,6 @@
    "source": [
     "# \"quick\" whole-sample reconstruction\n",
     "\n",
-    "pad = 30\n",
     "nthreads = len(os.sched_getaffinity(os.getpid()))\n",
     "\n",
     "whole_sample_recon = run_iradon(whole_sample_sino, ds.obincens, pad, shift, workers=nthreads, apply_halfmask=is_half_scan, mask_central_zingers=is_half_scan)"
@@ -446,8 +437,8 @@
     "axs[1].set_title(\"Binarised threshold\")\n",
     "axs[2].set_title(\"Convex hull\")\n",
     "\n",
-    "fig.supxlabel(\"<-- Y axis\")\n",
-    "fig.supylabel(\"Beam >\")\n",
+    "fig.supxlabel(\"<-- Sample Y axis\")\n",
+    "fig.supylabel(\"Sample x axis (Beam) >\")\n",
     "\n",
     "plt.show()"
    ]
@@ -592,8 +583,8 @@
     "axs[0].set_title(\"ssino\")\n",
     "axs[1].imshow(gs.recons[\"iradon\"], vmin=0, origin=\"lower\")\n",
     "axs[1].set_title(\"ID11 iradon\")\n",
-    "axs[1].set_xlabel(\"<-- Lab Y\")\n",
-    "axs[1].set_ylabel(\"Lab X\")\n",
+    "axs[1].set_xlabel(\"<-- Sample Y\")\n",
+    "axs[1].set_ylabel(\"Sample X\")\n",
     "\n",
     "plt.show()"
    ]
@@ -644,8 +635,8 @@
     "fig, a = plt.subplots(1,2,figsize=(10,5))\n",
     "rec = a[0].imshow(grainsinos[0].recons[\"iradon\"], vmin=0, origin=\"lower\")\n",
     "sin = a[1].imshow(grainsinos[0].ssino, aspect='auto')\n",
-    "a[0].set_xlabel(\"<-- Lab Y\")\n",
-    "a[0].set_ylabel(\"Lab X\")\n",
+    "a[0].set_xlabel(\"<-- Sample Y\")\n",
+    "a[0].set_ylabel(\"Sample X\")\n",
     "\n",
     "# Function to update the displayed image based on the selected frame\n",
     "def update_frame(i):\n",
@@ -948,8 +939,6 @@
     "                continue\n",
     "        \n",
     "        ds.phases = ds.get_phases_from_disk()\n",
-    "        if is_half_scan:\n",
-    "            ds.correct_bins_for_half_scan()\n",
     "        \n",
     "        print(\"Peaks\")\n",
     "        # cf_4d = ds.get_cf_4d_from_disk()\n",
@@ -968,9 +957,8 @@
     "        for grain_label, gs in enumerate(grainsinos):\n",
     "            gs.update_lab_position_from_peaks(cf_strong, grain_label)\n",
     "        \n",
-    "        sample_y0s = [gs.recon_y0 for gs in grainsinos]\n",
-    "        y0 = np.median(sample_y0s)\n",
-    "        shift = -y0/ds.ystep\n",
+    "        if is_half_scan:\n",
+    "            ds.correct_bins_for_half_scan(y0)\n",
     "        \n",
     "        print(\"Determining RGB colours\")\n",
     "        utils.get_rgbs_for_grains(grains)\n",

ImageD11/nbGui/S3DXRD/tomo_2_map_minor_phase.ipynb

@@ -139,18 +139,6 @@
     "is_half_scan = False"
    ]
   },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "metadata": {
-    "tags": []
-   },
-   "outputs": [],
-   "source": [
-    "if is_half_scan:\n",
-    "    ds.correct_bins_for_half_scan()"
-   ]
-  },
   {
    "cell_type": "code",
    "execution_count": null,
@@ -300,6 +288,16 @@
     "print(shift, y0, pad)"
    ]
   },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "if is_half_scan:\n",
+    "    ds.correct_bins_for_half_scan(y0)"
+   ]
+  },
   {
    "cell_type": "code",
    "execution_count": null,
@@ -994,8 +992,7 @@
     "        \n",
     "        ds.phases = ds.get_phases_from_disk()\n",
     "        \n",
-    "        if is_half_scan:\n",
-    "            ds.correct_bins_for_half_scan()\n",
+    "\n",
     "        \n",
     "        print(\"Importing peaks\")\n",
     "        cf_4d = ds.get_cf_4d_from_disk()\n",
@@ -1011,6 +1008,9 @@
     "        whole_sample_mask = major_phase_grainsinos[0].recon_mask\n",
     "        shift = major_phase_grainsinos[0].recon_shift\n",
     "        pad = major_phase_grainsinos[0].recon_pad\n",
+    "        y0 = major_phase_grainsinos[0].recon_y0\n",
+    "        if is_half_scan:\n",
+    "            ds.correct_bins_for_half_scan(y0)\n",
     "        \n",
     "        utils.assign_peaks_to_grains(grains, cf_strong, tol=peak_assign_tol)\n",
     "        for grain_label, g in enumerate(grains):\n",

ImageD11/nbGui/S3DXRD/tomo_3_refinement.ipynb

@@ -69,7 +69,6 @@
     "from ImageD11.sinograms.point_by_point import PBPRefine\n",
     "from ImageD11.peakselect import select_ring_peaks_by_intensity\n",
     "from ImageD11.sinograms import properties, roi_iradon\n",
-    "from ImageD11.sinograms import geometry\n",
     "from ImageD11.sinograms.sinogram import GrainSinogram, build_slice_arrays, write_slice_recon, read_slice_recon, write_h5, read_h5, write_pbp_strain\n",
     "from ImageD11.grain import grain\n",
     "from ImageD11 import cImageD11\n",
@@ -231,14 +230,13 @@
     "# set up a refinement manager object\n",
     "\n",
     "y0 = grainsinos[0].recon_y0\n",
-    "fpks = 0.9\n",
     "hkl_tol_origins = 0.05\n",
     "hkl_tol_refine = 0.1\n",
     "hkl_tol_refine_merged = 0.05\n",
     "ds_tol = 0.004\n",
     "ifrac = 7e-3\n",
     "\n",
-    "refine = PBPRefine(dset=ds, y0=y0, fpks=fpks, hkl_tol_origins=hkl_tol_origins, hkl_tol_refine=hkl_tol_refine, hkl_tol_refine_merged=hkl_tol_refine_merged, ds_tol=ds_tol, ifrac=ifrac, phase_name=phase_str)"
+    "refine = PBPRefine(dset=ds, y0=y0, hkl_tol_origins=hkl_tol_origins, hkl_tol_refine=hkl_tol_refine, hkl_tol_refine_merged=hkl_tol_refine_merged, ds_tol=ds_tol, ifrac=ifrac, phase_name=phase_str)"
    ]
   },
   {
@@ -442,7 +440,7 @@
     "        pmap = tensor_map.to_pbpmap(z_layer=0, default_npks=20, default_nuniq=20)\n",
     "        pmap.choose_best(1)\n",
     "\n",
-    "        refine = PBPRefine(dset=ds, y0=y0, fpks=fpks, hkl_tol_origins=hkl_tol_origins, hkl_tol_refine=hkl_tol_refine, hkl_tol_refine_merged=hkl_tol_refine_merged, ds_tol=ds_tol, ifrac=ifrac, phase_name=phase_str)\n",
+    "        refine = PBPRefine(dset=ds, y0=y0, hkl_tol_origins=hkl_tol_origins, hkl_tol_refine=hkl_tol_refine, hkl_tol_refine_merged=hkl_tol_refine_merged, ds_tol=ds_tol, ifrac=ifrac, phase_name=phase_str)\n",
     "        \n",
     "        refine.setmap(pmap)\n",
     "        refine.setpeaks(cf_2d)\n",

ImageD11/nbGui/S3DXRD/tomo_3_refinement_minor_phase.ipynb

@@ -69,7 +69,6 @@
     "from ImageD11.sinograms.point_by_point import PBPRefine\n",
     "from ImageD11.peakselect import select_ring_peaks_by_intensity\n",
     "from ImageD11.sinograms import properties, roi_iradon\n",
-    "from ImageD11.sinograms import geometry\n",
     "from ImageD11.sinograms.sinogram import GrainSinogram, build_slice_arrays, write_slice_recon, read_slice_recon, write_h5, read_h5, write_pbp_strain\n",
     "from ImageD11.grain import grain\n",
     "from ImageD11 import cImageD11\n",
@@ -319,7 +318,6 @@
     "# set up a refinement manager object\n",
     "\n",
     "y0 = grainsinos[0].recon_y0\n",
-    "fpks = 0.9\n",
     "hkl_tol_origins = 0.05\n",
     "hkl_tol_refine = 0.1\n",
     "hkl_tol_refine_merged = 0.05\n",
@@ -328,7 +326,7 @@
     "forref = [0, 2, 3, 4, 5, 6, 7, 8, 10, 12, 13]\n",
     "phase_str = minor_phase_str\n",
     "\n",
-    "refine = PBPRefine(dset=ds, y0=y0, fpks=fpks, hkl_tol_origins=hkl_tol_origins, hkl_tol_refine=hkl_tol_refine, hkl_tol_refine_merged=hkl_tol_refine_merged, ds_tol=ds_tol, ifrac=ifrac, phase_name=phase_str, forref=forref)"
+    "refine = PBPRefine(dset=ds, y0=y0, hkl_tol_origins=hkl_tol_origins, hkl_tol_refine=hkl_tol_refine, hkl_tol_refine_merged=hkl_tol_refine_merged, ds_tol=ds_tol, ifrac=ifrac, phase_name=phase_str, forref=forref)"
    ]
   },
   {
@@ -562,7 +560,7 @@
     "        pmap = tensor_map.to_pbpmap(z_layer=0, default_npks=20, default_nuniq=20)\n",
     "        pmap.choose_best(1)\n",
     "\n",
-    "        refine = PBPRefine(dset=ds, y0=y0, fpks=fpks, hkl_tol_origins=hkl_tol_origins, hkl_tol_refine=hkl_tol_refine, hkl_tol_refine_merged=hkl_tol_refine_merged, ds_tol=ds_tol, ifrac=ifrac, phase_name=minor_phase_str, forref=forref)\n",
+    "        refine = PBPRefine(dset=ds, y0=y0, hkl_tol_origins=hkl_tol_origins, hkl_tol_refine=hkl_tol_refine, hkl_tol_refine_merged=hkl_tol_refine_merged, ds_tol=ds_tol, ifrac=ifrac, phase_name=minor_phase_str, forref=forref)\n",
     "        \n",
     "        refine.own_filename = os.path.splitext(refine.own_filename)[0] + f'_{phase_str}.h5'\n",
     "        refine.icolf_filename = os.path.splitext(refine.icolf_filename)[0] + f'_{phase_str}.h5'\n",