%% PALA_InVivoULM_example.m : Post Processing - filtering, localization and tracking for in vivo data
% Simple script to perform ULM on in vivo data
% IQ are loaded, filtered and sent into ULM processing.
% Bubbles are detected, selected, then localized, and paired into tracks.
% The result is a list of interpolated tracks to reconstructed a full brain vascularization.
%
% This code can be easily adapted for any other in vivo datasets.
%
% Created by Arthur Chavignon 25/02/2020
%
% DATE 2020.12.17 - VERSION 1.1
% AUTHORS: Arthur Chavignon, Baptiste Heiles, Vincent Hingot. CNRS, Sorbonne Universite, INSERM.
% Laboratoire d'Imagerie Biomedicale, Team PPM. 15 rue de l'Ecole de Medecine, 75006, Paris
% Code Available under Creative Commons Attribution-NonCommercial-ShareAlike 4.0 International (see https://creativecommons.org/licenses/by-nc-sa/4.0/)
% ACADEMIC REFERENCES TO BE CITED
% Details of the code in the article by Heiles, Chavignon, Hingot, Lopez, Teston and Couture.
% Performance benchmarking of microbubble-localization algorithms for ultrasound localization microscopy, Nature Biomedical Engineering, 2021.
% General description of super-resolution in: Couture et al., Ultrasound localization microscopy and super-resolution: A state of the art, IEEE UFFC 2018
clear
% DEFINE THE ADDONS DIRECTORY ON YOUR COMPUTER
% cd ./..
% PALA_addons_folder = [cd]; % location of the addons folder
% DEFINE THE DATA DIRECTORY ON YOUR COMPUTER
% addpath(genpath(PALA_addons_folder))
%% Adapt parameters to your data
% A few parameters must be provided by the user depending of your input images (size of pixel, wavelength)
% These parameters will be copied and used later during the creation of the ULM structure.
% in this example, UF.TwFreq = 15MHz, UF.FrameRateUF = 1000Hz;
UF.TwFreq=14.8;
UF.FrameRateUF=1000;
% Here you put the size of your data
SizeOfBloc=[160128500];
% Here you put the size of pixel in the prefered unit. It can be um, mm, m, wavelength, or an arbitrary scale.
ScaleOfPixel=[11];% [pixel_size_z, pixel_size_x]
% In that example, the size of pixels is lambda x lambda. The localization process will be
% performs in wavelength. For velocity rendering, velocities in [wavelength/s] will be converted in [mm/s].
% The imaging frame rate is required for velocity calculation, and temporal filtering.
framerate=1000;% imaging framerate in [Hz]
% Number of blocs to process
Nbuffers=80;% number of bloc to process (used in the parfor)
% If pixel sizes are in wavelength, lambda must be provided for a velocity maps in mm/s,
lambda=1480/(UF.TwFreq*1e3);
%% ULM parameters
% this script can be run using different scaling, it can be wavelength, pixel size, mm, um.
% In this example, input pixel are isotropic and equal to lambda (pixelPitch_x = pixelPitch_y = lambda)
% All size defined later are expressed in lambda
res=10;% final ratio of localization rendering, it's approximately resolution factor of localization in scale(1) units.
% for a pixel size of 100um, we can assume that ULM algorithm provides precision 10 (res) times
% smaller than pixel size. Final rendering will be reconstructed on a 10x10um grid.
ULM=struct('numberOfParticles',90,... % Number of particles per frame. (30-100)
'size',[SizeOfBloc(1)SizeOfBloc(2)SizeOfBloc(3)],... % size of input data [nb_pixel_z nb_pixel_x nb_frame_per_bloc]
'scale',[ScaleOfPixel1/framerate],...% Scale [z x dt], size of pixel in the scaling unit. (here, pixsize = 1*lambda)
'res',res,... % Resolution factor. Typically 10 for final image rendering at lambda/10.
'SVD_cutoff',[5SizeOfBloc(3)],... % SVD filtering, to be adapted to your clutter/SNR levels
'max_linking_distance',2,... % Maximum linking distance between two frames to reject pairing, in pixels units (UF.scale(1)). (2-4 pixel).
'min_length',15,... % Minimum allowed length of the tracks in time. (5-20 frames)
'fwhm',[11]*3,... % Size [pixel] of the mask for localization. (3x3 for pixel at lambda, 5x5 at lambda/2). [fmwhz fmwhx]
'max_gap_closing',0,... % Allowed gap in microbubbles' pairing. (if you want to skip frames 0)
'interp_factor',1/res,... % Interpfactor (decimation of tracks)
