Required python packages that might not be installed:
In [16]:
import xray_vision.mpl_plotting.speckle as speckle_plot
import xray_vision.mpl_plotting as mpl_plot
# imports... all the imports
import numpy as np
import pandas as pd
import matplotlib as mpl
import matplotlib.ticker as mticks
import matplotlib.pyplot as plt
from mpl_toolkits.axes_grid1.axes_grid import ImageGrid
import skxray.core.roi as roi
import skxray.core.correlation as corr
import skxray.core.utils as core
import scipy.ndimage.measurements as ndim
In [17]:
interactive_mode = False
if interactive_mode:
%matplotlib notebook
else:
%matplotlib inline
backend = mpl.get_backend()
In [18]:
%run download.py
In [19]:
# we will provide a link to the data once it is approved
# for distribution by Dr. Fluerasu
data_dir = "Duke_data/"
duke_ndata = np.load(data_dir+"duke_img_1_5000.npy")
duke_dark = np.load(data_dir+"duke_dark.npy")
In [20]:
pixel_size = (0.1, 0.1) # in mm
image_center = (133, 143) # in pixels
# if we are 1.5 or later
if float('.'.join(mpl.__version__.split('.')[:2])) >= 1.5:
cmap = 'viridis'
else:
cmap = 'CMRmap'
In [21]:
# fake that this is 5 different points in the same sample
splits = [0, 500, 1500, 2000, 3500, duke_ndata.shape[0]]
split_data = []
ranges = []
for lower, upper in zip(splits, splits[1:]):
split_data.append(duke_ndata[lower:upper])
ranges.append('%s:%s' % (lower, upper))
data_labels = ['sample_point_%s' % (i+1) for i, region in
zip(range(len(split_data)), ranges)]
data_dict = {label: data for label, data in
zip(data_labels, split_data)}
In [22]:
# subtract the dark image from the dataset
data_dict = {label: data - duke_dark for label, data in
data_dict.items()}
In [23]:
# load the mask(s) and mask the data
mask1 = np.load(data_dir+"new_mask4.npy")
mask2 = np.load(data_dir+"Luxi_duke_mask.npy")
mask = mask1 + mask2
masked_data = {label: data * mask for label, data in data_dict.items()}
use the ImageGrid to make this plot look very nice because it handles a lot of annoying setup bits, such as:
In [24]:
fig = plt.figure(figsize=(10, 4))
from mpl_toolkits.axes_grid1.axes_grid import ImageGrid
grid = ImageGrid(fig, '111', nrows_ncols=(1, 3), share_all=True,
axes_pad=0.2, label_mode="L", cbar_location="bottom",
cbar_mode="single", cbar_size='3%')
# plot the images
im1 = grid[0].imshow(mask1, cmap=cmap)
im2 = grid[1].imshow(mask2, cmap=cmap)
im3 = grid[2].imshow(mask, cmap=cmap)
# instantiate the colorbar
grid.cbar_axes[0].colorbar(im3, ticks=mticks.MaxNLocator(integer=True))
# set the titles of the axes
grid[0].set_title('new_mask4.npy')
grid[1].set_title('Luxi_duke_mask.npy')
grid[2].set_title('combined mask')
if 'inline' in backend:
plt.show()
In [25]:
# load the mask(s) and mask the data
mask1 = ~mask1
mask2 = ~mask2
mask = ~mask
masked_data = {label: data * mask for label, data in data_dict.items()}
In [26]:
fig = plt.figure(figsize=(10, 4))
from mpl_toolkits.axes_grid1.axes_grid import ImageGrid
grid = ImageGrid(fig, '111', nrows_ncols=(1, 3), share_all=True,
axes_pad=0.2, label_mode="L", cbar_location="bottom",
cbar_mode="single", cbar_size='3%')
# plot the images
im1 = grid[0].imshow(mask1, cmap=cmap)
im2 = grid[1].imshow(mask2, cmap=cmap)
im3 = grid[2].imshow(mask, cmap=cmap)
# instantiate the colorbar
grid.cbar_axes[0].colorbar(im3, ticks=mticks.MaxNLocator(integer=True))
# set the titles of the axes
grid[0].set_title('new_mask4.npy')
grid[1].set_title('Luxi_duke_mask.npy')
grid[2].set_title('combined mask')
if 'inline' in backend:
plt.show()
In [27]:
# compute the average raw and masked images
avg_raw = np.average(np.vstack(data_dict.values()), axis=0)
avg_masked = np.average(np.vstack(masked_data.values()), axis=0)
In [28]:
fig = plt.figure(figsize=(15, 7))
grid = ImageGrid(fig, (1,1,1), nrows_ncols=(1, 2), share_all=True,
cbar_mode='each', axes_pad=.3)
im1 = grid[0].imshow(np.log(avg_raw), cmap=cmap)
im2 = grid[1].imshow(np.log(avg_masked), cmap=cmap)
grid[0].set_title("Duke Silica -- log of average")
grid[1].set_title("Duke Silica -- log of masked average")
grid.cbar_axes[0].colorbar(im1)
grid.cbar_axes[1].colorbar(im2)
# the following code is needed to make it easy to switch
# between '%matplotlib inline' and '%matplotlib notebook'
if 'inline' in backend:
plt.show()
In [29]:
# define the parameters for the ring regions of interest
starting_radius = 26
ring_width = 8
space_between_rings = 1
num_rings = 4
# image center was defined when we loaded the data
# get the inner/outer radii of the rings
ring_edges = roi.ring_edges(starting_radius, ring_width, space_between_rings, num_rings)
print('ring_edges are defined as [inner_radius, outer_radius]:\n%s' % ring_edges)
ring_centers = np.average(ring_edges, axis=1)
print('\nring_centers: %s' % ring_centers)
# rings is a labeled array where 0 is background and 1-4 are ROIS 1-4.
rings = roi.rings(edges=ring_edges, center=image_center, shape=avg_raw.shape)
# multiplying the rings with the image mask removes the bad information (beam stop, beam stop scatter,
# bad pixels, etc..) from the rings ROI. See the plot below
rings_mask = rings*mask
In [30]:
fig = plt.figure(figsize=(10, 5))
grid = ImageGrid(fig, '111', nrows_ncols=(1, 2), share_all=True, axes_pad=0.2,
label_mode="L", cbar_location="bottom", cbar_mode="single", cbar_size='2%')
im1 = mpl_plot.show_label_array(grid[0], rings, cmap=cmap)
im2 = mpl_plot.show_label_array(grid[1], rings_mask, cmap=cmap)
# turn on the colorbar
grid.cbar_axes[0].colorbar(im1, ticks=mticks.MaxNLocator(integer=True))
grid[0].set_title('Rings')
grid[1].set_title('Masked Rings')
if 'inline' in backend:
plt.show()
In [31]:
masked_data['sample_point_1'].shape
Out[31]:
In [32]:
intensity, index_list = roi.mean_intensity(masked_data['sample_point_1'], rings_mask)
In [33]:
intensity.shape
Out[33]:
In [34]:
roi_data = {}
for k, v in sorted(masked_data.items()):
intensity, index_list = roi.mean_intensity(v, rings_mask)
roi_data[k] = [arr for arr in intensity.T]
# roi_data[k] = intensity
# one line dict comprehension
# roi_data = {k: roi.mean_intensity(v, rings_mask) for k, v in masked_data.items()}
row_labels = ['roi_%s' % i for i in index_list]
df = pd.DataFrame(roi_data, row_labels)
In [35]:
df
Out[35]:
In [36]:
# fig, ax = plt.subplots(len(df.index))
fig, ax = plt.subplots(len(df.index), figsize=(15, 12), sharex=True)
arts = speckle_plot.mean_intensity_plotter(ax, df, cmap=cmap)
if 'inline' in backend:
plt.show()
In [37]:
arts
Out[37]:
If you guessed "a dataframe", then you would be correct! It is returning the matplotlib artists in the same layout as the data that was passed in.
What is a matplotlib artist? It is the thing that you can modify to change the appearance of the plot after it has been created. There is a "gotcha" here. If you are using matplotlib in "inline" mode (%matplotlib inline) then you will not be able to interactively modify a figure. Go back to the second cell code cell and make sure that "interactive" is set to True, if you want interactive plots. Then come back and play around with the next few cells. Otherwise, skip them.
In [38]:
arts.sample_point_1.roi_1.set_color('pink')
Let's pretend that there is an issue with the second sample (point_2)
Obviously there is not a problem, but if we want to use the helper function that removes a dataset
In [39]:
# make a copy of the dataframe so that we can mess it up for tutorial purposes
bad_df = df.copy()
In [40]:
bad_df['sample_point_2'] = bad_df['sample_point_2'] * 1.1
In [41]:
# fig, ax = plt.subplots(len(df.index))
fig, ax = plt.subplots(len(bad_df.index), figsize=(15, 12), sharex=True)
arts = speckle_plot.mean_intensity_plotter(ax, bad_df, cmap=cmap)
if 'inline' in backend:
plt.show()
In [42]:
dropped_column = df.drop('sample_point_2', axis=1)
# dropped_column = df.drop('point_2', axis=1, inplace=True)
In [43]:
dropped_column
Out[43]:
Note that there is no point_2 column any more. However, the original dataframe remains unchanged
In [44]:
df
Out[44]:
In [45]:
dropped_row = df.drop(df.index[[2]])
# dropped_row = df.drop(df.index[[2]], inplace=False)
In [46]:
dropped_row
Out[46]:
Note that there is no roi_3 anymore. However, the original dataframe remains unchanged
In [47]:
bad_data = ['sample_point_2']
# drop, by default, returns a copy of the dataframe. If you want to update it in-place, use 'inplace=True'
cleaned_df = bad_df.drop(bad_data, axis=1)
# bad_df.drop('point_2', axis=1, inplace=True)
In [48]:
cleaned_df
Out[48]:
In [49]:
roi_data = {roi: np.hstack(data) for roi, data in df.iterrows()}
In [50]:
roi_data
Out[50]:
In [51]:
fig, ax = plt.subplots(figsize=(15, 10))
speckle_plot.combine_intensity_plotter(ax, roi_data.values())
if 'inline' in backend:
plt.show()
In [52]:
for bad in bad_data:
del masked_data[bad]
In [53]:
masked_data.keys()
Out[53]:
In [54]:
average_images = {k: np.sum(data, axis=0) for k, data
in masked_data.items()}
In [55]:
averages = {}
for name, data in average_images.items():
averages[name] = roi.circular_average(data, calibrated_center=image_center,
threshold=0, nx=100, pixel_size=pixel_size)
In [56]:
fig, axes = plt.subplots(nrows=len(averages.keys()), ncols=2, figsize=(12, 20))
im_kw = {'vmax': 1, 'cmap': cmap}
line_kw = {'c': 'blue', 'marker': 'o'}
keys = sorted(list(averages.keys()))
for k, ax in zip(keys, axes):
bin_centers, ring_averages = averages[k]
im, line = speckle_plot.circular_average(ax, average_images[k], ring_averages, bin_centers,
im_kw=im_kw, line_kw=line_kw,
im_title='Average Image Data for %s' % k)
if 'inline' in backend:
plt.show()
In [57]:
roi_kymograph = roi.kymograph(masked_data['sample_point_1'], rings_mask, num=2)
In [58]:
fig, ax = plt.subplots(figsize=(14,10))
ax.set_xlabel('Pixel')
ax.set_ylabel('Frame')
ax.set_title('Kymograph')
plot_kw = {
'cmap': cmap,
'vmin': 0,
'vmax': 6,
'interpolation': 'nearest',
'aspect': 'auto',
}
im, cb = speckle_plot.kymograph(ax, roi_kymograph, **plot_kw)
# speckle_plot.roi_kymograph_plotter(ax, roi_kymograph)
if 'inline' in backend:
plt.show()
In [59]:
fig, ax = plt.subplots((len(average_images)), figsize=(12, 12))
for label, image in average_images.items():
roi_pixel_values, index = roi.roi_pixel_values(image, rings_mask)
speckle_plot.rois_as_lines(ax, roi_pixel_values,
labels=[label]*len(image))
if 'inline' in backend:
plt.show()
In [60]:
import skxray
In [61]:
skxray.__version__
Out[61]:
In [62]:
import xray_vision
In [63]:
xray_vision.__version__
Out[63]:
In [ ]: