From dfef1eadd9e3c02fa7a29a1eada2d659ce4bfd25 Mon Sep 17 00:00:00 2001 From: Ivan Usov Date: Wed, 25 Jan 2023 18:41:08 +0100 Subject: [PATCH] Add panel_plot_data --- pyzebra/app/app.py | 2 + pyzebra/app/panel_plot_data.py | 284 +++++++++++++++++++++++++++++++++ 2 files changed, 286 insertions(+) create mode 100644 pyzebra/app/panel_plot_data.py diff --git a/pyzebra/app/app.py b/pyzebra/app/app.py index 5be05b2..c2726ac 100644 --- a/pyzebra/app/app.py +++ b/pyzebra/app/app.py @@ -15,6 +15,7 @@ from pyzebra.app import ( panel_hdf_param_study, panel_hdf_viewer, panel_param_study, + panel_plot_data, panel_spind, ) @@ -68,6 +69,7 @@ doc.add_root( panel_hdf_viewer.create(), panel_hdf_anatric.create(), panel_ccl_prepare.create(), + panel_plot_data.create(), panel_ccl_integrate.create(), panel_ccl_compare.create(), panel_param_study.create(), diff --git a/pyzebra/app/panel_plot_data.py b/pyzebra/app/panel_plot_data.py new file mode 100644 index 0000000..b20c7da --- /dev/null +++ b/pyzebra/app/panel_plot_data.py @@ -0,0 +1,284 @@ +import base64 +import io + +import numpy as np +from bokeh.layouts import column, row +from bokeh.models import ( + Button, + CheckboxGroup, + ColumnDataSource, + DataRange1d, + Div, + FileInput, + LinearColorMapper, + NumericInput, + Panel, + Select, + Spacer, + TextInput, +) +from bokeh.plotting import figure +from scipy import interpolate + +import pyzebra +from pyzebra import app +from pyzebra.app.panel_hdf_viewer import calculate_hkl + + +def create(): + measured_data_div = Div(text="Measured data:") + measured_data = FileInput(accept=".hdf", multiple=True, width=200) + + def plot_file_callback(): + flag_ub = bool(redef_ub_cb.active) + flag_lattice = bool(redef_lattice_cb.active) + + # Define horizontal direction of plotting plane, vertical direction will be calculated + # automatically + x_dir = list(map(float, hkl_in_plane_x.value.split())) + + # Define direction orthogonal to plotting plane. Together with orth_cut, this parameter also + # defines the position of the cut, ie cut will be taken at orth_dir = [x,y,z]*orth_cut +- delta, + # where delta is max distance a data point can have from cut in rlu units + orth_dir = list(map(float, hkl_normal.value.split())) + + # Where should cut be along orthogonal direction (Mutliplication factor onto orth_dir) + orth_cut = 0.0 + + # Width of cut + delta = hkl_delta.value + + # Load data files + md_fnames = measured_data.filename + md_fdata = measured_data.value + + for ind, (fname, fdata) in enumerate(zip(md_fnames, md_fdata)): + # Read data + try: + det_data = pyzebra.read_detector_data(io.BytesIO(base64.b64decode(fdata))) + except: + print(f"Error loading {fname}") + return + + if ind == 0: + if not flag_ub: + redef_ub_ti.value = " ".join(map(str, det_data["ub"].ravel())) + if not flag_lattice: + redef_lattice_ti.value = " ".join(map(str, det_data["cell"])) + + num_slices = np.shape(det_data["counts"])[0] + + # Change parameter + if flag_ub: + det_data["ub"] = np.array(redef_ub_ti.value.strip().split()).reshape(3, 3) + + # Convert h k l for all images in file + h_temp = np.empty(np.shape(det_data["counts"])) + k_temp = np.empty(np.shape(det_data["counts"])) + l_temp = np.empty(np.shape(det_data["counts"])) + for i in range(num_slices): + h_temp[i], k_temp[i], l_temp[i] = calculate_hkl(det_data, i) + + # Append to matrix + if ind == 0: + h = h_temp + k = k_temp + l = l_temp + I_matrix = det_data["counts"] + else: + h = np.append(h, h_temp, axis=0) + k = np.append(k, k_temp, axis=0) + l = np.append(l, l_temp, axis=0) + I_matrix = np.append(I_matrix, det_data["counts"], axis=0) + + if flag_lattice: + lattice = np.array(redef_lattice_ti.value.strip().split()) + else: + lattice = det_data["cell"] + + # Define matrix for converting to cartesian coordinates and back + alpha = lattice[3] * np.pi / 180.0 + beta = lattice[4] * np.pi / 180.0 + gamma = lattice[5] * np.pi / 180.0 + + # reciprocal angle parameters + beta_star = np.arccos( + (np.cos(alpha) * np.cos(gamma) - np.cos(beta)) / (np.sin(alpha) * np.sin(gamma)) + ) + gamma_star = np.arccos( + (np.cos(alpha) * np.cos(beta) - np.cos(gamma)) / (np.sin(alpha) * np.sin(beta)) + ) + + # conversion matrix: + M = np.array( + [ + [1, 1 * np.cos(gamma_star), 1 * np.cos(beta_star)], + [0, 1 * np.sin(gamma_star), -np.sin(beta_star) * np.cos(alpha)], + [0, 0, 1 * np.sin(beta_star) * np.sin(alpha)], + ] + ) + + # Convert all hkls to cartesian + hkl = [[h[:, :, :], k[:, :, :], l[:, :, :]]] + hkl = np.transpose(hkl) + hkl_c = np.matmul(M, hkl) + + # Convert directions to cartesian: + x_c = np.matmul(M, x_dir) + o_c = np.matmul(M, orth_dir) + + # Calculate y-direction in plot (orthogonal to x-direction and out-of-plane direction) + y_c = np.cross(x_c, o_c) + + # Calculate distance of all points to plane + Q = np.array(o_c) * orth_cut + N = o_c / np.sqrt(np.sum(o_c**2)) + v = np.empty(np.shape(hkl_c)) + v[:, :, :, :, 0] = hkl_c[:, :, :, :, 0] - Q + dist = np.abs(np.dot(N, v)) + dist = np.squeeze(dist) + dist = np.transpose(dist) + + # Find points within acceptable distance of plane defined by o_c + ind = np.where(abs(dist) < delta) + + # Project points onto axes + x = np.dot(x_c / np.sqrt(np.sum(x_c**2)), hkl_c) + y = np.dot(y_c / np.sqrt(np.sum(y_c**2)), hkl_c) + + # take care of dimensions + x = np.squeeze(x) + x = np.transpose(x) + y = np.squeeze(y) + y = np.transpose(y) + + # Get slices: + x_slice = x[ind[0], ind[1], ind[2]] + y_slice = y[ind[0], ind[1], ind[2]] + I_slice = I_matrix[ind[0], ind[1], ind[2]] + + # Define meshgrid limits for plotting (upper and lower limit + how fine mesh should be) + min_x = np.min(x_slice) + max_x = np.max(x_slice) + min_y = np.min(y_slice) + max_y = np.max(y_slice) + delta_x = xrange_step_ni.value + delta_y = yrange_step_ni.value + + # Create interpolated mesh grid for plotting + grid_x, grid_y = np.mgrid[min_x:max_x:delta_x, min_y:max_y:delta_y] + I = interpolate.griddata((x_slice, y_slice), I_slice, (grid_x, grid_y)) + + # Update plot + display_min_ni.value = 0 + display_max_ni.value = np.max(I_slice) * 0.25 + image_source.data.update( + image=[I.T], x=[min_x], dw=[max_x - min_x], y=[min_y], dh=[max_y - min_y] + ) + + plot_file = Button(label="Plot selected file(s)", button_type="primary", width=200) + plot_file.on_click(plot_file_callback) + + plot = figure( + x_range=DataRange1d(), + y_range=DataRange1d(), + plot_height=450, + plot_width=600, + tools="pan,wheel_zoom,reset", + ) + plot.toolbar.logo = None + + color_mapper = LinearColorMapper(nan_color=(0, 0, 0, 0)) + image_source = ColumnDataSource(dict(image=[np.zeros((1, 1))], x=[0], y=[0], dw=[1], dh=[1])) + plot.image(source=image_source, color_mapper=color_mapper) + + hkl_div = Div(text="HKL:", margin=(5, 5, 0, 5)) + hkl_normal = TextInput(title="normal", value="0 0 1", width=70) + hkl_delta = NumericInput(title="delta", value=0.1, mode="float", width=70) + hkl_in_plane_x = TextInput(title="in-plane X", value="1 0 0", width=70) + + def redef_lattice_cb_callback(_attr, _old, new): + if new: + redef_lattice_ti.disabled = False + else: + redef_lattice_ti.disabled = True + + redef_lattice_cb = CheckboxGroup(labels=["Redefine lattice:"], width=110) + redef_lattice_cb.on_change("active", redef_lattice_cb_callback) + redef_lattice_ti = TextInput(width=490, disabled=True) + + def redef_ub_cb_callback(_attr, _old, new): + if new: + redef_ub_ti.disabled = False + else: + redef_ub_ti.disabled = True + + redef_ub_cb = CheckboxGroup(labels=["Redefine UB:"], width=110) + redef_ub_cb.on_change("active", redef_ub_cb_callback) + redef_ub_ti = TextInput(width=490, disabled=True) + + def colormap_select_callback(_attr, _old, new): + color_mapper.palette = new + + colormap_select = Select( + title="Colormap:", + options=[("Greys256", "greys"), ("Plasma256", "plasma"), ("Cividis256", "cividis")], + width=100, + ) + colormap_select.on_change("value", colormap_select_callback) + colormap_select.value = "Plasma256" + + def display_min_ni_callback(_attr, _old, new): + color_mapper.low = new + + display_min_ni = NumericInput(title="Intensity min:", value=0, mode="float", width=70) + display_min_ni.on_change("value", display_min_ni_callback) + + def display_max_ni_callback(_attr, _old, new): + color_mapper.high = new + + display_max_ni = NumericInput(title="max:", value=1, mode="float", width=70) + display_max_ni.on_change("value", display_max_ni_callback) + + # xrange_min_ni = NumericInput(title="x range min:", value=0, mode="float", width=70) + # xrange_max_ni = NumericInput(title="max:", value=1, mode="float", width=70) + xrange_step_ni = NumericInput(title="x mesh:", value=0.01, mode="float", width=70) + + # yrange_min_ni = NumericInput(title="y range min:", value=0, mode="float", width=70) + # yrange_max_ni = NumericInput(title="max:", value=1, mode="float", width=70) + yrange_step_ni = NumericInput(title="y mesh:", value=0.01, mode="float", width=70) + + range_layout = row( + # xrange_min_ni, + # xrange_max_ni, + # Spacer(width=10), + xrange_step_ni, + Spacer(width=50), + # yrange_min_ni, + # yrange_max_ni, + # Spacer(width=10), + yrange_step_ni, + ) + cm_layout = row(colormap_select, display_min_ni, display_max_ni) + column1_layout = column( + row(measured_data_div, measured_data, plot_file), + plot, + column( + hkl_div, + row( + hkl_normal, hkl_delta, Spacer(width=10), hkl_in_plane_x, Spacer(width=50), cm_layout + ), + row(column(Spacer(height=7), redef_lattice_cb), redef_lattice_ti), + row(column(Spacer(height=7), redef_ub_cb), redef_ub_ti), + range_layout, + ), + ) + column2_layout = app.PlotHKL().layout + + hdf_div = Div(text="HDF DATA") + ccl_div = Div(text="CCL DATA") + tab_layout = row( + column(hdf_div, column1_layout), Spacer(width=70), column(ccl_div, column2_layout) + ) + + return Panel(child=tab_layout, title="plot data")