Analiza danych spektralnej tomografii optycznej OCT z wykorzystaniem wielordzeniowych kart graficznych

Uniwersytet Mikołaja Kopernika w Toruniu Wydział Fizyki, Astronomii i Informatyki Stosowanej

Marcin Sylwestrzak

Analiza danych spektralnej tomografii optycznej OCT z wykorzystaniem wielordzeniowych kart graficznych

Praca doktorska napisana pod kierunkiem prof. dr. hab. Piotra Targowskiego

Toruń 2013

Pragnę podziękować prof. Piotrowi Targowskiemu, mojemu promotorowi, za inspirację i owocną współpracę w czasie studiów doktoranckich oraz za poświęcony mi czas na cenne dyskusje podczas realizacji projektu. Dziękuję również za wszelkie wysiłki organizacyjne, zapewnienie bardzo dobrych warunków pracy i prowadzenia badań oraz możliwość uczestniczenia i prezentacji pośrednich wyników prac na wielu konferencjach krajowych i zagranicznych.

Dodatkowo dziękuję całemu Zespołowi Fizyki Medycznej za bardzo owocną współpracę.

Osobne podziękowania kieruję do moich rodziców i

bliskich za cierpliwość i pełne wsparcie, na które zawsze

mogłem liczyć.

Spis treści

Wstęp ... 9

1. Spektralna tomografia optyczna ... 13

1.1. Podstawy fizyczne spektralnej tomografii optycznej ... 14

1.1.1. Ramię referencyjne ... 16

1.1.2. Ramię obiektowe ... 16

1.1.3. Rejestrowany sygnał ... 17

1.2. Obrazowanie strukturalne ... 18

1.2.1. Tworzenie A-skanów strukturalnych ... 18

1.2.2. Osiowa zdolnośd rozdzielcza metody ... 20

1.2.3. Zakres obrazowania ... 20

1.2.4. Opis procedur stosowanych w przetwarzaniu danych SOCT ... 22

1.2.4.1. Odjęcie tła ... 24

1.2.4.2. Transformacja



– k... 26

1.2.4.3. Numeryczna kompensacja dyspersji ... 29

1.2.4.4. Kształtowanie widma ... 29

1.2.4.5. Analiza częstotliwościowa widma interferencyjnego ... 32

1.2.4.6. Tworzenie optycznego A-skanu ... 33

1.3. Obrazowanie przepływów metodą dopplerowską ... 34

1.3.1. Faza w sygnale interferencyjnym ... 34

1.3.2. Metoda dopplerowska... 35

1.3.3. Analiza danych fazowych ... 37

1.4. Podsumowanie ... 38

2. Budowa tomografu i architektura oprogramowania ... 39

2.1. Wykorzystywany sprzęt i oprogramowanie ... 39

2.1.1. Charakterystyka podzespołów tomografu i zależności pomiędzy nimi ... 39

2.1.2. Wykorzystywane oprogramowanie przez rozwijaną aplikację ... 41

2.2. Opis projektów tworzących aplikację ... 42

2.3. Główne wątki aplikacji ... 44

Spis treści

2.3.1. Wątek akwizycyjny ... 44

2.3.2. Wątek obliczeniowo-wizualizacyjny ... 45

2.4. Protokoły akwizycji danych ... 46

2.5. Podsumowanie ... 50

3. Równoległe przetwarzanie danych na procesorach graficznych... 51

3.1. Rys historyczny ... 52

3.1.1. Skąd się wzięła wielordzeniowośd?... 52

3.1.2. Historia obliczeo na kartach graficznych ... 54

3.1.3. Zastosowania obliczeo na kartach graficznych ... 57

3.1.4. Procesory graficzne NVIDIA ... 59

3.2. Technologia CUDA ... 61

3.2.1. Podstawowe pojęcia ... 61

3.2.2. Architektura Fermi™ ... 63

3.2.3. Programowanie wielowątkowe w CUDA C ... 64

3.2.4. Model pamięci ... 67

3.2.4.1. Rejestry ... 68

3.2.4.2. Pamięd lokalna ... 68

3.2.4.3. Pamięd wspólna ... 68

3.2.4.4. Pamięd globalna ... 69

3.2.4.5. Pamięd stała ... 70

3.2.4.6. Pamięd tekstur ... 70

3.2.5. Narzędzia do optymalizacji kodu ... 70

3.3. Podsumowanie ... 71

4. Implementacja algorytmów przetwarzania danych OCT na GPU ... 73

4.1. Przegląd literatury dotyczącej wykorzystania procesorów GPU w OCT ... 73

4.1.1. Pierwsze zastosowania technologii CUDA w OCT ... 73

4.2. Przetwarzanie równoległe w procesorze graficznym ... 76

4.2.1. Testowanie rzeczywistej wydajności procesora graficznego ... 77

4.2.2. Implementacja kodu na procesor graficzny ... 80

4.2.3. Implementacja i optymalizacja algorytmów SOCT na procesor graficzny ... 82

4.2.3.1. Implementacja i optymalizacja analiz poprzedzających przekształcenie Fouriera .. 84

4.2.3.2. Wywołanie transformacji Fouriera ... 95

4.2.3.3. Implementacja procedur finalizujących obliczenia ... 96

4.2.4. Analiza wydajności zaimplementowanych algorytmów SOCT ... 104

4.3. Najważniejsze etapy optymalizacji ... 109

4.4. Podsumowanie ... 111

5. Wykorzystanie uzyskanej szybkości przetwarzania danych ... 113

5.1. Trójwymiarowe obrazowanie OCT w czasie rzeczywistym... 113

Spis treści

5.1.1. Podstawy wizualizacji w OpenGL

^

... 113

5.1.2. Wydajnośd teksturowania ... 115

5.1.3. Przetwarzanie i wizualizacja tomogramów OCT w czasie rzeczywistym ... 119

5.2. Automatyczne dostosowywanie protokołów skanowania ... 121

5.2.1. Adaptacja protokołów skanowania w metodzie dopplerowskiej ... 123

5.2.2. Adaptacja protokołów skanowania w metodzie redukcji szumu plamkowego ... 128

5.2.3. Adaptacja protokołów w badaniach oftalmologicznych ... 131

5.3. Podsumowanie ... 132

6. Podsumowanie rozprawy ... 133

Dodatki ... 137

A. Wyznaczanie krzywej transformacji  – k ... 137

B. Wyznaczanie wektora kompensacji dyspersji ... 139

C. Uzupełniający materiał ... 143

D. Lista publikacji autora ... 145

Literatura ... 147

Spis treści

Wstęp

Tomografia optyczna z użyciem światła częściowo spójnego (ang. Optical Coherence Tomgraphy) to nieinwazyjna technika pozwalająca na obrazowanie struktury wewnętrznej obiektów częściowo przeźroczystych dla światła. W wyniku pomiaru w sposób bezkontaktowy otrzymywane są dwuwymiarowe przekroje w głąb obiektu wzdłuż linii skanowania wiązką analizującą. Tak uzyskane tomogramy obrazują wszystkie nieciągłości struktury, które wiązka napotkała na swojej drodze. Dzięki nieinwazyjności metody tomografia optyczna jest obecnie powszechnie wykorzystywana w diagnostyce medycznej jako narzędzie do obrazowania struktur wewnętrznych tkanek i organów – przede wszystkim w okulistyce [1, 2], ale także w dermatologii [3, 4], gastrologii [5, 6] i innych [7, 8]. Poza zastosowaniami medycznymi tomografia optyczna wykorzystywana jest w badaniach materiałowych [9], w szczególności obiektów zabytkowych, takich jak obrazy olejne [10, 11]

czy szkła witrażowe [12], w których ze względu na wartość zabytkową badanych obiektów metody inwazyjne są niedopuszczalne. Badania te pozwalają właściwie ocenić stan zachowania zabytków oraz zaplanować, a także monitorować w czasie rzeczywistym niektóre zabiegi konserwatorskie [13, 14].

Wieloletni rozwój tomografii OCT wykazał, że zastosowanie dodatkowych analiz

numerycznych, specyficzne tryby rejestracji danych bądź wprowadzenie niewielkich

udoskonaleń do urządzeń pozwalają znacznie rozwinąć metodę i uzyskać dodatkowe

informacje. Jednym z takich rozszerzeń jest metoda dopplerowska (ang. Doppler OCT)

pozwalająca na wyznaczenie prędkości przepływu krwi w naczyniach krwionośnych oka [15,

16]. Innym rozwinięciem jest metoda wymagająca rozbudowy urządzenia w celu niezależnej

rejestracji dwóch składowych polaryzacyjnych światła. Dzięki naturalnej dwójłomności

niektórych tkanek tworzących oko, PSOCT (ang. Polarization Sensitive OCT) pomaga we

wczesnym wykrywaniu jaskry [17]. W przypadku pomiarów materiałowych metoda ta

pozwala analizować naprężenia wewnątrz badanych obiektów [18].

Wstęp

Tomografia OCT jest ciągle rozwijana przez wiele ośrodków badawczych na całym świecie. Obok głównego trendu dążącego do poprawy jakości obrazowania strukturalnego, znaczna część prac skierowana jest na nowe zastosowania tej techniki. Od kilku lat opracowywane są metody mające na celu uzyskanie informacji o funkcjonowaniu poszczególnych tkanek. Po osiągnięciu pewnego standardu w obrazowaniu prędkości przepływu krwi w naczyniach krwionośnych ludzkiego oka, wiele grup badawczych rozpoczęło prace służące obrazowaniu struktury naczyń w mózgu [19, 20] oraz w kierunku badania reakcji tkanki na pobudzanie światłem bądź impulsami napięciowymi [21]. Badania te w chwili obecnej prowadzone są na zwierzętach.

Jednym z ograniczeń tomografii OCT jest czasochłonna analiza numeryczna danych, która wraz z wizualizacją przekrojów na ekranie zwykle zajmuje więcej czasu niż akwizycja.

Obrazowanie w czasie rzeczywistym jest niezbędną funkcją oprogramowania sterującego tomografem. Pozwala ono na właściwe ustawienie pacjenta i przygotowanie go do pomiaru, a w przypadku badań materiałowych – odpowiednie umieszczenie badanego obiektu.

Realizacja tego zadania wymaga przede wszystkim szybkiego obrazowania badanego obszaru. W chwili obecnej na dobrze zoptymalizowanym oprogramowaniu pracującym wielowątkowo na procesorze głównym (CPU) uzyskuje się wydajność obrazowania wynoszącą ok. 10 klatek na sekundę dla dwóch przekrojów w rozdzielczości 500 × 1024 piksele. Film zarejestrowany podczas obrazowania w takim trybie siatkówki pacjenta dołączony został do pracy na nośniku CD (szczegółowy opis zawartości płyty znajduje się w dodatku 0). Jakość takiego podglądu ze względu na liczbę zarejestrowanych i przetwarzanych widm nie jest zbyt wysoka i pozwala jedynie na przybliżone ustawienie pacjenta przed wykonaniem właściwego pomiaru. W tej sytuacji uzyskanie obrazu drobnych zmian, istotnych dla diagnozowania, wymaga wykonania wielu pomiarów i jest w pewnym stopniu kwestią przypadku.

Od kilku lat bardzo dynamicznie rozwija się technologia przetwarzania danych na procesorach graficznych. Już w 2010 roku NVIDIA deklarowała wydajność wykonywania obliczeń zmiennoprzecinkowych wyższą o rząd wielkości niż osiągalna przez najszybsze procesory główne [22]. Zastosowanie równoległego przetwarzania danych na wielordzenio- wych procesorach graficznych otwiera wiele nowych możliwości w technice spektralnej tomografii optycznej, przede wszystkim zaś usuwa główne ograniczenie blokujące obrazowanie w czasie rzeczywistym – czasochłonną analizę danych.

Zamiarem podjętego projektu reformułowania algorytmów obliczeniowych spektralnej

tomografii OCT w celu równoległego przetwarzania danych na wielordzeniowych

procesorach graficznych było:

Wstęp

1. umożliwienie wysokojakościowego obrazowania w czasie rzeczywistym, tak aby uzyskiwane przekroje wyświetlane były na ekranie jak najszybciej, bez zbędnych opóźnień, które bardzo utrudniają ustawienie pacjenta bądź obiektu w badaniach materiałowych,

2. obrazowanie bez żadnych ustępstw dotyczących jakości (w celu przyśpieszenia obliczeń) prezentowanego na ekranie przekroju – prezentowany w czasie rzeczywistym tomogram ma być wyświetlony w najlepszej możliwej do osiągnięcia dla tej metody jakości,

3. trójwymiarowe obrazowanie w czasie rzeczywistym, dotychczas (w przypadku obliczeń na procesorze głównym) w ogóle nie osiągalne,

4. implementacja dodatkowych trybów obrazowania rozszerzających standardowe obrazowanie strukturalne w OCT o metodę dopplerowską (pozwalającą na obrazowanie przepływów naczyniach krwionośnych) oraz metodę redukcji szumu plamkowego,

5. implementacja algorytmów automatycznego dostosowywania protokołów skanowa- nia w metodzie:

 dopplerowskiej – w celu samoczynnego dostosowania protokołu do mierzonego przepływu, tak aby wykorzystywał on cały dostępny zakres dynamiczny mierzonych prędkości,

 redukcji szumu plamkowego – w celu rejestracji danych w taki sposób, aby uśrednianie pozwalało na zminimalizowanie szumu plamkowego, jednocześnie nie pogarszając zbytnio rozdzielczości poprzecznej obrazowania.

Wyżej przedstawione cele pozwoliły na sformułowanie następującej tezy niniejszej rozprawy: optymalne wykorzystanie typów pamięci i reformułowanie algorytmów w celu równoległego przetwarzania danych spektralnej tomografii optycznej na kartach graficznych pozwala na przyśpieszenie obliczeń, umożliwiające trójwymiarowe obrazowanie strukturalne oraz czynnościowe w czasie rzeczywistym z jednoczesnym adaptacyjnym modyfikowaniem protokołów skanowania.

Niniejszą tezę autor udowodnił, dostosowując algorytmy przetwarzania danych w

spektralnej tomografii optycznej do nowej obliczeniowej platformy sprzętowej, którą

stanowią procesory graficzne. Optymalna implementacja algorytmów na procesory graficzne

wymaga dokładnej znajomości architektury sprzętowej i umiejętności jej wykorzystania,

jednak sama wiedza jest tu niewystarczająca. Bardzo istotną kwestią jest doświadczenie, które

nabyć można jedynie poprzez wiele godzin testowania i analizowania wydajności

różnorodnych rozwiązań. Dostęp do sześciu odmiennych typów pamięci, możliwość

dowolnego podziału wykonywanych obliczeń pomiędzy miliony uruchamianych wątków oraz

Wstęp

ogólny podział zadań pomiędzy bloki wątków stanowi o niemal nieograniczonych możliwościach implementacyjnych. Z tego względu znalezienie wydajnych rozwiązań wymaga zarówno czasu, jak i cierpliwości. W dalszej części rozprawy przedstawiony zostanie przykład obrazujący, jak duży wpływ na wydajność obliczeń ma kolejność wykonywania instrukcji (wtedy gdy z punktu widzenia poprawności algorytmu jest ona dowolna).

W pierwszym rozdziale przedstawiona została ogólnie metoda spektralnej tomografii optycznej wraz z jej opisem matematycznym i fizycznym. Opisane zostały również metody analizy danych oraz ich wpływ na jakość otrzymywanego tomogramu. W drugim rozdziale została omówiona budowa tomografu, opisane podstawowe elementy i podzespoły wchodzące w jego skład oraz połączenia i zależności pomiędzy nimi. Rozdział ten zawiera również zwięzły opis całego oprogramowania przygotowanego przez autora w celu realizacji tego projektu w języku C/C++ oraz C#. W rozdziale trzecim znajduje się opis wykorzysty- wanej architektury sprzętowej i technologii równoległego przetwarzania danych na procesorach graficznych. Czwarty rozdział, objętościowo największy i jednocześnie najistotniejszy, przedstawia opis zastosowanych rozwiązań implementacyjnych oraz optymalizacyjnych wprowadzonych w rozwijanych algorytmach przetwarzania danych spektralnej tomografii optycznej. W piątym rozdziale przedstawiono rezultaty bazujące na szybkim przetwarzaniu danych za pomocą zaimplementowanych metod, a więc szybkie obrazowanie tomogramów dwuwymiarowych, wizualizacje trójwymiarowe uzyskane w czasie rzeczywistym oraz algorytmy automatycznej adaptacji protokołów skanowania.

1. Spektralna tomografia optyczna

Spektralna tomografia optyczna wykorzystuje interferencję promieniowania z zakresu podczerwieni, najczęściej z pierwszego (ok. 850 nm) bądź drugiego (ok. 1300 nm) okna telekomunikacyjnego. W tomografii OCT stosowane są źródła generujące promieniowanie o dużej szerokości spektralnej (nawet do 200 nm), jednocześnie cechujące się małą spójnością czasową, co zapewnia wysoką rozdzielczość osiowa (w głąb badanego obiektu) oraz dużą spójnością przestrzenną, gwarantującą wysoką czułość pomiaru. Wiązka światła o takich parametrach kierowana jest na obiekt, w który wnika i częściowo rozprasza się na niejednorodnościach jego struktury oraz odbija się od granic ośrodków o różnych współczynnikach załamania. Niewielka część promieniowania sondującego, która po odbiciu od obiektu zostaje skierowana przez soczewkę obrazującą do tomografu, zostaje zarejestrowana i poddana zostaje analizie numerycznej, dzięki której możliwe jest wyznaczenie wszystkich nieciągłości w obiekcie wzdłuż biegu promienia. W ten sposób otrzymywany jest rozkład współczynnika rozpraszania/odbicia światła wzdłuż kierunku wiązki skanującej, zwany – analogicznie jak w ultrasonografii – A-skanem (skanem amplitudowym). Po wyrażeniu amplitudy w skali poziomów szarości lub „fałszywych”

kolorów powstaje jedna pionowa linia wewnętrznego przekroju badanego obiektu. Dla zwięzłości, w niniejszej rozprawie określenie „A-skan” będzie stosowane zarówno do rozkładu amplitudy współczynnika rozpraszania, jak i do pojedynczej linii obrazu przekroju.

Podstawowym wynikiem pomiaru tomograficznego jest dwuwymiarowy widok przekroju badanego obiektu (B-skan – Rys. 1.1 ), składający się z wielu przyległych do siebie A-skanów.

Przedstawia on jego strukturę wewnętrzną, a uzyskiwany jest poprzez skanowanie obiektu wiązką.

Metoda tomografii optycznej dzieli się na dwie podstawowe odmiany: pierwsza, tzw.

czasowa (ang. Time domain OCT, TdOCT), została zaprezentowana już w 1991 roku [23], natomiast druga – fourierowska (ang. Fourier domain OCT, FdOCT) – pięć lat później [24].

Tą z kolei dzieli się na bardziej popularną spektralną tomografię optyczną (ang. Spectral

domain OCT, SdOCT, albo: SOCT) oraz tomografię optyczną z użyciem laserów strojonych

(ang. Swept Source OCT, SSOCT) [25]. Systemy pracujące w pierwszej odmianie (czasowej)

1 Spektralna tomografia optyczna

Rys. 1.1. Porównanie jakości obrazowania okolicy tarczy nerwu ludzkiego oka odmianą: a) czasową i b) spektralną tomografii optycznej OCT [28]

wymagają skanowania w głąb obiektu poprzez zastosowanie mechanicznie przesuwanego zwierciadła w ramieniu interferometru. Natężenie odbitego światła mierzone jest podczas ruchu zwierciadła za pomocą fotodiody. W ten sposób otrzymywane są kolejne A-skany.

Badanie tą metodą trwa stosunkowo długo, na tyle, że mimowolne ruchy pacjenta podczas badania oka wpływają negatywnie na jakość otrzymywanego tomogramu. W układach spektralnych zwierciadło referencyjne utrzymywane jest w stałym położeniu, a rejestracja światła zamiast na fotodiodzie odbywa się za pośrednictwem spektrometru. Z uzyskanego w ten sposób widma interferencyjnego generowany jest jeden A-skan przez zastosowanie transformacji Fouriera. Usunięcie ruchomych elementów z interferometru i rejestracja widm na liniowej kamerze o bardzo szybkim transferze danych ze spektrometru pozwoliły na około stukrotne przyśpieszenie akwizycji i zwiększenie osiowej rozdzielczości obrazowania z 10 µm do ok. 3 µm [26, 27]. Mimo istotnej niedogodności, jaką stanowi ograniczony zakres obrazowania (zwykle 2-3 mm), metoda ta niemal całkowicie wyparła odmianę czasową – TdOCT, oferując nieporównywalnie wyższą jakość i szybkość obrazowania (Rys. 1.1 [28]).

Druga z metod fourierowskich (SSOCT) wykorzystuje specjalny laser strojony generujący promieniowanie o długości fali zmieniającej się cyklicznie w czasie z bardzo dużą częstością. Dzięki temu możliwe jest zastąpienie rejestracji widma za pomocą kamery spektrometru przez próbkowanie w czasie na szybkiej fotodiodzie. Transformacja Fouriera tak zarejestrowanego sygnału również pozwala na uzyskanie A-skanów w sposób analogiczny jak w metodzie spektralnej. Podstawowymi przewagami tej odmiany tomografii nad SOCT są: większa szybkość obrazowania [29, 30], mniejsze spadki czułości przy wzroście głębokości obrazowania [31, 32] oraz większe zakresy obrazowania (do ok. 8 mm) w głąb obiektu [33].

1.1. Podstawy fizyczne spektralnej tomografii optycznej

W spektralnej tomografii optycznej informacja o strukturze wewnętrznej obiektu

uzyskiwana jest z natężenia sygnału interferencyjnego. Tomografy OCT projektowane są

1 Spektralna tomografia optyczna

najczęściej w konfiguracji interferometru Michelsona, jednak zdarzają się rozwiązania wykorzystujące interferometr Macha-Zehndera [34, 35]. Rysunek (Rys. 1.2) ilustruje zasadę działania tomografu SOCT w konfiguracji Michelsona. Procedura obrazowania przeprowadzana jest w następujący sposób: wiązka światła ze źródła kierowana jest na dzielnik wiązki, tam dzielona jest na dwie części, które propagują się w dwóch ramionach interferometru. W ramieniu referencyjnym wiązka odbija się od zwierciadła referencyjnego, natomiast badany obiekt umieszczany jest w ramieniu obiektowym. Światło odbija się częściowo od powierzchni obiektu, a częściowo wnika do jego wnętrza. Tam ulega rozproszeniu i częściowemu odbiciu. Część światła, które powraca do interferometru, niesie ze sobą informację o tym, z jakiej głębokości wewnątrz obiektu zostało ono odbite lub rozproszone. Obie wiązki powracają do dzielnika wiązki i interferując ze sobą tworzą widmo interferencyjne, które rejestrowane jest na spektrometrze. Tak zarejestrowany sygnał interferencyjny poddawany jest analizie numerycznej szczegółowo opisanej w dalszej części rozprawy. Umożliwia ona otrzymanie z każdego widma jednej pionowej linii tomogramu – A-skanu. Rejestracja wielu takich A-skanów obok siebie, realizowana poprzez skanowanie obiektu, pozwala na uzyskanie dwuwymiarowych przekrojów (B-skanów), a nawet na zarejestrowanie trójwymiarowego tomogramu składającego się z przyległych do siebie B-skanów.

Stosując opis formalny propagacji światła jako płaskiej fali elektromagnetycznej, natężenie pola elektrycznego zapisać możemy następującym równaniem:

)]

( exp[

) ( ) ,

( k t E

₀

k i t zk

E    , (1.1)

Rys. 1.2. Schemat budowy spektralnego tomografu OCT w konfiguracji interferometru Michelsona. Droga optyczna z

r

pokonywana jest przez światło pomiędzy dzielnikiem wiązki a zwierciadłem referencyjnym, z

w

jest nadwyżką drogi optycznej do poszczególnych warstw wewnątrz

obiektu ponad drogę z

r

.

1 Spektralna tomografia optyczna

gdzie E

0

(k) jest rozkładem spektralnym światła emitowanego ze źródła (stałym w czasie), k = 2π/λ jest liczbą falową w próżni,  jest częstością optyczną światła, a z i t są zmiennymi – odpowiednio przestrzenną i czasową. Dla fali rozchodzącej się w ośrodku o współczynniku załamania n, z jest drogą optyczną:



 n ds

z . (1.2)

1.1.1. Ramię referencyjne

Oznaczymy przez  amplitudowy współczynnik transmisji światła przez dzielnik wiązki na wprost, a przez R amplitudowy współczynnik odbicia światła od zwierciadła referencyjnego znajdującego się w odległości optycznej z

r

. Przy pierwszym przejściu (pod kątem prostym) przez dzielnik wiązki światło pomniejsza swoją amplitudę o 1   , a przy drugim (przechodząc w kierunku spektrometru na wprost) o  . Dla takich oznaczeń amplituda światła przechodzącego przez ramię referencyjne a następnie wychodzącego z dzielnika wiązki w kierunku spektrometru, opisana jest następującym równaniem:

)]

2 ( exp[

) 1 ( ) ( ) ,

( k t E

₀

k R i t z k

E

_r

     

_r

, (1.3)

0    1, zazwyczaj: 

= 0,5

0  R  1, zazwyczaj R ≈ 1.

Czynnik 2 w argumencie funkcji wykładniczej wynika z faktu, że światło pokonuje drogę pomiędzy dzielnikiem a zwierciadłem dwukrotnie (w obie strony). Równanie nie uwzględnia drogi przebytej przez światło wewnątrz dzielnika wiązki. Jest to świadome uproszczenie, wynikające z faktu, że na częstotliwość sygnału interferencyjnego ma wpływ jedynie różnica dróg optycznych, a można przyjąć, że wewnątrz dzielnika drogi te są identyczne dla obu wiązek.

1.1.2. Ramię obiektowe

Amplituda światła wchodzącego do ramienia obiektowego osłabiana jest na dzielniku wiązki o czynnik  w wyniku przejścia na wprost. Dalej światło ulega odbiciu od W warstw obiektu (których amplitudowe współczynniki odbicia oznaczymy przez r

_w

), a następnie, przechodząc przez dzielnik pod kątem prostym, osłabione zostaje o czynnik 1   . Ostatecznie amplituda promieniowania pola elektrycznego pochodzącego z ramienia obiektowego po przejściu przez dzielnik wiązki jest sumą amplitud fal świetlnych odbitych od poszczególnych warstw badanego obiektu:

)]

) (

2 ( exp[

) 1 ( ) ( ) , (

0

0

k r i t z z k

E t k

E

_{r w}

W

w w

ob

    







 _. ^(1.4)

1 Spektralna tomografia optyczna 1.1.3. Rejestrowany sygnał

W wyniku interferencji spójnego promieniowania pochodzącego z obu ramion interferometru wypadkowa amplituda pola elektrycznego fali świetlnej ma postać:

) , ( ) , ( ) ,

( k t E k t E k t

E 

_r



_ob

. (1.5)

Wielkością mierzoną na liniowej kamerze w spektrometrze jest natężenie światła, w związku z tym zarejestrowane widmo odpowiada widmowej gęstości mocy promieniowania:

. ) , ( ) , ( )

, ( ) , (

) , ( ) , ( )

, ( ) , (

) , ( ) , ( ) (

*

*

*

*



























^

t k E t k E t

k E t k E

t k E t k E t

k E t k E

t k E t k E k G

r ob

ob ob

ob r

r

r

(1.6)

gdzie  .  oznacza uśrednianie po czasie wynikające z faktu, że czas naświetlania widma na kamerze jest wielokrotnie dłuższy od okresu drgań fali świetlnej. W rezultacie rejestrowane widmo nie zależy od czasu.

Po podstawieniu wzorów (1.3) i (1.4) wyrażenie (1.6) można przedstawić w postaci:

 

 exp( 2 ( ) ) exp( 2 ( ) )  , )

1 ( ) (

) 2 exp(

) 2 exp(

) 1 ( ) (

) 1 ( ) (

) 1 ( ) ( ) (

' 0

' '

' 0

0 0

0 0

0







































W

w w w

w w w

w w

w

w w

W

w w W

w w

k z z i k

z z i r

r k

G

k z i k

z i r

R k

G

r k

G

R k

G k G

















(1.7)

gdzie G

₀

( k )  E

₀

( k ) E

₀^*

( k )  jest widmem (energetycznym) źródła. Pierwsze dwa składniki sumy opisują widmową gęstość mocy światła powracającego z ramienia referencyjnego i obiektowego. Kolejne człony to wyrazy opisujące modulację widma w wyniku interferencji.

Najważniejszym wyrazem z punktu widzenia obrazowania techniką OCT jest trzeci człon.

Zawiera on składowe powstałe w wyniku interferencji fali odbitej od zwierciadła referencyjnego z falami rozproszonymi lub odbitymi od kolejnych warstw wewnątrz obiektu.

Częstości tych składowych są równe (w dziedzinie liczb falowych) podwojonym różnicom

odpowiednich dróg optycznych 2z

w

. Tym samym niosą one informację o położeniach

kolejnych warstw obiektu względem zwierciadła odniesienia. Ostatni człon równania (1.7)

związany jest z wzajemną interferencją fal odbitych od różnych warstw wewnątrz obiektu,

jest to tzw. człon autokorelacyjny, tworzący artefakty na początku zakresu pomiarowego.

1 Spektralna tomografia optyczna

1.2. Obrazowanie strukturalne

Dotychczasowe rozważania, dotyczące podstaw fizycznych i zasady działania układów optycznych SOCT, zakończone zostały matematycznym opisem widmowej gęstości mocy światła rejestrowanej na kamerze spektrometru. W bieżącym podrozdziale przedstawione zostaną wszystkie operacje numeryczne niezbędne do uzyskania ostatecznego przekroju strukturalnego.

1.2.1. Tworzenie A-skanów strukturalnych

Z analizy postaci trzeciego składnika sumy w równaniu (1.7) wynika, że sygnał G(k) składa się z funkcji harmonicznych, których częstość jest proporcjonalna do różnicy dróg optycznych z

w

(Rys. 1.2, str. 15), czyli do różnicy odległości, które musi pokonać światło pomiędzy rozpraszającą warstwą w ramieniu obiektowym a zwierciadłem odniesienia w ramieniu referencyjnym. Z tego względu informacja o strukturze wewnętrznej obiektu uzyskiwana jest z zarejestrowanego na kamerze widma poprzez zastosowanie transformacji Fouriera:

 

 

 exp( 2 ( ) ) exp( 2 ( ) )  exp( ) . )

1 ( ) (

) exp(

) 2 exp(

) 2 exp(

) 1 ( ) (

) exp(

) 1 ( ) (

) exp(

) 1 ( ) (

) ( )

(

' 0

' '

' 0

0 0

0 0

0

dk izk k

z z i k

z z i r

r k

G

dk izk k

z i k

z i r

R k

G

dk izk r

k G

dk izk R

k G

k G FT z g

W

w w w

w w w

w w

w

w w

W

w w W

w w





































 

 

 







 





 





 























(1.8)

Uwzględnienie nieunormowanej funkcji spójności (z twierdzenia Wienera-Chinczyna) [36]:

 G k  G k izk dk FT

z 

^











 ( )

₀

( )

₀

( ) exp( ) , (1.9)

która zdefiniowana jest jako transformacja Fouriera widma źródła światła, pozwala na otrzymanie wzoru opisującego optyczny A-skan:

 

 

 

 







 

 













































 









 W

w w w

w w w

w w

w

w w

W

w w W

w w

z z z z

z z r r

z z z

z r R r z R z z

g

' 0

' '

'

0 0

)) (

2 ( )) (

2 (

) 2 ( ) 2 ( )

( ) ( ) 1 ( )

(   ^. (1.10)

Jest to funkcja zespolona, która jako transformata Fouriera funkcji rzeczywistej ma część

rzeczywistą symetryczną, a część urojoną antysymetryczną względem z = 0. Moduł

1 Spektralna tomografia optyczna powyższej funkcji zespolonej, po zlogarytmowaniu dla uwypuklenia słabych sygnałów, tworzy jeden optyczny A-skan strukturalny. Składając razem kolejne A-skany oraz kodując natężenia za pomocą skali fałszywych kolorów bądź szarości uzyskuje się B-skany. Każda składowa harmoniczna sygnału G(k) (1.7) w wyniku transformacji Fouriera przekształca się więc na impulsową funkcję spójności Γ(z) o stałej szerokości określonej przez postać funkcji G

0

(k). Funkcja ta określa więc osiową zdolność rozdzielczą tomografu (rozdział: 1.2.2.

Osiowa zdolność rozdzielcza metody).

Ponieważ wartości współczynników odbicia poszczególnych warstw obiektu są dużo mniejsze niż współczynnik odbicia zwierciadła r

w

<< R, tak więc dominującym składnikiem we wzorze (1.10) jest składnik pierwszy, który uzyskuje swoje maksimum dla z = 0 (Rys. 1.3).

Analogicznie jak we wzorze (1.7) jest on związany z natężeniem światła powracającego z ramienia referencyjnego. Drugi składnik, związany ze światłem odbitym od obiektu i powracającym przez ramię obiektowe interferometru, również osiąga swoje maksimum dla z = 0 i zazwyczaj wnosi mały wkład do amplitudy sygnału g(z). Konsekwencją postaci trzeciego składnika równania (1.7) jest tworzenie – w wyniku transformacji Fouriera – symetrycznego ciągu funkcji impulsowych znajdujących się po obu stronach osi odciętych, odpowiadających położeniom warstw wewnątrz badanego obiektu, w rezultacie tworzącego finalny przekrój przez strukturę obiektu. Czwarty wyraz odzwierciedla powstawanie wspomnianych już artefaktów autokorelacyjnych. Sygnały te powstają na początku zakresu pomiarowego, zwykle są znacznie słabsze niż sygnał tworzący tomogram, a często nawet niewidoczne na tle szumu.

) (z g



 



 





 W

w

rw

z R z

0

) ( ) ( ) 1 ( 





 



  







) 2 ( )

1 (

0

w W

w

w z z

r

 R



 



































 W

w w w

w w w

w w

w r z z z z z z

r

' 0

' '

' ( 2( )) ( 2( ))

) 1

( 



z

Rys. 1.3. Symulowane widmo amplitudowe po transformacji Fouriera widma interferencyjnego. Wykres obrazuje wpływ poszczególnych składników wzoru (1.10) na tworzony optyczny A-skan



 



  







) 2 ( )

1 (

0

w W

w

w z z

r

 R



1 Spektralna tomografia optyczna

1.2.2. Osiowa zdolność rozdzielcza metody

Zdolność rozdzielcza w obrazowaniu optycznym definiowana jest jako najmniejsza (geometryczna) odległość pomiędzy centrami rozpraszania możliwymi do zobrazowania jako odrębne elementy. Na podstawie wzorów (1.9) i (1.10) można sformułować wniosek, że osiowa zdolność rozdzielcza metody zależna jest od widma generowanego przez źródło światła. Tak rzeczywiście jest – szczegółowa analiza tego tematu była przedmiotem rozważań zawartych w rozprawie doktorskiej Iwony Gorczyńskiej [37]. Rozdzielczość osiowa w spektralnej tomografii optycznej najczęściej wyznaczana jest ze wzoru:

n

FWHM

z 



 

 

2

2

0

ln 2

1 , (1.11)

gdzie n – współczynnik załamania wewnątrz badanego obiektu, 

₀

– centralna długość fali i

∆

_FWHM

– szerokość spektralna w połowie maksimum widma mocy generowanego światła (FWHM, ang. Full Width at Half Maximum). Powyższy wzór jest słuszny dla źródeł światła o widmie gaussowskim. Warto zauważyć, że na mocy wzoru (1.11) osiowa zdolność rozdzielcza w ośrodku jest n krotnie lepsza niż w powietrzu.

Obecnie w OCT stosuje się źródła szerokopasmowe (np. wykorzystujące układ sprzężonych diod superluminescencyjnych) często o widmie w istotny sposób odbiegającym od kształtu funkcji Gaussa. W takich przypadkach bardziej poprawny jest wzór praktyczny [38]:

s s

s s

z n

min, max,

max,

2

min,

ln 4 1









 

  , (1.12)

gdzie 

min,s

i 

max,s

ograniczają zakres emisji źródła. W praktyce, m.in. ze względu na tłumienie światła w światłowodach, nawet dla źródeł światła o szerokości spektralnej zbliżającej się do 200 nm trudno jest uzyskać rozdzielczości poniżej 3 µm w powietrzu.

1.2.3. Zakres obrazowania

Zakres obrazowania w tomografii OCT to maksymalna osiowa odległość pomiędzy dwoma centrami rozpraszającymi, które mogą być zobrazowane na jednym A-skanie.

W SOCT odległość ta związana jest z maksymalnym opóźnieniem wiązki obiektowej względem referencyjnej, wynikającym z różnicy dróg optycznych w obu ramionach interferometru, które jeszcze może zostać zarejestrowane na spektrometrze jako najwyższa składowa harmoniczna.

Trzeci człon równania (1.7), opisujący tworzenie prążków interferencyjnych, można

zapisać za pomocą funkcji trygonometrycznych w następującej postaci:

1 Spektralna tomografia optyczna

...

) 2 cos(

2 )

1 ( ) ( ...

..

. ) (

0

0

 





 



k z r

R k

G k

G

_w

W

w



w

 . (1.13)

Zgodnie z wcześniejszymi oznaczeniami z

w

jest nadwyżką drogi optycznej do poszczególnych warstw wewnątrz obiektu. Przejście do odległości geometrycznych (∆s – geometryczna różnica dróg optycznych) wymaga uwzględnienia współczynnika załamania badanego ośrodka n. Tak więc jedna składowa oscylacji sygnału interferencyjnego może być opisana przez funkcję f(k):

) 2

cos(

) ( )

( k f

₀

k n k s

f   . (1.14)

Wraz ze wzrostem różnicy dróg w ramionach interferometru ∆s rośnie częstotliwość prążków interferencyjnych. Zgodnie z twierdzeniem Nyquista o próbkowaniu, maksymalna częstotliwość rejestrowanego przez spektrometr sygnału musi być taka, aby na jeden okres przypadały co najmniej dwie próbki ( Rys. 1.4 ), a więc dla maksymalnej różnicy dróg ∆s

max

, odległość∆k pomiędzy dwiema próbkami w dziedzinie liczby falowej k musi spełniać poniższy warunek:







_max

2 n k s _. ^(1.15)

Zakładając, że widmo jest rozłożone liniowo na linijce pikseli w dziedzinie długości fali  można przyjąć, że różnica pomiędzy długością fali padającej na sąsiednie piksele różni się o:

M

min

max



   ^

 , (1.16)

gdzie 

max

i 

min

są długościami fali padającymi na skrajne piksele kamery spektrometru.

Pomiędzy liczbą falową k a długością fali  zachodzi zależność odwrotna: k = 2/, tak więc najkrótsza fala 

min

padająca na linijkę pikseli opisywana jest przez liczbę falową k

max

. Odległości pomiędzy próbkami sygnału zarejestrowanymi na kamerze po przepróbkowaniu do dziedziny długości fali wynoszą:

 

 _



 2

₂

k (1.17)

) (k f

k

 k

Rys. 1.4. Warunek Nyquista próbkowania sygnału interferencyjnego dla

najwyższej częstotliwości oscylacji, odpowiadającej maksymalnej

różnicy dróg ∆z

max

1 Spektralna tomografia optyczna

Jak już zostało to wspomniane wcześniej, ∆ jest wartością stałą (1.16), dlatego odległość pomiędzy próbkami ∆k zależy jedynie od długości fali. Aby poprawnie przeprowadzić analizę numeryczną danych, kryterium Nyquista musi być spełnione dla wszystkich λ, w tym w najmniej korzystnym przypadku, a więc dla minimalnej długości fali 

min

. Wstawiając równanie (1.17) do (1.15) oraz uwzględniając (1.16), otrzymujemy:

min max 2 min

4  



 

 M

s

_MAX

n (1.18)

W powyższym wzorze dostrzec można pewną analogię do wzoru (1.11) definiującego rozdzielczość metody: zarówno zakres pomiarowy, jak i rozdzielczość są wprost proporcjonalne do kwadratu długości fali oraz odwrotnie proporcjonalne do szerokości widma (1.11) bądź zakresu spektralnego spektrometru (1.18). Oba te parametry są zwykle proporcjonalne do siebie, ponieważ spektrometry są tak projektowane, aby jak najefektywniej obrazować widmo na kamerze. Przyjmując to dość ogólne stwierdzenie, można sformułować ważny wniosek:

M z

z  



_max

 (1.19)

– zwiększanie rozdzielczości w SOCT przy zachowaniu tej samej liczby pikseli w detektorze prowadzi do zmniejszania zakresu pomiarowego.

1.2.4. Opis procedur stosowanych w przetwarzaniu danych SOCT

W przypadku detekcji fourierowskiej, z jaką mamy do czynienia w spektralnej tomografii optycznej, sygnał interferencyjny rejestrowany jest na liniowej kamerze w spektrometrze. Każde zebrane w ten sposób widmo wymaga wykonania szeregu procedur numerycznych, aby uzyskać z każdego z nich jeden A-skan. W wyniku rozwijania metody od ponad dziesięciu lat zarówno przez nasz zespół, jak i inne grupy badawcze na świecie, ustalił się standardowy zestaw procedur numerycznych wymaganych do uzyskania dobrej jakości obrazu tomograficznego. Składa się on z sześciu podstawowych kroków:

1. odjęcie tła,

2. transformacja – ,

3. numeryczna kompensacja dyspersji, 4. kształtowanie widma,

5. transformacja Fouriera,

6. prezentacja w skali logarytmicznej.

Wykonanie wszystkich tych kroków jest niezbędne w celu uzyskania tomogramów w

najlepszej jakości. Poniżej przedstawiony zostanie wpływ każdej z wyżej wymienionych

procedur na tworzony przekrój. Całkowicie nieczytelny tomogram uzyskany tylko poprzez

1 Spektralna tomografia optyczna

wykonanie transformacji Fouriera na zarejestrowanych widmach i odpowiednie wyświetlenie tych danych po zlogarytmowaniu sygnałów przedstawia rysunek (Rys. 1.5 a). W kolejnych krokach przedstawiony jest wpływ procedury odjęcia tła (Rys. 1.5 b) oraz przepróbkowania

Rys. 1.5. Wpływ poszczególnych procedur numerycznych na jakość tomogramu.

Przedstawiony przekrój pochodzi z badań obrazu olejnego odjęcie tła

1

kompensacja dyspersji 3

interpolacja

 – k 2

kształtowanie widma 4

transformacja Fouriera 5

20 log|FFT|

6

5, 6

1, 5, 6

1, 2, 5, 6

1, 2, 3, 5, 6

1, 2, 3, 4, 5, 6 Pełny zestaw procedur

numerycznych:

A-skan (a)

(b)

(d)

(e) (c) Zarejestrowane

widmo: Zastosowane Uzyskane tomogramy:

procedury

numeryczne:

1 Spektralna tomografia optyczna

widma z dziedziny długości fali  do dziedziny liczby falowej k (tzw. transformacja  – k, Rys. 1.5 c). Numeryczna kompensacja dyspersji (Rys. 1.5 d) pozwala na wyostrzenie obrazu w przypadku, gdy nie ma technicznej możliwości skompensowania jej sprzętowo (poprzez umieszczenie w obu ramionach interferometru elementów dyspersyjnych o jednakowych długościach). Pomaga ona również wyostrzyć tomogramy w przypadku pomiarów obiektów o dużej dyspersji bądź w przypadku wykonywania badań przez warstwę szkła. Numeryczne kształtowanie widma kompletuje zestaw procedur (Rys. 1.5 e), pozwalając na poprawę ostrości granic obrazowanych struktur poprzez minimalizację listków bocznych powstających w optycznym A-skanie w wyniku transformacji Fouriera sygnału o niegaussowskiej obwiedni. Operacja ta w niewielkim stopniu może pogorszyć rozdzielczość układu.

Dotychczasowy opis podstaw fizycznych tomografii optycznej przedstawiony był za pomocą funkcji ciągłych w dziedzinie czasu, długości fali bądź odległości. W praktyce rejestrowane widma są sygnałami dyskretnymi. Opis procedur numerycznych w dalszej części rozprawy zostanie przedstawiony z uwzględnieniem tego faktu; przyjęte też zostaną następujące oznaczenia:

m – numer próbki sygnału (w SOCT odpowiada numerowi piksela kamery), M – liczba próbek sygnału (w SOCT odpowiada liczbie pikseli kamery), g(m) – widmowa gęstość mocy zarejestrowana na kamerze spektrometru.

1.2.4.1. Odjęcie tła

Użyteczna informacja znajduje się jedynie w składowych zmiennych sygnału rejestrowanego na kamerze spektrometru, dlatego w przetwarzaniu danych SOCT w pierwszym kroku usuwana jest obwiednia widma źródła światła, a więc dwa pierwsze składniki równania (1.10). Ponadto w kroku tym zostają usunięte szumy koherentne pochodzące z optyki tomografu. Wszystkie elementy optyczne, znajdujące się wzdłuż biegu promienia w obu ramionach interferometru (przede wszystkim soczewki i kompensatory dyspersji), pomimo że posiadają warstwy antyrefleksyjne, odbijają znikome ilości światła z powrotem w kierunku interferometru. Zjawisko to dodatkowo wzmacniają drobinki zanieczyszczeń (np. kurzu), znajdujące się zarówno na soczewkach, jak i na skanerach.

Dopóki odległość pomiędzy takimi elementami jest większa niż zakres pomiarowy

(ok. 2 mm), zjawisko to nie wpływa na generowane obrazy. Jeżeli jednak dwie powierzchnie

odbijające umieszczone w jednym ramieniu interferometru znajdą się w odległości mniejszej

niż zakres pomiarowy bądź jeśli pomiędzy dwiema takimi powierzchniami znajdującymi się

w różnych ramionach interferometru różnica dróg optycznych będzie mniejsza niż zakres

pomiarowy, to światło odbite od nich będzie interferowało, a na widmach pojawi się

niewielka modulacja i w konsekwencji na obrazach tomograficznych zakłócenia te widoczne

będą w postaci poziomych linii. Przyczyna powstawania tego typu szumów jest identyczna

jak artefaktów związanych z wzajemną interferencją światła pochodzącego od wewnętrznych

1 Spektralna tomografia optyczna warstw badanych obiektów (ich opis matematyczny jest analogiczny do czwartego składnika równania (1.10)).

W celu usunięcia składowej stałej oraz opisanych powyżej szumów koherentnych pochodzących z optyki tomografu przeprowadzana jest następująca procedura: wiązka analizująca kierowana jest poza soczewkę obrazującą za pomocą skanerów galwanometrycznych, a na kamerze spektrometru rejestrowany jest sygnał (tzw. sygnał tła).

Zawiera on superpozycję światła powracającego z obu ramion spektrometru (oczywiście jest to przede wszystkim światło odbite od zwierciadła referencyjnego z możliwą niewielką modulacją wynikającą z refleksów od elementów optycznych w obu ramionach). Różnice pomiędzy sygnałem interferencyjnym zarejestrowanym podczas wykonywania pomiaru (stosunkowo duża modulacja) oraz sygnałem tła (znikoma modulacja) przedstawione zostały na wykresach (Rys. 1.6). Odjęcie sygnału tła od każdego zarejestrowanego widma pozwala na uzyskanie widm ze zminimalizowanymi zarówno szumami koherentnymi, jak i składową związaną z kształtem widma źródła światła. Sygnał taki, uzyskany w wyniku odjęcia sygnału tła (Rys. 1.6 b) od sygnału interferencyjnego z pomiaru (Rys. 1.6 a) przedstawiony został na wykresie (Rys. 1.8 a).

Rys. 1.6. Przykładowe widmo zarejestrowane na kamerze spektrometru podczas

pomiaru oraz sygnał zarejestrowany przy skierowaniu skanerów poza

obiekt, tzw. sygnał tła

1 Spektralna tomografia optyczna

W praktyce zauważalnie lepsze efekty (w postaci słabszych zakłóceń na przekrojach) są uzyskiwane, gdy w procedurze wykorzystuje się uśrednione widma tła zebrane podczas T akwizycji (np. T =100). Operację odjęcia tła możemy więc zapisać następującym równaniem:

) ( ) ( )

(

₀

1

m g m g m

g  

_tla

, (1.20)

gdzie:

m = 0, 1,..., M – jest numerem próbki zarejestrowanego sygnału,







^T

i tla i

tla

m g m

g

0

) ( )

( – jest uśrednionym sygnałem tła, g

_i^tla

(m) – jest sygnałem tła zebranym podczas i-tej akwizycji,

g

₀

(m) – jest „surowym” widmem zarejestrowanym podczas pomiaru, g

₁

(m) – jest widmem po wykonaniu operacji odjęcia tła.

Otrzymane za pomocą powyżej przedstawionych obliczeń widmo g

1

(m) stanowi sygnał wejściowy do następnego kroku w cyklu procedur numerycznych stosowanych do przetwarzania danych OCT.

1.2.4.2. Transformacja  – k

Ze wzoru (1.10) wynika, że oscylacje sygnału interferencyjnego opisane są za pomocą funkcji harmonicznych w dziedzinie liczby falowej k. W związku z tym jeśli w ramieniu obiektowym w zakresie pomiarowym urządzenia zostałoby umieszczone zwierciadło w takiej odległości, aby różnica dróg optycznych była niezerowa, urządzenie rejestrujące widmo interferencyjne w dziedzinie liczby falowej k zarejestruje widmo, cechujące się stałą częstotliwością oscylacji w całym zakresie spektralnym. Jednakże, w typowych spektrometrach z siatką dyfrakcyjną, projektowanych do detekcji sygnału interferencyjnego w spektralnej tomografii optycznej, widmo na kamerze rejestrowane jest w dziedzinie długości fali  [39]. Pomiędzy  a liczbą falową k zachodzi zależność odwrotna k = 2 /, tak więc rejestrowany sygnał jest sygnałem świergotowym (ang. chirp signal).

Celem operacji transformacji  – k jest takie przepróbkowanie widma zarejestrowanego w dziedzinie długości fali , aby dla danej różnicy dróg optycznych w interferometrze otrzymać sygnał interferencyjny w dziedzinie liczby falowej k, którego częstotliwość oscylacji będzie stała w całym zakresie spektralnym. Uzyskuje się to za pomocą funkcji m’=h(q) pozwalającej na znalezienie nowego indeksu piksela kamery m’

(m’

⁺

), odpowiadającego liczbie falowej k

q

: M q

k k k

k

_q ^M





 



^

1

0 1

0

, (1.21)

gdzie:

q = 0, 1,..., M-1 – numer próbki sygnału równoodległej w dziedzinie liczby falowej k, k

0

– liczba falowa światła padającego na pierwszy piksel kamery (m=0),

k

M-1

– liczba falowa światła padającego na ostatni piksel kamery (m=M-1).

1 Spektralna tomografia optyczna

Metoda wyznaczania tej krzywej przedstawiona została w dodatku A. Wynikowe numery pikseli m’ kamery nie są wyrażone liczbami naturalnymi, dlatego wymagane jest zastosowanie odpowiedniej interpolacji sygnału. Procedura przepróbkowania sygnału w OCT powinna charakteryzować się dwiema cechami: wysoką dokładnością odzwierciedlenia rzeczywistej wartości sygnału w punkcie dla zadanej wartości k oraz dużą szybkością działania. Niestety, te warunki są sobie wzajemnie przeciwstawne. Najszybszą metodą jest oczywiście interpolacja liniowa, w której wartość próbki w każdym punkcie obliczana jest jako średnia ważona z dwóch sąsiednich próbek:

   ( )   ( ) _{ } ( )    _{ } ( ) 1   _{ } ( )  

)

(

_{1 1 1}

2

q g h q h q h q g h q g h q

g       , (1.22)

gdzie:

q = 0, 1,..., M-1 – jest indeksem sygnału zależnym liniowo od liczby falowej k, g

₁

(m) – natężenie sygnału na m-tym pikselu kamery,

g

2

(q) – jest widmem po wykonaniu operacji transformacji  – k, h(q) – jest funkcją transformacji  – k,



^.

– jest operatorem zaokrąglenia do najbliższej mniejszej wartości całkowitej.

W spektralnej tomografii optycznej obok przedstawionej tu metody interpolacji liniowej stosowane są jeszcze dwie metody (wykorzystujące funkcje sklejane, ang. spline albo transformację Fouriera wraz z „doklejaniem zer”, ang. zero padding). Zostały one szczegółowo opisane w rozprawie doktorskiej Macieja Szkulmowskiego [40]. Ich największą zaletą jest lepsze odwzorowanie składowych interferencyjnych o największych częstotliwościach (zbliżających się do częstotliwości Nyquista), co w konsekwencji pozwala

Rys. 1.7. Krzywa h transformacji  – k

1 Spektralna tomografia optyczna

uzyskać nieco większy kontrast w B-skanie pod koniec zakresu pomiarowego. W typowych zastosowaniach techniki SOCT (do obrazowania stosunkowo cienkich struktur) różnica ta jest jednak minimalna i zazwyczaj do pominięcia. Kosztem zastosowania tych metod jest znaczne zwiększenie złożoności obliczeniowej. W tej sytuacji metody nieliniowe nie były rozwijane w ramach niniejszego projektu.

Konsekwencją niewykonania operacji transformacji  – k jest przeciek energii na kilka sąsiednich składowych częstotliwościowych w sygnale uzyskiwanym po wykonaniu transformacji Foureiera, co powoduje poszerzenie funkcji spójności, a więc pogorszenie rozdzielczości i rozmycie obrazu (Rys. 1.5 b). Różnice w widmach przed wykonaniem procedury przepróbkowania i po niej widoczne są jako przesunięcia próbek sygnału pomiędzy wykresami a i b (Rys. 1.8).

Rys. 1.8. Widmo interferencyjne po przetworzeniu kolejnymi procedurami:

a) odjęcia tła, b) transformacji  – k, c) kształtowania widma

1 Spektralna tomografia optyczna 1.2.4.3. Numeryczna kompensacja dyspersji

Dyspersja to zależność prędkości rozchodzenia się promieniowania elektromagnetycznego od jego częstości. Dla obrazowania OCT istotna jest dyspersja prędkości grupowej (GVD, ang. Group Velocity Dispersion). Zjawisko to, o ile nie jest poprawnie skompensowane, będzie powodować świergot sygnału interferencyjnego pochodzącego z pojedynczej głębokości i w konsekwencji rozmycie tomogramu [41].

Elementy optyczne niezbędne do budowy tomografu OCT (soczewki, światłowody) wykonane są najczęściej z różnego typu szkieł, a więc ośrodków dyspersyjnych. Dopóki w obu ramionach interferometru znajdują się elementy wykonane z identycznych szkieł i o takich samych grubościach, dyspersja jest skompensowana i nie powoduje żadnego pogorszenia jakości obrazu. Metodę te nie zawsze daje się zastosować do kompensacji dyspersji w materiałach wchodzących w skład badanego obiektu przez które obrazowane są jego głębiej leżące warstwy ponieważ materiały te często mają nieznane własności dyspersyjne. W tej sytuacji sprzętowe korygowanie dyspersji poprzez dodanie w ramieniu referencyjnym dodatkowego, odpowiednio dobranego szkła zazwyczaj nie prowadzi całkowitej kompensacji dyspersji. W przypadku niewielkiego niedopasowania należy zastosować numeryczną kompensację dyspersji. Sprowadza się ona do przemnożenia widma przez odpowiednio przygotowany wektor korekcyjny:

) ( ) ( )

(

₂

3

q g q D q

g   , (1.23)

gdzie:

g

₂

(q) – jest widmem po wykonaniu transformacji  – k, g

3

(q) – jest widmem po wykonaniu korekcji dyspersji, D(q) – jest wektorem korekcyjnym.

Wektor korekcyjny składa się z liczb zespolonych o jednostkowym module i odpowiednio dobranej fazie, tak więc dotychczas rzeczywisty sygnał g

2

(q), począwszy od tej operacji, staje się sygnałem zespolonym g

3

(q). Metoda tworzenia wektora korekcyjnego została przedstawiona w dodatku 2 (rozdział B). Zależy on tylko od jednego parametru, którego wartość należy dobrać empirycznie na podstawie prezentowanego na ekranie tomogramu. Dzięki szybkiemu procesowi przetwarzania danych na procesorze graficznym dobór parametru numerycznej kompensacji dyspersji można wykonać bardzo szybko, ponieważ każda zmiana wartości parametru widoczna jest niemal natychmiast jako wyostrzenie bądź rozmycie obrazu.

1.2.4.4. Kształtowanie widma

Źródła światła wykorzystywane w spektralnej tomografii optycznej to

szerokopasmowe diody superluminescencyjne lub lasery femtosekundowe generujące