Commit 10e41d78 authored by Félix Hartmann's avatar Félix Hartmann
Browse files

New floating window with spatiotemporal heatmaps of the curvature and

the radius.
parent 57c72f70
......@@ -78,7 +78,7 @@ else:
from numpy import (array, arange, zeros_like, zeros, ones, linspace, newaxis, NaN, isnan,
hstack, ma, exp, log, arctan2, sqrt, sin, cos, pi, mean,
timedelta64, argmin, argsort, nonzero, flatnonzero, diff, gradient,
rad2deg, deg2rad,
rad2deg, deg2rad, quantile,
)
from numpy.linalg import norm
......@@ -237,11 +237,6 @@ def convert_pkl_to_hdf(pklfile, hdf5file, display_message=False, display_progres
base_dir_path = os.path.dirname(pklfile) + '/'
if display_message:
interekt.display_message(
"""Ancien format de données: CONVERSION (rootstem_data.pkl -> """
"""interekt_data.h5)\nChargement des résultats""")
# On charge le fichier pkl
data_out = load_results(pklfile)
......@@ -1539,6 +1534,168 @@ def show_time_series(tige_id=None):
fig_time.show()
def show_heatmaps():
"""
Displays the spatiotemporal variations of geometric variables as heatmaps.
These variables are:
- the curvature C
- the radius R
- the product CR
"""
# Force la fermeture du menu popup dans tk
interekt.floatmenuisopen = False
# Retrieving the tige id from the hdf5 file
tige_hdf_id = get_hdf5_tigeid(hdf5file, cur_tige)
# Retrieving useful data as spacetime 2D-arrays
radii = tiges_data.diam[cur_tige]
s, angles = load_postprocessed_data(hdf5file, tige_hdf_id)[2:5:2]
# Scaling
length_unit, inv_length_unit = "pix", r"pix$^{-1}$"
scale_cmpix = h5store.get_pixelscale(hdf5file)
if scale_cmpix is not None:
length_unit, inv_length_unit = "cm", r"cm$^{-1}$"
s *= scale_cmpix
times, time_units, time_label, graduation, excluded_image, mask = get_time_data()
cur_time = times[cur_image]
# If the steady-state image has to be excluded from time series
# (for instance because it has been taken a very long time after other images)
if excluded_image is not None:
times = times[mask]
s = s[mask]
radii = radii[mask]
angles = angles[mask]
# Retrieve the detection step from the h5 file
detection_step = h5store.get_detection_step(hdf5file)
if detection_step is None:
detection_step = DETECTION_STEP_DEFAULT
# Size of the averaging window
W = int(4 * round(radii.mean() / detection_step))
# Smoothing and computation of the curvatures
smooth_s = zeros_like(angles) - 3000
C = zeros_like(angles) - 3000 # curvatures
R = zeros_like(angles) - 3000 # radii
cutoff_left, cutoff_right = 2, 2
for i in range(len(angles)):
imax = flatnonzero(angles[i].mask == False).max() + 1
angles_i = angles[i, cutoff_left:imax-cutoff_right]
s_i = s[i, cutoff_left:imax-cutoff_right]
R_i = radii[i, cutoff_left:imax-cutoff_right]
C_i, smooth_angles_i = get_tige_curvature(angles_i, s_i,
smoothing=True, window_width=W)
smooth_s[i, :len(s_i)] = s_i
C[i, :len(C_i)] = C_i
R[i, :len(R_i)] = R_i
smooth_s = ma.masked_equal(smooth_s, -3000)
smooth_s = ma.masked_invalid(smooth_s)
C = ma.masked_equal(C, -3000)
C = ma.masked_invalid(C)
R = ma.masked_equal(R, -3000)
R = ma.masked_invalid(R)
imax = flatnonzero(smooth_s[-1].mask == False).max() + 1
smooth_s = smooth_s[:,:imax]
C = C[:,:imax]
R = R[:,:imax]
CR = C * R
t_min, t_max = times.min(), times.max()
s_max = smooth_s.max()
C_min, C_max = C.min(), C.max()
R_min, R_max = R.min(), R.max()
# Create the figure with the tige name
tige_name = h5store.get_tiges_names(hdf5file, tige_hdf_id)
start_title = strings["Spatiotemporal heatmaps"] + " " + strings["for organ"]
fig_heatmaps = mpl.figure(start_title + " %s"%(tige_name), figsize=(12, 10))
grid = mpl.GridSpec(3, 1)
cmap = mpl.get_cmap('jet')
# Curvature heatmap
axC = fig_heatmaps.add_subplot(grid[0,:])
axC.set_xlim([-s_max/100, s_max + s_max/100])
axC.set_ylim([-t_max/50, t_max + t_max/50])
axC.set_xlabel('s (%s)'%length_unit)
axC.set_ylabel(time_label)
linear_max = quantile(abs(C), 0.99)
cmeshC = axC.pcolormesh(smooth_s[-1], times, C,
norm=SymLogNorm(linthresh=linear_max, linscale=4),
cmap=cmap, rasterized=True)
cbarC = mpl.colorbar(cmeshC, use_gridspec=True, pad=0.01, extend='both')
cbarC.ax.tick_params(labelsize=8)
cbarC.ax.get_yaxis().labelpad = 15
cbarC.ax.set_ylabel(strings["curvature"] + " (%s)"%inv_length_unit, rotation=270)
if excluded_image == cur_image:
# If the steady-state image has to be excluded from time series and is the current
# image, just create empty Line2D.
time_select, = axC.plot(axC.get_xlim(), [None, None], 'k--', lw=2)
else:
# Otherwise mark the current time as a horizontal dotted line.
time_select, = axC.plot(axC.get_xlim(), [cur_time]*2, 'k--', lw=2)
def on_press(event):
if event.ydata is not None and axC.contains(event)[0]:
# image closest to clicked point
image = ((times - event.ydata)**2).argmin()
time_select.set_ydata([times[image], times[image]])
# redraw figure
fig_heatmaps.canvas.draw()
# Change figure on main frame
interekt.plot_image(image, keep_zoom=True)
# Radius heatmap
axR = fig_heatmaps.add_subplot(grid[1,:])
axR.set_xlim([-s_max/100, s_max + s_max/100])
axR.set_ylim([-t_max/50, t_max + t_max/50])
axR.set_ylabel(time_label)
cmeshR = axR.pcolormesh(smooth_s[-1], times, R,
cmap=cmap, rasterized=True)
cbarR = mpl.colorbar(cmeshR, use_gridspec=True, pad=0.01, extend='both')
cbarR.ax.tick_params(labelsize=8)
cbarR.ax.get_yaxis().labelpad = 15
cbarR.ax.set_ylabel(strings["radius"] + " (%s)"%length_unit, rotation=270)
# Curvature × Radius heatmap
axCR = fig_heatmaps.add_subplot(grid[2,:])
axCR.set_xlim([-s_max/100, s_max + s_max/100])
axCR.set_ylim([-t_max/50, t_max + t_max/50])
axCR.set_xlabel(strings["s"] + ' (%s)'%length_unit)
axCR.set_ylabel(time_label)
linear_max = quantile(abs(CR), 0.99)
cmeshCR = axCR.pcolormesh(smooth_s[-1], times, CR,
norm=SymLogNorm(linthresh=linear_max, linscale=4),
cmap=cmap, rasterized=True)
cbarCR = mpl.colorbar(cmeshCR, use_gridspec=True, pad=0.01, extend='both')
cbarCR.ax.tick_params(labelsize=8)
cbarCR.ax.get_yaxis().labelpad = 15
cbarCR.ax.set_ylabel(strings["curvature"] + " × " + strings["radius"], rotation=270)
fig_heatmaps.canvas.mpl_connect('button_press_event', on_press)
grid.tight_layout(fig_heatmaps)
fig_heatmaps.show()
def show_angle_and_curvature():
# Force la fermeture du menu popup dans tk
interekt.floatmenuisopen = False
......@@ -1552,13 +1709,12 @@ def show_angle_and_curvature():
# Retrieving the tige id from the hdf5 file
tige_hdf_id = get_hdf5_tigeid(hdf5file, cur_tige)
# Retrieving useful data
# Retrieving useful data ('d' is for 'dummy')
d, d, s, d, angles, d, d, lines = load_postprocessed_data(hdf5file, tige_hdf_id)
# 'd' is for 'dummy'
# Create the figure with the tige name
tige_name = h5store.get_tiges_names(hdf5file, tige_hdf_id)
start_title = _("Angles and curvatures for organ")
start_title = _("Angles and curvatures") + " " + strings["for organ"]
fig_curvature = mpl.figure(start_title + " %s"%(tige_name), figsize=(12, 10))
# Scaling
......@@ -2026,11 +2182,11 @@ def show_growth_length():
smooth_angles = zeros_like(angles) - 3000
smooth_s = zeros_like(angles) - 3000
curvatures = zeros_like(angles) - 3000
cutoff = 15
cutoff_left, cutoff_right = 2, 2
for i in range(len(angles)):
imax = flatnonzero(angles[i].mask == False).max() + 1
angles_i = angles[i, :imax-cutoff]
s_i = s[i, :imax-cutoff]
angles_i = angles[i, cutoff_left:imax-cutoff_right]
s_i = s[i, cutoff_left:imax-cutoff_right]
curvatures_i, smooth_angles_i = get_tige_curvature(angles_i, s_i, smoothing=True,
window_width=W)
smooth_angles_i -= smooth_angles_i[0]
......@@ -2081,7 +2237,7 @@ def show_growth_length():
cbar.ax.tick_params(labelsize=8)
axCurvatureMap.set_xlim([-s_max/100, s_max + s_max/100])
axCurvatureMap.set_ylim([-t_max/50, t_max + t_max/50])
axCurvatureMap.set_xlabel('s (%s)'%length_unit)
axCurvatureMap.set_xlabel(strings["s"] + ' (%s)'%length_unit)
axCurvatureMap.set_ylabel(time_label)
# The growth length will be stored in growth_data.
......@@ -3676,39 +3832,33 @@ class Interekt:
#floatmenu.add_command(label="Inverser la base", command=_reverse_tige)
self.floatmenu.add_command(label=_("Time series"),
#label=_("Séries temporelles"),
command=show_time_series)
self.floatmenu.add_command(label=strings["Spatiotemporal heatmaps"],
command=show_heatmaps)
self.floatmenu.add_command(label=_("Angles and curvatures"),
#label=_("Angles et courbures"),
command=show_angle_and_curvature)
self.floatmenu.add_command(label=_("Estimate growth rate"),
#label=_("Estimer la vitesse de croissance"),
command=show_growth_rate)
self.floatmenu.add_command(label=_("Estimate growth length"),
#label=_("Estimer la longueur de croissance"),
command=show_growth_length)
self.floatmenu.add_command(label=_("Estimate beta"),
#label=_("Estimer beta"),
command=show_beta)
self.floatmenu.add_command(label=_("Estimate gamma"),
#label=_("Estimer gamma"),
command=show_gamma)
self.floatmenu.add_command(label=_("Estimate B"),
#label=_("Estimer B"),
command=show_B)
self.floatmenu.add_command(label=_("Settings"),
#label=_("Réglages"),
command=show_tige_options)
self.floatmenu.add_command(label=_("Suppress the base"),
#label=_("Supprimer la base"),
command=remove_tige)
self.floatmenuisopen = False
......@@ -4752,6 +4902,7 @@ def set_string_variables():
global strings
strings["organ"] = _("organ")
strings["Organ"] = _("organ")
strings["for organ"] = _("for organ")
strings["image"] = _("image")
strings["image number"] = _("image number")
strings["photo"] = _("photo")
......@@ -4768,6 +4919,7 @@ def set_string_variables():
strings["tip angle"] = _("tip angle")
strings["degree"] = _("(deg)")
strings["curvature"] = _("curvature")
strings["Spatiotemporal heatmaps"] = _("Spatiotemporal heatmaps")
strings["growth length"] = _("growth length")
strings["residuals"] = _("residuals")
strings["initial"] = _("initial")
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment