diff --git a/ImageD11/frelon_peaksearch.py b/ImageD11/frelon_peaksearch.py index 5cae1031..a62c78de 100755 --- a/ImageD11/frelon_peaksearch.py +++ b/ImageD11/frelon_peaksearch.py @@ -372,6 +372,8 @@ def process(ds, worker_args, ncpu=None): """ if ncpu is None: nthreads = max(1, ImageD11.cImageD11.cores_available() - 1) + else: + nthreads = int(ncpu) hname = ds.masterfile scan_name = ds.scans[0] frames_dset = scan_name + "/measurement/" + ds.detector diff --git a/ImageD11/nbGui/S3DXRD/pbp_1_indexing.ipynb b/ImageD11/nbGui/S3DXRD/pbp_1_indexing.ipynb index e58fc614..ac38b4d4 100755 --- a/ImageD11/nbGui/S3DXRD/pbp_1_indexing.ipynb +++ b/ImageD11/nbGui/S3DXRD/pbp_1_indexing.ipynb @@ -32,6 +32,12 @@ }, "outputs": [], "source": [ + "import os\n", + "\n", + "os.environ['OMP_NUM_THREADS'] = '1'\n", + "os.environ['OPENBLAS_NUM_THREADS'] = '1'\n", + "os.environ['MKL_NUM_THREADS'] = '1'\n", + "\n", "exec(open('/data/id11/nanoscope/install_ImageD11_from_git.py').read())\n", "PYTHONPATH = setup_ImageD11_from_git( ) # ( os.path.join( os.environ['HOME'],'Code'), 'ImageD11_git' )" ] diff --git a/ImageD11/nbGui/S3DXRD/run_astra_recon.py b/ImageD11/nbGui/S3DXRD/run_astra_recon.py index 9931a661..8fc54b18 100755 --- a/ImageD11/nbGui/S3DXRD/run_astra_recon.py +++ b/ImageD11/nbGui/S3DXRD/run_astra_recon.py @@ -1,16 +1,23 @@ def main(h5name, dsfile, group_name): + print('Loading dataset') # load the dataset ds = ImageD11.sinograms.dataset.load(dsfile) - + + print('Loading grain sinograms from disk') grainsinos = read_h5(h5name, ds, group_name=group_name) - - for gs in grainsinos: + + print('Reconstructing grain sinograms') + for inc, gs in enumerate(grainsinos): gs.recon(method="astra", astra_method="EM_CUDA") + print('Reconstructed ' + inc+1 + '/' + len(grainsinos)) # mask recon after running + print('Masking reconstructions') for gs in grainsinos: gs.recons['astra'] = np.where(gs.recon_mask, gs.recons['astra'], 0) - + + print('Reconstructions finished') + print('Writing grains to disk') write_h5(h5name, grainsinos, overwrite_grains=True, group_name=group_name) diff --git a/ImageD11/nbGui/S3DXRD/tomo_1_index.ipynb b/ImageD11/nbGui/S3DXRD/tomo_1_index.ipynb index 6c871b8a..4cbc0dd1 100755 --- a/ImageD11/nbGui/S3DXRD/tomo_1_index.ipynb +++ b/ImageD11/nbGui/S3DXRD/tomo_1_index.ipynb @@ -186,7 +186,7 @@ "\n", "ucell.makerings(cf_4d.ds.max())\n", "\n", - "fig, ax = plt.subplots(figsize=(16,9), layout='constrained')\n", + "fig, ax = plt.subplots(figsize=(10,5), layout='constrained')\n", "\n", "ax.scatter(cf_4d.ds, cf_4d.eta, s=1)\n", "ax.plot( ucell.ringds, [0,]*len(ucell.ringds), '|', ms=90, c=\"red\")\n", @@ -258,7 +258,7 @@ "source": [ "# now we can take a look at the intensities of the remaining peaks\n", "\n", - "fig, ax = plt.subplots(figsize=(16, 9), constrained_layout=True)\n", + "fig, ax = plt.subplots(figsize=(10, 5), constrained_layout=True)\n", "\n", "ax.plot(cf_4d.ds, cf_4d.sum_intensity,',', label='cf_4d')\n", "ax.plot(cf_strong.ds, cf_strong.sum_intensity,',', label='cf_strong')\n", @@ -324,7 +324,7 @@ "source": [ "# let's plot the assigned peaks\n", "\n", - "fig, ax = plt.subplots(layout='constrained', figsize=(16,9))\n", + "fig, ax = plt.subplots(layout='constrained', figsize=(10,5))\n", "\n", "# indexer.ra is the ring assignments\n", "\n", @@ -360,7 +360,7 @@ "# the sequence of minpks fractions the indexer will iterate through\n", "fracs = [0.9, 0.7]\n", "# the tolerance in g-vector angle\n", - "cosine_tol = np.cos(np.radians(90 - 0.25))\n", + "cosine_tol = np.cos(np.radians(90 - ds.ostep))\n", "# the max number of UBIs we can find per pair of rings\n", "max_grains = 1000\n", "\n", @@ -447,7 +447,7 @@ }, "outputs": [], "source": [ - "utils.plot_grain_sinograms(grains, cf_strong)" + "utils.plot_grain_sinograms(grains, cf_strong, min(len(grains), 25))" ] }, { diff --git a/ImageD11/nbGui/S3DXRD/tomo_2_map.ipynb b/ImageD11/nbGui/S3DXRD/tomo_2_map.ipynb index 2d3a202d..fa442451 100644 --- a/ImageD11/nbGui/S3DXRD/tomo_2_map.ipynb +++ b/ImageD11/nbGui/S3DXRD/tomo_2_map.ipynb @@ -26,6 +26,12 @@ }, "outputs": [], "source": [ + "import os\n", + "\n", + "os.environ['OMP_NUM_THREADS'] = '1'\n", + "os.environ['OPENBLAS_NUM_THREADS'] = '1'\n", + "os.environ['MKL_NUM_THREADS'] = '1'\n", + "\n", "exec(open('/data/id11/nanoscope/install_ImageD11_from_git.py').read())\n", "PYTHONPATH = setup_ImageD11_from_git( ) # ( os.path.join( os.environ['HOME'],'Code'), 'ImageD11_git' )" ] @@ -110,7 +116,9 @@ "outputs": [], "source": [ "# pick a phase\n", - "phase_str = 'Fe'" + "phase_str = 'Fe'\n", + "\n", + "ucell = ds.phases.unitcells[phase_str]" ] }, { @@ -178,6 +186,29 @@ "cf_strong.nrows" ] }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# now we can take a look at the intensities of the remaining peaks\n", + "\n", + "fig, ax = plt.subplots(figsize=(10, 5), constrained_layout=True)\n", + "\n", + "ax.plot(cf_4d.ds, cf_4d.sum_intensity,',', label='cf_4d')\n", + "ax.plot(cf_strong.ds, cf_strong.sum_intensity,',', label='cf_strong')\n", + "ax.plot(ucell.ringds, [1e4,]*len(ucell.ringds), '|', ms=90, c=\"red\")\n", + "\n", + "ax.semilogy()\n", + "\n", + "ax.set_xlabel(\"Dstar\")\n", + "ax.set_ylabel(\"Intensity\")\n", + "ax.legend()\n", + "\n", + "plt.show()" + ] + }, { "cell_type": "code", "execution_count": null, diff --git a/ImageD11/nbGui/S3DXRD/tomo_2_map_minor_phase.ipynb b/ImageD11/nbGui/S3DXRD/tomo_2_map_minor_phase.ipynb index f1c9180a..bc0bdde0 100755 --- a/ImageD11/nbGui/S3DXRD/tomo_2_map_minor_phase.ipynb +++ b/ImageD11/nbGui/S3DXRD/tomo_2_map_minor_phase.ipynb @@ -27,6 +27,12 @@ }, "outputs": [], "source": [ + "import os\n", + "\n", + "os.environ['OMP_NUM_THREADS'] = '1'\n", + "os.environ['OPENBLAS_NUM_THREADS'] = '1'\n", + "os.environ['MKL_NUM_THREADS'] = '1'\n", + "\n", "exec(open('/data/id11/nanoscope/install_ImageD11_from_git.py').read())\n", "PYTHONPATH = setup_ImageD11_from_git( ) # ( os.path.join( os.environ['HOME'],'Code'), 'ImageD11_git' )" ] diff --git a/ImageD11/nbGui/TDXRD/2_3DXRD_index.ipynb b/ImageD11/nbGui/TDXRD/2_3DXRD_index.ipynb index 42d712c8..e9844c87 100755 --- a/ImageD11/nbGui/TDXRD/2_3DXRD_index.ipynb +++ b/ImageD11/nbGui/TDXRD/2_3DXRD_index.ipynb @@ -306,7 +306,7 @@ "# the sequence of minpks fractions the indexer will iterate through\n", "fracs = [0.9, 0.75]\n", "# the tolerance in g-vector angle\n", - "cosine_tol = np.cos(np.radians(90 - 0.25))\n", + "cosine_tol = np.cos(np.radians(90 - ds.ostep))\n", "# the max number of UBIs we can find per pair of rings\n", "max_grains = 1000\n", "\n", @@ -342,6 +342,8 @@ "metadata": {}, "outputs": [], "source": [ + "# inverse pole figure of grain orientations\n", + "\n", "utils.plot_all_ipfs(grains)" ] }, @@ -434,7 +436,12 @@ "source": [ "# re-import our refined grains from the makemap procedure\n", "\n", - "grains2 = ImageD11.grain.read_grain_file(tmp_map_path)" + "grains2 = ImageD11.grain.read_grain_file(tmp_map_path)\n", + "\n", + "for g in grains2:\n", + " g.ref_unitcell = ucell\n", + "\n", + "utils.get_rgbs_for_grains(grains2)" ] }, { @@ -458,29 +465,9 @@ }, "outputs": [], "source": [ - "centre_plot = False\n", - "\n", - "fig = plt.figure(figsize=(12, 12))\n", - "ax = fig.add_subplot(projection='3d')\n", - "xx = [grain.translation[0] for grain in grains2]\n", - "yy = [grain.translation[1] for grain in grains2]\n", - "zz = [grain.translation[2] for grain in grains2]\n", - "# col = [utils.grain_to_rgb(grain) for grain in grains2] # IPF-Z colour instead\n", - "col = [float(grain.npks) for grain in grains2]\n", - "sizes = [0.01*(float(grain.intensity_info.split(\"mean = \")[1].split(\" , \")[0].replace(\"'\", \"\"))) for grain in grains2]\n", - "if centre_plot:\n", - " scatterplot = ax.scatter(xx-np.mean(xx), yy-np.mean(yy), zz, c=col, s=sizes)\n", - "else:\n", - " scatterplot = ax.scatter(xx, yy, zz, c=col, s=sizes)\n", - "ax.set_xlim(-200,200)\n", - "ax.set_ylim(-200,200)\n", - "ax.set_zlim(-200,200)\n", - "plt.colorbar(scatterplot)\n", - "ax.set_title(\"Grains coloured by n peaks\")\n", - "ax.set_xlabel(\"x\")\n", - "ax.set_ylabel(\"y\")\n", - "ax.set_zlabel(\"z\")\n", - "plt.show()" + "# 3D scatter plot of grain positions coloured by grain volume\n", + "\n", + "utils.plot_grain_positions(grains2, colour='npks', centre_plot=False, size_scaling=0.5)" ] }, { @@ -528,29 +515,9 @@ }, "outputs": [], "source": [ - "centre_plot = False\n", - "\n", - "fig = plt.figure(figsize=(12, 12))\n", - "ax = fig.add_subplot(projection='3d')\n", - "xx = [grain.translation[0] for grain in grains_filtered]\n", - "yy = [grain.translation[1] for grain in grains_filtered]\n", - "zz = [grain.translation[2] for grain in grains_filtered]\n", - "# col = [utils.grain_to_rgb(grain) for grain in grains_filtered] # IPF-Z colour instead\n", - "col = [float(grain.npks) for grain in grains_filtered]\n", - "sizes = [0.01*(float(grain.intensity_info.split(\"mean = \")[1].split(\" , \")[0].replace(\"'\", \"\"))) for grain in grains_filtered]\n", - "if centre_plot:\n", - " scatterplot = ax.scatter(xx-np.mean(xx), yy-np.mean(yy), zz, c=col, s=sizes)\n", - "else:\n", - " scatterplot = ax.scatter(xx, yy, zz, c=col, s=sizes)\n", - "ax.set_xlim(-200,200)\n", - "ax.set_ylim(-200,200)\n", - "ax.set_zlim(-200,200)\n", - "plt.colorbar(scatterplot)\n", - "ax.set_title(\"Grains coloured by n peaks\")\n", - "ax.set_xlabel(\"x\")\n", - "ax.set_ylabel(\"y\")\n", - "ax.set_zlabel(\"z\")\n", - "plt.show()" + "# 3D scatter plot of grain positions coloured by grain volume\n", + "\n", + "utils.plot_grain_positions(grains_filtered, colour='npks', centre_plot=False, size_scaling=0.5)" ] }, { @@ -584,6 +551,31 @@ "grains_filtered = ImageD11.grain.read_grain_file(map_path)" ] }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# compute IPF colours for grains\n", + "\n", + "for g in grains_filtered:\n", + " g.ref_unitcell = ucell\n", + "\n", + "utils.get_rgbs_for_grains(grains_filtered)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# inverse pole figure of grain orientations\n", + "\n", + "utils.plot_all_ipfs(grains_filtered)" + ] + }, { "cell_type": "code", "execution_count": null, @@ -592,9 +584,20 @@ }, "outputs": [], "source": [ - "# save grain data\n", + "# 3D scatter plot of grain positions coloured by grain volume\n", "\n", - "ds.save_grains_to_disk(grains_filtered, phase_name=phase_str)" + "utils.plot_grain_positions(grains_filtered, colour='npks', centre_plot=False, size_scaling=0.5)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# 3D scatter plot of grain positions coloured by inverse pole figure orientation\n", + "\n", + "utils.plot_grain_positions(grains_filtered, colour='z', centre_plot=False, size_scaling=0.5)" ] }, { @@ -605,29 +608,9 @@ }, "outputs": [], "source": [ - "centre_plot = False\n", - "\n", - "fig = plt.figure(figsize=(12, 12))\n", - "ax = fig.add_subplot(projection='3d')\n", - "xx = [grain.translation[0] for grain in grains_filtered]\n", - "yy = [grain.translation[1] for grain in grains_filtered]\n", - "zz = [grain.translation[2] for grain in grains_filtered]\n", - "# col = [utils.grain_to_rgb(grain) for grain in grains_filtered] # IPF-Z colour instead\n", - "col = [float(grain.npks) for grain in grains_filtered]\n", - "sizes = [0.01*(float(grain.intensity_info.split(\"mean = \")[1].split(\" , \")[0].replace(\"'\", \"\"))) for grain in grains_filtered]\n", - "if centre_plot:\n", - " scatterplot = ax.scatter(xx-np.mean(xx), yy-np.mean(yy), zz, c=col, s=sizes)\n", - "else:\n", - " scatterplot = ax.scatter(xx, yy, zz, c=col, s=sizes)\n", - "ax.set_xlim(-200,200)\n", - "ax.set_ylim(-200,200)\n", - "ax.set_zlim(-200,200)\n", - "plt.colorbar(scatterplot)\n", - "ax.set_title(\"Grains coloured by n peaks\")\n", - "ax.set_xlabel(\"x\")\n", - "ax.set_ylabel(\"y\")\n", - "ax.set_zlabel(\"z\")\n", - "plt.show()" + "# save grain data\n", + "\n", + "ds.save_grains_to_disk(grains_filtered, phase_name=phase_str)" ] }, { diff --git a/ImageD11/nbGui/TDXRD/2_3DXRD_index_friedel.ipynb b/ImageD11/nbGui/TDXRD/2_3DXRD_index_friedel.ipynb index 852f81eb..bebb5426 100755 --- a/ImageD11/nbGui/TDXRD/2_3DXRD_index_friedel.ipynb +++ b/ImageD11/nbGui/TDXRD/2_3DXRD_index_friedel.ipynb @@ -463,7 +463,7 @@ "# the sequence of minpks fractions the indexer will iterate through\n", "fracs = [0.5]\n", "# the tolerance in g-vector angle\n", - "cosine_tol = np.cos(np.radians(90 - 0.25))\n", + "cosine_tol = np.cos(np.radians(90 - ds.ostep))\n", "# the max number of UBIs we can find per pair of rings\n", "max_grains = 1000\n", "\n", @@ -522,9 +522,9 @@ "gridpars = {\n", " 'DSTOL' : 0.004,\n", " 'OMEGAFLOAT' : omega_slop,\n", - " 'COSTOL' : 0.002,\n", + " 'COSTOL' : cosine_tol,\n", " 'NPKS' : 10,\n", - " 'TOLSEQ' : [ 0.05, ],\n", + " 'TOLSEQ' : [hkl_tols_seq[-1],],\n", " 'SYMMETRY' : \"cubic\",\n", " 'RING1' : [1,5],\n", " 'RING2' : [1,5],\n", @@ -567,27 +567,9 @@ }, "outputs": [], "source": [ - "centre_plot = False\n", - "\n", - "fig = plt.figure(figsize=(12, 12))\n", - "ax = fig.add_subplot(projection='3d', proj_type=\"ortho\")\n", - "xx = [grain.translation[0] for grain in fittedgrains]\n", - "yy = [grain.translation[1] for grain in fittedgrains]\n", - "zz = [grain.translation[2] for grain in fittedgrains]\n", - "# col = [utils.grain_to_rgb(grain) for grain in grains_filtered] # IPF-Z colour instead\n", - "col = [float(grain.npks) for grain in fittedgrains]\n", - "sizes = [0.01*(float(grain.intensity_info.split(\"mean = \")[1].split(\" , \")[0].replace(\"'\", \"\"))) for grain in fittedgrains]\n", - "if centre_plot:\n", - " scatterplot = ax.scatter(xx-np.mean(xx), yy-np.mean(yy), zz, c=col, s=sizes)\n", - "else:\n", - " scatterplot = ax.scatter(xx, yy, zz, c=col, s=sizes)\n", - "ax.set_aspect(\"equal\")\n", - "plt.colorbar(scatterplot)\n", - "ax.set_title(\"Grains coloured by n peaks\")\n", - "ax.set_xlabel(\"x\")\n", - "ax.set_ylabel(\"y\")\n", - "ax.set_zlabel(\"z\")\n", - "plt.show()" + "# 3D scatter plot of grain positions coloured by grain volume\n", + "\n", + "utils.plot_grain_positions(fittedgrains, colour='npks', centre_plot=False, size_scaling=0.5)" ] }, { @@ -636,27 +618,9 @@ }, "outputs": [], "source": [ - "centre_plot = False\n", - "\n", - "fig = plt.figure(figsize=(12, 12))\n", - "ax = fig.add_subplot(projection='3d', proj_type=\"ortho\")\n", - "xx = [grain.translation[0] for grain in grains_filtered]\n", - "yy = [grain.translation[1] for grain in grains_filtered]\n", - "zz = [grain.translation[2] for grain in grains_filtered]\n", - "# col = [utils.grain_to_rgb(grain) for grain in grains_filtered] # IPF-Z colour instead\n", - "col = [float(grain.npks) for grain in grains_filtered]\n", - "sizes = [0.01*(float(grain.intensity_info.split(\"mean = \")[1].split(\" , \")[0].replace(\"'\", \"\"))) for grain in grains_filtered]\n", - "if centre_plot:\n", - " scatterplot = ax.scatter(xx-np.mean(xx), yy-np.mean(yy), zz, c=col, s=sizes)\n", - "else:\n", - " scatterplot = ax.scatter(xx, yy, zz, c=col, s=sizes)\n", - "ax.set_aspect(\"equal\")\n", - "plt.colorbar(scatterplot)\n", - "ax.set_title(\"Grains coloured by n peaks\")\n", - "ax.set_xlabel(\"x\")\n", - "ax.set_ylabel(\"y\")\n", - "ax.set_zlabel(\"z\")\n", - "plt.show()" + "# 3D scatter plot of grain positions coloured by grain volume\n", + "\n", + "utils.plot_grain_positions(grains_filtered, colour='npks', centre_plot=False, size_scaling=0.5)" ] }, { diff --git a/ImageD11/nbGui/TDXRD/2_3DXRD_index_grid.ipynb b/ImageD11/nbGui/TDXRD/2_3DXRD_index_grid.ipynb index a248c0e7..b55aad57 100755 --- a/ImageD11/nbGui/TDXRD/2_3DXRD_index_grid.ipynb +++ b/ImageD11/nbGui/TDXRD/2_3DXRD_index_grid.ipynb @@ -455,7 +455,7 @@ "gridpars = {\n", " 'DSTOL' : 0.004,\n", " 'OMEGAFLOAT' : omega_slop,\n", - " 'COSTOL' : np.cos(np.radians(90 - omega_slop)),\n", + " 'COSTOL' : np.cos(np.radians(90 - ds.ostep)),\n", " 'NPKS' : int(minpeaks),\n", " 'TOLSEQ' : makemap_tol_seq,\n", " 'SYMMETRY' : symmetry,\n", @@ -497,9 +497,26 @@ }, "outputs": [], "source": [ - "# re-import our refined grains from the makemap procedure\n", + "# import our initially indexed grains\n", + "\n", + "grains2 = ImageD11.grain.read_grain_file(map_path)\n", + "\n", + "for g in grains2:\n", + " g.ref_unitcell = ucell\n", + "\n", + "utils.get_rgbs_for_grains(grains2)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "80b7c444-1183-4a04-8c69-1f442da9db7f", + "metadata": {}, + "outputs": [], + "source": [ + "# inverse pole figure of grain orientations\n", "\n", - "grains2 = ImageD11.grain.read_grain_file(map_path)" + "utils.plot_all_ipfs(grains2)" ] }, { @@ -511,27 +528,21 @@ }, "outputs": [], "source": [ - "centre_plot = False\n", - "\n", - "fig = plt.figure(figsize=(12, 12))\n", - "ax = fig.add_subplot(projection='3d', proj_type=\"ortho\")\n", - "xx = [grain.translation[0] for grain in grains2]\n", - "yy = [grain.translation[1] for grain in grains2]\n", - "zz = [grain.translation[2] for grain in grains2]\n", - "# col = [utils.grain_to_rgb(grain) for grain in grains2] # IPF-Z colour instead\n", - "col = [float(grain.npks) for grain in grains2]\n", - "sizes = [0.01*(float(grain.intensity_info.split(\"mean = \")[1].split(\" , \")[0].replace(\"'\", \"\"))) for grain in grains2]\n", - "if centre_plot:\n", - " scatterplot = ax.scatter(xx-np.mean(xx), yy-np.mean(yy), zz, c=col, s=sizes)\n", - "else:\n", - " scatterplot = ax.scatter(xx, yy, zz, c=col, s=sizes)\n", - "plt.colorbar(scatterplot)\n", - "ax.set_title(\"Grains coloured by n peaks\")\n", - "ax.set_xlabel(\"x\")\n", - "ax.set_ylabel(\"y\")\n", - "ax.set_zlabel(\"z\")\n", - "ax.set_aspect(\"equal\")\n", - "plt.show()" + "# 3D scatter plot of grain positions coloured by grain volume\n", + "\n", + "utils.plot_grain_positions(grains2, colour='npks', centre_plot=False, size_scaling=0.5)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "a11efaf8-7246-49c5-bfc6-78c6654cf01a", + "metadata": {}, + "outputs": [], + "source": [ + "# 3D scatter plot of grain positions coloured by inverse pole figure orientation\n", + "\n", + "utils.plot_grain_positions(grains2, colour='z', centre_plot=False, size_scaling=0.5)" ] }, { @@ -590,27 +601,9 @@ }, "outputs": [], "source": [ - "centre_plot = False\n", - "\n", - "fig = plt.figure(figsize=(12, 12))\n", - "ax = fig.add_subplot(projection='3d', proj_type=\"ortho\")\n", - "xx = [grain.translation[0] for grain in grains3]\n", - "yy = [grain.translation[1] for grain in grains3]\n", - "zz = [grain.translation[2] for grain in grains3]\n", - "# col = [utils.grain_to_rgb(grain) for grain in grains2] # IPF-Z colour instead\n", - "col = [float(grain.npks) for grain in grains3]\n", - "sizes = [0.01*(float(grain.intensity_info.split(\"mean = \")[1].split(\" , \")[0].replace(\"'\", \"\"))) for grain in grains3]\n", - "if centre_plot:\n", - " scatterplot = ax.scatter(xx-np.mean(xx), yy-np.mean(yy), zz, c=col, s=sizes)\n", - "else:\n", - " scatterplot = ax.scatter(xx, yy, zz, c=col, s=sizes)\n", - "ax.set_aspect(\"equal\")\n", - "plt.colorbar(scatterplot)\n", - "ax.set_title(\"Grains coloured by n peaks\")\n", - "ax.set_xlabel(\"x\")\n", - "ax.set_ylabel(\"y\")\n", - "ax.set_zlabel(\"z\")\n", - "plt.show()" + "# 3D scatter plot of grain positions coloured by grain volume\n", + "\n", + "utils.plot_grain_positions(grains3, colour='npks', centre_plot=False, size_scaling=0.5)" ] }, { @@ -663,27 +656,9 @@ }, "outputs": [], "source": [ - "centre_plot = False\n", - "\n", - "fig = plt.figure(figsize=(12, 12))\n", - "ax = fig.add_subplot(projection='3d', proj_type=\"ortho\")\n", - "xx = [grain.translation[0] for grain in grains_filtered]\n", - "yy = [grain.translation[1] for grain in grains_filtered]\n", - "zz = [grain.translation[2] for grain in grains_filtered]\n", - "# col = [utils.grain_to_rgb(grain) for grain in grains_filtered] # IPF-Z colour instead\n", - "col = [float(grain.npks) for grain in grains_filtered]\n", - "sizes = [0.01*(float(grain.intensity_info.split(\"mean = \")[1].split(\" , \")[0].replace(\"'\", \"\"))) for grain in grains_filtered]\n", - "if centre_plot:\n", - " scatterplot = ax.scatter(xx-np.mean(xx), yy-np.mean(yy), zz, c=col, s=sizes)\n", - "else:\n", - " scatterplot = ax.scatter(xx, yy, zz, c=col, s=sizes)\n", - "ax.set_aspect(\"equal\")\n", - "plt.colorbar(scatterplot)\n", - "ax.set_title(\"Grains coloured by n peaks\")\n", - "ax.set_xlabel(\"x\")\n", - "ax.set_ylabel(\"y\")\n", - "ax.set_zlabel(\"z\")\n", - "plt.show()" + "# 3D scatter plot of grain positions coloured by grain volume\n", + "\n", + "utils.plot_grain_positions(grains_filtered, colour='npks', centre_plot=False, size_scaling=0.5)" ] }, { @@ -743,7 +718,12 @@ }, "outputs": [], "source": [ - "grains_final = ImageD11.grain.read_grain_file(new_filtered_map_path)" + "grains_final = ImageD11.grain.read_grain_file(new_filtered_map_path)\n", + "\n", + "for g in grains_final:\n", + " g.ref_unitcell = ucell\n", + "\n", + "utils.get_rgbs_for_grains(grains_final)" ] }, { @@ -755,27 +735,21 @@ }, "outputs": [], "source": [ - "centre_plot = False\n", - "\n", - "fig = plt.figure(figsize=(12, 12))\n", - "ax = fig.add_subplot(projection='3d', proj_type=\"ortho\")\n", - "xx = [grain.translation[0] for grain in grains_final]\n", - "yy = [grain.translation[1] for grain in grains_final]\n", - "zz = [grain.translation[2] for grain in grains_final]\n", - "# col = [utils.grain_to_rgb(grain) for grain in grains_filtered] # IPF-Z colour instead\n", - "col = [float(grain.npks) for grain in grains_final]\n", - "sizes = [0.01*(float(grain.intensity_info.split(\"mean = \")[1].split(\" , \")[0].replace(\"'\", \"\"))) for grain in grains_final]\n", - "if centre_plot:\n", - " scatterplot = ax.scatter(xx-np.mean(xx), yy-np.mean(yy), zz, c=col, s=sizes)\n", - "else:\n", - " scatterplot = ax.scatter(xx, yy, zz, c=col, s=sizes)\n", - "ax.set_aspect(\"equal\")\n", - "plt.colorbar(scatterplot)\n", - "ax.set_title(\"Grains coloured by n peaks\")\n", - "ax.set_xlabel(\"x\")\n", - "ax.set_ylabel(\"y\")\n", - "ax.set_zlabel(\"z\")\n", - "plt.show()" + "# 3D scatter plot of grain positions coloured by grain volume\n", + "\n", + "utils.plot_grain_positions(grains_final, colour='npks', centre_plot=False, size_scaling=0.5)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "b9178835-a3ce-4829-846f-b99a71a1300f", + "metadata": {}, + "outputs": [], + "source": [ + "# 3D scatter plot of grain positions coloured by grain volume\n", + "\n", + "utils.plot_grain_positions(grains_final, colour='z', centre_plot=False, size_scaling=0.5)" ] }, { @@ -845,7 +819,7 @@ "metadata": {}, "outputs": [], "source": [ - "ds.save_grains_to_disk(grains_filtered, phase_name=phase_str)" + "ds.save_grains_to_disk(grains_final, phase_name=phase_str)" ] }, { @@ -1009,7 +983,7 @@ " new_cf_3d_path = cf_3d_path + '.new'\n", " makemap_output = !makemap.py -p {oldparfile} -u {filtered_map_path} -U {new_filtered_map_path} -f {cf_3d_path} -s {symmetry} -t {makemap_tol_seq[-1]} --omega_slop={omega_slop} --no_sort\n", " grains_final = ImageD11.grain.read_grain_file(new_filtered_map_path)\n", - " ds.save_grains_to_disk(grains_filtered, phase_name=phase_str)\n", + " ds.save_grains_to_disk(grains_final, phase_name=phase_str)\n", "\n", " cf_3d = ImageD11.columnfile.columnfile(new_cf_3d_path)\n", " ImageD11.columnfile.colfile_to_hdf(cf_3d, ds.col3dfile, name='peaks')\n", diff --git a/ImageD11/nbGui/nb_utils.py b/ImageD11/nbGui/nb_utils.py index 11c883c2..24de1bb3 100644 --- a/ImageD11/nbGui/nb_utils.py +++ b/ImageD11/nbGui/nb_utils.py @@ -149,6 +149,7 @@ def prepare_mlem_bash(ds, grains, id11_code_path, n_simultaneous_jobs=50, cores_ # define memory needs and number of tasks for each array job #SBATCH --ntasks=1 #SBATCH --cpus-per-task={cores_per_task} +#SBATCH --mem-per-cpu=20G # date echo PYTHONPATH={id11_code_path} python3 {python_script_path} {grainsfile} $SLURM_ARRAY_TASK_ID {dsfile} {reconfile} {cores_per_task} > {log_path} 2>&1 @@ -314,23 +315,23 @@ def find_datasets_to_process(rawdata_path, skips_dict, dset_prefix, sample_list) return samples_dict - def save_ubi_map(ds, ubi_map, eps_map, misorientation_map, ipf_x_col_map, ipf_y_col_map, ipf_z_col_map): - with h5py.File(ds.pbpubifile, 'w') as hout: - grp = hout.create_group('arrays') - save_array(grp, 'ubi_map', ubi_map).attrs['description'] = 'Refined UBI values at each pixel' - save_array(grp, 'eps_map', eps_map).attrs['description'] = 'Strain matrices (sample ref) at each pixel' - save_array(grp, 'misorientation_map', misorientation_map).attrs[ - 'description'] = 'Misorientation to grain avg at each pixel' - ipfxdset = save_array(grp, 'ipf_x_col_map', ipf_x_col_map) - ipfxdset.attrs['description'] = 'IPF X color at each pixel' - ipfxdset.attrs['CLASS'] = 'IMAGE' - ipfydset = save_array(grp, 'ipf_y_col_map', ipf_y_col_map) - ipfydset.attrs['description'] = 'IPF Y color at each pixel' - ipfydset.attrs['CLASS'] = 'IMAGE' - ipfzdset = save_array(grp, 'ipf_z_col_map', ipf_z_col_map) - ipfzdset.attrs['description'] = 'IPF Z color at each pixel' - ipfzdset.attrs['CLASS'] = 'IMAGE' + raise ValueError('This function is deprecated" Use ImageD11.sinograms.tensor_map.TensorMap instead') + # with h5py.File(ds.pbpubifile, 'w') as hout: + # grp = hout.create_group('arrays') + # save_array(grp, 'ubi_map', ubi_map).attrs['description'] = 'Refined UBI values at each pixel' + # save_array(grp, 'eps_map', eps_map).attrs['description'] = 'Strain matrices (sample ref) at each pixel' + # save_array(grp, 'misorientation_map', misorientation_map).attrs[ + # 'description'] = 'Misorientation to grain avg at each pixel' + # ipfxdset = save_array(grp, 'ipf_x_col_map', ipf_x_col_map) + # ipfxdset.attrs['description'] = 'IPF X color at each pixel' + # ipfxdset.attrs['CLASS'] = 'IMAGE' + # ipfydset = save_array(grp, 'ipf_y_col_map', ipf_y_col_map) + # ipfydset.attrs['description'] = 'IPF Y color at each pixel' + # ipfydset.attrs['CLASS'] = 'IMAGE' + # ipfzdset = save_array(grp, 'ipf_z_col_map', ipf_z_col_map) + # ipfzdset.attrs['description'] = 'IPF Z color at each pixel' + # ipfzdset.attrs['CLASS'] = 'IMAGE' ### Sinogram stuff @@ -571,7 +572,50 @@ def plot_all_ipfs(grains): - +def plot_grain_positions(grains, colour='npks', centre_plot=False, size_scaling=0.5): + """ + colour: choose from 'npks' or one of 'x', 'y', 'z' for IPF scaling + centre_plot: choose whether to centre the plot horizontally (x and y) + size_scaling: we only know relative grain sizes, adjust this to scale the diameter of the points on the plot + """ + if colour.lower() not in ['npks', 'x', 'y', 'z']: + raise ValueError("colour should be one of ['npks', 'x', 'y', 'z']") + fig = plt.figure(figsize=(12, 12)) + ax = fig.add_subplot(projection='3d', proj_type="ortho") + xx = [grain.translation[0] for grain in grains] + yy = [grain.translation[1] for grain in grains] + zz = [grain.translation[2] for grain in grains] + if colour == 'npks': + ax.set_title("Grains coloured by number of peaks indexed") + col = [float(grain.npks) for grain in grains] + elif colour.lower() in ['x', 'y', 'z']: + rgbattr = 'rgb_' + colour.lower() + try: + col = [getattr(grain, rgbattr) for grain in grains] # IPF colour + except AttributeError: + # couldn't get the IPF attributes + # try to compute it first + # will still fail if we don't have reference unitcells + get_rgbs_for_grains(grains) + col = [getattr(grain, rgbattr) for grain in grains] # IPF colour + ax.set_title("Grains coloured by IPF", colour.lower()) + # sizes in MPL 3D scale the area of the plot + # intensity info is proportional to volume + # decrease to radius then scale to area with power(x, 2/3) + sizes = [size_scaling*np.power((float(grain.intensity_info.split("mean = ")[1].split(" , ")[0].replace("'", ""))), 2/3) for grain in grains] + if centre_plot: + scatterplot = ax.scatter(xx-np.mean(xx), yy-np.mean(yy), zz, c=col, s=sizes) + else: + scatterplot = ax.scatter(xx, yy, zz, c=col, s=sizes) + if colour == 'npks': + plt.colorbar(scatterplot) + + ax.set_aspect("equal") + ax.set_xlabel("x") + ax.set_ylabel("y") + ax.set_zlabel("z") + plt.show() + # backwards compatible do_index = ImageD11.indexing.do_index diff --git a/ImageD11/sinograms/tensor_map.py b/ImageD11/sinograms/tensor_map.py index 37769e0e..b37533f2 100644 --- a/ImageD11/sinograms/tensor_map.py +++ b/ImageD11/sinograms/tensor_map.py @@ -453,12 +453,6 @@ def strain_crystal_to_stress_crystal(strain_crystal, stiffness_tensor, B0, phase res[...] = stress_crystal -@numba.guvectorize([(numba.float64[:, :], numba.float64[:])], '(n,n)->()', nopython=True) -def sig_to_hydro(sig, res): - """Get hydrostatic stress scalar from stress tensor (frame invariant)""" - res[...] = np.sum(np.diag(sig)) / 3 - - @numba.guvectorize([(numba.float64[:, :], numba.float64[:])], '(n,n)->()', nopython=True) def sig_to_vm(sig, res): """Get von-Mises stress scalar from stress tensor""" @@ -820,6 +814,30 @@ def eps_crystal(self): eps_crystal_map = ubi_and_unitcell_to_eps_crystal(self.UBI, self.dzero_unitcell) self.add_map('eps_crystal', eps_crystal_map) return eps_crystal_map + + @property + def eps_hydro(self): + """ + The per-voxel hydrostatic strain tensor (frame invariant) + """ + if 'eps_hydro' in self.keys(): + return self.maps['eps_hydro'] + else: + eps_hydro_map = (self.eps_sample.trace(axis1=-2, axis2=-1)/3)[..., np.newaxis, np.newaxis] * np.eye(3) + self.add_map('eps_hydro', eps_hydro_map) + return eps_hydro_map + + @property + def eps_devia(self): + """ + The per-voxel deviatoric strain tensor + """ + if 'eps_devia' in self.keys(): + return self.maps['eps_devia'] + else: + eps_devia_map = self.eps_sample - self.eps_hydro + self.add_map('eps_devia', eps_devia_map) + return eps_devia_map # TODO - make multiphase - store C map for each voxel def get_stress(self, stiffness_tensor, phase_id): @@ -847,18 +865,50 @@ def get_stress(self, stiffness_tensor, phase_id): self.add_map('sig_crystal', stress_crystal_map) self.add_map('sig_sample', stress_sample_map) + + @property + def sig_sample(self): + """ + The per-voxel Biot stress tensor in the sample reference frame + """ + if 'sig_sample' in self.keys(): + return self.maps['sig_sample'] + else: + raise AttributeError('Stress not computed! Run self.get_stress() first') + + @property + def sig_crystal(self): + """ + The per-voxel Biot stress tensor in the crystal reference frame + """ + if 'sig_crystal' in self.keys(): + return self.maps['sig_crystal'] + else: + raise AttributeError('Stress not computed! Run self.get_stress() first') @property def sig_hydro(self): """ - The per-voxel hydrostatic stress scalar (frame invariant) + The per-voxel hydrostatic stress tensor (frame invariant) """ if 'sig_hydro' in self.keys(): return self.maps['sig_hydro'] else: - sig_hydro_map = sig_to_hydro(self.sig_crystal) + sig_hydro_map = (self.sig_sample.trace(axis1=-2, axis2=-1)/3)[..., np.newaxis, np.newaxis] * np.eye(3) self.add_map('sig_hydro', sig_hydro_map) return sig_hydro_map + + @property + def sig_devia(self): + """ + The per-voxel deviatoric stress tensor + """ + if 'sig_devia' in self.keys(): + return self.maps['sig_devia'] + else: + sig_devia_map = self.sig_sample - self.sig_hydro + self.add_map('sig_devia', sig_devia_map) + return sig_devia_map @property def sig_mises(self): @@ -1329,6 +1379,9 @@ def from_grainsinos(cls, grainsinos, method="iradon", use_gids=True, cutoff_leve phases[phase_inc] = phase except NameError: raise AttributeError("Some/all grains are missing reference unit cells! Can't continue") + + if not all([method in gs.recons for gs in grainsinos]): + raise AttributeError("Not all grainsinos have the method you specified! Check if your reconstructions worked") # work out the phase ID for each grain phase_ids = [phases_list.index(gs.grain.ref_unitcell) for gs in grainsinos]