Luận án Phương pháp sử dụng kỹ thuật nội suy để cải tiến DBIM

ầu được quan tâm do sự phát triển mạnh của phần cứng và phần mềm cho các thiết bị cùng với khả năng giải quyết những khuyết điểm còn tồn tại của B-mode. Kỹ thuật siêu âm cắt lớp được giới thiệu dưới đây cùng với những điểm mạnh và điểm yếu của nó. Chụp cắt lớp siêu âm sử dụng tán xạ ngược Khi một tia tới sóng âm gặp một môi trường không đồng nhất thì một phần năng lượng sẽ bị tán xạ theo mọi hướng. Bài toán chụp cắt lớp siêu âm bao gồm ước lượng sự phân bố của các tham số (tốc độ âm thanh, sự suy giảm âm, mật độ và những thứ khác) tán xạ cho một tập các giá trị đo của trường tán xạ bằng việc giải ngược các phương trình sóng. Vì thế, chụp cắt lớp siêu âm cho thấy định lượng thông tin của vật thể dưới sự khảo sát hay kiểm tra. Hiện tại mới chỉ có một vài hệ thống lâm sàng chụp siêu âm cắt lớp (utrasonic computerd tomography – UCT), hai trong số đó là CURE [5,6] và HUTT [7]. Tuy nhiên độ phân giải không gian và độ chính xác của các hệ thống này vẫn còn giới hạn vì bỏ qua vấn đề nhiễu xạ. Thiết bị thứ 3, máy scan TMS (Techniscal Medical Systems) [8] sử dụng cá thuật toán tán xạ ngược cho kết quả chính xác hơn. Tuy nhiên, tán xạ ngược âm gặp phải một số hạn chế như trong kết quả chụp cắt lớp Y-sinh không được thành công như các phương pháp tạo ảnh cắt lớp khác (cắt lớp X – quang, cắt lớp hạt nhân, và chụp cộng hưởng từ) thường được sử dụng cho chuẩn đoán y tế. Đầu tiên, phương pháp tán xạ ngược gặp phải vấn đề về hội tụ khi tái tạo lại đối tượng với “độ tương phản” lớn ( độ tương phản quyết định bởi tính chất của môi trường, biểu hiện bởi sự tán xạ âm thanh nhiều hay ít, chính là chênh lệch tốc độ truyền sóng giữa 2 môi trường). Rằng buộc này cho đến nay đã hạn chế những ứng dụng tán xạ ngược áp dụng cho việc tạo ảnh vùng ngực [9-11]. Thứ 2, số liệu tán xạ phải thu thập ở rất nhiều góc khác nhau từ 00 đến 3600 để thu được chất lượng chụp tốt. Đó cũng là lý do mà nghiên cứu chụp tán xạ ngược siêu âm lại tập trung vào tạo ảnh vùng ngực, để bao trùm được đầy đủ số liệu việc tạo ảnh ở tần số tương đối cao (lên đến 5 MHz). Trong trường hợp tạo ảnh vú, góc bao phủ đầy đủ thu được bằng cách cho vùng vú đó vào trong nước, cách này được sử dụng cho các cặp vợ chồng siêu âm khối u. Cuối cùng, hạn chế của chụp siêu âm cắt lớp sử dụng tán xạ ngược là tốc độ tính toán và chất lượng ảnh tái tạo. Phương pháp chụp cắt lớp sử dụng tán xạ ngược được đánh giá là cho kết quả chính xác và khả quan hơn các phương pháp chụp siêu âm trước đây nhưng vấn đề về tốc độ tính toán là một trở ngại lớn của phương pháp này, trong chuẩn đoán bệnh y học thì yêu cầu về tốc độ cũng như chất lượng cần được đảm bảo. Như vậy chụp cắt lớp siêu âm (thường được áp dụng cho các kiểm tra về ung thư vú) có nhiều ưu điểm trong siêu âm (cho chất lượng ảnh tốt, không độc hại bởi tác động ion) nhưng chất lượng ảnh chụp còn hạn chế, tốc độ tính toán chưa nhanh, vì thế luận văn này trình bày phương pháp nội suy kết hợp với xấp xỉ Born nhằm nâng cao chất lượng ảnh chụp và đảm bảo được tăng tốc độ tính toán trong chụp cắt lớp siêu âm sử dụng tán xạ ngược. Để xuất này cùng với những nghiên cứu khác trong bộ môn nằm trong hướng nghiên cứu chụp ảnh siêu âm cắt lớp tại Khoa ĐTVT [20-23]. Tổ chức luận văn Phần còn lại của luận văn này được tổ chức như sau: Chương 2 trình bày về các nguyên tắc làm việc bao gồm việc trình bày 2 phương pháp BIM (Born iterative method) và DBIM (Distorted born iterative method), so sánh 2 phương pháp. Chương 3 đưa ra phương pháp đề xuất để giải quyết vấn đề đặt ra. Chương 4 đưa ra những kết quả đã đạt được khi áp dụng phương pháp đề xuất cùng với những đánh giá và kết luận về những kết quả đã đạt được. CHƯƠNG 2: NGUYÊN LÝ HOẠT ĐỘNG Born Iterative Method (BIM) Hình 2.1 là sơ đồ cấu hình thu phát của hệ chụp siêu âm cắt lớp. Hình 2.1: Cấu hình hệ đo Việc thực hiện đo thực tế có thể làm theo 2 cách sau: Cách 1: Tất cả các máy phát và máy thu đều cố định trong suốt quá trình đo. Vật thể sẽ được xoay quanh trục trung tâm với 1 bước nhảy xác định. Nhận xét rằng một máy thu và Nr máy phát được đặt đối xứng nhau như hình 2 nhằm đảm bảo không bị hiện tượng dịch pha gây lỗi khi khôi phục ảnh [12] . Cách 2: Cố định vật thể, tại một vị trí máy phát xác định sẽ tiến hành đo trên Nr máy thu ở vị trí đối xứng. Trên thực tế chỉ cần một máy thu nhưng thực hiện Nr lần đo ứng với một vị trí máy phát. Sau đó khi dịch máy phát đi một góc thì Nr máy thu kia cũng tự động dịch chuyển một cách tương ứng. Vùng cần quan tâm (RIO – region of interest) bao gồm vật cần dựng ảnh. Vùng diện tích quan tâm này được chia thành NxN ô vuông (pixel) có kích thước là h. số lượng máy phát là Nt và máy thu là Nr. Với vùng tán xạ hình tròn như trong Hình 2.1, hàm muc tiêu (Object function) được tính bởi công thức (2.1). (2.1) Với c1 và c0 là tốc độ truyền sóng trong đối tượng và tốc độ truyền trong nước, f là tần số sóng siêu âm, ω là tần số góc (ω=2πf),R là bán kính của đối tượng. Phương trình truyền sóng được mô tả: ∇2+k02rpr=-Orpr (2.2) Với là số sóng tính trong môi trường B1, là tần số góc, tín hiệu áp suất tổng (là tổng của áp suất tới và áp suất gây bởi tán xạ), và là hàm mục tiêu. Giải phương trình (2.2) sẽ có nghiệm dạng tích phân tính theo hàm Green như sau: (2.3) Ở đó là sóng tán xạ, là sóng tới và G(.) là hàm Green. Chúng ta sử dụng phương pháp mômen để rời rạc hóa phương trình (2.3) [13] bằng cách chia lưới vùng diện tích chứa đối tượng (xem Hình 2.1). Sóng tới được biểu diễn dạng vecto kích thước N2×1 như sau: (2.4) Và áp suất tán xạ thu được có kích thước NtNr×1: (2.5) Ở đó là ma trận 1×N2 ứng với hệ số G0(r,r’) từ các pixel tới máy thu, là ma trận N2×N2 ứng với hệ số G0(r,r’) giữa các pixel, là ma trận đơn vị, và D(.) là toán tử chéo hóa [13]. Nếu có Nt máy phát và Nr máy thu thì tín hiệu áp suất tán xạ là vector NtNr×1 được viết lại: psc=BiDOp=MO (2.6) Với M=BDp là ma trận kích cỡ NtNr×N2. Những dữ liệu được xử lý sử dụng BIM để khôi phục lại mức độ tương phản tốc độ âm thanh. Bằng cách này chúng ta có thể xác định được có khối u tồn tại trong môi trường. BIM sử dụng xấp xỉ Born để tính toán trong các vòng lặp của bài toán phi tuyến tán xạ ngược. Nếu chúng ta sử dụng nhiều máy phát và máy thu hàm mục tiêu O có thể được ước lượng bằng cách sử dụng quy tắc Tikhonov. O=argminOp-MO22+γO22 (2.7) Với γ là tham số của quy tắc, tham số này cần được lựa chọn một cách cẩn thận vì nó ảnh hưởng đến sự ổn định của hệ thống [16]. Nếu γ mà lớn thì ảnh tái tạo sẽ khó nhưng nếu γ mà nhỏ thì việc tính toán sẽ phức tạp. Nhận xét rằng trong 2 phương trình (2.4) và (2.5) thì cả và đều chưa biết. Để giải quyết ta phải áp dụng xấp xỉ Born loại 1 cho điều kiện ban đầu trước khi bắt đầu vòng lặp: . Bên cạnh đo, tin hiệu sóng tới được cho là: (2.8) Với J0 là hàm Bessel bậc 0, là khoảng cách giữa máy phát và các điểm thứ k trong vùng chia lưới. Tín hiệu trong thực tế có thể đo được bằng cách lấy hiệu số của tín hiệu tại máy thu khi có đối tượng và khi không có đối tượng. Còn trong mô phỏng thì lại có thể được tính bằng phương trình (2.5) sử dụng hàm mục tiêu lý tưởng. Sau đó chúng ta sẽ sử dụng phương pháp BIM để khôi phục hàm mục tiêu bằng cách lặp hai phương trình (2.4) và (2.5) theo thuật toán BIM như sau: Thuật toán 1 The Born iterative method 1: Thiết lập và 2: Tính toán hai ma trận và 3: Lặp cho đến khi RRE < { 4: Tính sử dụng trong phương trình (2.4) 5: Tính RRE ứng với sử dụng (2.9) 6: Cập nhật hàm mục tiêu bằng cách giải (2.5). 7: n=n+1; } Ở đó sai ngưỡng sai số cho trước và RRE được định nghĩa như sau: (2.9) Distorted Born Iterative Method (DBIM) Sử dụng sơ đồ cấu hình hệ đo như trong Hình 2.1 phần BIM, bằng cách sử dụng DBIM để tái tạo lại độ tương phản âm thanh tán xạ để xác định khối u trong môi trường. Giải sử rằng có một không gian vô hạn chứa môi trường đồng nhất chẳng hạn là nước, số sóng là k0. Trong môi trường đó có vật với số sóng là kr phụ thuộc vào không gian trong vật. Phương trình truyền sóng của hệ thống có thể được cho như phương trình (2.2). Viết lại dưới dạng tích phân ta có: pr=pincr+psc(r) (2.10) pr=pincr+O(r')p(r')G0(r-r')dr' (2.11) Ở đó là sóng tán xạ, là sóng tới và G(.) là hàm Green. Or=kr2-k02 (2.12) là hàm mục tiêu cần dược khôi phục từ dữ liệu tán xạ. Bằng phương pháp moment (MoM) áp suất tổng có thể được tính: p=I-C.DOpinc (2.13) Áp suất tán xạ: psc=B.DO. p (2.14) Hai biến chưa biết là p và O trong công thức (2.13) và (2.14), trong trường hợp này áp dụng xấp xỉ Born loại 1 và theo (2.13), (2.14) ta có: Δpsc=B.Dp.ΔO=M.ΔO (2.15) Với M=B.Dp Với mỗi bộ phát và bộ thu, chúng ta có một ma trận M và một giá trị vô hướng Δpsc. Thấy rằng vector chưa biết O có N×N giá trị bằng với số pixel của RIO. Hàm mục tiêu (Object function) có thể được tính bằng cách lặp: On=O(n-1)+∆O(n-1) (2.16) Với On và O(n-1) là giá trị của hàm mục tiêu ở bước hiện tại và bước trước đó. ΔO có thể được tìm bằng quy tắc Tikhonov: ΔO=argmin∆O∆psct-Mt∆O22+γ∆O22 (2.17) Trong đó ∆psc là (NtNr×1) vector chứa giá trị sai khác giữa kết quả đo và kết quả tiên đoán tín hiệu siêu âm tán xạ; Mt là ma trận NtNr×N2 được tạo bởi NtNr phép đo Thuật toán 2: The distorted Born iterative method 1: Chọn giá trị khởi tạo On=O0 2: while(n<Nmax ) or( RRE < ), do { 3: Tính p, psc, C,và B tương ứng On sử dụng (2.13) và (2.14) 4: Tính ∆psc từ giá trị psc đo được và giá trị tiên đoán 5: Tính RRE tương ứng ΔO sử dụng công thức (2.18) 6: Tính giá trị On mới sử dụng (2.16) 7: n=n+1; } RRE=∆pscpsc,m (2.18) [13]. Bài toán ngược Để giải bài toán ngược khi có nhiễu ta phải sử dụng phương pháp “Nonlinear conjugate gradient method” (NCG) [13][14]. Vì thế ta có thuật toán để giải phương trình (2.17) như sau: Thuật toán 3: NCG method 1: Khởi tạo ΔO dưới dạng một vector 0 2: Khởi tạo b(0)=MtH.∆psct . 3: Khởi tạo x(0)=b(0) và r(0)=b(0) . 4: for n=1 đến giá trị lặp lớn nhất, do 5: qn=Mt.x(n-1) 6: α(n)=rn-1H.rn-1/(qnH.qn+γxn-1H.x(n-1)) 7: sn=MtH.qn 8: rn=rn-1-αn(sn+γxn-1) 9: β(n)=rnH.rn/rn-1H.rn-1 10: ΔOn=ΔOn-1+α(n)x(n-1) 11: x(n)=rn+β(n)x(n-1) 12: if r<tolerance, then 13: Break iterations 14: end if 15: end for Để lựa chọn tham số γ ta sử dụng công thức (2.19) [13]: γ=0.5σ02max⁡{10log2RRE,10-4} (2.19) Với σ02 được tính theo phương pháp lũy thừa lặp với xấp xỉ tỉ số Rayleigh [15]. Thuật toán 4: The power iteration method with Rayleigh quotient 1: Khởi tạo vector ngẫu nhiên đơn vị ω0 và σ0(0)2=0 2: for n=1 đến số vòng lặp lớn nhất, do 3: u(n)=Mt.ω(n-1) 4: σ0(n)2=unH.un/(ω(n-1)H.ωn-1) 5: if σ0n2-σ0n-12σ0n2<tolerance then 6: Break iterations 7: end if 8: r(n)=Mt.ω(n-1) 9: s(n)=MtH.r(n) 10: ω(n)=sn/(snH.sn) 11: end for Mô phỏng thuật toán BIM và DBIM Mô phỏng BIM Kịch bản mô phỏng 1: Tần số sóng siêu âm 1MHz Đường kính vùng tán xạ 1mm Vùng tán xạ được chia lưới 11×11 Số lượng máy phát Nt 9 Số lượng máy thu Nr 121 Chênh lệch tốc độ truyền sóng 2% Hình 2.2. là hàm mục tiêu lý tưởng được xây dựng bởi phương trình (2.1) Hình 2.2: Hàm mục tiêu lý tưởng Hình 2.3 và Hình 2.4 mô tả hàm mục tiêu cần khôi phục sau bước 1 và bước 5 (tức là sau vòng lặp thứ nhất và vòng lặp thứ 5) Hình 2.3: Hàm mục tiêu sau bước lặp đầu tiên Hình 2.4: Hàm mục tiêu sau bước lặp thứ 5 Hình 2.5 là đồ thị biểu thị sai số giữa kết quả khôi phục và hàm mục tiêu lý tưởng (sai số err) với err được tính như sau: err=i=1Nj=1NCij-CijCij (2.19) Trong đó Cij và Cij là hàm mục tiêu lý tưởng và ước lượng theo phương (i,j) Hình 2.5: Sai số qua từng bước lặp Kịch bản mô phỏng 2: Tần số sóng siêu âm 1MHz Đường kính vùng tán xạ 1mm Vùng tán xạ được chia lưới 20×20 Số lượng máy phát Nt 9 Số lượng máy thu Nr 121 Chênh lệch tốc độ truyền sóng 2% Hình 2.6 là hàm mục tiêu được xây dựng bởi phương trình (2.1) ứng với những tham số của kịch bản Hình 2.6: Hàm mục tiêu lý tưởng (N=17) Hình 2.7, Hình 2.8 và Hình 2.9 mô tả hàm mục tiêu cần khôi phục sau bước 1, bước 4 và bước 10 (tức là sau vòng lặp thứ nhất, vòng lặp thứ 4 và vòng lặp thứ 10) Hình 2.7: Hàm mục tiêu sau bước lặp thứ 1 Hình 2.8: Hàm mục tiêu sau bước lặp thứ 4 Hình 2.9: Hàm mục tiêu sau bước lặp thứ 5 Mô phỏng DBIM Kịch bản mô phỏng 3: Tần số sóng siêu âm 1MHz Đường kính vùng tán xạ 1mm Vùng tán xạ được chia lưới 20×20 Số lượng máy phát Nt 40 Số lượng máy thu Nr 20 Chênh lệch tốc độ truyền sóng 2% Hàm mục tiêu được xây dựng dựa vào (2.1) được cho ở đồ thị Hình 2.11 Hình 2.10: Hàm mục tiêu lý tưởng Hình 2.11, Hình 2.12 và Hình 2.13 mô tả hàm mục tiêu cần khôi phục sau bước 1, bước 4 và bước 5 (tức là sau vòng lặp thứ nhất, vòng lặp thứ 4 và vòng lặp thứ 5) Hình 2.11: Hàm mục tiêu sau bước lặp thứ 1 Hình 2.12: Hàm mục tiêu sau bước lặp thứ 4 Hình 2.13: Hàm mục tiêu sau bước lặp thứ 5 Hình 2.14 là đồ thị biểu thị sai số giữa kết quả khôi phục và hàm mục tiêu lý tưởng (sai số err) với err được tính như phương trình (2.19) của Kịch bản 2 và Kịch bản 3 Hình 2.14: Sai số qua từng bước lặp của BIM và DBIM So sánh và nhận xét Cả hai phương pháp Born Iterative Method (BIM) và Distorted Born Itertative Method (DBIM) được biết đến như là giải pháp tốt cho nhiễu xạ cắt lớp [17] và DBIM cho tốc độ hội tụ nhanh là giải pháp tốt nhất hiện nay. Như những kết quả ở kịch bản 1, kịch bản 2 và kịch bản 3 thì DBIM cho hội tụ nhanh hơn, đây cũng là phương pháp tốt nhất hiện nay cho siêu âm cắt lớp sử dụng tán xạ ngược. Vì thế luận văn này sẽ trình bày việc áp dụng đề xuất cho DBIM nhằm tăng chất lượng ảnh chụp và giảm thời gian tính toán. CHƯƠNG 3: PHƯƠNG PHÁP ĐỀ XUẤT Đề xuất Phương pháp đề xuất đưa ra gồm 3 phần: DBIM - Propose Gọi tổng số bước lặp của cả quá trình là Niter Step 0: Tìm số lần lặp tối ưu x trước khi nội suy, bước này xác định số lần lặp với ma trận có kích cỡ N1×N1 là bao nhiêu lần trong tổng số bước lặp để thu được ảnh có chất lượng tốt nhất. Step 1: Đầu tiên áp dụng khôi phục cho vùng lưới có kích cỡ N1×N1 chúng ta có thể dễ ràng có được sự hội tụ với chỉ 1 hay 2 lần lặp (là x lần được xác định ở Step 0). Kết quả thu được ở phần này là ảnh có mật độ trung bình của đối tượng. Step 2: Phần 2 áp dụng nội suy cho hàm mục tiêu kích cỡ N1×N1 thu được sau khi kết thúc quá trình lặp ở bên trên, ta thu được ma trận có kích thước N2×N2 (N2=2N1) Step 3: Cuối cùng sử dụng kết quả mà nội suy thu được mang trở lại DBIM (hoặc BIM) lặp (Niter - x)lần để tiếp tục quá trình khôi phục. Ở đây chúng ta sử dụng phương pháp nội suy gần nhất (Nearest neighbor) đây là phương pháp nội suy đơn giản và nhanh nhất. Có nhiều phương pháp nội suy như bilinear, bicubic, splinenhưng chúng ta chỉ tập trung vào phương pháp nội suy gần nhất vì yếu tố tiết kiệm thời gian tính toán [18]. Vậy ta có thuật toán để khôi phục sau: Thuật toán 5: Modified BIM 1: Thiết lập ON10 và p0=pN1inc 2: Tính toán hai ma trận và 3: Lặp cho đến khi (n<Nmax1) or (REEN1<εN1) { 4: Tính pN1 sử dụng ON1(n) trong phương trình (2.4) 5: Tính REEN1 ứng với ON1(n) sử dụng (2.9) 6: Cập nhật hàm mục tiêu ON1(n+1) bằng cách giải (2.5). 7: n=n+1; } 8: Nội suy ON1(n) để thu được ON2(0) 9: khởi tạo p0=pN2inc, n=0 10: Lặp cho đến khi (n<Nmax2) or (REEN1<εN2) { 11: Tính pN2 sử dụng ON2(n) trong phương trình (2.4) 12: Tính REEN2 ứng với ON2(n) sử dụng (2.9) 13: Cập nhật hàm mục tiêu ON2(n+1) bằng cách giải (2.5). 14: n=n+1; } ε là ngưỡng sai số cho trước quyết định bởi ồn nền (noise floor) [19]. Thuật toán 6: Modified DBIM 1: Chọn giá trị khởi tạo ON10 và p0=pN1inc 2: while(n<Nmax1 ) or (REEN1<εN1), do { 3: Tính p, psc, và BN1, CN1 tương ứng ON1(n) sử dụng (2.13) và (2.14) 4: Tính ∆psc từ giá trị psc đo được 5: Tính RRE tương ứng ΔON1(n) sử dụng công thức (2.18) 6: Tính giá trị ON1(n+1) mới sử dụng (2.16) 7: n=n+1; } 8: Nội suy ON1(n) để thu được ON2(0) 9: Khởi tạo p0=pN2inc, n=0 10: while(n<Nmax2 ) or (REEN2<εN2), do { 11: Tính p, psc, và BN2, CN2 tương ứng ON2(n) sử dụng (2.13) và (2.14) 12: Tính ∆psc từ giá trị psc đo được 13: Tính RRE tương ứng ΔON2(n) sử dụng công thức (2.18) 14: Tính giá trị ON2(n+1) mới sử dụng (2.16) 15: n=n+1; } ε là ngưỡng sai số cho trước quyết định bởi ồn nền (noise floor) [19]. Tìm giá trị x tối ưu. Theo như phương pháp đề xuất trong Chương 4 trước hết ta phải tìm số lần lặp x tối ưu cho ma trận N1×N1. Với tổng số bước lặp của cả quá trình là 4. Ta có thuật toán để so sánh sau: Thuật toán 7: 1: for x = 1 đến 3, do 2: DBIM – Propose. 3: Tính err theo công thức (2.19) 4: Vẽ đồ thị err ứng với từng giá trị của x 5: end for Như vậy sau khi thực hiện xong Thuật toán 7 ta có thể tìm được giá trị x tối ưu. Ta xét kịch bản sau: Kịch bản 4: Simulation parameters: Frequency = 1MHz N1 = 10, N2 = 20 Diameter of scatter area = 5*landa Percent of sound contrast 1% 5% Gaussian noise (SNR = 26 dB) Detector = 20, Transmiter = 40 Thực hiện Thuật toán 7 ta thu được bảng và đồ thị sau Bảng 4.1: err ứng với từng giá trị của x sau tổng số bước lặp là 4 (N1 = 10) x 1 2 3 err 0.0346 0.1245 0.2275 Hình 4.1: Đồ thị biểu diễn mối liên hệ giữa số lần lặp x và sai số error (N1 = 10) Kịch bản 5: Simulation parameters: Frequency = 1MHz N1 = 11, N2 = 22 Diameter of scatter area = 5*landa Percent of sound contrast 1% 5% Gaussian noise (SNR = 26 dB) Thực hiện Thuật toán 7 ta thu được bảng và đồ thị sau Bảng 4.2: err ứng với từng giá trị của x sau tổng số bước lặp là 4 (N1 = 11) x 1 2 3 err 0.0475 0.1063 0.2062 Hình 4.2: Đồ thị biểu diễn mối liên hệ giữa số lần lặp x và sai số error (N1 =11) Kịch bản 6: Simulation parameters: Frequency = 1MHz N1 = 12, N2 = 24 Diameter of scatter area = 5*landa Percent of sound contrast 1% 5% Gaussian noise (SNR = 26 dB) Thực hiện Thuật toán 7 ta thu được bảng và đồ thị sau Bảng 4.3: err ứng với từng giá trị của x sau tổng số bước lặp là 4 (N1 = 12) x 1 2 3 err 0.0180 0.0540 0.0997 Hình 4.3: Đồ thị biểu diễn mối liên hệ giữa số lần lặp x và sai số error (N1 = 12) Nhận xét: Qua các kịch bản 4, 5 và kịch bản 6 ta thấy đươc giá trị err tăng qua từng giá trị của x, như vậy với x=1 chính là số lần lặp tối ưu. Vậy chọn x=1 và sử dụng trong phương pháp đề xuất CHƯƠNG 4: KẾT QUẢ Trong chương này chúng ta sẽ đi vào thực hiện mô phỏng và so sánh giữa phương pháp đề xuất (DBIM – Propose) với phương pháp truyền thống (DBIM – Conventional) Mô phỏng và đánh giá về mặt chất lượng Kịch bản 7: Simulation parameters: Frequency = 1MHz N1 = 10, x = 1, N2 = 20 N = 20 Diameter of scatter area = 5*landa Percent of sound contrast 1% 5% Gaussian noise (SNR = 26 dB) Detector = 20, Transmiter = 40 Hình 4.4: Hàm mục tiêu lý tưởng (N = 20) Với DBIM thường ta có bảng kết quả err qua các bước lặp ở Bảng 4.4: Bảng 4.4: err của DBIM qua từng bước lặp (N = 20) Iter 1 2 3 4 err 0.5504 0.2769 0.1588 0.1059 Tổng thời gian tính toán của phương pháp này là 46.2 giây Với phương pháp đề xuất ta có err như ở Bảng 4.5: Ở đây ta không xét err ở bước lặp 1 vì thực hiện với ma trận N1×N1 Bảng 4.5: err của DBIM - Đề xuất qua từng bước lặp (N1 = 10,N2 = 20) Iter 1 2 3 4 err ---- 0.2374 0.1209 0.0346 Tổng thời gian tính toán của phương pháp này là 36.1 giây Hàm mục tiêu khôi phục qua các bước lặp của 2 phương pháp được cho ở Hình 4.5, Hình 4.6, Hình 4.7 và Hình 4.8 là hàm mục tiêu được khôi phục qua các bước lặp. a 10×10 pixel tomography b 20×20 pixel tomography Hình 4.5: Kết quả khôi phục sau bước lặp đầu tiên (N1 = 10, N = 20) a 20×20 pixel using interpolation b 20×20 pixel tomography Hình 4.6: Kết quả khôi phục sau bước lặp thứ 2 (N1 = 10, N = 20) a 20×20 pixel using interpolation b 20×20 pixel tomography Hình 4.7: Kết quả khôi phục sau bước lặp thứ 3 (N1 = 10, N = 20) a 20×20 pixel using interpolation b 20×20 pixel tomography Hình 4.8: Kết quả khôi phục sau bước lặp thứ 4 (N1 = 10, N = 20) Ta vẽ đồ thị so sánh sai số err của phương pháp đề xuất và phương pháp truyền thống, đồ thị được cho như ở Hình 4.9, đường màu đỏ là đường sai số err của phương pháp DBIM – Thường và đường màu đen nét đứt là của phương pháp DBIM – Đề xuất. Rõ ràng là đường màu đen sẽ cho kết quả khôi phục tốt hơn vì nó nằm dưới (sai số err là nhỏ hơn so với đường màu đỏ) Hình 4.9: Đồ thị so sánh err của DBIM Propose và DBIM conventional (N =20) Kịch bản 8: Simulation parameters: Frequency = 1MHz N1 = 11, x = 1, N2 = 22 N = 22 Diameter of scatter area = 5*landa Percent of sound contrast 1% 5% Gaussian noise (SNR = 26 dB) Detector = 22, Transmiter = 44 Hình 4.10. Hàm mục tiêu lý tưởng (N = 22) Với DBIM thường ta có bảng kết quả err qua các bước lặp ở Bảng 4.6: Bảng 4.6: err của DBIM qua từng bước lặp (N = 22) Iter 1 2 3 4 err 0.5648 0.2620 0.1504 0.1004 Tổng thời gian tính toán của phương pháp này là 71.7 giây Với phương pháp đề xuất ta có err ở Bảng 4.7: Bảng 4.7: err của DBIM - Đề xuất qua từng bước lặp (N1 = 11,N2 = 22) Iter 1 2 3 4 err ---- 0.2630 0.1305 0.0475 Tổng thời gian tính toán của phương pháp này là 57.4 giây Hàm mục tiêu khôi phục qua các bước lặp của 2 phương pháp được cho ở Hình 4.11, Hình 4.12, Hình 4.13 và Hình 4.14 là hàm mục tiêu được khôi phục qua các bước lặp. a 11×11 pixel tomography b 22×22 pixel tomography Hình 4.11: Kết quả khôi phục sau bước lặp thứ 1 (N1 = 11, N = 22) a 22×22 pixel using interpolation b 22×22 pixel tomography Hình 4.12: Kết quả khôi phục sau bước lặp thứ 2 (N1 = 11, N = 22) a 22×22 pixel using interpolation b 22×22 pixel tomography Hình 4.13: Kết quả khôi phục sau bước lặp thứ 3 (N1 = 11, N = 22) a 22×22 pixel using interpolation b 22×22 pixel tomography Hình 4.14: Kết quả khôi phục sau bước lặp thứ 4 (N1 = 11, N = 22) Ta vẽ đồ thị so sánh sai số err của phương pháp đề xuất và phương pháp truyền thống, đồ thị được cho ở Hình 4.15, đường màu đỏ là đường sai số err của phương pháp DBIM – Đề xuất và đường màu hồng là của phương pháp DBIM – Thường. Rõ ràng là đường màu đỏ sẽ cho kết quả khôi phục tốt hơn vì nó nằm dưới. Hình 4.15: Đồ thị so sánh err của DBIM Propose và DBIM conventional (N = 22) Về mặt hiệu suất thời gian thì phương pháp DBIM – Propose ở Kịch bản 8 có tổng thời gian tính toán là 71.7 giây, còn DBIM – Conventional là 57.4 giây. Như vậy phương pháp đề xuất cho tốc độ tính toán nhanh hơn 20% so với phương pháp thông thường. Ở Kịch bản 7 và Kịch bản 8 ta đã khảo sát sự khôi phục với vùng cần quan tâm RIO có độ phân giải tương đối thô, ở hai kịch bản tiếp theo này ta sẽ khảo sát với vùng được chia thành nhiều ô hơn. Kịch bản 9: Simulation parameters: Frequency = 1MHz N1 = 17, x = 1, N2 = 34 N = 34 Diameter of scatter area = 5*landa Percent of sound contrast 1%