Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

ImageD11_gui fixes and some scanning notebooks updates #226

Merged
merged 2 commits into from
Feb 26, 2024
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
218 changes: 41 additions & 177 deletions ImageD11/nbGui/S3DXRD/1_S3DXRD_index_z_slice_minor_phase.ipynb
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,7 @@
"source": [
"# Jupyter notebook based on ImageD11 to process scanning 3DXRD data\n",
"# Written by Haixing Fang, Jon Wright and James Ball\n",
"## Date: 19/02/2024"
"## Date: 21/02/2024"
]
},
{
Expand Down Expand Up @@ -53,7 +53,7 @@
"%matplotlib widget\n",
"from matplotlib import pyplot as plt\n",
"\n",
"import utils\n",
"import ImageD11.nbGui.nb_utils as utils\n",
"\n",
"import ImageD11.grain\n",
"import ImageD11.unitcell\n",
Expand All @@ -76,21 +76,19 @@
},
"outputs": [],
"source": [
"# OLD DATASETS\n",
"\n",
"# NOTE: For old datasets before the new directory layout structure, we don't distinguish between RAW_DATA and PROCESSED_DATA\n",
"# In this case, use this cell to specify where your experimental folder is, and do not run the cell below\n",
"# e.g /data/visitor/4752/id11/20210513\n",
"# e.g /data/visitor/ma4752/id11/20210513\n",
"\n",
"### USER: specify your experimental directory\n",
"\n",
"rawdata_path = \"/home/esrf/james1997a/Data/ma4752/id11/20210618\"\n",
"rawdata_path = \"/home/esrf/james1997a/Data/ihma439/id11/20231211/RAW_DATA\"\n",
"\n",
"!ls -lrt {rawdata_path}\n",
"\n",
"### USER: specify where you want your processed data to go\n",
"\n",
"processed_data_root_dir = \"/home/esrf/james1997a/Data/ma4752/id11/20240118/James\""
"processed_data_root_dir = \"/home/esrf/james1997a/Data/ihma439/id11/20231211/PROCESSED_DATA/James/20240221\""
]
},
{
Expand All @@ -103,8 +101,8 @@
"source": [
"# USER: pick a sample and a dataset you want to segment\n",
"\n",
"sample = \"MA4752_S4_2_XRD\"\n",
"dataset = \"DTL1z90\""
"sample = \"FeAu_0p5_tR_nscope\"\n",
"dataset = \"top_100um\""
]
},
{
Expand All @@ -119,8 +117,12 @@
"\n",
"dset_path = os.path.join(processed_data_root_dir, sample, f\"{sample}_{dataset}\", f\"{sample}_{dataset}_dataset.h5\")\n",
"\n",
"e2dx_path = \"/data/id11/nanoscope/Eiger/spatial_20210415_JW/e2dx.edf\"\n",
"e2dy_path = \"/data/id11/nanoscope/Eiger/spatial_20210415_JW/e2dy.edf\""
"# USER: specify the path to the parameter file\n",
"\n",
"par_path = os.path.join(processed_data_root_dir, 'Fe_refined.par')\n",
"\n",
"e2dx_path = os.path.join(processed_data_root_dir, '../../CeO2/e2dx_E-08-0173_20231127.edf')\n",
"e2dy_path = os.path.join(processed_data_root_dir, '../../CeO2/e2dy_E-08-0173_20231127.edf')"
]
},
{
Expand All @@ -139,17 +141,6 @@
"print(ds.shape)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"tags": []
},
"outputs": [],
"source": [
"par_path = 'nickel.par'"
]
},
{
"cell_type": "code",
"execution_count": null,
Expand Down Expand Up @@ -345,61 +336,32 @@
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"tags": []
},
"metadata": {},
"outputs": [],
"source": [
"# now we are indexing!\n",
"\n",
"# indexing will select all rings with a multiplicity below max_multiplity to search\n",
"max_multiplicity = 11\n",
"\n",
"min_counts_on_ring = 30\n",
"\n",
"n_peaks_expected = 0\n",
"rings = []\n",
"for i, dstar in enumerate(indexer.unitcell.ringds):\n",
" multiplicity = len(indexer.unitcell.ringhkls[indexer.unitcell.ringds[i]])\n",
" counts_on_this_ring = (indexer.ra == i).sum()\n",
" if counts_on_this_ring > min_counts_on_ring:\n",
" n_peaks_expected += multiplicity\n",
" if multiplicity < max_multiplicity:\n",
" rings.append((counts_on_this_ring, multiplicity, i))\n",
"\n",
"rings.sort()\n",
"\n",
"print(f\"{n_peaks_expected} peaks expected\")\n",
"print(f\"Trying these rings (counts, multiplicity, ring number): {rings}\")\n",
"\n",
"# USER: specify the HKL tolerances you want to use for indexing\n",
"hkl_tols_seq = [0.01, 0.02, 0.03]\n",
"\n",
"# USER: specify the fraction of the total expected peaks\n",
"fracs = [0.9, 0.8, 0.7]\n",
"\n",
"# ImageD11.cImageD11.cimaged11_omp_set_num_threads(1)\n",
"ImageD11.indexing.loglevel=3\n",
"\n",
"# indexer.uniqueness = 0.3\n",
"indexer.cosine_tol = np.cos(np.radians(90.25))\n",
"indexer.max_grains = 1000\n",
"\n",
"# iterate over HKL tolerances\n",
"for frac in fracs:\n",
" for tol in hkl_tols_seq:\n",
" indexer.minpks = n_peaks_expected*frac\n",
" indexer.hkl_tol = tol\n",
" \n",
" # iterate over rings\n",
" for i in range(len(rings)):\n",
" for j in range(i, len(rings)):\n",
" indexer.ring_1 = rings[i][2]\n",
" indexer.ring_2 = rings[j][2]\n",
" \n",
" indexer.find()\n",
" indexer.scorethem()\n",
" \n",
" print(frac, tol, len(indexer.ubis))"
"# the minimum number of peaks on a ring for a ring to be indexed on\n",
"min_counts_on_ring = 0\n",
"# the sequence of hkl tolerances the indexer will iterate through\n",
"hkl_tols_seq = [0.01, 0.02, 0.03, 0.04, 0.05, 0.06, 0.1]\n",
"# the sequence of minpks fractions the indexer will iterate through\n",
"fracs = [0.9, 0.8, 0.7, 0.6, 0.5]\n",
"# the tolerance in g-vector angle\n",
"cosine_tol = np.cos(np.radians(90.25))\n",
"# the max number of UBIs we can find per pair of rings\n",
"max_grains = 1000\n",
"\n",
"grains, indexer = utils.do_index(cf=cf_strong,\n",
" dstol=indexer.ds_tol,\n",
" max_mult=max_multiplicity,\n",
" min_ring_count=min_counts_on_ring,\n",
" hkl_tols=hkl_tols_seq,\n",
" fracs=fracs,\n",
" cosine_tol=cosine_tol,\n",
" max_grains=max_grains\n",
")"
]
},
{
Expand All @@ -410,9 +372,6 @@
},
"outputs": [],
"source": [
"# create grain objects\n",
"grains = [ImageD11.grain.grain(ubi, translation=np.array([0., 0., 0.])) for ubi in indexer.ubis]\n",
"\n",
"# set grain GIDs (useful if we ever delete a grain)\n",
"for i, g in enumerate(grains):\n",
" g.gid = i\n",
Expand Down Expand Up @@ -453,127 +412,32 @@
"\n",
"tol=0.05\n",
"\n",
"# column to store the grain labels\n",
"labels = np.zeros(cf_strong.nrows, 'i')\n",
"# get all g-vectors from columnfile\n",
"gv = np.transpose((cf_strong.gx,cf_strong.gy,cf_strong.gz)).astype(float)\n",
"# column to store drlv2 (error in hkl)\n",
"drlv2 = np.ones(cf_strong.nrows, 'd')\n",
"# iterate over all grains\n",
"print(f\"Scoring and assigning {len(grains)} grains\")\n",
"for g in tqdm(grains):\n",
" n = ImageD11.cImageD11.score_and_assign(g.ubi, gv, tol, drlv2, labels, g.gid)\n",
"\n",
"# add the labels column to the columnfile\n",
"cf_strong.addcolumn(labels, 'grain_id')\n",
"utils.assign_peaks_to_grains(grains, cf_strong, tol=tol)\n",
"\n",
"print(\"Storing peak data in grains\")\n",
"# iterate through all the grains\n",
"for g in tqdm(grains):\n",
" # store this grain's peak indices so we know which 4D peaks we used for indexing\n",
" g.mask_4d = cf_strong.grain_id == g.gid\n",
" g.peaks_4d = cf_strong.index[cf_strong.grain_id == g.gid]"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"tags": []
},
"metadata": {},
"outputs": [],
"source": [
"def plot_index_results(ind, colfile, title):\n",
" # Generate a histogram of |drlv| for a ubi matrix\n",
" indexer.histogram_drlv_fit()\n",
" indexer.fight_over_peaks()\n",
" \n",
" fig, axs = plt.subplots(3, 2, layout=\"constrained\", figsize=(9,12))\n",
" axs_flat = axs.ravel()\n",
" \n",
" # For each grain, plot the error in hkl vs the number of peaks with that error\n",
" \n",
" for grh in ind.histogram:\n",
" axs_flat[0].plot(ind.bins[1:-1], grh[:-1], \"-\")\n",
" \n",
" axs_flat[0].set(ylabel=\"number of peaks\",\n",
" xlabel=\"error in hkl (e.g. hkl versus integer)\",\n",
" title=title)\n",
" \n",
" # set a mask of all non-assigned g-vectors\n",
" \n",
" m = ind.ga == -1\n",
" \n",
" # plot the assigned g-vectors omega vs dty (sinograms)\n",
" \n",
" axs_flat[1].scatter(colfile.omega[~m],\n",
" colfile.dty[~m],\n",
" c=ind.ga[~m],\n",
" s=2,\n",
" cmap='tab20')\n",
" \n",
" axs_flat[1].set(title=f'Sinograms of {ind.ga.max()+1} grains',\n",
" xlabel='Omega/deg',\n",
" ylabel='dty/um')\n",
" \n",
" # Define weak peaks as all non-assigned peaks with intensity 1e-4 of max\n",
" cut = colfile.sum_intensity[m].max() * 1e-4\n",
" weak = colfile.sum_intensity[m] < cut\n",
" \n",
" # Plot unassigned peaks in omega vs dty\n",
" \n",
" axs_flat[2].scatter(colfile.omega[m][weak], colfile.dty[m][weak], s=2, label='weak')\n",
" axs_flat[2].scatter(colfile.omega[m][~weak], colfile.dty[m][~weak], s=2, label='not weak')\n",
" \n",
" axs_flat[2].set(title='Sinograms of unassigned peaks',\n",
" xlabel='Omega/deg',\n",
" ylabel='dty/um')\n",
" axs_flat[2].legend()\n",
" \n",
" # Plot d-star vs intensity for all assigned peaks\n",
" \n",
" axs_flat[3].scatter(colfile.ds[~m], colfile.sum_intensity[~m], s=2)\n",
" axs_flat[3].set(title='Intensity of all assigned peaks',\n",
" xlabel='d-star',\n",
" ylabel='Intensity',\n",
" yscale='log')\n",
" \n",
" # Plot d-star vs intensity for all unassigned peaks\n",
" \n",
" axs_flat[4].scatter(colfile.ds[m][weak], colfile.sum_intensity[m][weak], s=2, label='weak')\n",
" axs_flat[4].scatter(colfile.ds[m][~weak], colfile.sum_intensity[m][~weak], s=2, label='not weak')\n",
" \n",
" axs_flat[4].set(title='Intensity of all unassigned peaks',\n",
" xlabel='d-star',\n",
" ylabel='Intensity',\n",
" yscale='log')\n",
" axs_flat[4].legend()\n",
" \n",
" # Get the number of peaks per grain\n",
" \n",
" npks = [(ind.ga == i).sum() for i in range(len(ind.ubis))]\n",
" \n",
" # Plot histogram of number of peaks per grain\n",
" \n",
" axs_flat[5].hist(npks, bins=64)\n",
" axs_flat[5].set(title='Hist of peaks per grain',\n",
" xlabel='Number of peaks',\n",
" ylabel='Number of grains')\n",
" \n",
" for ax in axs_flat:\n",
" ax.set_box_aspect(0.7)\n",
" \n",
" plt.show()"
"utils.plot_index_results(indexer, cf_strong, 'First attempt')"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"tags": []
},
"metadata": {},
"outputs": [],
"source": [
"plot_index_results(indexer, cf_strong, 'First attempt')"
"utils.plot_grain_sinograms(grains, cf_strong, 25)"
]
},
{
Expand Down
Loading
Loading