Matteo Lenge

DEVELOPMENT AND VALIDATION OF INNOVATIVE ULTRASOUND FLOW IMAGING METHODS Candidate

Ultrasound is widely used for blood flow imaging because of the considerable advantages for the clinician, in terms of performance, costs, portability, and ease of use, and for the patient, in terms of safety and rapid checkup. The undesired limitations of conventional methods (1-D estimations and low frame-rate) are widely overtaken by new vector approaches, that offer detailed descriptions of the flow for a more accurate diagnosis of cardiovascular system diseases.

This PhD project concerns the development of novel methods for blood flow imaging. After studying the state-of-the-art in the field, a few approaches have been examined in depth up to their experimental validation, both in technical and clinical environments, on a powerful ultrasound research platform (ULA-OP).

Real-time novel vector methods implemented on ULA-OP were compared to standard Doppler methods in a clinical study. The results attest the benefits of the vector methods in terms of accuracy and repeatability.

Plane-wave transmissions were exploited to improve the transverse oscillation imaging method. Double oscillating fields were produced in large regions and exploited for the vectorial description of blood flow at high frame rates.

Blood flow maps were obtained by plane waves coupled to a novel velocity estimation algorithm operating in the frequency domain. The new method was demonstrated capable of high accuracy and reduced computational load by simulations and experiments (also in vivo). The investigation of blood flow inside the common carotid artery has revealed the hemodynamic details with unprecedented quality. A software solution implemented on a graphic processing unit (GPU) board was suggested and tested to reduce the computational time and support the clinical employment of the method.

Chapitre 4

Dans ce chapitre une méthode innovante pour la production de cartes bidimensionnelles du vecteur vitesse est présentée, apte à décrire le flux sanguin dans une section étendue de l'artère au cours de l'examen. La région d'intérêt est illuminée par un front d'onde plan transmis vers une xii R ÉSUM É direction préétablie. Pour chaque onde transmise, les échos rétrodiffusés sont recombinés grâce à une procédure hors ligne qui reconstruit l'image radiofréquence de toute la région d'intérêt. La différence de phase entre les images consécutives est estimée localement par un algorithme efficace qui travaille dans le domaine fréquenciel, dans le but de produire des cartes vectorielles avec une fréquence de quelques milliers d'images par seconde. L'évaluation des performances de la méthode a été obtenue grâce à de nombreuses simulations et études expérimentales, aussi bien in vitro que in vivo, par différentes configurations des paramètres mis en jeu.

Chapitre 5

La théorie des "oscillations transverses" vise à étendre en plusieurs dimensions les estimations du flux sanguin. La transmission d'ondes planes permet de générer des images à haute cadence d'imagerie. Dans ce chapitre, la transmission des ondes planes est combinée avec la théorie des oscillations transverses dans le but de générer des oscillations bidirectionnelles dans des régions étendues du champ de pression. Les performances de cette nouvelle modalité d'imagerie ont été évaluées grâce à une étude paramétrique. Son application pour les estimations de flux sanguin a été démontrée à travers d'autres simulations et expériences. Par rapport au chapitre précédent, la nouveauté de cette méthode concerne la stratégie de formation de l'image, alors que la carte vectorielle 2-D est obtenue à travers un estimateur déjà proposé précédemment.

Chapitre 6

L'algorithme d'estimation vectorielle dans le domaine fréquentiel, déjà décrit dans le chapitre 4, peut être implémenté de façon efficace sur des structures de calcul parallèle. Dans ce chapitre, une implémentation de la méthode sur une carte graphique (GPU) est décrite. Cette solution a été utilisée pour effectuer de nouvelles acquisitions in vivo. Les performances obtenues en termes de temps de calcul ont été confrontées avec celles obtenues sur un seul processeur (CPU). Enfin on présente des perspectives pour une version GPU optimisée de la méthode.

Sommario

Le tecniche ad ultrasuoni sono largamente utilizzate per descrivere il flusso sanguigno, in virtù dei numerosi vantaggi che offrono sia al clinico -in termini di prestazioni, rapporto costi-benefici, portabilità e relativa semplicità di utilizzo -sia al paziente, per l'innocuità delle indagini e la rapidità di diagnosi. Queste caratteristiche le rendono, per esempio, ideali a verificare la meccanica del cuore nel ciclo cardiaco, ed essenziali per il monitoraggio dello sviluppo del feto durante la gestazione.

Nell'ambito della diagnosi delle malattie cardiovascolari, i metodi Doppler ad ultrasuoni permettono di descrivere in modo dettagliato la dinamica del sangue all'interno dei vasi. L'analisi può riguardare un solo punto d'indagine (analisi spettrale) o un'intera regione bidimensionale (Color Flow Imaging). L'informazione sul flusso è solitamente combinata con quella morfologica, secondo modalità di presentazione duplex (B-Mode e analisi spettrale) o triplex (B-Mode, analisi spettrale e Color Flow Imaging). Tuttavia, queste strategie d'indagine presentano alcuni importanti limiti. In primis, può essere misurata soltanto la componente assiale della velocità del flusso e, per calcolarne il modulo, il clinico (di solito il medico angiologo) è obbligato a impostare manualmente l'angolo della direzione del flusso. Queste approssimazioni sono tollerate per le diagnosi all'interno di vasi estesi in larghezza e lunghezza, ma rendono difficile lo studio di sezioni che presentano biforcazioni e stenosi, dove il flusso può assumere percorsi impredicibili. In aggiunta, i metodi Doppler classici presuppongono un compromesso fra l'estensione della regione d'interesse e il frame-rate. I più moderni standard suggeriscono al clinico di utilizzare le modalità triplex per caratterizzare il flusso nella sua globalità, ma questa soluzione presenta lo svantaggio di restringere l'analisi quantitativa a un singolo volume campione, e il tasso di aggiornamento a poche decine di frame per secondo. Considerando l'elevata accelerazione che il cuore imprime al sangue nel picco sistolico -la velocità sanguigna all'interno dell'arteria carotide comune può variare da 0 a 1 m/s anche in meno di 20 ms -la risoluzione temporale della modalità triplex risulta insufficiente ad una descrizione del flusso esaustiva.

Questo progetto di ricerca ha riguardato lo sviluppo di metodi vettoriali vi SOMMARIO innovativi per blood flow imaging.

A partire da un'analisi dello stato dell'arte, sono state sperimentalmente validate, sia dal punto di vista tecnico sia in ambito clinico, alcune tecniche innovative di blood flow imaging "vettoriale" su regioni di dimensioni limitate, mediante un potente ecografo di ricerca: ULA-OP. Uno studio comparativo in ambito clinico ha evidenziato i vantaggi che presentano le metodiche Doppler vettoriali rispetto a quelle standard. È stato dimostrato che le nuove tecniche sono in grado di stimare in tempo reale il modulo della velocità senza dover ricorrere alla correzione manuale della direzione, e con benefici in termini di accuratezza e di riproducibilità.

Inoltre, è stata introdotta la trasmissione di onde piane per indagare la velocità vettoriale in regioni d'interesse estese. In quest'ambito sono stati sperimentati due approcci distinti. La velocità è stata infatti stimata sia attraverso un nuovo algoritmo di stima nel dominio della frequenza, sia sfruttando la teoria delle oscillazioni traversali del campo acustico ricevuto. L'accuratezza di tali metodi è stata verificata tramite simulazioni ed esperimenti. Sono state eseguite indagini in diverse regioni d'interesse, che hanno rivelato dettagli dell'emodinamica (ad esempio, il profilo del flusso al picco sistolico, o l'effetto dei movimenti delle valvole) con qualità elevata.

Nonostante che tali approcci siano in grado di ridurre significativamente il carico computazionale rispetto ad altri metodi di riferimento, l'implementazione in tempo reale rimane ancora proibitiva per i sistemi a ultrasuoni correntemente disponibili. Tuttavia, i metodi presentano elevate caratteristiche di parallelismo e, per l'algoritmo di calcolo nel dominio della frequenza, è stata implementata una soluzione su una scheda video (GPU) di ultima generazione, capace di ridurre significativamente il tempo computazionale.

Capitolo 1

La ricerca di nuovi strumenti per indagini di flusso basate sui metodi a ultrasuoni richiede l'impiego di strumenti e metodi allo stato dell'arte. In questo capitolo, sono descritti i fondamenti teorici della propagazione, della generazione e del rilevamento dell'onda acustica, e i vari approcci per la formazione dell'immagine anatomica. Vengono inoltre introdotte la teoria Doppler, alla base dei metodi classici a ultrasuoni utilizzati per le indagini di flusso, e le teorie più innovative in questo campo. Sono infine descritti gli strumenti hardware e software (l'ecografo di ricerca ULA-OP e il simulatore SIMAG) utilizzati nel corso dell'attività di ricerca, e le applicazioni cliniche più significative.

vii

Capitolo 2

Lo sviluppo di tecniche innovative per l'elaborazione del segnale ecografico comporta tre fasi distinte: dapprima l'idea, poi lo sviluppo e l'implementazione dell'algoritmo, infine la validazione tramite simulazioni ed esperimenti sia in phantom sia in vivo. Per semplificare l'attività di sviluppo è stata messa a punto una piattaforma software per il processing di dati ecografici, costituita da una collezione di blocchi funzionali implementati come classi in MATLAB. La piattaforma è stata utilizzata per applicazioni finalizzate a individuare nuovi marker diagnostici per le patologie a carico del sistema cardiovascolare. In questo capitolo sono presentati due studi basati su tale piattaforma, che riguardano la misura volumetrica del flusso sanguigno e la misura della dilatazione flusso-mediata delle pareti arteriose (flow-mediated dilation).

Capitolo 3

La misura della velocità di picco sistolico (PSV), ampiamente utilizzata per diagnosticare la stenosi dell'arteria carotide, è influenzata dalla dipendenza dall'angolo Doppler. Questo limite ha incoraggiato lo sviluppo di sistemi semi-automatici d'indagine "vettoriale" dei vasi. Tali metodiche, pur presentando evidenti vantaggi tecnici, non sono state tuttora sperimentate clinicamente in maniera approfondita. In questo capitolo viene presentato uno studio clinico in cui due metodiche vettoriali sono impiegate per la misura del PSV nelle arterie carotidi comune e interna, sia in volontari sani, sia in pazienti con differenti gradi di stenosi. Le misure sono messe in relazione con quelle ottenute dai metodi Doppler convenzionalmente utilizzati nella pratica clinica. I risultati di questo studio suggeriscono che l'introduzione di metodi vettoriali nelle macchine commerciali è matura per risolvere efficacemente i problemi di sovrastima tipici dei metodi standard che dipendono dall'angolo Doppler.

Capitolo 4

In questo capitolo è presentato un nuovo algoritmo per la produzione di mappe bidimensionali del vettore velocità, adatto a descrivere il flusso sanguigno in una sezione estesa dell'arteria sotto indagine. La regione d'interesse è illuminata da un fronte d'onda piano trasmesso in una direzione prefissata. Per ogni onda trasmessa, gli echi retrodiffusi sono ricombinati tramite una procedura che ricostruisce fuori linea l'immagine a radiofrequenza dell'intera regione d'interesse. La differenza di fase fra immagini consecutive è stimata localmente tramite un nuovo ed efficiente algoritmo che lavora nel dominio della frequenza, con lo scopo di produrre mappe vet-viii SOMMARIO toriali con una frequenza di migliaia di frame per secondo. La valutazione delle prestazioni del metodo è stata ottenuta tramite numerose simulazioni e prove sperimentali, sia in phantom sia in vivo, ottenute al variare dei parametri in gioco.

Capitolo 5

La teoria delle "oscillazioni trasversali" mira a estendere in più dimensioni le investigazioni del flusso sanguigno. La trasmissione di onde piane consente di generare immagini ad alto frame-rate. In questo capitolo, la trasmissione di onde piane è combinata con la teoria delle oscillazioni trasversali allo scopo di generare oscillazioni bidirezionali del campo acustico in regioni estese. Le prestazioni di questa nuova modalità di imaging sono state valutate tramite uno studio parametrico. La sua applicazione per stime di flusso sanguigno è stata dimostrata mediante ulteriori simulazioni ed esperimenti. Rispetto al precedente capitolo, la novità di questo metodo riguarda la strategia di formazione dell'immagine. La mappa vettoriale 2-D viene ottenuta tramite uno stimatore già precedentemente validato.

Capitolo 6

L'algoritmo di stima vettoriale nel dominio della frequenza (descritto nel capitolo 4) può essere implementato in maniera efficiente su architetture di calcolo parallelo. In questo capitolo, è descritta un'implementazione del metodo su una scheda grafica (GPU). Questa soluzione è stata utilizzata per elaborare nuove acquisizioni in vivo. Le prestazioni ottenute in termini di carico computazionale sono state messe a confronto con quelle ottenute tramite una soluzione su singolo processore (CPU). Infine, vengono discusse le prospettive per una implementazione del metodo ulteriormente ottimizzata.

Résumé

Les méthodes ultrasonores fournissent des images à haute résolution spatiale et temporelle avec une totale innocuité pour le patient et un diagnostic temps réel. Ces caractéristiques les rendent très performantes, par exemple, pour vérifier le suivi du cycle cardiaque, et essentiels pour le monitorage continu du développement du foetus pendant les différentes semaines de grossesse.

Dans le domaine des mesures de flux sanguin pour le diagnostic des maladies cardiovasculaires, les méthodes Doppler ultrasonores sont amplement utilisées pour décrire la dynamique du sang dans les vaisseaux. Normalement, elles sont employées en combinaison avec d'autres modalités dédiées à la description de l'anatomie (ex. : B-mode, en modalité duplex) et aussi du flux (Color Flow Mapping, en modalité triplex). Toutefois, ces stratégies d'imagerie présentent quelques limites importantes. Premièrement, on ne peut mesurer que la composante axiale de la vitesse du flux, et pour en calculer le module le clinicien (généralement l'angiologue) est obligé d'afficher manuellement l'angle entre la propagation de l'onde ultrasonore et de direction du flux. Ces approximations sont acceptables pour le diagnostic de vaisseaux ayant une forme linéaire simple, mais elles rendent difficile l'étude des zones qui présentent des bifurcations et des sténoses, où le flux peut prendre des parcours imprévisibles avec en particulier des turbulences. De plus, les méthodes Doppler classiques présupposent un compromis entre la taille de la région d'intérêt et la cadence d'imagerie. Les standard les plus modernes suggèrent au clinicien d'utiliser les modalités triplex pour caractériser le flux dans sa globalité, mais cette solution a le désavantage de restreindre l'analyse quantitative à un seul volume échantillonné, et un taux de mise à jour d'à peine une dizaine d'image par seconde. En considérant l'étendue élevée des accélérations que le coeur imprime au sang au moment du pic systolique -la vitesse sanguine à l'intérieur de l'artère carotidienne commune peut varier de 0 à 1 m/s en moins de 20 ms -la résolution temporelle de la modalité triplex est insuffisante pour la description complète du flux.

Dans ce projet de recherche, suite à un état de l'art des techniques

x R ÉSUM É existantes, des solutions innovatrices pour la mesure du flux sanguin ont été proposées. Une étude comparative clinique a mis en évidence les avantages que les méthodes Doppler vectorielles présentent par rapport aux méthodes standards. Il a été démontré que les nouvelles techniques sont capables d'estimer en temps réel la vitesse en deux dimensions (le vecteur 2-D de la vitesse), en évitant l'intervention manuelle pour indiquer la direction du flux et avec des bénéfices en termes de soin et de reproductibilité. En outre, la transmission d'ondes planes a été introduite pour étendre ces techniques à des régions d'intérêt étendues. La théorie des oscillations transversales a été employée avec la transmission d'ondes planes, afin d'obtenir le vecteur vitesse du flux sanguin. Un nouvel algorithme pour l'estimation de la vitesse dans le domaine fréquentiel a été développé et validé. Ces deux méthodes visent à dépasser le compromis demandé par les modalités triplex, puisque les cartes bidimensionnelles du vecteur vitesse à haute cadence sont obtenues dans des régions d'intérêt étendues. Les performances de ces méthodes ont été évaluées à travers des simulations et des expériences dans différentes région de l'anatomie du corps ont été faites et ont révélé des détails de l'hémodynamique (par exemple, le profil du flux au pic systolique, ou l'effet du mouvement des valves) avec une qualité sans précédent.

Même si de telles méthodes sont capables de réduire de façon significative la charge de calcul par rapport à d'autres méthodes de référence, l'implantation en temps réel reste encore prohibitive pour les systèmes ultrasonores couramment utilisés en clinique. Toutefois, les méthodes présentent des caractéristiques élevées de parallélisme et, pour la méthode de calcul dans le domaine fréquentiel, une solution sur une carte vidéo (GPU) de dernière génération a été implantée, dans le but de réduire le temps de calcul et soutenir l'expérimentation clinique de la méthode sur les systèmes de prochaine génération.

Chapitre 1

La recherche de nouveaux instruments pour des recherches de flux basées sur des méthodes ultrasonores demande l'emploi de différents instruments et méthodes. Dans ce chapitre sont décrits les fondements théoriques de la propagation, de la génération et la mesure de l'onde acoustique, et les différentes approches pour la formation de l'image anatomique. On introduit en outre la théorie Doppler à la base des méthodes classiques à ultrasons utilisées pour l'étude du flux, et les théories les plus novatrices en ce domaine. Enfin sont décrits les instruments hardware et logiciels (l'échographe de recherche ULA-OP et le simulateur SIMAG) utilisés au cours de l'activité de recherche, et les applications cliniques pour le diagnostic des maladies xi cardiovasculaires à charge de l'artère carotidienne.

Chapitre 2

Le développement de techniques innovatrices pour l'élaboration du signal échographique peut être divisé en trois phases distinctes : d'abord l'idée, puis le développement de l'algorithme et son implémentation, enfin la validation à travers simulations et expériences aussi bien in vitro que in vivo. Pour simplifier l'activité de développement une plateforme logicielle pour le traitement de données échographiques, constituée d'une collection de blocs fonctionnels a été implémentés sous formes classes en MATLAB. La plateforme a été utilisée pour des applications visant au développement de nouveaux marqueurs diagnostiques pour les pathologies affectant le système cardiovasculaire. Dans ce chapitre deux études de recherche sont présentées qui concernent la mesure volumétrique du flux sanguin et celle de la dilatation flux-médiate des parois artérielles.

Chapitre 3

La mesure de la vitesse de pic systolique (PSV), largement utilisée pour diagnostiquer la sténose de l'artère carotidienne, est influencée par la dépendance de l'angle Doppler. Cette limite encourage le développement de systèmes semi-automatiques pour l'étude des vaisseaux. Plusieurs méthodes capables de mesurer correctement la vitesse sanguine en plusieurs directions ont été définies. Même si ces méthodes présentent d'évidents avantages par rapport aux techniques déjà reconnus par la communauté scientifique, elles n'ont pas encore été expérimentées de manière approfondie. Dans ce chapitre une étude clinique est présentée, où deux méthodes vectorielles sont employées pour la mesure du PSV dans les artères carotidiennes commune et interne, aussi bien sur des sujets sains que sur des patients avec différents degrés de sténose. Les mesures sont mises en relation avec celles qui ont été obtenues par les méthodes Doppler conventionnellement utilisées dans la pratique clinique. Les résultats de cette étude suggèrent que l'introduction des méthodes vectorielles dans les machines commerciales pourrait résoudre efficacement les problèmes de surestimation typiques des méthodes standards qui dépendent de l'angle Doppler.

Abbreviations and Acronyms

ABUR

Introduction

Ultrasound imaging is featured by high resolution and total safety for the patients. Among other benefits, those are the essential characteristics for supervising the fetal development inside the mother's womb, or to verify the heart mechanics during the cardiac cycle. In blood flow imaging, Doppler ultrasound is widely used to check the vascular functionality. The clinician (often a specialist in angiology) uses duplex or triplex modalities, which allow explorations of both anatomical and functional aspects, i.e. showing the vessels and the blood flowing inside them. However, such safe and prompt strategies present some limitations. In primis, only the axial component of the velocity is estimated and, to calculate the velocity magnitude, the clinician must manually (and thus, approximately) adjust the flow direction. Such approximations are tolerated for straight and large vessels, but make difficult the investigation of, e.g., bifurcations and stenoses, where the flow may assume unpredictable patterns. Furthermore, classic Doppler measurements presume a tradeoff between the extent of the region under investigation and the frame-rate. The modern standards suggest the clinicians to use triplex modalities, i.e. the B-mode together with spectral Doppler and Color Flow Imaging. This solution affords the description of the complete hemodynamics, but restricts the velocity assessment in a single sample volume, and the rate in few tens of frames per second. Considering the blood acceleration, e.g. the blood velocity change in the common carotid artery from 0 to 1 m/s in less than 20 ms, the temporal resolution of triplex mode is insufficient for the complete description of flow.

Then, the use of plane waves for vector Doppler investigation over large region of interest has been introduced. The transverse oscillations theory was adapted to plane-wave imaging, in order to obtain vector blood flow estimations. Furthermore, a novel velocity estimation algorithm operating in the frequency domain was developed and validated. These methods aim at overcoming the tradeoff required by conventional triplex modes, because 2-D vector maps of the blood flow are obtained in extended regions of interest and at high frame-rate. The accuracy of the methods was demonstrated by simulations and experiments. Applications in many areas of the body were tested, and revealed hemodynamic details (e.g., the flow at systolic pulses, the effects of valve's movement) with unprecedented quality. Although significantly reducing the computational load with respect to other methods, the real-time implementation of these two algorithms remains prohibitive for the current ultrasound systems. Nevertheless, the methods are high parallelizable and, for the frequency-domain method, a software solution implemented on a GPU board was suggested and tested to reduce the computational time and support the clinical employment of the method on next generation systems.

Thesis Outline

The manuscript is organized as follows:

Chapter 1: The state-of-the-art in medical ultrasound imaging is described, including the ultrasound basics, some clinical applications and the hardware and software materials available for this project.

Chapter 2: A Matlab ® software platform constituted by a collection of classes is proposed, in order to allow a smart environment for the development of new flow imaging US techniques, especially for the early detection of cardiovascular diseases.

Chapter 3: Vector and spectral Doppler methods are compared in a clinical study oriented to the measurement of the peak systolic velocity (PSV) in the carotid artery; the comparison highlights the advantages of using vector "angle independent" techniques with respect to conventional methods. Chapter 4: A new efficient vector estimator operating in the frequency domain is here presented and validated by simulations and experiments; the study has produced interesting in-vivo results that describe the complex flow in the carotid and the radial arteries with unprecedented quality.

Introduction

Ultrasound Basics

Acoustics is a branch of classical mechanics that concerns the elastic propagation of mechanical waves, that are generated from an oscillating body and propagate through an elastic medium (gas, liquid, or solid). The vibration around the equilibrium position of the many particles of the medium globally gives rise to an elastic wave that propagates in the longitudinal and in the transverse directions. The propagation of acoustic waves from a place to another is merely a consequence of mechanics and the properties of the medium [START_REF] Feynman | The Feynman Lectures on Physics[END_REF]. Ultrasounds (US) are particular acoustic waves that vibrate with frequencies higher than 20 kHz, i.e. higher than the audible band of the human ear (Fig. 1.1).

Ultrasound Propagation

The Wave Equation

The wave equation is a second-order linear partial differential equation that describes the dynamics of waves. The scalar form of the ultrasound wave that propagates in the matter along the longitudinal (or axial) direction (z ), is:

∂ 2 T ∂z 2 = ρ β ∂ 2 T ∂t 2 (1.1)
where β [N/m 2 ] and ρ [kg/m 3 ] are respectively the elastic constant and the volumetric density of the medium, and T [N/m 2 ] is the stress. The solution of (1.1) brings to the well known plane wave function, expressed here in a complex form T (z, t) = T 0 e j(ωt±kz) , (1.2)

Ultrasound Basics

9

with the wave number (k), the pulsation (ω) and the wavelength (λ) are respectively defined as follows:

k = ω c = ω ρ β , ω = 2πf , λ = c f . (1.
3)

It can be noticed in (1.2) that the phase depends on the space (z) and the time (t). From (1.3), the speed of sound c is:

ρ β = 1 c 2 ⇔ c = β ρ . (1.4)
Equation (1.4) underlines the inherent dependency of wave mechanics from the physical and chemical characteristics of the medium, i.e. the coupling and the mechanical inertia of particles, which are macroscopically represented by the elasticity and the density. Although c is not homogeneous in the medium, in medical ultrasound applications it is usually calculated by using the mean density of the medium ρ 0 . The most common acoustic properties of biological media are reported in the Tab. 1.1. Among them, the acoustic impedance Z: Z = ρ 0 c [kg / (m 2 s)] summarizes how much ultrasound propagation is "impeded" by a particular medium.

Similarly to other waves, the ultrasound wave propagation is governed by several phenomena: the reflection, the refraction, the scattering and the attenuation.

Reflection and Refraction

The ultrasound wave propagates along a straight line if the medium is homogeneous, i.e. if the physical and chemical properties are independent from space and time. When the wave front proceeds from a medium (with Z = Z 1 ) to another medium (with Z = Z 2 ), the incident energy is partially reflected by the interface, while the remaining part is transmitted toward deeper regions. If the interface roughness is larger than the wavelength λ, the phenomenon follows the relationship of the geometrical optics, the Snell's law:

sin θ i sin θ t = c 1 c 2 (1.5)
where θ i and θ t are the angles of the incident and the transmitted waves, respectively, and c 1 and c 2 are the speeds of sound in the two media. At the interface, the reflected and the transmitted energies are proportional respectively to the reflection (R) and the transmission (T ) coefficients, [START_REF] Culjat | A review of tissue substitutes for ultrasound imaging[END_REF].

defined as follows:

R = (Z 2 cos θ i -Z 1 cos θ t ) 2 (Z 2 cos θ i + Z 1 cos θ t ) 2 , T = 4Z 2 Z 1 cos 2 θ i (Z 2 cos θ i + Z 1 cos θ t ) = 1 -|R| 2 (1.6)

Scattering

If the interface roughness is smaller or comparable to the incident wavelength λ, a part of the energy is transmitted and the other part is spread isotropically in all directions, according to the Rayleigh scattering [START_REF] Rayleigh | On the light from the sky, its polarization and colour[END_REF]. The phenomenon is represented by the back-scattering cross-section (σ) of the medium, defined as:

σ = S/I, (1.7)
where S is the total spread power and I the intensity of the incident wave. This parameter describes the effective power of the echo signal available for an imaging system [START_REF] Ramalli | Development of Novel Ultrasound Techniques for Imaging and Elastography : From Simulation to Real-Time Implementation[END_REF].

Ultrasound Basics 11 Attenuation

Reflection, refraction and scattering phenomena occur during the wave propagation and give rise to the dispersion of energy, i.e. the wave attenuation. Other accountable accomplices of that attenuation are the compression and expansion processes (thermo-elastic effect) and sliding movements of particles (viscous effect), that convert energy from mechanical to heat. Those effects are contemplated in the attenuation coefficient α [W/m 2 ] that has an inverse exponential relationship with the wave intensity I [dB / (cm MHz)]:

I = I 0 e -2α A z . (1.8)
The wave attenuation is usually related in dB/m:

α AdB = 10 W log I I 0 (1.9)
where I 0 is the initial intensity and W [m] is the tissue width or thickness.

Ultrasound Generation and Detection

The piezoelectric effect is a physical property of some materials that are typically used to generate and to detect the ultrasound waves. If a piezoelectric material is excited by an electric voltage, it modifies its dimensions and, vice versa, the alteration of its dimensions produces an electric voltage. The ultrasound transducer exploits this principle and employs an electronic circuit that interfaces the system with the piezoelectric sensor. The purpose is to convert the energy of an electrical wave into the energy of a pressure wave (i.e. the generation of the US wave), and vice versa for the detection of the back-scattered US echoes. For more details on design principles of transducers and probes, the reader is referred to [START_REF] Szabo | Diagnostic Ultrasound Imaging: Inside Out[END_REF].

Acoustic Beams

Single-element Transducer Beam

The acoustic field produced by a transducer in the surrounding space can be divided in two distinct zones (Fig. 1.2): the near field, or Fresnel zone, and the far field, or Fraunhofer zone. The source geometry identifies a border surface between the two zones at a depth z:

z = r 2 λ (1.10)
where r is the transducer radius and λ is the transmitted wavelength. The elements of the transducer generate waves that interfere in a constructive More important is that, in the near field, the intensity of the field oscillates very rapidly. The coordination system is also indicated: x for the lateral, y for the elevation, and z for the axial directions.

or a destructive manner, according to optics laws. The diffraction theory demonstrates that, in the Fraunhofer zone, the beam has approximately the same shape of the Fourier transform of the beam located on the active aperture [START_REF] Szabo | Diagnostic Ultrasound Imaging: Inside Out[END_REF]. For a finite rectangular aperture, the beam is a sinc-shaped function and it presents side lobes of lower intensity that flank the mainlobe.

Array Transducer Beam

In ultrasound imaging, the transducers are typically arranged in arrays, that are periodical structures composed by several elements, separated by a periodical distance (i.e., the pitch). For the array of transducers, the beam field is a combination of contributions by the entire rectangular aperture (the array factor ) and the finite single elements (element factor ). The array factor is described in the Fourier domain by the convolution between a Dirac-comb function and a Sinc function (i.e., the DFT of the rectangular aperture). The element factor is a Sinc function that takes account of the finite dimension of the single elements. The final field beam is the product between the array and element factors. Exploiting this basic principle, several array transducers have been developed to extend ultrasonic investigations towards many different anatomical areas.

Linear array transducers have the elements placed on a straight line to scan rectangular ROIs, such as vessel sections [Fig. 1.3(a)]. Convex array transducers have elements located on a convex surface, to scan wide circular sectors like the uterus in Fig. 1.3(b). Phased array transducers have little elements and boast high steering features. In echocardiographic applications [Fig. 1.3(c)] the beam is steered to cover a wide sectorial ROI. With similar performances in steering, but different applications, the transvaginal array transducers are used for gynecologic ultrasonography [Fig. 1.3(d)]. For more details on array transducers, the interested reader is referred to [START_REF] Szabo | Diagnostic Ultrasound Imaging: Inside Out[END_REF].

B-mode Imaging

The image of Fig. 1.3(a) was achieved by an imaging modality which is commonly available on commercial systems. The back-scattered echo signals are processed to obtain images which pixels having an amplitude proportional to the acoustic properties of the structures under investigation. The amplitude is associated to the brightness of a gray-scaled image, that's the reason why this modality is called B-mode (Brightness mode).

In transmission (TX), the elements of the transducer are excited by the electronic system, in order to generate short pulses (bursts) made by some radio-frequency sinusoidal cycles. Bursts are windowed using different weighting functions (Hamming, Hanning, etc.) to reduce the frequency spread, and periodically generated at Pulse Repetition Frequency rate (PRF). For each TX event, the beam is focalized in a limited area to concentrate the energy of the ultrasound waves. A concave acoustic length is superimposed on the array to provide the focalization on the elevation plane, i.e. the plane perpendicular to the focal plane. The focalization in the focal plane is obtained by the electronic focusing technique, in which each element of the transducer aperture is excited by delayed signals (Fig. 1.4). The focal point is selected in the azimuth plane by a steered beam. The steering lets different two-dimensional ROIs to be explored avoiding any mechanical movement, but introduces an electric cross-talk effect between neighbor elements that entails in an apparent widening of the transmitting elements [START_REF] Ramalli | Development of Novel Ultrasound Techniques for Imaging and Elastography : From Simulation to Real-Time Implementation[END_REF]. The active aperture (i.e., the ensemble of the probe elements involved in the transmission) is windowed by a set of weightings that modulates the elements amplitude (the so-called apodization), in order to limit the contribution of the grating lobes.

In most medical imaging techniques it is convenient to use the same probe to transmit and receive, and different active apertures can be defined. When the ultrasound wave is transmitted and reaches a point target with comparable dimensions, it is back-scattered with a spherical wavefront, according to scattering phenomenon described on p. 10. The temporal length of bursts is in the order of µs. The echoes are received by the transducer elements, converted to electric signals and, taking into account the Shannon theorem, sampled and stored by ultrasound system. To repeat the transmission-receive event, it must be considered that the ultrasound wave propagates forward and backward, with a velocity that is about 1500 m/s in the soft tissue. Considering a 50-depth region of interest, the maximum PRF is around 15 kHz and, consequently, the Pulse Repetition Interval time (T P RF = 1/P RF ) is even 100 times higher than the burst length. The time span between two TX events is conveniently employed for signal processing oriented to the quality improvement of the final images and to extract information. Considering the spherical wave-fronts scattered back by a single point, the beamforming consists in applying precise delay curves to the echo-signals, to select the samples due to the back-scattered contribution of the point, and to coherently sum the samples amplitudes [START_REF] Van Veen | Beamforming: a versatile approach to spatial filtering[END_REF]. This oper-Figure 1.4: Linear transducer array with lens focusing in elevation and electronically-controlled focusing in azimuth [START_REF] Wells | Ultrasonic imaging of the human body[END_REF].

ation, the so-called delay and sum beamforming, is applied at many depths in modern digital systems (dynamic focusing). The signals obtained by the beamforming processing are hereinafter referred to the radio-frequency signals (RF s ).

The radio-frequency signals are demodulated to obtain the base-band signals (IQ s ) in complex notation. This operation can be done by sampling directly the IQ s , or by applying the Hilbert transform [START_REF] Evans | Doppler Ultrasound: Physics, instrumentation and signal processing[END_REF]. Then the amplitude is calculated and rearranged to a matrix form IQ mn , where m and n are the indices along depths (axial direction) and lines (lateral direction), respectively. According to a gray-scaled map, each element (m, n) is displayed as a bright point. Each pixel of the final digital image assumes a value in brightness which is proportional to the back-scattered echoes produced by the volume that it represents.

To reconstruct all the lines of the image, several transmissions are used moving the active aperture along the transducer elements. Also multifocusing techniques are implemented, they exploit more transmissions with different depth of focalization. The main drawback is the drastic limitation of the frame-rate (F R). For example, in Fig. 1.3(d), 2 focal points are used to scan a 16 cm depth ROI and the final frame-rate (36 Hz) is just enough to display.

Medical Ultrasound Imaging High Frame-Rate Imaging

New imaging techniques has been proposed in order to reconstruct images of broad regions at high frame-rate (HFR) [START_REF] Lu | 2D and 3D high frame rate imaging with limited diffraction beams[END_REF]. One of the most prolific HFR methods is the ultrafast compound imaging, firstly proposed by Tanter et al. [START_REF] Tanter | Ultrafast compound imaging for 2-D motion vector estimation: application to transient elastography[END_REF], developed by Montaldo et al. [START_REF] Montaldo | Coherent plane-wave compounding for very high frame rate ultrasonography and transient elastography[END_REF] and extended by Bercoff et al. [START_REF] Bercoff | Ultrafast compound Doppler imaging: providing full blood flow characterization[END_REF] and Macé et al. [START_REF] Mace | Functional ultrasound imaging of the brain: theory and basic principles[END_REF]. Instead of using the conventional focused transmissions, plane-waves inclined at many steering angles are transmitted. Each transmission produces a low-resolution image of a broad ROI. To obtain the final image at high resolution, the back-scattered echoes are coherently recombines (i.e. the images are compounded).

Considering the coordination system in Fig. 1.5(a), the plane pulsedwave is transmitted into the medium and the echo signals RF ( p , t) are received by the P elements ( p : p = 1, 2, ..., P ) of the probe. The traveling time from the transducer to the point (x, z) and back to the single element p of the transducer is: (1.12)

τ ( p , x, z) = (z + z 2 + (x -p ) 2 )/c. ( 1
In this way, one entire frame is produced from a single acoustic pulse with a parallel processing. The lack of focusing in transmission causes a low contrast in the final image, because the image resolution is obtained only in receive by adding coherently the contribution of each scatter. The transmission with different steering angles and the coherent compounding of the back-scattered echoes is an efficient way to re-synthesize the transmit focusing, as deeply explained in [START_REF] Montaldo | Coherent plane-wave compounding for very high frame rate ultrasonography and transient elastography[END_REF]. Results (Tab. I in [START_REF] Montaldo | Coherent plane-wave compounding for very high frame rate ultrasonography and transient elastography[END_REF]) had demonstrated that a compounded image, obtained by 12 plane waves at different sector angles and with 128-element active aperture, has a resolution in a wide ROI of about 1 mm, so it is comparable to that of a standard mono-focus B-mode image in the focal point, with the advantage that the compounded image is achieved by a ten times higher frame-rate.

The Synthetic Aperture (SA) method [START_REF] Misaridis | Use of modulated excitation signals in medical ultrasound. part III: High frame rate imaging[END_REF] presents strong conceptual analogies with the coherent plane-wave compounding. In the SA modality, a small number of transducer elements are activated to transmit a broad beam, and an optimal focusing on receive is obtained by receiving from all the elements. Each individual image has a relatively low resolution but is acquired in a very short time, so that more images can be combined to generate high-resolution images.

Both the SA and ultrafast methods can improve the focusing power and the image contrast shortening the data-acquisition period, but have the drawback that the processing for receive beamforming is more computationally demanding [START_REF] So | Medical ultrasound imaging: To GPU or not to GPU?[END_REF]. Nonetheless, their heavy computational load and high parallelism are the suitable conditions for the implementation on a Graphics Processing Unit (GPU) board. In particular, the ultrafast imaging has been already employed in many applications [START_REF] Tanter | Ultrafast imaging in biomedical ultrasound[END_REF], and its basic principles have been the rich soil for the development the algorithms for blood flow imaging presented in Chapters 4 (p. 63) and 5 (p. 89). The main prospects of such high powerful GPU boards in US imaging are delineated in Chapter 6 (p. 117).

Imaging System Resolutions

The spatial resolution of an imaging system is in general defined as the smallest separation of two point sources necessary for the sources to be resolved [START_REF] Webb | Introduction to Biomedical Imaging[END_REF]. In the longitudinal direction, the axial resolution is the minimum distance between two targets that belongs to the plane perpendicular to the focal plane and parallel to the transducer plane. The axial resolution is related to the sample volume, which is the duration δτ of the transmitted Medical Ultrasound Imaging pulse:

δτ = N c λ = cN c f 0 . (1.13)
The axial resolution is directly proportional to the frequency in transmission f 0 and inversely proportional to the number N c of transmitted cycles. The equation (1.13) compels the following tradeoff [START_REF] Coppini | Bioimmagini[END_REF]: higher f 0 and lower N c would give a better resolution, but the more is the f 0 , the more is the attenuation of the ultrasound wave (p. 11); similarly, with lower N c , the transmitted energy is lower and the signal-to-noise ratio (S/N) drops down. In the azimuth direction, the lateral resolution is determined by the minimum distance of the two targets placed on planes orthogonal to the transducer one. It depends on the transmission frequency, on the transducer geometry (e.g., the pitch), and especially on the focalization scheme. For example, in focalized imaging, the lateral resolution is optimal in the focal point but rapidly decreases in the adjacent zones.

Blood Flow Imaging

Ultrasounds have been used for medical imaging since the 1950's [START_REF] Satomura | Ultrasonic Doppler method for the inspection of cardiac functions[END_REF]. Its successfully diffusion is due to the safety for the patient (no ionizing radiations are used), joined to the high diagnostic capabilities and the relatively low cost of the systems [START_REF] Szabo | Diagnostic Ultrasound Imaging: Inside Out[END_REF]. Last but not least, most techniques are implemented in smart real-time applications that give high benefit for the clinical employment [START_REF] Moore | Point-of-care ultrasonography[END_REF].

The medical techniques of Figs. 1.3(a) and (c) have been achieved by combining both the anatomical and functional imaging methods, in order to give the clinician an exhaustive characterization of the anatomy and the physiology of the districts under investigation. For example, Grant et al. [START_REF] Grant | Carotid artery stenosis: Gray-scale and Doppler US diagnosis[END_REF] recommend that all internal carotid artery (ICA) examinations should be performed with gray-scale B-mode, spectral Doppler US and color Doppler (CFI), the so-called triplex mode. Nowadays such requirements are more than fulfilled with modern systems, but severe tradeoffs on imaging mode quality and/or frame rate are necessary, because such methods require complex sequencing, high-energy ultrasound transmission, and high processing power [START_REF] Tanter | Ultrafast imaging in biomedical ultrasound[END_REF]. The combination of the continual improvements in electronics, a better understanding of the interaction of ultrasound with tissues and more sophisticated signal processing algorithms, will lead to imaging systems of increased complexity and diagnostic capabilities [START_REF] Szabo | Diagnostic Ultrasound Imaging: Inside Out[END_REF]. In the following, the main techniques for conventional blood flow imaging are introduced, and the most recent technical solutions are then listed.

Doppler Imaging

In more than 60 years, ultrasound systems and techniques have reached so high performances that nowadays they are widely used for the diagnosis of cardiovascular diseases [START_REF] Evans | Doppler Ultrasound: Physics, instrumentation and signal processing[END_REF]. The analysis of the blood flow patterns can demonstrate or even exclude the presence of such pathologies (e.g., stenosis, blood clots and heart valve defects), enabling the clinician to grade the severity of the disease and to choose the most appropriate treatment [START_REF] Thrush | Peripheral Vascular Ultrasound: How, Why, And When[END_REF]. Blood flow measurements are typically achieved by spectral ultrasound Doppler methods, focused to detect blood flow velocity in a single sample volume or along an entire line [START_REF] Tortoli | An FFT-based flow profiler for high-resolution in vivo investigations[END_REF], and by colour flow imaging, which is capable to investigate bi-dimensional regions of interests (ROIs), with frame-rate in tens of Hz [START_REF] Kasai | Real-time twodimensional blood flow imaging using an autocorrelation technique[END_REF]. Doppler techniques take name from their discoverer, the Austrian physicist and mathematician Christian Johann Doppler, who tracked the lines of the Doppler effect theory in 1842. In pulsed-wave ultrasound, a target that moves along an ultrasound field with a velocity v produces a modulation of the transmitted frequency (f 0 ) which is revealable in receive. That frequency is called Doppler frequency f d :

f d = 2f 0 v c cos ϑ, (1.14)
in which ϑ is the beam-to-flow angle, i.e. the angle between the ultrasound beam and the direction of the target velocity. To track the blood movement, the scattering properties of blood cells are exploited. In spectral investigations, bursts are periodically transmitted inside a volume located in the vessel (i.e., the sample volume), and the so-called pulsed-wave Doppler technique can be described by a simplified model [START_REF] Jensen | Estimation of Blood Velocities Using Ultrasound: A Signal Processing Approach[END_REF], that represents the back-scattered echo signal s i by the following equation:

s i = A i sin 2π f 0 2v cos ϑ c iT P RF + ϕ 0 (1.15)
where i refers to the i th pulse of the periodical sequence repeated every T P RF time, A i is the signal amplitude and ϕ 0 a phase shift due to the wave propagation. Note that the term in brackets of (1.15) is the Doppler frequency in (1.14). This model shows how the phase shift is used in pulsedwave Doppler systems to detect the blood velocity. If the Shannon theorem is not respected (i.e., f d > P RF/2), the phase wraps around and the aliasing phenomenon occurs. In other words, the maximum velocity detectable by the Doppler system without ambiguities is

v max = c • P RF 4f 0 cos ϑ . (1.16)
The Doppler spectrum is affected also by the broadening effects, due to the pulse transit-time and to the different velocities of the groups of scatterers that are in movement inside the sample volume [START_REF] Newhouse | Analysis of transit time effects on Doppler flow measurement[END_REF]. Clinical applications of spectral Doppler methods usually investigate one sample volume by a single line of view, like in Fig. 1.3(a), and display the spectrogram which shows the evolution of the blood velocity at a single depth. Extended applications (i.e., the Multigate Spectral Doppler technique [START_REF] Tortoli | An FFT-based flow profiler for high-resolution in vivo investigations[END_REF]) provide the entire spectral profile inside the vessel lumen, detecting the spectral Doppler signal of several gates. The spectral processing of the Doppler signals and its extensions are widely used in this work (see §2.1.2 and §3.2). Alternative approaches are the autocorrelation and the cross-correlation estimators [START_REF] Kasai | Real-time twodimensional blood flow imaging using an autocorrelation technique[END_REF]. The former estimates the mean phase shift between pulse emissions in the time domain; the latter tracks the position of scatterers between consecutive emissions. For a full description of these techniques, the reader is referred to [START_REF] Evans | Doppler Ultrasound: Physics, instrumentation and signal processing[END_REF].

Color Flow Imaging

The Color Flow Imaging (CFI) exploits the autocorrelation estimator [START_REF] Kasai | Real-time twodimensional blood flow imaging using an autocorrelation technique[END_REF] to represent the blood movement in extended ROIs by superimposing a color-coded velocity map onto a B-mode image. Such modality has greatly improved cardiovascular examinations: the velocity color-map makes vessel detection easier and gives a more global information on blood movement with respect to spectral stand-alone investigations [START_REF] Evans | Doppler Ultrasound: Physics, instrumentation and signal processing[END_REF]. However, such improvements hide the drawback that CFI has an extremely low frame-rate, in the order of tens of Hz. This is due to the fact that blood cells have from 10 to 100 times lower scattering coefficients with respect of the signal from vessel boundaries (i.e. the clutter ), and the related signals are 40-60 dB lower than the clutter [START_REF] Jensen | Stationary echo canceling in velocity estimation by timedomain cross-correlation[END_REF]. Its rejection is demanded to highpass filters, but several samples (8 ÷ 16) are necessary to increase the S/N and to obtain one velocity estimation, with the inevitably decrease of the temporal resolution.

Power Doppler. An alternative approach to CFI is called Power Doppler [START_REF] Rubin | Power Doppler US: a potentially useful alternative to mean frequency-based color Doppler US[END_REF]. It is the measurement of the total Doppler power from each sample volume which is displayed as a range of colors [START_REF] Evans | Colour flow and motion imaging[END_REF]. The power of the Doppler signal is governed in first approximation by the volume of moving blood rather than its velocity [START_REF] Shung | The effects of hematocrit, shear rate, and turbulence on ultrasonic Doppler spectrum from blood[END_REF], because the phase information is removed. In the power mode, the angular cosine effects are reduced and the aliasing at high velocities has little effects on the displayed power [START_REF] Szabo | Diagnostic Ultrasound Imaging: Inside Out[END_REF]. The advantages are the increment of the sensitivity, which is exploitable for the detection of small blood vessels, and a clearer visualization of complex vascular anatomy. The drawback is the total suppression of the velocity information, which is crucial in many circumstances.

Novel Methods for Blood Flow Imaging

All the conventional Doppler techniques have three main limitations: 1) the estimation of the only axial component of the velocity, in one single sample volume (spectral Doppler) or in an extended ROI (CFI); 2) the limited velocity span due to aliasing; 3) the low frame-rate of the B-mode imaging, which is much lower in CFI investigations, due to the time-consuming line-by-line acquisition.

In this Section, after a brief introduction on limitations of the conventional approaches, some of the imaging techniques that tempt to overcome these limitations are introduced: the multi-beam approaches, the speckletracking, the transverse oscillations, the synthetic aperture, the plane-wave, and the ultrafast Doppler.

Doppler methods detect only the axial blood velocity component (i.e. the one in the US propagation direction) in a sample volume. This is not ruinous for simple vessels, because a blood velocity parallel to the boundaries can be assumed, and the direction can be manually assessed. For the exhaustive investigations of regions that include curves, stenosis or bifurcations, where the blood flow assumes complex patterns, the automatic assessment of the modulus and the angle (i.e. the vector ) of blood velocity is essential. Moreover, if the blood velocity was assessed in a 2-D or 3-D ROI (i.e. not only restricted to a sample volume), a larger portion of the vessel would be investigated, and more information would be available for the description of the flow in its whole.

Different techniques have been proposed to extend the measurement in more dimensions. The popular multi-beam approach [START_REF] Fox | Multiple crossed-beam ultrasound Doppler velocimetry[END_REF], [START_REF] Dunmire | Cross-beam vector Doppler ultrasound for angleindependent velocity measurements[END_REF] trigonometrically combines two or more Doppler measurements obtained from different transducers, to estimate the two components (2-D) of the velocity vector in the scan plane, or even all three spatial components (3-D) of the vector [START_REF] Bonnefous | Measurement of the complete (3D) velocity vector of blood flows[END_REF]. It has also been proposed to exploit the spectral broadening effect to get the direct information about the magnitude or the direction of the transverse velocity component in a sample volume [START_REF] Newhouse | Three-dimensional vector flow estimation using two transducers and spectral width[END_REF], [START_REF] Tortoli | Transverse Doppler spectral analysis for a correct interpretation of flow sonograms[END_REF] and, more recently, for real-time 2-D vector Doppler investigations [START_REF] Tortoli | An automatic angle tracking procedure for feasible vector Doppler blood velocity measurements[END_REF].

Other studies have been dedicated to the investigation of a two dimensional ROI, to produce corresponding vector flow maps [START_REF] Marion | Blood flow evaluation in high-frequency, 40 MHz imaging: A comparative study of four vector velocity estimation methods[END_REF]. The road has been paved by the speckle tracking methods [START_REF] Trahey | Angle independent ultrasonic detection of blood flow[END_REF], [START_REF] Bohs | Speckle tracking for multi-dimensional flow estimation[END_REF], which follow the scatterers movement by tracking the corresponding speckle over a given region. Such speckle patterns, or "signatures", remain relatively constant as the blood or tissue moves, and can thereby be tracked using pattern matching techniques. In other words, it determines the local displacement vectors of speckle patterns in successive digitalized ultrasound images (RF or demodulated) using a template matching algorithm. The location of the best match defines the vector of motion for the particular kernel. Bigger kernels increment the S/N but reduce accuracy and also spatial velocity resolution decrease. Given the time between acquisition of the kernel and search regions, the 2-D velocity vector is calculated. When this is done over many small kernels throughout the region of interest, a vector velocity map is generated. The displacement is typically achieved through the estimation of the 2-D cross-correlation peak, or the minimization of suitable cost functions -e.g. the sum of absolute or square differences -in the time-domain.

More recently, Jensen et al. [START_REF] Jensen | A new method for estimation of velocity vectors[END_REF] have proposed to include a special dynamic apodization function in the receive beamforming to create a transverse oscillation in the pulse-echo field. This oscillation is exploited to calculate 2-D or even 3-D [START_REF] Pihl | Measuring 3D velocity vectors using the transverse oscillation method[END_REF] vector flow maps by an autocorrelation estimator. The 2-D method was off-line validated on healthy volunteers [START_REF] Pedersen | Comparison of real-time in vivo spectral and vector velocity estimation[END_REF] and real-time implemented on a commercial ultrasound system (Vector Flow Imaging, BK Medical ApS, Denmark).

2-D vector flow maps have also been obtained using a synthetic aperture approach [START_REF] Jensen | Estimation of velocity vectors in synthetic aperture ultrasound imaging[END_REF]. In that technique, the signals received from multiple spherical-wave transmissions are beam-formed according to a double scheme that gives the angle and the magnitude of the velocity vector by the maximum of the cross-correlation function. To increase the achievable frame-rate, plane-waves (PWs) methods have been introduced by the Jensen's research group.

Udesen et al. [START_REF] Udesen | High frame-rate blood vector velocity imaging using plane waves: Simulations and preliminary experiments[END_REF] have recently demonstrated the feasibility of obtaining 2-D vector maps at high frame-rate, by plane-wave excitations coupled with 2-D cross-correlation based displacement estimator. Impressive in-vivo results [START_REF] Hansen | Fast blood vector velocity imaging using ultrasound: In-vivo examples of complex blood flow in the vascular system[END_REF] have been achieved using this method and the RAS-MUS research platform, however the 2-D cross-correlation estimator is too time-consuming, and keeps far from the real-time implementation. More recently, Ricci et al. [START_REF] Ricci | Real-time vector velocity assessment through multigate Doppler and plane waves[END_REF] have implemented dual-beam vector Doppler imaging based on the transmission of two plane waves with opposite steering angles. Ekroll et al. [START_REF] Ekroll | Simultaneous quantification of flow and tissue velocities based on multi-angle plane wave imaging[END_REF] uses the same plane-wave acquisition scheme of [START_REF] Ricci | Real-time vector velocity assessment through multigate Doppler and plane waves[END_REF] to reconstruct 2D-vector maps of wide regions, for both the arterial flow and tissue velocities.

Different approaches have been pursued by Bercoff et al. [START_REF] Bercoff | Ultrafast compound Doppler imaging: providing full blood flow characterization[END_REF]. Their so-called ultrafast compound Doppler imaging uses the ultrafast imaging method [START_REF] Montaldo | Coherent plane-wave compounding for very high frame rate ultrasonography and transient elastography[END_REF], previously described [ §1. 1.4], to implement a faster CFI mode. This modality solves the frame-rate limitations. All blood flow dynamics can be acquired during a limited number of cardiac cycles, so that the 1.3 ULA-OP exam duration is reduced and the operator can retrospectively select pulsed Doppler windows to obtain a more refined analysis of complex blood flows. This technique has the drawback of not considering the angle correction and has been recently implemented on a commercial scanner (Aixplorer, Supersonic Imagine, Aix en Provence, France).

In two following Chapters, two novel methods based on plane-wave imaging are presented. Chapter 4 draws on an estimator for the 2-D vector assessment of the blood flow [START_REF] Lenge | High-frame-rate 2-D vector blood flow imaging in the frequency domain[END_REF], that operates in the frequency domain to reduce the computational load. A GPU-implementation of this estimator, aiming at the further reduction of the computational time, has been presented in Chapter 6. In Chapter 5, an efficient approach based on planewave imaging and transverse oscillations is proposed [START_REF] Lenge | A new method for 2D-vector blood flow imaging based on unconventional beamforming techniques[END_REF], [START_REF] Lenge | Plane-wave transverse-oscillation for high frame rate 2-D vector blood flow imaging[END_REF]. Both these solutions aim at overcoming the main limitations of conventional Doppler methods, estimating the 2-D vector in wide ROIs and at high frame-rate, with an high accuracy, demonstrated by results of simulations, in-vitro, and in-vivo experiments. In the following sections, the hardware and the software materials used for these works are outlined, and the most significant clinical application of these methods (the diagnosis of the carotid artery disease) is introduced.

ULA-OP

The ULtrasound Advanced-Open Platform (ULA-OP1 , [START_REF] Boni | A reconfigurable and programmable FPGAbased system for nonstandard ultrasound methods[END_REF][START_REF] Tortoli | ULA-OP: an advanced open platform for ultrasound research[END_REF]) is a compact and complete ultrasound system designed at University of Florence for research purposes. It is connected to a PC through USB 2.0. A specific software running in the PC is used to set the system parameters and display in real-time the results of the processing. ULA-OP consists in a rack where a backplane connects 2 main boards: an analog board (AB) and a digital board (DB). An expansion connector is also available on the backplane to connect optional boards like the Acquisition Board for Ultrasound Research (i.e., ABUR). The AB includes the RF front-end (transmitters, multiplexers and low-noise receivers) while the DB hosts the digital devices in charge of numerical beamforming and signal/image processing.

ULA-OP integrates a fully programmable transmission (TX) section including 64 independent arbitrary wave generators [START_REF] Ricci | Multichannel FPGA-based arbitrary waveform generator for medical ultrasound[END_REF], which can be routed to 192-element probes. The receiving (RX) section (Fig. 1.6) acquires the signals from 64 elements of the probe selected through the switch matrix. These signals are amplified and converted to digital at 50 MHz with 12-bit format. Five Field Programmable Gate Arrays (FPGAs) performs the beamforming and the coherent demodulation in real-time. The baseband signal is then processed by software modules running on an embedded digital signal processor (DSP), to enable custom and/or standard (e.g. B-, Doppler-, CFI-) modes. The TX and RX sections can be fully reprogrammed before the start of each T P RF . Several computations can concurrently run in the DSP. The DSP formats the results in graphical frames that are sent to the PC through the USB channel for being displayed in the corresponding windows. ULA-OP includes 1.25 GB of DDR (i.e. Double Data Rate) memory where data can be stored from every point of the RX chain, to be subsequently downloaded to a PC. The 256 MB connected to DSP are typically dedicated to the acquisition of beamformed data or demodulated IQ data. Such an amount of memory is sufficient to store several seconds of RX data. On the other hand, the 1 GB of memory connected to the FPGAs can be dedicated to the acquisition of pre-beamforming data, but, for a typical acquisition, less than 1s of data can be stored.

The Acquisition Board for Ultrasound Research. The architecture of the Acquisition Board for Ultrasound Research (ABUR) board is presented in [START_REF] Boni | Multi-channel raw-data acquisition for ultrasound research[END_REF]. The 64 RX analog signals coming from the connector are sampled at 50 MHz by 12-bit analog-to-digital (AD) converters, like in the ULA-OP system. Each AD device has 8 channels and output the results in LVDS bit-streams that are directly sent to 2 high-density FPGAs (Arria GX, Altera). Here, the signals are de-serialized, packed and stored, at a data rate of 57.6 Gbit/s, in 4 dual in-line commercial DDR2 memory modules. The capacity of each module is made equal to 9 GB by suitably exploiting the 9-bit format (8 + 1 error correction bit) of the memory cells. A third FPGA controls the data flow according to the user requests. Both the acquisition length within each T P RF and the total number of PRIs to 1.4 SIMAG 25 be acquired can be programmed. On user request the board starts data sampling and storage, triggering the acquisition task, which saves signal samples to the memories according to the programmed settings.

Similarly to the ULA-OP system, the ABUR connects to a PC through a dedicated USB 3.0 channel. A dedicated PC software has been created to manage the acquisition board without interfering with the concurring ULA-OP real-time application. The software presents an user interface that allows to set the relevant acquisition parameters. In particular, the starting depth, the number of samples for each T P RF and the total number of PRIs to be saved, can be set by the user. Once the samples are acquired in the board internal memory, the software allows the user to choose the size and offset of the data that has to be downloaded and saved on a file on the host PC. Since downloading 36 GB of data can be really time-consuming, the software and the FPGA firmware have been duly optimized for maximum performance. A transfer rate of 150 MB/s was measured when the ABUR is connected to a recent PC.

The ABUR connected to ULA-OP. The acquisition board is hosted on the expansion connector of the ULA-OP backplane. The ULA-OP AB board, that interfaces to the probe, produces the 64 amplified analog signals. These signals are routed to the backplane and distributed in parallel to the ULA-OP DB, which performs the real-time tasks, and the ABUR, which saves the data as described in the previous section. Through the backplane, the ABUR has access to the system clock and the synchronization signals as well.

The phantom and the in-vivo experiments of the novel ultrasound imaging techniques presented in the following Chapters have been achieved by ULA-OP coupled with ABUR.

SIMAG

SIMAG is an ultrasound system simulator conceived for research and didactic purposes [START_REF] Ramalli | Development of Novel Ultrasound Techniques for Imaging and Elastography : From Simulation to Real-Time Implementation[END_REF]. It is based on Field II ® software [START_REF] Jensen | FIELD: A program for simulating ultrasound systems[END_REF][START_REF] Jensen | Calculation of pressure fields from arbitrarily shaped, apodized, and excited ultrasound transducers[END_REF], the ultra- sound simulator worldwide used for the design of ultrasound systems. Field II ® exploits the theory of linear systems: when the transducer is excited by a Dirac delta function, the ultrasound field in a point is got by the spatial impulse response. The images of the transmitted beam (one-way field) and received beam (two-way field) are obtained, respectively, by the convolution between the excitation signal and the impulse response in transmission and receive. In SIMAG, the user selects the scanning modality in which the transducer, the emitted pulse and beamforming schemes in transmission and receive are settable by an user-friendly interface (e.g., see Fig. 1.7). The standard tools for the development of ultrasound signal processing algorithms are implemented, including signal filtering, blood flow phantoms, and the classical B-mode representations. The major strength of SIMAG is that it has been conceived for MATLAB ® (The MathWorks Inc., Natick, MA) and based on an object-oriented paradigm. This feature facilitates the development of new signal processing algorithms for ultrasound applications. Simulations of the blood flow imaging algorithms presented in Chapters 3, 4 and 5 have been achieved by SIMAG.

Medical Imaging for Carotid Artery Disease Diagnosis

Stroke (i.e., in Latin ictus) is one of the most important cause of death in the world, with the cardiovascular diseases and the neoplasms. Between

Medical Imaging for Carotid Artery Disease Diagnosis

27

10% and 12% of deaths for a year are due to stroke, which represents also the major cause of disability and dementia [START_REF] Organisation | Global atlas on cardiovascular disease prevention and control[END_REF]. It is defined as the loss of brain function due to a disturbance in the blood supply to the brain.

The cardiovascular system consists of the heart and the branched arrangement of blood vessels that transports blood inside the human organism, including the brain. The blood vessels are classified into the arteries, that carry oxygenated blood from the heart to tissues, the capillaries, responsible of the tissue perfusion, and the veins that bring back the blood to the central system (lungs and heart). The vessels anatomy reflects its functionality [START_REF]Anatomia e fisiologia[END_REF]. The walls of huge vessels are organized in stratified layers. More externally, the adventitia is formed by fibrous tissues that maintain the vessel structure. The media is organized in two laminar elastic tissues that give elastic properties to the arteries. More internally, the intima is directly in contact with blood and is covered by endothelial tissue. One of the most important vessels is the carotid artery, which task is bringing the blood from the aortic arch to the brain and the face. The carotid artery follows the neck (common carotid artery, CCA) and bifurcates in the internal carotid artery (ICA), that supplies the brain, and the external carotid artery, which nourishes other portions of the head such as face, scalp, skull, and meninges (Fig. 1.8).

The mathematical model of the blood flow in a vessel is governed by fluid dynamics. According to the Poiseuille's law [START_REF] Waite | Biofluid Mechanics in Cardiovascular Systems[END_REF], the theoretical parabolic laminar flow in a cylindrical pipe can be described by the formula:

v(r) = v 0 1 - r R 2 , ( 1.17)
where R is the pipe inner radius and v 0 is the velocity peak in the centre of the tube. Note that, according to ϑ, the 2-D theoretical velocity profile can be expressed in a Cartesian coordinate system by its axial (v z ) and lateral (v x ) components. To obtain the total flow passing through a cross-section of the tube, the differential flow is integrated between r = 0 and r = R. This model is a simplified description of the flow dynamics inside the vessel and is further used in this work. More complex blood simulators have been developed in [START_REF] Ekroll | Ultrasound imaging of blood flow based on high frame rate acquisition and adaptive signal processing[END_REF].

The sclerosis of the arterial walls (i.e., the atherosclerosis) is mainly caused by extra-circulatory factors, such as old age, high blood pressure, previous stroke or transient ischemic attack (TIA), diabetes, high cholesterol, tobacco smoking and atrial fibrillation. Stroke is frequently due to the rupture of plaques formed inside the intima lesions by the invasion by macrophage foam cells. Not all plaques are vulnerable, sometimes they become calcified and stiffen the arteries. The development, the erosion and the rupture of the plaques are extremely linked to the blood flow patterns, because turbulence created by curves, stenosis and bifurcations produce frictional forces on the wall surfaces. For that reason, the majority of carotid artery disease develops at the carotid bifurcation [START_REF] Thrush | Peripheral Vascular Ultrasound: How, Why, And When[END_REF], and in the most cases they are asymptomatic before the obstruction of the 70 ÷ 80% of the artery lumen.

Computed tomography (CT) angiography was the gold standard to investigate stenosis in the carotid arteries [START_REF] Ekroll | Ultrasound imaging of blood flow based on high frame rate acquisition and adaptive signal processing[END_REF]. However, because of the invasiveness of its intrinsic dangerous procedures, that consists of the insertion in the artery of a contrast agent by a catheter, increasing de facto the risk of a stroke event, CT angiography has been gradually replaced by other noninvasive techniques with competitive diagnostic accuracies [START_REF] Edelman | The history of MR imaging as seen through the pages of radiology[END_REF].

Magnetic resonance angiography (MRA) is largely used for blood flow imaging of many arteries in the body, including carotid arteries and in particular cerebral vessels, for which it is considered as the gold standard [START_REF] Hansen | In vivo validation of a blood vector velocity estimator with MR angiography[END_REF]. Many techniques have been developed to exploit MR images of vessels, which provide detailed image of the arterial walls and plaques instead of only the lumen, as provided by computed tomography (CT) [START_REF] Ekroll | Ultrasound imaging of blood flow based on high frame rate acquisition and adaptive signal processing[END_REF]. Nevertheless, the most effective methods need a contrast agent injection to improve enough the S/N. Moreover, MRA requires high resources (MR 1.5 Medical Imaging for Carotid Artery Disease Diagnosis 29 scanner, reserved areas, practiced clinical staff) and complex procedures to be executed. If compared to ultrasound techniques, MRA blood flow images has lower temporal resolution, because several heart cycles have to be averaged. These observations explain why Doppler measurement are so widespread for preventive investigations in medical check-ups, while the MRA is considered as 2 nd level in-depth analysis. Nowadays, Doppler ultrasound has been adopted as the more advisable tool for the preventive screening and diagnosis of stenotic arteries. Quantitative evaluations of the degree of stenosis are obtained measuring the peak systolic velocity (PSV). For example, Grant et al. [START_REF] Grant | Carotid artery stenosis: Gray-scale and Doppler US diagnosis[END_REF] diagnose a stenosis of maximum the 50% when ICA PSV is less than 125 cm/s and plaque or intimal thickening is visible. The main problem of this kind of approaches is the variability of Doppler measurements. Aside from the angle correction, which is even more critical for stenotic arteries, different equipments, procedures and operators introduce high variability in the velocity estimations [START_REF] Ekroll | Ultrasound imaging of blood flow based on high frame rate acquisition and adaptive signal processing[END_REF]. Such limitations strongly encourage the development of new automatic non-invasive imaging techniques. developed. In this Chapter, the signal processing platform and two applications for volume flow measurement [START_REF] Ricci | Volume flow assessment through simultaneous B-mode and multigate Doppler[END_REF] 1 and the flow-mediated dilation [START_REF] Ramalli | An integrated system for the evaluation of flow mediated dilation[END_REF] 2 are presented.

Chapter 2

Signal Processing Platform for Ultrasound Imaging

Summary. The development of novel signal processing techniques for blood flow imaging is structured in three consecutive steps: first the idea, second the algorithm and its implementation, and third the validation by in-silico, in-vitro or in-vivo experiments. In order to facilitate the development of innovative algorithms for US imaging and flow imaging, a group of modules for the off-line signal processing of echographic signals has been

Signal Processing Platform

The signal processing platform (SPP) has been organized as a collection of functional blocks implemented as MATLAB ® classes. Each block is responsible of a specific task and works in cascade with the other blocks, i.e. the output of an upstream block forms the input for a downstream block. Fig. 2.1 shows the flow chart and the links among the main functional blocks.

Data

The Data class works as an interface between SPP and ULA-OP [ §1.3], converting the data stored in a specific file into a MATLAB ® two dimen- sional matrix. Since the acquired file could be very large and not fit the computer memory, the class is implemented to read a reduced set of data each time it is requested. The last read data, i.e. the Data class output, is the input of both B-mode and Profile classes.

B-mode Computing

The B-mode class organizes the data in lines and frames, interpolates them along either axial and lateral direction, and filters them in time and in 2-D spatial domain. Furthermore, the B-mode class applies a logarithmic scale compression and converts the results in a color scale. In a few words, the class implements the B-mode processing [ §1. 1.4]. The B-mode result is further processed by the Diameter class which extracts the vessel diameter and the position of both near and far walls [START_REF] Gemignani | A system for real-time measurement of the brachial artery diameter in b-mode ultrasound images[END_REF]. The wall-tracking algorithm here employed locates the points of the vessel walls, exploiting a detector based on the First Order Absolute Moment (FOAM) statistical operator. The detector performs a recursive calculation of the mass center point of the gray level variability, that occurs in a predefined circular area which must include the wall. This procedure quickly converges towards the edge following the minimum distance direction. The initial circular area is selected by the operator in the first B-mode frame, then the algorithm automatically tracks the edges in the subsequent frames. The FOAM algorithm was selected because of its good capability of tracking the inner region of the lumen.

Profile Computing

The Profile class processes the last read data similarly to the ULA-OP real-time software [START_REF] Tortoli | ULA-OP: an advanced open platform for ultrasound research[END_REF], previously described in §1.2.1. Multigate Spectral Doppler (MSD) profiles, containing the information about the distribution The WSR (wall shear rate) block has two inputs, i.e. the outputs coming from Speed and Diameter. In particular, a polynomial least square fit is applied on the velocity points and the resulting velocity profile is used to evaluate the speed gradient with respect to the radius (shear rate). The positions of the walls estimated by the Diameter block are used to extract the corresponding WSR.

The Flow block exploits the Speed and Diameter outputs to assess the volume flow. If v z (r, t) is the axial velocity profile measured along the vessel diameter at instant t and depth r, and -R(t) and R(t) are the wall positions tracked by the FOAM algorithm at instant t, the instantaneous volume flow Q(t) is obtained as:

Q(t) = π +R(t) -R(t) v(r, t)|r|dr. (2.1)
which represents the integral of the velocity profile v(r, t) around the symmetry axis located at the diameter center (r = 0) [ 1.5].

Volume Flow Assessment

Introduction

The typical approach for blood volume flow assessment consists in three steps: measuring the blood velocity in a single sample volume (SV) placed in the middle of the vessel, estimating the average vessel diameter, and applying a flow model (e.g. the quasi-static Poiseuille parabolic flow distribution [START_REF] Doucette | Validation of a doppler guide wire for intravascular measurement of coronary artery flow velocity[END_REF], in (1.17)) for inferring the final estimate. This simple approach is affected by poor accuracy and repeatability for two main reasons: it cannot account for the complex flow configurations usually present in the artery and, moreover, it ignores the wall distension which, in important regions like the common carotid artery (CCA), could easily exceed the 10% of the diameter.

Better results can be achieved by exploiting more sophisticated models. For example, the Womersley model [START_REF] Womersley | An elastic tube theory of pulse transmission and oscillatory flow in mammalian arteries[END_REF] accounts for blood pulsatility, but its application is limited to straight tubes with rigid walls [START_REF] Leguy | Model-based assessment of dynamic arterial blood volume flow from ultrasound measurements[END_REF], thus neglecting the flow distortions due to curvatures. MSD methods allow the measurement of the actual velocity profile across the artery with high temporal resolution [START_REF] Krams | Effect of vessel curvature on Doppler derived velocity profiles and fluid flow[END_REF]. This technique can be exploited to provide an accuracy improvement that is larger when the real profile differs from the ideal shapes predicted by the aforementioned mathematical models.

This Section draws on the work [START_REF] Ricci | Volume flow assessment through simultaneous B-mode and multigate Doppler[END_REF] 3 and concerns the implementation of a volume flow measurement method based on the simultaneous assessment of the velocity profile, obtained from MSD data, and of the wall positions, obtained from B-mode data. The method represents an evolution of [START_REF] Ricci | Volumetric blood flow assessment through multigate spectral Doppler[END_REF] which was based on a dual-beam custom board [START_REF] Ricci | A programmable real-time system for development and test of new ultrasound investigation methods[END_REF] with no imaging capabilities. The implemented methods exploit the SPP with an ad hoc graphical user interface (Fig. 2.2). The proposed method was tested on a flow phantom using the ULA-OP research system [START_REF] Boni | A reconfigurable and programmable FPGAbased system for nonstandard ultrasound methods[END_REF], and finally validated by Ricci et al. in [START_REF] Ricci | Accuracy and reproducibility of a novel dynamic volume flow measurement method[END_REF].

Methods

ULA-OP and SPP.

The ULA-OP [ §1.3] was programmed to transmit alternate B-mode and pulsed-wave bursts through a linear array probe placed in longitudinal position over the vessel. The velocity distribution was obtained by applying the Speed chain (Data, Profile, Denoising, Frequency Extraction, Speed ), that processed the slow-time data gathered from 512 sample volumes aligned along the pulsed-wave line, covering the whole vessel at a suitable Doppler angle. Simultaneously, the wall positions were tracked by the Diameter chain (Data, B-mode, Diameter ), that exploited the gray level gradient-based algorithm [START_REF] Gemignani | A system for real-time measurement of the brachial artery diameter in b-mode ultrasound images[END_REF] of each B-mode frame to locate the points where the pulsed-wave line crosses the walls (P 1 and P 2 in Fig. 2.2). Finally, the Flow estimated the volume flow, assuming a cylindrical flow distribution.

Phantom set-up. The tests were carried out at Lund University (Sweden) by an experimental set-up. A programmable progressive cavity pump [START_REF] Eriksson | A computercontrolled arbitrary flow wave form generator for physiological studies[END_REF], controlled by a PC, forced a blood mimicking fluid to run from a reservoir into a tube immersed in a water tank. Suitable fittings in the measuring site allowed assembling pipes of different materials and dimensions. In the reported experiments, rigid plastic and flexible silicon pipes of 8 mm diameter were used. The probe was held at the top of the measuring site and was aligned in longitudinal position parallel to the pipe axis.

The blood mimicking fluid was prepared by dissolving 3.4 g of Orgasol ® (Arkema Inc. Philadelphia, PA) in 2 l of demineralised water [START_REF] Ramnarine | Validation of a new blood-mimicking fluid for use in Doppler flow test objects[END_REF]. The Orgasol ® is constituted by 10 µm polyamide spherical particles, capable of producing an US back-scattering waves similar to that generated by human blood erythrocytes. To support an uniform particle distribution, 1.8 g of surface-active agent was also added to the fluid. Finally, to guarantee per-Signal Processing Platform for Ultrasound Imaging fect homogeneity, the solution was mixed for several hours by a magnetic stirrer before each experimental session.

Experiments

The volume-flow method was tested with more than 1700 experiments carried out in steady, sinusoidal and pulsatile flow condition inside rigid and flexible pipes. For each flow/pipe combination (Tab. 2.1), different steering angles and flow rates were used. For example, with flexible pipe and sinusoidal flow trend, the nominal volume flow was changed in the range (60 ÷ 630) ml/min in 30 ml/min step (19 different rates) and the steering angle in the range (10 ÷ 22) • in 2 • step (7 different angles). For each configuration, 3 measurements of several seconds were performed, for a total of 19 × 7 × 3 = 399 measurements.

For each measurement, about 50 MB of IQ data were saved in both B-and PW-mode and post-processed by SPP. B-mode images were recalculated on 512 × 192 points, obtaining a frame sequence at 50 Hz rate. The edge detector was then applied to the B-mode sequence, and the two points where the investigation line intersects the walls were finely located and tracked. Pulsed-wave data were processed by a 30 Hz wall filter and 128-point FFT with 50% overlap. Velocity estimates were produced at 62.5 Hz rate. Finally, after interpolating the wall position measurements to match the rate of blood velocity, the instantaneous blood volume flow was calculated.

Results

The measurements were analyzed to evaluate the accuracy and reproducibility of the volume-flow method. In particular, the accuracy was assessed comparing the measurement Q m toward the reference value Q g pro-

Experiments
Err % = Q m -Q g Q m • 100.
(2.2) Fig. 2.3 reports the histogram of the errors calculated from all measurements. A Gaussian curve is fitted to the histogram for comparison. The mean error (Err % ) and standard deviation (S d ), listed in Tab. 2.2 together with 95% confidence bounds, are -3.7% and ±5.3%, respectively. The value of the coefficient of determination R 2 of the fitting confirms that the percentage error population is well represented by the normal distribution.

In Fig. 2.4, the estimated flow is compared to the reference. The vertical axis reports the percentage error, Err % , while the horizontal axis reports the nominal flow. The central horizontal dashed line shows the mean value (-3.7%), while the upper and lower horizontal dashed lines extend for S d (±5.3%).

The reproducibility was firstly qualitatively evaluated by comparing 3 measurements obtained for each flow/pump/pipe configuration. For example, Fig. 2.5 reports the measurements obtained for a 420 ml/min nominal flow in pulsed conditions at 12 • steering angle. The acquisitions were aligned at the first peak for comparison.

The repeatability was then quantified by means of the coefficient of variability (CV ), which is obtained by the standard deviation (S d ) of the repeated measurements divided by their average value. From the 1743 mea- surements, i.e. 581 triplets, we calculated the CV s populations obtaining a 1.7% mean value and 1.7% standard deviation.

Discussion and Conclusion

The measurement of volumetric blood flow requires simultaneous knowledge of the absolute velocity distribution over a cross-sectional area of the investigated vessel as well as the estimate of the area. These pieces of information are difficult to be obtained by means of standard US investigations. The volume-flow method overcomes some of the limitations by simultaneously detecting the actual velocity profile and the wall positions. This method, tested on a very large set of controlled conditions and with more than 1700 experiments, features a slight underestimation (-3.7%), probably due to the effect of the finite size of the sample volume on pulsedwave measurements [START_REF] Gemignani | A system for real-time measurement of the brachial artery diameter in b-mode ultrasound images[END_REF]. The relative error is normal distributed (Fig. 2.3) with a 5.0% S d , which does not depend on volume flow rate (Fig. 2.4). The reproducibility is very good also in challenging conditions: the example of Fig. 2.5 was obtained in pulsed condition with high flow slew-rates, which produced a distension of the silicon pipe larger than 10%. The good reproducibility is also confirmed by the low CV calculated over the whole experimental data. The proposed method does not require particular computational power nor relies on special transmission/reception strategies, so that it can be easily integrated in current echographic systems. Moreover its accuracy could be further improved by exploiting a vector Doppler methods [START_REF] Ricci | Accuracy and reproducibility of a novel dual-beam vector Doppler method[END_REF] for velocity profile assessment, which has already been integrated in the ULA-OP system [START_REF] Tortoli | An automatic angle tracking procedure for feasible vector Doppler blood velocity measurements[END_REF].

Flow-Mediated Dilation

Introduction

Flow-mediated dilation (FMD) is widely used to assess the endothelial function by ultrasound. Blood flow in the brachial artery is restricted by a cuff for about 5 minutes. When the cuff is removed, the wall shear stress (i.e. the product of wall shear rate, WSR, by blood viscosity) increases, thus stimulating the release of a vasodilator, the nitric oxide, from the endothelial cells into the smooth muscle. By measuring the diameter change due to reactive hyperemia, the possible endothelial dysfunctions can be detected [START_REF] Celermajer | Endothelium-dependent dilation in the systemic arteries of asymptomatic subjects relates to coronary risk factors and their interaction[END_REF].

Although FMD is thought to be evoked by shear stress, no direct measurements of the stimulus (WSR) yielding the FMD have been reported until 2011 [START_REF] Tortoli | Simultaneous ultrasound assessment of brachial artery shear stimulus and flow-mediated dilation during reactive hyperemia[END_REF]. This can be mainly attributed to the fact that the measurement of wall shear rate, i.e. the blood velocity gradient near the arterial walls, requires the simultaneous estimation of velocity in multiple sample volumes located at different distances from the probe. The computational complexity of such Multigate Spectral Doppler measurement adds to the inherent complexity of automatic diameter change measurement. The latter has been traditionally obtained from the analysis of B-mode echographic images. For this purpose, gradient-based processing [START_REF] Beux | Automatic evaluation of arterial diameter variation from vascular echographic images[END_REF][START_REF] Cinthio | A method for arterial diameter change measurements using ultrasonic B-mode data[END_REF], artificial neural networks [START_REF] Newey | Online artery diameter measurement in ultrasound images using artificial neural networks[END_REF] or first order absolute central moment estimations [START_REF] Demi | The first absolute central moment in lowlevel image processing[END_REF] have been proposed. Such methods, however, have not been so far integrated together with a WSR measurement method in a single equipment dedicated to FMD.

This section draws on the work [START_REF] Ramalli | An integrated system for the evaluation of flow mediated dilation[END_REF] 4 and concerns the description of an integrated system capable of estimating either the stimulus (WSR change) and the effect (diameter change) in FMD investigations. To reach this goal, the hardware of the ULA-OP has been enriched with a specific acquisition system capable of storing the large amount of raw echo-data produced during this long examination (the ABUR, §1.3). Furthermore, a graphical user interface for SPP has been developed to analyze the data in order to evaluate the endothelial functionality (Fig. 2.6).

Methods

The ULA-OP [ §1.3] was programmed to transmit alternate B-mode and pulsed-wave bursts through a linear array probe placed in longitudinal position over the vessel. As in [START_REF] Ricci | Volume flow assessment through simultaneous B-mode and multigate Doppler[END_REF], the diameter were measured by exploiting the Diameter chain (Data, B-mode, Diameter ). At the same time, the velocity profile was obtained by the Speed chain (Data, Profile, Denoising, Frequency Extraction, Speed ), that processes by a 256-point FFTs the data gathered from 512 sample volumes. Finally, the WSR block estimates the corresponding wall shear rate.

Results and Discussion

The ULA-OP and the SPP have been exploited for the acquisitions of healthy volunteers and patients. Fig. 2.7 shows the output of the main processing blocks, as they are directly plotted by the graphic methods of the developed classes. In particular, Fig. Furthermore, in order to facilitate the use of the processing platform, a user friendly graphical interface has been developed (Fig. 2.6). This allows the user to change all the parameters, and to check the instantaneous values of both diameter and WSR over the entire exam. The user can also check whether the processing chain is correctly working by looking at the B-mode image, the MSD profile and at the speed estimation and fitting.

In Fig. 2.8 an example of a FMD processing result is reported. The upper plot shows the WSR averaged between the near and far wall WSR values, the middle plot reports the diameter changes and the lower one shows the blood pressure. The graph is splitted in two time intervals, the first corresponding to the acquisition before the cuff restriction (baseline) and the second, 5 minutes long, starting one minute before releasing the cuff. The red lines report the values averaged on each cardiac cycle. These graphs allow extracting the main parameters needed to investigate the endothelial function, i.e. the WSR peaks, the area under the curve of WSR, the percentage diameter distension and the delay occurring between WSR peak and diameter peak.

Conclusion

In this work, an integrated system for FMD application, consisting of a hardware and a software part, has been presented. The hardware of an existing scanner (ULA-OP) has been upgraded to acquire the large amount of data, which is needed for FMD exams. By SPP, the WSR is directly estimated from the velocity profile gradient, avoiding assumptions concerning the flow nature. The approach followed in the design of the SPP blocks, i.e. based on object-oriented programming and on a specific task for each block, allows changing or updating each part of the processing chain. For example, different algorithms [START_REF] Gemignani | A system for real-time measurement of the brachial artery diameter in b-mode ultrasound images[END_REF][START_REF] Cinthio | A method for arterial diameter change measurements using ultrasonic B-mode data[END_REF] for arterial diameter detection have been tested simply developing two classes differing only in one method. In the same way the flow speed can be evaluated starting from the mean spectral frequency or by extracting the maximum spectral frequency. The integrated FMD system is currently adopted in a multicenter study as part of the SUMMIT Project granted by UE.

Conclusions

In this Chapter, a platform for the processing of ultrasound data have been presented. The platform is able to process the B-mode and the Doppler data, in order to extract meaningful characteristics of the blood dynamics, such as the velocity, the flow, the vessel diameter and the wall shear-rate. Such measurements increase the knowledge about the blood flow, and are additional information to understand the etiology of the cardiovascular diseases.

Chapter 3

Clinical Comparison of Blood Velocity Measurements by Standard and Vector Doppler Approaches

Summary. Carotid artery peak systolic velocity (PSV) measurements, although severely affected by the angle dependency, are widely used for stenosis assessment. Many vector Doppler (VD) methods have been so far proposed, but their clinical usefulness was not demonstrated yet. Blood PSVs in the common and internal carotid arteries of both healthy volunteers and patients with stenosis have here been measured by two VD methods and the results compared with those produced by conventional spectral Doppler. The results of this study (included in [START_REF] Lenge | Blood velocity measurement in healthy and diseased carotid arteries by vector Doppler techniques[END_REF] 1 and [START_REF] Tortoli | Comparison of carotid artery blood velocity measurements by vector and standard Doppler approaches[END_REF] 2 ) suggest that the introduction of vector Doppler methods in commercial machines may be considered mature and capable of overcoming the angle-dependent over-estimation typical of the standard spectral Doppler approach.

Introduction

Stroke is one of the major causes of death [START_REF] Koton | Stroke incidence and mortality trends in US communities, 1987 to 2011[END_REF][START_REF]for Health Statistics, Compressed Mortality File, 1999-2009 ([END_REF] and of serious physical and cognitive long term disability in adults [START_REF] Brault | Prevalence and most common causes of disability among adults -United States 2005[END_REF]. The ischemic stroke is more frequent [START_REF] Goldstein | Council for High Blood Pressure Research, and Council on Peripheral Vascular Disease, and Interdisciplinary Council on Quality of Care and Outcomes Research[END_REF] and, in most cases, can be caused by emboli or artery thrombosis. Atherosclerotic plaques have an important etiologic role, indeed carotid artery stenosis causes approximately (10 -20)% of ischemic strokes [START_REF] Kistler | Carotid endarterectomy revisited[END_REF]. Early identification and characterization of atherosclerotic plaques is thus fundamental for risk stratification in the patient. The stenosis can be assessed by angiographic (CT/RM) and/or Doppler ultrasound methods. Magnetic resonance is usually considered more accurate than ultrasound [START_REF] Nicolaides | Ultrasound and Carotid Bifurcation Atherosclerosis[END_REF], although recent ex-vivo studies have shown similar accuracies for the two approaches [START_REF] Pan | Assessment of carotid artery stenosis by ultrasonography, conventional angiography, and magnetic resonance angiography: correlation with ex vivo measurement of plaque stenosis[END_REF]. However, the cost-effectiveness, portability, safety and repeatability of duplex scanning make it the test of choice in the vast majority of patients [START_REF] Nicolaides | Ultrasound and Carotid Bifurcation Atherosclerosis[END_REF].

Blood velocity is one of the most significant parameters to identify carotid stenosis, especially if that is over 70%. The spectral Doppler method, in particular, can be used to measure the peak systolic velocity (PSV) and end-diastolic velocity (EDV) close to the stenosed region. These parameters correlate with the stenosis rate estimated using the reference technique, i.e. the angiographic approach [START_REF] Grant | Carotid artery stenosis: Gray-scale and Doppler US diagnosis[END_REF][START_REF] Heijenbrok-Kal | Optimal peak systolic velocity threshold at duplex US for determining the need for carotid endarterectomy: A decision analytic approach[END_REF][START_REF] Von Reutern | Neurosonology Research Group of the World Federation of Neurology[END_REF].

However, the reliability of Doppler velocity measurements is severely affected by their intrinsic operator-and machine-dependency. For example, sample volume (SV) placement [START_REF] Mynard | Effect of velocity profile skewing on blood velocity and volume flow waveforms derived from maximum Doppler spectral velocity[END_REF] and gain setting are fundamental for the PSV measure, while the inter /intra operator variability can cause significant measurement differences [START_REF] Lui | Human factors as a source of error in peak Doppler velocity measurement[END_REF][START_REF] Von Reutern | Measuring the degree of internal carotid artery stenosis[END_REF]. One of the major limitations of the Doppler approach is surely represented by its dependency on the beam-to-flow angle, ϑ. Indeed, the Doppler equation (1.14) states that the detected Doppler frequency, f d , depends on the target velocity, v, as well as on cos ϑ:

v = f d 2 f 0 c cos ϑ . (3.1)
An error of few degrees in the estimate of ϑ leads to a significant error in the velocity estimate due to the nonlinear contribution, especially at large angles. For this reason, the sonographers typically try to maintain the Doppler angle as close to 60 • as possible in clinical practices [START_REF] Evans | Doppler Ultrasound: Physics, instrumentation and signal processing[END_REF][START_REF] Beach | Standardized ultrasound evaluation of carotid stenosis for clinical trials: University of washington ultrasound reading center[END_REF]. The estimation of the velocity magnitude thus needs the manual adjustment of the flow direction by the sonographer, who typically assumes it is parallel to the vessel walls. The operator aims at setting a suitable cursor parallel to the vessel wall, and the machine consequently calculates ϑ. In this way, he implicitly introduces an approximation that is effective only for rectilinear sections of big arteries, and risks to be inappropriate for stenosed arteries which are the most important to be diagnosed. Most carotids have a significant velocity component orthogonal to the probe axis and, most importantly, this method assumes that a single Doppler angle is involved, while the velocity vector subtends a full range of angles at the transducer aperture, thus producing the well known geometrical broadening effect [START_REF] Newhouse | The dependence of ultrasound Doppler bandwidth on beam geometry[END_REF][START_REF] Daigle | Overestimation of velocity and frequency values by multi element linear array Doppler[END_REF][START_REF] Thrush | Intrinsic spectral broadening: a potential cause of misdiagnosis of carotid artery disease[END_REF].

The estimated angle is typically used, like in (3.1), to convert all detected Doppler frequencies, f d , to velocities, and the vertical axis of spectrograms is thus directly indicated in velocity units. In particular, in the spectral Doppler method, the value corresponding to the spectrogram's peak is assumed to be coincident with the PSV. The negligence of geometrical broadening has here deleterious effects, since the machine assumes that the maximum Doppler frequency is received at the centre of the aperture whereas it is received at the edge of the aperture [START_REF] Newhouse | The dependence of ultrasound Doppler bandwidth on beam geometry[END_REF][START_REF] Hoskins | A comparison of single-and dual-beam methods for maximum velocity estimation[END_REF]. Although a significant PSV overestimation is thus obtained, the method is still implemented in most commercial US machines, and widely adopted in Doppler investigations [START_REF] Trial Collaborative Group | General anaesthesia versus local anaesthesia for carotid surgery (GALA): a multicentre, randomised controlled trial[END_REF].

Although many efforts have been spent by clinicians to make the best use of the traditional Doppler technology, the residual variability due to machine-and operator-reasons suggests adopting innovative methods that overcome the Doppler angle dependency. The multi-beam [START_REF] Dunmire | Cross-beam vector Doppler ultrasound for angleindependent velocity measurements[END_REF][START_REF] Tortoli | An automatic angle tracking procedure for feasible vector Doppler blood velocity measurements[END_REF], the transverse oscillations [START_REF] Jensen | A new method for estimation of velocity vectors[END_REF], speckle tracking [START_REF] Bohs | Speckle tracking for multi-dimensional flow estimation[END_REF], and plane-wave [START_REF] Bercoff | Ultrafast compound Doppler imaging: providing full blood flow characterization[END_REF][START_REF] Udesen | High frame-rate blood vector velocity imaging using plane waves: Simulations and preliminary experiments[END_REF][START_REF] Ricci | Real-time vector velocity assessment through multigate Doppler and plane waves[END_REF][START_REF] Ekroll | Simultaneous quantification of flow and tissue velocities based on multi-angle plane wave imaging[END_REF][START_REF] Lenge | High-frame-rate 2-D vector blood flow imaging in the frequency domain[END_REF][START_REF] Lenge | A new method for 2D-vector blood flow imaging based on unconventional beamforming techniques[END_REF] techniques have gained increasing interest [START_REF] Pedersen | Comparison of real-time in vivo spectral and vector velocity estimation[END_REF][START_REF] Ricci | Accuracy and reproducibility of a novel dual-beam vector Doppler method[END_REF][START_REF] Hansen | In vivo comparison of three ultrasound vector velocity techniques to MR phase contrast angiography[END_REF].

Among the dual beam approaches, the one capable of tracking the Doppler angle by using the transverse Doppler principle [START_REF] Tortoli | Transverse Doppler spectral analysis for a correct interpretation of flow sonograms[END_REF] has recently been implemented in the ULA-OP platform [START_REF] Tortoli | An automatic angle tracking procedure for feasible vector Doppler blood velocity measurements[END_REF] to allow repeatable PSV measurements in the common carotid arteries (CCAs) of 13 healthy subjects. In [START_REF] Pedersen | Comparison of real-time in vivo spectral and vector velocity estimation[END_REF], a modified commercial scanner integrating the transverse oscillation vector technique has been used to estimate multiple parameters (including the PSV and the EDV) in the CCAs of 16 healthy volunteers. In [START_REF] Ekroll | Simultaneous quantification of flow and tissue velocities based on multi-angle plane wave imaging[END_REF], a dual plane wave vector Doppler technique was tested also in patients with carotid artery disease, showing images of blood flow through stenoses with secondary flow patterns. A similar approach was used by [START_REF] Ricci | Real-time vector velocity assessment through multigate Doppler and plane waves[END_REF] to estimate, in real time, the (multigate) velocity vectors along one line having a direction coincident with the axis of the used linear array probe. Promising results from the carotid bulb of healthy volunteers were presented.

This Chapter draws on the works [START_REF] Lenge | Blood velocity measurement in healthy and diseased carotid arteries by vector Doppler techniques[END_REF] and [START_REF] Tortoli | Comparison of carotid artery blood velocity measurements by vector and standard Doppler approaches[END_REF] 3 , and concerns the PSVs measurements comparison in the CCAs and internal carotid arteries (ICAs) of both healthy volunteers and patients with different degrees of stenosis, Standard and Vector Doppler Approaches by means of three different approaches. First, the vessels were investigated by a commercial machine using the conventional spectral Doppler method. Then the ULA-OP research platform was used to perform PSV measurements, at the same locations, by means of the aforementioned angle tracking [START_REF] Tortoli | An automatic angle tracking procedure for feasible vector Doppler blood velocity measurements[END_REF] and plane wave [START_REF] Ricci | Real-time vector velocity assessment through multigate Doppler and plane waves[END_REF] vector Doppler (VD) methods. Average PSVs collected from 15 subjects were pairwise compared by the Bland-Altman plot and the linear regression analysis. Two expert sonographers were also recruited to evaluate the inter -operator and intra-operator repeatability of each Doppler method on an additional group of 8 healthy volunteers.

Materials and Methods

Method 1: Conventional Spectral Estimation

The majority of commercial scanners employs the classic spectral method (Method 1), in which the PSV estimation in a sample volume (SV) is obtained by the conversion of the maximum Doppler frequency (f max ) to velocity, according to the Doppler equation (3.1). In Fig. 3.1, the result of a CCA investigation is depicted.

Method 2: Angle Tracking Vector Doppler

The angle tracking method [START_REF] Tortoli | An automatic angle tracking procedure for feasible vector Doppler blood velocity measurements[END_REF][START_REF] Ricci | Accuracy and reproducibility of a novel dual-beam vector Doppler method[END_REF] uses two sub-apertures of a linear array probe to transmit, with independent steering angles, two focused beams (hereinafter referred to as "tracking" and "measurement" beams, respectively) which intersect in the SV of interest [Fig. 3.2(a)]. By using the symmetry properties of transverse Doppler spectra [START_REF] Newhouse | Ultrasound Doppler probing of flows transverse with respect to beam axis[END_REF], it is possible to automatically maintain the tracking beam at 90 • to the flow. In this way, the flow direction is continuously identified, and the knowledge of the beamto-flow angle, ϑ, can be used to estimate the velocity by the measurement beam. In this paper, the maximum frequencies of the Doppler spectra produced from the investigated SV have been converted to velocities by the following equation [START_REF] Tortoli | Improved blood velocity estimation using the maximum Doppler frequency[END_REF]:

v max = f max 2f 0 c (cos ϑ + k sin ϑ) , ( 3.2)
where k is a factor of proportionality depending on the transducer geometry. Notice that the angle is here used for evaluating both the frequency associated to the maximum velocity by the Doppler equation (first term in the denominator) and the intrinsic spectral broadening effects [START_REF] Evans | Doppler Ultrasound: Physics, instrumentation and signal processing[END_REF] (second term in the denominator).

Method 3: Plane Wave Vector Doppler

Steered plane waves are transmitted by two sub-apertures of a linear array probe. The left and right sub-apertures are symmetrically located with respect to the probe axis and use opposite steering angles, so that a wide region of interest (ROI) is insonified [START_REF] Ricci | Real-time vector velocity assessment through multigate Doppler and plane waves[END_REF]. The back-scattered echoes received by the same apertures are independently beamformed along the probe axis [START_REF] Jensen | Estimation of velocity vectors in synthetic aperture ultrasound imaging[END_REF]. The beamformed data are processed by the Welch method [START_REF] Evans | Doppler Ultrasound: Physics, instrumentation and signal processing[END_REF], to calculate the mean frequencies f d l , f d r of the Doppler power spectral densities related to the "left" and the "right" aperture, respectively. Such mean frequencies are then trigonometrically combined to produce the velocity magnitude and direction [START_REF] Dunmire | Cross-beam vector Doppler ultrasound for angleindependent velocity measurements[END_REF]. The plane wave VD method is thus capable of estimating, in real time, the velocity vector in many small SVs along the beamforming line [Fig. 3.2(b)].

Experimental Setup

A commercial scanner (HDI5000, Philips, Amsterdam, The Netherlands) operating in classic duplex mode with a 5 MHz linear array probe

Materials and Methods

51

was used to measure the PSV in carotid arteries by the conventional spectral estimation method. For each investigated artery, a screenshot like that in Fig. 3.1 was produced. The velocity value indicated in the lower left corner of the picture corresponds to the average of all the peak velocities (typically 3) extracted from the visible spectrogram.

Methods 2 and 3 were implemented on the ULA-OP research scanner ( §1.3) coupled with the LA533 192-element linear array probe (Esaote S.p.a., Florence, Italy). For both methods, the ULA-OP sequencer was programmed to alternate one block of 96 B-mode transmissions and one block of 256 vector Doppler transmissions, i.e. 128 for each sub-aperture, 128 being the ensemble size. In the latter case, the FPGA-based arbitrary waveform generators excited the elements of each sub-aperture at rate P RF/2, by means of 5-cycle tone bursts at 5.5 MHz. The received echo-signals (converted to digital at 50 MHz with 12-bit resolution) were dynamically beamformed by the FPGAs and further processed by the DSP to extract the velocity values. The refresh rate for imaging and Doppler data was typically around 30 Hz.

ULA-OP is supervised by a user-interface on a host PC connected by a USB 2.0 port. By the interface, the operator sets the applicative module, thus defines the corresponding transmit-receive (TX-RX) strategy, and selects the data storing modality. In both Angle Tracking and Plane Wave VD modalities, the user-interface was organized in two sections (Fig. 3.2) including the B-mode image, that highlights the selected SV (top) and the corresponding velocity magnitude profile (bottom). Although the plane wave mode allows calculating multiple vectors along the investigation line, the velocity magnitude profile is here related to the SV selected by the operator [identified by a circle in Fig. 3.2(b)]. In both vector methods, the PSV, calculated as the average of the last 3 peak systolic velocities, is shown in the upper-left corner of the velocity section.

Measurement Protocol for Methods Comparison

This study was approved by the local ethics committee and informed consent was received from each participant. 15 patients with significant carotid artery stenosis were recruited. Before the examination, each patient stayed in a comfortable supine position until the heart achieved a stable rate. The subject maintained the relaxed position while an expert sonographer performed the examinations.

Each investigation started from the left carotid artery and was concluded with the right carotid artery. For each subject, the operator put the ultrasound probe on the volunteer's neck, so that a longitudinal scan of the vessel was obtained, and the Doppler line(s) suitably crossed the vessel of Standard and Vector Doppler Approaches interest. The SV was set at least 2 cm proximal to the bifurcation apex in CCAs and 2 cm distal to the bifurcation apex, or in the point of maximum velocity, for ICAs.

First, the sonographer explored the ROI with the commercial machine. When the section of the vessel was identified, the physician set the angle cursor and manually fixed its inclination in parallel with the vessel walls. After adjusting the probe inclination and the beam steering to get the Doppler angle as close to 60 • as possible, the sonographer activated the spectral analysis. The value indicated in the spectrogram window was recorded as PSV.

The Plane Wave VD investigation was then started by taking care, based on the B-mode images, that the investigated ROI, and in particular the SV selected for PSV measurement, were the same as in Method 1.

Once the flow and the image were sufficiently stabilized, the acquisition was stopped. The last value indicated in the velocity window was assumed as PSV for the vessel under investigation. Maintaining the probe in the same position over the same ROI, the measurement was repeated by the Angle Tracking VD application.

Repeatability Measurements Protocol

An additional group of 8 healthy volunteers was recruited to evaluate, by two expert sonographers, the inter -operator and intra-operator repeatability performance of the three Doppler methods. The probe positioning procedure was similar to that described above. For each subject, the first sonographer investigated the right and the left CCA by using the commercial echograph. Each measurement was repeated for three times by the same operator by alternating the left and right carotids, for a total of 6 exams. Then, the same operation sequence was repeated by the second sonographer. Once terminated, the machine was substituted by ULA-OP and the procedure was repeated for the two vector methods, with the same investigation modality detailed in the previous paragraph. The complete examination of a volunteer lasted about 30 minutes, producing a total of 36 measurements.

Data Processing

The methods performance was compared by a pairwise Bland-Altman representation of the averaged PSV values [START_REF] Bland | Statistical methods for assessing agreement between two methods of clinical measurement[END_REF] and by a linear regression analysis.

The intra-operator repeatability was evaluated by comparing the 3 measurements available for each CCA. In particular 16 triplets were available for each method (8 volunteers x 2 sides) and operator. The coefficient of

Results

53

variability (CV ), i.e. the standard deviation normalized with respect to the mean value, was calculated for each triplet. A group of 16 CV s was obtained for each method and operator. The mean value CV s for each population was considered as the final result.

The inter -operator repeatability was assessed by matching the measurements produced by the 2 sonographers for the same method and carotid. In this case the 3 measurements available for each method/carotid/operator were first processed by excluding the most diverging value. The 2 remaining measurements were averaged. This produced 16 couples of values for each method, where each couple was constituted by the measurements obtained by the 2 different operators in the same condition. The difference normalized with respect to the mean value from each couple was calculated and the results constituted 3 groups, one for each method, populated by 16 elements. These groups were finally analyzed by extracting the mean, the standard deviation, the 25 th and 75 th percentiles, and the outliers.

Results

When, for a specific vessel, the sonogram quality obtained with any of the methods was not considered sufficient to perform a reliable PSV measurement, all of the results related to that vessel were discarded. An example is when the patient is affected by an atrial fibrillation, which leads to floating values of the PSV. After eliminating the outliers, full data from 22 CCAs and 15 ICAs of the patients remained available for the analysis. Since all data from the group of healthy volunteers were of good quality, none of them was discarded.

Methods Comparison

The panels in Fig. 3.3 report the comparison between the 3 possible method pairs according to Bland and Altman [START_REF] Bland | Statistical methods for assessing agreement between two methods of clinical measurement[END_REF]. For each measurement pair referred to the same patient, the mean value is reported on the horizontal axis and the difference on the vertical axis. Each mark presents either a CCA (black circle) or an ICA (red cross) measurement couple, respectively. The horizontal lines correspond to the average difference (continuous line) and ±(1.96S d ) (standard deviation) range (dotted lines). These values are separately listed in Tab. [START_REF] Culjat | A review of tissue substitutes for ultrasound imaging[END_REF]

Intra-operator Repeatability Measurements

The Tab. 3.4 reports the mean of the coefficient of variability CV s measured for evaluation of the intra-operator repeatability. The CV s , calculated as detailed in the previous section, are distinguished by method and operator and expressed in percent.

Intra-operator Repeatability [CV (%)].

Method 1 Method 2 Method 3 CCA 5.9 6 .0 6 .9 ICA 6.2 5 .3 6 .2

Inter -operator Repeatability Measurements

The assessment of inter -operator data from the 3 methods is performed by comparing the measurements carried out by the 2 sonographers in the same conditions. As detailed in the previous section, the measurements have been processed in three populations of 16 values. Each value repre- sents the difference between the measurements operated by 2 sonographers normalized with respect to the mean, and expressed in percent. Fig. 3.5 reports, for each method, the 25 th and 75 th percentile in the population

Results

55

Discussion

57

(box top and bottom edges), the median value (horizontal segments), the maximum and minimum (whiskers), the outliers (crosses). The mean and standard deviation of the population are reported in Fig. 3

Discussion

In this study, three different methods for the measurement of the blood peak velocity have been compared through examinations on volunteers and patients. The first method is the one currently used in daily clinical practice, the others are innovative methods with automatic angle compensation that have been implemented in the ULA-OP research platform. The 3 methods exploit completely different TX-RX strategies and velocity estimation techniques. Methods 1 and 2 are based on focused TX and single SV analysis, while Method 3 transmits plane waves and uses a multigate Standard and Vector Doppler Approaches approach. In Method 1 the angle is manually compensated, Method 2 exploits a specific reference/measurement dual-beam approach, Method 3 is based on triangulation. Furthermore, Methods 1 and 3 calculate the blood velocity by directly converting the maximum and the mean frequency, respectively, of the spectrum by the Doppler equation, while in Method 2 the conversion is applied to the maximum frequency after spectral broadening compensation.

In spite of such differences, the comparison between measurements pairs reported in Fig. 3.3 highlights an evident asymmetry: the agreement between the velocities measured by the VD methods (bottom panel) is much higher than that obtained when a VD method is compared to the Method 1 (top and middle panels). This result is quantified in Tab. 3.1. The mean velocity difference, lower than 3.5 cm/s when the VD methods are compared (see rightmost column), rises to more than 20 cm/s when Method 1 is compared to Methods 2 or 3 (see left and central column). These evaluations are valid for both ICA and CCA investigations. These data highlight a significant overestimation of Method 1 with respect to both Methods 2 and 3. The averaged PSVs are in fact in the range (80 -90) cm/s for Method 1 and as low as (50 -60) cm/s for Methods 2 and 3 (Tab. 3.2).

The same tendencies are confirmed by the linear regression analysis. The statistical results of Fig. 3.4(a) and (b) show poor agreement between the spectral and the vector methods. On the other hand, a better agreement is observed between the vector methods in Fig. 3.4(c), and confirmed by the values in Tab. 3.3. Compared to the CCAs, ICA investigations typically present lower signal to noise ratio due to the higher inclinations and deeper placement of the SV, and the results are further deteriorated. This is clearly shown in Fig. 3.4(b), in which the significant overestimation in the ICA measurements (symbolized as crosses) produces the drastic reduction of the coefficient a (see the first row in Tab. 3.3), and the agreement of results between the spectral and the vector method is compromised.

These results are the direct consequence of the fact that in the standard spectral Doppler method the maximum frequencies are converted to velocity by the Doppler equation (3.1). The effects of intrinsic spectral broadening are thus disregarded [START_REF] Yang | Assessment of spectral Doppler in preclinical ultrasound using a small-size rotating phantom[END_REF]. It is worth reminding that, since the broadening depends on the sin ϑ, such effects are different at different angles, i.e. the broadening is more consistent at high Doppler angles rather than at low angles [START_REF] Logason | The importance of Doppler angle of insonation on differentiation between 50-69% and 70-99% carotid artery stenosis[END_REF][START_REF] Steel | Angledependence and reproducibility of dual-beam vector Doppler ultrasound in the common carotid arteries of normal volunteers[END_REF].

In Tab. 3.3, the agreement between vector measurements (Method 3 vs. 2) is very high, presenting a regression coefficient (a) of 0.83 and an intercept value (coefficient b) of 7.00 cm/s (see rightmost column). This is strikingly low considering the lack of practice of the clinical operator, who was more familiar with the conventional spectral methods, implemented on the commercial machines employed in clinical routine.

The low values obtained for R 2 (Tab. 3.3, last row) suggest that, for the better comprehension of the clinical employment of vector methods, a wider clinical investigation is required, and that is the purpose of the analysis for repeatability (both intra-and inter -operator).

The reproducibility measurements confirm the results reported by Steel et al. [START_REF] Steel | Angledependence and reproducibility of dual-beam vector Doppler ultrasound in the common carotid arteries of normal volunteers[END_REF], who concluded that dual-beam vector Doppler systems can measure blood velocity with similar reproducibility to that of standard systems. The intra-operator repeatability tests have given similar CV s , i.e. between 5% and 7% with all methods for both operators (Tab. 3.4). Such low values are not surprising for the spectral Doppler approach, which was so familiar for both operators that they could easily reproduce the same conditions in consecutive measurements on the same subject/vessel. A similar Doppler angle (and possible related error) may thus be expected in consecutive measurements. However, from this point of view, obtaining similar CV s with the VD methods looks remarkable, since the operators had no previous experience with such methods.

As far as the inter -operator measurements are concerned, we did not find significant differences between the measurements performed by the 2 sonographers for all of the 3 methods (-[2-3]% according to values listed in Tab. 3.5). The measurements dispersion is about 11.4% for Methods 1 and 3 and 7.9% for Method 2 (see last row in Tab. 3.5 and the box extension in Fig. 3.5). The extremes of the measurement difference are around ±20% in Method 1, ±15% in Method 2 and ±10% (with an outlier at 30%) in Method 3 (see whiskers in Fig. 3.5). These results may be compared to those obtained by [START_REF] Corriveau | Interobserver variability of carotid Doppler peak velocity measurements among technologists in an ICAVL-accredited vascular laboratory[END_REF], who reported an inter -operator variability of (25 -40)% for PSV measurements. However, it should be underlined that their tests referred to the ICA of patients with stenosis, i.e. a more difficult measurement condition.

The use of different machines for the application of different methods might be considered a limitation of the current study. This was partly unavoidable, since the vector Doppler methods are new and their implementation was possible only on a research platform like ULA-OP. On the other hand, for Method 1, we preferred letting the operators work on machines that were already familiar to them, so that the results were not affected by a lack of practice. The fact that the operators needed some training, although limited, to use the new VD methods in a non standard machine, adds further value to the obtained results.

The in-vivo application of VD techniques has not highlighted specific critical situations. The relatively large number of discarded data mentioned in the Results section was mainly due to the conditions of the patient (irregular heartbeat or carotid occlusion) and not to the used method. In

Clinical Comparison of Blood Velocity Measurements by

Standard and Vector Doppler Approaches about 10% of the cases, aliasing was observed in the VD methods, because they use a P RF/2 TX-RX rate. However, this limitation could be easily overcome by adopting more sophisticated TX-RX strategies in which different sub-apertures are used in TX and RX phases, respectively [START_REF] Dunmire | Cross-beam vector Doppler ultrasound for angleindependent velocity measurements[END_REF].

ICA exams, of course, have turned out to be more critical, because of the higher difficulties (i.e. poorer signal-to-noise ratio) connected to scanning a deeper ROI, the reduced diameter and larger inclination with respect to the CCA. In severely stenosed ICAs (Fig. 3.6), in particular, higher velocities were measured (range: [86 -115] cm/s with VD, and [110 -201 cm/s with the spectral Doppler method), with inherent higher aliasing risk. Notwithstanding such difficulties, a significant number (15) out of the 30 ICA PSV measurements could be included in this study.

Conclusion

In this work, a comparison between the spectral Doppler and two novel vector Doppler methods has been presented. The results of this study suggest that the introduction of vector Doppler methods in commercial machines may finally be considered mature. The high measurement correlation obtained with the two VD techniques strengthens their suitability to perform accurate PSV measurements in carotid arteries.

In clinical trials, the spectral Doppler is considered an excellent tool for the detection of many cardiovascular diseases, especially if coupled to the B-mode (duplex mode) and Color Flow Imaging (triplex mode). The introduction of VD methods may finally represent the occasion to remove the further angle dependency correlated to spectral broadening [START_REF] Hoskins | A comparison of single-and dual-beam methods for maximum velocity estimation[END_REF]. As shown in this paper, this may be obtained calculating the mean frequency, provided that the SV is small enough to avoid a significant bias due to velocity gradient broadening. As an alternative closer to the spectral Doppler method, the velocity can be extracted from the maximum spectral frequency, by considering a suitable, angle-dependent, correction factor.

However, the incoming of new advanced methods allows for more thorough and accurate description of blood flow patterns in extended regions, e.g. stenosis and bifurcations. These results promote the inclusion of vector Doppler methods in commercial systems. Such applications are candidate to be complementary tools for vessel investigations, aimed at understanding the origin of stenosis and thereby paving the way to innovative solutions for the preventive diagnoses of cardiovascular diseases. [START_REF] Lenge | High-frame-rate 2-D vector blood flow imaging in the frequency domain[END_REF]

2-D Vector Blood Flow Imaging in the Frequency Domain

Summary. This Chapter introduces a novel method to produce bidimensional maps of 2-D velocity vectors. The region of interest is illuminated by plane waves transmitted at the pulse repetition frequency in a fixed direction. For each transmitted plane wave, the back-scattered echoes are recombined offline to produce the radio-frequency image of the ROI. The local 2-D phase shifts between consecutive speckle images are efficiently estimated in the frequency domain, to produce vector maps up to 15 kHz PRF. Simulations and in-vitro steady-flow experiments, with different setup parameters, and in pulsatile flow conditions have been conducted to thoroughly evaluate the methods performance, and presented in

Introduction

As previously illustrated in §1.2, spectral Doppler methods are typically used to estimate the axial component (1-D) of blood velocity in a single sample volume or along an entire line [START_REF] Evans | Colour flow and motion imaging[END_REF]. The estimate can be extended to the transverse component (2-D) by means of multi-beam approaches [START_REF] Dunmire | Cross-beam vector Doppler ultrasound for angleindependent velocity measurements[END_REF][START_REF] Tortoli | An automatic angle tracking procedure for feasible vector Doppler blood velocity measurements[END_REF], which combine more Doppler measurements. Speckle-tracking methods [START_REF] Trahey | Angle independent ultrasonic detection of blood flow[END_REF], employing a 2-D cross-correlation, follow the scatterers movements in 2-D. They have been recently coupled to high frame-rate imaging methods to obtain extended 2-D vector maps [START_REF] Udesen | High frame-rate blood vector velocity imaging using plane waves: Simulations and preliminary experiments[END_REF]. This technique has the strong potential to extend both the region of interest and the frame-rate, in order to map fast flows in wider vessel sections. In [START_REF] Ekroll | Simultaneous quantification of flow and tissue velocities based on multi-angle plane wave imaging[END_REF] and [START_REF] Bassi | Real-time vector velocity profile measurement based on plane waves transmission[END_REF], an alternative dual-beam Doppler vector imaging method has recently been presented, based on the transmission of two plane-waves with opposite steering angles.

This Chapter draws on [START_REF] Lenge | Frequency-domain high frame-rate 2D vector flow imaging[END_REF] and [START_REF] Lenge | High-frame-rate 2-D vector blood flow imaging in the frequency domain[END_REF] 1 and concerns a novel method to produce bi-dimensional maps of 2-D velocity vectors. The method employs PW transmission to illuminate the ROI in a single shot. No compounding [START_REF] Montaldo | Coherent plane-wave compounding for very high frame rate ultrasonography and transient elastography[END_REF] is adopted in the reconstruction of the RF images, to keep the Doppler pulse repetition frequency (PRF) as high as possible. The phase shifts caused by blood movements over the ROI are estimated in the frequency domain, as previously proposed for 1-D elastography investigations [START_REF] Ramalli | Frequency-domain-based strain estimation and high-frame-rate imaging for quasi-static elastography[END_REF]. The preliminary method validation presented in [START_REF] Lenge | Frequency-domain high frame-rate 2D vector flow imaging[END_REF] has been thoroughly confirmed in [START_REF] Lenge | High-frame-rate 2-D vector blood flow imaging in the frequency domain[END_REF], by in-vitro and in-vivo experiments. It has been also estimated that the computation of the frequency domain algorithm is more than 50 times faster than the computation of the reference 2-D cross-correlation algorithm. An efficient implementation of the method on a GPU board has been joined ( §6).

Flow Imaging Method

Plane Wave Imaging

RF images are obtained by transmitting plane waves along a single direction at PRF rate with a linear-array probe, and beamforming (in parallel) the RF echo data received by the active array elements. By avoiding any compounding [START_REF] Montaldo | Coherent plane-wave compounding for very high frame rate ultrasonography and transient elastography[END_REF] [ §1. 1.4], the Doppler sampling frequency is kept co- incident with the PRF. In the following, it is assumed that the RF image is reconstructed from the echo signals, sampled at f s rate, obtained at interline distance (δ x ) equal to the probe pitch. To remove the tissue clutter before any Doppler processing, each RF frame is pre-filtered by subtracting the average of the previous 20 RF frames from itself [START_REF] Udesen | High frame-rate blood vector velocity imaging using plane waves: Simulations and preliminary experiments[END_REF].

2-D Displacement Estimation in the Frequency Domain

Let us consider a rectangular matrix of M B ×N B samples, selected from an RF frame as illustrated in Fig. 4.1. Hereafter, the RF data in the matrix, weighed by a Hanning window in both directions, will be called a "block". The pixels contributing to the block are spaced by δ z = c/(2f s ) in the axial direction, and by the distance, δ x , between adjacent reconstructed lines in the lateral direction. A block will be described by the matrix s mn , where m and n are the indices in the axial and lateral directions, respectively.

If the scatterers inside the ROI are moving, as happens in blood, the blocks, s mn and smn , taken at consecutive pulse repetition intervals, t and t, respectively, are thus different. If the scatterers, during the time span τ = t -t, have covered short Δz and Δx displacements in the axial and lateral directions, respectively, smn can be considered as equivalent to a 2-D-shifted version of s mn .

As indicated in the block diagram in Fig. 4.2, the scatterers' displacement can be estimated starting from the 1-D and 2-D discrete Fourier

S f m n = M B -1 m=0 s mn e -i( 2π M B )f m m (4.1) S f m f n = N B -1 n=0 M B -1 m=0 s mn e -i[( 2π M B )f m m+( 2π N B )f n n] (4.2)
where f m and f n are the spatial frequencies in the axial and lateral directions, respectively. The well-known shift-property of the DFT correlates the spectra of two consecutive blocks through

Sf m n ∼ = Sf m n+Δx = S f m n e i( 2π M B )f m Δz (4.3) Sf m f n = S f m f n e i( 2π M B )f m Δz+i( 2π N B )f n Δx . (4.4)
Note that (4.3) assumes that the lateral displacements are not high enough to determine the passage of scatterers from one block to the next. Because the minimum lateral distance between two (partially overlapped) adjacent blocks is equivalent to the distance between two lines, δx, the maximum lateral displacement should be ≤ δx/2. For example, for 15 kHz frame rate and δx = 245 µm, the lateral velocity should be lower than 1.83 m/s. From (4.3) and (4.4), the displacements along z and x can be computed through the frequency-domain phase shifts, estimated as

Δφ f m n = arg[ Sf m n ] -arg[S f m n ] = 2π M B f m Δz (4.5) Δφ f m f n = arg[ Sf m f n ] -arg[S f m n ] = 2π M B f m Δz + 2π N B f n Δx. ( 4

Flow Imaging Method 67

Δz can be extracted from (4.5) as

Δz = M B 2πf m Δφ f m n . (4.7)
Substituting Δz in (4.6) and solving for Δx, the lateral displacement is obtained as

Δx = N B 2πf n (Δφ f m f n -Δφ f m n ). (4.8)
In (4.7) and (4.8), the displacements are expressed in sample units. In principle, one pair of frequencies (f m ,f n ) would be sufficient to calculate the phase shifts and the relative bi-dimensional displacement. To improve the robustness in operative conditions, it is however convenient to average a number (n f ) of estimates obtained from different pairs of frequencies. As further discussed later, one possible choice is selecting n f equally spaced frequencies, covering a z -axis bandwidth B z around the central frequency of the transmitted signal (f 0 ), and 2n f equally spaced lateral frequencies, n f over the lateral positive band B x = (0, +1/2δx) and n f over the lateral negative band B x = (-1/2δx, 0). The mean displacements are calculated as follows:

Δz = 1 2n f M B N B n=1 2n f m=1 M B 2πf m Δφ f m n (4.9) Δx = 1 2n 2 f n f n=1 2n f m=1 N B 2πf n (Δφ f m f n -Δφ f m n ). (4.10)
Finally, the mean components of the velocity vector in the block during the interval τ are v z = Δzδz/τ (4.11)

v x = Δxδx/τ. (4.
12)

The complete 2-D vector map is obtained by applying the described method to all blocks, which may be partially overlapped. In this paper, a 90% overlap will be considered. To smooth the velocity profiles, a slowtime low-pass temporal filter is applied on the velocities estimated from consecutive frames. In the following, the cut-off frequency of such filter will be fixed to 300 Hz, high enough to preserve the pulsating nature of flow in the main human arteries.

Method Validation

The proposed method has been tested through simulations and experiments. In both cases, a linear-array probe with 110% bandwidth centered at 8 MHz (pitch 245 µm) has been considered. As illustrated in Fig. 4.1, the probe was placed in front of a cylindrical pipe having 4 mm internal radius, 1-mm-thick walls, and axis at distance z 0 from the center of the array. The tube axis was tilted by an angle ϑ with respect to the beam axis. 64 elements of the probe were simultaneously excited by five-cycle sinusoidal bursts at f 0 central frequency, weighted by a Hanning window, to produce plane waves of 15.4 mm width (i.e., the width of the active aperture). In receive mode, the RF sampling frequency (f s ) was set to 50 MHz, corresponding here to δz = 14.8 µm for a sound speed of c = 1480 m/s. The velocity was estimated on 32 uniformly distributed lines, covering a 7.7-mm wide ROI, corresponding to the central part of the region insonified by plane waves. This choice is suggested by the poor contrast typically obtained along the peripheral lines with this type of transmission [START_REF] Montaldo | Coherent plane-wave compounding for very high frame rate ultrasonography and transient elastography[END_REF]. Nevertheless, the estimation region turns out to be sufficiently wide to permit a suitable evaluation of blood flow dynamics.

Simulations

The Field II ® software [START_REF] Jensen | FIELD: A program for simulating ultrasound systems[END_REF][START_REF] Jensen | Calculation of pressure fields from arbitrarily shaped, apodized, and excited ultrasound transducers[END_REF] has been used to simulate a flow within a 30-mm-long vessel, located at z 0 depth, with beam-to-flow angle ϑ. Inside the vessel, a parabolic flow moves with (v 0 ) peak velocity. The scatterer density was fixed at 25 scatterers/mm 3 , and the back-scattering intensity was 20 dB lower than that of the walls. In Tab. 4.1 the bold values indicate the so-called reference setup (RS). Each parameter was individually changed to verify its effect on the method performance. The speed of sound was considered 1480 m/s, and 520 RF images were reconstructed for each simulation.

A pulsatile flow simulation was also implemented, reproducing a physiological arterial velocity profile during a cardiac cycle (1.0 s long) at RS conditions. In this case, the vessel barycenter that was set to 20 mm and the central streamline reached a peak velocity of only 20 cm/s in the systolic phase, so the simulation time can be limited by derating the PRF to 3 kHz.

Experiments

The in-vitro set-up included a blood-mimicking fluid flowing from an upper reservoir into a cylindrical tube immersed in a water tank. The tube was made by Rilsan ® , a plastic material characterized by an attenuation of data were acquired at 5 kHz frame rate and covered one cardiac cycle.

Performance Metrics

In steady flow conditions, by using the peak velocities specified in the previous section, a steady parabolic flow was guaranteed, so that the axial and lateral velocity components, v z and v x , could be predicted for each point of the ROI. The algorithm performance was evaluated by comparing such values with the values, vz and vx , estimated in simulations and experiments. For the axial component, the relative bias (B v z ) and the standard deviation (σ v z ) are defined as follows:

B v z = 1 N l N l j=1 1 v 0 z N d N d i=1 1 N f N f -1 k=1 |v z (z i , x j , p k ) -v z (z i , x j )| (4.13) σ v z = 1 N l N l j=1 1 v 0 z 2 N d N d i=1 1 N f N f -1 k=1 (v z (z i , x j , p k ) -v z (z i , x j )|) 2 , (4.14)
where vz (z i , x j , p k ) is the estimated axial velocity in the point of coordinates (z i , x j ) of the frame p k , v z (z i , x j ) is the theoretical axial velocity component, v 0 z is the axial velocity component of the peak velocity, N l is the number of estimation lines [START_REF] Fox | Multiple crossed-beam ultrasound Doppler velocimetry[END_REF], N d is the number of axial samples in each profile and N f is the total number of frames (499). Note that N d depends on the vessel radius (8 mm), the beam-to-flow angle (ϑ) and the sampling frequency (50 MHz). B v x and σ v x are evaluated in equivalent way for the lateral velocity components.

Results

71

Results

Simulation Results

Static-Flow Simulations:

The performance sensitivity to the variation of parameters listed in Tab. 4.1 has been evaluated. In each test, all RS values were kept unaltered, except the value corresponding to the parameter under consideration. First, Fig. 4.4(a) shows the performance metrics obtained at different central frequencies (f 0 ) and frame-rates (F R). The errors obtained when the beam-to-flow angle was swept from 90 • (i.e. the tube parallel to the probe) to 30 • with 10 • steps, are reported in Fig. 4.4(b). In Fig. 4.4(c), the block barycenter depth z 0 is varied from 15 to 35 mm, and in Fig. 4.4(d) the flow peak speed v 0 is changed in 30, 60, and 100 cm/s.

The bias and the standard deviation of both velocity components have been estimated also by comparing rectangular (not apodized), Tukey (50% tapering), and dynamic sinc (with constant receive f-number = 2) apodization windows. As shown in Fig. 4.4(e), better results, in terms of lateral estimation are achieved without apodization, whereas the Tukey and the sinc apodizations improve the performance for axial estimation. Fig. 4.4(f) presents the performance obtained with different numbers of estimation frequencies, n f . The parameter which appears to be more affected by this number is the bias of the lateral velocity component (B v x ), but extending n f beyond 5 does not improve the results because the difference between adjacent frequencies becomes too low to produce uncorrelated estimates.

The possible effects of z -axis bandwidth changes have been evaluated, but the results in Fig. 4.4(g) show that no major improvements are obtained for bandwidths wider than 0.5 MHz.

Figs. [START_REF] Rayleigh | On the light from the sky, its polarization and colour[END_REF].4(h) and 4.4(i) show the effects of varying the block size. The axial estimation is not heavily influenced by the axial block size but the lateral velocity estimation is considerably affected by the lateral block size. Indeed, the simulations indicate that the wider the block is, the better is the performance, and lateral dimensions lower than 2 mm should not be used (Fig. 4.14 in §4.5).

Pulsatile-Flow Simulations:

The accompanying movie [START_REF] Lenge | Sequence of 2-D vector maps obtained with the pulsatile-flow simulation[END_REF] shows the sequence of 2-D vector maps obtained with the pulsatile-flow simulation.

Experimental Results

In-Vitro Experiments: With reference to the parameters reported in Tab. 4.2, the beam-to-flow angle was limited within the range 90 • to 70 • to avoid working close to the Rilsan ® tube critical angle, at which total reflection occurs. This point is further discussed in §4.5. As shown in Fig. 4.7(a), the errors, and the axial bias (B v z ) in particular, noticeably increase by reducing the angle. The distance between the probe and the tube (barycenter depth) affects the performance, too, especially at higher depths where the lateral bias increases because of the drop of the lateral resolution [Fig. 4.7(b)]. The bias remains close to 20% when the peak velocity is varied [Fig. 4.7(c)] except at 5 cm/s, where a 33% bias on the lateral estimation is obtained, resulting from the clutter filtering, which removes the slowest velocity components. Fig. 4.7(d) shows that a receive apodization slightly affects the performance, reducing the bias of the axial velocity against a small increase of the lateral bias. Fig. 4.7(e) presents the dependence by the number of frequencies, n f , whereas Fig. 4.7(f) shows the effects of the bandwidth, B z , over which the axial frequencies are selected. Finally, Figs. 4.7(g) and 4.7(h) illustrate the performance for different axial and lateral block sizes. It is worth noting that in all test conditions the bias is due to a velocity underestimation in both lateral and axial directions. This result is similar to that obtained in conventional pulsed Doppler systems, when the mean frequency is extracted from the Doppler spectrum, which is explained by the spread of velocities present in a sample volume [START_REF] Evans | Doppler Ultrasound: Physics, instrumentation and signal processing[END_REF]. Considering the performance metrics, the most accurate is that of lateral velocity at the systolic peak, V x (T 1 ), showing B v x ± σ v x = (2.7 ± 2.4)%, and the worst is that of V x (T 3 ) which presents B v x ± σ v x = (9.6 ± 9.3)%. Note that these results are related with the peak velocity component (v x 0 and v z 0 ).

A simple in-vitro experiment has been made to illustrate the method accuracy in estimating the velocity vector angle. The beam-to-flow angle was slowly changed between 80 • to 100 • . Since flow was maintained laminar, the true Doppler angle was considered the one deducted from the B-mode images. In Fig. 4.8, the estimated beam-to-flow angle in one block at the center of the pipe (black circles) is compared with the true Doppler angle (red asterisks).

In-Vivo Experiments:

The in-vivo accompanying movie [START_REF] Lenge | Blood velocities in the CCA and in the jugular vein recorded during a cardiac cycle from one 28-year-old healthy volunteer[END_REF] shows the blood velocities in the CCA and in the jugular vein recorded during a cardiac cycle from one 28-year-old healthy volunteer. Fig. 4.9 illustrates the blood velocity components at the center of the CCA whereas Fig. 4.10 compares the spatial profiles simultaneously measured in the two vessels at three different instants: close to the peak systole (T 1 ), in late systole (T 2 ) and in diastole (T 3 ). beam-to-flow angle near 80 • . In Fig. 4.12, two other CCAs are depicted and, in Fig. 4.13, two screenshots of the femoral vessels, obtained in different time instants, are presented.

Discussion

A novel method for the 2-D blood motion estimation in extended bidimensional ROIs has been presented. The plane wave imaging method has already been thoroughly evaluated by other authors, who estimated the related point spread function and compared it to that obtained with conventional (focused) B-mode imaging (see, e.g., [START_REF] Montaldo | Coherent plane-wave compounding for very high frame rate ultrasonography and transient elastography[END_REF] and [START_REF] Udesen | High frame-rate blood vector velocity imaging using plane waves: Simulations and preliminary experiments[END_REF]). PWs have already been used by other groups for Doppler imaging applications [START_REF] Bercoff | Ultrafast compound Doppler imaging: providing full blood flow characterization[END_REF][START_REF] Ekroll | Simultaneous quantification of flow and tissue velocities based on multi-angle plane wave imaging[END_REF]. With respect to previous approaches, the method presented here is characterized by not using any compounding, so that each velocity estimate is based on the transmission of only two consecutive PWs. Although this approach unavoidably restricts the attainable image quality and measurement accuracy, it is crucial to maintain the Nyquist limit equal to P RF/2. Furthermore, the subsequent slow-time low-pass filtering applied to the estimated velocities is intended to smooth the velocity profile in each pixel but does not cancel the fastest physiological velocity contributions, as it would happen if filtering were directly applied to the raw image data.

The results shown in Figs. 4.3,4.5,4.6,4.8 and Figs. 4.9, 4.10, 4.11, 4.12, 4.13, as well as the accompanying movies [START_REF] Lenge | Sequence of 2-D vector maps obtained with the pulsatile-flow simulation[END_REF] and [START_REF] Lenge | Blood velocities in the CCA and in the jugular vein recorded during a cardiac cycle from one 28-year-old healthy volunteer[END_REF], indicate that the method is capable of estimating the velocity components in pulsatile flow conditions. In particular, the in-vivo results of Figs. 4.9, 4.10, 4.11, 4.12, 4.13 show that physiological details can be clearly detected. An example is represented by the M-shaped profile at late systole (T 2 ), which is known to be due to a high velocity gradient coupled to the large CCA curvature [START_REF] Tortoli | Interaction between secondary velocities, flow pulsation and vessel morphology in the common carotid artery[END_REF]. Fig. 4.11(b) demonstrates that the method also works nicely in vessels with reduced lumen. Fig. 4.11(a) indicates some blooming effect, probably highlighted by the transverse vessel wall orientation to the incident plane wave. To limit this effect in the final images, a simple mask has been added to avoid any velocity calculation when the B-mode data intensity is above a given threshold.

Quantitative data have been extracted from the simulations of laminar flow in the different conditions reported in Table 4.1. The simulation results (Fig. 4.4) indicate that in general, although the velocities are calculated from only two consecutive frames, accuracies better than 10% and precisions better than 5% can be obtained. It may be noted that even though, according to (4.5) and (4.6), the error on the axial movement is expected to influence the error in the lateral movement, the effect is not dramatic, because B v x and B v z are quite close to each other. The obtained accuracies appear to be substantial, considering that they have been obtained by limiting the clutter filter to a simple echo canceller, as in [START_REF] Udesen | High frame-rate blood vector velocity imaging using plane waves: Simulations and preliminary experiments[END_REF], and slow-time averaging to a mean temporal filter over 20 consecutive estimates. Such constraints have been suggested by the need to limit the computational load in view of a possible real-time implementation.

The simulations clearly show that the method performance is quite good in all cases where aliasing is not produced. The aliasing effects are highlighted by varying both the central frequency and the frame rate [Fig. 4.4(a)]. At constant reference velocity and beam-flow angle, higher aliased estimations are obtained at higher central frequencies or lower frame-rate: increasing the frame-rate from 10 to 15 kHz and decreasing the central frequency from 8 MHz to 6 MHz, the lateral bias has been reduced from 26.1% to 6.2%. The previous solutions are effective if the beam-to-flow angle remains higher than 60 • [Fig. 4.4(b)]: in this case, the method presents performance with errors lower than 10%. Reducing the Doppler angle, the aliasing limit becomes remarkable and the performance falls down, especially for the lateral velocity where the bias rises to 40% at 30 • . Such a high bias value is especially due to the aliasing limit in the axial direction, which causes the incorrect tracking of the fastest movements, and to the speckle de-correlation. Changing the barycenter depth [Fig. 4.4(c)], the performance are stable, although a deterioration in the lateral bias (B v x ) is notable at higher depths, due to the known reduction of the lateral resolution that affects the plane wave imaging. Varying the peak velocity [Fig. 4.4(d)], the performance are almost stable maintaining the error bias lower than 7%.

The main critical issue suggested by the simulation results is related to the lateral block size, which, according to Fig. 4.4(i), should not be kept below 2 mm. This value should not be considered too high, especially considering the Hanning windowing which is superimposed to each block. The lateral block size constraint is particularly relevant in in vivo acquisitions, in which the blood movement is typically parallel to the probe surface (e.g., in CCA investigations). The experimental results concerning the axial [Fig. 4.7(g)] and the lateral [Fig. 4.7(h)] block sizes confirm that increasing the number of samples in the lateral direction improves the performance, likely because it partially compensates the poor lateral resolution (δx), which is more than 10 times lower than the axial resolution (δz). Opting for a wider block is also convenient to balance the decrease of lateral contrast resulting from the intrinsic limits of the plane wave imaging reconstructions. Furthermore, by looking at the distribution of velocity profiles in Fig. 4.14, it can be observed that the profiles extracted from the boundaries of the image are degraded with respect to the central ones. Using wider blocks, the estimations get closer to the expected profiles and performance improves.

Experimental results present poorer performance than simulations, as expected considering the effects of attenuation and wall reflections as well as the difficulty of reproducing well controlled flow conditions (i.e., probe-tube alignment, flowmeter accuracy, etc.). In particular, both the axial and the lateral bias tend to increase at decreasing Doppler angles [Fig. 4.7(a)] [START_REF] Swillens | Two-dimensional blood velocity estimation with ultrasound: speckle tracking versus crossed-beam vector Doppler based on flow simulations in a carotid bifurcation model[END_REF]. This is not due to aliasing, because the maximum velocity here was only 16.6 cm/s. The aliasing effects have been previously considered, where it was shown that this problem can be tackled by tuning the frame rate and the central transmission frequency (p. 82). In the experiments here reported, the performance limitation can be attributed to the use of plane waves in combination with a tube of Rilsan ® material, which is characterized by a critical angle larger than 60 • . When the plane wave collides with the water-Rilsan ® interface at low beam-to-flow angles (e.g., 70 • ), most of the power is reflected from the first wall interface instead of being transmitted in the blood-mimicking fluid. The back-scattered echoes have a very low energy, and the poor S/N degrades, in particular, the axial estimations at the boundary of the ROI. The in-vivo investigations seem to not be affected by this problem, confirming that it is probably specific to the plastic material used in the experiments.

2-D Vector Blood Flow Imaging in the Frequency Domain

The performance of our method depends on depth [Fig. 4.7(b)]. At lower depths, the artifacts introduced by the grating lobes cause a reduction of the image contrast, influencing the axial bias. By increasing the depth, the contrast increases in favor of the axial estimation, but lateral errors increase because of the drop in lateral resolution.

Using larger numbers of frequencies (n f ) [Fig. 4.4(f) and Fig. 4.7(e)] improves the performance, but not dramatically when n f > 3 because averaging more estimations obtained by correlated speckle patterns becomes ineffective. Hence, to avoid fruitlessly increasing the computational load, n f in the range 3 to 5 represents a good trade-off. The aforementioned performance is not easily comparable with the performance reported by other authors, typically using different experimental conditions and adopting different metrics. Nevertheless, it is worth trying to compare the results of this study with the results obtained by Udesen et al. [START_REF] Udesen | High frame-rate blood vector velocity imaging using plane waves: Simulations and preliminary experiments[END_REF] in steady-flow measurements. In this paper, the bias and standard deviation have been evaluated by comparing the individual estimated velocity profiles with the expected profiles. In [START_REF] Udesen | High frame-rate blood vector velocity imaging using plane waves: Simulations and preliminary experiments[END_REF], the accuracy and precision of the averaged velocity profile produced by 480 frames turned out to be B v z ± σ v z = (0.7 ± 0.5)% and B v x ± σ v x = (6.4 ± 4.6)% for the axial and lateral velocity components, respectively. If the same metrics were adopted in this paper, comparable values would be obtained:

B v z ± σ v z = (2.0 ± 0.9)% and B v x ± σ v x = (5.8 ± 2.7)%.
Comparing the performance with the other vector methods presented in Chapter 3 [START_REF] Tortoli | An automatic angle tracking procedure for feasible vector Doppler blood velocity measurements[END_REF][START_REF] Ricci | Real-time vector velocity assessment through multigate Doppler and plane waves[END_REF], similar results are obtained. For example, the Multi-Gate Vector Doppler method obtains a 0.6 cm/s bias on 0.50 m/s peak velocity, corresponding to a 12% relative bias. However, the frequencydomain technique has the powerful advantage of calculating the maximum velocity in an extended section of the vessel, as indicated by the maximum value in the color-bars of Figs. 4.11, 4.12, and 4.13.

Computational Load Performances:

The proposed method's complexity has been evaluated in comparison with the block-matching method based on 2-D cross-correlation presented in [START_REF] Udesen | High frame-rate blood vector velocity imaging using plane waves: Simulations and preliminary experiments[END_REF], used as reference. Because the phases can be extracted from a look-up table, the computational cost is assessed by considering the needed number of sums (N sum ) and multiplications (N mult ). The complexity of sums or of multiplications can be considered equivalent in current digital signal processor (DSP) implementations. Accordingly, we have considered the total number of operations (N = N sum + N mult ) needed to obtain a velocity vector from a M B × N B block with n f = 5. The processing cost of the entire ROI inside a frame can be obtained by multiplying this value by the number of possible blocks defined inside the ROI. Referring to the method's formula in §4.2.2, N can be estimated as

N ∼ = 20M B N B + 200N B , (4.15)
where the velocity is obtained by a temporal average on 20 frames. To reach the sub-pixel accuracy, the 2-D cross-correlation algorithm entails interpolation procedures of the RF image, which can be partially alleviated by suitable optimizations [START_REF] Giunta | Fine estimators of two-dimensional parameters and application to spatial shift estimation[END_REF]. Supposing an interpolation factor of 4 in the lateral direction, as used in [START_REF] Udesen | High frame-rate blood vector velocity imaging using plane waves: Simulations and preliminary experiments[END_REF], and maintaining the same axial sampling frequency (f s ), the cross-correlation cost is

N ∼ = M 2 B N 2 B + 6M B N B . ( 4.16)
Replacing in (4.15) and (4.16) the RS block size (i.e., in samples, M B = 102, N B = 12), the computational demand of the proposed algorithm is decreased by a factor of 56 with respect to the reference algorithm. Of course, for both methods, further optimization could be conceived. Note that this factor increases up to 117 if n f is set to 3, as suggested by the

Conclusion and Perspectives

In this Chapter, a new method for blood motion estimation in extended bi-dimensional regions of interest has been proposed. The potentialities and the limits of the method have been investigated in terms of performance metrics (bias and standard deviation) and computational load. Plane wave imaging without any compounding has been proved as a suitable solution to insonify extended ROIs while exploiting the maximum Doppler PRF, although with unavoidable consequences in terms of image quality. PWs have been coupled with the Fourier-domain algorithm to create bi-dimensional flow maps. The results show that, although no statistical optimization of the estimator has been made, the method is capable of providing vector images with estimation errors lower than 10% and 20% in most simulations and experiments, respectively. Furthermore, the computational load of this algorithm has been demonstrated to be lighter than the reference block-matching method based on 2-D cross-correlation. Finally, the pulsatile-flow simulation and the CCA acquisition designate the method as a possible imaging tool for clinical investigation of cardiovascular diseases.

Future perspectives of the method may concern the statistical optimization. The weighted average, used to convert the phase shifts into the displacements in (4.7) and (4.8), introduces an alteration in the phase error distribution, that after the averaging is no longer uniform. An increase of performance is thus possible, developing the statistical mathematical model to obtain the optimal estimator, with the minimal variance [START_REF] Kay | Fundamentals of Statistical Signal Processing[END_REF]. Optimizations in terms of computational load and execution time, that are further described in Chapter 6, open the door towards the clinical tests.

Introduction

As already described in §1.2, one of the most important limitations of conventional US blood flow imaging is the estimation of only one component of the blood velocity vector (i.e. along the axial direction) [START_REF] Evans | Colour flow and motion imaging[END_REF]. With the purpose of extending the estimation in all the other directions, one of the most prolific techniques proposed in literature is the transverse-oscillation theory (TO). In general, while the estimator presented in Chapter 4 is more affine to speckle-tracking algorithms3 , the TO approach recalls the Doppler technique and, as suggested by D. H. Evans in [START_REF] Evans | Colour flow and motion imaging[END_REF], it can be enumerated among the multi-beam approaches (previously introduced in 1.2.3). The basic principle of TO approach takes inspiration from the Doppler principle. The following formula is adapted from (1.14) and defines the Doppler frequency f d :

f d = 2| v| cos ϑ c f 0 , ( 5.1)
in which | v| cos ϑ is the target velocity along the beam direction, f 0 is the transmitted frequency and 2/c is the reciprocal of the US velocity in the tissue, considering the forward and back propagation. In other words, a scatterer moving inside the ultrasound field gives rise to a sampled signal that is a function of the spatial field, the sampling time and the scatterer's velocity [START_REF] Jensen | A new method for estimation of velocity vectors[END_REF], and introduces a modification in the point spread function of the acoustic field. In Doppler systems, the Doppler frequency is revealed in the received signals by a phase modulation, that depends on the scatter velocity and the central frequency of the axial oscillation. The phase shift is then used to detect the blood velocity, as described in 1.2.1. Introducing a controlled oscillation in the transverse direction, the estimation can be extended to the transverse velocity according to the same principle sketched above. The received signal can be adapted from 1.15 as [START_REF] Jensen | Estimation of Blood Velocities Using Ultrasound: A Signal Processing Approach[END_REF][START_REF] Jensen | A new estimator for vector velocity estimation[END_REF]:

r i = cos[2πf d iT P RF + ϕ 0 ] cos 2π v x d x iT P RF + ϕ 0 (5.2)
where d x is the lateral modulation period and v x is the transverse velocity. The frequency of the transverse modulation is the term enclosed in brackets, i.e. f x = v x /d x . By a following spatial quadrature or a Hilbert transform, the one-sided spectrum can be obtained, and the sign of the velocity (positive or negative) detected.

The idea of introducing a transverse oscillation for multi-dimensional velocity estimations by US has two ancestors at the same time, i.e. [START_REF] Jensen | A new method for estimation of velocity vectors[END_REF] and [START_REF] Aderson | Multi-dimensional velocity estimation with ultrasound using spatial quadrature[END_REF]. J. J. Jensen [START_REF] Jensen | A new method for estimation of velocity vectors[END_REF] and his group had made the effort for the extensive development and improvement of the method in those years. Indeed, after the thorough explanation of the theory in [START_REF] Jensen | A new method for estimation of velocity vectors[END_REF], the TO method was exploited for blood [START_REF] Jensen | A new estimator for vector velocity estimation[END_REF] and tissue [START_REF] Liebgott | PSF dedicated to estimation of displacement vectors for tissue elasticity imaging with ultrasound[END_REF] motion estimations, and deeply investigated in [START_REF] Udesen | Investigation of transverse oscillation method[END_REF]. Then, the theory was validated in-vivo by a comparison with the MR phase contrast angiography [START_REF] Hansen | In vivo validation of a blood vector velocity estimator with MR angiography[END_REF], including also other ultrasound vector velocity techniques [START_REF] Hansen | In vivo comparison of three ultrasound vector velocity techniques to MR phase contrast angiography[END_REF]. Subsequently, TO imaging has been implemented in real-time on a commercial ultrasound system (Vector Flow Imaging, BK Medical ApS, Herlev, Denmark), adapted for phased-array investigations [START_REF] Pihl | Phased-array vector velocity estimation using transverse oscillations[END_REF] and, more recently, extended for 3-D vector estimations in a sample-volume ROI [START_REF] Pihl | A transverse oscillation approach for estimation of three-dimensional velocity vectors, part I: concept and simulation study[END_REF].

To generate the transverse oscillations in the pulse-echo field, many approaches are possible [START_REF] Liebgott | Transverse oscillations for tissue motion estimation[END_REF]. In principle, for TO imaging, only two beamformers are needed in receive, and that constitutes an advantage if compared with cross-correlation approaches [START_REF] Bohs | Speckle tracking for multi-dimensional flow estimation[END_REF][START_REF] Udesen | High frame-rate blood vector velocity imaging using plane waves: Simulations and preliminary experiments[END_REF][START_REF] Udesen | Investigation of transverse oscillation method[END_REF]. However, this may not be completely true right now, with the modern technologies of parallel computing. Moreover, TO approaches are not yet proposed to overcome the limitations in terms of frame-rate, because a focused line-by-line approach is still adopted in transmission in the same way as the conventional methods.

In this work, plane waves are combined with transverse oscillations, to propose a new approach for blood flow imaging. Plane-wave imaging is capable to procure high frame-rate descriptions of the fast blood flows in large ROIs. The images are obtained by a plane-wave transmission, and the back-scattered echoes are coherently summed by introducing a dynamic apodization that generates, in the pulse-echo field of the entire RF image, controlled oscillations in both the axial and the transverse directions. Such double oscillating images are finally processed by a phase-based estimator, previously proposed for elastographic applications [START_REF] Basarab | Phase-based block matching applied to motion estimation with unconventional beamforming strategies[END_REF], to obtain the 2Dvector maps that describe the blood flow in extended ROIs. The method was preliminary proposed in [START_REF] Lenge | A new method for 2D-vector blood flow imaging based on unconventional beamforming techniques[END_REF] and [START_REF] Lenge | Plane-wave transverse-oscillation for high frame rate 2-D vector blood flow imaging[END_REF] 4 , and here validated by further simulations and experiments in several conditions.

Plane-Wave Transverse-Oscillation for High Frame Rate 2-D

Vector Blood Flow Imaging

Methods

Plane-Waves Transverse-Oscillations Imaging

Conventional TO approaches use focusing of the emitted field in transmission, and a delay-and-sum beamforming on echo signals to obtain the dynamic beamforming in receive. One possible solution to generate the transverse oscillations in the pulse-echo field is applying, during the beamforming in transmission or reception, a special dynamic apodization, according to a line-by-line procedure. In order to get double-oscillating fields more quickly with respect to the conventional TO approaches (also 100 times quicker), plane waves are transmitted in the region of interest, and the transverse oscillations are generated during the beamforming in receive, in the frame of Fraunhoffer approximation.

Plane waves are transmitted at high frame-rate, by simultaneously exciting a subset of adjacent elements of a linear array transducer. The same elements are used to receive the echoes, back-scattered by the region of interest (ROI) under investigation. In mathematical form, each point of the i th RF image (in the following, RF i ) is obtained by a delay and sum beamforming:

RF i (x, z) = P j=1 w j RF j (τ (ξ j , x, z)) , (5.3)
where j indicates the probe element, ξ j is the position of the transducer element, (x, z) are the axial and lateral coordinates of the point, w j is the weighting function of the apodization in receive, RF j indicates the received radio-frequency signal, P is the active aperture in receive, and τ is the traveling time from the point (x, z) to the transducer element (ξ j ). The formula (5.3) is the same used for the receive beamforming in conventional modalities, where the back-scattered echoes are modulated only along the axial direction and the pulse-echo field is sensible only to axial movements. The main difference is that, for each single transmission, the entire RF image is reconstructed, while only one line is reconstructed in conventional approaches.

The transverse oscillations are suitably obtained by adopting a dynamic beamforming in receive, employing for w j the weights of a double Gaussian apodization function:

w j (ξ j , x, z) = 1 2 (e -π( ξ j -x-x 0 σ 0 ) 2 + e -π( ξ j -x+x 0 σ 0 ) 2 ), (5.4)
where x 0 is the position of the Gaussian peaks on the transducer, and σ 0 is their width at half maximum. Note that w j depends on z because,

Methods

93

according to the TO theory, x 0 has to be adapted as a function of depth (z) to obtain an oscillation with a constant transverse wavelength (λ x ):

x 0 = λ z z λ x = c f z z λ x , ( 5.5)
where λ z = c/f z is the axial wavelength, c the speed of sound and f z the burst central frequency (i.e., the axial frequency). The apodization function is defined on a convenient subset of elements and, as notable from (5.4), it is laterally centered on x. In general, if L is the width of each mainlobe of the Gaussian apodization, one beamforming line located at x is obtained by the beamforming of the 2L signals received by the elements in the neighborhood of x (i.e. L in each direction). An exemplification of the imaging technique with L = 16 is depicted in Fig. 5.1.

Velocity Estimator

Let us consider two double-oscillating images obtained at consecutive instants, RF i = r(x, z) and RF i+1 = s(x, z), and the lateral (Δ x ) and the axial (Δ z ) displacements of the scatterers inside the ROI. From the TO theory (5.2), the double-oscillating images can be approximated by the following mathematical model:

r(x, z) = cos (2πf x x) cos (2πf z z) (5.6) s(x, z) = cos (2πf x (x -d x )) cos (2πf z (z -d z )). (5.7)
where z = 1, ..., M and x = 1, ..., N are respectively the axial and the lateral coordinates in the M × N rectangular ROI, and f x = 1/λ x is the lateral frequency. The estimation of the 2-D displacement is obtained by calculating the difference of the phase between analytical images according to a block-shift approach [127]5 . First, the 2-D discrete Fourier transform (DFT) is calculated for both images

R(f x , f z ) = F 2D {r(x, z)} S(f x , f z ) = F 2D {s(x, z)} .
By calculating the inverse 2-D discrete Fourier transform in only the single quadrants I and II, the analytical images [r I (x, z), rII (x, z)] and [s I (x, z), Vector Blood Flow Imaging sII (x, z)] are thus obtained, and the phase images [φ I r (x, z), φ II r (x, z)], [φ I s (x, z), φ II s (x, z)] are then extracted. Let us define a block as a rectangular matrix of M B × N B samples, selected in the analytical image (Fig. 5.1). Thanks to the a priori signal model with the double oscillations ((5.6) and (5.7)), the phase differences between two blocks of consecutive analytical images are:

ψ I mn = 2πf x Δx mn + 2πf z Δz mn (5.8)
ψ II mn = -2πf x Δx mn + 2πf z Δz mn , (

where m = 1, ..., M B and n = 1, ..., N B are the indices of the block samples in the axial and the lateral directions, respectively. From the linear phase differences (5.8) and (5.9), the least square method is adopted to obtain the displacement estimation in the block minimizing the variance [START_REF] Basarab | Twodimensional least-squares estimation for motion tracking in ultrasound elastography[END_REF]:

Δ x = 1 4πf x M B m=1 N B n=1 ψ I mn -ψ II mn (5.10) Δ z = 1 4πf z M B m=1 N B n=1 ψ I mn + ψ II mn .
(5.11)

The velocity vector of the block is obtained multiplying the 2-D displacements for the frame-rate (F R):

v = (Δ x x + Δ z ẑ)F R = vx x + vz ẑ.
(5.12)

An extended vector map is obtained repeating the estimation on all the blocks defined into the ROI. In this paper, a lateral overlap of 90% and an axial overlap of 80% are considered. Finally, velocities are filtered by a mean temporal filter on 75 frames, and no spatial filtering is applied.

Method Validation

In order to test the performances of the method, Field II [START_REF] Jensen | FIELD: A program for simulating ultrasound systems[END_REF][START_REF] Jensen | Calculation of pressure fields from arbitrarily shaped, apodized, and excited ultrasound transducers[END_REF] has been used to simulate pulse-echo field (PEF) and blood flow images. 64 elements of the linear probe (110% bandwidth @ 8 MHz, 245 µm pitch) have been excited by N c sinusoidal cycles at f z axial frequency, with a PRF of 15 kHz. Transverse oscillations at λ x lateral wavelength have been generated apodizing the 64 channels in receive with a double Gaussian function defined on 32 elements (L = 16 elements for each mainlobe), and an RF i extended on 32 lateral lines has been obtained (Fig. 5.1).

Pulse-echo Fields

A cube of scatterers of dimension (20×4×40) mm 3 in lateral, elevation, and axial directions, respectively, has been considered. The scatterers have been randomly generated with an uniform intensity level and density of 25/mm 3 . 17 speckle images were obtained in the different setup conditions, listed in Tab. 5.1. As displayed in Fig. 5.2, for the speckle image a grid of 100 rectangular blocks (M P = 20 in the axial, N P = 5 in the lateral directions), uniformly distributed between 7.5 and 45 mm depths, has been defined. The block size has been adapted for each simulation setting a lateral size of 2λ x and an axial size of 8λ z . The ultrasound velocity c was set to 1480 m/s.

For each block, the mean lateral frequency fx and standard deviation Vector Blood Flow Imaging σ f x were estimated according to [START_REF] Pihl | Phased-array vector velocity estimation using transverse oscillations[END_REF]:

fx = f s /2 -f s /2 f sx /2 -f sx /2 f x G(f x , f z ) 2 df x df z f s /2 -f s /2 f sx /2 -f sx /2 G(f x , f z ) 2 df x df z (5.13)
and

σ f x = f s /2 -f s /2 f sx /2 -f sx /2 (f x -fx ) 2 G(f x , f z ) 2 df x df z f s /2 -f s /2 f sx /2 -f sx /2 G(f x , f z ) 2 df x df z , ( 5.14)
where

G(f x , f z ) = F 2D {g(x, z)}
is the 2-D frequency spectrum of each block g(x, z), f s is the temporal sampling frequency and f sx is the spatial sampling frequency. By using λx = 1/ fx , the performance metrics were averaged on the M P × N P blocks of the entire speckle image and were expressed in terms of the lateral wavelength (the bias B λ x and the standard deviation σλ x ):

B λ x = 1 λ x M P m=1 N P n=1 λx mn -λ x , σλ x = 1 M P N P M P m=1 N P n=1 σ λ x , ( 5.15)
where λx mn and λ x are respectively the block and the designed (theoretical) wavelengths. The performance for the axial frequency (B f z ± σf z ) were similarly calculated. Note that the values for the lateral wavelength, here expressed in mm, correspond to [START_REF] Rayleigh | On the light from the sky, its polarization and colour[END_REF][START_REF] Ramalli | Development of Novel Ultrasound Techniques for Imaging and Elastography : From Simulation to Real-Time Implementation[END_REF][START_REF] Szabo | Diagnostic Ultrasound Imaging: Inside Out[END_REF][START_REF] Wells | Ultrasonic imaging of the human body[END_REF][START_REF] Van Veen | Beamforming: a versatile approach to spatial filtering[END_REF] samples.

Pulse-echo

Blood Flow Simulations

The scatterers move inside a vessel phantom located at z 0 mm depth, inclined of ϑ with respect to the beam direction, with a 4 mm internal radius (Fig. 5.1). The scatterers have the same density as the PEF simulations (25/mm 3 ), and those inside the vessel move with a parabolic steady flow, with v 0 velocity peak. The vessel walls has 20 dB higher intensity than the blood scatterers. Double-oscillating-field images of the vessel phantom have been obtained by Field II adopting 15 kHz of frame-rate. To remove the back-scattered contributions of the vessel boards, a stationary echo-canceling filter [START_REF] Udesen | Investigation of transverse oscillation method[END_REF] has been used on 20 consecutive radio-frequency images. The lateral interline distance (δx) was 245 µm (i.e. the probe pitch, δx = 1/f sx ), and the axial pixel size (δz = 14.8 µm) was derived from the speed of sound (1480 m/s) and the sampling frequency in receive (f s = 50 MHz). Finally, the phase-based algorithm has been employed adopting a block size of M B samples in the axial and N B samples in the lateral directions. Simulations have been carried out to analyze the blood flow measurement performance by adopting a parametric study. Each pa-This solution realizes an impedance matching between water and Rilsan ® , with the purpose of limiting the intensity of the plane-wave reflections at the interface 6 . The experiment was executed in reference setup conditions (bolded in Tab. 5.2), except that the beam-to-flow angle was 80 • and, to avoid turbulence, the maximum velocity was set to 10 cm/s. The method beamforming has been accomplished by an off-line processing. For this experiment, a bigger size for the axial block was chosen for S/N compensation, and the block sizes were changed according to Tab. 5.3. The performances metrics are the same as simulations presented in the previous Section and used for the frequency-domain estimator ( §4.3.3). Two in-vivo investigations were conducted for the evaluation of the blood flow in the right common carotid artery (CCA) and in the left brachial artery of a 36-year-old healthy volunteer. The ULA-OP was programmed by the same configuration, only the frame-rate was increased to 7 kHz with the purpose of avoiding aliasing during the systolic peak.

Experimental Setup

Results

Pulse-echo Field Results

The method capability to obtain the designed oscillations has been assessed by the analysis of the pulse-echo fields (PEFs). The study is obtained according to a parametric study, in which each parameter has been singularly varied in the range of values listed in Tab. 5.1.

Blood Flow Simulations Results

Simulations of the blood flow have been obtained with the parametric study indicated in Tab. 5.2. Fig. 5.7(a) depicts the velocity profiles along 5 of the 32 lateral lines used for the performance metrics estimation ((4.13) and (4.14)) at RS conditions, and the averaged profile is shown in Fig. 5.7(b). The parameters have been changed around the RS setup, and the results are showed as following: the no. of cycles (N c ), in Fig. 5.8; the axial frequency (f z ), in Fig. 5.9; the lateral wavelength (λ x ), in Fig. 5.10; the beam-to-flow angle (ϑ), in Fig. 5.11; the barycentre depth (z 0 ), in Fig. 5.12; the peak velocity (v 0 ), in Fig. 5.13; the axial (M B δ z ) and the lateral (N B δ x ) block sizes, respectively in Figs. 5.14 and 5.15. For each case, the trend of the performance metrics and the velocity profiles of particular configurations are shown. In Fig. 5.16, a snapshot of the blood flow simulation at RS conditions is also shown.

Experimental Results

A preliminary in-vitro experiment was conducted for the evaluation of the performance in more realistic S/N conditions. The images were processed by the velocity estimator changing the block size, as indicated by Tab. 5.3, and results are described by Fig. 5.17 and Fig. 5.18, where a snapshot is shown. The blood dynamics into the vessels are described during 1 second, to cover an entire cardiac cycles.

Discussion

In this Chapter, a new method for 2-D vector blood flow imaging, based on plane waves and transverse oscillations, has been presented. The method exploits the plane-wave imaging to acquire images at high frame-rate. The transverse oscillations theory is used to obtain the double-oscillating field in a broad region. A 2-D estimator, based on the phase shift of the analytic images [START_REF] Basarab | Phase-based block matching applied to motion estimation with unconventional beamforming strategies[END_REF], has been adapted to process the RF signals following a block approach, in order to obtain 2-D vector maps of the blood flow. The results of the pulse-echo field analysis (Figs. 5.4-5.6) and of the blood flow simulations (Figs. 5.8-5.10) obtained by changing the PW-TO parameters (the no. of cycles, the axial frequency and the lateral wavelength) are in agreement, and confirms the preliminary studies [START_REF] Lenge | A new method for 2D-vector blood flow imaging based on unconventional beamforming techniques[END_REF].

The distribution of the PEF bias (Fig. 5.3) is almost uniform for the axial frequency, while the lateral wavelength presents a bias which is higher than the 10% for depths lower than 10 mm and higher than 40 mm. In the former case, the reason of such performance is that, at so low depths, the wave front is not planar and the Fraunhofer assumption (that is fundamental for the TO theory [START_REF] Jensen | A new method for estimation of velocity vectors[END_REF]) is not valid. In the latter case, the bias is due to the cut of the apodization function, that occurs at the highest depths [beyond 35 mm, see Fig. 5.4(b)]. The effect is the reduction of the distance between peaks (2x 0 , see Fig. 5.1), and as a consequence the lateral wavelength increases, according to (5.5). Indeed, the lateral velocity profile presents a deformed shape due to the S/N reduction [Fig. 5.9(b) and Fig.

5.12(b)].

The decrease of the central frequency (f z ) is also confirmed by Fig. is adequate, for the phase-based estimator, to assess the blood flow along each of the estimation lines (also the more external in the ROI, see Fig. 5.7). This is also confirmed by the frame snapshot (Fig. 5.16), obtained at RS conditions. Fig. 5.11 shows the performance of the blood flow simulations obtained varying the beam-to-flow angle (ϑ). Good results are obtained for the axial velocity, that presents errors lower than 15% up to 40 • angles. In the lateral velocity, biases and standard deviations remains lower than 15% up to 60 • , but beyond this limit the trend rapidly increases. This is clearly visible in Figs. [START_REF] Ramalli | Development of Novel Ultrasound Techniques for Imaging and Elastography : From Simulation to Real-Time Implementation[END_REF].11(g) and (h), where no parabolic lateral velocity profile can be seen, especially at higher depths. These results are obtained also in other blood flow methods. In [START_REF] Udesen | Investigation of transverse oscillation method[END_REF], the investigations on the TO methods are provided by using wider main-lobes of the TO peaks, focused beams in transmission and 32 pulses per velocity estimate, but the performance present almost the similar values. Also the speckle tracking method, proposed by Bohs et al. [START_REF] Bohs | Speckle tracking for multi-dimensional flow estimation[END_REF] that used focusing imaging, and Udesen et al. approach [START_REF] Udesen | High frame-rate blood vector velocity imaging using plane waves: Simulations and preliminary experiments[END_REF], based on plane-wave imaging with Barker coding, are affected by this problem, and the same is for the frequency-domain estimator Vector Blood Flow Imaging proposed in [START_REF] Lenge | High-frame-rate 2-D vector blood flow imaging in the frequency domain[END_REF]. As explained in [START_REF] Bohs | Speckle tracking for multi-dimensional flow estimation[END_REF], the drop of performance with higher beam-to-flow angles in not due to the transmission strategy, but it is caused by the huge mismatch in the spatial frequencies between the lateral and the axial dimensions, that results in the decorrelation of the speckle patterns, especially for the high velocity gradients and near to the walls.

If lower velocities have to be scanned (Fig. 5.13), the method is enough sensible to provide good results. Exploiting the machine PRF is thus possible to introduce compounding techniques or inter-leaved imaging TX-RX strategies or modalities.

Standard deviations of simulations are not affected by varying the block size [Figs. 5.14 and 5.15]. The axial biases slightly increase using bigger blocks, while the lateral biases present a minimum error for an axial block size of 0.95 mm (64 samples) and a lateral block size of 2.94 mm (12 samples).

The simulation study aimed at defining the most suitable setup for experiments. For in-vitro tests, the axial block size has been incremented with respect to simulations (64 samples instead of 32), in order to increase the signal energy used for the estimation in the block. This choice has been endorsed by results of Figs. 5.19 and 5.20. The method provided errors lower than the 15%, in agreement with other methods [START_REF] Udesen | High frame-rate blood vector velocity imaging using plane waves: Simulations and preliminary experiments[END_REF][START_REF] Lenge | High-frame-rate 2-D vector blood flow imaging in the frequency domain[END_REF]. With respect to other works that make use of Rilsan ® pipes [START_REF] Lenge | High-frame-rate 2-D vector blood flow imaging in the frequency domain[END_REF], using the Agar matching lets to reduce the side-effects of the reflections, and better performance have been obtained especially regarding the axial velocity (biases lower than 8%). The in-vitro results 5.17 have supported further in-vivo investigations.

The method was capable to clearly describe the flow inside big arteries, such as the CCA of Fig. 5.19, as well as in restricted arteries like the brachial artery 5.20, that shows a diameter of about 2.5 mm. Such results demonstrate that the method is sufficiently sensible for the assessment of flows in stenosed arteries, that are associated to pathological conditions. Vector Blood Flow Imaging Comparison with the frequency-domain approach. The same capabilities were demonstrated by the frequency-domain algorithm (FD), which was capable of describing fast flows inside both big arteries (like CCAs) and vessels with reduced lumina [like brachial and radial arteries, see also Fig. 4.11(b), p. 74]. A more direct comparison highlights that the two methods present some analogies. Indeed, both use plane waves in TX, to insonify broad regions at high frame-rate, and both are based on a block comparison approach on a 2-D mesh (the block-shift approach, as previously named), although the shift estimations are done in different domains (the frequency-domain and the analytic domain).

In terms of accuracies, comparing the results of blood flow simulations (Figs. 4.4 and 5.11-5.15) and experiments (Figs. 4.7 and 5.17), both methods present similar performance for the lateral velocities. The TO method obtains better performance for the axial velocity, especially in experiments, which can be due to the Agar-based phantom introduced to reduce reflections. However, for the FD method employs 5 cycles in the TX burst, while in the TO methods 8 cycles were adopted from other previous studies [START_REF] Udesen | Investigation of transverse oscillation method[END_REF]. The difference between the two methods can be ascribed to their different maturity. The transverse-oscillation theory and the phase-based estimator were more deeply investigated and validated, while the frequency-domain method needs further improvements and optimizations. About the computational loads, the two methods are very similar considering the number of calculus needed for the vector estimation in a block. Under the same hypothesis of the FD (p. 4.5), that sums and multiplications have similar complexity, and assuming that the computational load of the fft function is N log N , the final complexity is:

N ∼ = 3N d N l (log(N d N l ) + 6M B N B ), (5.16)
with N d the number of depths and N l the lateral lines. If a block size of M B × N B = (64× 8) samples is considered for both, the computational load is very similar in terms of number of operations, with the same performance discussed above. Another term of comparison is the number of beamformers that each method needs to be implemented. To reconstruct N l lines with both the RF and the demodulated (IQ) images, the FD method needs N l beamformers Vector Blood Flow Imaging and the same RF images used by the estimator can be demodulated. The double-oscillating images used in the TO estimator are not suitable to obtain a B-mode representation, and the TO scheme based on the apodization in receive hinders the efficient compounding of tilted plane waves (ultrafast imaging), unless different techniques to obtain the transverse oscillations are adopted.

In future works, the transverse oscillations can be generated by filtering techniques [START_REF] Liebgott | Fourier domain beamforming for transverseoscillations[END_REF]. This novel approach has been adopted by Salles et al. [START_REF] Salles | Plane wave transverse oscillation (PWTO): An ultra-fast transverse oscillation imaging mode performed in the fourier domain for 2D motion estimation of the carotid artery[END_REF], for the evaluation of the 2-D wall motion in the CCA. This method introduces the transverse oscillations avoiding the use of apodization in receive. Moreover, the lateral wavelength can be decided after the image reconstruction and also ultrafast investigations can be included in the method, on the condition that aliasing limits are respected. In perspective, both methods can be merged to characterize the dynamics of both the blood and the vessel walls [START_REF] Cinthio | Longitudinal movements and resulting shear strain of the arterial wall[END_REF]. Finally, more advanced anti-clutter strategies, higher number of elements in TX/RX, and the more efficient computing techniques based on parallel computing are other possible steps to improve the method.

Conclusion

In this Chapter, a new method for US blood flow imaging, based on the plane-wave transmissions and the transverse oscillation theory, is presented. The method performs double-oscillating images that are suitably exploited by a phase-based estimator, to achieve 2-D vector maps of simulated blood flows. The performance of the method are described by simulated pulseecho fields and blood flow images, obtained by changing the specific parameters of the US system and the blood flow vessel. The results outline that the PW-TO method is capable of producing the designed oscillations in the RF pulse-echo fields with accuracy and precisions lower respectively than the 4% and the 6% in standard conditions. The phase-based algorithm provides 2-D vector maps with errors lower than 10% in most of conditions. Preliminary experiments conducted in-vitro and in-vivo confirmed the simulations results and introduce the method among the possible clinical solutions for the non-invasive diagnose of cardiovascular diseases. The v x in colored in red, v z in blue. The B-mode was obtained by demodulating the directly TO data. The movement traced at the lowest depths, outside the pipe, is due to the relative movement between the Agar and the Rilsan ® pipe, that were not reciprocally fixed.

Conclusion 113

Introduction

Many efforts have been made to improve the quality of ultrasonic images, and to extend the fields of application. In ultrafast imaging, the RF image is obtained through a parallel processing in receive, i.e. by coherently adding the echoes coming from each same sample volume in the region of interest (ROI) (equation (1.12)). It was already demonstrated by Montaldo et al. that ultrafast imaging is a very efficient way to perform ultrasound imaging (see [START_REF] Montaldo | Coherent plane-wave compounding for very high frame rate ultrasonography and transient elastography[END_REF] and p. 16). Another approach addressed to increase the image quality was pursued by Jensen's group by developing the synthetic aperture (SA) method [START_REF] Misaridis | Use of modulated excitation signals in medical ultrasound. part III: High frame rate imaging[END_REF]. As introduced in §1.1.4, both SA and ultrafast methods are too time-consuming to be processed in real-time on current commercial US systems, and high levels of computational power are required from the new architectures which aim at hosting that kind of algorithms.

New advanced US systems (e.g., Aixplorer ® [Supersonic Imagine, Aix en Provence, France], SARUS [START_REF] Jensen | SARUS: A synthetic aperture realtime ultrasound system[END_REF] and ULA-OP [START_REF] Boni | A reconfigurable and programmable FPGAbased system for nonstandard ultrasound methods[END_REF]) have been developed with different basic principles and several key applications. SARUS (Synthetic Aperture Real-time Ultrasound System, [START_REF] Jensen | SARUS: A synthetic aperture realtime ultrasound system[END_REF]) is a research system for acquiring and processing synthetic aperture data. Up to 1024 channels are used to transmit custom excitation signals and to acquire the echosignals from the same number of transducer elements. Sampled data are stored in 2 GB RAM and processed by five FPGAs, and both raw and beamformed data can be stored in the system for further storage, processing, and parallel transferring by using the 1-Gbit Ethernet (3.2 GB/s theoretical rate).

ULA-OP (ULtrasound Advanced Open Platform, §1.3) has been used in this PhD project. It is a fully programmable research platform, designed by field-programmable gate arrays (FPGAs) and one digital signal processor (DSP), and capable, among other features, of the acquisition and the storing of large quantities of raw RF data (up to 36 GB, sampled at 50 MSPS, acquired from 64 elements) that can be efficiently transferred to PC (USB 3.0, up to 180 GB/s) for advanced post-processing procedures [ §1.3].

Other US systems (such as the Aixplorer ® [START_REF] Bercoff | Ultrafast ultrasound imaging[END_REF] and, more recently, the Prodigy ® [S-Sharp Corporation, New Taipei City, Taipei]) are based on graphics processing unit (GPU) architectures. Current GPU boards are equipped by thousands of cores, for computing performance up to 1.7 Tflop/s (2014 Nvidia GeForce technology, Nvidia Corporation, Santa Clara, CA) and data transfer on high-speed buses at ∼8 GB/s (PCIe 3.0 x16 protocol).

In blood flow studies, such advanced systems and imaging modalities may be exploited to overcome the limitations of the conventional Doppler in triplex modalities, which are the most recommended to clinicians for the extra-cranial circulation assessment [START_REF] Thrush | Peripheral Vascular Ultrasound: How, Why, And When[END_REF]. The speckle-tracking estimators for the 2-D vector investigations in extended ROIs, so far proposed by Bohs et al. [START_REF] Bohs | Speckle tracking for multi-dimensional flow estimation[END_REF], were adapted to plane-wave imaging by Udesen et al. [START_REF] Udesen | High frame-rate blood vector velocity imaging using plane waves: Simulations and preliminary experiments[END_REF], in order to obtain high frame-rates. Despite requiring heavy computational loads, the suggestive in-vivo images obtained by this approach [START_REF] Hansen | Fast blood vector velocity imaging using ultrasound: In-vivo examples of complex blood flow in the vascular system[END_REF] have inspired the development of more advanced techniques for blood flow imaging. The basic principle proposed by Montaldo et al. (called ultrafast imaging, [START_REF] Montaldo | Coherent plane-wave compounding for very high frame rate ultrasonography and transient elastography[END_REF]) was later applied in improved solutions (in terms of temporal resolution and sensitivity) for 1-D Doppler investigations of blood flow in wide ROIs [START_REF] Bercoff | Ultrafast compound Doppler imaging: providing full blood flow characterization[END_REF], and lately adapted to pioneering applications for cerebral blood flow imaging [START_REF] Mace | Functional ultrasound imaging of the brain: theory and basic principles[END_REF]. Configurations of multi-angle plane waves have also been proposed for flow velocity estimation [START_REF] Ricci | Real-time vector velocity assessment through multigate Doppler and plane waves[END_REF], as well as for tissue motion estimation [START_REF] Ekroll | Simultaneous quantification of flow and tissue velocities based on multi-angle plane wave imaging[END_REF].

More recently, plane-wave imaging has been coupled to an efficient estimator working in the frequency domain (FD) [START_REF] Lenge | High-frame-rate 2-D vector blood flow imaging in the frequency domain[END_REF]. This novel approach has been developed and validated, and a complete description is presented in Chapter 4. In summary, the estimator provides extended 2-D vector maps of blood flow from high frame-rate RF images, obtained by plane-wave imaging and coherent reconstruction beamforming (no compounding). The 2-D vector map is obtained by splitting the entire ROI into a 2-D mesh of rectangular blocks, and by calculating the velocity comparing the phase shift between relative blocks directly in the frequency domain (i.e., performing the 2D-DFT of each block and calculating its shift of the phase).

The peculiarity of this method is that it provides extended ROI investigations (like CFI) at full-PRF frame-rate (like in spectral Doppler, p. [START_REF] Coppini | Bioimmagini[END_REF]). Moreover, it solves the limitation common to both CFI and SD, i.e. the estimation in the axial dimension only.

A disadvantage of the method is the heavy computational load, which prevents a possible real-time implementation. With respect to the classic 2-D cross-correlation algorithm [START_REF] Udesen | High frame-rate blood vector velocity imaging using plane waves: Simulations and preliminary experiments[END_REF], the FD reaches similar performance in terms of accuracy, while slashing the computational load more than 100 times. Also with this drastic reduction of complexity, the implementation of the method on CPU remains too time-consuming to allow a profitable clinical application. Since the FD method is based on a block approach, many blocks can be processed in parallel, and the entire algorithm can be parted in several sections according to a pipeline execution. Such features are ideal for the implementation on a GPU board.

In this Chapter, the main traits of GPU boards are reviewed. Then, a GPU implementation of the HFR-FD method is depicted. Further optimizations and speed-up results in terms of computation time are discussed. Finally clinical in-vivo investigation and future perspective are delineated. Frequency-Domain

Materials and Methods

The Graphics Processing Unit

The GPU is an advanced electronic co-processor specialized in executing intense data-parallel computations, and featured by thousands of cores. It can be seen as an ensemble of computing units, each equipped by multi-core processors, which share a local memory and a common pipeline. Compared to that of a CPU, the GPU architecture has a much higher number of cores (Fig. 6.1), while the control unit is very reduced, and the branch prediction 1 is limited with respect to CPUs. The GPU threads are lighter, and need a low level of context-switch. In a few words, the GPU devotes more transistors to Data Processing, and less to the control unit [START_REF] Nvidia | CUDA C Programming Guide[END_REF]. The chip is hosted on a GPU board (the device), connected to the motherboard of the PC (the host), by means of a PCIe interface. The board is equipped by a high-speed GDDR RAM 2 , accessible by all of the cores, and by other types of memory (shared, constant, texture, cache registers) organized in hierarchical levels with specific properties and functions.

The term GPU was coined by Nvidia in the August of 1999, when the Nvidia GeForce 256 was marketed as the world's first 'GPU', or Graphics Processing Unit, a single-chip processor with integrated transform, lighting, triangle setup/clipping, and rendering engines that is capable of processing a minimum of 10 million polygons per second. In 2002, ATI Technologies 1 In computer architecture, a branch predictor is a digital circuit that tries to guess which way a branch (e.g. an if-then-else structure) will go before this is known for sure. The purpose of the branch predictor is to improve the flow in the instruction pipeline. Branch predictors play a critical role in achieving high effective performance in many modern pipelined microprocessor architectures such as x86 [136]. 2 The Random Access Memory (RAM) is a form of computer data storage and allows data items to be read and written in roughly the same amount of time regardless of the order in which data items are accessed [138]. In the GPU, the RAM is shared by all the cores and writing/reading events can be concurrent to processing operations. GDDR means Graphics Double Data Rate.

released the Radeon 9700 defined as visual processing unit or VPU. The two competitors have designed GPU devices with very different basic principles for the access to the memory subsystem [START_REF] Zhang | Architecture comparisons between nvidia and ATI GPUs: Computation parallelism and data communications[END_REF].

Modern GPU are suitable not only for 3-D computer graphics, but also in some scientific fields, such as 3-D computational fluid dynamics, bioinformatics, computational modeling, computer vision, or financial modeling, image processing, and, in particular, medical image processing [START_REF] Eklund | Medical image processing on the GPU -past, present and future[END_REF]. Other types of applications, such as optimized problems or data compression, need a powerful control unit to adopt instruction prefetch techniques, that reduce the wait states and run better on CPU.

General-purpose Computing on GPU

GPUs were conceived for the fast rendering of 3-D computer graphics, especially in computer gaming industry. Recent developments in GPUs include the support for programmable shaders 3 , which can manipulate vertexes and textures with many of the operations supported by CPUs (oversampling and interpolation techniques to reduce aliasing), and very highprecision color spaces. Since most of these computations involve matrix and vectors operations, engineers and scientists have increasingly investigated the use of GPUs for non-graphical calculations [142]. Using of GPU for computing tasks has been defined "general-purpose computing on GPU" (GP-GPU).

To exploit these features, the programmer should implement a Turin machine 4 equivalent to the GPU architecture. The GPU is re-programmable using the graphic libraries (i.e., DirectX and OpenGL) that support the shader model. In other words, by means of shaders, the programmer access to the memory on the GPU device, data are then submitted to the "rendering" and later stored into the memory accessible in output. Many projects have been activated to allow CPU-based applications to directly access the 3 In the field of computer graphics, a shader is a computer program that is used to do shading: the production of appropriate levels of color within an image, or, in the modern era, also to produce special effects or do video post-processing. Shaders calculate rendering effects on graphics hardware with a high degree of flexibility. Most shaders are coded for a graphics processing unit (GPU), though this is not a strict requirement. Shading languages are usually used to program the programmable GPU rendering pipeline, which has mostly superseded the fixed-function pipeline that allowed only common geometry transformation and pixel-shading functions; with shaders, customized effects can be used. The position, hue, saturation, brightness, and contrast of all pixels, vertexes, or textures used to construct a final image can be altered on the fly, using algorithms defined in the shader, and can be modified by external variables or textures introduced by the program calling the shader (vertex and pixel shaders) [141] 4 A Turing machine is an abstract model that defines a machine able to execute algorithms. Frequency-Domain resources of a GPU for more general purpose computing without the limitations of using a graphics API. Among the most popular projects there are CUDA (Compute Unified Device Architecture) and OpenCL (Open Computing Language), that allow specified kernel functions defined in a normal C program to run on the GPU's stream processors. In this way, a C program can execute serial code on the CPU and operate on large matrices in parallel on the GPUs.

The CUDA platform is the earliest widely adopted programming model for GPU computing. It was introduced by Nvidia Corporation in 2006 with the purpose of solving many complex computational problems in a more efficient way than on a CPU. CUDA was born with a software environment that allows developers to use high-level programming languages (e.g, C/C++ and FORTRAN). It consists in a hardware and software coprocessing architecture for parallel computing based on a scalable programming model [START_REF] Nickolls | The GPU computing era[END_REF]. The GPU works as a co-processor for the CPU, with its own RAM and capable to execute thousands of cores in parallel [START_REF] Nvidia | CUDA C Programming Guide[END_REF]. While the serial programming is demanded to the CPU, the parts of a program that are data-parallel and computational-intensive are demanded to the GPU. In such model, a kernel function is defined and identically executed by all the threads on the device.

In 2009, the Khronos Group TM defined the OpenCL, an open standard for parallel programming of heterogeneous systems. A multitude of architectures including GPU, CPU, DSP and FPGA are compatible with OpenCL, and the purpose is to increase the portability of the code towards many other devices. The slogan "program once, run everywhere" has been received by the most important manufacturer of programmable devices (Intel, AMD, Nvidia, Altera, Texas Instruments and ARM) that made available specific software development kits (SDKs) for the OpenCL programming of their products.

Among the possible solutions, the CUDA programming language makes easier programming the GPU for general-purpose calculations, and it was chosen for the implementation of the frequency-domain algorithm (FD). First, the code portability given by OpenCL is not directly translated into performance portability, meaning that a code developed for a specific architecture can reach better performance with respect to a more portable solution developed in OpenCL standard and run on the aforementioned architecture. This is reasonable, because usually top performance are reached with assembly instructions, and the SDKs for OpenCL have to abstract from the fastest assembly instructions used in the original SDK. The issue must consider also the level of optimization necessary to the application, and the time necessary to obtain such optimization, that depends on programmer's experiences and skills. In general, the gap in performance 6.2 Materials and Methods 123 between CUDA and OpenCL was initially wide, but currently it is diminishing as OpenCL compilers improve [START_REF] Eklund | Medical image processing on the GPU -past, present and future[END_REF]. Therefore, sooner or later, it is likely that the tradeoff becomes convenient, considering the continuous releases of renewed architectures and the time needed by code porting.

The OpenCL framework is broadly supported and it is the most widely used by developers for GPGPU platform. However, CUDA has currently some advantages. Nvidia continuously supplies optimized libraries (e.g., CUBLAS, CUFFT, CURAND, etc.) for many different purposes, and Mathworks encapsulated many commonly used algorithms into GPU accelerated versions collected in the parallel computing toolbox of Matlab ® . These features allow researchers the fast prototyping of novel algorithms, without going into details of GPU programming [START_REF] So | Medical ultrasound imaging: To GPU or not to GPU?[END_REF][START_REF] Eklund | Medical image processing on the GPU -past, present and future[END_REF].

Considering the purpose of this study, that consists in a preliminary test of the GPU capability to speed-up the FD algorithm, a powerful consumer Nvidia GPU board and the CUDA model have been adopted to run a derated Matlab version of the FD algorithm. In the following Section, the proposed solution is outlined in details.

GPU Implementation of the Frequency-Domain Algorithm

Acquisition Procedure

The data acquisition was based on the ULA-OP system ( §1.3), coupled with the LA533 linear probe (110% bandwidth centered at 8 MHz) manufactured by Esaote S.p.A. (Florence, Italy). Plane waves were generated at 7 kHz PRF, by simultaneously exciting 64 elements of the linear probe with 5 sinusoidal cycles at 6 MHz frequency.

The external and internal carotid arteries, the subclavian artery and the renal vein of 30-year-old healthy volunteer were investigated by an expert sonographer. After 5 minutes of rest, the volunteer assumed a supine position. The sonographer performed a preliminary scan with ULA-OP in duplex mode. Once the ROI was identified, and maintaining the probe position unaltered, the ULA-OP transmission was switched to plane wave mode. In receive, the 64 RF signals (50 MHz frequency with a 12-bit dynamic) were A/D converted and stored in the high capacity memory board. 2400 RF samples were taken for each of the 64 channels, corresponding to a depth of 35 mm, to cover a ROI of about (35 x 15.8) mm 2 . For 3 s long acquisitions at 7 kHz, about 5 GB data were saved on the memory board. Finally, the data were transferred from the memory board to the PC, in approximately 30 s by the USB 3.0 connection (about 150 MB/s transfer rate, §1.3). Frequency-Domain

Data Processing

The processing started with the coherent reconstruction of the RF signals. The beamforming was executed by a non-optimized Matlab ® script, that employs about 1 hour of processing to reconstruct 21000 RF images.

The FD method was executed on a PC supported by the Intel ® i7-3770 Processor (3.4 GHz) and the GPU board Nvidia GeForce GTX TITAN Black. The latter is a top-level consumer graphic card featured by 15 multi-processors, with 192 processor cores each (for a total of 2880 cores), and 6 GB of GDDR5 memory with very fast internal bus (up to 144 GB/s bandwidth). The GPU board was installed on the PC motherboard (PCIe 3.0, up to 16 GB/s bandwidth). The processing performed on each couple of consecutive RF images mainly consists of FIR filtering, 3-point DFT, phase extraction and displacement estimation. Such operations are demanded to the GPU device according to a pipeline architecture. A Matlab ® version of the algorithm cyclically transfers an ensemble (N = 1000) of RF images to the GPU device, while the processing runs in parallel on the thousands of cores. Such hiding technique was adopted to reduce the long-latency of data transfers through the PCIe bus. However, the parallel computing toolbox of Matlab ® does not allow the perfect control on the workflow.

The method was also optimized for the CPU execution, in order to make a comparison between the two implementations.

For both implementations, a derated version of the FD method presented in Chapter 4 was adopted, in order to reduce the computation time without loosing in performance. A block size of 64×8 samples was adopted, and the vector map was processed on a mesh of 240 × 56 blocks for each RF image. The map was obtained adopting an overlap of 61% and 66% in the axial and lateral directions, respectively. Such values consisted in The scripts were finally integrated in the SIMAG software ( §1.4).

Results

2-D Vector Maps

The flows in the renal vein, in the external and internal carotid arteries and in the subclavian artery are reported in Figs. 6.2, 6.3-6.4, 6.5-6.6, respectively. For the arteries, the four screen-shots represent the flow at different instants of the cardiac cycle, in particular the end diastolic, systolic, end systolic and start diastolic. Such instants are represented by a red point superimposed on a schematic electrocardiogram (ECG). The timing between the blood-flow event and the ECG (and also other cardiovascular events) can be illustrated by the well-known Wiggers diagram (Fig. 6.7). It can be appreciated that the method works with both surface and deep vessels. For the renal vein shown in Fig. 6.2, the two frames show different flows that are caused by asynchronous contractions of the skeletal muscles.

Computational Time Results

The computation times needed to process the RF data by the CPU and the GPU board, respectively, have been evaluated. For each second of acquired RF-data, 71.8 s were necessary for the CPU execution, while the GPU execution needed 11.5 s, with a resulting 6.2x speed-up factor. In the estimation, the data transfer from the RAM to the GDRAM was taken into account. Frequency-Domain which are normal values for healthy conditions [START_REF] Grant | Carotid artery stenosis: Gray-scale and Doppler US diagnosis[END_REF]. The brain is located on the left. the arterial flow and the heart sound (adapted from [START_REF]Wiggers diagram[END_REF]).

Results

Discussion

In this Chapter, a GPU implementation of the FD algorithm was evaluated. We adopted a derated version of the FD method, implemented on a top-level consumer GPU board by Matlab functions. The method was applied for the investigation of different vessels of the body and produced realistic descriptions of the blood flow also in reduced lumen vessels. Encouraging results were obtained in terms of speed-up, and an optimized C version of the algorithm is expected to produce the further reduction of the computational time. Considering that the GPU technology is interested by a rapid increase of performance, aiming at answering the constant demand for greater realism in computer games [START_REF] Nickolls | The GPU computing era[END_REF], even better results are expected with the next generation boards.

In this context, the beamforming complexity was not considered. However, ultrafast and synthetic aperture imaging have already been implemented in efficient solutions for GPU computing [START_REF] Yiu | GPU-based beamformer: Fast realization of plane wave compounding and synthetic aperture imaging[END_REF]. These methods increase the frame-rate, the focusing power or the image contrast, by shortening the data-acquisition period and, consequently, increasing the processing time for the image formation.

As suggested in [START_REF] Bercoff | Ultrafast compound Doppler imaging: providing full blood flow characterization[END_REF], compounding of tilted plane waves could be proficiently exploited to solve the limitation of the decreased quality due to plane wave imaging, and to allow empowered triplex modalities. GPU-based solu-

Conclusion and Perspectives 131

tions has been already proposed [START_REF] Yiu | GPU-based beamformer: Fast realization of plane wave compounding and synthetic aperture imaging[END_REF] for high frame-rate imaging. However, it is well known that the data transfer is the real performance bottleneck in GPU systems [START_REF] Eklund | Medical image processing on the GPU -past, present and future[END_REF][START_REF] Geist | Paving the roadmap to exascale[END_REF]. Hence, the choice between GPU-based and more hardware-oriented, or even hybrid, architectures, for such beamforming algorithms [START_REF] So | Medical ultrasound imaging: To GPU or not to GPU?[END_REF] remains an open question. Especially taking into account the increasing number of elements, involved in potential 3-D ultrasound investigation allowed by 2-D probes, the question remains challenging.

Conclusion and Perspectives

This Chapter concerns the implementation of the frequency-domain algorithm on a GPU board. New clinical in-vivo results are showed and a temporal evaluation of the method is provided by CPU and a GPU board on a desktop PC. The results highlight that a significant speed-up is obtainable by adopting a commercially available GPU board, demonstrating that the method is suitable for a GPU implementation. Work is in progress to develop a C version of the code. Enhanced in-vivo clinical investigation should be possible, in which slow-motion movies of HFR blood flow images are made available within a short time after the data acquisition.

Conclusions

Doppler methods are widely used for the analysis of blood flow patterns, to demonstrate or exclude the presence of cardiovascular diseases. However, these methods are still affected by some limitations preventing the achievement of accuracy and reproducibility levels needed for reliable diagnosis. The main limitation is that the estimation is performed in the axial direction only. The estimation of the velocity magnitude needs the manual adjustment of the flow direction by the clinician, who typically assumes it is parallel to the vessel walls. In this way, he introduces an approximation that is effective only for rectilinear sections of big arteries, and risks to be inappropriate for stenosed arteries which are the most important to be diagnosed. Furthermore, the modern clinical standards suggest the use of the triplex modalities, which allow the velocity assessment only in one sample volume selected by the clinician, or over wide regions but with low temporal resolutions.

This PhD project was focused on the investigation of Doppler methods addressed to overcome the above limitations, and two major goals have been achieved. First, the maturity of vector methods and their suitability to increase the reproducibility of Doppler exams were demonstrated. Second, novel approaches to extend the analysis to wide regions of interest at high frame rate were proposed.

The common factor of these novel approaches is the transmission of plane waves, which represents a breakthrough modality for ultrasound imaging in toto. Plane-wave imaging allows the investigation of the region of interests at full PRF. The related loss of image quality is conveniently solved by the ultrafast imaging techniques, which efficiently recombine the low quality images obtained from multiple tilted plane-waves in order to obtain a high resolute image. The plane waves have paved the way to many different applications. In this PhD work, plane-wave transmissions have been introduced to extend the 2-D vector velocity estimate to a wide region of interest. Such bi-dimensional vector maps reveal the hemodynamic details of arterial pulsatile flow with unprecedented quality, enabling more effective diagnosis of cardiovascular diseases.

Plane-wave imaging makes the conventional estimators (spectral analysis, auto-correlation, cross-correlation) unsatisfactory for processing at high frame-rate, and introduces additional challenges for post-processing algorithms. The estimators adopted or introduced in this PhD work reduce the computational load with respect to other conventional approaches, but their real-time implementation still requires computing capabilities that are not compatible with the transportable echographic systems so far developed. Nevertheless, the "'real-time"' feature is still among the strong points of ultrasound techniques, together with safety and non-invasiveness for the patient.

In this context, the GPU-based solutions are of increasing interest, and their capability of producing significant speed-up factors has been demonstrated also for the displacement estimator proposed in this Thesis. The GPU technology will surely have impact on the next generation of ultrasound systems. In few years, new GPU boards will be endowed of computational power sufficient to make ultrafast imaging suitable for real-time operation. They could positively affect also the development of 3-D imaging systems based on matrix arrays. The latter ones are proposing new challenges related to the control of high number of elements, and innovative architectural solutions are needed. At least for research purposes, but also for the clinical employment, hybrid systems combining hardware (e.g. FPGAs and DSPs) and software (e.g., GPU) solutions could be ideal to face the most difficult tasks emerging in this field. In a realistic perspective, the system will have to beamform the data and to display the anatomic grayscale image in real-time, with superimposed the color-coded velocity map in extended ROIs and the vector in at least one point selected by the operator. A post-processing on the beamformed data could provide the vector map all over the ROI.

  Conference on Acoustics, Speech and Signal Processing (ICASSP), Firenze (Italy), 2014, pp. 5125-5129. M. Lenge, D. Righi, S. Ricci, H. Liebgott, and P. Tortoli, "Blood velocity measurement in healthy and diseased carotid arteries by vector Doppler techniques," in 2014 IEEE International Ultrasonics Symposium (IUS), Chicago (IL, USA), 2014, pp. 345-348. S. Ricci, M. Cinthio, M. Lenge, R. Matera, J. Albinsson, and P. Tortoli, "Volume flow assessment through simultaneous B-mode and Multigate Doppler," in 2012 IEEE International Ultrasonics Symposium (IUS), Dresden (Germany), 2012, pp. 1588-1591. A. Ramalli, L. Bassi, M. Lenge, C. Palombo, K. Aizawa, and P. Tortoli, "An integrated system for the evaluation of flow mediated dilation," in 2014 IEEE International Conference on Acoustic, Speech and Signal Processing (ICASSP), Firenze (Italy), 2014, pp. 5182-5185. E. Boni, A. Cellai, A. Ramalli, M. Lenge, and S. Ricci, "Multichannel Raw-Data Acquisition for Ultrasound Research," in 2014 17th Euromicro Conference on Digital System Design (DSD), 2014, pp. 647-650.

Figure 1 . 1 :

 11 Figure 1.1: Sound spectrum of mechanical waves.

Figure 1 . 2 :

 12 Figure 1.2: Scheme of the transducer beam field: the near field, where the beam presents a quite constant width, and the far field, where the beam diverges [5].

Figure 1 . 3 :

 13 Figure 1.3: Examples of medical ultrasound imaging applications: (a) The Common Carotid Artery scanned by a Color Flow Imaging (CFI) technique (see p.

  20 for more details); (b) an obstetric ultrasonography of a fetus at 12 weeks lying down the mother's uterus: it is observable the head (the bright white circle), that faces to the left, the right hand and the left foot; (c) a CFI echocardiography: the left atrium (LA), the right ventricle (RV), the aorta (Ao), the left ventricle (LV) of the heart are visible during a mitral regurgitation (MR); (d) a transvaginal ultrasonography, that displays the embryo at 4 weeks from the conception.

. 11 )

 11 Each point (x, z) of the image is obtained by delaying the echo signals RF ( p , t) by the time τ ( p , x, z) and summing them in the direction of the array p [Fig. 1.5(b)]: s(x, z) = P p=1 RF ( p , τ ( p , x, z)).

Figure 1 . 5 :

 15 Figure 1.5: Plane-wave imaging schemes: (a) transmission and (b) coherent reconstruction in receive.

Figure 1 . 6 :

 16 Figure 1.6: Architecture of the ULA-OP data processing chain.

Figure 1 . 7 :

 17 Figure 1.7: TX/RX settings interface, by which the user can set the excitation functions, the filtering operations, the beamforming in transmission (plane-wave imaging) and in reception (transverse oscillations), and see the relative plots. In this particular case, the setup for plane-wave transverse oscillations method (Chapter 5) is adopted: by the "Signal and Demodulation" panel, 8 cycles sinusoidal bursts at 6 MHz are set to shape a plane wave in transmission, selected by the "TX/RX beamforming" panel; in reception, a dynamic focusing using the "TO apodization" is selected, and the relative plots are shown.

Figure 1 . 8 :

 18 Figure 1.8: Superficial dissection of the right side of the neck: the Common, the External and the Internal Carotid Arteries [61].

Figure 2 . 1 :

 21 Figure 2.1: Sketch of the main blocks of the signal processing platform.

Figure 2 . 2 :

 22 Figure 2.2: The SPP interface for the volume flow assessment application. On top, from left to right, the B-mode, the Profile and the Speed plots. On bottom, the control panel and the trends of the Flow (upper) and the Diameter (lower).

Figure 2 . 3 :

 23 Figure 2.3: Relative error distribution (dark bars) compared to Gaussian (red-line) fitted to relative error values.

Figure 2 . 4 :

 24 Figure 2.4: Comparison between the volume flow measurements of the proposed method and the reference values.

Figure 2 . 5 :

 25 Figure 2.5: Three overlapped measurements obtained with pulsed flow, 12 • steering angle, 420 ml/min nominal flow.

  2.7(a) reports the B-mode image of a brachial artery logarithmically compressed and represented in black and white color scale. The horizontal red lines on the image are the outputs of the Diameter class and represent the estimated position of near and

Figure 2 . 6 :

 26 Figure 2.6: The SPP interface for FMD application.

Figure 2 . 7 :

 27 Figure 2.7: Output of the main processing blocks, how they are plotted by the graphic methods of the developed classes. (a) B-mode and Diameter extraction; (b) Multigate Spectral Doppler profile; (c) denoised profile and mean frequency extraction (white line); (d) speed estimation (white line) and fitting (red line); shear rate considering the unfitted speed (white line) and the fitted one (red line).

Figure 2 . 8 :

 28 Figure 2.8: A FMD processing result. The average WSR between near and far wall (top), the estimated diameter (middle) and the blood pressure (bottom). In each case the red line corresponds to the result of a low pass filter.

Figure 3 . 1 :

 31 Figure 3.1: Typical display of the conventional spectral estimation in duplex mode. The B-mode image of the investigated region, with superimposed Doppler line, angle cursor, and SV, is shown on top. The spectrogram corresponding to the SV is displayed on the bottom and the average PSV, calculated from the spectrogram peaks, is shown in the bottom-left corner. Note that the SV here is long about 2 mm.

Figure 3 . 2 :

 32 Figure 3.2: ULA-OP user-interfaces. The same left CCA investigated as in Fig. 3.1 was here scanned by the sonographer with ULA-OP and processed according to Method 2 (a) and Method 3 (b). (a) Angle-tracking VD: once the SV is located in the ROI, the tracking (yellow) and the measuring (cyan) lines, which intersect in the SV, are superimposed on the B-mode (top-left). (b) Plane-wave VD: the velocity is estimated along the entire vertical line and shown through red or blue arrows, according to the direction; the SV selected by the operator is emphasized by the yellow ellipse. In both cases, the velocity magnitude profile in the SV is shown at the bottom, together with the average PSV and the heartbeat. Here the SV is long about 0.7 mm.

Figure 3 . 3 :

 33 Figure 3.3: Comparison between PSV values measured by Methods 1-2 (top), 1-3 (centre), 2-3 (bottom). The red crosses and black circles represent ICA and CCA measurements, respectively. The continuous and dotted lines report the mean velocity difference and the ±1.96S d range, centered around the mean.

Figure 3 . 4 :

 34 Figure 3.4: PSV values comparison between Methods: 1-2 (a), 1-3 (b), 3-2 (c). The circles and the crosses represent the CCA and the ICA measurements, respectively. The black dotted line is the theoretical trend, and the red dotted line is the interpolation curve of the measurements values.

Figure 3 . 5 :

 35 Figure 3.5: Inter -operator analysis for the three tested methods. For each method, the box represents the 25 th and 75 th percentiles, the central line is the median value, the whiskers extend among the boundaries, the cross shows an outlier.

Figure 3 . 6 :

 36 Figure 3.6: Results of the investigation of a stenosed ICA (70% occlusion) by means of: a) Spectral Doppler; b) Angle Tracking; c) Plane Wave VD methods.

1

  These works have been supported by the Italian Ministry of Education, University and Research (PRIN 2010-2011), the European Fund for Regional Development for the 2007-2013 programming period (POR FESR 2007-2013 CReO ASSO and MIMAUS Projects), by the Franco-Italian University (Galileo Project 2011-2012 n.26075WL), and by the European Union's Seventh Framework Programme (FP7/2007-2013) for the Innovative Medicine Initiative under grant agreement number IMI/115006 (the SUMMIT consortium).

Figure 4 . 1 :

 41 Figure 4.1: The setup used in simulations and experiments: in the tube, the rectangular region producing the data for an estimation block (s mn ) is highlighted.

66 2 -

 2 D Vector Blood Flow Imaging in the Frequency Domain transforms (DFTs) of s mn , which are defined as follows:

. 6 )Figure 4 . 2 :

 642 Figure 4.2: Block diagram illustrating the main steps of the frequency-domain algorithm.

Fig. 4 .

 4 3 reports the simulated and estimated axial and lateral velocities on the vessel axis. In Fig.4.5, the axial and lateral spatial velocity profiles at the checkpoints (T 1 , T 2 , T 3 in Fig.4.3) are reported. A sample frame captured at the systolic peak is shown in Fig.4.6.

Figure 4 . 3 :

 43 Figure 4.3: Pulsatile-flow simulation axial (black) and lateral (gray) velocities estimated in the centre of the vessel. T 1 , T 2 and T 3 are the time checkpoints used in Fig. 4.5 to show the corresponding spatial velocity profiles.

Figure 4 . 4 :

 44 Figure 4.4: Simulations bias and standard deviation of vz (black) and vx (gray) estimations. The varying parameters are reported in Tab. 4.1.

74 2 -

 2 D Vector Blood Flow Imaging in the Frequency Domain

Figure 4 . 5 :

 45 Figure 4.5: Instantaneous axial (a) and lateral (b) velocity profiles of the pulsatileflow simulation at the T 1 , T 2 and T 3 time checkpoints indicated in Fig. 4.3. Solid lines indicate the estimated velocities and dashed lines are the simulated ones.

Fig. 4 .

 4 11(a) shows a frame of the movie captured at systole. Another example is shown in Fig.4.11(b), obtained from a radial artery with diameter of about 2 mm insonified at a

Figure 4 . 6 :

 46 Figure 4.6: Vector maps captured at systole (checkpoint T 1 ) in the pulsatile flow simulation. Each estimated local velocity magnitude is coded according to the colorbar on the left, while only a down-sampled set of arrows is displayed.

Figure 4 . 7 :

 47 Figure 4.7: In-vitro performance metrics (bias and standard deviation) on vz (black) and vx (gray) estimations. The varying parameters are reported in Tab.

4. 2 .

 2

4. 5 77 Figure 4 . 8 :

 57748 Figure 4.8: In-vitro beam-to-flow angle trend recorded during from a sample volume located at the center of the pipe. Since the vessel was manually rotated, the angular acceleration was not constant. The estimated angle is compared with the angle of the vessel walls extracted from the B-mode images.

Figure 4 . 9 :

 49 Figure 4.9: In-vivo axial (black) and lateral (gray) velocities in a sample volume located at the centre of the CCA. T 1 , T 2 and T 3 are the time checkpoints at 100, 175, 400 ms.

Figure 4 . 10 :

 410 Figure 4.10: (a) Axial and (b) lateral spatial velocity profiles as a function of time checkpoints of Fig. 4.9. From 11.5 to 15 mm depths, the profiles of the jugular vein presents a parabolic shape. For depths between 15 and 22 mm, the profiles related to the carotid artery assumes different shapes, according to the cardiac cycle. Note that the ordinates of (a) and (b) are differently scaled.

Figure 4 . 11 :

 411 Figure 4.11: In (a), vector map in the CCA captured at the systolic peak (checkpoint T 1 ). The opposite flow (towards the heart) in the external jugular vein can be observed on the top. In (b), vector map in the radial artery captured during diastole (the hand is in the left).

80 2 -

 2 D Vector Blood Flow Imaging in the Frequency Domain

Figure 4 . 12 :

 412 Figure 4.12: In (a), vector map in a CCA captured during the M-shaped profile. The flow stands more on the external wall. The same behavior, but in the internal wall, is depicted in (b), that shows another CCA captured just before the systolic peak. In both, the brain is on the right.

4. 5 81 Figure 4 . 13 :

 581413 Figure 4.13: In-vivo B-mode images (left) with overlapped 2-D vector maps (right) of femoral vessel flow: the femoral artery is centered at about 8 mm, above the vein bifurcation. In (a), the blood map describes the artery flow in diastole, when blood flows across the opened valve towards heart; (b) describes the artery flow toward the foot in the systolic peak, when the venous flow is blocked by the closed valve.

82 2 -

 2 D Vector Blood Flow Imaging in the Frequency Domain

Figure 4 . 14 :

 414 Figure 4.14: Simulated lateral averaged velocity profiles obtained with the RS, calculated along different estimation lines. Each figure reports the parabolic profiles with different lateral block sizes (N x = N B δx = 1, 2, 3, 4, 5 mm), where the numbers 17, 24, and 32 correspond to the lateral position of the considered estimation lines (32 is the most lateral).

86 2 -

 2 D Vector Blood Flow Imaging in the Frequency Domain results in Fig. 4.4(f) and Fig. 4.7(e).

Figure 5 . 1 :

 51 Figure 5.1: The PW-TO imaging method. On top, three examples of apodization functions defined on 32 RX channels. To obtain the TO beamforming for the point A = (-3.4, 10.0) mm, colored in red, the Gaussian main-lobes (16-elements width) of the apodization (the red dash-dotted lines) are centered on the element j = 17 and are defined on the adjacent 32 RF channels (i.e. j = 1 ÷ 32). Similarly, the solid blue line corresponds to the apodization function for B = (0.0, 25.0) mm, and the green dotted line corresponds to C = (+3.4, 40.0) mm.

Figure 5 . 2 :

 52 Figure 5.2: PEF grid. In (a), the speckle image with superimposed the grid of blocks for the PEF evaluation. In (b) and (c), the RF and the TO signals are respectively shown, both referred to a vertical and a horizontal line selected in the block highlighted in red. In (d), the four quadrants (I, II, III, IV ) of the 2D-DFT of the block are reported and where the axial (f z ) and the lateral (f x = 1/λ x ) frequencies are indicated. Note that f x = 0.68 mm -1 is λ x = 1/f x = 1.47 mm.

Fig. 5 .

 5 [START_REF] Culjat | A review of tissue substitutes for ultrasound imaging[END_REF] shows the spatial distribution of the biases of the axial frequency (B f z ) and the lateral wavelength (B λ x ), measured in the blocks of Fig.5.2(a). The performances obtained varying the number of cycles (N c ), Vector Blood Flow Imaging the axial frequency (f z ) and the lateral wavelength (λ x ) are respectively reported in Figs.5.4, 5.5 and 5.6, where also some selected maps of the apodization functions are reported.

Figure 5 . 3 :

 53 Figure 5.3: PEFs biases: spatial distribution at RS of (a) the axial frequency (B f z ) and (b) the lateral wavelength (B λ x ). The averaged performance metrics are: B f z ± σ f z =(0.3 ± 1.5)% and B λ x ± σ λ x =(-3.3 ± 5.7)%.

5. 4 101 Figure 5 . 4 :

 410154 Figure 5.4: No. of cycles in transmission (N c ) effect on PEFs: (a) PEFs performance metrics (B f z ± σ f z ) and (B λ x ± σ λ x ); (b) the apodization along the central line (j = 32) obtained in RS (N c = 8). The filled markers in (a) represent the RS.

Figure 5 . 5 :

 55 Figure 5.5: Axial frequency (f z ) effect on PEFs: (a) PEFs performance metrics; (b) the apodization along the central line (j = 32) obtained with f z = 4 MHz. The filled markers in (a) represent the RS.

Figure 5 . 6 :

 56 Figure 5.6: Lateral wavelength (λ x ) effect on PEFs: (a) PEFs performance metrics; (b) The apodization along the central line (j = 32) obtained with λ x = 1.96 mm. The filled markers in (a) represent the RS.

Figs. 5 .

 5 19 and 5.20 show two snapshots of the in-vivo investigations conducted in the CCA [Fig.5.19] and in the brachial artery [Fig.5.20].

Figure 5 . 7 :

 57 Figure 5.7: (a) Lateral (gray) and axial (black) velocity profiles of the blood flow simulation (BFS) obtained at RS. Solid lines indicate the estimated velocities (v) and dashed lines are the simulated ones (v). The numbers 5, 11, 17, 23 and 29 are a subset of the 32 lateral lines where the estimation is performed. (b) The lateral and the axial velocity profiles obtained averaging all the 32 profiles, as indicated in (4.13) and (4.14). The solid lines are the estimated velocities, dashed lines are the simulated ones and dot-dashed lines are the ±3 standard deviation. The averaged performance metrics are: B v z ± σ v z =(3.5 ± 1.4)% and B v x ± σ v x =(8.5 ± 4.3)%.

  5.5(a) where, according to (5.5), the cutting of the apodization function at 25 mm depth is due to the setting of the axial frequency to low values [see Fig. 5.5(b), where f z =4 MHz].

Figure 5 . 9 :

 59 Figure 5.9: Axial frequency (f z ) effect on BFS: (a) BFS performance metrics for the axial (black) and the lateral (gray) velocity profiles; (b) velocity profiles for f z = 4 MHz. In (a) the filled markers represent the RS. In (b), solid lines are the estimated velocities, dashed lines are the simulated ones and dot-dashed lines are the ±3 standard deviation.

Figure 5 . 10 :

 510 Figure 5.10: Lateral wavelength (λ x ) effect on BFS: (a) BFS performance metrics for the axial (black) and the lateral (gray) velocity profiles; (b) velocity profiles for λ x = 1.96 mm. In (a) the filled markers represent the RS. In (b), solid lines are the estimated velocities, dashed lines are the simulated ones and dot-dashed lines are the ±3 standard deviation.

Figure 5 . 11 :

 511 Figure 5.11: Beam-to-flow angle (ϑ) effect on BFS: (a) BFS performance metrics for the axial (black) and the lateral (gray) velocity profiles; (b)-(h) velocity profiles for ϑ = [90, 80, 70, 60, 50, 40, 30] • . In (a), the filled markers represent the RS. In (b)-(h), solid lines are the estimated velocities, dashed lines are the simulated ones and dot-dashed lines are the ±3 standard deviation. The absolute error for ϑ = 90 • (b) is (0.3 ± 0.1). In titles of (b)-(h), the symbol "deg" is used for • .

Figure 5 . 12 :

 512 Figure 5.12: Barycentre depth (z 0 ) effect on BFS: (a) BFS performance metrics for the axial (black) and the lateral (gray) velocity profiles; (b) velocity profiles for z 0 = 35 mm. In (a) the filled markers represent the RS. In (b), solid lines are the estimated velocities, dashed lines are the simulated ones and dot-dashed lines are the ±3 standard deviation.

Figure 5 . 13 :

 513 Figure 5.13: Peak velocity (v 0 ) effect on BFS: (a) BFS performance metrics for the axial (black) and the lateral (gray) velocity profiles; (b) velocity profiles for v 0 = 30 cm/s. In (a) the filled markers represent the RS. In (b), solid lines are the estimated velocities, dashed lines are the simulated ones and dot-dashed lines are the ±3 standard deviation.

Figure 5 . 14 :

 514 Figure 5.14: Axial block size (M B ) effect on BFS: (a) BFS performance metrics for the axial (black) and the lateral (gray) velocity profiles; (b) velocity profiles for M B = 160 samples. In (a) the filled markers represent the RS. In (b), solid lines are the estimated velocities, dashed lines are the simulated ones and dot-dashed lines are the ±3 standard deviation.

Figure 5 . 15 :

 515 Figure 5.15: Lateral block size (N B ) effect on BFS: (a) BFS performance metrics for the axial (black) and the lateral (gray) velocity profiles; (b) velocity profiles for N B = 20 samples. In (a) the filled markers represent the RS. In (b), solid lines are the estimated velocities, dashed lines are the simulated ones and dot-dashed lines are the ±3 standard deviation.

Figure 5 . 16 :

 516 Figure 5.16: Snapshot of the simulated blood flow at RS conditions.

Figure 5 . 17 :

 517 Figure 5.17: Axial (a) and lateral (b) block sizes: performance metrics for the in-vitro experiment. The filled markers represent the RS.

Figure 5 . 18 :

 518 Figure 5.18: Snapshot of the in-vitro experiment. On top, a snapshot shows the flow inside the Rilsan ® tube enveloped by the Agar-based material. On top, the frame with superimposed the 2-D vector map, the color map and the green vertical cursor related to the velocity profiles on the right. The instantaneous velocity profiles are indicated by continuous lines. The expected velocity profiles are in dotted lines. On bottom, the velocity in the point indicated by the cursors, obtained during the simulation (in continuous lines) and the expected value (in dotted lines).

Figure 6 . 1 :

 61 Figure 6.1: Comparing the CPU and the GPU architectures. The GPU has more transistors but less powerful control units and cache memories [137].

6. 3 125 Figure 6 . 2 :

 312562 Figure 6.2: A section of a renal artery. As suggested by its Latin name, it is responsible for the supplying of the kidney and takes blood from the abdominal aorta. The kidney is located on the left.

Figure 6 . 3 :

 63 Figure 6.3: The external (on top) and the internal (on bottom) carotid arteries taken at 2 cm distal to the bifurcation apex, at end diastolic (left) and systolic (right) instants. The velocities are lower than 60 cm/s also for the systolic peak,

Figure 6 . 4 :

 64 Figure 6.4: The external (on top) and the internal (on bottom) carotid arteries taken at 2 cm distal to the bifurcation apex, at end systolic (left) and start diastolic (right) instants.

Figure 6 . 5 : 129 Figure 6 . 6 :Figure 6 . 7 :

 651296667 Figure 6.5: A section of the subclavian artery at end diastolic (left) and systolic (right) instants. The subclavian artery takes blood from the aortic arch to supply the arm and some parts of the head and the thorax. Such structures are located on the left.

Table 1 .

 1

1: Acoustic properties of the most common biological tissues

Table 2 .

 2 Setup.

	Flow/pipe	Flow rate	Total
	condition	[ml/min]	experiments
	Rigid/Constant	60 ÷ 750, step: 30	504
	Rigid/Sinusoidal	60 ÷ 630, step: 30	420
	Flexible/Sinusoidal 60 ÷ 600, step: 30	399
	Flexible/Pulsed	350 ÷ 560, step: 30 420
	Total experiments		1743

1: Experimental conditions. For each experiment, the steering angle was changed in (10 ÷ 22) • with step: 2 • . The total number takes account of this.

2.

2 Volume Flow Assessment 37 duced

  by the calibrated pump, and then extracting the percentage error Err

% :

Signal Processing Platform for Ultrasound Imaging

	Fitting Results.
		Value	95% confidence
			bounds
	Mean error (Err % )	-3.7% (-4.0 ÷ -3.5)%
	Standard deviation (S d )	±5.3% (±5.0 ÷ ±5.6)%
	Coefficient of determination (R 2 ) 0.96	-
	Total experiments		1743

Table 2 . 2 :

 22 Fitting results.

Standard and Vector Doppler Approaches

  .1, while Tab. 3.2 reports the mean velocity magnitudes measured by each of the 3 methods. Tab. 3.3 shows the regression coefficients (a,b) and the coefficient of determination (R 2 ), while Fig. 3.4 shows the relationship between the measurements of each comparison and the fitting results.

	Bland-Altman comparison of PSVs [cm/s].
		Method 1 vs. 2 Method 1 vs. 3 Method 2 vs. 3
	CCA	25.1 ± 13.2	21.9 ± 15.4	3.2 ± 8.6
	ICA	37.9 ± 22.9	37.5 ± 21.5	0.4 ± 13.1
	Table 3.1: Mean difference of PSVs measured in the pairwise comparison, accord-
	ing to Bland-Altman.
			Mean PSVs [cm/s].
			Method 1 Method 2 Method 3
		CCA	80.6	55.5	58.7
		ICA	89.0	51.0	51.5

Table 3 . 2 :

 32 Mean PSVs measured in methods 1, 2, 3.

		Linear Regression Analysis.
		Method 1 vs. 2 Method 1 vs. 3 Method 3 vs. 2
	a	0.46	0.42	0.83
	b [cm/s]	14.84	19.98	7.00
	R 2	0.63	0.59	0.62

Table 3 .

 3

3: Linear model: Y = aX + b. R 2 is the coefficient of determination.

Table 3 . 4 :

 34 Comparison between measurements by the same operator.

  .5 and detailed in Tab. 3.5 for each method.

	Inter -operator repeatability [CV (%)].
		Method 1 Method 2	Method 3
	CCA	-2.7	-2.6	-2.5
	ICA	11.3	7.9	1 1 .4

Table 3 . 5 :

 35 Comparison between measurements of the two operator.

  1 and [111] 2 .

	1 [49] M. Lenge, A. Ramalli, E. Boni, H. Liebgott, C. Cachard, and P. Tortoli, "High-
	frame-rate 2-D vector blood flow imaging in the frequency domain," published in IEEE
	TUFFC, pp. 1504-1514, 2014. Author's contributions: conceiving and developing the
	method; performing experiments; processing the data; validating the method; writing
	the paper.

2

[111] M. Lenge, A. Ramalli, E. Boni, A. Cellai, H. Liebgott, C. Cachard, and P. Tortoli, "Frequency-domain high frame-rate 2D vector flow imaging," presented at 2013 IEEE IUS, Prague (Czech Republic), pp. 643-646, 2013. Author's contributions: conceiving and developing the method; performing experiments; processing the data; validating the method; writing the paper.

Table 4 . 2 :

 42 Experimental setup. Bold values indicate reference setup.

Table 5 . 1 :

 51 The reference setup (RS) is referred to the parameters in bold type.

	Fields Setup

Table 5 . 3 :

 53 Bold values indicate the reference setup (RS). The axial block sizewas set to 0.95 mm instead of the 0.47 mm used for simulations (Tab. 5.2).

	Velocity estimator parameters Values
	Axial block size (M B δ z )	[0.47, 0.95, 1.42, 1.89, 2.37] mm
	Lateral block size (N B δ x )	[0.98, 1.96, 2.94, 3.92, 4.90] mm
	Note that the values for the axial block size, here expressed in mm, correspond to
	M B = [32, 64, 96, 128, 160] samples. The lateral block size correspond to N B =
	[4, 8, 12, 16, 20] samples.

The description of the ULA-OP system draws on [1] E. Boni, A. Cellai, A. Ramalli, M. Lenge, and S. Ricci, "Multi-channel Raw-Data Acquisition for Ultrasound Research," presented

at 2014 17th Euromicro Conference on Digital System Design (DSD), Verona (Italy), pp. 647-650, 2014. Author's contributions: performing experiments; processing the data.

This section draws on [1] E. Boni, A. Cellai, A. Ramalli, M. Lenge, and S. Ricci, "Multi-channel Raw-Data Acquisition for Ultrasound Research," presented at

17th Euromicro Conference on Digital System Design (DSD), Verona (Italy), pp. 647-650, 2014. Author's contributions: performing experiments; processing the data.

[START_REF] Ricci | Volume flow assessment through simultaneous B-mode and multigate Doppler[END_REF] S. Ricci, M. Cinthio, M. Lenge, R. Matera, J. Albinsson, and P. Tortoli, "Volume flow assessment through simultaneous B-mode and multigate Doppler," presented at

2012 IEEE IUS, Dresden (Germany), pp. 1588-1591, 2012. Author's contributions: developing the methods; processing the data. 2 [65] A. Ramalli, L. Bassi, M. Lenge, C. Palombo, K. Aizawa, and P. Tortoli, "An integrated system for the evaluation of flow mediated dilation," presented at 2014 IEEE ICASSP, Firenze (Italy), pp. 5182-5185, 2014. Author's contributions: conceiving and developing the methods.

The work has been supported by the Swedish Foundation for International Cooperation in Research and Higher Education (STINT) and the Swedish research council.

This work has been supported by the European Fund for Regional Development for the 2007-2013 programming period (POR FES 2007-2013 CReO SAMURAI Project) and by the European Union's Seventh Framework Programme (FP7/2007-2013) for the Innovative Medicine Initiative under grant agreement number IMI/115006 (the SUMMIT consortium).

[START_REF] Lenge | Blood velocity measurement in healthy and diseased carotid arteries by vector Doppler techniques[END_REF] M. Lenge, D. Righi, S. Ricci, H. Liebgott, and P. Tortoli, "Blood velocity measurement in healthy and diseased carotid arteries by vector Doppler techniques," presented at

IEEE IUS, Chicago (IL, USA), pp. 345-348, 2014. Author's contributions: developing the methods; performing experiments; processing the data; writing the paper. 2 [84] P. Tortoli, M. Lenge, D. Righi, G. Ciuti, H. Liebgott, and S. Ricci, "Comparison of carotid artery blood velocity measurements by vector and standard Doppler approaches," in press on

UMB, 2015. Author's contributions: developing the methods; performing experiments; processing the data; writing the paper.

These works has been supported by the European Fund for Regional Development for the 2007-2013 programming period (POR FESR 2007-2013 CReO, ASSO Project).

[START_REF] Lenge | A new method for 2D-vector blood flow imaging based on unconventional beamforming techniques[END_REF] M. Lenge, A. Ramalli, A. Cellai, H. Liebgott, C. Cachard, and P. Tortoli, "A new method for

2D-vector blood flow imaging based on unconventional beamforming techniques," presented at 2014 IEEE ICASSP, Firenze (Italy), pp. 5162-5166, 2014. Author's contributions: conceiving and developing the method; performing experiments; processing the data; validating the method; writing the paper.[START_REF] Feynman | The Feynman Lectures on Physics[END_REF] [51] M. Lenge, A. Ramalli, P. Tortoli, C. Cachard, and H. Liebgott, "Plane-wave transverse oscillation for high frame rate 2-D vector blood flow imaging," submitted for publication in IEEE TUFFC. Author's contributions: conceiving and developing the method; performing experiments; processing the data; validating the method; writing the paper.

Speckle-tracking algorithms have been proposed to determine the motion in 2-D dimensions[START_REF] Bohs | Speckle tracking for multi-dimensional flow estimation[END_REF]. They perform a frame-to-frame analysis of the movement of the speckle patterns in B-mode pulse-echo images[START_REF] Evans | Colour flow and motion imaging[END_REF], by a block-matching technique borrowed from image processing.

The work has been supported by the Programme Avenir Lyon St. Etienne -année 2014 (PALSE no. 40/CA/2013), the Franco-Italian University (2011-2012 n.26075WL and 2013-2014 n.30224WH), the CeLyA and PRIMES LabEx, and the European Fund for Regional Development for the 2007-2013 programming period (POR FESR 2007-2013 CReO, ASSO Project).

In[START_REF] Basarab | Phase-based block matching applied to motion estimation with unconventional beamforming strategies[END_REF], as in other speckle-tracking approaches[START_REF] Bohs | Speckle tracking for multi-dimensional flow estimation[END_REF], the movement is estimated by searching the speckle of a little kernel in a big region of interest. In video processing, this approach is called block-matching. Here, two blocks of the same size are selected from consecutive images and are compared to calculate the phase shift, and the approach is called block-shift. The two blocks are defined on the same indexes of the grid.

The introduction of an Agar-based phantom has been adopted in consequence of the results obtained in the validation of the frequency-domain algorithm, presented in §4.4.2. See also §4.5 for further details.

Simulations Setup Blood vessel parameters

Values Beam-to-flow angle (ϑ) [START_REF] Nicolaides | Ultrasound and Carotid Bifurcation Atherosclerosis[END_REF][START_REF] Cinthio | A method for arterial diameter change measurements using ultrasonic B-mode data[END_REF][START_REF] Krams | Effect of vessel curvature on Doppler derived velocity profiles and fluid flow[END_REF][START_REF] Ekroll | Ultrasound imaging of blood flow based on high frame rate acquisition and adaptive signal processing[END_REF][START_REF] Lenge | A new method for 2D-vector blood flow imaging based on unconventional beamforming techniques[END_REF][START_REF] Bohs | Speckle tracking for multi-dimensional flow estimation[END_REF][START_REF] Evans | Colour flow and motion imaging[END_REF] • Barycenter depth (z 0 ) [START_REF] Misaridis | Use of modulated excitation signals in medical ultrasound. part III: High frame rate imaging[END_REF][START_REF] Kasai | Real-time twodimensional blood flow imaging using an autocorrelation technique[END_REF][START_REF] Newhouse | Three-dimensional vector flow estimation using two transducers and spectral width[END_REF] mm Peak velocity (v 0 )

[ 3 0 , 60, 100] cm/s Method parameters Values Central frequency (f 0 ) [START_REF] Szabo | Diagnostic Ultrasound Imaging: Inside Out[END_REF][START_REF] Van Veen | Beamforming: a versatile approach to spatial filtering[END_REF] MHz Frame-rate (F R) [ about 2 dB/(cm MHz) and a sound propagation velocity of 2600 m/s [START_REF] Ricci | Accuracy and reproducibility of a novel dual-beam vector Doppler method[END_REF].

The fluid consists of 3 g of 10 µm polyamide spherical particles (Orgasol ® , Arkema Inc., Philadelphia, PA) suspended in 2 dm 3 of demineralized water, with a resulting density of 2.6 × 1010 scatterers/mm 3 . The flow rate was controlled by a bypass flow-meter (Asameter D6, ASA S.r.l., Milan, Italy), having ±2.5% accuracy. The tube had the same geometry used in simulations, with a 70 cm extension before the investigated section to guarantee a parabolic flow. To avoid any turbulence inside the tube, the peak velocity was set to 16.6 cm/s, compatible with a 2.5 kHz PRF. The experimental results were obtained through the ULA-OP research system ( §1.3). This was equipped with the LA533 probe (Esaote S.p.A., Florence, Italy) having the characteristics reported in §4.3. The ULA-OP arbitrary waveform generators [START_REF] Ricci | Multichannel FPGA-based arbitrary waveform generator for medical ultrasound[END_REF] produced the same Hanning-weighted bursts used in simulations. The RF echo data received at 2.5 kHz PRF in plane wave mode by each active probe element were stored in the highcapacity (36 GB) memory board, and later processed offline. The parameters used in the experiments are reported in Tab. 4.2, and the central frequency f 0 was set to 6 MHz. One experiment has been conducted in RS conditions, manually changing the beam-to-flow angle between 80 • to 100 • during one acquisition.

Complementary in-vivo test were also made on the common carotid artery (CCA), the radial artery and the femoral vessels of a healthy volunteer. In this case, ULA-OP was programmed to interleave standard Bmode and plane wave transmissions at 10 kHz PRF. A real-time reference B-mode image was useful to check the correct probe-vessel alignment. The The method performance has been evaluated by the relative bias and the standard deviation of the estimated velocity (v) compared to the theoretical parabolic profiles (v) calculated on 32 estimation lines (Fig. 5.1), uniformly distributed in the center of the ROI. Performance metrics of the lateral and the axial velocity components are the same used for the frequency-domain estimator ( §4.3.3).

Experiments

RF images (3 kHz frame-rate) were achieved by the ULA-OP system [START_REF] Boni | A reconfigurable and programmable FPGAbased system for nonstandard ultrasound methods[END_REF] and stored in the high capacity (36 GB) memory board [START_REF] Boni | A high performance board for acquisition of 64-channel ultrasound RF data[END_REF] installed in ULA-OP. The system was equipped with a linear array probe (LA533, Esaote S.p.A., Florence, Italy) located in front of a cylindrical Rilsan ® pipe of 4-mm internal radius. An electronic-controlled pump was pushing a blood mimicking fluid (Orgasol ® /water= 1.5 g/dm 3 ) into the tube.

The phantom was immersed in a water tank and another cubical phantom, made by an Agar-based tissue mimicking material [START_REF] Teirlinck | Development of an example flow test object and comparison of five of these test objects, constructed in various laboratories[END_REF], enveloped it.

Plane-Wave Transverse-Oscillation for High Frame Rate 2-D

Vector Blood Flow Imaging Similar performance is obtained by reducing the lateral wavelength to λ x = 0.98 mm [5.6(a)], respecting the (5.5), and the bias of the lateral velocity overcomes the 12% [Fig. 5.10(a)]. With higher lateral wavelengths (e.g., λ x = 1.96 mm), the TO apodization reaches higher depths [up to 45 mm without cuttings, as Fig. 5.6(b) shows]. However, since the number of samples of the apodization is lower, the transverse oscillations are less defined and the velocity profile assumes a pointed shape, especially for the highest velocities [Fig. 5.10(b)].

The performance is not heavily influenced by the no. of cycles used in transmission (Fig. 5.4), as confirmed by the Udesen et al. results [START_REF] Udesen | Investigation of transverse oscillation method[END_REF]. If N c is reduced down to 2 [Fig. 5.4(a)], the PEFs performance gets worse, the lateral velocity profile is heavily flattened and the velocity performance drops down (Fig. 5.8).

Concerning blood flow simulations obtained varying the blood vessel parameters, the reference setup suggested by previous tests [START_REF] Lenge | A new method for 2D-vector blood flow imaging based on unconventional beamforming techniques[END_REF] has proved to be suitable for the correct evaluation of the blood flow. The apodization at RS entirely covers the range of depths under investigations [up to 35 mm, see Fig. Chapter 6

GPU-Based High Frame-Rate Blood Flow Imaging in the Frequency-Domain

Summary. A GPU computing solution for the frequency domain algorithm (described in Chapter 4) is here proposed. The implementation is tested for processing new in-vivo acquisitions, and the performance, in terms of computational time, are compared with that of a CPU solution. Finally, the perspectives of an optimized GPU version of the approach are discussed.

Author's Publications

Journal Papers