Budapesti Műszaki Főiskola Neumann János Informatikai Kar Szoftvertechnológia Intézet
TUDOMÁNYOS DIÁKKÖRI DOLGOZAT
VASTAGBÉL DAGANATOK DIAGNOSZTIKÁJÁNAK AUTOMATIZÁLÁSA
Szerzők:
Bándi István mérnők informatikus szak, IV. évf. Reményi Attila mérnők informatikus szak, IV. évf.
Konzulensek:
Kozlovszky Miklós ügyvivő szakértő Vámossy Zoltán főiskolai docens
Budapest, 2009.
1 Tartalomjegyzék 1
Tartalomjegyzék ............................................................................................................. 1
2
Bevezetés ........................................................................................................................ 4
3
4
2.1
Problémaleírás ........................................................................................................ 4
2.2
Célkitűzés................................................................................................................ 4
Biológia háttér ................................................................................................................ 6 3.1
A vastagbél ............................................................................................................. 6
3.2
A vastagbél nyálkahártya ........................................................................................ 6
3.3
A vastagbél betegségei [5] ...................................................................................... 8
3.3.1
Gyulladásos betegségek ................................................................................... 8
3.3.2
A diszplázia ................................................................................................... 10
Informatikai háttér ........................................................................................................ 12 4.1
Hasonló programok ............................................................................................... 12
4.2
Képfeldolgozó keretrendszerek ............................................................................ 12
4.2.1
Intel IPP ......................................................................................................... 12
4.2.2
OpenCv .......................................................................................................... 12
4.3
Masszívan párhuzamos programozás alternatívái ................................................ 13
4.4
CUDA (Compute Unifed Device Architecture) [7] .............................................. 15
4.4.1
CUDA alapjai ................................................................................................ 15
4.4.2
Programozási modell ..................................................................................... 15
4.4.3
Kernelek és szálak ......................................................................................... 16
4.4.4
Memória hierarchia........................................................................................ 18
4.4.5
CUDA API .................................................................................................... 19
5
Probléma megoldásának megközelítése ....................................................................... 20
6
Program architektúrája ................................................................................................. 21 6.1
Digitális metszet kezelés ....................................................................................... 23
6.1.1 6.2 7
A digitális metszetek kezelése a gyakorlatban .............................................. 24
Cuda használata .Net-ben ..................................................................................... 27
Sejtmagdetektálás ......................................................................................................... 29 7.1
A probléma ismertetése ........................................................................................ 29
7.2
Konvertálás HSV modellbe [23]........................................................................... 29
7.3
Az algoritmus ........................................................................................................ 31 Oldal 1
7.4
Ultimate Erode ...................................................................................................... 31
7.5
Sejtmagok paramétereinek mérése ....................................................................... 34
8
7.5.1
Kitöltő algoritmus .......................................................................................... 34
7.5.2
Körvonal detektáló ........................................................................................ 35
7.5.3
Átmérők megállapítása [21] .......................................................................... 36
Mirigydetektálás ........................................................................................................... 38 8.1
A probléma megismerése: .................................................................................... 38
8.2
Az algoritmus ........................................................................................................ 38
8.2.1 9
Háló készítése [21] ........................................................................................ 39
Hámdetektálás .............................................................................................................. 41 9.1
A probléma megismerése ...................................................................................... 41
9.2
A detektáló algoritmus .......................................................................................... 41
9.3
Sejtmagok differenciálása ..................................................................................... 42
10
Képfeldolgozás gyorsítása GPGPU segítségével ..................................................... 43
10.1 Alapok, használat .................................................................................................. 43 10.1.1
Futásidejű API ............................................................................................... 43
10.1.2
Meghajtó API ................................................................................................ 46
10.1.3
A szegmentáló kernel .................................................................................... 47
10.1.4
A teljesítménymérő alkalmazás ..................................................................... 48
10.1.5
A zajszűrő kernel ........................................................................................... 50
10.2 Globális memória .................................................................................................. 52 10.2.1
Beépített struktúra használata ........................................................................ 54
10.2.2
Négy bájt olvasása, három bájt írása ............................................................. 55
10.2.3
Textúra alkalmazása ...................................................................................... 56
10.2.4
12 bájt ............................................................................................................ 57
10.3 Osztott memória .................................................................................................... 58 10.3.1
Redundancia megszüntetése .......................................................................... 59
10.3.2
Adatlehívás sorrendjének megváltoztatása .................................................... 59
10.4 Optimalizálás számításra ...................................................................................... 60 10.5 Egyéb algoritmusok .............................................................................................. 60 11
Paraméterek mérése .................................................................................................. 62
11.1 Megvalósított paraméterek ismertetetése .............................................................. 62 12
Tesztelés ................................................................................................................... 63 Oldal 2
13
12.1.1
Sejtmagok detektálása ................................................................................... 63
12.1.2
Mirigy detektálása ......................................................................................... 65
Értékelés ................................................................................................................... 68
13.1 GPGPU ................................................................................................................. 68 13.1.1
Kernelek összehasonlítása ............................................................................. 68
13.1.2
GPU – CPU összehasonlítása ........................................................................ 69
13.1.3
G80 – GT200 összehasonlítás ....................................................................... 72
13.2 Továbbfejlesztési lehetőségek .............................................................................. 73 13.2.1
GPGPU .......................................................................................................... 74
14
Ábrajegyzék .............................................................................................................. 75
15
Irodalomjegyzék ....................................................................................................... 78
Oldal 3
2
Bevezetés
2.1 Problémaleírás Napjainkban a szövetminták mikroszkopikus vizsgálata túlnyomórészt manuálisan történik. A kutatónak, vagy orvosnak a mikroszkóp fölé kell hajolnia, egy kényelmetlen testhelyzetben, szemével a mikroszkópban végtelenbe kell fókuszálnia, majd a jegyzeteire tekintve közelre, ez a szemnek és az agynak egyaránt megterhelő. Egy minta manuális kielemzése, akár több 10 percet is igénybe vehet, így egy teljes mintasor elemzése fárasztó, valamint sok időt és energiát emészt fel. Arról nem is beszélve, hogy az orvos is ember, így pár órás munka után, a fáradtsági faktor miatt megnő a valószínűsége, hogy a mintából hibás diagnózist állít fel. A probléma egy másik olvasata a gazdaságosság. A költségek csökkentése érdekében a kórházak, kutatóintézetek sok esetben csökkentik a laboratóriumok, ill. a kutatók, orvosok számát. Ez a módszer a költségeket látszólagosan csökkenti, mert a túlterheltség miatt megnő a hibaarány. Egy laboratóriumban sokféle mintát elemeznek, ezek vizsgálata különböző módszert igényel, mi a vastagbél daganatok diagnosztikájára koncentrálunk. Ezúton szeretnénk köszönetet nyilvánítani a SOTE Sejtanalitika Laborjának, és a 3DHistech Kft-nek, különösen Molnár Bélának, Ficsor Leventének, Valcz Gábornak, Szabó Dánielnek és Jónás Viktornak a dolgozat és a szoftver létrejöttében nyújtott segítségükért.
2.2 Célkitűzés A problémára megoldást nyújt egy szövettani immundiagnosztikai kiértékelő és speciálisan a vastagbél daganatok differenciáldiagnosztikáját, prognosztikáját támogató szoftver. A program feladata a HE (Haematoxylin and Eosin) festést kapott digitalizált minták számítógépes feldolgozása a 3DHistech Kft. által több év óta fejlesztett MIRAX szoftver termékcsaláddal együttműködve. A feldolgozás során a szoftver a morfológiai paraméterek alapján azonosítja az egyes szövetkomponenseket: felszíni hámot, mirigyeket, sejtmagokat. A szöveti képek feldolgozását követően diagnosztikai elemzés készül a szegmentált szövettípusokra a diagnosztikai paraméterek alapján, célzottan az ép, gyulladásos, hyperplasztikus és daganatos elváltozások azonosítására. Célunk: a szövetkomponensek felismerését és szegmentációját támogató algoritmusok implementálása, a szoftver használatát segítő grafikus felhasználói felület tervezése, az indikátorként alkalmazható paraméterek definiálása, és kvantitatív értékeik meghatározása, melyek segítségével a szoftver a későbbiekben döntéstámogatást biztosíthatna az orvos felhasználók számára, az egyes paraméterek mérését megvalósító algoritmusok implementálása, a mért paraméterek kiértékelését végző eljárások kidolgozása, Oldal 4
az ép, gyulladásos, daganatos szöveti struktúrák megkülönböztetésére. A szoftver képfeldolgozó algoritmusait, a gyorsabb futás érdekében, az nVidia által megalkotott CUDA technológiát felhasználva implementáljuk. Ez a technológia kiválóan alkalmas a számítás intenzív, jól párhuzamosítható feladatok, így a képfeldolgozás felgyorsítására. A CUDA-val kihasználható a számítógépekben lévő, eddig leginkább megjelenítésre használt GPU-k számítási kapacitása. Ezzel a program futása gyorsabb, és (a CPU tehermentesítése miatt) használata kényelmesebb lesz. A szoftver kevésbé időkritikus részeit és a grafikus felületet C#-ban, a .Net keretrendszert felhasználva valósítjuk meg.
Oldal 5
3 Biológia háttér 3.1 A vastagbél A vastagbél az emésztőrendszer utolsó szakasza, átmérője átlagosan 5-6 cm, hossza 150180 cm [14]. Bélrendszer feladatai a táplálék emésztése, felvétele, felszívása és az emészthetetlen salakanyag külvilágra juttatása. A vastagbélrák az egyik leggyakoribb daganatfajta, Magyarországon évente 7000 – 8000 új megbetegedést diagnosztizálna, amivel a 2. leggyakoribb daganat a tüdőrák után. A vastagbél szöveti szerkezetét az alábbi ábra szemlélteti.
1. ábra Vastagbél szöveti szerkezete [18]
Mivel a vastagbél daganat a nyálkahártyából – mucosa – fejlődik ki, ennek a felépítését részletesebben is megvizsgáljuk.
3.2 A vastagbél nyálkahártya A vastagbél nyálkahártyája az alábbi rétegekből épül fel [18]:
hámsejt réteg
propria réteg: a nyálkahártya kötőszöveti rétege [16]
lamina muscolaris mucosae: nyálkahártya izomrétege
Oldal 6
2. ábra Vastagbél nyálkahártyájának szerkezete [18]
A hámréteg jellemzője, hogy a sejtek közötti távolság kicsi, szorosan egymás mellet helyezkednek el, megnyúlt hasáb alakúak, és sejtmagjaik egyvonalban vannak [5]. A nyálkahártyába betüremkedő „képződmények” mirigyek, melyek feladata, hogy nyák termelésével síkosítsák a vastagbél falát. A mirigyek alsó részén egy osztódó sejtréteg, a germinatív vagy generatív réteg található.
3. ábra Hámréteg, és egy mirigy **
Az egészséges nyálkahártya tulajdonságai [18]:
a nyálkahártya felszíne sima Oldal 7
a mirigyek sűrűsége 6 - 8 db / mm
a mirigyek mérete normális, a nyálkahártyáig érnek, irányultságuk párhuzamos, és nem ágaznak el
csak a mirigy alsó részén található osztódó sejtréteg
sejtjeik száma és felépítése normális
a kötőszövetben és az izomrétegben nincs elváltozás
3.3 A vastagbél betegségei [5] Egy egészséges vastagbél nyálkahártya sejt daganatos elváltozását több ok is előidézheti. Az elváltozás hatására a sejt bizonyos tulajdonságai megváltoznak, így képes lesz magát kivonni a saját, és a környező sejtek növekedési kontrolja alól. Ezek a sejtek új populációt hoznak létre, és a túlnövik a környezetüket, ún. polip vagy hiperplázia jön létre. A hiperplázia azt jelenti, hogy a sejtek normálisnál ugyan gyorsabban szaporodnak, de felépítésük nem változik. Ha a folyamat itt nem áll meg, diszplázia alakul ki, amikor már a sejtek szerkezetében is változások történnek. Később ebből alakulhat ki a rák. 3.3.1 Gyulladásos betegségek
A vastagbél gyulladásos betegségei a rák kialakulásának szempontjából fontosak, hosszabb ideig tartó gyulladás esetén megnő a rák kockázata. A gyulladás a hámsejtek szerkezetének megváltozásával jár, ami elősegíti a mutációt, és így áttételesen a rák kialakulását. A gyulladt nyálkahártya az alábbiakban különbözik az egészségestől:
a nyálkahártya felszíne megváltozik
a mirigyek sűrűsége csökken
a mirigyek párhuzamos, egyenes lefutása megváltozik
A nyálkahártya felületének simasága megszűnik, bizonyos helyeken kitüremkedések, máshol hámsejthiány, bemaródás alakul ki. A mirigyek sűrűsége az 7 - 8 db / mm-ről 4-5 db / mm-re csökken, a távolság közöttük megnő. Nemcsak a mirigyek helyzete, de a szerkezetük is drasztikusan megváltozik:
a mirigyek párhuzamossága megszűnik, átmérőjük megváltozik (4. ábra),
a mirigyek egyenes lefutása megszűnik, elágazások jönnek létre (5. ábra A mirigyekben elágazások jelennek meg),
az izomréteg és a mirigy alapja közötti távolság megnő, és mirigyenként különbözik (6. ábra Az izomréteg és a mirigyek alapja közti távolság megnő, és mirigyenként különbözik).
A gyulladás során a kötőszövetben – propira réteg – is elváltozások jelennek meg, a sejtek száma megnő, egészséges egyenletes eloszlásuk megváltozik. Ez minden vastagbél gyulladás esetén megjelenik, a fent felsorolt morfológiaváltozások nem feltétlenül!
Oldal 8
4. ábra A mirigyek párhuzamossága megszűnik
5. ábra A mirigyekben elágazások jelennek meg
Oldal 9
6. ábra Az izomréteg és a mirigyek alapja közti távolság megnő, és mirigyenként különbözik
3.3.2 A diszplázia
A diszplázia által okozott elváltozások túlnyomórészt megegyeznek a gyulladás által okozottakkal, de a diszplázia a gyulladástól függetlenül is kialakulhat. A diszpláziából közvetlenül is kialakulhat rákos daganat. A diszplázia az alábbi elváltozásokat okozza (a gyulladásnál már megismerteken kívül):
a hámréteg sejtmagjai megnyúlnak, sűrűn helyezkednek el nincsenek egy vonalban,
a sejtmagok megnagyobbodnak, alakjuk megváltozik,
a sejtmagok jobban festődnek,
a mirigyek generatív sejtrétege megvastagszik, és a sejt felsőbb részein is megjelenhet (7. ábra A generatív réteg a sejt felső részén is megjelenik).
Oldal 10
7. ábra A generatív réteg a sejt felső részén is megjelenik
Oldal 11
4 Informatikai háttér 4.1 Hasonló programok Mivel a feladat napjaink egyik megoldatlan problémája, így nehéz hasonló programokat találni. A HistoQuant nevű program célja részben megegyezik a mi feladatunkkal. Ez a program is tud sejtmagokat detektálni és ezek bizonyos jellemzőit is meg tudják mérni. A program azonban nincs ingyenes vagy próba verziója, így nem tudtam pontos információt szerezni milyen hatásfokkal képes detektálni.
4.2
Képfeldolgozó keretrendszerek
Napjainkban már nagyon sok képfeldolgozással foglalkozó nyílt és zárt forráskódú keretrendszer elérhető. Ezek újraimplementálása fölösleges és nagyon időigényes. A feladat megoldásához felhasználjuk az OpenCv-t és az Intel IPP-t. 4.2.1 Intel IPP
Integrated Performance Primitives az Intel egy többplatformos zárt forráskódú szoftverkönyvtár, amely kiemelten optimalizált funkciói fejlett alkalmazások kidolgozását teszik lehetővé. Az Intel Signal Processing Library (SPL), Image Processing Library (IPL), Intel JPEG Library (IJL) és az Intel Recognition Primitives Library (RPL) ezek funkcióit egybevonva született meg az IPP. Az IPP sok területhez ad fejlesztési eszközt. Ilyen területet a képfeldolgozás, gépilátás, hangfelismerés, adat tömörítés, kriptográfia, hangfeldolgozás, video kódolás stb. Mi azonban csak az első kettőt fogjuk használni. Az IPP kihasználja a legújabb processzorok nyújtotta technológiákat, mint a MMX, Streaming SIMD Extensions (SSE), Streaming SIMD Extensions 2 (SSE2), Streaming SIMD Extensions 3 (SSE3) technológia. Ajánlott az Intel architektúra használata, de a fent említetett technológiák nagy része megtalálható már az AMD processzorokban is. 4.2.2
OpenCv
Az OpenCV (Open Source Computer Vision) egy nyílt forráskódú keretrendszer. 1999ben kezdte el fejleszteni az Intel C, C++ nyelvben és 2006-ban adták ki az első verziót. Az OpenCV több platformon fut Windows, Linux, Mac OS X. Az Intel IPP-vel szoros kapcsolatban áll. Az OpenCV önmagában is ki tudja használni a többprocesszoros gépeket, de amennyiben az IPP telepítve van a számítógépre, ki tudja használni az IPP további előnyeit is. Az OpenCV több fő részből áll. CXCore tartalmazza az alapvető struktúrákat és ezeken értelmezett alapvető műveleteket. CVCam tartalmazza a kamerával kapcsolatos függvényeket. HighGUI modul az grafikus felülethez az eszközöket. Machine Learning modulban találhatók meg a neurális hálók és hasonló főleg osztályozással foglalkozó függvények. CV modul pedig a képfeldolgozás, képelemzéssel mozgásdetektálással foglalkozó függvényeket gyűjti össze. Mi főleg a CV előnyeit fogjuk kihasználni.
Oldal 12
4.3 Masszívan párhuzamos programozás alternatívái A masszívan párhuzamos architektúrák az asztali számítógépekben napjainkban olcsón elérhetőek. Három cég foglalkozik ilyen processzorok fejlesztésével: az nVidia, az Intel és az AMD. Ezen sorok írásakor azonban az Intel megoldása még nem elérhető, és kevés információ jelent meg róla. Az AMD és az nVidia is árul játékosoknak szánt videokártyák mellett professzionális felhasználásra szánt adatpárhuzamos gyorsítókat is, melyek jellemzője, hogy a központi egység teljesen ugyanaz, mint a játékosoknak szánt kártyákon, azonban általában több memóriával kaphatóak és nem található rajtuk video kimenet. A fejlesztőkörnyezeteket tekintve tágabb a választék, hiszen az egyes gyártók specifikus rendszerei mellett számos gyártófüggetlen megoldás is létezik. Az AMD megoldásai Az AMD első, asztali számítógépekben is elérhető, egységes árnyalókkal rendelkező architektúrája a 2007-ben piacra dobott R600-as volt. A cég a professzionális felhasználásra szánt kártyáknak a FireStrem nevet adta. Az R600-ra építve a következő években sorra jelentek meg a továbbfejlesztések: 2007-ben az RV670, 2008-ban pedig az RV770. [19] Az R600 mellé kiadott fejlesztőkörnyezet a Stanford Egyetemen kifejlesztett BrookGPU fordító és futtatókörnyezet egy (az AMD hardverekre optimalizált) változatát, a Brook+-t használja. A fejlesztőkörnyezet biztosít egy alacsonyabb szintű hozzáférést is a GPU erőforrásaihoz a CTM (Close To Metal) révén. [20] Az nVidia megoldásai [8] Az nVidia első kereskedelmi forgalomban kapható architektúrája mely programozható árnyalókkal rendelkezett, a 2006-ban bejelentett G80-as volt, és ezzel indult útnak a professzionális felhasználásra szánt Tesla termékvonal is. A következő lépés a 2008 elején forgalomba került G90-es volt, melyet a GT200 követett 2008 júniusában. A G80-nal egy időben jelent meg a hozzá tartozó fejlesztő környezet a CUDA (Compute Unified Device Architecture), melynek célja, hogy a hagyományos C nyelvet a lehető legkevesebb kiegészítéssel tegye használhatóvá az nVidia architektúrák programozására. Gyártófüggetlen megoldások A BrookGPU-t, ahogy már említettem, a Stanford Egyetemen fejlesztették ki, azzal a céllal, hogy hatékonyabbá és egyszerűbbé tegyék a GPU-ra történő fejlesztést. Ahogy a neve is mutatja, az ANSI C nyelv egy változatát, a Brook-ot használja, melyet kifejezetten adatpárhuzamos problémák egyszerűbb megközelítésére fejlesztettek ki. A fordító a Brook-ban írt kódot OpenGl, DirectX vagy az ATI féle CTM utasításokká fordítja, és azt futtatja az aktuálisan használt GPU-n. Tartalmaz emulátort is, mely a hibakeresést egyszerűsíti. Nagy előnye, hogy nyílt forráskódú, ingyenesen elérhető. [2] A RapidMind nem használ külön nyelvet, hanem a hagyományos C++ kódot elemezve próbálja a kódot párhuzamos architektúrákra optimalizálni. Képes kihasználni a sokmagos processzorokat, adatpárhuzamos gyorsítókat, sőt a Sony, Toshiba és az IBM által kifejlesztett Cell processzort is. Kutatási célokra ingyenesen elérhető. [4] A DirectCompute API-t a Microsoft fejleszti a DirectX11 részeként, mely ezen sorok írásakor még nem elérhető. Támogatni fogja a DirectX10 kompatibilis GPU-kat is. [10] Az OpenCL specifikációját 2008 végén hozta nyilvánosságra a Khronos Group. Céljai között szerepel, hogy hasonlóan a RapidMind-hoz, sokmagos processzorokra, GPU-kra és a Cell-re egy nyelvet és fejlesztőkörnyezetet használva lehessen hatékonyan fejleszteni. Az Oldal 13
OpenCL teljesen nyílt platform, bármelyik gyártó implementálhatja architektúrájára, így az egész iparágnak jelenthet egységes megoldást. [17]
a
saját
Megemlítendő még a PeakStream, mely a RapidMind-hoz hasonló rendszert hozott létre [3]. Érdekességképpen megjegyezném, hogy sem a PeakStream, sem a RapidMind nem önálló, előbbit a Google vásárolta fel 2007. június 5.-én [1], utóbbit az Intel 2009. augusztus 19.-én [4]. Véleményem szerint az ideális választás az OpenCL lenne, hiszen teljesen ingyenes és független, a Khronos Group-on keresztül mégis szinte minden jelentős informatikai cég, kivétel a Microsoft, támogatja. Azonban ez egy nagyon új szabvány, a projekt indulásakor még csak a specifikáció volt elérhető, ezért más alternatíva után kellett néznem. Választásom végül a CUDA-ra esett, ennek legfontosabb okai, hogy teljesen ingyen elérhető és programozási modelljének alapjai megegyeznek az OpenCL specifikációjában leírtakkal, így később egyszerűbb lesz az átállás.
Oldal 14
4.4 CUDA (Compute Unifed Device Architecture) [7] 4.4.1 CUDA alapjai
Az nVidia által fejlesztett CUDA egy fejlesztőkörnyezet, melynek segítségével kizárólag a cég saját architektúráit lehet programozni, a már elterjedt C nyelv kicsit módosított változatával. Egy CUDA kompatibilis GPU-nak több száz magja van, és ezek több ezer párhuzamos szálat képesek futtatni. A fejlesztés során az elsődleges szempont volt, hogy a C nyelvet a lehető legkevesebb módosítással alkalmassá tegyék a GPGPU-k programozására. Első verziója 2006 novemberében jelent meg, az első végleges, 1.0-ás verzió pedig 2007 júniusában került a fejlesztőkhöz. A programozás szempontjából a legfontosabb része a fordító. A fordítási folyamat kicsit különbözik a már megszokottól, hiszen a „.cu” kiterjesztésű forráskódot tartalmazó állomány két féle kódot is tartalmaz: a központi processzoron futtatandó ún. hoszt kódot, és a grafikai egységen végrehajtandó ún. eszköz (device) kódot. A fordító első lépésként szétválasztja a két féle kódot, és a CPU kódot átadja egy hagyományos C fordítónak, míg a GPU-n végrehajtandó kóddal maga foglalkozik. A fordítás során a forrás először egy köztes kódra fordul, ez a PTX, melynek szerepe hasonló a .Net által használt MSIL-hez: elfedi a tényleges használt architektúrák közötti különbséget, és garantálja, hogy ugyanaz a kód az összes architektúrán futtatható lesz. Fontos megjegyezni, hogy bár a PTX-ről natív kódra történő fordítás során a fordító optimalizálja a kódot, az egyes architektúrák valamelyest más követelményeket támasztanak a hatékony kóddal szemben. Egy CUDA kompatibilis processzor alapegysége a feldolgozás szempontjából a „Streaming Processor”, a továbbiakban SP. Nyolc darab SP van egy „Streaming Multiprocessor”-ba szervezve, erre a későbbiekben mint multiprocesszor vagy MP fogok hivatkozni. A hardver részletes ismertetése ennek a dolgozatnak nem célja, ezért csak a programozás szempontjából lényeges fogalmakat és tulajdonságok magyarázatára térek ki. Minden multiprocesszor tartalmaz egy ún. SIMT egységet, melynek feladata a szálak létrehozása, ütemezése és futtatása. Ez az egység a szálakat ún. „warp”-okba csoportosítja, egy warp 32 szálat tartalmaz. Egy warp az aktuálisan a multiprocesszor által futtatott szálakat jelenti. Az egy MP – egy SIMT egység felépítésből következik, hogy adott időpillanatban az SP-k, azonos utasítást tudnak csak végrehajtani, nincs minden SP-nek külön PC-je (Program Counter: a következő végrehajtandó utasításra mutat). Emiatt előfordulhat, hogy az amúgy párhuzamosan, függetlenül végrehajtható szálak között az elágazások miatt soros függőség alakul ki, ez lassítja a program futását. 4.4.2 Programozási modell
Programozási modell vagy paradigma alatt egy módszertant értünk, mely alkalmazásával az egyes problémákat képesek vagyunk számítógéppel leírni, ill. megoldani. A C nyelv a strukturált programozási paradigmát valósítja meg, mely Neumann féle soros feldolgozásra épít. Ezzel szemben a CUDA a párhuzamos programozási modellt implementálja. A két paradigma jelentősen különbözik, de a C nyelv néhány kiegészítéssel alkalmassá tehető a párhuzamos programozásra. A CUDA a C nyelvet az alábbi elemekkel egészíti ki:
függvény típusokkal, amikkel meghatározhatjuk, hogy egy függvény hol fusson, és honnan legyen hívható (CPU vagy GPU)
Oldal 15
változó típusokkal, amikkel meghatározhatjuk egy változó helyét a GPU memóriájában
a kernel hívási direktívájával
négy beépített változóval, amikkel meghatározhatók a grid és a blokk dimenziói, valamint a blokkok és a szálak indexei.
Ezzel elérhetjük, hogy a párhuzamosság csak a kernelek hívásakor jelenjen meg, a kernelhívások között, és a kerneleken belül is strukturált kódot írhatunk. A CUDA által alkalmazott párhuzamos programozási modell lényege, hogy egy szál a feldolgozni kívánt adatstruktúra egy vagy néhány elemét dolgozza fel, pl. egy kép esetében egy vagy néhány pixelt. A CUDA kódot futtatni képes GPU-k jellemzően több tucat, vagy akár több száz szálat képesek párhuzamosan futtatni. Az nVidia által használt terminológiában a hoszt a központi vezérlő egységet, míg a „device” – továbbiakban eszköz – a videokártyát vagy adatpárhuzamos gyorsítót jelenti. A kettő között a kapcsolatot a mai rendszerekben a PCI-Express 2.0 szabványú csatoló biztosítja, melynek elméleti sávszélessége 8 GB/s. A hoszt és az eszköz aszinkron működésű, egy kernelhívás után a vezérlés azonnal visszakerül a CPU-hoz, ezért nem blokkolja a program futását.
8. ábra Kernelhívás
4.4.3 Kernelek és szálak
A kernel egy olyan függvény, amit a CPU hív, de a GPU hajt végre, hívásakor meg kell adnunk, hogy hány párhuzamos szálon szeretnénk egyszerre futtatni. A kernelre az alábbi megszorítások érvényesek:
csak a GPU memóriáját érheti el
nem lehet visszatérési értéke (void) Oldal 16
nem lehet rekurzív
a statikus változók használata nem megengedett
A GPU szálak ún. „könnyű” szálak, létrehozásukhoz, és a szálak közti váltásokhoz szükséges idő jóval kevesebb, mint a CPU-n futatott szálak esetén, és kevésbé erőforrás igényesek. Egy kernel hívása az alábbi módon történik: kernelnnév<<
>>(paraméterek). A gridDim-mel és a blockDim-mel a szálak hierarchiáját határozhatjuk meg, a memSize-zal pedig az osztott memória méretét; ez opcionális. A szálak hierarchiáját az alábbi ábra szemlélteti:
9. ábra Szálhierarchia [6]
A szálakat ún. blokkokba rendezzük. Minden blokk ugyanannyi szálat tartalmaz, ezek számát és dimenzióit adjuk meg a blockDim-mel a kernel hívásakor. Egy blokkban a szálakat egy, kettő vagy három dimenzióba rendezhetjük, egy blokk a mai GPU-k esetén maximum 512 szálat tartalmazhat. Minden blokknak külön osztott memóriája van, így csak az azonos blokkban lévő szálak képesek ezen keresztül együttműködni. A kernelek egy gridben futnak, és egy GPU egy időben csak egy kernelt futtathat. Egy MP egy időben egy blokkot futtat, de tartozhat hozzá több blokk is, ilyenkor sorban futtatja őket. A blokkok és a blokkon belüli szálak száma kritikus a futási idő szempontjából, de számukat mindig az aktuálisan feldolgozni kívánt adat struktúrájához kell igazítani, és nem a hardverhez. A hardver osztja be, hogy melyik blokkot, melyik maghoz rendeli.
10. ábra Blokkok kiosztása a magokhoz különböző GPU-knál [6]
Oldal 17
4.4.4 Memória hierarchia
A CUDA-ban hat különböző memória területet különböztetünk meg, melyek elhelyezkedését és egymáshoz való kapcsolatukat az alábbi ábra szemlélteti.
11. ábra Az egyes memória fajták kapcsolata, elhelyezkedése [12]
A szálak által leggyorsabban a végrehajtó egységekkel közös félvezetőlapkán található regiszterek, és osztott memória érhető el, ezek késleltetése mindössze néhány ciklus. A regiszterek használata triviális, minden kernelben definiált változó az aktuális szál privát használatú regisztereibe kerül. Az osztott memória már sokkal nehezebben megfogható, hiszen a központi processzora írt programok esetében ezzel a fogalommal nem találkozhatunk. Ha mégis analógiát szeretnénk keresni, akkor a működése és használhatósága miatt a CPU-ban található gyorsító tárakhoz hasonlíthatjuk, azzal a lényeges különbséggel, hogy a GPU-ban a használata nem automatikus, a programozó dolga, hogy kihasználja az előnyeit. Logikailag az osztott memória a blokkokhoz tartozik, fizikailag viszont multiprocesszoronként van egy. Hatékony kihasználásának módjára a későbbiekben még visszatérek. A többi négy memóriafajta közös tulajdonsága, hogy fizikailag a GPU-val egy NYÁK-ra szerelt DRAM-ban találhatóak, emiatt a késleltetés nagyjából százszor akkora, mint a regiszterek és az osztott memória esetében. Ne feledjük azonban, hogy a videokártyák által biztosított sávszélesség többszöröse a központi egységnél megszokottnak. Memória órajel
Szélesség
Sávszélesség
nVidia GT200 [9]
1100 MHz
512 bit
141 GB / s
Intel Nehalem [15]
1066 MHz
192 bit
25,6 GB /s
1. táblázat GPU - CPU sávszélesség
A négy memória típus közül egyedül a lokális szálszintű, a többit az összes szál kezelheti. A lokális memória a regisztertér kiterjesztésének tekinthető, tipikusan akkor kerül alkalmazásra, amikor a kernelben egy tömböt definiálunk, amit aztán dinamikusan címzünk. Ennek oka az, hogy a szálakhoz tartozó privát regiszterek nem indexelhetőek. Oldal 18
Emiatt csak az olyan tömbök kerülhetnek közvetlenül a regiszterekbe, amelyek indexei már a fordítási időben eldönthetőek. A globális memóriát az összes szál írhatja és olvashatja, de gyorsító tár nem tartozik hozzá, emiatt a kezelésénél különösen ügyelni kell hatékonyságra. A konstans és a textúra memóriához ugyan tartozik egy gyorsító tár a GPU-ban, amin keresztül sokkal hatékonyabban elérhetők, de alkalmazhatóságukat behatárolja, hogy csak olvashatóak. Az utóbbi három memória típus kezelése – allokálás, felszabadítás – a hoszt kódon keresztül történik. Az egyes memóriák tulajdonságait az alábbi táblázat foglalja össze. Elhelyezkedés Gyorsító tár
Elérhetőség
Láthatóság
Élettartam
Regiszter
integrált
nem értelmezett
olvasható írható
/ 1 szál
szál
Lokális
külső
nincs
olvasható írható
/ 1 szál
szál
Osztott
integrált
nem értelmezett
olvasható írható
/ 1 blokk
blokk
Globális
külső
nincs
olvasható írható
/ összes szál + lefoglalástól CPU felszabadításig
Konstans
külső
van
csak olvasható
összes szál + lefoglalástól CPU felszabadításig
Textúra
külső
van
csak olvasható
összes szál + lefoglalástól CPU felszabadításig
2. táblázat Memória fajták összefoglalása [12]
4.4.5 CUDA API
A CUDA kétféle API-t biztosít a GPU programozására: a runtime, azaz futásidejű API-t, és a driver, azaz meghajtó API-t. A futásidejű API-t használva az eszköz kód egy állományba kerül a hoszt kóddal, míg a meghajtó API-t használva a kernelek egy külön fájlba kerülnek, majd a fordítás után ezeket be lehet tölteni a programba, és a paraméterek átadása után futtatni. Az utóbbit használva a kód jóval bonyolultabb lesz, de nagyobb programok esetén a használata elkerülhetetlen. Mindkettő működését és használatát a későbbiekben részletesen tárgyalom.
Oldal 19
5 Probléma megoldásának megközelítése Képfeldolgozó algoritmusokkal szeretnénk kinyerni az információt. A szöveti minta hatékony detektálása és kielemzése napjaink egyik megoldatlan problémája. Mivel a teljese probléma megoldása több éves kutatást igényel érdemes a programot úgy tervezni, hogy az később kis munka befektetésével módosítható, továbbfejleszthető legyen. A teljes problémát megpróbáljuk kisebb feladatokra osztani melyek minimálisan kapcsolódnak egymáshoz. Minden egyes feladathoz megfeleltetünk a programban egy modult. Mivel a minták sok szempontból (méret, vágási sík, festődés erősége, festődés szín stb.) nagyon eltérőek lehetnek, minden algoritmust úgy kell fejleszteni, hogy minél több típusú mintát tudjon detektálni, valamint az egyes kulcsfontosságú problémák (sejtmagdetektálás, mirigydetektálás, hámdetektálás, egyes bonyolultabb paraméterek mérése) esetén több algoritmust implementálunk. Amennyiben egy mintán nem megfelelő eredményt produkál az alapérdemezett algoritmus lehetőség van egy másik a mintához jobban alkalmazkodó algoritmus használatára. A digitalizált szöveteket úgynevezett digitális slide-okon tároljuk, ezek képezik a programunk bementét. Problémát okozhat a program fejlesztése során, hogy ezek a digitális slide-ok igen nagyméretűek lehetnek. Ezt a problémát több síkon kezeljük. A digitális slide betöltése után lehetőséget adunk a felhasználónak hogy kijelölje azt a területet, ami fontos az vizsgálat szempontjából, így redukálódik a vizsgálat terület. Az összes algoritmust fel kell készíteni, hogy végrehajtás előtt megvizsgálja a bemeneti adat méretét és amennyiben az túl nagy, képes legyen azt felbontani kisebb, a végrehajtás szempontjából hatékonyabb méretűekre bontani és ezeket egymás után végre tudja hajtani. A digitális slide-ok nagy mérete nagy futási időt igényel. A probléma nagy adatpárhuzamosságot tartalmaz, így ha kihasználjuk az adatpárhuzamosságot GPU segítségével a futási időt az eredeti töredékére csökkenthetjük.
Oldal 20
6 Program architektúrája A program architektúrájának megtervezéséhez tudni kell, hogyan végzik el ezt a feladatot az orvosok napjainkban. A minta elkészítése után az orvos behelyezi a mintát egy mikroszkópba. Minta festés során sérülhet, bizonyos esetekben ez lehetetlené teszi az minta pontos kiértékelését. Az első feladta ennek megállapítása. Meg kell vizsgálni a festődési szintjét. Amennyiben nem festődik kellőképpen áttetsző lesz a minta, ha túlságosan megfestődik az egész minta egyszínű lesz, így nem lehet elkülöníteni az egyes szöveti komponenseket. Amennyiben megfelelő a festődés szintje, meg kell vizsgálni a minta kimetszése során mennyire sérült a minta. Kisebb vágási hibák nagyon gyakran előfordulnak, így ezeket tudni kell kezelni későbbiekben a programnak is. A festődésen és a vágási hibákon kívül néha előfordulhat, hogy az igen vékony minta gyűrődik. A kisebb gyűrődést tartalmazó mintákat ki lehet értékelni, csupán a gyűrődési területet ki kell hagyni. Amennyiben az orvos meggyőződött arról, hogy megfelelő a minta megvizsgálja a szövetkomponenseket. Megvizsgálja a sejtmagok, mirigyek és a hám elhelyezkedését, mennyiségét, és alakjukat. Ezek az adatok orvos számára kellő adatot jelentenek, hogy meg tudja határozni a szövet állapotát. A program a hagyományos kiértékelési munkafolyamatra épül. A program első feladata a digitális minta betöltése. A programnak tudnia kell kezelni az alapvető kép kiterjesztéseket és a digitális slide-okat. Ezen feladatok elvégzésére valamint annak megvizsgálása, hogy a minta megfelelően ki lehet értékelni létrehozunk egy LOAD nevű modult. Megnyitás után a programnak meg kell vizsgálnia, hogy a minta megfelelően kiértékelhető, amennyiben lehetséges fel kell ismernie és meg kell találni az egyes szöveti komponenseket (sejtmag, mirigy, hám) a későbbi kiértékeléshez. A szöveti komponensek pontos detektálása egy nagyon bonyolult probléma így ez újabb önálló modult igényel (DETECTION). Feltéve hogy sikeresen detektáltuk a szöveti komponenseket meg kell határozni azok tulajdonságait. Ezt a feladatot az PARAMETER MEASUREMENT nevű modul fogja ellátni. Mint azt a hagyományos módszer tette, a kiszámolt paraméterek és ezekből készített statisztikai adatok alapján döntést kell hozni a minta állapotát illetően. A döntéshozó modul neve PARAMETER EVALUATION.
12. ábra Rendszerterv fejlődése I
Oldal 21
Ahhoz, hogy ez egy használható program legyen, létre kell hozni a GUI-t és a kimenetet kezelő modult (SAVE). A modulok végrehajtása soros mivel minden modul támaszkodik az azt megelőző modul eredményeire, de bizonyos esetekben nem csupán azokra. A paraméter mérésre szükségünk lesz a bementi képre is, így a LOAD és a PARAMETER MEASUEMENT modulok között direkt kommunikáció van.
13. ábra Rendszerterv fejlődése II
A program gyorsításához használni fogunk képfeldolgozó keretrendszereket, valamit Cuda –nyelvben megírt GPU függvényeket. Mivel a program C#-nyelvben készül és a képfeldolgozó keretrendszerek és GPU kód nem C#-ban íródott ezért ezekhez szükségünk lesz wrapperekre. A program legvégső állapotában teljesen automatizált módon fog működni, ez lehetővé teszi, hogy nagy adatmennyiséget önállóan fel tudjon dolgozni. A nagy adatmennyiség maga után vonja, hogy a programnak tudnia kell kezelni adatbázisokat.
Oldal 22
14. ábra Rendszerterv fejlődése III
6.1 Digitális metszet kezelés A célkitűzéseink között magas prioritással szerepelt, hogy a programunk képes legyen a 3DHistech Kft. által fejlesztett digitális metszet formátumot kezelni. Ezzel elérhetjük, hogy a felhasználó külön program használata nélkül, közvetlenül a metszeten jelölje ki azt a területet, amin aztán szeretné a vizsgálatot lefuttatni. A cégtől kaptunk egy függvénygyűjteményt, az ún. SlideAC-t, mely segítségével a digitális metszetek egyszerűen kezelhetőek. A digitális metszet működési elve a képpiramisra épül: a szkennelés során a digitális mikroszkóp végigpásztázza a metszetet, majd ez elkészült képekből, az átfedések kezelése után elkészül a képpiramis legalsó, azaz legnagyobb felbontású szintje.
Oldal 23
15. ábra Képpiramis [13]
Annak érdekében, hogy az alacsonyabb nagyításon ne kelljen túl nagy képeket kezelni, a legnagyobb felbontású kép fölé még további kilenc készül, és mindegyik szint fele akkora nagyításban készül el, mint alatta lévő. Az így elkészült képpiramis összesen 33%-kal foglal több helyet, mint ha csak a legnagyobb felbontású kép kerülne eltárolásra, de az alacsonyabb felbontású képek jelentősen gyorsítják mind a feldolgozást, mind a megjelenítést.
1. egyenlet A képpiramis által maximálisan felhasznált tárterület
A képpiramis minden egyes szintje kisebb képekből áll össze, és mindegyik ilyen kis képelem 256 x 256 pixel méretű a nagyítástól függetlenül. A SlideAC nem enged képpont szintű hozzáférést a metszethez, hanem a kép pixel koordináta rendszere fölé definiál egy újat, amelynek alapegysége egy képelem. Ennek következménye, hogy a lekérdezett kép szélessége és magassága is 256 többszöröse lehet csak. 6.1.1 A digitális metszetek kezelése a gyakorlatban
A digitális metszetek kezelését az alábbi szempontokat figyelembe véve terveztem meg: a navigálás a metszet legkisebb nagyítású képén történjen a használható nagyítási szintek egyezzenek meg a képpiramis szintjeivel a navigálás képelemenként történjen, és a megjelenített képen lehessen kijelölni azt a területet, amin később a feldolgozást szeretnénk lefuttatni a megjelenített kép mérete 512 x 512 pixel legyen a navigálás során a távlati képen jelenjen meg az aktuálisan mutatott kép befoglaló négyzete. Így a feladat két részre bontható. Egyrészt a navigáláskor kell az aktuális nagyításnak, és az egér mozgásának megfelelő képet megjeleníteni, másrészt a megjelenített képen a kijelölt területet a feldolgozáshoz szükséges nagyítási szinten lekérdezni. A SlideAC-hez Oldal 24
kapott dokumentáció alapján a digitális metszetről az alábbi függvény segítségével lehet képet lekérdezni: HRESULT GetImage([in] long x1, [in] long y1, [in] long x2, [in] long y2, [in] long Magnification, [in] long edChn, [in] long GreenChn, [in] long BlueChn, [in, out] IMrxBitmapImage* pPicture). A navigáció szempontjából a legfontosabb paraméterek az (x1, y1), (x2, y2) koordináták és a „Magnification”, azaz a nagyítási szint. Az (x1, y1) számokkal lehet megadni a lekérdezendő kép bal felső, míg az (x2, y2) számokkal a jobb alsó képelem koordinátáit. A nullás nagyítási szint felel meg az szkennelés során felvett képek nagyítási szintjének, a kilences pedig a piramis legfelső szintjét jelenti.
16. ábra Navigáció
A távlati képen egy pixel megfelel a legnagyobb nagyításon 2 x 2 képelemnek. Így a befoglaló négyzet mérete az alábbi képlet segítségével számítható:
2. egyenlet Slide kezelés: befoglaló négyzet mérete
A felső egyenletben az „m” a nagyítási szintet, az „a” a négyzet oldalainak hosszát jelöli pixelben. A következő lépésként a lekérdezendő kép bal felső képelemének koordinátáit kell kiszámolni, ezt az alábbi képlettel lehet megtenni:
3. egyenlet Bal felső képelem koordinátáinak számítása
Az egyenletben (tX, tY) a képelem koordinátáit, (X, Y) az egér pozícióját, az „m” pedig a nagyítási szintet jelöli. Tehát már megvan a bal felső képelem koordinátája és Oldal 25
természetesen a nagyítási szint, azaz már csak a jobb alsó képelem hiányzik. Ennek értéke (tX + 1, tY + 1), mivel a feltételek között az szerepelt, hogy a megjelenítendő kép 512 x 512 pixel méretű legyen. A következő lépésként a megjelenített képen a felhasználó által kijelölt területet kell a SlideAC segítségével lekérdezni.
17. ábra Képrészlet kijelölése feldolgozásra
Ehhez szükség van a megjelenített kép nagyítási szintjére, a bal felső képelemének koordinátáira, melyet az előző lépésben meghatároztunk, illetve az új nagyítási szintre valamint a kijelölt téglalap bal felső és jobb alsó sarokpontjának koordinátáira. Első lépésként meg kell határozni a megjelenített kép bal felső képelemének koordinátáit az új nagyítási szinten:
4. egyenlet Bal felső képelem koordinátáinak kiszámítása az új nagyítási szinten
A fenti egyenletben a „diff” a régi és az új nagyítási szint különbségének abszolút értéke, a (bX, bY) koordináták pedig a bal felső képelem koordinátái az új nagyítási szinten, a továbbiakban erre báziskoordinátákként hivatkozok. A következő lépésként a megjelenített kép koordináta rendszerében kell kiszámolni az új nagyítási szintnek megfelelő képelem szélességet. Például azonos nagyítási szinten ez az érték 256, hiszen a lekérdezett kép két képelemből épült fel, eggyel kisebb nagyítási szinten 128 és így tovább, általánosan:
5. egyenlet Képelem szélessége az új nagyítási szinten
Az utolsó lépésként meg kell határozni azoknak a képelemeknek a koordinátái, melyek tartalmazzák a felrajzolt téglalap bal felső, és jobb alsó sarokpontjait, a következő módon: Oldal 26
6. egyenlet Befoglaló téglalap sarkainak koordinátái
Az egyenletben (pX, pY) koordináták rendre a téglalapot meghatározó sarokpontok. Az így megkapott értékek azonban csak a bázishoz képest jelölik ki a lekérdezendő képet, ezért még el kell tolni őket a báziskoordinátákkal.
6.2 Cuda használata .Net-ben Ahogy azt már leírtam a CUDA alapvetően a C nyelvet támogatja, a projekt fő nyelve viszont a C#, ezért valamilyen módon meg kell oldani, hogy a menedzselt .Net kódból is meg lehessen hívni a GPU-ra írt függvényeket. Triviális megoldásként adódik, hogy a külön C-ben írt alkalmazásokban letesztelt kerneleket, a hoszt kóddal együtt egy közös gyűjteménybe pakoljuk, majd egy dll-t fordítsunk belőle, és ezt importáljuk a fő programba. A nem menedzselt kódban írt függvények hívása a menedzselt C#-ból ugyan megoldható, de tapasztalataim alapján bonyolult, nehézkes és a hibakeresés így sem a kernelekben, sem a hoszt kódban nem megoldható. Előnyösebb a GASS CUDA.Net keretrendszer használata, mely a meghajtó API-t valósítja meg .Net alatt, így legalább a hoszt kód C#-ba kerül, ami leegyszerűsíti a hibakeresést. Továbbra is szükség van egy külön fájlra, de ez már csak a kerneleket tartalmazza, vagy a köztes PTX nyelvre fordítva, vagy bináris formában. A PTX azért előnyösebb, mert így garantált, hogy bármelyik architektúrán lefut a kód. A CUDA.Net használata ugyan jóval egyszerűbb mint ha a C-ben írt dll függvényeit használnánk, de még így is több lépcsőfokot kell végigjárni. Először is meg kell írni, és le kell tesztelni magát a kernelt, mivel a kernelben hibát keresni csak a futásidejű API-t használva lehet, ezért célszerű minden kernelhez egy külön kis alkalmazást írni. Ha a kernel készen van, akkor át kell másolni egy külön állományba, amit aztán a megfelelően felparaméterezett CUDA fordítóval le kell fordítani. Az eredmény egy csak PTX, vagy GPU bináris kódot tartalmazó fájl, amit aztán már fel lehet használni a meghajtó API segítségével a .Net-ben. Ezt a folyamatot a későbbiekben még részteltesen tárgyalom.
Oldal 27
18. ábra CUDA használata .Net-ben
Oldal 28
7 Sejtmagdetektálás 7.1 A probléma ismertetése Első lépésként meg kell ismerni magát a sejtmagot. A sejtmagok területe egy ép vastagbélszövetről készült mintán átlagosan 25-35 um2, alakjuk rendkívül változatos lehet, de általánosságban elmondhatjuk, hogy lekerekített, alig szegmentált. Számunkra fontosabb tulajdonság, hogy a sejtmagok a minta festése során jobban festődnek, mint a környezetük, sötétebbek lesznek.
19. ábra Sejtmagok
Célunk, hogy a binarizálás során, egy olyan maszkot készítsünk, ami lefedi a sejtmagok teljes területét. Ezt megkönnyíti, ha szokásos RGB szintér helyett a képet HSV modellbe konvertáljuk.
7.2 Konvertálás HSV modellbe [23] A HSV az RGB színtér gyakori reprezentációja, amely sokkal precízebben írja le a színkapcsolatokat, mit az RGB. A HSV és az RGB egyaránt három dimenzióval rendelkeznek, de egészen eltérő módon határoznak meg egy-egy színt.
Oldal 29
20. ábra Az RGB modell ábrázolása [22]
Az RGB rendszert a három alapszín illetve ezek angol megfelelői alkotják (Red, Green, Blue). Ezen három színből bármilyen szín kikeverhető. A három dimenzió - x, y, z tengelyek - mindegyikéhez hozzárendelve a három alapszín egyikét megkapjuk az RGB színteret, melyben az origó megfelel a szín (fény) hiánynak, az x, y, z koordinátája pedig az adott szín piros zöld, kék összetevőinek. Az RGB modell mondható a legelterjedtebb színábrázolási módszernek az informatikában. HSV a színeket egy cilinderen értelmezi. RGB térnek megfelelően itt is minden egyes szín egy pontot jelöl meg a cilinder terében. A szín árnyalata meghatározza a cilinder körlapján a szöget, a telítettség a centrális ponttól való távolságot, a világosság pedig a cilinderen belüli magasságot határozza meg.
21. ábra HSV modell két ábrázolása [22]
Mivel a HSV egy egyszerű transzformációja az RGB színkódolási rendszernek, ezért létezik egy egyértelmű leképzés a két rendszer között.
Oldal 30
Kihasználva ezt az egyértelmű leképzést felépíthetünk egy LUT (Look Up Table) – t, mely nagymértékben gyorsítani fogja az átkódolást. HSV rendszer alkalmazása nagymértékben megkönnyíti azon paraméterek meghatározását, amelyekkel maszkolhatjuk a sejtmagokat, hiszen mint azt a mintán láthatjuk, a sejtmagok főként lilás színnel rendelkeznek, színtelítettségük viszonylag magas, világosságuk alacsonyabb a környezeténél.
7.3 Az algoritmus Sajnos a képek, amiken az egyes sejtmagokat, ill. később a mirigyeket, hámréteget, detektálni szeretnénk, zajosak. A mi esetünkben a zajok forrása elég sokrétű lehet, zaj keletkezik az egyes minták szkennelése folyamán, és a jpg tömörítés során is. Szkennelés során előfordulhat, hogy a lencse vagy a tárgylemez szennyeződéseket tartalmaz, sőt, az egyes mikroszkópok fénytörési tulajdonságai is eltérhetnek. Ez utóbbi nem pontszerű, pixel szintű zajokat eredményez, hanem az egész minta színtelítettségét eltolja. Ez a későbbiek folyamán nehezíteni fogja a binarizálásnál egy általános küszöb érték meghatározását. A zajok miatt az első lépésként egy zajszűrést kell elvégeznünk. Zajszűrés esetén fennállhat annak a veszélye, hogy képen számunkra hasznos információkat is elmosunk. A mi esetünkben azonban ilyen káros mellékhatásról nem beszélhetünk, hiszen célunk nem egy éles határra rendelkező objektum éleinek detektálása, hanem egy maszk készítése. A zajszűréssel elérhetjük, hogy ez a maszk kevésbé lesz szegmentált, és reményeink szerint jobban lefedi a valódi sejtmagokat. A zajszűrés után következhet magának a maszknak az elkészítése. A maszk alatt az alábbi tulajdonságokkal rendelkező képet értjük:
mérete ugyanakkora, mint az eredeti kép mérete,
ahol az eredeti képen sejtmagot valószínűsítünk, fehér színű,
egyéb helyeken fekete.
Ahhoz, hogy a képet a szegmentálni tudjuk, meg kell határozni a fentebb ismertetett HSV modell egyes csatornáinak küszöbértékeit. Ha a küszöbértékek megvannak, már el tudjuk készíteni a maszkot. Zajszűrés ellenére az elkészített maszk nem független a zajtól ezért végső lépéként a maszkot is kezelnünk kell. A maszk utólagos kezelésére a legjobb módszer az Open. Ezzel a lépéssel szétválaszthatjuk az esetlegesen összetapadt sejtmagokat és eltüntethetjük a maszkon található kisebb zajokat is.
7.4 Ultimate Erode Gyakran előfordul, hogy a sejtmagok egymáshoz közel helyezkednek el. Ilyen helyek például a mirigy fala valamint a hámréteg. Ezeken a helyeken a sejtmagok szorosan egymás mellett helyezkednek el. Oldal 31
22. ábra Összetapadt sejtmagok
Ilyen esetben ezek több sejtmag maszkja összeérhet. Ezeket külön kell választani a pontos detektáláshoz.
23. ábra Összetapadt sejtmagok binarizálása
Fentebb már említettem, hogy az összeragadt sejtmagokat Open algoritmussal kezeljük. Azonban a fenti képen jól látható, hogy egyetlen Open algoritmus nem nyújt megoldást a problémára. Ezt a feladatot az Ultimate Erode nevű algoritmust fejlesztettem ki. Az algoritmus a binárizált képet vár bemenetként. Jelen esetben feketével jelöltem a két sejtmag összetapadt maszkját (24. ábra/a). Oldal 32
24. ábra Ultimate Erode
Az algoritmus megvizsgálja és feljegyzik a maszk külső pixeleit. Azt a pixelt nevezzük külső pixelnek, amelynek legalább egy nem a maszkhoz tartozó pixellel szomszédos nyolcas szomszédságot alkalmazása esetén.(24. ábra/b) Ez a lépés megegyezik egy Eróziós lépéssel [24], azzal a kivétellel, hogy jelen esetben el kell tárolni a külső pixeleket. A külső pixeleket eltávolítjuk a maszkon, majd újra végrehajtjuk az előző lépést mindaddig, amíg a maszk teljesen feldolgozásra kerül. Amennyiben több sejtmag maszkja összeér. A csatlakozási szakasz valamennyivel keskenyebb lesz, mint maga a sejtmag átmérője. Ezt használja ki az algoritmus. Mint az a képen is látszik a folytonos méret csökkentés során a maszk egy idő után kettészakad. Így amennyiben tovább folytatjuk az összetapadt sejtmagmaszkok centrumát kapjuk meg. Mivel a minta csupán egy igen vékony metszet az eredeti szövetből a sejtmagok területe igen eltérő lehet. Attól függ, hogy a gömbszerű sejtmagot épp hol metszették el. Ennek tudatában nem elegendő csupán a legutolsó pillanatban megnézni a megmaradt centrumokat, hiszen így a kisebb sejtmagok maszkjai már eltűntek. Az eredményen maximumkeresést kell végezni. A fenti képen két maximum gócot fogunk találni. Ezek a gócok megfelelnek egy-egy sejtmag centrumának. A baloldalon elhelyezkedő maximumgóc egy darab pixelből áll így a centrumpont megegyezik ennek koordinátájával. A jobboldalon elhelyezkedő maximumgóc 3 pixelből áll ebben az esetben a centrumpontot ezek koordinátáinak átlagolásával határozzuk meg. A minta bonyolultsága miatt a binarizálása során előfordul, hogy olyan részeket is lefed a maszk, amik nem sejtmagot jelölnek. Az algoritmus végén ezért végigjárjuk a centrumpontokat és megvizsgáljuk, hogy azok hányadik ciklusban kerültek végleges feldolgozásra. A ciklusszámból következni lehet az eredeti maszk méretére. Amennyiben a ciklusszám nagy nem lehet egyetlen sejtmag, vagy több összetapadt sejtmag vagy valamilyen binarizálási hiba. Minkét esetben ezeket a programnak fel kell ismernie és közölnie kell a felhasználóval. Amennyiben a ciklusszám túl kicsi nem jelölhet sejtmagot így ezeket el kell dobni.
Oldal 33
7.5 Sejtmagok paramétereinek mérése A sejtmagokat lefedő maszk elkészült, a következő feladatunk, hogy ezt felhasználva a sejtmaghoz tartozó paramétereket meghatározzuk. A meghatározandó paraméterek közé tarozik a terület, kerület, ill. a sejt morfológiai jellemzői. 7.5.1 Kitöltő algoritmus
A kitöltő algoritmust használjuk a sejtmag terültének megmérésére, valamint a sejteket befoglaló téglalap meghatározására. Az algoritmus működése egyszerű. Megvizsgáljuk az adott pontot, és ha az egy megadott színtartományba esik akkor átfestjük, valamint minden még nem átfestett szomszédjára végrehajtjuk ezen lépéseket. Forráskódba könnyen átírható rekurzió használatával azonban nem a legszerencsésebb megoldás mivel a rekurziós módon megírt algoritmusok nagy mennyiségű adatok esetén sokkal lassabb futási időt eredményeznek, mint az iteratív megvalósítás és sok esetben veremtúlcsordulást is okozhatnak. Az iteratív megvalósításhoz két tárolóra van szükség. Az első a vizsgálatra várakozó elemeket, a második a már megvizsgált elemeket tartalmazza. A gyors futás érdekében ezen tárolóknak a következő funkciókat kell megvalósítaniuk: o gyors elemhozzáadás o gyors elem eltávolítás o gyors tartalmazás lekérdezés. Célunk egy ponthalmaz beazonosítása, amelyhez első megközelítésben FIFO (First In First Out) puffert alkalmaznánk az első halmaz esetén, de jobban átgondolva a problémát rájöhetünk, hogy a beazonosítás sorrendje nem hordoz számunkra hasznos információt. Egyszerűség kedvéért nevezzük el az elsőt in, a másodikat out tárolónak.
25. ábra A kitöltő algoritmus működési elve
Kezdetben az in tároló csak a kiindulási pontot tartalmazza, az out tároló üres. Az algoritmus első lépésben vesz az in-ből egy elemet. Megvizsgálja, hogy annak a pontnak a színe az adott tartományon belűre esik. Amennyiben teljesült a feltétel, a szomszédos pontok közül azokat, amelyek sem az in, sem az out tárolóban nem szerepelnek betesszük az in halmazba. Utolsó lépésben a megvizsgált elemet minden esetben áthelyezzük az out halmazba, itt végezzük a paraméterek mérését. A területmérést úgy végezzük, hogy megszámoljuk az áthelyezett pixeleket, befoglaló téglalaphoz eltároljuk a legkisebb és Oldal 34
legnagyobb x és y koordinátákat, és meghatározzuk a súlypontot. Az algoritmus addig fut ameddig az in halmaz üressé nem válik.
26. ábra A nyolcas és négyes szomszédság
Szomszédos pontoknak tekinthetjük a nyolcas vagy a négyes szomszédságot. Négyes szomszédság használata ajánlott, hiszen ez is teljesen bejárja a kétdimenziós teret és maximálisan 4 szer vizsgál meg minden egyes pontot ellenben a nyolcas szomszédság használata esetén. 7.5.2 Körvonal detektáló
Szeretnénk a sejtmagok kerületét is megmérni. Erre a legjobb módszer, ha egyszerűen körbejárjuk a maszkon a sejtet határoló pixeleket és megszámoljuk őket. Ehhez egy körvonal detektáló algoritmust használunk. Az algoritmus kiinduló pontja a megmérendő objektum bal felső pontja, azaz a legkisebb y koordinátájú sor legkisebb x koordinátájú oszlopában lévő pixel. A körbejárás iránya megegyezik az óramutató járásának irányával. Első lépésként definiálunk egy 3x3-as maszkot, melynek a középpontja az aktuális pont, továbbá a maszkon egy kitüntetett pontot, a forráspontot. A forráspont jelöli azt a pontot ahonnan érkeztünk, azaz az algoritmus előző iterációjában megtalált határpixelt.
27. ábra A körvonal detektáló maszk **
Az első lépésben forráspontunk még nincs, az algoritmus az első lehetséges következő kontúrpontból indul, azaz az első vizsgált pont az (1, 0) lesz. Azért ez az első lehetséges pont, mert a kezdőpont meghatározásánál kizártuk, hogy a kezdőpont fölött, ill. tőle balra objektum pontja helyezkedjen el. Az algoritmus lényege, hogy a forráspontból kiindulva az aktuális vizsgálandó pixelről eldönti, hogy az objektumhoz tartozik-e vagy nem. Ha igen, akkor feljegyzi a talált pontot, és a következő iteráció középpontja az aktuálisan megtalált pont lesz, forráspontja pedig az aktuális középpont. Ha nem, akkor folytatja a vizsgálatot a Oldal 35
következő pixelen, addig, amíg objektumot nem talál. A „következő pixel” jelen esetben azt a pontot jelenti, ami az óramutató járásával megegyező irányban az aktuálisan vizsgált pixelt követi. Például az (1,0) pont rákövetkezője ebben az értelemben az (1,1). Nézzünk meg egy példát az algoritmus működésére.
28. ábra A mirigydetektáló algoritmus működése **
A 2. ábra egy a körvonal detektáló algoritmus működését mutatja be. A 2/a ábra mutatja a kiinduló helyzetet, a fekete négyzetek alkotják az objektumot, a fehérek a hátteret. Az első dolgunk, hogy meghatározzuk az algoritmus kiindulópontját, ez a 2/b ábrán piros színnel jelzett pont lesz. Ezután meg kell találnunk az óramutató járásának irányában a következő objektumot alkotó pixelt. Mivel az első lépésben még nincs forráspontunk, ezért ahogy már említettem, az első lehetséges objektumot alkotó pontból indulunk ki, ez az aktuális középponthoz képest az (1,0) pont. A kiinduló pont nem része az objektumnak, ezért a következő is meg kell vizsgálnunk. Ha az óramutató járásának irányában továbblépünk, a következő pont már az objektum része lesz, ezért fejlegyezzük a kontúrpontok közé, és a következő iterációban ez a pont lesz a középpont, a mostani középpont pedig a forráspont. Ezt az állapotot a 2/c ábra szemlélteti. Ebben az állapotban már definiált a forráspont, ezért innen indítjuk a belső ciklust. Ebben az esetben a következő objektumot alkotó pont a (1,0) relatív koordinátájú pont lesz. Ezt fejlegyezzük a kontúrpontok közé, és a keresést tovább folytatjuk. Az eljárás akkor áll le, amikor a megtalált pont a kiindulópont lesz. 7.5.3 Átmérők megállapítása [21]
Feladatunk egy olyan algoritmus megalkotása mely meghatározza a sejtmag legkisebb és legnagyobb átmérőjét. Mivel a sejtmagok általában ovális vagy kör alakúak a sejtmag átlói jó közelítéssel egybeesnek az őt befoglaló téglalap átlóival. Első lépésben meghatározzuk a sejtmagot befoglaló téglalapot a sejtmag legnagyobb és legkisebb x, y koordinátája alapján. Ezt követően meghatározzuk a téglalap átlóit. Itt fontos megjegyezni, hogy az átlók számát úgy kell meghatározni, hogy annak száma minimális legyen, de még jó közelítést adjon. Amennyiben az átlók száma nagy, a program futása lassú lesz, ha alacsonyra állítjuk pontatlan eredményt fogunk kapni. Miután meghatároztuk az átlókat, meg kell határozni minden egyes átló esetén a ki és belépési pontot. Ennek leggyorsabb módja, ha a téglalap falától haladunk az átló mentén a középpont felé. Az ábrán jól látszik, hogy így sokkal kisebb távolságot kell megvizsgálnunk, mint ha a középponttól indulnánk. Oldal 36
29. ábra Sejtmagmaszk befoglaló téglalapban [21]
Az így kapott átlónkénti két pont távolsága megadja az átlók hosszát, ezeken már csak egy egyszerű minimum és maximumkeresést kell végrehajtanunk. Optimalizációt jelent, ha a befoglaló téglalap átlóit nem szinusz és koszinusz függvényekkel határozzuk meg, hanem a téglalap falán adott egységenként lépkedünk.
Oldal 37
8 Mirigydetektálás 8.1 A probléma megismerése A mirigyek háromféle alakot is felvehetnek, attól függően, hogy hogyan metszették el őket. Ha a hossztengelyükre merőlegesen vagy nagyobb szöget bezárva metszik el őket, akkor kör vagy ellipszis, ha a hossztengelyükkel párhuzamosan, akkor nyitott félcső alakúak lesznek. Az sem garantált, hogy az elliptikus alakúak közepe „lyukas”. Könnyen előfordulhatnak olyan metszetek, ahol a mirigy közepén sejtmagok találhatók. Mi az egyszerűség kedvéért, csak a kör, ill. ellipszis alakú mirigyekkel foglalkozunk, a félcső alakúakat hámrétegnek vesszük, és külön foglalkozunk velük. A hámréteg detektálásról később lesz szó.
30. ábra Mirigyek, a hossztengelyükre merőlegesen elvágva
Egy tulajdonság azonban minden mirigyre jellemző, sűrű sejtmagsor határolja őket. Ezt használjuk ki a mirigyek detektálásánál, felhasználjuk a már korábban elkészített sejtmagtérképünket. A mirigyek pozícióját, elhelyezkedését az egyes sejtmagok egymáshoz viszonyított helyzetéből próbáljuk meghatározni.
8.2 Az algoritmus A teljes mintát ábrázoló képből indulunk ki. Ennek a képnek a maszkja úgy fog előállni, hogy a már korábban meghatározott sejtmagok helyére egy pár pixeles foltot rajzolunk. Az elhelyezett pixelekből felépítjük a későbbiekben (Háló készítés című fejezetben) ismertetett algoritmust. Így egy olyan maszkot kapunk ahol az üresen hagyott foltok helyén mirigyek lehetnek. Azt hogy ott ténylegesen mirigyek vannak-e, a méretük alapján tudjuk eldönteni.
Oldal 38
Ehhez végigmegyünk a maszk összes pontján, és ahol mirigyet sejtünk, azaz fekete folt van, meghívjuk a már korábban ismertetett kitöltő algoritmust, ami átfesti a megfelelő pixeleket, így azokat később már nem fogjuk vizsgálni. A kitöltő algoritmus ráadásul visszaadja a kitöltött terület nagyságát is, amit fel tudunk használni arra, hogy eldöntsük mirigyről van-e szó. Ezt azért tudjuk megtenni, mert tapasztalat alapján meghatározhatunk egy minimum területet, ami alatt már nem fogadjuk el a kitöltött területet mirigynek. Az algoritmus bizonyos esetekben olyan területeket is kiválaszt, amik nem mirigyek. Ezek tipikusan a mirigyek közötti olyan területek ahol a sejtmagsűrűség alacsony. Ezeket a területek ki tudjuk szűrni, ha minden esetben, megvizsgáljuk az alakzat kompaktáságát
7. egyenlet Kompaktság [25]
A kompaktsági érték az olyan területek esetében, amik valóban mirigyeket fednek el kicsi lesz, mivel kör vagy valamilyen ovális alakzathoz hasonlítanak. Azonban mirigyek közötti alacsony sejtmagsűrűségű területet lefedő maszkok esetében ez az érték magas lesz. Végezetül egy dilatációt kell végrehajtanunk, mert a megkapott maszk csak a határoló sejtmagok súlypontjáig ér. A mirigyeket lefedő maszk elkészült, az utolsó lépés az, hogy a ugyanazokat paramétereket ugyanúgy megmérjük, mint a sejtmagoknál.
a
8.2.1 Háló készítése [21]
Háló készítéséhez felhasználjuk a már korábban meghatározott sejtmagtérképet. Első lépésben készítünk egy új képet, melyen kezdetben csak a sejtmagok szerepelnek egységnyi kiterjedéssel. A sejtmagok közepéből csigavonalas bejárást indítunk adott távolságig. Amennyiben a bejárás során találunk egy fekete pixelt, másik sejtmagot, egy vonalat húzunk a két pont közé, és ha elértük a maximális kapcsolat számot leállítjuk az algoritmust. Ezeket a lépéseket minden egyes sejtmagnál végrehajtjuk.
Oldal 39
31. ábra Háló készítése [21]
A fenti ábrán piros színnel jelöltük az aktuálisan vizsgált sejtmagot, szürkével azokat a pixeleket, amelyeket a csigavonalas bejárás már vizsgált, feketével pedig a környezetében található sejtmagokat. A csigavonalas bejárásra azért van szükség, hogy a legközelebbi sejtmagokat találjuk meg először. A lépések leírásában jól látható hogy két értéktől függ a végeredmény. Függ a csigavonalas bejárás maximális távolságától és a maximálisan összeköthető sejtmagok számától. Ennek megfelelő meghatározása fontos, de mivel tapasztati úton történik, a munkafolyamat későbbi fázisában kerül sor. Fő feladatunk, hogy az elkészült háló segítségével meghatározzuk a mirigy maszkját. Korábbiakban meghatároztuk a mirigy közepét igaz nem túl pontosan, de a célnak jelenleg megfelel. Sejtmag közepéből indítjuk a korábbiakban már ismertetett kitöltő algoritmust. Optimális esetben így jó közelítést kapunk a sejtmag maszkjára. Azonban ha jobban megfigyeljük az ábrát, láthatjuk, hogy vannak olyan mirigyek melyeknek a közepén is található sejtmag. Ilyen esetekben a kitöltő algoritmus a valódi mirigymaszk csupán egy részét adja eredményül. Az ilyen hibákat kiküszöbölhetjük, ha a kitöltési algoritmust 45 fokonként, növekvő távolságban többször is elindítjuk ugyanarra a mirigyre. Amennyiben ez már egy bejárt terület nem vesztünk túl sok processzor időt, ha azonban egy még nem azonosított terület, akkor további részét találtuk meg a mirigynek.
Oldal 40
9 Hámdetektálás 9.1 A probléma megismerése A felszíni hám jellemzője, hogy a minta szélén található, sejtmagjai egy sorban helyezkednek el, köztük távolság kicsi, gyakran szinte egybeolvadnak. A sejtmagokat egyenletes plazma réteg veszi körül, mert a sejtek hasáb alakúak.
32. ábra Hám
9.2 A detektáló algoritmus Mivel a hám nagyobb kiterjedésű, ezért az egész mintát kell vizsgálni. A hám detektálásánál felhasználjuk a csak rá jellemző tulajdonságait. Hámot alkotó sejtmag csak a minta szélén helyezkedhet el, ezért első lépésként a minta körvonalát határozzuk meg. Ehhez a mirigydetektálásnál ismertetett módon elkészítjük a sejtmag-hálót, majd az egész háló által körbehatárolt területet kitöltjük. A következő lépésként felhasználjuk a már korábban ismertetett körvonal detektáló algoritmust, egy kis kiegészítéssel. Jelen esetben a körvonallal együtt el akarjuk tárolni azt is, hogy kontúrpont hám része-e, vagy nem. Ehhez a kontúrpont struktúrát kiegészítjük egy logikai változóval. Kezdetben minden kontúrpontnál igaz értékre állítjuk ezt a változót, csak a későbbi vizsgálatok során fogjuk tudni eldönteni, hogy valójában hám része-e a kontúrpont. A pontosításhoz két eljárás szükséges, mindkettő a hám egy-egy alaptulajdonságára épít. Mindkét eljárás visszatérési értéke logikai típusú, igaz, ha hámról szó, hamis. ha nem. A két eljárás visszatérési értékének „és” kapcsolata fogja megadni, hogy ténylegesen hámról van-e szó. Az első eljárás a hámnak arra a tulajdonságára épít, hogy a sejtmagok sűrűn helyezkednek el. Megpróbáljuk kizárni azokat a területeket, ahol a sejtmagsűrűség alacsony. Ilyen területek tipikusan olyan helyeken fordulhatnak elő, amelyek a metszés során megsérültek. Az eljárás első lépése, hogy a kitöltött sejtmaghálóra felfestjük a magokat, a korábban meghatározott koordináták alapján. Minden koordináta helyére egy több pixeles foltot rajzolunk. A következő lépés az, hogy a kontúrpontokon végighaladva megpróbáljuk kizárni a vizsgálatból az alacsony magsűrűséggel rendelkező területeket. Oldal 41
Ez azonban még mindig nem elég a pontos detektáláshoz, hiszen a minta szakított részén is előfordulhatnak nagyobb magsűrűséggel rendelkező területek. Ezeket a második eljárás során próbáljuk meg kizárni.
33. ábra A minta szakított oldala
A második eljárás során azt használjuk ki, hogy a hámot alkotó sejtek hasáb alakja miatt, a sejtmagokat egyenletes plazma réteg veszi körül. Első lépésként binarizáljuk az eredeti mintát ábrázoló képet, úgy hogy az elkészült maszk a minta teljes egészét lefedje. Ezt fogjuk összehasonlítani kitöltött hálós képpel. Ahol a kettő határa egybeesik, nem lehet hám, mert hiányzik a plazmaréteg. A két eljárás metszetével már valószínűleg nagy biztonsággal meg tudjuk határozni a hámréteg elhelyezkedését.
9.3 Sejtmagok differenciálása Három szöveti szerkezetben található sejtmagokat kell elkülönítenünk. Sejtmagok találhatók a mirigyekben, a hámszövetben és a mirigyek közti állományban, melyet más néven lamina propria-nak nevezünk. A sejtmagkereső algoritmusunk ellenben nem tud különbséget tenni a különböző struktúrákhoz tartozó sejtmagtípusok között, ezért volt fontos, hogy beazonosítsuk a szövetben található különböző struktúrákat. Mivel beazonosítottuk a mirigyekhez és a hámhoz tartotó sejtmagokat, egyértelművé válik, hogy az eddig be nem sorolt sejtmagok a lamina propria, azaz kötőszöveti réteghez tartoznak.
Oldal 42
10 Képfeldolgozás gyorsítása GPGPU segítségével Ebben a fejezetben bemutatom a kétféle API használatát, és megvizsgálom, hogy a fejlesztőkörnyezethez mellékelt dokumentációkban [7][12] vázolt optimalizálási eljárásokat hogyan lehet alkalmazni képfeldolgozó algoritmusok gyorsítására.
10.1 Alapok, használat Ahogy már említettem, a CUDA két féle API-t biztosít a GPU erőforrásainak kihasználásához. Mindkettőnek megvan a saját felhasználási területe: a futásidejű API-t egyszerűbb használni, de csak kisebb programokban alkalmazható, míg a meghajtó API-t használva akár más nyelvekből is meg lehet hívni a GPU-ra írt függvényeket, de használata bonyolultabb. A használt API-tól függetlenül, egy GPU-t is kihasználó program szerkezete ugyanaz.
34. ábra Egy CUDA-t használó program felépítése
Minden CUDA-t használó program három fő részre osztható a GPU szempontjából: 1. kernel indításának előkészítése, kernel indítása, 2. kernel futtatása, 3. eredmények kezelése. Az első részben elsősorban a kernel bemenő paramétereinek beállítása történik: itt foglaljuk le a szükséges helyet a GPU memóriájában, majd másoljuk át a CPU memóriájából az adatokat. A kernel indításakor pedig a szálhierarchiát adjuk meg. Mivel a grafikus egység párhuzamosan dolgozik a központi processzorral, amint fel szeretnénk használni a kimenő adatokat, meg kell hívnunk egy szinkronizáló utasítást, ami addig blokkolja a program futását, amíg az összes szál le nem futott. Csak ezután lehetünk biztosak a kimenetek helyességében. A szinkronizálás után a kimenő adatokat visszamásolhatjuk a CPU memóriájába, majd a már nem használt memória területeket fel kell szabadítani. A következőkben egy egyszerű programon (egy RGB kép binarizálásán) mutatom be, hogy a két API segítségével hogyan lehet a GPU-ra programot írni. 10.1.1 Futásidejű API
A feladat tehát egy RGB kódolású kép szegmentálása valamilyen küszöbértékkel. A küszöbölés során az egyszerűség kedvéért, minden pixel intenzitásainak átlagát hasonlítom Oldal 43
egy határértékhez. Természetesen a legelső lépés magának a képnek a megnyitása, azaz a háttértárról operatív tárba történő betöltése. Ehhez az OpenCV könytárat használom, ami a betöltés után egy egyszerűen kezelhető struktúrába menti el a kép legfontosabb adait. Ezekből most az alábbiakra van szükségem:
char *imageData: egy mutató, mely a képi adat kezdetére mutat a CPU memóriájában, int width: a kép szélessége pixelben, int height: a kép magassága pixelben, int widthStep: a kép egy sorának lefoglalt memória, mértékegysége egy bájt.
A „widthStep”-re azért van szükség, mert az egy sor számára lefoglalt memória nem feltétlenül egyenlő a kép szélessége szorzva az egy képpont által lefoglalt területtel: a hatékonyabb memória elérés miatt, az OpenCV néggyel osztható számmá egészíti ki. Először a kép számára le kell foglalni a szükséges memória helyet a globális memóriában: cudaMalloc((void**) &d_imageData, widthStep * height). A függvény lefoglal a második paraméterben megadott mennyiségű memóriát, és visszaad egy, a globális memóriában értelmezett mutatót, mely ennek a memória területnek az elejére mutat. Ezt a pointert később a GPU-n futó szálak fogják használni, a központi processzor címterében nem értelmezett. Most fel kell töltenünk adattal a lefoglalt memória területet: cudaMemcpy(d_imageData, h_imageData,widthStep*height, cudaMemcpyHostToDevice). Ez az a függvény, amely a csatolósínt felhasználva a CPU operatív tárjából, a grafikus egység memóriájába másolja az adatokat. Első paramétere az előző lépésben megkapott mutató, a cél pointer, a második paramétert a kép megnyitásakor kaptuk meg az OpenCVtől, ez a forrás pointer. A harmadik paraméter a másolandó mennyiség, a negyedik pedig egy felsorolás típus, mely a másolás irányát jelöli ki. A következő lépés magának a kernelnek az indítása. Ehhez először meg kell határozni a szálhierarchiát, azaz a grid, és a blokkok dimenzióit. A blokk méretét a hatékonyabb feldolgozás miatt érdemes úgy megválasztani, hogy az egy blokkban lévő szálak száma az egy multiprocesszor által párhuzamosan futtatni képes szálak számának, azaz a warp méretének többszöröse legyen. Ez a mai architektúrákon 32. A grid dimenzióit a feldolgozni kívánt adatstruktúra, és a már megválasztott blokk méretéhez kell igazítani, úgy hogy biztosítva legyen az összes adategység feldolgozása. A példára visszatérve: ha a kép szélessége maradék nélkül osztható a választott blokk szélességgel, akkor a grid szélessége az osztás eredménye, ha nem, akkor az osztás eredménye eggyel növelve:
8. egyenlet Grid dimenziók számítása
A képletben a „gWidth” a grid, a „width” a kép, a „bWidth” pedig a blokk szélességét jelenti. Természetesen ugyanígy számolhatjuk ki a szükséges grid magasságot is rendre a magasságokat behelyettesítve. Ezzel a módszerrel biztosíthatjuk, hogy akkor is megfelelő Oldal 44
számú szálat indítsunk, ha a kép szélessége nem osztható a blokk szélességével. Ez viszont magával hozza, hogy lesznek szálak, amelyek a kép határán kívül lesznek, azaz „semmi dolguk” nem lesz. Ezt a kernelben kezelni kell. Ez pazarlásnak tűnhet, de ebben az esetben csak a jobb szélső, és az alsó blokkokban lesznek olyan szálak, amik semmit sem csinálnak. Míg ha rossz blokkméretet választunk, akkor minden blokkban lesz olyan warp, amelyik nincs kitöltve.
35. ábra Egy képre feszített szálhierarchia
Miután kiszámoltuk a blokk és a grid méreteit, indíthatjuk is a kernelt a már ismertetett módon: binarize<<>>(d_imageData, thres, widthStep, height); A paraméterek többsége már ismert, a „thres” a szegmentálás küszöbértéke. A meghatározott blokk és grid méretekből következik, hogy minden szálra egy képpont feldolgozása jut. Ez első megközelítésnek megfelel, a későbbiekben viszont látni fogjuk, hogy ez nem jelent optimális memóriaelérést. Mivel a kernel párhuzamosan fut a hoszt kóddal, a kernelhívás után a vezérlés egyből visszaadódik CPU-hoz. Ezért mielőtt fel szeretnénk használni a függvény kimenetét, meg kell hívnunk a szinkronizáló utasítást: cudaThreadSynchronize(); Ez a sor addig blokkolja a program futását, amíg a grafikus egység dolgozik. Ezután már biztosak lehetünk a kimenetek helyességében, visszamásolhatjuk a CPU memóriájába.
Oldal 45
Ehhez a már ismert másoló függvényt használjuk, csak a cél és a forrás helyet cserél, és az irány megfordul: cudaMemcpy(h_imageData,d_imageData,widthStep*height, cudaMemcpyDeviceToHost); Ez a példaprogram az egyszerűség kedvéért a már megnyitott képet írja felül a szegmentált képpel, de természetesen új helyre is visszamásolhatjuk az adatokat. Ezután már nem maradt más hátra, mint a GPU memóriájában felszabadítani a lefoglalt területet. 10.1.2 Meghajtó API
Ebben a fejezetben nem a C nyelvben használható meghajtó API-t mutatom be, hanem annak a C#-ban is használható megfelelőjét a CUDA.Net-et. A kettő logikai felépítése ugyanaz, a kulcsszavak és a használat módja kicsit eltér egymástól, hiszen a CUDA.Net kihasználja a C# objektum orientált felépítését. A példaprogram ugyanaz, mint amit fentebb ismertettem. Ahogy már korábban utaltam rá, a meghajtó API használatakor a kernel kódja fizikailag elkülönül a hoszt kódtól, és a fordítása is külön folyamat. Ahhoz hogy a megírt kerneleket használni tudjuk, át kell másolni őket egy külön fájlba, majd a fordítót úgy kell paraméterezni, hogy az eredmény egy bináris vagy PTX utasításokat tartalmazó állomány legyen. A könnyebb újrahasznosíthatóság érdekében célszerű a megírt fordító parancsot egy batch fájlba tenni. A fordítót az alábbi módon kell paraméterezni: nvcc -ccbin "C:\Program Files\Microsoft Visual Studio 9.0\VC\bin" binarize.cu –ptx Az első paraméterként egy C fordítót tartalmazó mappát kell átadni, bár ez ebben az esetben nincs használatban, hiszen az állomány csak CUDA kódot tartalmaz. A következő paraméter a lefordítandó fájl neve, az utolsó pedig azt adja meg, hogy mi legyen a fordítás eredménye. Ide a PTX helyett a CUBIN kulcsszót írva nem köztes kód lesz az eredmény, hanem bináris. A fordító folyamat kimeneteként kapott, jelen esetben, binarize.ptx vagy binarize.cubin fájlok már felhasználhatóak a meghajtó API-val. Az első dolgunk egy CUDA objektum létrehozása, melyen keresztül tudjuk az API függvényeit elérni. Az objektum tartalmaz egy modult, amibe be lehet tölteni egy, a lefordított kerneleket tartalmazó fájlt: cuda.LoadModule(".\\binarize.ptx");
A következő lépésként létre kell hozni egy függvény objektumot: CUfunction binarise = cuda.GetModuleFunction("binarize");
A paraméterként átadott sztring a kernel neve. Természetesen, ha egy modul több kernelt tartalmaz, akkor több ilyen függvény objektum is létrehozható. Mivel mutatók a menedzselt C# kódban nem használhatóak, a GPU globális memóriájába mutató pointernek egy külön struktúrát hoztak létre: CUdeviceptr d_imageData = new CUdeviceptr();
Ahogy a futásidejű API-nál, itt is a memóriakezelés következik. Ehhez a CUDADriver statikus osztály metódusait lehet használni: CUDADriver.cuMemAlloc(ref d_imageData, widthStep * height); CUDADriver.cuMemcpyHtoD(d_imageData, h_imageData, widthStep * height);
Ahhoz hogy elindíthassuk a betöltött kernelt, már csak a paramétereit és a szálhierarchiát kell beállítanunk. A paraméterek átadása bonyolultabb, mint az előző esetben: Oldal 46
int ofset = 0; cuda.SetParameter(binarise, ofset, (uint)d_imageData.Pointer); ofset += Marshal.SizeOf(d_imageData.Pointer); cuda.SetParameter(binarise, ofset, thres); ofset += Marshal.SizeOf(thres); cuda.SetParameter(binarise, ofset, widthStep); ofset += Marshal.SizeOf(widthStep); cuda.SetParameter(binarise, ofset, height); ofset += Marshal.SizeOf(height); cuda.SetParameterSize(binarise, (uint)(ofset));
A „cuda” objektum SetParameter függvényének első paramétere a betöltött függvény objektum, amelynek a paramétereit be szeretnénk állítani. A második paraméter egy ofszet, vagy eltolás, ami az aktuális beállítandó paraméter előtti paraméterek hossza bájtban. Az első paraméternél ez nulla, a következőnél az első paraméter hossza, a harmadiknál az első kettő paraméter hosszának összege és így tovább. Végül be kell állítani az összes átadott paraméter hosszainak összegét. Már csak a szálhierarchia megadása van hátra. Először a blokkok méretét kell beállítani: cuda.SetFunctionBlockShape(binarise,blockSize.Width,blockSize.Height, 1);
A grid méretét pedig az indításkor is megadhatjuk: cuda.Launch(binarise, gridWidth, gridHeight));
A CUDA.Net-ben nincs külön szinkronizáló utasítás, a fentebb is használt függvény futása akkor ér véget, amikor az összes szál lefutott. Ha aszinkron módon szeretnénk kernelt indítani, akkor a Launch helyett a LaunchAsync függvényt használhatjuk. Ha a GPU végzett, az adatot visszamásolhatjuk: CUDADriver.cuMemcpyDtoH(h_imageData, d_imageData, memSize);
Nem szabad elfelejteni, hogy fel kell szabadítani a globális memóriában lefoglalt területet. 10.1.3 A szegmentáló kernel
Ugyan a kernel kódja a két API-t használva ugyanaz, a szignatúrája némileg eltér egymástól. A futásidejű API-t használva az alábbi szignatúrát kell használni: __global__ static void binarize(unsigned char* imageData, int thres, int widthStep, int height). Míg a meghajtó API-t használva, a fordítási folyamat eltérése miatt ez így modosul: extern "C" __global__ void binarize(unsigned char* imageData, int thres, int widthStep, int height). A „__global__” kulcsszóval lehet a fordítónak jelezni, hogy ez a függvény egy kernel: a hoszt hívja, de az eszközön fut. A kernelben az első dolgunk annak a meghatározása, hogy az aktuális szál melyik adategység, jelen esetben pixel, feldolgozásáért felel, azaz kellenek az aktuális képpont X és Y koordinátái: unsigned int X = (blockIdx.x * blockDim.x + threadIdx.x) * 3; unsigned int Y = blockIdx.y * blockDim.y + threadIdx.y;
Ehhez fel kell használni néhány beépített változót:
a blockDim az indításkor átadott blokkdimenziókat tartalmazza, Oldal 47
a blockIdx az aktuális szál blokkjának az indexét, a threadIdx pedig az aktuális szál indexét a blokkon belül. Ezek a változók alapértelmezetten, minden kernelben elérhetőek, helyezkednek el, értéket pedig a szál indításakor a hardvertől kapnak. azért kell hárommal megszorozni, mert egy szál három bájtnyi adatot képre mutató pointer bájt típusú. Az így kiszámolt értékek a szál regisztereibe kerülnek.
a regiszterekben Az X koordinátát dolgoz fel, míg a privát használatú
Miután megvannak a koordináták, már el tudjuk dönteni, hogy ez a szál dolgoz-e fel adatot vagy nem, ha igen, akkor le kell kérdeznünk a pixel RGB értékeit: unsigned char *actElem = (imageData + Y * widthStep) + X; int b = actElem[0]; int g = actElem[1]; int r = actElem[2];
Az actElem az aktuális szál által feldolgozandó pixel első értékére mutat a globális memóriában, indexelésénél figyelni kell arra, hogy az RGB intenzitások BGR sorrendben vannak eltárolva. Minden egyes ilyen utasítás egy load parancsot jelent a memóriavezérlő felé, és mivel a globális memóriának nincs gyorsító tára a chipen, az összes a DRAM-ot terheli. Olyan ez, mintha a CPU-n futó kód összes load utasítása cache hibát adna. Paraméterként csak egy küszöbértéket kapott a kernel, ezért a színes képet, mint szürkeárnyalatos szegmentáljuk, ehhez viszont ki kell számítani az RGB intenzitások átlagát. A következő lépésként a feltétel vizsgálat után visszaírjuk a megfelelő értékeket: if (rgb < thres) { actElem[0] = 0;
actElem[1] = 0;
actElem[2] = 0;
} else { actElem[0] = 255; actElem[1] = 255; actElem[2] = 255; }
Itt ugyanaz mondható el, mint az olvasásnál, csak fordított iránnyal. Minden egyes írás közvetlenül a globális memóriába megy. 10.1.4 A teljesítménymérő alkalmazás
A CUDA fejlesztőkörnyezet tartalmaz egy teljesítménymérő alkalmazást is, amellyel a GPU-n futó kód paramétereit lehet megmérni, elemezni. A program által mért adatok azonban nem feltétlenül felelnek meg a valóságnak. Néhány paramétert warponként, néhányat multiprocesszoronként mér, és ugyanolyan feltételek mellett is előfordulhat, hogy futásonként más eredményt kapunk [11]. Ez azonban nem jelenti a teljesítménymérő alkalmazás használhatatlanságát. Az egyes értékek abszolút vonatkozásban ugyan semmit sem jelentek, a különböző kernelek összehasonlítására viszont alkalmasak. Emiatt nem sorolom fel az összes mutatót, amit mérni tud, inkább a fentebb ismertetett példakernel elemzésén keresztül mutatom be a fontosabbakat, a későbbiek során pedig csak a különbségekre hívom fel a figyelmet. Az összes tesztet, ha csak másként nem jelölöm, egy 512 × 512 méretű képen végzem, egy GeForce GTX295 típusú videokártya egyik Oldal 48
processzorát felhasználva. A tesztrendszer adatai részletesen a 13.1.2. fejezetben találhatóak meg. Futási idő CPU GPU
282,88
Kernel
102,112
GPU CPU
412,512
GridX
GridY
BlockX
BlockY
32
32
16
16
3. táblázat Szegmentáló kernel mérési erdményei I.
A táblázat első oszlopban a GPU-n mért futási idő szerepel, mértékegysége egy mikro szekundum. A következő oszlopok a kernelnek átadott szálhierarchiát mutatják. Látható, hogy egy blokkban 16 × 16 = 256 db szál található, amit így a SIMT egység 256 / 32 = 8 db warpra fog osztani. A blokk méretét, ha csak erre külön nem utalok, a későbbi tesztelések során nem fogom változtatni. A blokk kiterjedése és a kép mérete meghatározza, hogy a grid méretének mekkorának kell lennie: 512 / 16 = 32 db blokkra van szükségünk mindkét irányban, hogy teljesen lefedjük a képet.
Kernel
Regiszterek Kihasználtság
MP / Feltételes blokk utasítás
Elágazó f.u.
Utasítások száma
5
35
123
16347
1
1891
4. táblázat Szegmentáló kernel mérési eredményei II.
Az első oszlop az egy kernel által felhasznált regiszterek számát mutatja. A kihasználtság, angolul occupancy, az aktuálisan aktív warpok aránya a multiprocesszorokban az aktív warpok maximális számához viszonyítva. Az aktív warpok maximális száma architektúránként különböző lehet, jelen esetben 32. A kihasználtságot három dolog befolyásolja: a blokkméret és a regiszterek száma és a felhasznált osztott memória mérete. Ahhoz hogy egy warp aktív lehessen, azaz készen álljon a futtatásra, a benne lévő szálak által használt erőforrásoknak készen kell állnia: pl. szálhierarchia adatok (threadIdx) legyenek betöltve a regiszterekbe, ill. legyen elég hely az osztott memóriában stb. Ha tehát egy szál túl sok regisztert használ, vagy egy blokk túl nagy osztott memóriát foglal le, akkor az aktív warpok száma kisebb lesz a maximálisnál. Fele akkora kihasználtság azonban nem feltétlenül jelent fele akkora teljesítményt: ha a kernel a futása közben egy nagy késleltetésű utasításhoz érkezik, pl. globális load utasítás, akkor a SIMT egység egy következő aktív warp futtatására fog váltani. Teljesítménycsökkenés csak akkor következik be, ha az összes aktív warp futása blokkolt, pl. adatra vár, azaz ha a warpok váltogatásával nem sikerült áthidalni a nagy késleltetésű utasítást. Látható, hogy egy ilyen egyszerű kernel egy modern GPU-n nem befolyásolja a kihasználtságot. Az MP / blokk oszlop, angolul „sm cta launched”, az egész kernel futása alatt az egy multiprocesszor által feldolgozott blokkok számát jelenti. Mivel az általam használt GPUban 30 multiprocesszor van, blokkból pedig jelen esetben 32 × 32 = 1024 van, így egy MPnek legalább 34 db blokkot kell feldolgoznia, de néhánynak 35-öt. A következő oszlop a feltételes elágazó utasítások számát mutatja. Önmagában ennek nincs jelentősége a futásidő szempontjából, csak akkor, ha az egy warpon belüli szálaknál a feltétel vizsgálat
Oldal 49
különböző eredmény ad. Erre a jelenségre már utaltam korábban. Ilyenkor a különböző „utat” bejáró szálak sorosan futnak le. Az utolsó oszlop a kernel által összesen végrehajtott utasítások számát mutatja, ez is inkább érdekesség, gyakorlati haszna nincs.
Kernel
Globális Load 32B
Globális Store 32B
Globális Load kérés
Globális Store kérés
9888
12006
840
1104
5. táblázat Szegmentáló kernel mérési eredményei III.
Ezek a számok a kernel globális memória kezelését jellemzik. A számok darabot jelentenek. Az itt megjelent eredményeket most nem elemzem, mert a fentebb említett okok miatt csak két kernel mért értékeinek összehasonlításának van értelme. 10.1.5 A zajszűrő kernel
A továbbiakban az optimalizálásokat egy binomiális szűrővel közelített Gauss konvolúción mutatom be. A képen végzett konvolúció, egy maszk, vagy ablak alapú művelet, amelynél a kimenet adott pixelének intenzitás értékét a bemeneti kép ugyanannak a pixelének valamekkora környezete határozza meg.
36. ábra Konvolúció [26]
A művelet lépései: 1. a maszkot ráhelyezzük a bementi kép minden pixelére, úgy hogy a maszk középpontja a feldolgozandó pixelt fedje; 2. a maszk értékeit szorozzuk az adott érték által lefedett pixel intenzitásával; Oldal 50
3. a kimeneti képen a maszk középpontjának megfelelő pixel intenzitása a szorzatok összege lesz. A maszkot, vagy ablakot szokás még kernelnek is nevezni, azonban ezt a kifejezést én csak a GPU-n futó függvény megnevezésére fogom használni. Az ablak méretét, és az együtthatók értékeit aszerint választhatjuk meg, hogy mit szeretnénk a képpel elérni. Jelen esetben a cél a sejtmagdetektálás első lépésének gyorsítása, azaz egy zajszűrő, simító művelet elvégzése a képen az alábbi maszkkal: 1
2
1
2
4
2
1
2
1
37. ábra Konvolúciós maszk
Mivel alap esetben a maszk súlyainak összege egy kell hogy legyen, a lineáris kombináció eredményét még el kell osztani 16-tal. Látható, hogy mind a szorzások, mind az osztás kettő hatványával történik, így a szorzó, osztó műveletek helyett a sokkal hatékonyabb shiftelést lehet alkalmazni. Egy ilyen konvolúció megvalósítása a központi processzoron dióhéjban úgy néz ki, hogy egy dupla ciklussal végiglépkedünk a kép pixelein, majd a kisebb maszk esetén közvetlen címzést alkalmazva, nagyobb maszknál egy újabb dupla ciklussal elvégezzük a lineáris kombináció kiszámítását. Ez a grafikus egységet alkalmazva úgy módosul, hogy a külső dupla ciklust elhagyjuk, és minden pixelre egy szálat indítunk. A szálakon belül ugyanazokat a műveleteket kell elvégeznünk, mint a CPU-nál a dupla cikluson belül. A teljes kernel kód bemásolásától terjedelmi okok miatt most, is a későbbiekben is eltekintek, csak a jellemző különbségekre hívom fel a figyelmet. A konvolúciós művelet jellegéből fakadóan, az eredményeket nem írhatjuk vissza forrás képre, hiszen az befolyásolná a többi eredményt. Emiatt a hoszt kód annyival bővül, hogy le kell foglalni az eredmény kép számára is a szükséges helyet. A kernelnek az indításkor átadott paraméterek értelemszerűen módosulnak: kell még egy pointer az új intenzitások elmentéséhez, és a küszöbértékre itt már nincs szükség. Az aktuális szál által feldolgozandó képpont koordinátáinak kiszámítása megegyezik az előző példában látottakkal. A számítás során viszont le kell kérdezni a szomszédos pixelek intenzitását is, ami jelentősen megnöveli az adatforgalmat. Mivel a maszk elég kicsi, ehhez a ciklusok elkerülése végett közvetlen indexelést használok: if (Y > 0) { actElem sumB += sumG += sumR += }
= (Y - 1) * widthStep + (X); srcImageData[actElem] << 1; srcImageData[actElem + 1] << 1; srcImageData[actElem + 2] << 1;
Ez a néhány sor az aktuális pixel fölötti képpont intenzitását veszi figyelembe kétszeres szorzóval. Az „srcImageData” a forrás képi adat kezdetére mutató pointer. A címzés módja némileg módosult: ahelyett hogy minden pixelhez külön mutatót rendelnék, a kép kezdetére mutató pointert indexelem a megfelelő értékekkel. A közvetlen címzés miatt ilyen blokkból összesen kilenc van, különböző belépési feltételekkel és szorzókkal. A Oldal 51
végén aztán a sumB, simG és sumR értékeit osztom 16-tal, és a visszaírom a cél kép megfelelő helyére. Most vizsgáljuk meg ezt a kernelt is a teljesítménymérő alkalmazással. Az első táblázatban csak a kernel futási ideje módosult, ez jelen esetben 425,344 µs. A grid dimenziói sem változtak, ugyanannyi szálat kell indítani, mint az előbbi esetben.
Kernel
Regiszterek Kihasználtság
MP / blokk
Feltételes utasítás
Elágazó f.u.
Utasítások száma
12
35
8132
112
68610
1
6. táblázat Naiv zajszűrő kernel mérési eredményei I.
Ez a táblázat még nem magyarázza a futási idő növekedését: a kihasználtság a több felhasznált regiszter ellenére nem csökkent, az egy warpon belül ténylegesen elágazó utasításokból sincs több.
Kernel
Globális Load 32B
Globális Store 32B
Globális Load kérés
Globális Store kérés
102735
9888
7560
840
7. táblázat Naiv zajszűrő kernel mérési eredményei II.
Ebből a táblázatból kiolvasható, hogy a globális load utasítások száma jelentősen megnőtt. Ez abból következik, hogy egy szálnak egy pixelnyi adat helyett kilenc pixel intenzitás értékeit kell betöltenie. Elméleti úton számolva ez kilencszeres növekedést jelentene a forgalomban, azonban a mért értékekben több mint tízszeres a növekedés. Ezt az eltérést mérési hibának tulajdonítom.
10.2 Globális memória A fejlesztőkörnyezethez mellékelt, a GPU-k programozásáról szóló dokumentumokban az optimalizálásról szóló részben első helyen szerepel a globális memória elérések hatékonyabbá tétele. Különösen igaz ez az olyan memória intenzív műveletekre, mint a képfeldolgozás. Az algoritmusok a jellemző szűk keresztmetszet szerint lehetnek számítás, vagy memória intenzívek. A GPU-ra történő optimalizálás során fontos már az elején meghatározni, hogy az aktuális algoritmus melyik csoportba esik, mert a két esetben némileg más feltételeknek kell megfelelni. Korábban említettem, hogy kihasználtság csökkenése nem feltétlenül eredményezi a teljesítmény csökkenését. Ez azonban csak a memória intenzív műveletekre igaz, ahol a nagy késleltetésű load/store utasításokat át lehet hidalni a warpok cserélgetésével. A számítás intenzív műveleteknél sokkal kritikusabb a kihasználtság magasan tartása. Az egyes képkezelő műveletek viszont memória intenzívek, különösen a példának használt binomiális szűrő, amelyben a pixelkoordinátáknál használt fixpontos szorzáson és összeadáson kívül csak shiftelés található. Ellenben ezeket a műveleteken nagy mennyiségű adaton kell elvégezni, tehát az algoritmus futási idejét a memória elérés hatékonyabbá tételével lehet a legnagyobb mértékben csökkenteni. A kérdés már csak az, hogy mit értünk hatékony memória elérésen? A válasz függ a használt architektúrától: a régebbi processzorok valamelyest szigorúbb feltételeket
Oldal 52
szabnak. Én a tesztben használt grafikus egységre szeretnék optimalizálni, így ezt mutatom be részletesebben. A kérdés megválaszolásához először meg kell érteni, hogy a globális memória hogyan épül fel. A globális memóriát az elérés szempontjából, logikai szegmensekre bonthatjuk. A memóriavezérlő egyszerre csak egy szegmensnyi adatot tud mozgatni. Egy szegmens háromféle méretű lehet: 32, 64 vagy 128 bájtos. A másik fontos dolog a memória kérések ütemezése. Említettem már, hogy a SIMT egység a blokkokat 32 szálat tartalmazó csoportokra, warpokra bontja szét. A multiprocesszor nyolc darab feldolgozó egysége azonban egyszerre csak 16 szálat tud kezelni, ezért a warpokat középen elfelezi. Ez a 16 szál, a fél warp, ha nincs elágazás, mindig ugyanazt az utasítást hajtja végre, ezeket tekinthetjük aktív szálaknak. Ha ez egy globális load/store utasítás, akkor a 16 kérés egyszerre továbbítódik a memóriavezérlő felé.
38. ábra Szegmensek I. [12]
Azt, hogy éppen mekkora a szegmens méret az egyes szálak által kért bájtok mennyisége határozza meg: egy bájt esetén a szegmens 32 bájtos, két bájt esetén 64, négy vagy nyolc bájt esetén 128. Az egyszerre beérkező 16 kérés kiszolgálásának folyama: 1. a vezérlő megkeresi a memóriában azt a szegmenst, amelyből az aktív szálak közül a legkisebb számozású kér adatot, 2. a többi aktív szál közül megkeresi azokat a szálakat, amelyek ugyanebből a szegmensből kérnek, és betölti a szükséges adatokat, 3. azokat a szálakat, amelyek megkapták az adatot, inaktívnak jelöli, és visszaugrik az első pontra. A fentebb ismertetett kernelben egy szál egy bájtot kér, hiszen mindig az aktuális pixel egy csatornájának intenzitását kérdezi le. Ebből következik, hogy a hatékony kihasználás érdekében az egyszerre továbbított 16 kérésnek egy 32 bájtos szegmensbe kell esnie. Vizsgáljuk meg, hogy mi történik egy képpontnyi adat, azaz három bájt lekérdezésekor: sumB += srcImageData[actElem]; sumG += srcImageData[actElem + 1]; sumR += srcImageData[actElem + 2];
Oldal 53
A kék csatorna lekérdezésekor minden képpont első bájtjára van szükségünk. Ez 16 szál esetén pontosan 16 bájtnyi adat, azonban a kék csatornák nem egymás mellett helyezkednek el a memóriában. Emiatt legalább 16 × 3 = 48 bájt adatot kell mozgatni, de mivel a memóriavezérlő csak szegmenseket tud kezelni, ez legalább két 32 bájtos olvasási műveletet jelent. A zöld és a piros csatornák ugyanígy működnek, így 16 képpontnyi adat, azaz 48 bájt helyett 3 × 64 = 192 bájtnyi adat mozog. Ez a memória átbocsátóképesség pazarlását jelenti. A hatékony memória elérés alatt tehát a teljes sávszélesség minél hatékonyabb kihasználást értem. A cél az, hogy egy lekérdezett, vagy visszaírt szegmensnyi adat teljes tartalma ki legyen használva. 10.2.1 Beépített struktúra használata
Első megközelítésben tehát egy szál, egy képpont feldolgozásáért felel. Triviális, hogy gyorsíthatnánk azzal, ha az egyes csatornákat nem külön-külön kérdeznénk le, hanem egy struktúrát felhasználva egyben. Ehhez még csak nem is kell külön struktúrát írni, a CUDAban rengeteg előre megírt struktúra van az ilyen esetek kezelése. Szinte minden alap adattípushoz létezik egy-, kettő-, három- és négyelemű vektor típus. Ez a módosítás a hoszt kódot nem érinti, ugyanannyi szálat kell indítani. A kernel kódjában az egy képpont adatainak lekérdezése, és a visszaírás módosul: uint3 sum = make_uint3(0, 0, 0); uchar3 actElem; if (X > 0 && Y > 0) {//top left actElem = src[(Y - 1) * width + (X - 1)]; sum.x += actElem.x; sum.y += actElem.y; sum.z += actElem.z; }
A visszaírás során, ugyanígy uchar3 típusú változót használok. Ezzel a lépéssel tehát három sor helyett egy sorban kérdezi le a kernel a három bájtot. Vizsgáljuk meg, hogy ez milyen változásokat hozott a globális memória elérések számában:
Kernel
Futási idő
Globális Load 32B
Globális Store 32B
Globális Load kérés
Globális Store kérés
416,96
102753
9792
7344
816
8. táblázat Mérési eredmények beépített struktúra használatával
Az értékek majdnem teljesen ugyanazok, mint az előző példában, a hiba betudható mérési hibának. Ennek az okára nem találtam egyértelmű leírást. Az biztos, hogy a fejlesztőkörnyezethez mellékelt dokumentációk nem tartalmaznak leírást arra vonatkozóan, hogy hogyan alakul a memóriaelérés hárombájtos lekérdezések esetén. Magyarázat csak egy-, kettő- vagy négybájtos esetre van. Ebből, és mért paraméterek hasonlóságából arra következtetek, hogy ilyen esetben hardver szálanként három egybájtos kérést továbbít a memóriavezérlőhöz.
Oldal 54
10.2.2 Négy bájt olvasása, három bájt írása
A hárombájtnyi adat írása, olvasása tehát nem vált be. Következő ötletként felmerült, hogy az egy képpontnyi adat helyett mozgassunk 4 bájtot. Ez magával hozza azt, hogy az egy szál által lekérdezett adat „átnyúlik” a szomszédos pixelre is, azaz lekérdezésenként egy bájttal több adat kerül mozgatásra a szükségesnél. További hátránya, hogy írásnál továbbra is meg kell elégednünk három bájttal, a szomszédos értékeket nem írhatjuk felül. Előnye viszont, hogy így az olvasások, bár több adatot mozgatnak, hatékonyabban tudnak teljesülni. Ez a megoldás továbbra sem hoz változást a hoszt kódban, és a visszaírás sem módosul, egyedül a globális memória olvasását kell módosítani: uchar4 actElem; uchar3* src3 = (uchar3*)srcImageData; uchar4* src4 = (uchar4*)srcImageData; if (X > 0 && Y > 0) {//top left actElemCount = (Y - 1) * width + (X - 1); src4 = (uchar4*)(&src3[actElemCount]); actElem = src4[0]; sum.x += actElem.x; sum.y += actElem.y; sum.z += actElem.z;
} A lekérdezést némileg megbonyolította a változtatás. A könnyebb átláthatóság kedvéért bevezettem egy új változót: az „actElemCount” az aktuálisan lekérdezendő pixel sorszámát jelenti. Első lépésként a háromelemű vektor segítségével meghatározom, hogy az adott sorszámú elem milyen címen található a memóriában. A konkrét lekérdezéshez viszont már a négyelemű vektort használom, és a negyedik intenzitás értéket egyszerűen figyelmen kívül hagyom. Lássuk, ez milyen változást hozott: Futási idő Kernel
224,16
GL 32B
GL 64B
GL 128B
GS 32B
GL kérés
GS kérés
19760
7482
3645
9888
3360
840
9. táblázat Mérési eredmények: négy bájt olvasása, három írása
A táblázatban megjelent két új oszlop, ill. bevezettem néhány rövidítést, azért hogy elférjenek egy táblázatban. A GL a globális load, a GS pedig a globális store műveleteket jelenti. Mivel most már ténylegesen négy bájt adatot olvas egy-egy szál, ahogy már láttuk, az memóriavezérlőhöz egyszerre beérkező 16 kérésnek egy 128 bájtos szegmensbe kell esni. A memóriavezérlő azonban képes az elérést optimalizálni: ha a 16 kérés egy 64 bájtos szegmense esik, akkor csak 64 bájtot mozgat (39. ábra, a), ha viszont átcsúszik két 64 bájtos szegmensbe, amik egy 128 bájtos területet alkotnak akkor egy 128 bájtos tranzakció lesz az eredmény (39. ábra, b).
Oldal 55
39. ábra Szegmensek II. [12]
Ha azonban a két 64 bájtos szegmens két 128 bájtoshoz tartozik, akkor egy 64 bájtos és egy 32 bájtos adatcsomagot kérdez le a memóriavezérlő (40. ábra).
40. ábra Szegmensek III. [12]
A probléma a kernellel az, hogy ugyan négybájtos vektorokat olvasnak a szálak, de ezeket a memória minden harmadik bájtjáról, ezért fordulhat elő, hogy nem egy 64 bájtos szegmensbe esnek a kérések. Mindazonáltal az olvasást így is sikerült hatékonyabbá tenni, ami meglátszik mind a kérések, mind a továbbított szegmensek mennyiségén. Az írás továbbra sem változott, a lekérdezések gyorsításával viszont majdnem kétszeres gyorsítást sikerült elérni. 10.2.3 Textúra alkalmazása
A mai grafikus egységeket, ahogy a nevük is mutatja, elsősorban 3D-s képalkotásra tervezték. Ezért nem meglepő, hogy külön hardveregységek találhatók bennük, amik kifejezetten a képekkel, képalkotással kapcsolatos műveleteket gyorsítják. Ilyen a már említett textúra memóriatér, és a hozzá tartozó gyorsító tár is. A textúra memóriatér csak logikailag különbözik a globálistól, fizikailag ugyanúgy a gyorsító- vagy videokártyán található DRAM-ból lehet elérni. Működésének lényege, hogy a már memóriában lévő adatra egy textúra referenciát hozunk létre, amelyen keresztül a chipbe integrált gyorsító tár hatékonyan tudja menedzselni a memória busz forgalmát. Hátránya, hogy írásra nem lehet használni, így az eredmények elmentését továbbra is bájtonként kell megoldani. Használata már befolyásolja a hoszt kódot is. Első lépésként deklarálnunk kell magát a referenciát: texture texRef;
Ezután létre kell hozni egy ún. CUDA tömböt, amibe majd a képi adatot fogjuk másolni, mert a referencia csak egy ilyen típusú tömbre mutathat: cudaArray *d_imageArray;
Oldal 56
cudaChannelFormatDesc U8ChannelDesc = cudaCreateChannelDesc(); cudaMallocArray(&d_imageArray, &U8ChannelDesc, byteWidth, height);
Az allokáló függvény utolsó két paraméterével adhatjuk meg a lefoglalni kívánt terület szélességét és magasságát, a „byteWidth” a kép a szélességét jelöli bájtban. A következő lépésben az eddigiekben alkalmazott másoló függvény egy változatával tudjuk a tömböt feltölteni adattal: cudaMemcpy2DToArray(d_imageArray, 0, 0, byteWidth, height, cudaMemcpyHostToDevice);
h_imageData,
widthStep,
A kernel indítása előtt az utolsó feladatunk, hogy a textúra referenciát hozzá kössük a tömbhöz: cudaBindTextureToArray(texRef, d_imageArray);
A kernel lefutása után a már fentebb ismertetett teendőkön kívül még szét kell kapcsolni a textúrát a tömbtől, és nem szabad elfelejtkezni a tömb felszabadításról sem: cudaUnbindTexture(texRef); cudaFree(&d_imageArray);
A kernelt is módosítani kell, a textúra gyorsító tár eléréséhez, egy külön függvény használatos: float xTex = (float)X; float yTex = (float)Y; if (X > 0 && Y > 0) {//topleft sumB += tex2D(texRef, xTex - 3.0f, yTex - 1.0f); sumG += tex2D(texRef, xTex - 2.0f, yTex - 1.0f); sumR += tex2D(texRef, xTex - 1.0f, yTex - 1.0f); }
A függvény első paramétere maga a referencia, a második és a harmadik pedig az X és Y koordináták lebegőpontos típusban. A visszaírás most sem módosult.
Kernel
Futási idő
Regiszterek Kihasználtság GL
GS 32B
GL kérés
GS kérés
202,24
18
9888
0
840
0,75
0
10. táblázat Mérési eredmények textúrával
A 10. táblázat tartalmazza a mért paramétereket. A textúra használata esetén a globális memória kérések száma nulla, azaz a teljesítménymérő alkalmazás nem tudja mérni a memória busz forgalmát, ha az a gyorsító táron keresztül történik. Ebből következik, hogy a 32, 64 vagy 128 bájtos adatcsomagok száma is nulla. Látható, hogy a visszaírás nem módosult. A szálak által használt regiszterek száma viszont annyira megemelkedett, hogy ez már befolyásolta a kihasználtságot: visszaesett 75%-ra, azaz a maximális 32 warpból csak 24 lehet aktív. A kernel futási ideje azonban így is csökkent valamelyest. 10.2.4 12 bájt
Az eddigi optimalizálások során csak az olvasásra koncentráltam. Ahhoz, hogy az írás hatékonyabb legyen, a textúra nem használható, hiszen az írás művelet nem megengedett. Azt kellene elérni, hogy a szálanként egyszerre lekérdezett és a visszaírt adat is egy-, kettőOldal 57
vagy négybájtos legyen. Ennek gátat szab, hogy egy RGB pixel három bájtot foglal el. Az ötlet az, hogy egy szál annyi adattal dolgozzon, amennyinél már megoldható, hogy négybájtosával lehessen olvasni és írni is. A megoldás a három és a négy legkisebb közös többszöröse a 12. Ha egy szál négy pixelt dolgoz fel, akkor a 12 bájt adat lekérdezhető és visszaírható háromszor négybájtos csomagokban. A hoszt kódban ez csak annyi változást hoz, hogy X irányban negyedannyi szálat kell indítani, a kernelt viszont jelentősen megbonyolítja, és több regiszterre is lesz szükség.
Kernel
Futási idő
GridX
GridY
Regiszterek Kihasználtság
166,944
8
32
24
0,5
11. táblázat Mérési eredmények 12 bájt használatával I.
Látható, hogy valóban több regisztert használ egy kernel, ami most már jelentősen lecsökkentette a kihasználtságot. Ezzel együtt a futási idő tovább csökkent, ami bizonyíték arra az állításra, hogy a memória intenzív műveleteknél a kihasználtság magasan tartása nem kritikus, sokkal fontosabb a memóriaelérés optimalizálása. GL 32 B Kernel
912
GL 64 B 4988
GL 128 B
GS64 B
GS 128 B
GL kérés
GS kérés
7482
1248
1248
1080
216
12. táblázat Mérési eredmények 12 bájt használatával II.
Az eredményekből kiolvasható, hogy ez a módosítás nem csak az írást, de az olvasást is hatékonyabbá tette. Az író műveleteknél teljesen megszűntek a 32 bájtos tranzakciók, és ezzel párhuzamosan a kérések száma is jelentősen leesett. A load műveleteknél a kérések száma kb. harmadára esett, és a számarányukat tekintve is túlsúlyba kerültek a 128 bájtos átvitelek.
10.3 Osztott memória Az eddigi optimalizálások során a globális memória összpontosítottam, és figyelmen kívül hagytam a chipre integrált osztott memóriát. Ahogy már szó volt róla, az osztott memória logikailag egy blokkhoz, fizikailag egy multiprocesszorhoz tartozik, azaz az egy blokkban lévő szálak tudnak rajta keresztül kommunikálni. Késleltetése egy-két ciklus, mérete minden MP-ben 16 KB, azaz összesen 480 KB található a lapkán. Használata némileg átalakítja kernel szerkezetét.
41. ábra A kernel szerkezete osztott memória használata esetén
Az eddigiekben a lekérdezés egybe esett a számolással. Az osztott memória használata ezzel szemben két jól elkülöníthető részre osztja a kernelt: a szinkronizálásig csak olvassuk az adatokat, utána pedig a feldolgozás után visszaírjuk. A kernelben használható szinkronizáló utasítás addig blokkolja a szál futását, amíg az azonos blokkban található szálak mindegyike el nem éri ezt a pontot. Erre azért van szükség, hogy a feldolgozás Oldal 58
során biztosak lehessünk abban, hogy az osztott memóriában minden szükséged adat rendelkezésre áll. Az eddig kipróbált módszerek megegyeztek abban, hogy minden szál kilencszer annyi adatot kérdezett le, mint amennyit visszaírt. Ez rengeteg redundáns adatmozgatást jelent, hiszen egy adategységet így kilencszer kell kiolvasni. Az osztott memóriát tehát fel lehet arra a használni, hogy minden szál – kivétel a blokk szélén lévő szálak – csak egy egységnyi adatot kérdezzen le, mert a feldolgozás során a szomszédos pixelek által lekérdezett adatokat is fel tudják használni. További előnye, hogy mivel az olvasás elkülönül a feldolgozástól, lehetőségünk van más sorrendben lehívni az adatokat, így is hatékonyabbá téve globális memória elérését. Használata a hoszt kódot csak annyiban változtatja meg, hogy a kernel indításánál meg kell adni a blokkonként felhasználandó osztott memória mennyiségét. Azt, hogy pontosan mekkora legyen ez a mennyiség, a blokkméretből kiszámítható:
9. egyenlet Osztott memória méretének számítása
.A számításnál figyelembe kell venni a blokkok határait is. Így a 16 × 16 szálat tartalmazó blokknak 3600 bájt osztott memóriára van szüksége. Ha egy változót az osztott memóriában szeretnénk elhelyezni, akkor ezt a „__shared__” kulcsszóval tudjuk megtenni: __shared__ uchar4 smem[50][18];
10.3.1 Redundancia megszüntetése
Először szüntessük meg a redundáns adatolvasásokat. Ennek első lépése, hogy minden szál lehívja a saját adatát. Ez 3 × 4 bájt olvasás szálanként: int x3 = threadIdx.x * 3; actElemCount = (Y) * width smem[x3 + 1][threadIdx.y + smem[x3 + 2][threadIdx.y + smem[x3 + 3][threadIdx.y +
+ (X); 1] = src[actElemCount]; 1] = src[actElemCount + 1]; 1] = src[actElemCount + 2];
Ezután kezelni kell a blokk széleit is, mert a szélső szálaknak olyan adatokra is szüksége van, amelyek blokk határain kívül esnek. Ha ez megvan, akkor a szinkronizáló utasítás segítségével blokkolni kell a program futását addig, minden szál be nem fejezte az olvasást. A feldolgozás és a visszaírás nem változik. GL 32 B Kernel
807
GL 64 B
GL 128 B
GS64 B
GS 128 B
GL kérés
GS kérés
1398
1398
1248
1248
416
216
13. táblázat Mérési eredmények: osztott memória I.
Látható, hogy az olvasások száma tovább esett, különösen feltűnő a 64 és 128 bájtos műveletek számának csökkenése. A futási valamelyest gyorsult: 159,584 µs. 10.3.2 Adatlehívás sorrendjének megváltoztatása
Az olvasás és a feldolgozás időbeni elkülönülése lehetőséget biztosít a globális memória még hatékonyabb kihasználására. Azzal, hogy megnövekedett az egy szál által feldolgozott adatmennyiség, sikerült előrelépést elérni az olvasásban és az írásban egyaránt. Ez azonban Oldal 59
újra előhozta azt a problémát, hogy a memóriavezérlőhöz egyszerre beérkező 16 utasítás által kért adatok, nem egy egységes memória területet foglalnak le. Az osztott memória segítségével lehetőség nyílik arra, hogy az adatokat ne olyan sorrendben hívjuk le, amilyen sorrendben szükség van rá: unsigned int base = Y * width + (blockIdx.x * 3) * blockW; smem[1 + threadIdx.x][1 + threadIdx.y] = src[base + threadIdx.x]; smem[17 + threadIdx.x][1 + threadIdx.y] = src[base + 16 + threadIdx.x]; smem[33 + threadIdx.x][1 + threadIdx.y] = src[base + 32 + threadIdx.x];
Ehhez meg kell határozni egy bázis értéket, melynek segítségével a blokk adott sorának az első elemét lehet megcímezni. Ezután a bázishoz a szál blokkon belüli azonosítójának X komponensét hozzáadva az első 16 lekérdezés sorba rendezhető. A továbbiakban bázist kell 16-osával eltolni. Az osztott memóriát X irányban szintén a szál azonosítójával, megfelelően eltolva kell címezni. Vizsgáljuk meg az eredményeket.
Kernel
GL 32 B
GL 64 B
GL128 B
GS64 B
GS 128 B
GL kérés
GS kérés
807
1398
0
1248
1248
416
216
14. táblázat Mérési eredmények: osztott memória II.
Látható, hogy ugyanannyi kérés érkezett a memóriavezérlőhöz, de nem volt szükség arra, hogy egy egész 128 bájtos szegmens mozgatásra kerüljön. A 32 bájtos olvasás a blokk széleinek kezelése miatt maradt meg. A futási időt ez már drasztikusan nem változtatta: 154,592 µs.
10.4 Optimalizálás számításra Bár a képkezelő műveletek alapvetően memória intenzívek, nagyobb késleltetésű aritmetikai műveletek – osztás, szorzás, logaritmus stb. – is előfurdalnak. Az alapvetőbb műveletekhez a CUDA biztosít beépített függvényeket, melyekkel (bizonyos korlátozások mellett) gyorsíthatók ezek a számolások. A legtipikusabb ilyen a fixpontos szorzás. Erre minden kernelnek szüksége van, ahhoz hogy meghatározza, melyik adategységet kell feldolgoznia. A fixpontos szorzás áteresztőképessége 2 utasítás / ciklus / MP. Ha a 32 bites egész típusú változóból elég 24 bit, akkor használhatjuk az _umul24() függvényt, amellyel az áteresztőképesség 8 utasítás / ciklus / multiprocesszorra módosul. Ugyanígy léteznek függvények egyszeres pontosságú lebegőpontos műveletekhez, szögfüggvényekhez, gyökvonáshoz stb.
10.5 Egyéb algoritmusok A binomiális szűrővel történő zajszűrés csak az első lépése a sejtmagdetektálásnak. Ahhoz hogy a GPU minél hatékonyabban működjön, érdemes minél több algoritmust implementálni rá, így csökkentve a lassú csatolósín használatát. A sejtmagdetektálás következő lépései az RGB-HSV konverzió és a binarizálás. A detektáló folyamat hatékonyságának növelése érdekében célszerű, a binomiális szűrő eredményét a GPU memóriában hagyni, és a következő lépéseket ezen elvégezni. A folyamatot az alábbi szemlélteti.
Oldal 60
42. ábra Sejtmagdetektálás GPU-val
A zajszűrés utáni két lépést össze lehet vonni egy kernelbe, így a kapott HSV értékek elmentésének kihagyásával globális memória forgalmát lehet csökkenteni. A színterek közötti konverzió miatt ez a kernel több lebegőpontos utasítást is tartalmaz, ami lehetőséget biztosít a fentebb említett függvények hatásának vizsgálatára. A két kernel fontosabb mért paramétereit az alábbi táblázat tartalmazza: Futási idő Regiszterek
Kihasználtság
Utasítások
Kernel1
154,144
24
0,5
54 871
Kernel2
971,936
45
0,25
271 959
15. táblázat Lebegőpontos függvények használatának eredménye
Az első kernelben hagyományos lebegőpontos műveleteket alkalmaztam, míg a másodikban ezeket kicseréltem az aritmetikai függvényekre. Az eredmény meglepő, a futási idő több hatszorosra emelkedett. Ennek legfőbb oka, hogy a megnövekedett regiszterhasználat miatt a kihasználtság a kritikus szint alá csökkent, de érdekes még a végrehajtott utasítások kb. ötszörös emelkedése is. További hátránya a kernelnek, hogy a hagyományos műveleti jelek hiánya átláthatatlanná teszi a kódot.
Oldal 61
11 Paraméterek mérése Miután detektáltuk a sejtmagokat, mirigyeket és hámrétegeket meg kell állapítani ezen szöveti struktúrák tulajdonságait, amelyek alapján meg tudjuk különböztetni egymástól a beteg mintát az éptől. Mivel nem ismerünk olyan tulajdonságot melynek tudatában egyértelműen meg lehet határozni a szövet állapotát, több paramétert kell vizsgálni. Minden paraméter estében meghatározunk egy minimum és maximális értéket, melyeket nem léphet át a kapott eredmény egészséges minta esetében. Amennyiben valamelyik érték kilép ebből a tartományból, úgy azt a programnak jeleznie kell az orvosnak, aki meg tudja állapítani, hogy az eltérés tekingető normálisnak, vagy a minta valóban beteg. Mivel a beteg minták elég változatosak lehetnek és mivel bizonyos elváltozások még épp minta esetén is előfordulhatnak, érdemes minél több paramétert figyelembe venni. Valószínű, hogy beteg minta esetén több paraméter érték is eltérő lesz az átlagtól.
11.1 Megvalósított paraméterek ismertetetése Bizonyos paramétereket meg kell mérnünk detektálás fázisban. Ilyen paraméter a terület és bizonyos esetekben a kompaktság is. Minden szövet struktúrának saját tulajdonságai vannak így a paramétereket szintenként kell meghatározni. Sejt szinten már a detektálásnál megmértük a sejt területét a Kitöltő algoritmus, kerületét meghatároztuk a Körvonal detektáló című fejezetben leírtak szerint, valamint a legkisebb és legnagyobb átlókat kiszámítottuk az Átmérők megállapítása [21] című fejezet szerint. Ezek után érdemes lenne meghatározni minden sejtmag befoglaló téglalapját (bounding box). Befoglaló téglalap meghatározásakor a korábban már rögzített külső pontokra X és Y koordináták szerint egy maximum és egy minimumkeresést kell végrehajtani. Így megkapjuk a befoglaló téglalap bal felső és jobb alsó sarkát. A sejtmagok gömbszerűek, így ha ezt bármilyen síkban metsszük a kompaktsági ért kékek csak kis mértékben térhetnek el. Amennyiben bizonyos területen nagyon eltérő kompaktságú sejtmagokat találunk, zajra vagy valamilyen hibára következtethetünk. Kompaktság mérése a fentebb leírt módon történik. Sejtmag szintem ki kell továbbá számolni néhány orvosi szempontból fontos statisztikai értéket is. Meg kell határozni az átlagos sejtnagyságot, ennek szórását, felszín és sejtmag arány, kerület sejtmag arány. Ezen mérőszámok meghatározása egyszerű mivel kiszámításukhoz már minden információt ismerünk. Mirigy detektálása után szintén rendelkezésünkre áll néhány paraméter. A mirigy detektálás során felépítjük a sejtmagtérképet majd ebből a hálót. A háló után minden különálló területet megvizsgálunk kitöltő algoritmussal. A mirigyeknek így ismerjük a területét. A detektálási folyamat végén minden mirigynek kiszámoljuk a kompaktágát, így ezt is rögzíthetjük, mint újabb paraméter. Mirigy esetében is ki kell számolni a befoglaló téglalapot és a legkisebb és legnagyobb átlókat a fentebb ismertetett módon. Hám szinten meghatározzuk a detektált hám területét, megvizsgáljuk hogyan aránylik a teljes területhez képpest, valamint megnézzük ilyen arányban fedi le a teljes minta kerületét. Oldal 62
12 Tesztelés A program célja a szöveti komponensek detektálása. Ennek megfelelően a tesztelési fejezetben a szöveti komponensek detektálásának hatásfokát fogjuk vizsgálni. 12.1.1 Sejtmagok detektálása
Mivel a sejtmagokat szemmel kell megszámolni és leellenőrizni valamint összeszámolni érdemes kis méretűk képet vizsgálni.
43. ábra kis tesztkép
RGB binarizálást alkalmazva az program 428 sejtmagot számolt meg. Szemmel megszámolva 313 sejtmagot találtam. Ezekből az adatokból kiszámolhatjuk, hogy 36%-kal több sejtmagot talált a program minta amennyi van a képen. A két eredmény között 115 sejtmagnyi különbség van ezek közül 62 darab sejtmagot detektált a mirigyek belsejében. Itt nincsenek sejtmagok, azonban a belső szerkezetük megzavarja az algoritmust. A maradék 53 tévesen detektált sejt az Ultimate Erode algoritmus miatt jelentkezik. A nagyobb sejteket úgy hajlamos kettébontani így ezeket a sejteket kétszeresen detektálja.
44. ábra sejtmagdetektálás RGB binarizálással
HSV detektálás estében 313 sejtmag helyett 308 sejtmagot detektált. Itt azonban meg kell jegyezni, hogy az Ultimate Erode itt is produkál duplázódást. A képen 19 sejtmagot találtam ahol kétszer találta meg a sejtmagot. Az algoritmus 24 sejtmagot nem talált meg. A nem detektált sejtmagok nagy része a mirigyek falában helyezkednek. Az algoritmus a nagymértékű összetapadás miatt nem találja meg. A HSV használatával ki tudtuk küszöbölni, hogy a mirigy közepében is detektáljon sejtmagokat.
Oldal 63
45. ábra sejtmagdetektálás HSV binarizálással
RGB automatikus detektálás esetén 313 sejtmag helyett 536 sejtmagot detektált. A 223 sejtmaggal talált többet, mint amennyi valójában van a képen. Ebből 85 tévesen megjelölt terület van a sejtmag közepében a maradék 138 téves detektálást az Ultimate Erode algoritmus generál.
46. ábra sejtmagdetektálás automatizált RGB binarizálással
Algoritmusok
Sejtmagok száma
Talált sejtmagok
Tévesen detektált sejtmagok száma
RGB Threshold
313
428 (136%)
62 (20%)
HSV Threshold
313
308 (98%)
19 (6%)
RGB Automatikus
313
536 (171%)
85 (27%)
Algoritmusok
Sejtmagok száma
Talált sejtmagok
Tévesen detektált sejtmagok száma
RGB Threshold
122
150 (123%)
6 (5%)
HSV Threshold
122
140 (115%)
3 (2%)
RGB Automatikus
122
260 (213%)
32 (26%)
Oldal 64
Algoritmusok
Sejtmagok száma
Talált sejtmagok
Tévesen detektált sejtmagok száma
RGB Threshold
158
181 (114%)
13 (8%)
HSV Threshold
158
188 (118%)
11 (7%)
RGB Automatikus
158
244 (154%)
28 (18%)
Algoritmusok
Sejtmagok száma
Talált sejtmagok
Tévesen detektált sejtmagok száma
RGB Threshold
166
200 (120%)
7 (4%)
HSV Threshold
166
182 (109%)
4 (2%)
RGB Automatikus
166
284 (171%)
9 (5%)
Algoritmusok
Sejtmagok száma
Talált sejtmagok
Tévesen detektált sejtmagok száma
RGB Threshold
284
388 (136%)
27 (10%)
HSV Threshold
284
368 (129%)
24 (8%)
RGB Automatikus
284
481 (169%)
39 (14%)
Algoritmus
Talált sejtmagok átlaga
számának Tévesen detektált sejtmagok számának átlaga
RGB Threshold
126%
9%
HSV Threshold
114%
5%
RGB Automatikus
176%
18%
12.1.2 Mirigy detektálása
A mirigyek detektálása a sejtmagtérképre épül. A sejtmagtérképből felépíti a hálót és megkeresi azokat a helyeket, ahol mirigy lehet. Minden lehetséges területnél két paraméter alapján döntünk méret és kompaktsági faktor. Amennyiben ezeknek megfelel a vizsgált maszk, feljegyezzük, mint mirigy.
Oldal 65
47. ábra tesztkép mirigydetektáláshoz
48. ábra detektált mirigyek a tesztképen
A fenti képen látható a program által detektál mirigyeket. A képen 109 mirigy található a program 114 mirigyet talált, de ebből 12 hibás találat volt és 7 mirigyet a program nem tudott detektálni. A program nem tudja detektálni a mirigyeket, ha azok túl kicsik, ha roncsolódtak, illetve ha azok csak részben láthatóak a képen. A téves detektálás abban az esetbej jelentkezik, ha a mirigyek között nem túl nagy a távolság, de ez a szabad tér terület megfelelhet egy átlagos mirigy területével és a sejtsűrűség alacsony.
Oldal 66
Kép
Átlag:
Mirigyek száma
Program által Hibásan detektált detektált mirigyek mirigyek száma száma
Nem detektált mirigyek száma
109
114 (105%)
12 (11%)
7 (6%)
53
52 (98%)
3 (6%)
4 (8%)
30
42 (140%)
13 (43%)
1 (3%)
114%
20%
6%
Oldal 67
13 Értékelés Tesztelés fejezetben meghatároztuk milyen hatásfokkal sikerült detektálnunk a sejteket és a mirigyeket. Sejtmagok detektálásánál a talált sejtmagok száma mindig több mint a valós sejtmagok száma. Ez azért lehetséges, mert az összeragadt sejtmagok szétválasztására írt Utlimate Erode nevű algoritmus egyszerű, nem összeragadt sejtmagokat is szétválaszt. Az algoritmus hátránya hogy nagyon érzékeny a zajra. Amennyiben a sejtmag maszkja csak 12 pixellel vékonyabb valahol már összetapadt sejtmagként értelmezi és szétbontja. A probléma részben kezelhető, ha minden sejtmagot megvizsgálunk, és amelynek a mérete nem sokkal nagyobb mint egyetlen átlagos sejtmag, akkor nem boltjuk fel. Az implementált sejtmagdetektáló algoritmusok 5% -18% között detektálnak hibásan sejtmagot. A hibásan detektált sejtmagok nagy része a mirigyekben található. A mirigyekben nem található sejtmag, de a belső szerkezete bizonyos esetben megtéveszti az algoritmust. Csökkenthetjük a hibásan detektált sejtmagok számát, ha a mirigydetektálás után töröljük azokat a sejtmagokat, amik egy mirigyben találhatóak. Átlagban a mirigydetektálásnál is többet kapunk eredményül, mint ahány mirigyünk valójában van a mintán. Jelen esetben azonban nem duplázódás okozza a problémát, hanem hibásan detektált mirigyek. Tesztek szerint a detektált mirigyek 30% hibásan detektált mirigy. A hibás detektálás a mirigyek közötti nagyobb kiterjedésű kötőszöveti réteg okozza. Ezeken a helyeken, ha a sejtmagsűrűség alacsony előfordulhat, hogy a sejtmagok valamilyen kör vagy ovális alakzatra hasonló alakzatban helyezkednek el. Ezeket a hibásan detektált mirigyeket ki tudjuk szűrni, ha minden detektált mirigyre készítünk egy hisztogramot és amennyiben a fehérhez közeli színekből sokat találunk, valószínűleg mirigyről van szó, különben a kötőszövet egy részét találata meg az algoritmus. Jelenleg a program futási ideje több minta esetén órák is lehetnek. A programban bizonyos algoritmus optimalizálásával illetve további algoritmusok GPU-ra implementálásával csökkenthetjük.
13.1 GPGPU 13.1.1 Kernelek összehasonlítása
Első lépésként foglaljuk össze különböző kernelek futási idejét. Típus
Naiv
Futási idő 425.344
Struktúra 4BO,3BÍ Textúra
12Byte
Osztott1
Osztott2
416.96
166.944
159.584
154.592
224.16
202.24
16. táblázat Különböző kernelek futási ideje
Oldal 68
49. ábra Különböző kernelek futási ideje
Látható, hogy a legnagyobb ugrást, kb. kétszerest, a naiv megközelítés elhagyása jelenti. További jelentősebb előrelépést a 12 bájtos feldolgozás bevezetése hoz, ez kb. 20%-kal gyorsabb, mint a textúra alapú feldolgozás. Az osztott memória használata azonban nem okozott jelentős gyorsulást. Mindent összevetve a leglassabb és a leggyorsabb megoldás között kb. 2,8-szoros különbség van. 13.1.2 GPU – CPU összehasonlítása
Az összehasonlítás során arra a kérdésre keresem a választ, hogy megéri-e használni a GPU-t képfeldolgozásra. A teszteléshez használt rendszer főbb jellemzőit az alábbi táblázat tartalmazza. Típus
Magok száma
Órajel
Mem. órajel
Mem. szélesség
CPU
Core2Duo
2
1866 Mhz
667 Mhz
128 bit
GPU
GTX295
240
1242 Mhz
999 Mhz
448 bit
17. táblázat Tesztrendszer főbb adatai
Bár a központi egység kétmagos, a feldolgozás során csak az egyik mag kerül felhasználásra, hasonlóan a videokártyához, ami két processzort tartalmaz, de csak az egyiket használom. A tesztelő algoritmus a sejtmagdetektálás (7. fejezet) első három lépése: a zajszűrés, az RGB-HSV konverzió és a binarizálás. A GPU-nál a leggyorsabb megvalósítást alkalmazom, míg a CPU-nál az OpenCV függvényeit használom. A tesztelés során az algoritmusok összesített futási idejét mérem különböző méretű képeken, a megjelenített eredmények legalább tíz futás átlagát jelentik. A grafikus egységnél mért idők tartalmazzák az összes hoszt függvény futási idejét is a memóriafoglalástól a felszabadításig.
Oldal 69
Képméret 256
512
1024
2048
4096
8192
10240
CPU
4.6
15.6
63.49
244
976
4371
7216
GPU
6.255
7.269
12.292
27.3
96.2
366
577
18. táblázat GPU és CPU összehasonlítása
50. ábra GPU és CPU összehasonlítása
Az összes tesztkép négyzetes alapú, a méret pixelben értendő. A megadott futási idők mértékegysége most is, és a későbbiekben is egy milliszekundum. A legkisebb és a legnagyobb kép mérete közötti jelentős különbség miatt, a jobb átláthatóság végett az időtartamok megjelenítéséhez logaritmikus skálát választottam. Látható, hogy a kisebb képeknél, az időigényes GPU-CPU másolási műveletek miatt, még nem éri meg a grafikus egységet használni. A GPU a legkisebb képet még lassabban dolgozza fel, a későbbiek folyamán azonban alkalmazásának előnyei egyre jobban megmutatkoznak. Vizsgáljuk meg képméretenként a futási idők arányát. Képméret
256
512
CPU/GPU 0.735412 2.1461
1024
2048
4096
8192
10240
5.165148 8.937729 10.14553 11.94262 12.50607
19. táblázat GPU és CPU futási idejének aránya
Oldal 70
51. ábra GPU és CPU futási idejének aránya
Az 50. ábra adataiból kiolvasható, hogy a CPU futási ideje végig arányos a képmérettel: négyszer akkora kép feldolgozásához kb. négyszer annyi idő szükséges. Ezzel szemben a GPU-n mért futási idők az első néhány képméretnél nem növekszenek a méretnövekedésnek megfelelően. Egy 512 × 512-es kép feldolgozásához csak kb. 16%-kal több idő szükséges, mint egy negyedakkora képhez. Ennek megfelelően a GPU előnye a kisebb méreteknél jelentősen növekszik, a 51. ábra grafikonja meredeken emelkedik, míg a nagyobb képeknél lassul a növekedés. Mivel sokkal nagyobb kép nem fér el az általam alkalmazott GPU memóriájában csak következtetni tudok a grafikon alakulásából, hogy az arány kb. 13 körül állandósul. A következőkben arra keresem a választ, hogy mi lehet annak az oka, hogy kisebb felbontásoknál kisebb különbségek jönnek ki. Ehhez vizsgáljuk meg, hogy a kép méretének növekedése hogyan befolyásolja a kernelek és az adatmásolás futási idejét. 256
512
1024
2048
4096
8192
10240
Kernel
0.116672
0.302624
1.076416
4.02016
16.17629
63.1225
96.4061
Másolás
0.1784
0.704256
3.18572
14.30096
66.8092
260.397
388.891
20. táblázat Kernel és másolás futási ideje
52. ábra Kernel és másolás futási ideje
Oldal 71
A megjelenített értékek a két kernel és másolás futási idejének összegét jelentik. A grafikonokról leolvasható, hogy a másolás a teljes tartományban arányos a képmérettel, míg a kernelek futási ideje a kisebb felbontásokon nem növekszik együtt feldolgozott adatmennyiséggel. Ebből azt a következtetést lehet levonni, hogy a GPU nagyobb adatmennyiség feldolgozás esetén tud hatékonyan működni. Ezt bizonyítja az a grafikon, amelyik a kernelek futási idejét a teljes futás százalékában fejezi ki (53. ábra).
53. ábra Kernel aránya
A grafikonról leolvasható, hogy a GPU nagyjából a 2048 × 2048-as felbontásnál éri el a maximális hatékonyságot, ezután a kernel futási idejének aránya 20% körül stabilizálódik. Minden összevetve elmondható, hogy bár már kisebb (512 × 512-es) képeket is gyorsabban lehet feldolgozni a grafikus egységgel, a GPU előnyei nagyobb felbontásban jelentkeznek igazán, ahol kb. 12-szeres különbséget sikerült elérni. 13.1.3 G80 – GT200 összehasonlítás
Érdekes lehet még két GPU architektúra összehasonlítása: az első grafikus egységet, amely képes a CUDA kód futtatására, hasonlítom össze az nVidia legfrissebb processzorával. A tesztelés során arra a kérdésre keresem a választ, hogy a fentebb vázolt optimalizáló technikák hogyan hatnak a futási időre egy régebbi architektúrán. A régebbi tesztkártya egy 8800GT, az újabb a fentebb is használt GTX295, a két típus paramétereit az alábbi táblázat tartalmazza. Magok száma
Mag órajel
Mem. órajel
Mem. szélesség
8800GT
112
1500 Mhz
900 Mhz
256 bit
GTX295
240
1242 Mhz
999 Mhz
448 bit
21. táblázat G80, GT200 paraméterek
A kernelek közül az utolsót nem tudta futtatni a régebbi processzor. A táblázatban látható értékek mértékegysége egy mikro szekundum, a grafikonon feltüntettem a GT200-as architektúrán mért futási időket is. Oldal 72
Típus
Naiv
Futási idő 5475.58
Strukt.
4BO,3BÍ Textúra
12Byte
Osztott1
Osztott2
5448.22
2831.17
813.12
489.952
-
694.464
22. táblázat Különböző kernelek futási ideje G80-as architektúrán
54. ábra Különböző kernelek futási ideje G80-as architektúrán
A kapott értékeket vizsgálva feltűnik, hogy a régebbi architektúrán sokkal nagyobb a különbség a naív és az optimalizált megvalósítások között. Ennek oka a globális memória elérés fejlődése. Az újabb processzort használva a nem hatékony kernelek kevésbé rosszak, mint a régebbi egységen. Ez okozza azt, hogy a G80-as architektúrán a leglassabb és leggyorsabb kernel közötti különbség kb. 11-szeres, az újabb kártyán kapott 2,8-szoros eltéréssel szemben, ill. hogy a textúra alkalmazása hatékonyabb, mint a 12 bájtos megvalósítás. A leggyorsabb kerneleket figyelembe véve a két processzor közötti különbség kb. 3-szoros.
13.2 Továbbfejlesztési lehetőségek Sejtmagok detektálásának javítása miatt az értékelésben leírt okok miatt érdemes lenne megvalósítani az Ultimate Erode algoritmust úgy, hogy vegye figyelembe a sejtmagok méretét. Kizárólag abban az esetben bontsa szét a sejtmagokat ha a maszk mérete jelentősen nagyobb mint egy átlagos sejtmag. További javítás elérése miatt, mirigydetektálás után felül kell vizsgálni a sejtmagokat. Amelyek valamelyik mirigy területén találhatóak, azokat el kell dobni. Mirigydetektálás javítása miatt módosítani kellene a mirigydetektáló algoritmust, hogy a detektálás után számolja ki minden mirigy hisztogramját és amennyiben ezek nem hasonlítanak egy átlagos mirigy hisztogramjára valószínűleg hibásan detektált mirigyről, kötőszöveti réteg egy részéről van szó. A minta állatának meghatározását segítheti, ha több paramétert implementálunk és mérünk. Ilyen paraméterek lehetnének a hám legkisebb és legnagyobb átlója és a mirigyek közötti legnagyobb és átlagos távolság. Oldal 73
Sok minta esetén problémát okozhat a futási idő. A futási időt csökkenthetjük, ha több algoritmust implementálunk GPU-ra is. 13.2.1 GPGPU
A grafikus egységet tekintve két fő továbbfejlesztési irány jelentkezik: a minél szélesebb körű alkalmazás, és további optimalizálások. Az bebizonyosodott, hogy a GPU-t hatékonyan lehet képfeldolgozó műveletek gyorsítására használni, tehát érdemes minél több algoritmust implementálni rá. Ilyen lehet pl. a sejtmagdetektálás következő lépése, az „ultimate erode” (7.4 fejezet) algoritmus megvalósítása. A jelenleg alkalmazott mirigy- és hámdetektáló algoritmusok viszont alapvetően soros természetűek, ezért vagy új megközelítést kell kitalálni, vagy a meglévőeket kell újragondolni, ahhoz hogy hatékonyan implementálhatóak legyen a grafikus egységen. A további optimalizálások terén három továbblépési irányt jelölnék meg: az osztott memória használatával megjelent ún. bankütközés jelenségének feloldása, több GPU kihasználása, egyéb, képfeldolgozásban használt optimalizációs eljárás bevezetése: pl. Gauss szűrő szeparálhatóságának kihasználása.
Oldal 74
14 Ábrajegyzék 1. ábra Vastagbél szöveti szerkezete [18] .............................................................................. 6 2. ábra Vastagbél nyálkahártyájának szerkezete [18] ............................................................ 7 3. ábra Hámréteg, és egy mirigy ** ....................................................................................... 7 4. ábra A mirigyek párhuzamossága megszűnik ................................................................... 9 5. ábra A mirigyekben elágazások jelennek meg .................................................................. 9 6. ábra Az izomréteg és a mirigyek alapja közti távolság megnő, és mirigyenként különbözik ........................................................................................................................... 10 7. ábra A generatív réteg a sejt felső részén is megjelenik .................................................. 11 8. ábra Kernelhívás .............................................................................................................. 16 9. ábra Szálhierarchia [6] ..................................................................................................... 17 10. ábra Blokkok kiosztása a magokhoz különböző GPU-knál [6]..................................... 17 11. ábra Az egyes memória fajták kapcsolata, elhelyezkedése [12] ................................... 18 1. táblázat GPU - CPU sávszélesség ................................................................................... 18 2. táblázat Memória fajták összefoglalása [12] ................................................................... 19 12. ábra Rendszerterv fejlődése I ........................................................................................ 21 13. ábra Rendszerterv fejlődése II ....................................................................................... 22 14. ábra Rendszerterv fejlődése III ...................................................................................... 23 15. ábra Képpiramis [13] ..................................................................................................... 24 1. egyenlet A képpiramis által maximálisan felhasznált tárterület ...................................... 24 16. ábra Navigáció ............................................................................................................... 25 2. egyenlet Slide kezelés: befoglaló négyzet mérete ........................................................... 25 3. egyenlet Bal felső képelem koordinátáinak számítása .................................................... 25 17. ábra Képrészlet kijelölése feldolgozásra ....................................................................... 26 4. egyenlet Bal felső képelem koordinátáinak kiszámítása az új nagyítási szinten............. 26 5. egyenlet Képelem szélessége az új nagyítási szinten ...................................................... 26 6. egyenlet Befoglaló téglalap sarkainak koordinátái .......................................................... 27 18. ábra CUDA használata .Net-ben ................................................................................... 28 19. ábra Sejtmagok .............................................................................................................. 29 20. ábra Az RGB modell ábrázolása [22] ............................................................................ 30 Oldal 75
21. ábra HSV modell két ábrázolása [22] ............................................................................ 30 22. ábra Összetapadt sejtmagok........................................................................................... 32 23. ábra Összetapadt sejtmagok binarizálása ...................................................................... 32 24. ábra Ultimate Erode ....................................................................................................... 33 25. ábra A kitöltő algoritmus működési elve....................................................................... 34 26. ábra A nyolcas és négyes szomszédság ......................................................................... 35 27. ábra A körvonal detektáló maszk ** ............................................................................. 35 28. ábra A mirigydetektáló algoritmus működése ** .......................................................... 36 29. ábra Sejtmagmaszk befoglaló téglalapban [21] ............................................................. 37 30. ábra Mirigyek, a hossztengelyükre merőlegesen elvágva ............................................. 38 7. egyenlet Kompaktság [25] ............................................................................................... 39 31. ábra Háló készítése [21] ................................................................................................ 40 32. ábra Hám ....................................................................................................................... 41 33. ábra A minta szakított oldala ......................................................................................... 42 34. ábra Egy CUDA-t használó program felépítése ............................................................ 43 8. egyenlet Grid dimenziók számítása ................................................................................. 44 35. ábra Egy képre feszített szálhierarchia .......................................................................... 45 3. táblázat Szegmentáló kernel mérési erdményei I. ........................................................... 49 4. táblázat Szegmentáló kernel mérési eredményei II. ........................................................ 49 5. táblázat Szegmentáló kernel mérési eredményei III. ....................................................... 50 36. ábra Konvolúció [26] .................................................................................................... 50 37. ábra Konvolúciós maszk................................................................................................ 51 6. táblázat Naiv zajszűrő kernel mérési eredményei I. ........................................................ 52 7. táblázat Naiv zajszűrő kernel mérési eredményei II. ....................................................... 52 38. ábra Szegmensek I. [12] ................................................................................................ 53 8. táblázat Mérési eredmények beépített struktúra használatával ....................................... 54 9. táblázat Mérési eredmények: négy bájt olvasása, három írása ........................................ 55 39. ábra Szegmensek II. [12] ............................................................................................... 56 40. ábra Szegmensek III. [12].............................................................................................. 56 10. táblázat Mérési eredmények textúrával ......................................................................... 57 11. táblázat Mérési eredmények 12 bájt használatával I. .................................................... 58 Oldal 76
12. táblázat Mérési eredmények 12 bájt használatával II. ................................................... 58 41. ábra A kernel szerkezete osztott memória használata esetén ........................................ 58 9. egyenlet Osztott memória méretének számítása .............................................................. 59 13. táblázat Mérési eredmények: osztott memória I. ........................................................... 59 14. táblázat Mérési eredmények: osztott memória II. ......................................................... 60 42. ábra Sejtmagdetektálás GPU-val ................................................................................... 61 15. táblázat Lebegőpontos függvények használatának eredménye ..................................... 61 43. ábra kis tesztkép ............................................................................................................ 63 44. ábra sejtmagdetektálás RGB binarizálással ................................................................... 63 45. ábra sejtmagdetektálás HSV binarizálással ................................................................... 64 46. ábra sejtmagdetektálás automatizált RGB binarizálással .............................................. 64 47. ábra tesztkép mirigydetektáláshoz ................................................................................. 66 48. ábra detektált mirigyek a tesztképen ............................................................................. 66 16. táblázat Különböző kernelek futási ideje....................................................................... 68 49. ábra Különböző kernelek futási ideje ............................................................................ 69 17. táblázat Tesztrendszer főbb adatai ................................................................................. 69 18. táblázat GPU és CPU összehasonlítása ......................................................................... 70 50. ábra GPU és CPU összehasonlítása ............................................................................... 70 19. táblázat GPU és CPU futási idejének aránya ................................................................ 70 51. ábra GPU és CPU futási idejének aránya ...................................................................... 71 20. táblázat Kernel és másolás futási ideje .......................................................................... 71 52. ábra Kernel és másolás futási ideje................................................................................ 71 53. ábra Kernel aránya ......................................................................................................... 72 21. táblázat G80, GT200 paraméterek ................................................................................. 72 22. táblázat Különböző kernelek futási ideje G80-as architektúrán .................................... 73 54. ábra Különböző kernelek futási ideje G80-as architektúrán ......................................... 73
Oldal 77
15 Irodalomjegyzék 1. Vance, Ashlee. Google shivs server crowd with PeakStream buy. The Register. [Online] 2007. 06 05. http://www.theregister.co.uk/2007/06/05/google_buys_peakstream/. 2. Universty, Stanford. BrookGPU. http://graphics.stanford.edu/projects/brookgpu/.
[Online]
2009.
3. Stokes, Jon. PeakStream unveils multicore and CPU/GPU programming solution. Ars Technica. [Online] 2006. 09 18. http://arstechnica.com/hardware/news/2006/09/7763.ars. 4. Rapidmind. RAPIDMIND. [Online] 2009. http://www.rapidmind.com/. 5. Pál, Röhlich. Szövettan. Budapest : Egyetemi Tankönyv, 1999. 6. NVIDIA. Volume1: Introduction to CUDA programming. CUDA ZONE - CUDA Education. [Online] 2008. http://www.nvidia.com/object/cuda_education.html. 7. —. Programming Guide 2.3. CUDA ZONE - CUDA Documentation. [Online] 2009. július 1. http://www.nvidia.com/object/cuda_develop.html. 8. —. NVIDIA Products. http://www.nvidia.com/page/products.html.
[Online]
NVIDIA,
2009.
9. —. NVIDIA GeForce GTX 280 / 260. NVIDIA Products. [Online] 2009. http://www.nvidia.co.uk/object/geforce_gtx_280_uk.html. 10. —. DirectCompute Support on NVIDIA’s CUDA Architecture GPUs . [Online] NVIDIA, 2009. http://www.nvidia.com/object/directcompute.html. 11. —. CUDA Documentation. CUDA http://www.nvidia.com/object/cuda_develop.html.
-
Zone.
[Online]
2009.
12. —. Best Practises Guide 2.3. CUDA ZONE - CUDA Documentation. [Online] 2009. július. http://www.nvidia.com/object/cuda_develop.html. 13. Nikos Drakos, Ross Moore. Edge detection with image pyramid. [Online] 2009. http://fourier.eng.hmc.edu/e161/lectures/canny/node2.html. 14. Márta, Máriáss. Vastagbélrák. [Online] http://www.hazipatika.com/topics/daganatos/denc/Vastagbelrak?did=295.
2009.
15. Intel. Intel Core i7 Processor Overview. Intel Products. [Online] 2009. http://www.intel.com/products/processor/corei7/index.htm. 16. György, Fehér. A háziállatok funkcionális anatómiája. Budapest : Mezőgazda Kiadó, 2004. Oldal 78
17. Group, Khronos. OpenCL. [Online] 2009. http://www.khronos.org/opencl/. 18. Gábor, Valcz. Az ép és beteg vastagbél szöveti felépítése. biotech.nik.bmf.hu. [Online] 2008. október. [Hivatkozva: 2009. március 10.] http://biotech.nik.bmf.hu/web/uploads/slides/DKM/2007_11_05_Az_ep_es_a_beteg_vasta gbel_szovettani_szerkezete.pdf. 19. AMD. ATI Graphics for Desktop PCs. [Online] AMD, http://www.amd.com/us/products/desktop/graphics/Pages/desktop-graphics.aspx.
2009.
20. —. AMD Professional Graphics. [Online] AMD, 2009. http://www.amd.com/us/products/workstation/graphics/Pages/workstation-graphics.aspx. 21. Levente, Ficsór. Colonbiopsziás hisztológiai metszetek számítógéppel automatizált kiértékelése. Szegedi Tudományegyetem : ismeretlen szerző, 2003. 22. Mundhenk, T. Nathan. HSV And H2SV Color Space. [Online] University of Southern California, 2009. http://ilab.usc.edu/wiki/index.php/HSV_And_H2SV_Color_Space. 23. Gonzalez, Rafael C. és Woods, Richard E. Digital Image Processing, 2nd ed. New Jersey : Prentice Hall Press, 2002, old.: 295. 24. —. Digital Image Processing, 2nd ed. New Jersey : Prentice Hall Press, 2002, old.: 525. 25. —. Digital Image Processing, 2nd ed. New Jersey : Prentice Hall Press, 2002. old.: 661. 26. Apple. vImage Programming Guide. Mac Dev Center. [Online] Apple, 2009. http://developer.apple.com/mac/library/documentation/Performance/Conceptual/vImage/In troduction/Introduction.html.
Oldal 79
