changing pcolormesh plotting + fixing counter printing + adding positions...

changing pcolormesh plotting + fixing counter printing + adding positions extracting to load_tracking

doc/source/_static/example/load_tracking.py

@@ -3,10 +3,61 @@ import sys
 import argparse
 
 import numpy as np
+import pandas as pd
 from skimage.io import imread
+from skimage.measure import regionprops_table
 import tifffile as tifff
 
 
+# number of images
+NUM_IMAGES = 195
+
+def load_image(idx: int, PATH):
+    """
+    Freely inspired from napari tutorial: https://napari.org/tutorials/applications/cell_tracking.html
+    Load result of manual tracking as an image. 
+
+    Parameters
+    ----------
+    idx : int
+        Index of the image to load.
+
+    PATH : str
+        Path to image folder
+
+    Returns
+    -------
+    image : np.ndarray
+       The image specified by the index, idx
+    """
+    filename = osp.join(PATH, '01_GT/TRA', f'man_track{idx:0>3}.tif')
+    return imread(filename)
+
+
+def regionprops_plus_time(idx, stack):
+    """
+    Freely inspired from napari tutorial: https://napari.org/tutorials/applications/cell_tracking.html
+    Return the unique track label, centroid and time for each track vertex.
+
+    Parameters
+    ----------
+    idx : int
+        Index of the image to calculate the centroids and track labels.
+    stack : numpy.ndarray
+        stack of images containing the results of manual tracking
+
+
+    Returns
+    -------
+    data_df : pd.DataFrame
+       The dataframe of track data for one time step (specified by idx).
+    """
+
+    props = regionprops_table(stack[idx, ...], properties=('label', 'centroid'))
+    props['frame'] = np.full(props['label'].shape, idx)
+    return pd.DataFrame(props)
+
+
 def parse_args(args=None):
     """
     Parse arguments for main()
@@ -16,6 +67,16 @@ def parse_args(args=None):
 
     parser.add_argument('data_dir',
                         help='path of the data directory')
+    
+    parser.add_argument('-i', '--no_image',
+                        action="store_true",
+                        default=False,
+                        help='do not extract image')
+
+    parser.add_argument('-p', '--positions',
+                        action="store_true",
+                        default=False,
+                        help='extract positions')
 
     parsed_args = parser.parse_args(args)
 
@@ -28,26 +89,50 @@ def main(args=None):
     Freely inspired from napari tutorial: https://napari.org/tutorials/applications/cell_tracking.html
     """
 
-    # number of images
-    NUM_IMAGES = 195
-
     # load and check argument
     args = sys.argv[1:] if args is None else args
     parsed_args = parse_args(args)
 
     data_dir = osp.realpath(parsed_args.data_dir)
+    dont_get_image = parsed_args.no_image
+    get_positions = parsed_args.positions
 
+    # check data_dir
     if not osp.exists(data_dir):
         raise Exception("ERROR: the passed data directory does not exist. Aborting...")
 
     if not osp.isdir(data_dir):
         raise Exception("ERROR: the passed data directory is not a directory. Aborting...")
 
-    # extract images
-    stack = np.asarray([imread(osp.join(data_dir, '01', f't{i:0>3}.tif')) for i in range(NUM_IMAGES)])
+    if not dont_get_image:
+        print("extracting images...")
+        # extract images
+        stack = np.asarray([imread(osp.join(data_dir, '01', f't{i:0>3}.tif')) for i in range(NUM_IMAGES)])
+
+        #save_image
+        tifff.imsave(osp.join(data_dir,'stack.tif'), stack)
+
+    if get_positions:
+        print("extracting positions...")
+        stack = np.asarray([load_image(i,data_dir) for i in range(NUM_IMAGES)])
+        data_list = []
+        for idx in range(NUM_IMAGES): 
+            sys.stdout.write("\033[K")
+            print('image {}/{}'.format(idx + 1,NUM_IMAGES), end ='\r', flush = True)
+            data_list.append(regionprops_plus_time(idx,stack))
+        data_df_raw = pd.concat(data_list).reset_index(drop=True)
+
+        # sort the data lexicographically by track_id and time
+        data_df = data_df_raw.sort_values(['label', 'frame'], ignore_index=True)
+
+        # create the final data array: track_id, T, Z, Y, X
+        data = data_df.loc[:, ['label', 'frame', 'centroid-0', 'centroid-1', 'centroid-2']]
+
+        # rename columns
+        data.columns = ['track', 'frame', 'z', 'y', 'x']
 
-    #save_image
-    tifff.imsave(osp.join(data_dir,'stack.tif'), stack)
+        # save
+        data.to_csv(osp.join(data_dir,'positions.csv'), index=False)
 
 
 if __name__ == '__main__':

doc/source/user_guide/quickstart.rst

@@ -229,6 +229,8 @@ To run the script, open a terminal and run: ::
     cd <path_to_script_folder>
     python load_tracking.py <path_to_dataset_folder>
 
+You can also generate the positions by adding the flag :code:`-p`, it will generate the *positions.csv* file present in the archive. 
+
 Warning: if you try to open the generated tiff file with `Fiji <https://fiji.sc/>`_, you will see that the t and z dimensions are not separated. You will have to run "stack to hyperstack" with z=35 and t=195. But this is only if you want to see the file in Fiji, you don't need to do this for **Track Analyzer**!
 
 

track_analyzer/calculate.py

@@ -720,8 +720,8 @@ def average_local_vel(df, kernel_radius, dim=3):
     new_df = pd.DataFrame()
     coordinates = ['vx', 'vy', 'vz'] if dim == 3 else ['vx', 'vy']
     for i in df['frame'].unique():
-        sys.stdout.write('\033[2K\033[1G')
-        print('\rcomputing local average on frame ' + str(i), end='\r')
+        sys.stdout.write("\033[K")
+        print('\rcomputing local average on frame ' + str(i), flush=True, end='\r')
 
         df_frame = frame_groups.get_group(i)
         for coordinate in coordinates:

track_analyzer/plotting.py

@@ -27,6 +27,7 @@ import sys
 import os.path as osp
 import datetime
 import itertools
+import copy
 
 import matplotlib.pyplot as plt
 import matplotlib as mpl
@@ -105,8 +106,8 @@ def plot_traj(df, frame, data_dir, groups=None, image={'image_fn': None, 't_dim'
     The trajectory path can be removed to keep only the dots if traj_parameters['traj'] is False.
     It can be plotted in 3D with plot3D, elevation and angle set the 3D view
     """
-    sys.stdout.write('\033[2K\033[1G')
-    print('plotting frame {}'.format(int(frame)), end='\r')
+    sys.stdout.write("\033[K")  # go back to previous line
+    print('plotting frame {}'.format(int(frame)), flush=True, end='\r')
 
     # get config parameters
     plot_config = make_plot_config() if plot_config is None else plot_config
@@ -281,8 +282,8 @@ def plot_scalar_field(data_dir, df, data, field, frame, image={'image_fn': None,
     """Plot scalar field as colormap. The data needs to be generated before. """
 
     if not dont_print_count:
-        sys.stdout.write('\033[2K\033[1G')
-        print('plotting {} {}'.format(field, int(frame)), end='\r')
+        sys.stdout.write("\033[K")  # go back to previous line
+        print('plotting {} {}'.format(field, int(frame)), flush=True, end='\r')
 
     if plot_dir is None:
         plot_dir = osp.join(data_dir, field)
@@ -298,7 +299,7 @@ def plot_scalar_field(data_dir, df, data, field, frame, image={'image_fn': None,
     info = tpr.get_info(data_dir)
     image_size = [info['image_width'], info['image_height']]
     invert_yaxis = plot_config['invert_yaxis']
-    cmap = plt.get_cmap(cmap)
+    cmap = copy.copy(plt.get_cmap(cmap))  # make a shallow copy of cmap as modifying Matplotlib colormaps is deprecated
 
     # extract data
     X = data[frame]['X']
@@ -326,7 +327,10 @@ def plot_scalar_field(data_dir, df, data, field, frame, image={'image_fn': None,
     val_masked = np.ma.array(val, mask=np.isnan(val))
     [vmin, vmax] = [val_masked.min(), val_masked.max()] if vlim is None else vlim
     cmap.set_bad('w', alpha=0)  # set NAN transparent
-    C = ax.pcolormesh(X, Y, val_masked, cmap=cmap, alpha=0.5, vmin=vmin, vmax=vmax)
+
+    # shading=nearest so color value is centered on grid points
+    # for more info on pcolormesh behavior, see https://matplotlib.org/stable/gallery/images_contours_and_fields/pcolormesh_grids.html#sphx-glr-gallery-images-contours-and-fields-pcolormesh-grids-py
+    C = ax.pcolormesh(X, Y, val_masked, cmap=cmap, alpha=0.5, vmin=vmin, vmax=vmax, shading='nearest')  
 
     if show_axis:
         ax.grid(False)
@@ -351,8 +355,8 @@ def plot_vector_field(data_dir, df, data, field, frame, plot_on_field=None, dim=
     """ Plot vector field"""
 
     if not dont_print_count:
-        sys.stdout.write('\033[2K\033[1G')
-        print('plotting {} {}'.format(field, int(frame)), end='\r')
+        sys.stdout.write("\033[K")  # go back to previous line
+        print('plotting {} {}'.format(field, int(frame)), flush=True, end='\r')
 
     if plot_dir is None:
         plot_dir = osp.join(data_dir, field)
@@ -444,8 +448,8 @@ def plot_Voronoi(data_dir, df, frame, data, show_local_area=True,
     """
 
     if not dont_print_count:
-        sys.stdout.write('\033[2K\033[1G')
-        print('plotting {} {}'.format('voronoi', int(frame)), end='\r')
+        sys.stdout.write("\033[K")  # go back to previous line
+        print('plotting {} {}'.format('voronoi', int(frame)), flush=True, end='\r')
 
     if plot_dir is None:
         plot_dir = osp.join(data_dir, 'voronoi')
@@ -1277,7 +1281,8 @@ def plot_total_traj(data_dir, df, dim=3, plot_dir=None, plot_fn=None, plot_confi
         ylim = ax.get_ylim()
         ax.set_ylim(ylim[1], ylim[0])
 
-    filename = osp.join(plot_dir,
-                        datetime.datetime.now().strftime("%Y%m%d_%H%M%S") + '.png') if plot_fn is None else plot_fn
-    fig.savefig(filename, dpi=plot_config['dpi'], bbox_inches='tight')
+    fig.tight_layout()
+
+    filename = osp.join(plot_dir,'total_traj'+plot_config['format']) if plot_fn is None else plot_fn
+    fig.savefig(filename, dpi=plot_config['dpi'])
     plt.close(fig)