====== Cvičení 6 : Ultrazvuk ======

V rámci dnešního cvičení si projdeme několik příkladů na ultrazvuk ''(ultrasound, US)'' a také se podíváme na analýzu nasnímaných US sekvencí.


==== Zadání ====

  * Spočtěte hvězdičkou označené příklady ze [[#Příklady US|souboru příkladů]] [1b]
  * [[#Extrakce tepové frekvence]]: Odhadněte z nasnímaných záznamů tepovou frekvenci subjektu 
    * pomocí FFT analýzy  [2b]
    * pomocí autokorelace [2b]  
  * [[#Extrakce - Bonus]]: Detekujte a zobrazte pixely, jejichž časový profil má vysokou podobnost s detekovanou tepovou frekvencí. [1b]

==== Příklady US ====


Popis ultrazvuku a jeho charakteristiky najdete v [[http://cw.felk.cvut.cz/lib/exe/fetch.php/courses/a6m33zsl/ultrasound.pdf|přednáškách]].

  * šíření zvuku v látce $c = \frac{1}{\sqrt{\rho K}}$
  * akustická impedance $Z_t = \rho c$
  * odrazivost zvoleného přechodu $R = \frac{Z_k - Z_t}{Z_k + Z_t}$


** Šíření materiálem **

 Mějme materiál o hustotě $1200 \text{kg/m}^{3}$ a ultrazvukové vlny šířící se v něm rychlostí $1540 \text{m/s}$. Určete objemový modul pružnosti pro tento materiál.

** Šíření materiálem II. **

Pokud se materiálem o hustotě $1000\,\text{kg}/\text{m}^3$ šíří zvukové vlny rychlostí $1500\,\text{m/s}$, jaký je objemový modul pružnosti tohoto materiálu? Jaká je specifická akustická impedance tohoto materiálu?


** Odrazivost tkáně I. **

Mejme tkáň o akustické impedanci $800$kRayl a kost o akustické impedanci  ́$6$MRayl. Jaká energie se na rozhraní této tkáně s kostı odrazí zpět?


** [*] Frekvence snímání **

Mějme ultrazvuk, kterým snímáme pacienta. Rychlost ultrazvuku je 1540m/s a hloubka použitelného skanu je 40cm. Určete maximální snímací frekvenci s jakou můžeme pacienta skenovat, tak bychom zachovali hloubku snímání. Uvažujme, že jeden snímek je tvořen z 200 paprsků.

** [*] Odrazivost tkáně II. **

Mějme tkáň hustotě $900 \text{kg/m}^{3}$, ultrazvukové vlny šířící se v ní rychlostí $1540 \text{m/s}$ a kost o akustické impedanci $6.75 \cdot 10^{6}$Rayl. Určete odrazivost na rozhraní této tkáně.

/*
** [*] Útlum v materiálu **

Ultrazvuk prochází materiálem, kde jeho průchodem amplituda po 70cm klesne na polovinu. Určete lineární koeficient útlumu a amplitudu vlnění po průchodu 15cm tímto materiálem. 
*/
==== Extrakce tepové frekvence ====

Srdeční tep budeme zjišťovat pomocí analýzy rozídlných frekvencí (přes časovou osu) obsažených ve výřezu obsahujícím tepající krkavicí a sice na sekvenci nasnímané [[http://cmp.felk.cvut.cz/~herinjan/docs/raw_long.avi|podélně]] a na sekvenci nasnímané [[http://cmp.felk.cvut.cz/~herinjan/docs/raw_tra.avi|příčně]]. Pro práci s video-soubory v MATLABu použijeme 
 nástroj [[http://www.mathworks.com/help/matlab/ref/videoreader.html?s_tid=doc_ta|VideoReader]].


Pomocí frekvenční analýzy budeme sledovat časově závislý signál $j(t)$ střední intenzity výřezu. Ve zprávě uveďte vypočtené tepové frekvence a v postupu zmíněné grafy.

  - Extrahujte vhodný výřez $P$ (''imcrop'') a spočtěte pro každý frame $t$ průměrnou intenzitu $j(t) = \frac{1}{|P|}\sum_P I(x,y,t)$ a vykreslete funkci $j(t)$ do grafu. {{ :courses:a6m33zsl:us_viz_crop_compl.png?nolink&800 |}}
  - Normalizujte signál $j(t)$ odečtením střední hodnoty (přes klouzavý interval o délce $w$, např. $w=32$.
  - Spočtěte dominantní frekvenci signálu $j(t)$ (tepová frekvence) pomocí následujících metod a srovnejte jejich výsledek:
    * Fourierova transformace 
      - Pomocí funkce [[http://www.mathworks.com/help/matlab/ref/fft.html|fft]] spočtěte spektrum signálu $J(f) = |\mathcal{F}(j)|$
      - Za pomocí funkce [[http://www.mathworks.com/help/matlab/ref/fftshift.html|fftshift]] převeďte $J(f)$ na vector o rozsahu $-\frac{f_r}{2} \leq f \leq \frac{f_r}{2}$, kde $f_r$ je //frame rate// vstupního videa. 
      - Najděte dominantní frekvenci (odpovídající tepu) a vyjádřete jí v úderech za minutu.
      - Vykreslete $J(t)$ a použije tepovou frekvenci na ose $x$. {{ :courses:a6m33zsl:us_jf.png?nolink&400 |}}
 
    * Autokorelace signálu 
      - Je-li potřeba, vyhlaďte signál před použitím autokorelace; například pomocí konvoluce Gaussovským kernelem (funkce [[http://www.mathworks.com/help/matlab/ref/conv.html|conv]]). Gaussián se pro připomenutí spočte pomocí $$g(x) = \frac{1}{\sqrt{2 \pi} ~ \sigma} e^{-\frac{x^2}{2 \sigma^2}}$$
      - Spočtěte autokorelaci filtrovaného signálu (funkce [[http://www.mathworks.com/help/signal/ref/xcorr.html|xcorr]])
      - Z periody autokorelačního signálu (získáme s pomocí funkce [[http://www.mathworks.com/help/signal/ref/findpeaks.html|findpeaks]]) spočtěte tepovou frekvenci. {{ :courses:a6m33zsl:us_xcorr.png?nolink&400 |}}


==== Extrakce - Bonus ====
 Najděte ve vstupní sekvenci pixely, jejichž časový profil má vysokou podobnost s detekovanou srdeční frekvencí $f_h$. Tj. hledáme pro pixely $x$, čas $t$ a frekvenci $f_h$ parametry modelu: $s(x,t) = a(x) $sin$(2\pi f_h t) + b(x)$cos$(2\pi f_h t) + c(x)$ 
    * Sestavíme soustavu lineárních rovnic \\ <latex>
\begin{bmatrix}
  s(x,1) \\
  \vdots \\
  s(x,N)
 \end{bmatrix} =
 \begin{bmatrix}
  \textrm{sin}(2\pi f_h t_1) & \textrm{cos}(2\pi f_h t_1) & 1 \\
  \vdots  & \vdots  & \vdots \\
  \textrm{sin}(2\pi f_h t_N) & \textrm{cos}(2\pi f_h t_N) & 1
 \end{bmatrix}
 \begin{bmatrix}
  a \\
  b \\
  c
 \end{bmatrix} \quad (=) \quad S = A.V
</latex>\\
     * pomocí operace \ ([[http://www.mathworks.com/help/matlab/ref/mldivide.html|mldivide]]) a získáme $V = [a b c]$
     * Pro každý pixel spočteme $a(x)^2 + b(x)^2$.
     * Zobrazte pixely, které mají se signálem srdečního tepu vysokou podobnost, $a(x)^2 + b(x)^2 > \tau$.