Warning
This page is located in archive.

Cvičení 4 : Počítačová tomografie (CT), přímá projekce

Na dnešním cvičení se budeme zabývat principem měření v počítačové tomografii (CT). CT pracuje při měření objektu (v kartézských souřadnicích x, y) pomocí sady projekcí pod různými úhly \theta. To odpovídá Radonově transformaci vstupního obrázku a dnešním úkolem je naprogramovat tuto transformaci pomocí vlastní funkce myRadon. Při rekonstrukci dat CT z měřeného signálu se pak používá inverzní Radonova transformace, kterou se budeme zabývat příště.

Zadání

  1. Najděte analyticky či graficky maximální absolutní hodnotu p_m = \max |p| a q_m = \max |q|, kde p=x\cos\theta+y\sin\theta, q=-x\sin\theta+y\cos\theta a (x,y)\in [-x_m,x_m] \times [ -y_m, y_m]. [0.5b.]
  2. Diskretizujte rovnici (1). [0.5b.]
  3. Naprogramujte funkci imgRadon=myRadon(img,theta) která vrátí přímou Radonovu transformaci imgRadon [(2p_m+1) \times k] obrázku image [m\times n] pro úhly theta [1\times k]. Úhly zadávejte ve stupních. Dodržte šablonu funkce. [3b.]
  4. Pomocí vámi napsané funkce transformujte Shepp-Loganův fantom, který si můžete vytvořit funkcí phantom, parametrem funkce je velikost fantomu. Do reportu vložte výsledek spočítaný na fantomu s velikostí 256px pro úhly 0,1\ldots 179^\circ. Nezapomeňte na správná měřítka a popisy os. [1b.]

Podrobný popis

Stručná ilustrace souřadných systémů pro danou projekci pod úhlem \theta je na následujícím obrázku.

 Radonova transformace pro projekci p, theta

Souřadná soustava (p, q) vznikne otočením soustavy (x, y) o úhel \theta. Pro transformaci tedy platí


 \begin{equation*}
p = x \cos(\theta) + y \sin(\theta) \qquad \qquad q = -x \sin(\theta) + y \cos(\theta)
\end{equation*}

a pro inverzní transformaci pak platí

 \begin{equation*}
x = p \cos(\theta) - q \sin (\theta) \qquad \qquad y = p \sin(\theta) + q \cos(\theta)
\end{equation*}

Výsledná Radonova projekce ze souřadného systému obrázku (x,y) do Radonova prostoru (\theta,p) je dána integrálem

    J(\theta,p) = \int_{-\infty}^\infty f\left ( p\cdot\cos(\theta) - q\cdot\sin(\theta) ~ , ~ p\cdot\sin(\theta) + q\cdot\cos(\theta) \right ) dq

Možný postup

  1. Je dobré nastavit souřadnice obrázku tak, aby bod [0, 0] ležel uprostřed obrázku
  2. Určete minimální a maximální hodnoty pro p, q (1. krok zadání)
  3. Zvolte diskretizaci pro p, q v rozsahu  [\min p, \max p] , respektive  [\min q, \max q]
  4. Pro všechny úhly \theta a všechny diskretizační body p_h, q_h spočtěte J(\theta,p_h)
    1. pro uvažované body p_h, q_h spočtěte souřadnice x, y při otočení o \theta
    2. intenzitu obrázku v daném bodě, který nebude ve většině případů ležet přímo ve středu pixelu, dostaneme interpolací (funkce interp2
    3. spočtěte numerickou hodnotu J(\theta,p_h)
courses/a6m33zsl/lab04_ct_fwd.txt · Last modified: 2018/03/17 10:54 by herinjan