TECHNICKÁ UNIVERZITA V LIBERCI
Fakulta mechatroniky a mezioborových inženýrských studií
Studijní program: N 2612 – Elektrotechnika a informatika
Studijní obor: 3902T005 – Automatické řízení a inženýrská informatika
Tvorba modelů kostí zpracováním obrazových dat CT a MR
Creation of bones models by treating image CT and MR dates
Diplomová práce
Autor: Bc. Ondřej Štípek
Vedoucí práce: Ing. Lukáš Matela, Ph.D.
Konzultant: Ing. Josef Novák, Ph.D.
V Liberci dne 15.5.2008
TECHNICKÁ UNIVERZITA V LIBERCI
Fakulta mechatroniky a mezioborových inženýrských studií
Ústav řízení systémů a spolehlivosti Akademický rok: 2007/08
ZADÁNÍ DIPLOMOVÉ PRÁCE
Jméno a příjmení: Bc. Ondřej Štípek
studijní program: N 2612 - Elektrotechnika a informatika
obor: 3902T005 - Automatické řízení a inženýrská informatika
Vedoucí ústavu Vám ve smyslu zákona o vysokých školách č.111/1998 Sb. určuje tuto diplomovou práci:
Název tématu:
Tvorba 3D modelů kostí zpracováním obrazových dat CT a MR
Zásady pro vypracování:
1. Seznamte se s technologií pořizování snímků z CT a MR 2. Seznamte se s algoritmy zpracování 2D a 3D obrazových dat
3. Navrhněte způsob jak z CT a MR snímků vyseparovat tvary kostí od měkkých tkání 4. Navrhněte způsob jak připravit segmentované tvary kostí do podoby vhodné pro další
zpracování užitím MKP
5. Realizujte aplikaci s uživatelským rozhraním pro tvorbu 3D modelů
Prohlášení
Byl jsem seznámen s tím, že na mou diplomovou práci se plně vztahuje zákon č. 121/2000 o právu autorském, zejména § 60 (školní dílo).
Beru na vědomí, že TUL má právo na uzavření licenční smlouvy o užití mé diplomové práce a prohlašuji, že s o u h l a s í m s případným užitím mé diplomové práce (prodej, zapůjčení apod.).
Jsem si vědom toho, že užít své diplomové práce či poskytnout licenci k jejímu využití mohu jen se souhlasem TUL, která má právo ode mne požadovat přiměřený příspěvek na úhradu nákladů, vynaložených univerzitou na vytvoření díla (až do jejich skutečné výše).
Diplomovou práci jsem vypracoval samostatně s použitím uvedené literatury a na základě konzultací s vedoucím diplomové práce a konzultantem.
V Liberci dne 15.5.2008 ...
podpis diplomanta
Poděkování
Chtěl bych poděkovat vedoucímu mé diplomové práce, Ing. Lukáši Matelovi, Ph.D., za cenné rady, připomínky a výklad potřebné teorie při vypracování diplomové práce, Krajské nemocnici Liberec, a.s., za poskytnutí snímků pro testování mých postupů a radiologickým asistentům Františku Šebkovi a Filipu Rasochovi za informace o metodách získávající snímky z lidského těla. V neposlední řadě bych rád poděkoval rodičům za vytvoření skvělých studijních podmínek.
Ondřej Štípek
Anotace:
Cílem diplomové práce je tvorba 3D modelů kostí zpracováním obrazových dat CT a MR. V jejím úvodu je popsán současný stav řešení problému. Další její část je zaměřena na seznámení se s technologií pořizování CT a MR snímků, s jejím vývojem, členěním a možnostmi využití. Práce rovněž hovoří o digitalizovaném obrazu, zejména o formátu DICOM, ve kterém jsou snímky z CT a MR uloženy, a základních metodách zpracování obrazových dat. Jsou v ní popsány metody segmentace, pomocí nichž je možné rozlišit kosti od měkké tkáně. Jde o metody adaptivního prahování Chow a Kaneko, Otsu, lokálního prahování a algoritmus Niblack. Následuje popis postupu tvorby 3D modelů kostí a způsob, jak tyto segmentované modely připravit pro jejich další zpracováni metodou konečných prvků. Pro poloautomatickou tvorbu 3D modelů kostí bylo v prostředí Matlab® vytvořeno jednoduché grafické uživatelské rozhraní (GUI) umožňující prohlížení snímků, provádění jejich segmentace a ukládání výstupů do souborů.
Klíčová slova: CT a MR snímky, DICOM, segmentace kostí, adaptivní prahování, 3D modely kostí.
Annotation:
The aim of this thesis is creation of the 3D bones models by processing CT and MR data. The theoretical part consists of description current way of a data processing.
There is also a part considering technology of CT and MR charts enrolling, their evaluation, classification and their possibility of utilization. The digitalizing of a screen, mainly DICOM format, is also mentioned. CT and MR charts are saved in this format.
Theoretical part ends with a characterization of basic methods of screen data processing. There are described methods of segmentation used for distinction of bones and tissue. Such as adaptive thresholding Chow a Kaneko, Otsu, local thresholding and algorithm Niblack. A theoretical part is followed by description of the 3D bones models creation and their set up for subsequent processing via Finite element method.
For a semi-automatic creation of the 3D bones models was prepared a simple graphical user interface (GUI). The graphical user interface allows viewing of charts, performing their segmentation and saving outputs to files.
Keywords: CT and MR charts, DICOM, bones segmentation, adaptive thresholding, 3D bones models.
Obsah
1 ÚVOD ... 9
2 ZÍSKÁVÁNÍ SNÍMKŮ LIDSKÉHO TĚLA... 11
2.1 RENTGEN...11
2.2 POČÍTAČOVÁ TOMOGRAFIE (COMPUTED TOMOGRAPHY,CT) ...13
2.2.1 Kategorie CT systémů ...16
2.2.2 Spirální počítačová tomografie (Spirální CT)...18
2.2.3 Artefakty na CT snímcích ...19
2.3 MAGNETICKÁ REZONANCE (MR,MRI) ...21
2.4 DALŠÍ METODY ZÍSKÁVÁNÍ SNÍMKŮ...23
2.4.1 Ultrazvuk...23
2.4.2 Skenování řezů lidského těla...23
3 DIGITALIZOVANÝ OBRAZ A POSTUPY JEHO ZPRACOVÁNÍ.... 24
3.1 DIGITALIZOVANÝ OBRAZ...24
3.2 CT A MR SNÍMKY...25
3.3 FORMÁT DICOM ...25
3.3.1 Zobrazení DICOM snímků ...27
3.4 HISTOGRAM...28
3.5 FILTRACE OBRAZU...29
3.5.1 Vyhlazování obrazu ...30
3.5.2 Hranové detektory ...30
3.6 MATEMATICKÁ MORFOLOGIE...30
3.7 3D METODY ZPRACOVÁNÍ OBRAZU...32
3.7.1 3D filtrace...32
3.7.2 Vyplnění děr 3D matice ...33
3.7.3 Odstranění hraničních objektů ve 3D ...33
4 METODY SEGMENTACE CT SNÍMKŮ ... 34
4.1 PROBLÉMY PŘI SEGMENTACI CT SNÍMKŮ...35
4.2 PRAHOVÁNÍ...35
4.3 ADAPTIVNÍ PRAHOVÁNÍ...37
4.3.1 Metoda Chow a Kaneko ...37
4.3.2 Metoda Otsu (metoda optimálního prahu) ...37
4.3.3 Lokální prahování ...39
4.3.4 Algoritmus Niblack ...42
4.4 POROVNÁNÍ METOD ADAPTIVNÍHO PRAHOVÁNÍ...43
5 PRAKTICKÁ REALIZACE TVORBY 3D MODELŮ KOSTÍ ... 46
5.1 NAČTENÍ CT SNÍMKŮ...47
5.2 SEGMENTACE...47
5.2.1 Segmentace kostí ...47
5.2.2 Automatické určení prahu...49
5.3 VYTVOŘENÉ GRAFICKÉ UŽIVATELSKÉ ROZHRANÍ (GUI) ...50
5.4 NÁVRH DALŠÍHO ZPRACOVÁNÍ VYTVOŘENÝCH 3D MATIC...51
5.5 VIZUALIZACE VYTVOŘENÝCH 3D MATIC...51
6 ZÁVĚR ... 53
LITERATURA... 55
PŘÍLOHA A... 57
A.1 SOFTWARE PRO PROHLÍŽENÍ SNÍMKŮ VE FORMÁTU DICOM ...57
A.2 ZDROJE DAT...57
A.3 POROVNÁNÍ METOD ADAPTIVNÍHO PRAHOVÁNÍ...58
A.3.1 CT snímek hlavy (čelist a krční páteř)...58
A.3.2 CT snímek páteře ...60
A.4 SPIRÁLNÍ CT ZAŘÍZENÍ Z KRAJSKÉ NEMOCNICE LIBEREC, A.S...61
PŘÍLOHA B... 62
B.1 UŽIVATELSKÝ MANUÁL APLIKACE PRO SEGMENTACI KOSTÍ...62
B.1.1 Ovládání uživatelského rozhraní...63
B.1.2 Postup realizace segmentace...64
PŘÍLOHA C... 67
C.1 SEZNAM ELEKTRONICKÝCH PŘÍLOH NA CD...67
Seznam použitých zkratek a symbolů:
2D Dvourozměrný prostor 3D Třírozměrný prostor
MKP Metoda konečných prvků
CT (Computed Tomography) Počítačová tomografie H.u. (Hounsfield unit) Hounsfieldovo číslo
RTG Rentgen MR, MRI (Magnetic Resonance Imaging) Magnetická rezonance
OCR (Optical Character Recognition) Optické rozpoznávání znaků GUI (Graphic User Interface) Grafické uživatelské rozhraní DICOM (Digital Imaging and
Communications in Medicine)
Standard digitálního zobrazení a komunikace v medicíně
LUT tabulka (LookUp Table) Transformační tabulka jasů pixelů px (Pixel) Obrazový prvek (Bod digitálního obrazu) ω Pravděpodobnost výskytu
μ Střední hodnota
σ2 Rozptyl (Střední kvadratická odchylka)
IGES (The International Graphics Exchange Specification)
Výměnný souborový formát pro 3D počítačem podporované navrhování CATIA (Computer-Graphics Aided Three
Dimensional Interactive Application)
Software pro 3D počítačové konstruování
ANSYS
Programový balík založený na MKP pro řešení fyzikálních i geometrických (lineárních i nelineárních) úloh mnoha různých kategorií
1 Úvod
Díky rozvoji výpočetní techniky a jejímu zavedením do medicíny došlo v 70. letech 20. století k výraznému pokroku v oblasti zobrazovacích metod v medicíně.
Byly objeveny zobrazovací metody výpočetní tomografie a magnetická rezonance. Tyto metody jsou vhodné pro zobrazení různých orgánů v těle. Výpočetní tomografie (CT) vychází z principu rentgenu, jsou zde tedy dobře viditelné kosti, chrupavky a kovové předměty uvnitř pacientova těla, které jsou zobrazeny na snímku bíle nebo světle šedě.
Tmavě šedými barvami jsou zobrazeny tkáně jako je mozek a střeva. Černou barvou jsou zobrazeny plíce a tělní dutiny.
Magnetická rezonance (MR) je založena na vytvoření silného magnetického pole, které přístroj vytváří a ovlivňuje pohyb vodíkových iontů v těle. Vychází z toho, že vody je v těle nadpoloviční většina. Molekula vody obsahující vodík je tak skoro v každé buňce. Jelikož každá tkáň má jiný obsah vody, jsou na snímku zobrazeny orgány různými jasy šedé barvy. Tato metoda se využívá zejména při vyšetření cévních či nádorových onemocněních mozku, při poranění páteře a výhřezech meziobratlových plotének. Kosti, které obsahují minimum vody, nejsou vhodné pro sledování touto zobrazovací metodou.
Na rozdíl od rentgenového snímku, kde je barva v jednom místě snímku určena z průmětu tkání mezi rentgenkou a filmem, nedochází u výpočetní tomografie a magnetické rezonance k překrývání tkání a orgánu v těle. U CT a MR snímků je barva každého bodu určena z jednoho bodu v těle. Tělo je jakoby rozřezáno na tenké plátky, ze kterých jsou určeny výsledné jasy šedé barvy. Ve skutečnosti se vyšetření pomocí těchto metod provádí bez mechanického poškození. Jde tedy o neinvazivní metody. Výstupem CT a MR zařízení není pouze jeden snímek, ale sada snímků, které jsou rovnoběžné a s konstantní vzdáleností od sebe.
Pokud budeme snímky analyzovat a umístíme výsledky do 3D matice, můžeme tak vytvářet 3D modely různých částí těla. Ve své diplomové práci se zabývám problémem vytváření 3D matic nasegmentovaných kostí. Nejde o zobrazování 3D modelů kostí. Zjednodušený proces řešení mé diplomové práce je zobrazen na obrázku 1.1.
Téma diplomové práce vzešlo od konzultanta Ing. Josefa Nováka, Ph.D., na jehož pracovišti se plánuje řešení a realizace náhrad páteřních obratlů z vytvořených 3D matic pomocí metody konečných prvků (MKP).
Obrázek 1.1: Proces řešeného problému
Vytvářením 3D modelů z CT a MR snímků se zabývají týmy lidí, odborníků a výzkumné firmy. Jde tedy o problém, který není do dnešní doby uzavřen. Většinou se specializují pouze na vybrané části těla z důvodů rozmanitosti ve tvarech kostí a různorodosti měkkých tkání. Důvodem realizování této práce je také vysoká pořizovací cena programu, který by tyto 3D modely vytvářel. Většinou jde o speciální programy, které jsou dodávané s CT či MR zařízeními pro medicínské účely. Často jsou tyto systémy uzavřené a dodaný software tak nelze využívat samostatně.
2 Získávání snímků lidského těla
Snímky lidského těla se využívají v medicíně a praktiky získávání těchto snímků jsou velice známé. Díky vývoji výpočetní techniky se začaly postupem času uplatňovat teoretické znalosti získávání snímků lidského těla a docházelo a nadále dochází k vývoji nových zařízení. Stručně zde uvedu nejpoužívanější způsoby a metody získávání snímků lidského těla.
2.1 Rentgen
Dne 8. listopadu roku 1895 začal provádět experimenty s elektrickými výboji ve vakuové trubici Wilhelm Conrad Röntgen a objevil záření, které nazval paprsky X.
Jelikož šlo o doposud neznámé záření, bylo pojmenováno po svém objeviteli. Röntgen dostal za tento objev v roce 1901 jako první Nobelovu cenu za fyziku. Princip rentgenu je vidět na obrázku 2.1, kde zdroj záření (rentgenka) je zaměřena na vybranou část těla a procházející paprsky jsou zaostřeny na fotografický materiál (film), na kterém po vyvolání vznikne negativní obraz.
Obrázek 2.1: Princip rentgenu
Rentgenové záření je forma elektromagnetického záření o vlnových délkách 10 nanometrů až 0,1 nanometrů [5]. Vzniká prudkým zbrzděním velmi rychle letících elektronů o látku s vysokým atomovým číslem. Hlavní vlastnosti záření jsou:
a) Ionizační schopnosti (narušení elektrické rovnováhy atomu), fotochemický efekt, luminiscenční efekt (na vrstvě luminoforu dochází k fluorescenci).
b) Různé látky pohlcují rentgenové záření v různé míře (závisí zejména na protonovém čísle prvků). Čím je vyšší protonové číslo objektu (např. kosti, které obsahují vápník), tím je rentgenové záření více pohlcováno.
c) Pohlcení rentgenového záření závisí na tloušťce exponovaného materiálu. Čím je větší tloušťka materiálu, tím více se záření pohltí.
U atomů s nízkým protonovým číslem (tzv. lehké atomy) může dojít ke Comptonovu rozptylu. Při tomto jevu dochází k typu absorpce, kdy foton rentgenového záření narazí na elektron, který uvolní z jeho dráhy. Foton ztratí část své energie, změní směr pohybu a pokračuje dál jako rozptýlené záření o větší vlnové délce [15]. Ke Comptonovu rozptylu pak dochází zejména v tukové tkání a při vysoké kilovoltáži zařízení [4]. Na snímku tím dochází ke snížení kontrastu.
Obrázek 2.2: Rentgenový snímek pánve [5]
Na obrázku 2.2 je vidět rozmanitost jasů jednotlivých tkání snímané oblasti lidského těla (pánve). Největší jas vykazují kosti, což je velice vhodné pro hledání kostí v lidském těle, diagnózu zlomenin kostí, atd. Jelikož jde pouze o průmět do jedné roviny, kde nejsme schopni určit přesnou polohu zkoumaných objektů, protože dochází k jejich překrývání, není tato metoda pořizování snímků vhodná pro tvorbu modelů kostí využitím analýzy obrazu.
2.2 Počítačová tomografie (Computed Tomography, CT)
Tato metoda (někdy i nazývána výpočetní tomografie) odstraňuje nedostatky klasického rentgenu, kde docházelo k překrývání jednotlivých zkoumaných objektů na snímku. Jelikož se jedná o výpočetně náročný proces, její zavedení bylo možné až se zavedením počítačů do medicíny.
Objevitelem počítačové tomografie byl v roce 1971 Brit Sir Godfrey Newbold Hounsfield, který byl následně u realizace prvního počítačového tomografu ve výzkumných laboratořích firmy EMI v Londýně. Nezávisle stejný objev učinil i Američan Allan McLeod Cormack. Oba dva dostali v roce 1979 Nobelovu za fyziologii a medicínu [16].
Název počítačového tomografu pochází z řeckého slova tome, které znamená řez. Pacient leží při vyšetření na vyšetřovacím stole, který je zasunut do přístroje.
Po kruhové ose obíhá zařízení složené z rentgenky a detektorů kolem pacienta (obrázek 2.3). Pacient je tak prozářen z nejrůznějších úhlů v jedné rovině. Za pomocí počítače, který je nedílnou součástí tomografu, je zrekonstruován plošný řez části pacientova těla. Každý řez je vytvořen matematickou rekonstrukcí předmětu ze znalosti průmětů předmětu do různých směrů. To spočívá ve vyřešení velké soustavy rovnic.
Jednotlivé řezy jsou rozděleny do matic malých objemových elementů (voxelů) se čtvercovou základnou a výsledek řešení soustavy rovnic spočívá v přiřazení každému voxelu konstantní hodnotu útlumu [7].
Obrázek 2.3: Princip počítačového tomografu [7]
Matematické odvození rekonstrukce prostorového obrazu na základě znalosti jeho projekcí popsal ve své práci roku 1917 Johann Radon. K řešení tzv. obrazové rekonstrukce se používá Radonova transformace (viz. obrázek 2.4 a rovnice 2.2.1, 2.2.2 a 2.2.3[2]), která z obrazu vytváří projekce, a zpětná Radonova transformace, která rekonstruuje z projekcí obraz.
Obrázek 2.4: Princip Radonovy transformace[2]
(2.2.1) (2.2.2) (2.2.3) Absorpční vlastnosti tkáně jsou ohodnoceny Hounsfieldovým číslem (H.u.), které je v rozsahu od -1000 H.u. pro vzduch (který má na výsledném snímku nejmenší jas) až do +1000 H.u.[7]. pro kosti (které mají na snímku nejvyšší jas). Každý voxel na snímku je pak zobrazen rozdílným jasem. Jelikož je rozmezí Hounsfieldových čísel tak veliké, je výsledný kontrast mnohanásobně vyšší než na klasickém rentgenovém snímku.
Touto metodou získáváme mnoho obrazů řezů lidského těla. Jednotlivé jasy jsou dány útlumem tkáně v daném místě. Získaný snímek tak obsahuje pouze informace o jednotlivém bodu v těle. Snímků se většinou pořizuje více (cca několik stovek) v rovnoběžných rovinách s konstantní vzdáleností řezů lidského těla (obrázek 2.5). Vzdálenost mezi jednotlivými řezy se pohybuje od 0,4mm až 2mm.
θ
θ sin
cos + ⋅
⋅
=x y
r
θ
θ cos
sin + ⋅
⋅
−
= x y
s
∞∫
∞
− ⋅ − ⋅ ⋅ + ⋅
= f r s r s ds
r
Rθ( ') ( cosθ sinθ, sinθ cosθ)
Obrázek 2.5: Získané CT snímky řezů pacienta v rovnoběžných rovinách
Jelikož je zde velice dobře rozpoznatelná kost, protože vykazuje vysoký jas, a CT snímky se zde vytváří v mnoha rovnoběžných rovinách, je tato metoda vhodná pro vytváření modelů kostí. Na obrázku 2.6 je zobrazen CT snímek hlavy.
Obrázek 2.6: CT snímek hlavy
2.2.1 Kategorie CT systémů
Systémy výpočetní tomografie lze rozdělit dle vývoje do několika kategorií [7].
U CT systémů první generace jsou rentgenka a jeden až dva detektory umístěny proti sobě. Dochází k jejich současnému lineárnímu posunu a k společné rotaci (obrázek 2.7). U těchto systémů dochází k snímání dvou tomografických vrstev současně.
Obrázek 2.7: CT systém první generace
CT systémy druhé generace mají již více detektorů a svazek záření má vějířovitý tvar. Lineární posun rentgenky a detektorů je oproti předchozí generaci menší. Obsahují již sendvičový a lamelový kolimátor, který usměrňuje rentgenové záření do určitého místa, a je tak snížena radiační dávka pacienta, která byla dříve vlivem radiačních stínů vyšší (obrázek 2.8). Jelikož jsou detektory seskupeny ve tvaru vějíře, nelze použít Radonovu transformaci (paprsky již nejsou paralelní). Je zde použita Fan-beam[7] (vějířová) rekonstrukce.
Obrázek 2.8: CT systémy druhé generace
U CT systémů třetí generace je již odstraněn lineární posun rentgenky a detektory jsou umístěny ve vějířovém svazku (obrázek 2.9). Rotuje tedy pouze rentgenka, která je pulzně buzená, a stovky detektorů. Jde asi o nevíce používaný druh CT systémů.
Obrázek 2.9: CT systémy třetí generace
U CT systémů čtvrté generace rotuje pouze rentgenka a detektory umístěné po celé délce kruhu jsou stabilní (obrázek 2.10). Může se jednat i o inverzní způsob uspořádání, kdy rotuje svazek detektorů a rentgenky umístěné po celé délce kruhu jsou stabilní.
Obrázek 2.10: CT systémy čtvrté generace
Dělení CT systémů do kategorií je spíše symbolické, většina výrobců samozřejmě do systému přidá i jiné vlastní metody a technologie získávání CT snímků.
2.2.2 Spirální počítačová tomografie (Spirální CT)
K velkým změnám došlo v roce 1985 [7], kdy bylo zavedeno spirální CT.
Elektrická energie se zde přenáší bez použití kabelů z nehybného zdroje energie do rotačního prstence. Díky tomu dochází k trvalé rotací RTG lampy kolem vyšetřovacího stolu, který kontinuálně zajíždí do vyšetřovacího tunelu přístroje (otvoru v gantry).
Výsledkem obou pohybů je spojité snímání ve tvaru spirály kolem vyšetřovaného objektu (obrázek 2.11), proto se přístroj nazývá spirální (v angličtině Spiral CT nebo Helical CT).
Obrázek 2.11: Princip spirálního CT
Můžeme například získat snímky plic za 20 až 30 sekund při zadržení dechu pacienta[8]. Odpadá tak získávání mnoha snímků jednotlivých řezů, které můžou být znehodnoceny vlivem dechu či pohybu vyšetřovaného pacienta. Vyšetření je mnohem rychlejší než u předešlých CT systémů. Spirální CT získává obrazová data z celého objemu vyšetřované oblasti najednou (ne pouze z jednotlivých vrstev). Tento objem dat je následně v počítači rekonstruován do libovolného počtu axiálních CT snímků, trojrozměrného zpracování krevních cest (aorty, atd.) a 3D CT obrazů (MPR – multiplanární rekonstrukce, MIP – projekce maximální intenzity a vytváření 3D objektů [9]), které například usnadní představení komplexních zlomenin.
Nové „více snímkové“ („multi-slice“) CT systémy jsou již vyvinuty tak, že mohou získat až 64 řezů dat v spirálním režimu CT a některé dosahují rychlostí až 120 ot/min.
Tyto systémy mohou sebrat až 128krát více dat než předchozí spirálové CT systémy, které se otáčely 60 ot/min. a sbíraly pouze jeden řez dat. Mnoha snímkové CT systémy tak umožní rychlejší vyšetření. Rychlejším snímáním je docílena eliminace artefaktů vzniklá pohybem a dýcháním pacienta.
Nejmodernější mnohosnímkové spirální CT systémy mohou sebrat 64 řezů za otáčku s rychlostí posunu vyšetřovacího stolu 87mm/s [10].
Na obrázku 2.12 je zobrazeno spirální CT zařízení z Krajské nemocnice Liberec, a.s.. Další fotografie částí CT zařízení jsou zobrazeny v příloze A. 4. Jde o CT zařízení firmy Philips, typ Brilliance™ 16 (16 řezu za otáčku).
Obrázek 2.12: CT zařízení Philips Brilliance™ 16 z Krajské nemocnice Liberec, a.s.
Spirální CT je hlavní zařízení pro zobrazování hrudi, plic a kostí kvůli jeho schopnosti spojitého a rychlého získávání snímků ve vysokém rozlišení.
2.2.3 Artefakty na CT snímcích
Při snímání pomocí CT může díky principu vzniknout mnoho artefaktů a také při vadném či špatně obsluhovaném počítačového tomografu může dojít k určitým defektům.
Jeden z podstatných artefaktů může vzniknout při snímání části těla, kde se nachází kovové části (náhrady lebečních kostí, zubní plomby, atd.). Na obrázku 2.13 je CT snímek hlavy, kde jsou vidět „stíny“ vzniklé v důsledku přítomnosti zubních plomb.
Obrázek 2.13: Artefakt způsobený přítomností zubní plomby
Další artefakty mohou vzniknout v důsledku pohybu pacienta nebo i jeho dýcháním. Takovéto snímky jsou pro získávání modelů kostí nevhodné. Těmto artefaktům lze předcházet nebo je alespoň zmírnit rychlostí pořizováním snímků a také použitím spirálního CT.
2.3 Magnetická rezonance (MR, MRI)
Magnetická rezonance vůči ostatním zobrazovacím metodám má větší přesnost při zobrazení většiny vnitřních orgánů. Orgány, jako jsou nervy či mozek bylo možné zobrazovat až pomocí této metody. Je ovšem nejnáročnější na pořizovací a provozní náklady, ale i délku vyšetření pacienta. Pacient, podobně jako při CT vyšetření, leží na vyšetřovacím stole, který je zasunut do otvoru přístroje (gantry).
Jde o neinvazivní metodu, která využívá magnetických vlastností jader některých prvků a rezonančních jevů, k nimž dochází ve velkém magnetickém poli a při elektromagnetickém vlnění s vysokou frekvencí. Nemá tedy žádné prokazatelné vedlejší účinky na vyšetřovaný organismus. Někdy bývá tato metoda nazývána i jako jaderná tomografie.
Magnetickou rezonanci vyvíjeli od roku 1973 Brit Sir Peter Mansfield a Američan Paul Christian Lauterbur. Oba dva dostali v roce 2003 Nobelovu za fyziologii a medicínu [16]. První snímek pořízený pomocí této metody vznikl v roce 1977.
Princip magnetické rezonance využívá toho, že atomová jádra s lichým počtem protonů či neutronů mají nenulový vlastní spin, neboli rotaci kolem vlastní osy.
V lidském těle je ze všech prvků nejvíce zastoupen vodík, jehož jádra se chovají jako miniaturní magnety. V podmínkách tepelné rovnováhy a v nepřítomnosti magnetického pole bude směr i orientace jejich magnetických momentů náhodný a magnetický moment je nulový. Pokud je ovšem jádro umístěno v konstantním magnetickém poli, dochází k rotaci ve směru působení tohoto pole. Při přerušení magnetického pole se jádro vrátí do své klidové polohy. Přidáním druhého kolmo působícího vysokofrekvenčního magnetického pole (transverzálního) k zevnímu magnetickému poli dojde k excitaci jader [19]. Toto magnetické pole bývá o velikosti v rozmezí 0,5 až 3 Tesla.
Aby došlo k rotačnímu pohybu jádra, musí frekvence magnetického pole odpovídat jeho rezonanční frekvenci, která se někdy též nazývá vlastní frekvence.
Volbou velikostí obou magnetických polí lze velmi přesně určit, která jádra budou v rezonanci. Při rezonanci je magnetický moment jádra překlopen o 90° a osa pak rotuje podle transverzálního pole. Velikost naklonění je úměrná amplitudě transverzálního pole a době jeho působení. Po vypnutí tohoto pole přecházejí jádra z excitovaného stavu do rovnovážné polohy dané stálým magnetickým polem. Při této změně vyzařují získanou energii o shodné frekvenci, kterou bylo vytvořeno
vysokofrekvenční magnetické pole. Amplituda vlnění odpovídá počtu jader ve studovaném objektu. Tato energie se indukuje v cívce. Velikost změřeného napětí je závislá na poloze a typu vyšetřované tkáně. Na základě změřeného napětí a dalších parametrů je pomocí složitých postupů frekvenční analýzy za pomocí počítače vytvořena matice číselných údajů, které reprezentují protonovou skladbu voxelů jako funkci jejich polohy. Aby mohl být výsledný snímek zobrazen, jsou k těmto číselným údajům přiřazeny různé stupně šedi.
Přístroj je vizuálně velice podobný CT. Je zde také veliký prstenec, který ovšem ve svých útrobách ukrývá veliký magnet (supravodivou cívku) a množství radiofrekvenčních cívek, které vytváří pulzní vysokofrekvenční magnetické pole a také přijímají vysílané rezonanční frekvence. Používají se gradienty magnetického pole různých směrů, díky nimž získáváme průměty spinových hustot v těchto směrech [19].
Rozlišovací schopnost je až 0,1×0,1×0,2 mm na voxel [11]. Váha nejmodernějších přístrojů je přibližně 6000 kg[11]. Délka vyšetření je závislá na vyšetřované oblasti, požadavků lékaře a typu použitého přístroje. Vyšetření trvají většinou v rozmezí 10 až 30 minut[5].
Obrazy řezů lidského těla jsou stejně jako u CT složeny z rovnoběžných řezů, které lze promítat v různých rovinách. Jednotlivé jasy jsou dány získanými parametry vyšetřované tkáně, kterých je mnohem více než u CT, kde měříme pouze útlum tkáně.
Na výsledných snímcích je velice špatně rozpoznatelná kost, což vyplývá z principu magnetické rezonance, protože kost obsahuje velice málo atomů vodíku. Má tedy nepatrný jas, který se téměř shoduje se vzduchem. Na obrázku 2.14 je zobrazen MR snímek hlavy. Tato metoda není tedy vhodná pro vytváření modelů kostí, ale spíše pro zkoumání měkkých tkání jako je mozek, nervy, svaly, atd.
Obrázek 2.14: MR Snímek hlavy [19]
2.4 Další metody získávání snímků
2.4.1 Ultrazvuk
Ultrazvuk je metoda, která využívá frekvencí vyšších, než dokážeme slyšet, tedy přibližně 20kHz. Někdy se též nazývá sonografie. Ve srovnání s MR nebo CT jde o relativně levnou a přenosnou technologii. Nemá také žádné vedlejší účinky na vyšetřovaného pacienta. Proto bývá velice často používána.
Ultrazvukové vlny procházejí tělem a jsou odráženy nebo pohlcovány jednotlivými orgány, které jsou umístěny v různé hloubce a dle síly odrazu či pohlcení mají různý jas. Odražené vlny jsou převedeny ve formě jasově modulovaného obrazu na monitor. Získaný obraz je podobný CT obrazu. Vlnění je vyšetřovací sondou neustále vysíláno a přijímáno. Na monitoru je tak vidět pohyb orgánů. Proto se tato metoda využívá hlavně při vyšetřování plodu v období těhotenství. Na Obrázku 2.15 je zobrazen ultrazvukový snímek lidského plodu.
Obrázek 2.15: Ultrazvukový snímek lidského plodu [17]
Kvalita snímků je podstatně horší než u CT a MR, proto tato metoda není vhodná pro vytváření modelů kostí.
2.4.2 Skenování řezů lidského těla
Jde o destruktivní metodu, která se nevyužívá v lékařství. Nejprve se provede hluboké zmražení lidského těla, které je následně rozřezáno na přibližně 1mm tlusté plátky. Poté jsou plátky skenovány do počítače [12].
Získané obrazy mají vysokou kvalitu i rozlišení. Proto by byly pro tvorbu modelů kostí nejvhodnější. Jelikož jde o destruktivní metodu, nenalézá v medicíně uplatnění.
3 Digitalizovaný obraz a postupy jeho zpracování
V této kapitole se budu věnovat digitalizovanému obrazu, zejména pak formátu DICOM, v kterém jsou uloženy snímky z CT a MR. Základní metody zpracování obrazu zde uvedu pouze pro kompletnost, vzhledem k jejich známosti a detailnímu popisu v [1]. Více se budu zabývat lokálním prahování a jeho použitím u CT snímků.
3.1 Digitalizovaný obraz
Ve většině elektronických zařízeních jakými jsou mobilní telefony a digitální fotoaparáty, může být uložen digitalizovaný obraz. Omezíme se na obraz, který má pouze jasové stupně šedi. Čidla pro vstup obrazové funkce jsou většinou zdrojem spojitého signálu. Aby mohla být obrazová funkce dále zpracována v elektronickém zařízení, musí dojít k její digitalizaci. Jelikož je obraz dvourozměrný, digitalizace spočívá ve vzorkování spojitého obrazového signálu v matici o rozměru M×N bodů.
Při vzorkování musí být dodržena Shannonova věta o vzorkování. U obrazu by měl být interval vzorkování zvolen tak, aby byl menší nebo rovný polovině rozměru nejmenších detailů v obraze[1]. Dochází také ke kvantování spojité jasové úrovně obrazové funkce každého vzorku do K intervalů. Většinou se používá kvantování 8 bitů na obrazový element a bývá rozděleno do stejných intervalů. 8 bitů odpovídá 256 hodnotám jasů šedi (28 = 256 různých hodnot). Výsledná jasová funkce v digitalizovaném obrazu má tak pouze celočíselné hodnoty. Čím větší jsou rozměry M a N výsledné matice a hodnoty kvantování, tím je dokonaleji původní spojitý obrazový signál aproximován.
Navzorkované body matice se v počítačové grafice nazývají pixely, jde o zkratku z anglických slov Picture Element. Počet sloupců matice udává horizontální rozlišení a počet řádků vertikální rozlišení.
Operace s obrazem je velice podobná práci s maticemi. Mřížka matice se nejvíce používá čtvercová (Obrázek 3.1 a), i když je příčinou problémů se spojitostí oblastí. Snadno se realizuje a odpovídá pravoúhlé matici. Odstranění problému se spojitostí lze dosáhnout použitím hexagonální mřížky, která vykazuje pravidelnost vzhledem k šesti okolí bodům každého bodu a jednotlivé body jsou od sebe stejně vzdálené (Obrázek 3.1 c). Její problém je v samotné realizaci a není také vhodná pro provádění některých operací.
Obrázek 3.1: a) Čtvercová matice jasů b) Interpretace matice jako obraz c) Hexagonální mřížka
3.2 CT a MR snímky
Jelikož pro vytváření 3D modelů kostí jsou vhodnější CT snímky, budu se v této části zabývat pouze jejich vlastnostmi.
CT snímky jsou již z principu šedotónové. Mezi hlavní parametry CT snímku patří rozlišení, které bývá většinou 512×512 pixelů. Při tomto rozlišení jsou od sebe pixely vzdáleny 0,1 až 1 mm. Vzdálenost závisí zejména na zvoleném zvětšení při snímání. Při vyšetření dojde nejdříve k tzv. hrubému snímání, kdy se provede sejmutí skoro celého těla pacienta a následně se označí oblast zájmu pro detailnější snímkování. Dalším parametrem je barevná hloubka, neboli počet jasových úrovní šedé barvy, které se v obraze vyskytují. CT snímky nejvíce užívají maximálně 4096 jasů šedi. Pixel je tak uložen na 12 bitech. Maximální hodnota jasu se v obraze vyskytuje pouze při výskytu hutných kostí a nejčastěji při výskytu kovových předmětů.
Jak jsem zjistil při načítání snímků v prostředí Matlab®, jednotlivé pixely zabírají v paměti 16 bitů z čehož je využito maximálně zmíněných 12 bitů.
Důležitým parametrem pro tvorů 3D modelů kostí je také vzdálenost mezi jednotlivými CT snímky, resp. vzdálenost mezi jednotlivými řezy. Obecně můžeme říci, že čím menší je vzdálenost mezi řezy, tím můžeme docílit přesnějších 3D modelů kostí. Tato vzdálenost se pohybuje od 0,1 mm až po 1 cm a je především závislá na účelu a velikosti zkoumaného objektu v těle pacienta. Tato vzdálenost bývá větší než vzdálenost mezi pixely na CT snímku.
3.3 Formát DICOM
Objevením počítačové tomografie, různých dalších digitálních diagnostických zobrazovacích metod v 70. letech 20. století a s rostoucím použitím počítačů v medicíně vznikala potřeba standardizovat metody pro přenos digitálního obrazu a přidruženou informaci mezi zařízeními vyráběnými různými výrobci.
a) b) c)
Americká vysoká škola radiologie (ACR) a Národní sdružení elektrotechnických výrobců (NEMA) vytvořili v roce 1983 výbor pro vyvinutí standardu [18]. Mezi hlavními body byla podpora komunikace digitální obrazové informace bez ohledu na výrobce zařízení, dále pak usnadnění vývoje a rozšíření archivace digitálního obrazu. Ale také vývoj komunikačních systémů, které mají rozhraní pro komunikaci s dalšími nemocničními zařízením a umožnění vzniku diagnostických databází, které mohou být přístupné z široké škály zařízení.
V roce 1985 vznikla první verze standardu DICOM (Digital Imaging and Communications in Medicine), který umožnil digitální komunikaci mezi diagnostickým a léčebným vybavením a systémy různých výrobců. Je tak zajištěna kompatibilita všech radiologických zařízeních napříč všemi zeměpisnými regiony. Pracovní stanice, CT a MR snímací zařízení, zařízení pro digitalizaci obrazu, sdílené archivy, tiskárny, atd. mohou spolu komunikovat pomocí standardu DICOM napříč systémem sítí.
Následkem toho jsou CT a MR snímky vytvořeny a předány rychleji, lékaři tak mohou stanovit diagnózu a terapii dříve.
DICOM standard usnadňuje součinnost lékařského zobrazovacího zařízení specifikující [18]:
• Síťovou komunikaci – sadu protokolů používaných zařízeními pro síťovou komunikaci.
• Syntaxi a sémantiku příkazů a přidružené informace pro výměnu dat.
• Prostředky pro komunikaci – prostředky pro uložení služeb (ukládání na vyměnitelná média, atd.), formát souborů a adresářovou strukturu.
• Informace potvrzující shodu se standardem, která musí být dodána s realizovaným zařízením od dodavatele.
DICOM standard nespecifikuje[18]:
• Nárok shody se standardem v provedení detailů vlastností každé součásti.
• Nárok na shodu se standardem všech vlastností a funkcí zrealizovaného systému integrovaném ve skupině zařízení.
• Shodu provedení se standardem pomocí testování nebo jinými ověřovacími procedurami.
Pro potřeby mé diplomové práce jsem se seznámil pouze s prostředky pro uložení služeb, zejména formátu souborů. Formát souborů DICOM umožňuje ukládat dva typy dat, podobně jako u běžné digitální fotografie. Jedná se o textová data, tzv. metadata, kde jsou uloženy informace o pacientovi, identifikační číslo pacienta (rodné číslo), poloha pacienta, dále také informace o snímacím zařízení, vzdálenosti snímků, číslo snímku, jeho rozlišení, atd. Druhým typem jsou obrazová data, která obsahují jeden či více snímků v jednom souboru. Jelikož v prostředí Matlab® jsou již naprogramovány funkce pro načítání z DICOM souborů, jejich zápis a získávání metadata informací, nemusel jsem tento standard příliš podrobně studovat. Kompletní popis standardu DICOM je přiložen na CD (Příloha C.1).
Snímky, které jsem získal pro mou diplomovou práci z Krajské nemocnice Liberec, a.s. byly uloženy ve formátu DICOM. Každý snímek byl uložen v samostatném souboru, měl 12bitovou barevnou hloubku a jednalo se o šedotónové snímky s rozlišením 512×512 pixelů. Soubor se snímkem neměl příponu, obsahoval pouze název. V sadě snímků bylo od 100 do 550 řezů pacienta. Z důvodu ochrany osobních údajů byly informace o pacientech před předáním odstraněny.
3.3.1 Zobrazení DICOM snímků
Pro zobrazení snímků jsem použil program MxLiteView verze 1.24 od Philips Medical Systéme Inc., který byl přiložen na CD u získaných snímků. Tento program umožňuje načtení snímků, jejich prohlížení, zobrazení metadata informací, měření rozměrů objektů na snímku, export souborů do obrazového formátu BMP, atd. Hodnoty jasů snímků v tomto programu se pohybovaly od –1024 až po 3071, celkem tedy 4096 hodnot jasů. V prostředí Matlab® jsem pracoval s hodnoty jasů od 0 do 4096, záporné hodnoty neexistují. Záporné hodnoty jasů v programu MxLiteView mají zřejmě souvislost s Hounsfieldovým číslem, jelikož jde o program pro účely medicíny.
Aby bylo možné CT snímek zobrazit na obrazovce počítače, musíme ho transformovat z 12bitové do 8bitové hloubky (256 jasů šedi). To provedeme vytvořením tzv. LUT tabulky (LookUp Table), pomocí které transformujeme hodnoty jasů jednotlivých pixelů. Výsledný obrázek má již 8bitovou hloubku. Převodní LUT tabulka má dva řádky. V prvním řádku je hodnota jasu pixelu v původním 12bitovém snímku a v druhém řádku je odpovídající hodnota jasu pixelu ve výstupním 8bitovém snímku.
Dochází tak ke ztrátě informace, ale pro zobrazení snímku je to dostačující.
Pro medicínské účely se používají různě převodní LUT tabulky. Například pro zobrazení plic, mozku, kostí, atd. Jako parametr převodu je zde zadána šířka a střed
okna. Tím jsou dány dvě hodnoty jasů, neboli prahy. Pixely mající hodnotu jasu nižší než první práh jsou černé a pixely mající hodnotu jasu větší než druhý práh budou na výsledném snímku bílé. Hodnoty jasů mezi prahy budou lineárně rozloženy do 0 až 255 jasů šedi. Pro zobrazení kostí je v programech pro medicínské účely nastavena šířka okna na 1800 hodnot jasů a střed okna na 600. jas. Odpovídající prahy jsou tedy –300 a 1500. Pro zobrazení plic je šířka okna 1600 hodnot jasů a střed okna –600. jas. Prahy jsou v tomto případě –1024 a 200. Dále existují parametry pro zobrazení páteře, mozku, tuku atd. Grafické zobrazení LUT tabulky pro zobrazení kostí a plic je zobrazeno na obrázku 3.2.
Obrázek 3.2: Grafické zobrazení LUT a) pro zobrazení kostí b) pro zobrazení plic v programu MxLiteView verze 1.24
3.4 Histogram
Histogram je graf, který nám umožní získat představu o rozdělení jasových úrovní v obraze (Obrázek 3.3). Jedná se o grafické zobrazení vektoru jasu hf (zi), i = 0, …, L-1, který má počet složek rovný počtu jasových úrovní L. Hodnota každé složky odpovídá četnosti bodů příslušného jasu v obraze [1].
Obrázek 3.3: a) Histogram b) Vyhlazený histogram mediánem pro velikost okolí 15
a) b)
a) b)
Využívá se v metodách při zpracování obrazu, např. při segmentaci. Histogram digitalizovaného obrazu obsahuje řadu lokálních maxim a minim, díky nimž vznikají problémy při jeho analýze. Proto se většinou vyhlazuje filtrem, který potlačí jeho lokální poruchy. Nejjednodušší typ filtru pro vyhlazení je klouzavý průměr či medián.
3.5 Filtrace obrazu
Filtrace obrazu je transformace, která se často užívá při předzpracování obrazu.
Filtrace převádí hodnoty jasu vstupního obrazu na jiné hodnoty výstupního obrazu. Při této operaci se k výpočtu nové hodnoty jasu pixelu využívají malá okolí O právě zpracovávaného pixelu. Cílem filtrace je zvýraznit nebo potlačit některé jeho vlastnosti.
Podle cíle rozdělíme lokální předzpracování do dvou skupin:
• Vyhlazování, které vede k potlačení vyšších frekvencí obrazové funkce.
Výsledkem je odstranění náhodného šumu v obraze.
• Gradientní operace, které naopak zvýrazňují vyšší frekvence obrazové funkce. Výsledkem je zvýraznění hran v obraze, ale i šumu, což bývá nežádoucí.
Vidíme, že vyhlazování a gradientní operace jsou (alespoň ve své lineární podobě) v protikladu. Řešení tohoto rozporu umožňují některé nelineární metody, které vyhlazují a zároveň jsou šetrné k hranám a detailům v obraze. Lineární operace počítají hodnotu ve výstupním obraze g(x,y) jako lineární kombinací hodnot vstupního obrazu f v malém okolí reprezentativního pixelu (x,y)[1]. Příspěvek jednotlivých pixelů v okolí O je vážen koeficienty h podle vztahu (3.5.1).
(3.5.1) Na vztah (3.5.1) lze pohlédnout jako na diskrétní konvoluci s konvolučním jádrem h. Jádru h se také říká konvoluční maska.
Musíme si uvědomit, že v průběhu předzpracování obrazu nezískáme žádnou novou informaci. Předzpracováním můžeme jen některou informaci potlačit nebo zvýraznit.
Filtrů pro vyhlazování nebo detekci hran a čar existuje velké množství, uvedu zde pouze ukázku některých z nich.
∑ ∑
∈⋅
−
−
=
) , (
) , ( ) , ( )
, (
n
m O
n m f n y m x h y
x g
3.5.1 Vyhlazování obrazu
Mezi základní metody vyhlazování obrazu patří prosté průměrování, kde je nový jas každého bodu vypočten aritmetickým průměrem původních jasů ve zvoleném okolí.
Nevýhodou této metody je rozmazávání hran v obraze.
Další jednoduchou metodou je filtrace mediánem. U této metody se nový jas každého bodu vypočte určením mediánu ve svém lokálním okolí. Medián se určí tak, že uspořádáme vzestupně hodnoty jasu v lokálním okolí a vybereme prvek, který je uprostřed této posloupnosti. U této metody nedochází k tak velkému rozmazávání hran jako u metody prostého průměrování a dobře potlačuje impulsní šum.
Mezi známé metody vyhlazování patří také metoda rotující masky. U této metody rotuje malá maska kolem bodu zájmu a v každé masce je vypočten rozptyl jasů. Za homogenní okolí bodu zájmu se vybere maska s nejmenším rozptylem a nový jas se určí jako aritmetický průměr původních jasů této masky. Tato metoda nerozmazává hrany a má mírně ostřící charakter.
3.5.2 Hranové detektory
Hrana je vektorová veličina a je určena velikostí a směrem[3]. Velikost a směr vychází z gradientu obrazové funkce v obrazovém elementu. Hrana indikuje body v obrazu, ve kterých dochází ke změnám obrazové funkce.
Hranové operátory můžeme rozdělit do dvou skupin dle způsobu hledání hran:
• U první skupiny je derivace obrazové funkce aproximována pomocí diferencí realizovaných diskrétní konvolucí.
• U druhé skupiny se hrany hledají v místech, kde druhá derivace obrazové funkce prochází nulou.
Do první skupiny můžeme zahrnout Robertův, Sobelův, Laplaceův, Robinsonův a Kirschův operátor a operátor Prewittové. Do druhé skupiny můžeme zahrnout Marrovu teorii a Cannyho hranový detektor.
3.6 Matematická morfologie
Matematická morfologie vychází z vlastností bodových množin a s jejíchž pomocí lze zpracovávat či předzpracovávat binární obrazy, které mají pouze dvě barvy, ale také obrazy s více úrovněmi jasu. Metody matematické morfologie nacházejí své uplatnění především při odstraňování šumu, zjednodušení tvaru a zdůraznění struktury
objektů, jako je např. ztenčování a zesilování, ale také pro popis objektů číselnými charakteristikami, jako je např. obvod a plocha objektu.
Hlavní výhodou těchto metod je jejich krátký čas na zpracování a snadná realizovatelnost. Uvedu zde pouze základní operace matematické morfologie, vzhledem k jejich detailnímu popisu v [1].
Dilatace ⊕ skládá body dvou množin pomocí vektorového součtu (Minkovského množinový součet). Dilatace X⊕B je bodovou množinou všech možných vektorových součtů pro dvojice pixelů, vždy jeden z množiny X a jeden z množiny B. Dilatace je definována dle předpisu (3.6.1).
(3.6.1).
Eroze \ skládá body dvou množin pomocí vektorového rozdílu (Minkovského množinový rozdíl). Eroze je duální operací k dilataci, ale ani dilatace, ani eroze nejsou invertovatelné. Eroze \ skládá dvě množiny dle předpisu (3.6.2).
(3.6.2)
Otevření ○ jde o operaci, při které je eroze následovaná dilataci. Otevření množiny X strukturním elementem B se označuje X○B a je definováno dle vzorce (3.6.3). Tato operace oddělí objekty spojené a odstraní detaily v obraze, které jsou menší než strukturní element.
(3.6.3)
Uzavření ● jde o operaci, při které je dilatace následovaná erozí. Uzavření množiny X strukturním elementem B se označuje X●B a je definováno dle vzorce (3.6.4).
Uzavření spojí objekty, které jsou blízko u sebe, zaplní malé díry a zaplní úzké zálivy.
(3.6.4)
Top-hat transformace, nebo-li vrchní část klobouku. Jedná se o množinový rozdíl mezi otevřením a původním obrazem (3.6.5). Části v obraze, které se nevejdou do strukturního elementu B použitého pro otevření, se odstraní. Po odečtení otevřeného obrazu od původního zůstanou právě jen odstraněné části obrazu.
(3.6.5)
{
d E d x b x X b B}
B
X⊕ = ∈ 2: = + , ∈ , ∈
X\B=
{
d∈E2:d+b∈X pro∀b∈B}
X ○B (= X\B)⊕B
X ●B=(X⊕B)\B
X X
Y = −( ○ )B
Bottom-hat transformace, nebo-li spodní část klobouku. Jedná se o množinový rozdíl mezi původním obrazem a uzavřením (3.6.6). Tato operace zobrazí malé díry a mezery mezi objekty.
(3.6.6) 2D vyplnění děr, je morfologická operace, která se provádí na binárních i šedotónových obrazech. Pokud na binárním obrazu obklopují bílé pixely kolem dokola nějaké černé pixely. Pak tyto černé pixely se stanou bílými. Černé pixely představují tzv. díru a po provedení této funkce dojde k její zaplnění. Důležitá je volba velikosti okolí reprezentativního pixelu, pro kterou považujeme oblast za uzavřenou. Nejčastěji se používá 4-okolí a 8-okolí.
Odstranění hraničních objektů ve 2D, je morfologická operace, která se také provádí na binárních i šedotónových obrazech. Tato funkce odstraní takové objekty z obrazu, které se dotýkají jeho okraje. I zde je důležitá volba velikosti okolí reprezentativního pixelu, pro kterou považujeme oblast za uzavřenou.
3.7 3D metody zpracování obrazu
3D metody zpracování obrazu vychází z některých 2D metod zpracování obrazu. Tyto metody, jak již z názvu vyplývá, musí být aplikovány na 3D data, která reprezentují 3D modely objektů (např. srdce, kostí, těla, atd.). 3D metody mají podstatně větší časovou a výpočetní náročnost než 2D metody. Je proto důležité, aby byly správně a optimálně naprogramovány s ohledem na rychlost výpočtu.
3.7.1 3D filtrace
Trojrozměrná filtrace je jedna z metod, která je pouze rozšířením 2D metody.
Platí zde stejné vlastnosti jako ve 2D. Filtrace převádí hodnoty jasu vstupních dat na jiné hodnoty výstupního obrazu. Na rozdíl od 2D filtrace se k výpočtu nové hodnoty jasu využívá 3D okolí O právě zpracovaného pixelu. I cíle filtrace jsou shodné jako ve 2D. Podle cíle je dělíme na vyhlazování, pro potlačení např. šumu v obraze, a gradientní operace, které zvýrazňují např. hrany v obraze.
Zaměříme se pouze na lineární operace, u kterých se hodnota ve výstupním obraze g(x,y,z) počítá jako lineární kombinace hodnot vstupního obrazu f v malém okolí reprezentativního pixelu (x,y,z). Příspěvek jednotlivých pixelů v okolí O je vážen koeficienty h podle vztahu (3.7.1).
(3.7.1) X
Y (= ●B) −X
∑ ∑
∈⋅
−
−
−
=
) , , (
) , , ( ) , , ( )
, , (
n m
k O
n m k f n z m y k x h z
y x g
Na vztah (3.7.1) lze pohlédnout jako na diskrétní konvoluci s konvolučním jádrem h. Konvoluční jádro h je na rozdíl od 2D filtrace trojrozměrná matice. I zde můžeme použít vyhlazovací filtry prostým průměrováním, mediánem, atd.
3.7.2 Vyplnění děr 3D matice
Jednou z užitečných funkcí ve 3D zpracování obrazu je vyplnění děr v 3D matici. Této funkce může být využito například u matice, která obsahuje pouze binární informaci (jedničky a nuly). Pokud jedničky kolem dokola obklopují nějaké nuly, pak tyto nuly budou přepsány na jedničky. Nebo-li tím dojde k tzv. vyplnění díry. Závisí na volbě velikosti okolí reprezentativního pixelu, pro které považujeme oblast za uzavřenou. U 3D morfologických operací se používá nejčastěji 6-okolí, 18-okolí nebo 26-okolí.
3.7.3 Odstranění hraničních objektů ve 3D
Tato funkce odstraní objekty, které se dotýkají okraje 3D matice. Zde je také nutná volba velikosti okolí pro uzavřenost hranice. Tato funkce je velice užitečná například pro tvorbu a další zpracování 3D matice. Pokud bychom měli matici, která obsahuje například jeden model páteřního obratle a dva neúplné modely páteřních obratlů. Předpokládáme, že neúplné modely jsou na okraji matice a pro další zpracování jsou nepoužitelné. Při použití funkce odstranění hraničních objektů by došlo k odstranění těchto neúplných obratů a v matici by tak zůstal jenom jeden úplný model páteřního obratle.
V mých postupech jsem se této funkce pokoušel využít, ale narazil jsem na problém, kdy objekty byly spojeny na více místech několika málo pixely a došlo tak k odstranění všech objektů, což je nežádoucí.
4 Metody segmentace CT snímků
Segmentace je jedna z nejnáročnějších metod analýzy obrazu, která vede k nalezení různých objektů v obraze (např. kostí, svalů, osob, zvířat nebo také poznávací značky automobilu). Výsledkem segmentace je rozdělení do částí, které jsou oblastí zájmu v dalším zpracování. Tuto analýzu provádí lékař, který si prohlédne snímek a následně určí jednotlivé tkáně, orgány či jiné objekty zájmu.
Díky současným výkonům počítačů mohou být aplikovány i složité algoritmy segmentace a následná 3D vizualizace segmentovaných oblastí zájmu. Tato vizualizace je velice vhodná pro představení komplexních objektů, které si již lékař nedokáže tak jednoduše představit. Cílem je navrhnout algoritmus, který dokáže rozlišit oblast zájmu podobně jako člověk. Jde o rozdělení obrazových bodů snímku na body zájmové oblasti a body ostatních oblastí. Segmentace je hlavním bodem při tvorbě 3D modelů kostí.
Segmentaci můžeme rozdělit dle několika hledisek. Podle míry zásahu uživatele jí dělíme na automatickou segmentaci bez zásahu uživatele, poloautomatickou segmentaci, kde uživatel může ovlivňovat proces segmentace změnou parametrů či drobných úprav výsledků, a manuální segmentaci, ve které uživatel segmentaci ovlivňuje zcela sám pomocí svého analyzování snímku. Dále můžeme segmentaci rozdělit dle typu zpracovávaných dat na 2D segmentaci, kde se segmentace provádí analýzou 2D snímku, a 3D segmentaci, kde jsou analyzována data 3D snímku, například u CT nebo MR sada 2D snímků za sebou.
Výsledek segmentace je pak reprezentován většinou hranicí objektu, která je velice vhodná pro klasifikaci a kontrolu správnosti segmentované oblasti uživatelem, nebo také barvením objektu, která je pak vhodná pro tvorbu 3D modelů (obrázek 4.1).
Tato kapitola patří tematicky do kapitoly 3 (Digitální obraz a postupy jeho zpracování). Z mého hlediska se jedná o významnou část práce, ve které se snažím přinést nové metody v segmentaci kostí. Jelikož je mnoho algoritmů segmentace, v této kapitole uvedu jejich stručný přehled. Zejména zde uvedu metody adaptivního prahování, které jsem nastudoval, implementoval a experimentálně testoval v prostředí Matlab®. Některé z metod jsem testoval i v prostředím National Instruments™ Vision Assistant verze 8.5.
Obrázek 4.1: Formy znázornění segmentace
4.1 Problémy při segmentaci CT snímků
Jedním z problémů při segmentaci je případný výskyt různých stínů v snímku, které vznikají pohybem pacienta nebo výskytem kovových objektů. Takto znehodnocený snímek je velice složité nasegmentovat, proto je vhodné tyto snímky vyloučit již před segmentací. Případně zvolit jinou metodu jejich získávání, či ji upravit.
Například moderní CT dokáží získávat snímky pod různými úhly nakloněním prstence (gantry). Tím je možné se vyhnout určitým oblastem s výskytem kovových objektů.
Čehož se využívá při pořizování snímků hlavy. Vyhneme se oblasti úst, v které mohou být zubní plomby.
Hlavním problémem při segmentaci kostí je jejich velká rozmanitost ve tvarech, velikostech, polohách a intenzitách jasu kostní tkáně.
4.2 Prahování
Prahování je nejednoduší segmentační metodou, která při segmentaci využívá statistickou analýzu obrazových dat. Mnoho objektů se velmi často liší svoji intenzitou jasu od intenzity jasu pozadí. Rozdělení úrovní jasů v obraze tak může pomoci určit práh pro převod obrazu do binární reprezentace, neboli oddělit objekty od pozadí.
Tato transformace vstupního obrazu f na výstupní binární obraz g je dána vztahem,
(4.2.1)
⎩⎨
⎧
<
= ≥
T f(x,y)
T f(x,y) g(x,y)
pro 0
pro 1
kde T je předem daná konstanta, která se nazývá práh. Ve výsledném obraze mají pixely náležící zkoumanému objektu hodnotu 1 a pixely náležící pozadí mají hodnotu 0.
Touto metodou postupně testujeme všechny pixely v obraze.
Pro správný výsledek je velice důležitá optimální volba prahu, kterou lze určovat pokusně nebo např. pomocí analýzy histogramu. Na obrázku 4.2 je zvolen práh mezi dvěma vrcholy histogramu, kde jeden vrchol může reprezentovat pixely pozadí a druhý pixely objektu.
Obrázek 4.2: Zvolený práh pomocí analýzy histogramu
Rozšíření této metody je prahování s více prahy. Výsledkem tedy není pouze binární obraz, ale obraz s vyšším počtem jasových úrovní, které je dáno počtem prahů (N prahům odpovídá N+1 jasových úrovní výstupního obrazu). Modifikací prahování je poloprahování, které se užívá pro vizuální hodnocení výsledků člověkem (Obrázek. 4.3). Získáme tak obraz, v kterém je odstraněné pozadí, ale objektům byla zachována jejich hodnota jasů.
Obrázek 4.3: a) originální snímek z CT b) výsledný obraz vzniklý poloprahování
a) b)
Práh, který je určen z celého obrazu, se nazývá globální. Většinou ovšem nemůžeme použít stejný práh na celé ploše obrazu. Může to být způsobeno změnami hodnot jasů objektů i pozadí, nerovnoměrností osvětlení či nestejnými vlastnostmi snímacího zařízení v celé ploše obrazu. Možným odstraněním těchto problémů je adaptivní prahování.
4.3 Adaptivní prahování
U adaptivního prahování je obraz rozdělen do několika oblastí. Hodnota prahu je určována pro jednotlivé oblasti obrazu a je funkcí jejich lokálních parametrů. Obecně můžeme říci, že práh je počítán pro každý pixel originálního obrazu. Výsledný obraz je dán vztahem 4.3.1. Tato metoda může lépe vyhovovat při změně světelných podmínek v obrazu (např. vzniku stínů, intenzivního osvětlení, atd.).
Definice prahové hodnoty T na vstupním obrazu f může být dána vztahem, (4.3.1) kde p(x,y) je určitá lokální vlastnost daného bodu. Výpočet vlastnosti je obvykle počítán v okolí bodu, aby bylo možné vzít v úvahu vliv osvětlení a šum. Vlastnost může být například průměrná hodnota nebo medián v předdefinovaném okolí, jehož je daný bod středem.
Pro hledání adaptivního prahu se využívá několik postupů, mezi které patří metoda Chow a Kaneko, algoritmus Otsu, „lokální prahování“ a algoritmus Niblack.
Předpoklad těchto metod je, že malé obrazové oblasti budou mít přibližně rovnoměrné osvětlení. Výsledkem je lepší nalezení prahu.
4.3.1 Metoda Chow a Kaneko
Metoda od C. K. Chow a T. Kaneka dělí obraz do pole překrývajících se podobrazů. Pro každý podobraz je analýzou histogramu nalezen optimální práh. Práh pro jednotlivé pixely je vypočten interpolací nalezených prahů jednotlivých podobrazů.
Nevýhoda této metody je výpočetní náročnost, a proto není vhodná pro aplikace pracující v reálném čase. Složité je také vytvoření optimálního algoritmu pro automatické určování prahu, které je vzhledem k rozličnosti velikostí a tvarů kostí velmi složité. Tato metoda je spíše vhodná pro segmentaci textu.
4.3.2 Metoda Otsu (metoda optimálního prahu)
Tato metoda vychází z relativního histogramu a základem jsou statistické výpočty. Předpokládáme existenci dvou tříd C0, která reprezentuje body objektů, a C1,
)]
, ( ), , ( [
,y) T f x y p x y
T(x =
která reprezentuje body pozadí. Tyto třídy jsou rozděleny prahem T, který dělí body na patřící do třídy C0, které jsou v intervalu 〈0,T–1〉, a třídy C1, které jsou v intervalu
〈T,N–1〉. Hodnota N odpovídá počtu jasových úrovní obrazu.
Pravděpodobnost výskytu bodu ze třídy C0 je dána vztahem[14]
(4.3.2) a adekvátně označíme i pravděpodobnost výskytu bodu ze třídy C1 vztahem (4.3.3).
(4.3.3) V obou vzorcích je pi i-tá hodnota z relativního histogramu. Celková pravděpodobnost výskytu bodu obou tříd je dána jejich součtem.
Střední hodnota pro třídu C0 a C1 je dána vztahy [14]
(4.3.4) kde ω(T) je pravděpodobnost výskytu bodu do jasu T-1, μ(T) je střední hodnota bodů do jasu (T-1) a μcelk.je celková střední hodnota celého originálního obrazu, která je dána vztahem (4.3.5).
(4.3.5) Pro třídu C0 a C1 dále definujeme rozptyl, který je dán vztahy (4.3.6).
(4.3.6) Nyní zavedeme tzv. mezitřídový rozptyl σB2 (4.3.7)[14], díky kterému se pomocí nalezení maxima určí hodnota prahu histogramu.
(4.3.7) Dosazením do vzorce za μcelk, vztah (4.3.5) můžeme rovnici upravit do tvaru (4.3.8).
(4.3.8) Dalšími úpravami při použití vztahů (4.3.2), (4.3.3) a (4.3.4) dostaneme výsledný vzorec ve tvaru (4.4.9) [14].
(4.3.9)
∑
−=
=
= 1
0
0 ( )
T i
i T
p ω
ω
∑
−=
−
=
= 1
1 1 ( )
N T i
i T
p ω
ω
∑
∑
−=
−
= −
= −
=
=
= 1
1 1
1
0 0
0 1 ( )
)
; ( ) (
)
( N
T i
celk T i
i i
T T ip
T T ip
ω μ μ
μ ω ω μ μ ω
∑
−=
⋅ +
⋅
=
= 1
0
1 1 0 0 N
i i
celk ip ω μ ω μ
μ
∑
∑
−=
−
=
⋅
= −
⋅
= − 1
1 1 2 2
1 1
0 0
2 2 0
0
)
; ( )
( N
T i T i
i
i i p
p i
ω σ μ
ω σ μ
(
0)
2 1(
1)
20 2
celk celk
B ω μ μ ω μ μ
σ = ⋅ − + ⋅ −
(
1 0)
21 0
2 ω ω μ μ
σB = ⋅ ⋅ −
)]
( 1 [ ) (
)]
( ) (
[ 2
2
T T
T
celk T
B ω ω
μ ω
σ μ
−
⋅
−
= ⋅
Hledaný práh T dostaneme dosazením do výsledného výrazu (4.3.9) pro všechna T∈〈0, N-1〉 a nalezením T, pro které bude výraz nabývat maximální hodnoty.
Nalezený práh se nachází v místě histogramu, kde vzdálenost středních hodnot tříd C0
a C1 je maximální a dochází k optimálnímu oddělení objektů a pozadí.
Pokud použijeme pouze dvě třídy, pro objekty a pozadí, není tato metoda vhodná pro segmentaci kostí z CT snímků. Pomocí dvou tříd by došlo k vyseparování pouze pozadí a lidského těla na CT snímku. Rozšířením této Otsu metody je tzv. Multi Otsu metoda, u které můžeme definovat vyšší počet tříd než dvě. Což je pro segmentaci kostí již vhodnější. Porovnání je zobrazeno na obrázku 4.4. Použijeme-li například tři třídy, dojde k vyseparování pozadí, měkké tkáně a kostí. Výsledný obraz není již binární. Počet jasů ve výsledném obraze odpovídá počtu tříd. Jelikož kosti a měkké tkáně mají velký rozsah jasů, není automatické určení prahů pro segmentaci kostí příliš přesné. Zůstávají zde nespojité hranice kostí.
Obrázek 4.4: a) snímek z CT b) nasegmentovaný snímek pro dvě třídy c) nasegmentovaný snímek pro tří třídy
Tuto metodu vymyslel N. Otsu v roce 1979, po kterém je pojmenována.
V některých publikacích je tato metoda nazývána klastrování (Clustering) nebo InterVariance[13].
4.3.3 Lokální prahování
Dalším způsobem jak určit lokální práh je statisticky vyhodnocovat hodnoty intenzity jasu okolí každého bodu. Volba nejvíce vhodné statistiky závisí velkou měrou na povaze vstupního obrazu. Například pokud obraz obsahuje velký rozdíl osvětlení, pak použití průměru může být velice účinné pro odstranění tohoto nežádoucího vlivu pro segmentaci.
a) b) c)
