Skip to content

Commit

Permalink
Merge pull request #242 from jadball/scanbooks
Browse files Browse the repository at this point in the history 
Optionally scale sinogram rows by ring current
  • Loading branch information
@jadball
jadball authored Mar 6, 2024
2 parents c0674b2 + b145eff commit fe74de9
Show file tree
Hide file tree
Showing 3 changed files with 100 additions and 13 deletions.
43 changes: 32 additions & 11 deletions ImageD11/nbGui/S3DXRD/2_S3DXRD_sinograms_map.ipynb
View file
Original file line number Diff line number Diff line change
Expand Up @@ -141,6 +141,7 @@
" order = np.lexsort((np.arange(npks), sinoangles))\n",
" sinoangles = sinoangles[order]\n",
" ssino = sino[order].T\n",
" \n",
" return sinoangles, ssino, hits[order].T\n",
"\n",
"def do_sinos(g, hkltol=0.25):\n",
Expand Down Expand Up @@ -409,6 +410,17 @@
"ds = ImageD11.sinograms.dataset.load(dset_path)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# determine ring currents for sinogram row-by-row intensity correction\n",
"\n",
"utils.get_ring_current_per_scan(ds)"
]
},
{
"cell_type": "code",
"execution_count": null,
Expand Down Expand Up @@ -710,17 +722,6 @@
"utils.plot_grain_sinograms(grains, cf_strong)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"tags": []
},
"outputs": [],
"source": [
"np.sum(grains[0].mask_4d)"
]
},
{
"cell_type": "code",
"execution_count": null,
Expand Down Expand Up @@ -842,6 +843,18 @@
" pass"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# we can optionally correct the grain sinograms by scaling each row by the ring current:\n",
"\n",
"for grain in tqdm(grains):\n",
" utils.correct_sinogram_rows_with_ring_current(grain, ds)"
]
},
{
"cell_type": "code",
"execution_count": null,
Expand Down Expand Up @@ -1241,6 +1254,7 @@
"cf_strong_dstol = 0.01\n",
"\n",
"is_half_scan = False\n",
"correct_sinos_with_ring_current = True\n",
"\n",
"peak_assign_tol = 0.25\n",
"\n",
Expand Down Expand Up @@ -1281,6 +1295,9 @@
" if \"slice_recon\" in hin.keys():\n",
" print(f\"Already reconstructed {dataset} in {sample}, skipping\")\n",
" continue\n",
" \n",
" # determine ring currents for sinogram row-by-row intensity correction\n",
" utils.get_ring_current_per_scan(ds)\n",
" \n",
" cf_4d = ImageD11.columnfile.columnfile(ds.col4dfile)\n",
" cf_4d.parameters.loadparameters(par_path)\n",
Expand Down Expand Up @@ -1364,6 +1381,10 @@
" for i in tqdm(pool.map(do_sinos, grains), total=len(grains)):\n",
" pass\n",
" \n",
" if correct_sinos_with_ring_current:\n",
" for grain in tqdm(grains):\n",
" utils.correct_sinogram_rows_with_ring_current(grain, ds)\n",
" \n",
" print(\"Running iradon\")\n",
" \n",
" run_this_iradon = partial(run_iradon_id11, pad=pad, y0=y0, sample_mask=whole_sample_mask, workers=1, apply_halfmask=is_half_scan, mask_central_zingers=is_half_scan)\n",
Expand Down
51 changes: 49 additions & 2 deletions ImageD11/nbGui/S3DXRD/2_S3DXRD_sinograms_map_minor_phase.ipynb
View file
Original file line number Diff line number Diff line change
Expand Up @@ -161,6 +161,7 @@
" order = np.lexsort((np.arange(npks), sinoangles))\n",
" sinoangles = sinoangles[order]\n",
" ssino = sino[order].T\n",
" \n",
" return sinoangles, ssino, hits[order].T\n",
"\n",
"def do_sinos(g, hkltol=0.25):\n",
Expand Down Expand Up @@ -234,8 +235,23 @@
" blobs_sorted = sorted(blobs, key=lambda x: x[2], reverse=True)\n",
" try:\n",
" largest_blob = blobs_sorted[0]\n",
" grain.x_blob = largest_blob[1]\n",
" grain.y_blob = largest_blob[0]\n",
" \n",
" # we now have the blob position in recon space\n",
" # we need to go back to microns\n",
" \n",
" # first axis (vertical) is x\n",
" # second axis (horizontal) is y\n",
" \n",
" x_recon_space = largest_blob[0]\n",
" y_recon_space = largest_blob[1]\n",
" \n",
" # centre of the recon image is centre of space\n",
" \n",
" x_microns = (x_recon_space - grain.recon.shape[0]//2) * ds.ystep + y0\n",
" y_microns = -(y_recon_space - grain.recon.shape[1]//2) * ds.ystep + y0\n",
" \n",
" grain.x_blob = x_microns\n",
" grain.y_blob = y_microns\n",
" except IndexError:\n",
" # didn't find any blobs\n",
" # for small grains like these, if we didn't find a blob, normally indicates recon is bad\n",
Expand Down Expand Up @@ -377,6 +393,17 @@
"ds = ImageD11.sinograms.dataset.load(dset_path)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# determine ring currents for sinogram row-by-row intensity correction\n",
"\n",
"utils.get_ring_current_per_scan(ds)"
]
},
{
"cell_type": "code",
"execution_count": null,
Expand Down Expand Up @@ -714,6 +741,18 @@
" pass"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# we can optionally correct the grain sinograms by scaling each row by the ring current:\n",
"\n",
"for grain in tqdm(grains):\n",
" utils.correct_sinogram_rows_with_ring_current(grain, ds)"
]
},
{
"cell_type": "code",
"execution_count": null,
Expand Down Expand Up @@ -1137,6 +1176,7 @@
"cf_strong_dstol = 0.005\n",
"\n",
"is_half_scan = False\n",
"correct_sinos_with_ring_current = True\n",
"\n",
"peak_assign_tol = 0.25\n",
"\n",
Expand Down Expand Up @@ -1172,6 +1212,9 @@
" if \"slice_recon\" in hin.keys():\n",
" print(f\"Already reconstructed {dataset} in {sample}, skipping\")\n",
" continue\n",
" \n",
" # determine ring currents for sinogram row-by-row intensity correction\n",
" utils.get_ring_current_per_scan(ds)\n",
" \n",
" cf_4d = ImageD11.columnfile.columnfile(ds.col4dfile)\n",
" cf_4d.parameters.loadparameters(par_path)\n",
Expand Down Expand Up @@ -1236,6 +1279,10 @@
" for i in tqdm(pool.map(do_sinos, grains), total=len(grains)):\n",
" pass\n",
" \n",
" if correct_sinos_with_ring_current:\n",
" for grain in tqdm(grains):\n",
" utils.correct_sinogram_rows_with_ring_current(grain, ds)\n",
" \n",
" print(\"Running iradon\")\n",
" \n",
" run_this_iradon = partial(run_iradon_id11, pad=pad, y0=y0, sample_mask=whole_sample_mask, workers=1, apply_halfmask=is_half_scan, mask_central_zingers=is_half_scan)\n",
Expand Down
19 changes: 19 additions & 0 deletions ImageD11/nbGui/nb_utils.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2,6 +2,7 @@
import subprocess
import time

import h5py
import numba
import numpy as np
from matplotlib import pyplot as plt
Expand All @@ -20,6 +21,24 @@
from scipy.optimize import curve_fit



def correct_sinogram_rows_with_ring_current(grain, ds):
grain.ssino = grain.ssino / ds.ring_currents_per_scan_scaled[:, None]


def get_ring_current_per_scan(ds):
"""Gets the ring current for each scan (i.e rotation/y-step)
Stores it inside ds.ring_currents_per_scan and a scaled version inside ds.ring_currents_per_scan_scaled"""
if not hasattr(ds, "ring_currents_per_scan"):
ring_currents = []
with h5py.File(ds.masterfile, "r") as h5in:
for scan in ds.scans:
ring_current = float(h5in[scan]["instrument/machine/current"][()])
ring_currents.append(ring_current)

ds.ring_currents_per_scan = np.array(ring_currents)
ds.ring_currents_per_scan_scaled = np.array(ring_currents/np.max(ring_currents))

def correct_pixel(pixel, spline_file):
sr, fr = pixel
sc, fc = ImageD11.blobcorrector.correctorclass(spline_file).correct(sr, fr)
Expand Down

0 comments on commit fe74de9

Please sign in to comment.