USG/PWI
Spis treści
Obrazowanie falą płaską
W obrazowaniu fala płaską (PWI - Plane Wave Imaging) wykres opóźnień nadawczych względem przetworników ma postać prostej - przetworniki odpalane są w równych odstępach czasu, aby uformować wiązkę o (w przybliżeniu) płaskim froncie falowym. W przeciwieństwie do rozważanej wcześniej metody klasycznej, wiązka taka jest nieogniskowana - stąd oczekujemy podobnej ostrości / rozdzielczości na całej głębokości pomiarowej. O froncie falowym w PWI możemy myśleć jak o froncie fali pochodzącej ze źródła punktowego (ognisko wirtualne) oddalonego istotnie daleko za aperturą.
O ile nie ograniczają nas możliwości techniczne sprzętu (np. mała moc obliczeniowa), w PWI staramy się nadawać i odbierać pełną aperturą. Punktem wyjścia jest sytuacja w której wszystkie przetworniki nadają jednocześnie. Obrót wykresu opóźnień nadawczych odpowiada przesunięciu ogniska wirtualnego.
Parametry potrzebne do dalszej pracy:
f0=5.5e6 # Częstotliwość nadawcza przetworników [Hz]
fs=50e6 # Częstotliwość próbkowania [Hz]
pitch = 0.00021 # Deklarowana odległość między środkami przetworników nadawczo-odbiorczych
na = 3 # Liczba nadań
theta=[-10,0,10] # kąty dla kolejnych nadań
NT=192 # Liczba przetworników w pełnej aperturze
Ntr=192 # Pełna subapertura nadawcza
Surowe dane RF
Do dyspozycji mamy ponownie dwa pliki z pomiarów na fantomie nitkowym (usg2_nitki.npy) i cystowym (usg2_cysty.npy). W każdym pliku znajduje się trójwymiarowa tablica [math]NT\times N \times na [/math]
gdzie [math]N[/math] odpowiada czasowi rejestracji (maksymalnej głębokości obrazowanej) danych z pojedynczego nadania.
Kluczowe dla dalszych obliczeń będzie wyznaczenie chwili zerowej, czyli momentu rozpoczęcia akwizycji danych. Na różnych urządzeniach przyjmowane są różne konwencje - czasami akwizycja danych rozpoczyna się w momencie, gdy rozpoczyna się nadawanie; czasami akwizycja rozpoczyna się jeszcze przed nadawaniem.
Zadanie: proszę podejrzeć surowe dane RF dla kilku kolejnych nadań i spróbować określić (w przybliżeniu) w jakim momencie rozpoczyna się akwizycja danych:
a) przed rozpoczęciem nadawania
b) w momencie rozpoczęcia nadawania
c) w chwili gdy nadała połowa przetworników
Rekonstrukcja obrazu z pojedynczego strzału
W przypadku klasycznej rekonstrukcji dokonywaliśmy najpierw przesunięcia w czasie próbek RF, a następnie przesunięty (syntezowany) obraz interpolowaliśmy do siatki odpowiadającej rzeczywistym proporcjom obrazowanej struktury. W przypadku PWI wygodniej będzie nam najpierw określić siatkę, a następnie dla kolejnych punktów na tej siatce liczyć wartość natężenia. Poniższy kod wygeneruje dwie tablice w których przechowywane są współrzędne poziome i pionowe w [m] kolejnych punktów siatki. Siatka jest generowana przy założeniu, że punkt (0,0) znajduje się w środku apertury.
maxWidth=0.02 # Połowa szerokości obrazowania [m]
l2mbd = c/(fs*2)
ar=int(2*maxWidth/l2mbd)+1 #Określamy gęstość poziomą siatki zgodnie z twierdzeniem ???
maxdepth=0.0610304 #maksymalna głębokość [m]
mindepth=0.00298 # minimalna głębokość [m]
br=(maxdepth-mindepth)/((maxWidth*2)/ar)
X0 = np.linspace(-maxWidth,maxWidth,ar)
R0 = np.linspace(mindepth,maxdepth,br)
Grid_div_x,Grid_div_z = np.meshgrid(X0,R0,indexing='xy') #tablice ze współrzędnymi
Podobnie jak wcześniej wartość natężenia w punkcie (x,y) będziemy liczyć jako sumę natężeń sygnałów z każdego przetwornika odbiorczego odpowiednio opóźnionych:
[math]S(x,y)=\sum^{NT-1}_{k=0} s_k(t_{tr}(x,y)+t_k(x,y)) [/math]
gdzie [math]s_k(i)[/math] to sygnał z danego przetwornika odbiorczego; [math]t_{tr}(x,y)[/math] czas przejścia fali do punktu obrazowania; [math]t_k(x,y)[/math] czas przejścia fali odbitej z punktu do przetwornika [math]k[/math].
Dla ustalonego przetwornika [math]k[/math] i punktu (x,y) opóźnienie liczone jest jako czas przejścia fali od początku nadania do punktu i spowrotem do danego przetwornika odbiorczego.
W przypadku PWI możemy założyć, że
[math]t_{tr}(x,y) = (x\cdot sin(\theta)+y\cdot cos(\theta))/c + t_{start} + c [/math]
gdzie [math]c[/math] jest pewnym przesunięciem stałym dla danego nadania - wyznaczonym przez położenie chwili zerowej (moment rozpoczęcia akwizycji).
Do wykonania
- Stworzyć funkcję, która dla zadanej siatki generować będzie wartości natężeń dla pojedynczego strzału z danego kąta. Podobnie jak wcześniej, interesuje nas obraz w skali decybelowej z ograniczoną dynamiką.
- Porównać obraz syntezowany dla zadanej wcześniej siatki z obrazem syntezowanym dla siatki o dwukrotnie większych odległościach między punktami
Złożenie (compounding) obrazów z nadań pod różnymi kątami
Jednym z prostych sposobów poprawienia obrazu jest złożenie ze sobą (np. przez średnią arytmetyczną) kilku obrazów z różnych kątów nadawczych.