Tài liệu Biến hình và xử lý ảnh - Ghép hình: Báo cáo cuối kỳ
Môn: Biến hình và xử lý ảnh.
Giảng viên: thạc sĩ Phạm Thế Bảo.
Tên đề tài: Ghép hình
Nhóm sinh viên thực hiện:
Huỳnh Nguyễn Duy Nhân Mssv: 0511027
Nguyễn Quang Thái Mssv: 0511323
1.Cơ sở lý thuyết:
Xử lý ảnh trong miền tần số
1.1 Biến đổi Fourier 2 chiều rời rạc :
Đặt f(x,y) ký hiệu cho ảnh M*N với x = 0,1,2,…,M-1 và y = 0,1,2,…,N-1.Trong không gian hai chiều , biến đổi Fourier rời rạc (DFT) của f , ký hiệu bởi F(u,v) , được cho bởi phương trình sau:
Cho u = 0,1,2,…M-1 và v = 0,1,2,..,N-1. Chúng ta có thể mở rộng hàm mũ thành dạng sin hay cosine với biến u và v trong việc xác định tần số của nó (x và y được tính tổng). Miền tần số đơn giản là hệ tọa độ được dựa theo F(u,v) với u và v là các biến tần số. Nó tương tự như là miền không gian trong chương trước, theo đó hệ tọa độ được dựa theo f(x,y) với x, y là biến không gian. Vùng hình chữ nhật M*N đươc định nghĩa bởi u= 0,1,2,……,M-1 và v = 0,1,2,……,N-1 thường được xem như là hình chữ nhật tần số. Rõ ràng ...
33 trang |
Chia sẻ: hunglv | Lượt xem: 1544 | Lượt tải: 2
Bạn đang xem trước 20 trang mẫu tài liệu Biến hình và xử lý ảnh - Ghép hình, để tải tài liệu gốc về máy bạn click vào nút DOWNLOAD ở trên
Báo cáo cuối kỳ
Môn: Biến hình và xử lý ảnh.
Giảng viên: thạc sĩ Phạm Thế Bảo.
Tên đề tài: Ghép hình
Nhóm sinh viên thực hiện:
Huỳnh Nguyễn Duy Nhân Mssv: 0511027
Nguyễn Quang Thái Mssv: 0511323
1.Cơ sở lý thuyết:
Xử lý ảnh trong miền tần số
1.1 Biến đổi Fourier 2 chiều rời rạc :
Đặt f(x,y) ký hiệu cho ảnh M*N với x = 0,1,2,…,M-1 và y = 0,1,2,…,N-1.Trong không gian hai chiều , biến đổi Fourier rời rạc (DFT) của f , ký hiệu bởi F(u,v) , được cho bởi phương trình sau:
Cho u = 0,1,2,…M-1 và v = 0,1,2,..,N-1. Chúng ta có thể mở rộng hàm mũ thành dạng sin hay cosine với biến u và v trong việc xác định tần số của nó (x và y được tính tổng). Miền tần số đơn giản là hệ tọa độ được dựa theo F(u,v) với u và v là các biến tần số. Nó tương tự như là miền không gian trong chương trước, theo đó hệ tọa độ được dựa theo f(x,y) với x, y là biến không gian. Vùng hình chữ nhật M*N đươc định nghĩa bởi u= 0,1,2,……,M-1 và v = 0,1,2,……,N-1 thường được xem như là hình chữ nhật tần số. Rõ ràng hơn, hình chữ nhật tần số cùng kích thước với ảnh đầu vào.
Đảo lại, biến đổi Fourier rời rạc cho bởi:
Với x = 0,1,2,……,M-1 và y = 0,1,2,…….N-1. Vì vậy, cho F(u,v), chúng ta có thể thu được f(x,y) trở lại bằng biến đổi ngược DFT. Giá trị của F(u,v) trong phương trình này đôi khi được suy ra như hệ số của khai triển Fourier.
Trong một vài công thức của DFT, 1/MN đứng trước phép biến đổi và trong một số khác nó nằm ở phép biến đổi ngược. Để thống nhất việc cài đặt trong Matlab của biến đổi Fourier, chúng ta có thể giả sử trong suốt cuốn sách này 1/MN đứng trước phép biến đổi ngược như đã nêu trong phươn trình. Bởi vì dãy chỉ số bắt đầu trong Matlab bắt đầu ở 1, không phải là 0, F(1,1) và f(1,1) trong Matlab tương ứng với giá trị toán học F(0,0) và f(0,0) trong biến đổi và biến đổi ngược của nó.
Giá trị của biến đổi tại gốc tọa độ của miền tần số [ví dụ F[0,0]] được gọi là thành phần dc của thành phần biến đổi Fourier. Thuật ngữ này xuất phát từ khoa học điện tử, trong đó ‘dc’ ký hiệu cho ‘direct current’ (hiện tại của tần số 0). Chúng ta có thể dễ dàng chứng minh rằng F(0,0) bằng MN lần giá trị trung bình của f(x,y).
Ngay cả khi f(x,y) là số thực, thì biến đổi của nó trong dạng tổng quát vẫn là số phức. Phương pháp chính của việc phân tích là tính toán phổ của nó [ ví dụ độ lớn của F(u,v)] và biểu diễn nó như là một ảnh. Đặt R(u,v) và I(u,v) biểu diễn thành phần thực và ảo của F(u,v), phổ Fourier được định nghĩa như sau :
Góc pha của biến đổi được định nghĩa là :
Hai hàm có trước dùng để biểu diễn F(u,v) trong mối liên quan biểu diên cực của giá trị số phức :
Năng lượng phổ được định nghĩa như là hình vuông của độ lớn :
Để dễ hình dung chúng ta có thể chọn một trong hai cách biểu diễn là |F(u,v)| hay P(u,v).
Nếu f(x,y) là số thực, biến đổi Fourier là đối xứng liên hợp qua gốc tọc độ :
Điêu này cũng có nghĩa là phổ Fourier cũng đối xứng qua gốc tọa độ :
Chúng ta cũng có thể kiểm chứng rằng :
Nói cách khác, DFT là vô hạn tuần hoàn theo cả u và v, với chu kỳ xác định trước bởi M và N. Tuần hoàn cũng là một đặc trưng của biến đổi ngược DFT :
Có nghĩa là môt ảnh thu được bởi phép biến đổi ngược Fourier cũng vô hạn tuần hoàn. Điều này có thể gây ra nhầm lân vì không phải ảnh nào thu được từ biến đổi fourier cũng đều tuần hoàn. Chúng ta có thể hiểu nó đon thuần là một tính chất toán học của DFT và phép biến đổi ngược của nó. Chúng ta cũng nên nhớ rằng cài đặt tính toán của DFT chỉ tính trên một chu kỳ, vì vậy chúng ta làm việc với kích thước MN.
Hàm tuần hoàn rất quan trọng khi chúng ta xem xét dữ liệu DFT liên quan đến chu kỳ biến đổi như thế nào. Ví dụ, hình 1.1(a) biểu diễn phổ của biến đổi một chiều, F(u). Trong trường hợp này, biểu thức tuần hoàn F(u) = F(u+M), từ đó ta có |F(u)| = |F(u+M)|. Đồng thời do tính đối xứng nên |F(u)| = |F(-u)|. Tính chất tuần hoàn chỉ ra rằng F(u) có chu kỳ chiều dài M, và tính đối xứng chỉ ra rằng độ lớn của biến đổi được trung tâm hóa tại gốc tọa độ, như trong hình 1.1(a).
Hình 1.1 :
(a) Phổ Fourier minh họa mỗi nửa chu kỳ trong khoảng [0, M-1].
(b) Trung tâm của phổ trong cùng một khoảng thu được bằng việc nhân f(x) và (-1)x trước khi tính toán biến đổi Fourier.
Hình này và những lời chú thích trước đó minh họa rằng độ lớn của biến đổi giá trị từ M/2 tới M-1 là sự lăp lại của giá trị trong một nửa chu kỳ bên trái gốc tọa độ. Bởi vì biến đổi DFT một chiều chỉ được cài đặt cho M điểm (ví dụ cho giá trị của u trong khoảng [0, M-1]), nên nó tính toán theo biến đổi một chiều của mỗi nửa chu kỳ trong khoảng [0, M-1]. Không khó để chứng minh rằng một chu kỳ thích hợp là chu kỳ thu được bằng việc nhân f(x) với (-1)x trước khi tính toán biến đổi. Nó có tác dụng chuyển gốc tọa độ của biến đổi tới điểm u = M/2, như hình 1.1(b). Bây giờ giá trị của phổ tại u = 0 trong hình 1.1(b) tương ứng với |F(-M/2)| trong hình 1.1(a). Tương tự giá trị tại |F(M/2)| và |F(M-1)| trong hình 1.1(b) tương ứng với |F(0)| và |F(M/2-1)| trong hình 1.1(a).
Điều này củng xảy ra tương tự trong hàm biến đổi hai chiều. Tính toán biến đổi DFT hai chiều bây giờ tương ứng với điểm biến đổi trong khoảng hình chữ nhật như trong hình 1.2(a), vùng hình chữ nhật đổ bóng chỉ giá trị của F(u,v) thu được bằng cài đặt của phương trình biến đổi Fourier hai chiều được định nghĩa trong phần đầu tiên của chương này. Vùng hình chữ nhật đứt nét là sự lặp lại, giống như hình 1.1(a). Phân tích một cách trực quan ta thấy rằng đó chính là sự di chuyển phổ từ gốc tọa độ đến trung tâm của vùng tần số. Việc này được thực hiện bằng cách nhân f(x,y) với (-1)x+y trước khi tính toán biến đổi Fourier hai chiều. Chu kỳ khi ấy sẽ được biểu diễn rõ ràng hơn như trong hình 1.2(b). Như đã bàn luận ở phần trước cho hàm 1 chiều, giá trị của phổ tại (M/2, N/2) trong hình 1.2(b) là cùng giá trị với phổ tại (0,0) trong hình 1.2(a), và giá trị tại (0,0) trong hình 1.2(b) thì cùng giá trị với giá trị tại (-M/2, - N/2) trong hình 1.2(a). Tương tự, giá trị tại (M-1, N-1) trong hình 1.2(b) là cùng giá trị với giá trị tại (M/2-1, N/2-1) trong hình 1.2(a).
Hình 1.2 (a) M*N phổ Fourier (phần đổ bóng), chỉ ra mỗi ¼ chu kỳ chứa trong phổ dữ liệu. (b) Phổ thu được bằng việc nhân f(x,y) với (-1)x+y trước khi tính toán biến đổi Fourier. Chỉ có một chu kỳ được biểu diễn đổ bóng. bởi vì đây là dữ liệu sẽ thu được bằng việc cài đặt của phương trình cho F(u,v).
Những bàn luận ở phần trên trong việc trung tâm hóa bằng cách nhân f(x,y) và (-1)x+y là một khái niệm quan trọng. Khi làm việc trong MATLAB, tính toán không bao gồm việc nhân (-1)x+y và sau đó sắp xếp dữ liệu lại bằng hàm fftshift. Hàm này và cách sử dụng nó sẽ được bàn trong phần kế tiếp.
1.2 Tính toán và hình dung biến đổi DFT hai chiều trong MATLAB:
Biến đổi DFT và biến đổi ngược của nó thu được trong thực nghiệm bằng thuật toán biến đổi Fourier nhanh (FFT). FFT của ảnh kích thước M*N f thu được trong hộp công cụ với hàm fft2 với cú pháp như sau:
F = fft2(f).
Hàm này trả về biến đổi fourier cũng có kích thước M*N, với dữ liệu được sắp xếp như trong hình 1.2(a); đó là bởi vì gốc tọa độ được đặt ở góc trên bên trái, và với 4 góc phần tư gặp nhau tại tâm của vùng tần số.
Như giải thích trong phần 1.3.1, sẽ cần thiết khi ảnh đầu vào được thêm vào những số 0 khi biến đổi Fourier được sử dụng như hàm lọc. Trong trường hợp này, cú pháp trở thành
F = fft2(f, P, Q)
Với cú pháp này, fft2 thêm vào đầu vào những số 0 được yêu cầu để hàm kết quả có kích thước P*Q.
Phổ Fourier thu được bằng cách sử dụng hàm abs:
S = abs(F)
Hàm tính toán độ lớn (căn bậc hai của tổng bình phương phần thực và phần ảo) của mỗi phần tử trong mảng.
Phân tích trực quan của phổ bẳng việc biểu diễn nó như ảnh là một khía cạnh quan trọng trong miền tần số. Như minh họa xem như chúng ta có ảnh mẫu f, trong hình 1.3(a). Chúng ta tính toán biển đổi Fourier của nó và biểu diễn phổ theo những bước sau:
>> F = fft2(f);
>> S= abs(F);
>> imshow(S, [ ])
Hình 1.3(b) minh họa kết quả. 4 điểm sáng màu trắng trong góc của ảnh là bởi vì tính chất đã được đề cập trong phần trước.
Hàm fftshift có thể sử dụng để di chuyển gốc tọa độ của biến đổi tới trung tâm của vùng tần số. Cú pháp như sau:
Fc = fftshift(F)
Trong đó F là biến đổi fft2 và Fc là biến đổi trung tâm hóa. Hàm fftshift thực thi bằng việc trao đổi các góc phần tư của F. Ví dụ, nếu a = [1 2;3 4], fftshift(a) = [4 3; 2 1]. Khi áp dụng biến đổi này sau tính toán fft2 thì kết quả cho ra cũng giống như nhân (-1)x+y trước khi tính toán. Tuy nhiên điều đó không có nghĩa là hai biến đổi này có thể thay thế lẫn nhau. Cụ thể hơn, nếu ta ký hiệu F là biến đổi fourier thì F[(-1)x+yf(x,y)] bằng fftshift(fft2(t)) nhưng lại không bằng fft2(fftshift(f)).
Trong ví dụ hiện tại chúng ta thêm:
>> Fc = fftshift(F);
>> imshow(abs(Fc),[ ]).
Tương ứng với ảnh 1.3(c). Kết quả của việc trọng tâm hóa được minh chứng trong ảnh này
Mặc dầu sự chuyển đổi được hoàn thành như mong muốn, nhưng sự sắp xếp động của giá trị trong phổ thì quá lớn (0 tới 24000) trội hơn so với 8 bit biểu diễn giá trị độ sáng. Như đã bàn luận trong mụch 3.2.2, khó khăn này có thể giải quyết bằng biến đổi log. Ta thực hiện những câu lệnh sau:
>> S2 = log(1+abs(Fc));
>> imshow(S2, [ ])
Kết quả trong hình 1.3(d).Sự gia tăng những chi tiết nhìn thấy được được minh chứng trong ảnh.
Hình 1.3
a b
c d
(a) Ảnh đơn giản
(b) Phổ fourier
(c) Trung tâm hóa phổ
(d) Phổ dễ nhìn hơn thông qua biến đổi log.
Hàm ifftshift nghịch đảo lại việc trông tâm hóa. Cú pháp là:
F = ifftshift(Fc)
Hàm này chuyển đổi một hàm có trung tâm tại tâm hình chữ nhật thành một hàm có trung tâm tại góc trên bên trái hình chữ nhật. Chúng ta sẽ sử dụng tính chất này trong mụch 1.4.
Chúng ta nên nhớ rằng trung tâm của miền tần số là tại (M/2,N/2) nếu biến u và v chạy từ 0 đến M-1 và N-1. Ví dụ, tâm của miền tần số có kích thước 8*8 là tại điểm (4,4), chính xác là điểm thứ 5 dọc theo trục tọa độ tính từ (0,0). Trong MATLAB khi biến chạy từ 1 tới M và 1 tới N thì trọng tâm nằm tại [M/2 +1, N/2 +1]. Trong trường hợp này trung tâm của chúng ta nằm tại điểm (5,5), tính từ (1,1).
Nếu M, N là lẻ, trung tâm trong MATLAB được tính bằng cách làm tròn M/2, N/2 tới số nguyên gần nhất nhỏ hơn nó. Phần còn lại giống như những phân tích trong phần trên. Ví dụ, trung tâm của vùng 7*7 là (3,3) nếu tính từ (0,0) và (4,4) nếu tính từ (1,1). Trong cả hai trường hợp trên, trung tâm là điểm thứ 4 tại gốc tọa độ. Nếu chỉ một chiều là lẻ thì tọa độ trung tâm theo chiều ấy cũng đượclàm tròn tương tự. Sử dụng hàm floor trong MATLAB (nhớ rằng gốc tọa độ đặt tại (1,1)), trung tâm của miền tần số được tính như sau :
[floor(M/2) +1, floor(N/2) +1]
Trung tâm tính theo biểu thức này được dùng cho cả giá trị chẵn và lẻ của M, N.
Cuối cùng, chúng ta cần biết rằng biến đổi ngược Fourier được tính toán dựa theo hàm ifft2, có cú pháp:
>> f = ifft2(F)
Trong đó F là biến đổi Fourier của f là ảnh kết quả. Nếu F là số thực thì kết quả cũng là số thực. Tuy nhiên trong thực nghiệm thì đầu ra của ifft2 thường có phần ảo rất nhỏ (nhờ phép làm tròn ) là đặc trưng của tính toán. Vì vậy người ta đi tính phần thưc của kết quả sau khi tính hàm ngược để thu được ảnh bao gồm giá trị thực. Ta thực hiện như sau :
>> f = real(ifft2(F));
Như trong trường hợp trước, hàm này thay thế ifft2(F, P, Q), là lệnh thêm 0 vào F để có kích thước P*Q trước khi tính hàm ngược. Tuy nhiên lựa chọn này không được dùng trong sách này.
1.3 Lọc trong miền tần số:
Lọc trong miền tần số là khái niệm đơn giản. Trong phần này chúng tôi đưa ra tóm tắt tổng quan những khái niễm xoay quanh lọc trong miền tần số và cài đặt nó trong MATLAB.
1.3.1 Khái niệm căn bản:
Khái niệm căn bản cho cả lọc tuyến tính trong cả miền không gian và tần số là lý thuyết tích chập. Ta có thể biểu diễn như sau:
Hay ngược lại:
Ở đây, biểu tượng “*” chỉ tích chập giữa hai hàm, các biểu thức ở hai phía tạo thành những cặp biến đổi Fourier. Ví dụ, trong biểu thức thứ nhất chỉ ra rằng, tích chập của hai hàm không gian có thể tính bằng bến đổi ngược Fourier của tích từng biến đổi Fourier của hai hàm. Ngược lại, ta có thể tính tích biến đổi Fourier của hai bằng việc lấy tính tích chập của hai hàm đó. Ta cũng có suy diễn tương tự cho biểu thức thứ hai.
Trong phần lọc, chúng ta quan tâm trước hết đến 2 biểu thức này. Lọc trong miền không gian chúng ta tính tích chập của ảnh f(x,y) với mặt nạ lọc h(x,y). Theo lý thuyết tích chập, chúng ta có thể thu được cùng một kết quả trong miền tần số bằng việc nhân F(u,v) với H(u,v), (H(u,v) là biến đổi Fourier của lọc không gian). Người ta xem H(u,v) như hàm lọc chuyển đổi.
Ý tưởng lọc trong miền tần số là tìm hàm lọc chuyển đổi hiệu chỉnh F(u,v) theo một cách đặt biệt nào đó. Ví dụ, bộ lọc trong hình 1.4(a) có hàm chuyển đổi có dạng như sau:
Hình 1.4
Hàm chuyển đổi của hàm lọc băng thông thấp
Định dạng sử dụng cho bộ lọc DFT. Chú ý rằng có rất nhiều hàm lọc miền tân số.
Khi hàm lọc nhân với F(u,v), nó làm mỏng đi thành phần tần số cao của F(u,v), trong khi những thành phần tần số thấp không thay đổi. Bộ lọc với những đặc trưng như vậy gọi là bộ lọc băng thông thấp. Như bàn luận trong mụch 1.5.2, kết quả của lọc băng thông thấp là làm ảnh mờ đi (mịn hóa ảnh). Hình 1.4(b) chỉ ra một bộ lọc tương tự sau khi nó được xử lý với fftshift. Đây là định dạng bộ lọc được sử dụng thương xuyên trong sách khi làm việc với miền tần số trong đó biến đổi Fourier của đầu vào chưa được trọng tâm hóa.
Dựa theo lý thuyết tích chập, chúng ta biết rằng để thu được ảnh lọc tương ứng trong miền không gian chúng ta đơn giản là tính biến đổi ngược của biến đổi Fourier của tích H(u,v)F(u,v). Điều quan trọng cần nhớ là tiến trình này cũng đồng nhất với việc sử dụng tích chập trong miền không gian với mặt nạ h(x,y), là hàm biến đổi Fourier ngược của H(u,v). Trong thực nghiệm, tích chập tổng quát trong miền không gian là sử dụng một mặt nạ nhỏ để cố thu được của miền tần số tương ứng của nó.
Như đã lưu ý trong mụch 1.1, ảnh và biến đổi của nó tự động được xử lý tuần hoàn nếu chúng ta chọn DFT để cài đặt hàm lọc. Không khó để nhận ra rằng hàm tích chập tuần hoàn có thể là nguyên nhân ảnh hưởng tới sự tác động giữa những chu kỳ đó liên quan tới khoảng thời gian của những phần khác 0 của hàm. Sự tác động này gọi là lỗi vòng tròn khép kín (wraparound error), có thể tránh việc thêm vào hàm những số 0.
Giả sử rằng hàm f(x,y) và h(x,y) có kích thước A*B và C*D. Chúng ta định dạng hai hàm mở rộng, đều có kích thước P*Q bằng việc thêm 0 vào hàm f và g. Chúng ta có thể chỉ ra được lỗi vòng tròn khép kín có thể tránh được bằng cách chọn:
P>= A+C-1.
và
Q>= B+D-1.
Phần lớn trong chương này đều làm việc với hàm có cùng kích thước M*N, trong trường hợp này chúng ta thêm giá trị : P>=2M-1 và Q>=2N-1.
Hàm này gọi là hàm paddedsize, tính toán giá trị chẵn cực tiểu của P và Q thỏa yêu cầu phương trình trước đó. Ta có thể lựa chọn thêm vào ảnh để định dạng ảnh là hình vuông (có kích thước chiều dài và rộng bằng nhau). Việc thực thi thuật toán FFT phụ thuộc mạnh mẽ vào số thừa số nguyên tố trong P và Q. Thuật toán tổng quát sẽ nhanh hơn nếu P và Q là số chính phương và chậm đi khi P và Q là số nguyên tố. Trong thực nghiệm thì chúng ta nên sử dụng ảnh vuông và bộ lọc để công việc lọc xảy ra trên cùng cả hai hướng của ảnh. Hàm paddedsize cung cấp sự linh hoạt cho ta thực hiện lời khuyên trên thông qua việc lựa chọn tham số.
Trong hàm paddedsize, vecto AB, CD, và PQ có các phần tử [A B], [C D], [P Q], trong đó những giá trị đã được định nghĩa trước đó.
function PQ = paddedsize(AB, CD, PARAM)
% PADDEDSIZE tính toán kích thước cần thiết thêm vào cho FFT dựa theo hàm lọc.
% PQ = PADDEDSIZE(AB), trong đó AB là vecto gồm hai phần tử, tính PQ = 2*AB.
%
% PQ = PADDEDSIZE(AB, ‘PWR2’) tính toán vecto PQ để PQ(1) = PQ(2) = 2^nextpow2(2*m), trong đó m là max(AB).
%
% PQ = PADDEDSIZE(AB, CD), trong đó AB và CD là vecto gồm hai phần tử, tính toán PQ. Những phần tử của PQ là số nguyên chẵn nhỏ nhất lớn hơn hoặc bằng AB + CD -1
%
% PQ = PADDEDSIZE(AB,CD,’PWR2’) tính toán vecto PQ để PQ(1) = PQ(2) =
2^nextpow2(2*m), trong đó m là max([AB CD]).
if nargin = = 1
PQ = 2*AB;
esleif nargin = = 2 & -ischar(CD)
PQ = AB + CD -1;
PQ = 2* ceil(PQ/2);
esleif nargin = = 2
m = max(AB); % chiều cực đại.
% Tìm kiếm số mũ gần nhấn lớn hơn 2 lần m.
P = 2^nextpow2(2*m);
PQ = [P, P];
esleif nargin = = 3
m = max([AB CD]); % chiều cực đại
P = 2^nextpow2(2*m);
PQ = [P, P];
else
error(‘đối số nhập vào là không đúng.’)
end
Với PQ tính được từ hàm paddedsize, chúng ta sử dụng lệnh fft2 để tính như sau:
F = fft2(f, PQ(1), PQ(2))
Câu lệnh thực hiện đơn giản là thêm đủ số 0 cho f để ảnh có kích thước là PQ(1)*PQ(2), và sau đó tính FFT như đã mô tả trước đó. Chú ý rằng khi chúng ta sử thủ tục trên thì hàm lọc trong miền tần số cũng phải có kích thước là PQ(1)*PQ(2).
Ảnh f trong hình 1.5(a) được sử dụng trong ví dụ để minh họa hàm lọc sử dụng và không sử dụng thủ tục thêm vào. Trong phần bàn luận dưới đây, chúng ta sử dụng lpfilter để sinh ra lọc băng thông thấp Gauss [tương tự hình 1.4(b)] với giá trị đặc biệt sigma(sig). Hàm này được bàn luận chi tiết trong mụch 1.5.2, nhưng cú pháp của nó dễ hiểu vì vậy chúng ta chỉ sử dụng nó ở đây còn phần giải thích để lại sau.
Những lệnh sau biểu diễn công việc lọc mà không dùng thủ tục thêm:
>> [M, N] = size(f);
>> F = fft2(f);
>> sig = 10;
>> H = lpfilter(‘gaussian’, M, N, sig);
>> G = H.*F;
>> g = real(ifft2(G));
>> imshow(g, [ ])
Trong hình 1.5(b) chỉ ra ảnh g. Như mong muốn, ảnh bị mờ đi nhưng cạnh nằm ngang thì không. Lý do có thể được giải thích trong hình 1.6(a). Hình 1.6(a) cho thấy chu kỳ trong tính toán DFT. Nhưng đường thẳng màu trắng nhỏ trong ảnh là cho dễ nhìn, nó không phải là một phần dữ liệu. Đường thẳng đứt nét dùng để chỉ định một cách tùy ý ảnh có kích thước M*N được xử lý bởi fft2. Tùy ý ở đây có nghĩa là ta chọn một chu kỳ bất kỳ trong dãy chu kỳ vô hạn ở trên. Tưởng tượng về sự hoạt động của bộ lọc mờ trong chuỗi chu kỳ vô hạn, chúng ta có thể nhận thấy rõ ràng là bộ lọc đi từ đỉnh của hình chữ nhật đứt nét, nó sẽ bao quanh từng vùng của ảnh cho tới phần cuối cùng của chu kỳ bên phải nó. Vì vậy vùng sáng và vùng tối khi ở chung một bộ lọc thì cho ra kết quả là độ xám trung bình, làm mờ ảnh đầu ra. Hình 1.5(b) cho ra kết quả đó. Nói cách khác, khi bộ lọc ở vùng sáng của vùng hình chư nhật đứt nét, nó sẽ bắt gặp vùng đồng nhất, khi trung bình của một vùng không đổi(vùng có các giá trị bằng nhau) thì vùng đó không bị làm mờ. Những phần khác của hình 1.5(b) được giải thích tương tự theo cách đó.
Hình 1.5
a b c
(a) Mẫu ảnh có kích thước 256*256
(b) Ảnh được lọc bằng bộ lọc băng thông thấp trong miền tần số mà không dùng kỹ thuật thêm vào.
(c) Ảnh được lọc bằng bộ lọc băng thông thấp với kỹ thuật thêm vào.
Hình 1.6:
a
b
Chuỗi chu kỳ vô hạn của ảnh trong hình 1.5(a). Vùng chữ nhật có nét đứt biểu diễn dữ liệu được xử lý bằng fft2.
Cùng một chu kỳ như trên nhưng sau khi được thêm 0. Đường thẳng mỏng màu trắng trong cả hai ảnh dùng để cho ta quan sát cho tiện, nó không phải là một phần dữ liệu
Ta xem bộ lọc có áp dụng kỹ thuật:
>> PQ = paddedsize(size(f));
>> Fp = fft2(f, PQ(1), PQ(2), 2*sig);
>> Hp = lpfilter(‘gaussian’, PQ(1), PQ(2), 2*sig);
>> Gp = Hp*Fp;
>> gp = real(ifft2(Gp));
>> gpc = gp(1:size(f,1), 1:size(f,2));
>> imshow(gp, [ ]);
Ở đó chúng ta dùng 2*sig bởi vì kích thước bộ lọc này gấp 2 lần bộ lọc không dùng kỹ thuật thêm.
Hình 1.7 chỉ ra kết quả đầy đủ khi dùng kỹ thuật thêm được trích xuất từ hàm ifft2 sau khi lọc. Ảnh này có kích thước 512*512 pixel.
1.3.2 Những bước cơ bản trong bộ lọc DFT:
Trong bàn luận trước đó chúng ta có thể tổng kết thành những hàm sau trong Matlab, ở đó f là ảnh được lọc, g là kết quả, với giả sử hàm lọc H(u,v) có cùng kích thước với ảnh được thêm số 0:
Dùng hàm paddedsize tính tham số thêm vào:
PQ = paddedsize(size(f));
Tính toán biến đổi Fourier:
F = fft2(f, PQ(1), PQ(2));
3. Sinh ra hàm H có kích thước PQ(1)*PQ(2) dùng bất kỳ phương pháp này trong phần còn lại của chương này. Hàm lọc được định dạng như trong hình 1.4(b). Nếu nó chưa được trung tâm hóa thì ta phải dùng H = fftshift(H) trước khi dùng bộ lọc.
4. Tiến hành nhân:
G = H.*F;
5. Tính phần thực của biến đổi ngược FFT của G:
g = real(ifft2(G));
6. Trả kích thước hình chữ nhật bên trái phía trên (hình 1.7) về màn kích thước gốc:
g = g(1:size(f,1), 1:size(f,2));
Phương thức lọc được tổng kết trong hình 1.8(a). Những bước xử lý (như việc xác định kích thước ảnh, tính tham số thêm vào, sinh ra bộ lọc)có thể được gói gọn trong phương thức. Hậu xử lý bao gồm việc tính toán phần thực dựa theo kết quả, thu lại kích thước ảnh gốc, chuyển về lớp uint8 hay uint16 để lưu trữ.
Hình 1.8
f(x,y)
Ảnh đầu ra
g(x,y)
Ảnh lọc
Tiền xử lý
Biến đổi ngược Fourier
Hàm lọc H(u,v)
Biến đổi Fourier
Tiền xử lý
Hàm lọc H(u,v) trong hình 1.8 nhân cả phần thực và phần ảo của F(u,v). Nếu H(u,v) là thực, pha của kết quả không đổi, chúng ta có thể kiểm chứng trong phương trình pha (hình 1.1), khi tiến hành nhân hàm lọc cho cả phần thực và phần ảo của ảnh thì tan của tỷ số giữa chúng không đổi, dẫn đến pha không đổi. Việc lọc tính theo cách trên gọi là hàm lọc hiệu chỉnh pha zero(zero phase shift filters). Nó chỉ là một loại trong hàm lọc tuyến tính trong chương này.
Chúng ta biết rằng trong lý thuyết tuyến tính, dưới một điều kiện vừa phải, việc thêm vào hệ thống tuyến tính một xung lực(impulse) là đặc trung cho hệ thống bởi khi chúng ta làm việc với một dữ liệu rời rạc hữu hạn như trong cuốn sách này, kết quả trả ra cũng tuyến tính. Nếu hệ thống tuyến tính chỉ là hàm lọc không gian, chúng ta có thể hoàn toàn xác định bộ lọc nhờ quan sát kết quả. Hàm lọc như thế gọi là bộ lọc đáp lại hữu hạn(finite impulse response) FIR. Tất cả hàmm lọc tuyến tính trong sách này đều là hàm lọc FIR.
1.3.3 Hàm M cho phép lọc trong miền tần số:
Chuỗi bước lọc miêu tả ở trên được dùng xuyên suốt trong chương này và một phần kế trong chương tiếp theo, vì vậy chúng ta cần thiết có sẵn hàm M chấp nhận ảnh đầu vào và hàm lọc, thực hiện tất cả chi tiết lọc bên trong, đầu ra của bộ lọc, thu lại ảnh gốc.Hàm đó được thực hiện như sau:
function g = dftfilt(f,H)
%DFTFILT biểu diễn phép lọc trong miền không gian.
% G = DFTFILT (F,H) lọc F trong miền tần số sử dụng hàm lọc H. Đầu ra, G, là ảnh được lọc, có cùng kích thước với F. DFTFILT tự động thêm 0 vào f cho nó cùng kích thước với H. Hàm paddedsize có thể được sử dụng để xác định kích thước thích hợp của H.
%
% DFTFILT giả sử rằng F và H là số thực, chưa được trọng tâm hóa, hàm lọc đối xứng vòng.
% Tính fft của đầu vào đã được thêm 0.
F = fft2(f, size(H,1), size(H,2));
%Biểu diễn bộ lọc
g = real(ifft2(H.*F));
% thu ảnh lại kích thước cũ
g = g(1:size(f,1), 1:size(f,2));
Kỹ thuật tạo ra những miền tần số được bàn luận trong phần dưới đây.
1.4 Tính toán hàm lọc trong miền tần số từ hàm lọc trong miền không gian:
Trong trường hợp tổng quát thì hàm lọc trong miền không gian tính toán hiệu quả hơn miền tần số khi miền tần số nhỏ. Định nghĩa nhỏ là một câu hỏi phức tạp phụ thuộc vào những nhân tố như máy móc và thuật toán, kích thước truyền tải, dữ liệu được xử lý tốt như thế nào, hay những nhân tố khác vượt ra khỏi những lĩnh vực trong cuốn sách này. Một so sánh được thực hiện bởi Brigham [1988] sử dụng hàm một chiều để chứng minh rằng hàm lọc sử dụng thuật toán fft có thể nhanh hơn cài đặt trong miền không gian khi hàm thực hiện lệnh trên 32 điểm. Vì vậy chúng ta cần biết cách chuyển hàm lọc trong miền không gian tương đương với hàm lọc trong miền tần số để có một so sánh có ý nghĩa giữa hai cách tiếp cận.
Một cách tiếp cận hiển nhiên cho việc tạo ra hàm lọc trong miền tần số H tương ứng với miền không gian h cho trước là đặt H = fft2(h, PQ(1), PQ(2)), trong đó kết quả của vecto PQ phụ thuộc vào kích thước của ảnh mà chúng ta muốn lọc, như đã bàn trong phần trước. Tuy nhiên, chúng ta quan tâm trong phần này gồm hai chủ đề: (1) làm thế nào chuyển hàm lọc miền không gian thành hàm lọc miền tần số tương đương; (2) làm thế nào để so sánh kết quả hàm lọc miền không gian sử dụng hàm imfilter, và hàm lọc miền tần số sử dụng kỷ thuật như đã đề cập ở trên. Bởi imfilter sử dụng sự tương quan giữa gốc tọa độ của hàm lọc được xem như là trung tâm của nó, thì cần một lượng dữ liệu tiền xử lý nhất định được yêu cầu để làm hai cách tiếp cận tương đương với nhau. Hôp công cụ cung cấp hàm freqz2, thực hiện chính xác điều này và xuất ra hàm lọc tương ứng trong miền không gian.
Hàm freqz2 tính toán tần số tương ứng của hàm lọc FIR, như đã đề cập ở phần cuối của mụch 1.3.2, ta chỉ xem xét hàm lọc tuyến tính trong sách này. Kết quả cho ra một hàm lọc mong muốn trong miền tần số. Cú pháp của lệnh như sau:
H = freqz2(h, R, C)
Trong đó h là hàm lọc không gian hai chiều và H là hàm lọc hai chiều trong miền tần số tương ứng. Ở đây, R là số dòng, C là số cột ta muốn ở H . Nếu freqz2 được viết mà không có đầu ra thì giá trị tuyệt đối của H được biểu diễn trong màn hình Matlab dưới dạng 3-D. Cơ chế xoay quanh việc dùng ham freqz2 dễ dàng hiểu trong ví dụ sau .
Ví dụ 1.2 : so sánh hàm quá trình lọc trong miền không gian và miền tần số.
Xét ảnh f có kích thước 600*600 như trong ảnh 1.9(a). Chúng ta tạo ra hàm lọc miền tần số H, tương ứng với hàm lọc sobel trong miền không gian,hàm lọc sobel tang cường cho những cạnh thẳng đứng . Chúng ta sau đó so sánh kết quả của việc lọc f trong miền không gian với mặt nạ sobel và kết quả biểu diễn trong miền tần số tương ứng. Trong thực nghiệm, với hàm lọc có kích thước nhỏ như Sobel sẽ được cài đặt trực tiếp trong miền không gian. Tuy nhiên ta dùng hàm lọc này cho much5 đích minh họa bởi vì hệ số của nó đơn giản và kết quả trực quan dễ so sánh. Những hàm lọc không gian lớn hơn cũng thực hiện theo cách đó.
Hình 1.9 :
a b
(a) ảnh xám
(b) phổ Fourier của nó.
Hình 1.9(b) là ảnh của phổ Fourier của f, nó được tính toán như sau :
>> F = fft2(f) ;
>> S = fftshift(log(1+ abs(F))) ;
>> S = gscale(S);
>> imshow(S);
Tiếp theo chúng ta tạo ra hàm lọc không gian dùng hàm fspecial:
h = fspecial(‘sobel’);
h =
1 0 -1
2 0 -2
1 0 -1
Để xem sơ đồ của hàm lọc không gian tương ứng chúng ta gõ lệnh:
>> freqz2(h)
Hình 1.10(a) cho thấy kết quả (kỹ thuật tính toán sơ đồ được bàn luận trong mụch 1.5.3), với cạnh bị giấu đi. Hàm lọc được tính toán như sau:
>> PQ = paddedsize(size(f));
>> H = freqz2(h, PQ(1), PQ(2));
>> H1 = ifftshift(H);
Trong đó hàm ifftshift cần thiết cho việc sắp xếp lại dữ liệu để gốc tọa độ nằm ở góc trên bên trái vùng tần số. Hình 1.10(b) chi ra sơ đồ của abs(H1). Hình 1.10(c) và (d) chỉ ra giá trị tuyệt đối của H và H1 trong ảnh được thực hiện bằng câu lệnh:
>> imshow(abs(H), [ ] );
>> figure, imshow(abs(H1), [ ] );
Hình 1.10
a b
c d
Tiếp theo chúng ta tạo ra ảnh lọc. Trong miền không gian ta sử dụng
>> gs = imfilter(double(f), h);
Thêm đường biên của ảnh là 0 như mặc định. Ảnh lọc trong miền tần số thu được bằng tính toán trong miền tần số được xử lý như sau:
gf = ftfilt(f, H1);
Hình 1.11(a) và (b) cho ra kết quả bằng câu lệnh:
>> imshow(gs, [ ]);
>> figure, imshow(gf, [ ] )
Tổng mức xám trong ảnh phụ thuộc vào yếu tố là gs và gf có giá trị âm hay không. Điều này làm ảnh hưởng tới giá trị trung bình của ảnh. Giá trị trung bình sau đó được định tỷ lệ trong lệnh imshow. Như miêu tả trong mụch 6.6.1 và 10.1.3 mặt nạ Sobel, h, sinh ra ở trên dùng để xác định cạnh thẳng đứng của ảnh dùng giá trị tuyệt đối của kết quả. Vì vậy chúng ta còn phải chỉ ra giá trị tuyệt đối của ảnh được tính toán. Hình 1.11(c) và (d) chỉ ra ảnh thu được bằng câu lệnh sau:
>> figure, imshow(abs(gs), [])
>> figure, imshow(abs(gf), [])
Cạnh có thể được nhìn thấy rõ hơn bằng việc thiết lập ảnh nhị phân theo ngưỡng:
>> figure, imshow(abs(gs) > 0.2*abs(max(gs(:))))
>> figure, imshow(abs(gf) > 0.2*abs(max(gf(:))))
Trong đó việc nhân với 0.2 được chọn tùy ý để chỉ ra rằng chỉ những cạnh có độ lớn lớn hơn 20% giá trị cực đại của gf và gs mới được biểu diễn. Hình 1.12(a) và (b) chỉ ra kết quả
Ảnh thu được sử dụng hàm lọctrong miền không gian và miền tần số cho tất cả cách mụch tiêu đồng nhất, chúng ta tiến hành tính toán sự khác biệt như sau:
>> d = abs(gs - gf);
Cực đại và cực tiểu của khác biệt:
>> max(d(:))
ans =
5.4015e-012
>> min(d(:))
ans =
0
1.5 Sinh ra hàm lọc một cách trực tiếp từ miền không gian
Trong phần này chúng ta minh họa cách cài đặt hàm lọc trực tiếp trong miền không gian. Chúng ta quan tâm tới những hàm lọc đối xứng vòng được đặc trưng như những hàm khác nhau của khoảng cách từ gốc tọa độ của biến đổi. Những hàm M được phát triển để cài đặt những hàm lọc này là cơ sở để dễ dàng mở rộng và cài đặt những hàm khác trong cùng một ngữ cảnh. Chúng ta bắt đầu bằng việc cài đặt một vài hàm lọc mịn (băng thông thấp). Sau đó chúng ta chỉ ra làm thế nào để dùng wireframe và khả năng thiết lập bề mặt của Matlab, những chức năng này làm tăng khả năng hình dung về hàm lọc. Chúng ta cũng có những bàn luận vắn tắt về hàm lọc sắc nét(băng thông cao).
1.5.1 Thiết lập dãy lưới trong việc cài đặt hàm lọc trong miền tần số:
Vấn đề chủ yếu của hàm M trong phần bàn luận dưới đây là tính toán hàm khoảng cách từ bất kỳ điểm nào tới một điểm đặc biệt trong hình chữ nhật tần số. Bởi vì việc tính toán FFT trong Matlab giả sử rằng gốc tọa độ nằm ở góc trên bên trái của hình chữ nhật tần số, mà khoảng cách tính toán của chúng ta lại liên quan tới điểm này. Dữ liệu có thể được sắp xếp lại cho mụch đích quan sát bằng cách dùng hàm fftshift (chuyển giá trị của gốc tọa độ về trung tâm của hình chữ nhật tần số).
Những hàm M sau đây chúng tôi gọi là dftuv, cung cấp dãy lưới cho việc tính toán khoảng cách và những ứng dụng tương tự. (xem thêm trong mụch 2.10.1 dể hiểu về giải thích của hàm meshgrid được sử dụng trong đoạn code này). Dãy lưới được sinh ra bởi dftuv theo một thứ tự phục vụ cho việc xử lý fft2 và ifft2, vì vậy không có sự yêu cầu sắp xếp lại dữ liệu trong trường hợp này.
function [U,V] = dftuv(M,N)
%DFTUV tính toán ma trận lưới tần số
% [U,V] = DFTUV(M,N) tính toán ma trận lưới tần số U và V, U và V rất hữu dụng cho việc tính toán hàm lọc trong miền tần số, hàm lọc này sau đó có thể dùng trong hàm DFTFILT. U và V đều có kích thước M*N.
%
% Thiết lập trật tự các biến:
u = 0: (M-1);
v = 0: (N-1);
% Tính toán chỉ số cho việc dùng lưới.
idx – find(u>M/2);
u(idx) = u(idx) –M;
idy = find(v> N/2);
v(idy) = v(idy) –N;
% Tính toán dãy lưới
[V,U] = meshgrid(v,u);
Như minh họa trên, nhưng câu lệnh sau đây tính toán bình phương khoảng cách từ mỗi điểm trong hình chữ nhật kích thước 8*5 tới gốc tọa độ hình chữ nhật :
>> [U,V] = dftuv(8,5) ;
>> D = U.^2 + V.^2 ;
D =
0 1 4 4 1
1 2 5 5 2
4 5 8 8 5
9 10 13 13 10
16 17 20 20 17
9 10 13 13 10
4 5 8 8 5
1 2 5 5 2
Chú ý rằng khoảng cách là 0 tại góc trên bên trái và những khoảng cách lớn hơn là tại tâm của hình chữ nhật tần số, theo những định dạng cở sở trong hình 1.2(a). Chúng ta có thể dùng hàm fftshift để tính toán khoảng cách liên quan tới tâm của hình chữ tần số.
>> fftshift(D)
ans =
20 17 16 17 20
13 10 9 10 13
8 5 4 5 8
5 2 1 2 5
4 1 0 1 4
5 2 1 2 5
8 5 4 5 8
13 10 9 10 13
Khoảng cách bây giờ là 0 tại tọa độ (5,3), và dãy đối xứn qua điểm này.
1.5.2 Hàm lọc băng thông thấp trong miền tần số :
Hàm lọ băng thông thấp theo vành (ILPF) :
1 nếu D(u,v) <= Do
H(u,v) =
0 nếu D(u,v) > Do
Trong đó Do là một giá trị đặc biệt không âm và D(u,v) là khoảng cách từ điểm (u,v) tới tâm hàm lọc. Quỹ tích những điểm D(u,v) = Do là hình tròn. Nhớ rằng hàm lọc H nhân với biến đổi Fourier của ảnh, chúng ta có thể thấy rằng hàm lọc theo vành bỏ đi (bằng cách nhân với 0) tất cả những thành phần của F bên ngoài hình tròn và không đổi tất cả những thành phần nằm trong hình tròn (nhân với 1). Mặc dù hàm lọc này không nhận ra được bằng cách dùng thành phần điện tử, nhưng nó có thể được thiết lập nhờ hàm chuyển. Tính chất của hàm lọc theo vành thường hữu dụng trong việc giải thích những khái niệm như vòng tròn khép kín.
Hàm lọc băng thấp butterworth (BLPF) thứ tự n, với khoảng cách Do tính từ gốc tọa độ, có dạng hàm chuyển như sau :
Không giống như ILPF, hàm chuyển BLPF không có bị gián đoạn tại Do. Với bộ lọc có hàm chuyển mịn, thông thường chúng ta định nghĩa quỹ tích ngưỡng tần số mà tại đó H(u,v) bị giảm tới một giá trị thập phân đặt biệt dựa theo giá trị cực đại của nó. Trong phương trình trước thì H(u,v) = 0.5 (giảm 50% từ giá trị cực đại của nó là 1) khi ấy D(u,v) = Do.
Hàm chuyển của bộ lọc băng thông thấp Gauss (GLPF) được cho bởi :
Trong đó là độ lệch chuẩn. Bằng cách đặt = Do, chúng ta thu được biểu thức sau :
Khi D(u,v) = Do hàm lọc giảm xuống 0.607 giá trị cực đại của nó là 1.
Ví dụ 1.4 : hàm lọc băng thông thấp
Như minh họa chúng ta áp dụng bộ lọc băng thông thấp Gauss cho ảnh có kích thước 500*500 pixel, trong ảnh 1.13(a). Chúng ta dùng giá trị Do bằng 5% chiều dài ảnh đã áp dụng kỹ thuật thêm vào.Với những bàn luận trong mụch 1.3.2 ta có
Hình 1.13 :
a b
c d
(a) Ảnh gốc
(b) Bộ lọc băng thông thấp Gauss
(c) Phổ của (a)
(d) Ảnh được xử lý.
Chúng ta có thể quan sát hàm lọc như trong hình 1.13(b) bằng cách :
>> figure, imshow(fftshift(H), [])
Tương tự phổ cũng được biểu diễn như sau :
>> figure, imshow(log(1+abs(fftshift(F))), [ ]);
Cuối cùng là ảnh 1.13(d):
>> figure, imshow(g, [ ])
Như mong muốn ảnh đã được làm mờ đi.
Hàm tiếp theo tạo ra hàm chuyển đổi của tất cả các bộ lọc băng thông thấp được bàn đến trong chương này.
function [H,D] = lpfilter(type, M, N, D0, n)
% lpfilter tính toán hàm lọc băng thông thấp trong miền tần số.
% H = lpfilter(type, M, N, D0, n) tạo hàm chuyển của bộ lọc băng thông thấp, H, của một loại đặc trưng nào đó (kích thước M*N). Để quan sát bộ lọc như ảnh hay sơ đồ chúng ta có thể tiến hành trung tâm hóa H bằng cách H = fftshift(H)
%
% Những giá trị của TYPE, D0, n có thể nhận là:
% ‘ideal’ hàm lọc theo vành với ngưỡng tần số là D0, n không cần được cung cấp, D0 phải dương.
% ‘btw’ hàm lọc băng thông thấp butterworth thứ tự n, và ngưỡng D0. Giá trị mặc định cho n là 1.0. D0 phải dương
% ‘gaussian’ hàm lọc băng thông thấp gauss với ngưỡng (độ lệch chuẩn) D0. n không cần được cung cấp. D0 phải dương.
%
% Dùng hàm dftuv để thiết lập dãy lưới cần thiết cho tính toán khoảng cách theo yêu cầu.
[U,V] = dftuv(M, N);
% Tính khoảng cách D(u,v)
D = sqrt(U^.2 + V^.2);
% Bắt đầu với việc tính toán hàm lọc.
switch type
case ‘ideal’
H = double(D<= D0);
case ‘btw’
if nargin = = 4
n = 1;
end
H = 1./(1 + (D./D0).^(2*n));
case ‘gaussian’
H = exp(-(D.^2));
otherwise
error(‘loại hàm lọc không tồn tại’);
end
Hàm lpfilter được dùng trong mụch 1.6 như là cơ sở cho việc tạo ra hàm lọc băng thông cao.
1.5.3 Wireframe và thiết lập bề mặt.:
Sơ đồ của hàm một biến được giới thiệu trong mụch 3.3.1 Trong phần bàn luận tiếp theo, chúng tôi giới thiệu wireframe và sơ đồ mặt phảng 3-D, rất hữu ích trong việc quan sát hàm lọc chuyển đổi hai chiều. Cách dễ nhất để vẽ sơ đồ wireframe của hàm 2 chiều cho trước, H, là dùng hàm mesh:
mesh(H)
Hàm này vẽ một wireframe cho x = 1:M và y = 1:N , trong đó [M, N] = size(H). Wireframe thường không được chấp nhận trong trường hợp M và N lớn và chúng ta vẽ mọi điểm thứ k theo cú pháp sau:
mesh(H(1:k:end, 1:k:end))
Dưa quy luật đó, chỉ có từ 40 đến 60 đoạn con được chia dọc theo mỗi trục cung cấp sự cân bằng giữa độ phân giải và hình thức xuất hiện của wireframe.
Trong Matlab có thêm chức năng thực thi hàm mesh theo màu dựa vào câu lệnh
colormap([0 0 0])
Tạo ra một wireframe màu đen (chúng ta sẽ bàn luận hàm colormap trong chương 6). Trong Matlab người ta cũng bỏ qua lưới và 2 trục tọa độ trong quá trình thực thi hàm mesh. Ta có thể dùng những lệnh sau:
grid off
axis off
Ta có thể kích hoạt hai lưới và 2 trục tọa độ trở lại bằng cách thay thế off là on trong hai câu lệnh. Cuối cùng việc xác định vị trí của người quan sát được điều khiển bằng câu lệnh view có cú pháp như sau:
view(az, el)
Như trong hình 1.14, az, el là góc phương vị và góc nâng (tính theo độ). Những mũi tên chỉ theo hướng dương. Giá trị mặc định là az = -37, el = 30, tại đó vị trí người quan sát đặt tại góc phần tư được định nghĩa bởi trục –x, -y, và nhìn vào góc phần tư được định nghĩ bởi trục x,y như trong hình 1.14.
Để xác định vị trí hình học hiện tại, chúng ta dùng:
>> [az, el] = view;
Để thiết lập góc nhìn tới giá trị mặc định chúng ta dùng hàm
>> view(3)
Góc nhìn có thể được hiệu chỉnh theo tương tác bằng cách click lên nút Rotate 3D trong hộp công cụ của cửa sổ, sau đó ta tiến hành click và kéo hình ảnh trong cửa sổ giao diện.
Như trong bàn luận của chương 6, chúng ta có thể cụ thể hóa vị trí người nhìn theo hệ tọa độ Decas (x,y,z), củng tương tự như khi ta làm việc với dữ liệu RGB. Tuy nhiên, Với mụch đích quan sát tổng quát, phương pháp chủ yếu xoay quanh 2 tham số và mang tính trực giác nhiều hơn.
Xét bộ lọc băng thông thấp Gauss tương tự như trong ví dụ 1.4:
>> H = fftshift(lpfilter(‘gaussian’,500, 500,50));
Hình 1.15(a) cho thấy wireframe giới thiệu trong lệnh:
>> mesh(H(1:10:500, 1:10:500))
>> axis ([0 50 0 50 0 1])
Trong đó lệnh axis đã được miêu tả trong mụch 3.3.1, ngoại trừ việc nó chứa sự sắp xếp thứ 3 dành cho trục z.
Hình 1.15
a b
c d
(a) Sơ đồ thu được nhờ hàm mesh
(b) Trục và lưới bị mất
(c) Sự khác biệt liên quan tới góc nhìn nhờ dùng hàm view
(d) Góc nhìn khác
Như đã chú ý trong phần trước của chương này, wireframe dùng hàm mặc định, chuyển từ màu xanh ở phía dưới dần thành màu đỏ lên tới đỉnh. Chúng ta có thể chuyển sơ đồ thành màu đen và bỏ đi trục tọa độ và lưới bằng cách:
>> colormap([0 0 0])
>> axis off
>> grid off
Hình 1.15(a) chỉ ra kết quả, hình 1.15(c) được thực hiện bằng câu lệnh:
>> view(-25,30)
Nó giúp ta quan sát hơi lệch qua bên phải mộ tí và bỏ luôn hai trục tọa độ. Cuối cùng hình 1.15(d) chỉ ra kết quả khi ta đặt góc nâng là 0:
>> view(-25,0)
Ví dụ này đã chỉ ra sức mạnh của việc vẽ sơ đồ của đon giản mesh.
Đôi khi chúng ta muốn vẽ một hàm như một mặt phẳng thay vì wireframe . Hàm surf có thể giúp ta làm điều đó với cú pháp:
surf(H)
Hàm này tạo ra môt sơ đồ như hàm mesh, ngoại trừ một việc rằng 4 bên của trong sơ đồ của hàm mesh được tô đầy màu (cò gọi là đổ bóng trên bề mặt). Để chuyển màu thành màu xám chúng ta dùng lệnh:
colormap(gray)
Các hàm axis, grid, và view tương tư như cách dùng hàm mesh. Ví dụ, trong hình 1.16(a) là kết quả của chuỗi câu lệnh:
>> H = fftshift(lpfilter(‘gaussian’, 500, 500, 50));
>> surf(H(1:10:500, 1:10:500))
>> axis([0 50 0 50 0 1])
>> colormap(gray)
>> grid off; axis off
Đổ bóng bề mặt có thể được làm mịn và các đường thẳng lưới có thể bỏ đi bằng nội suy bằng hàm lệnh
shading interp
Đánh câu lệnh này tại dấu nháy tạo ra hình 1.16(b)
Hình 1.16
a b
(a) Sơ đồ thu được bằng hàm surf
(b) Kết quả dùng hàm shading interp.
Khi mụch tiêu là vẽ sơ đồ hàm giải tích hai biến, chúng ta dùng hàm meshgrid để tạo ra giá trị tọa độ từ đó chúng ta tạo ra một ma trận rời rạc có thể dùng trong hàm mesh hay surf. Ví dụ, với hàm sau:
Từ -2 tới 2 và độ tăng là 0.1cho cả x và y, chúng ta viết như sau:
>> [Y,X] = meshgrid(-2:0.1:2, -2:0.1:2);
>> Z = X*.exp(-X^2-Y.^2);
Sau đó dùng hàm mesh(Z) và surf(Z) như trước.
.
1.6 Làm sắt nét hóa trong hàm lọc miền tần số:
Hàm lọc băng thông thấp làm mờ ảnh thì trái lại hàm lọc băng thông cao làm sắt nét ảnh hơn bằng cách bỏ đi những tần số thấp và giử lại những tần số cao của biến đổi Fourier. Trong phần này chúng ta xem xét một vài cách tiếp cận trong lọc băng thông cao.
1.6.1 Cơ bản cảu việc lọc băng thông cao:
Cho hàm lọc Hlp(u,v) của hàm lọc băng thông thấp, chúng ta tính hàm chuyển đổi của lọc băng thông cao bằng mối quan hệ đơn giản:
Hhp(u,v) = 1- Hlp(u,v)
Vì vậy hàm lpfilter phát triển trong phần trước có thể dùng như cơ sở cho hàm bộ băng thông cao trong phần này như sau:
function H = hpfilter(type, M, N, D0, n)
% hpfilter tính hàm lọc băng thông cao trong miền tần số.
% H = hpfilter(type, M, N, D0, n) tạo hàm chuyển cho một loại bộ lọc băng thông cao cụ thể là H có kích thước M*N
% Những giá trị mà TYPE có thể nhận là:
%
% ‘ideal’ hàm lọc băng thông cao theo vành với ngưỡng tần số D0, n không cần được cung cấp. D0 phải dương.
%
% ‘btw’ hàm lọc băng thông cao butterworth thứ tự n, và ngưỡng D0. Giá trị mặc định cho n là 1.0.D0 phải dương.
%
% ‘gaussian’ hàm lọc băng thông cao gaussian với độ lệch chuẩn D0. n không cần được cung cấp. D0 phải dương.
%
% Hàm chuyển Hhp của bộ lọc băng thông cao là 1 – Hlp, trong đó Hlp là hàm chuyển %tương ứng của bộ lọc băng thông thấp. Vì vậy, chúng ta cần dùng hàm lpfilter để tạo ra %bộ lọc băng thông cao.
if nargin = = 4
n = 1; % giá trị mặc định của n
end
% tạo ra hàm lọc băng thông cao
Hlp = lpfilter(type, M, N, D0, n);
H = 1 – Hlp;
Ví dụ 1.6:
Trong hình 1.17 chỉ ra sơ đồ và ảnh của bộ lọc băng thông cao theo vành, lọc butterworth và Gauss. Sơ đồ trong hình 1.17(a) tạo ra theo câu lệnh:
>> H = fftshift(hpfilter(‘ideal’, 500, 500, 50));
>> mesh(H(1:10:500, 1:10:500));
>> axis([0 50 0 50 0 1])
Hình 1.17
a b c bộ lọc băng thông cao theo vành, butterworth, Gauss
d e f hình ảnh của bộ lọc tương ứng
>> colormap([0 0 0])
>> axis off
>> grid off
Hình ảnh tương ứng trong 1.17(d) được tạo ra bằng câu lệnh
>> figure, imshow(H, [])
Trong đó đường biên màu đen được đặt chồng lên trong ảnh để phát họa biên của nó. Những bộ lọc còn lại cũng được làm tương tự. Trong ví dụ này thì bộ lọc băng thông cao butterworth có thứ tự là 2.
Hình 1.18(a) cũng kiểm tra trên cùng một mẫu như trong ví dụ 1.13(a). Hình 1.18(b) là kết quả của việc áp dụng bộ lọc băng thông cao Gauss và được tính toán dùng những câu lệnh sau:
>> PQ = paddedsize(size(f));
>> D0 = 0.05*PQ(1);
>> H = hpfilter(‘gaussian’, PQ(1), PQ(2), D0);
>> g = dftfilt(f,H);
>> figure, imshow(g, [])
Hình 1.18:
a b
(a) ảnh gốc
(b) Kết quả của lọc băng thông cao với Gauss.
Như trong hình 1.18(b), cạnh và những chỗ chuyển sắc nét trong ảnh được tăng cường. Tuy nhiên, bởi vì giá trị trung bình của ảnh là F(0,0), hàm lọc băng thông cao có giá trị ở gốc tọa độ của biến đổi Fourier khác 0 nên ảnh sẽ bị mất rất nhiều cảnh quan trọng so với tổng thể trong ảnh gốc. Vấn đề này được giải quyết ở mụch tiếp theo.
1.6.2 Sự nhấn mạnh tần số trong bộ lọc tần số cao.
Như đã đề cập trong mụch 1.7, bộ lọc băng thông cao có giá trị dc khác 0. Vì vậy chúng ta phải giảm giá trị trung bình của ảnh tới 0. Một phương án tiếp cận là thêm vào phần bù cho bộ lọc băng thông cao.Khi phần bù kết hợp với việc nhân bộ lọc bởi một hằng số lớn hơn 1, phương pháp này gọi là sự nhấn mạnh trong bộ lọc tần số cao. Việc nhân như vậy cũng làm gia tăng tần số thấp, nhưng tần số thấp ảnh hưởng tới sự tăng cường ít hơn tần số cao, khi ấy cần có sự so sánh giữa phần bù và hằng số nhân tùy theo mụch đích chọn lựa. Phương thức nhấn mạnh tần số cao có hàm chuyển là
Hhfe(u,v) = a + b*Hhp(u,v)
Trong đó a là phần bù, b là hằng số nhân, và Hhp(u,v) là hàm chuyển của bộ lọc băng thông cao.
Ví dụ 1.8:
Hình 1.19(a) chỉ ra hình chụp tia X của ngực, f, người chụp tia X không thể nhìn bằng thấu kính bởi vì ảnh kết quả thường hơi mờ. Mụch tiêu của ví dụ này là làm sắt nét hóa hình 1.19(a). Bởi vì màu ảnh xám trong trường hợp cụ thể này có xung hướng là màu tối, vì vậy chúng tôi muốn dùng ví dụ này để cho thấy làm thế nào xử lý trong miền không gian có thể sử dụng để bổ xung cho bộ lọc miền tần số.
Hình 1.19
a b
c d
(a) Ảnh gốc
(b) Kết quả bộ loc băng thông cao.
(c) Kết quả của nhấn mạnh tần số cao.
(d) Ảnh (c) sau khi cân bằng histogram.
Hình 1.19(b) chỉ ra kết quả của việc lọc hình 1.19(a) với bộ lọc băng thông cao butterworth thứ tự 2, giá trị D0 bằng 5% trục thẳng đứng của ảnh đã được thêm vào. Bộ lọc băng thông cao không quá mức nhạy cảm với giá trị của D0, do vậy nên bán kính không nhỏ tới mức những tần số gần gốc tọa độ của biến đổi bị bỏ qua. Như mong muốn thì kết quả, kết quả lọc đã bắt đầu xuất hiện những nét đặc biệt, nhưng nó lại chỉ ra những thành phần chính trong ảnh một cách mờ nhạt . Lợi thế của sự nhấn mạnh tần số trong bộ lọc tần số cao (trong trường hợp này là a = 0.5 và b = 2.0) được chỉ ra trong hình 1.19(c), trong đó tổng mức xám phụ thuộc vào những thành phần tần số yếu được giữ lại. Chuỗi câu lệnh được sử dụng để tạo ra ảnh được xử lý trong hình 1.19, trong đó f ký hiệu cho ảnh đầu vào [câu lệnh cuối cùng tạo ra hình 1.19(d)]:
>> PQ = paddedsize(size(f));
>> D0 = 0.05*PQ(1);
>> HBW = hpfilter(‘btw’, PQ(1), PQ(2), D0, 2);
>> H = 0.5 +2*HBW;
>> gbw = dftfilt(f, HBW);
>> gbw = gscale(gbw);
>> ghf = dftfilt(f,H);
>> ghf = gscale(ghf);
>> ghe = histeq(ghf, 256);
Như trong hình 1.19(d), cân bằng hitogram được tương thích với phương pháp lọc băng thông cao để mang lại hiệu quả tăng cường ảnh tốt hơn cho ảnh. Chú ý rằng sự rõ ràng của cấu trúc xương và những chi tiết khác không được quan sát thấy trong 3 ảnh trước đó. Ảnh được nâng cao chất lượng cuối cùng có mang một chút nhiễu, nhưng đó là điều hết sức thông thường trong ảnh tia X vì khi ấy mức xám được mở rộng. Kết quả thu được là sự kết hợp của việc nhấn mạnh tần số và cân bằng histogram tốt hơn nhiều khi áp dụng đơn lẻ một trong hai phương pháp.
2.Ghép hình:
Từ các hình ban đầu như sau:
Để có hình ảnh kết quả:
Chúng ta tiến hành các bước như sau:
Tiền hành cắt các đối tượng cần thao tác trên cả 3 hình bằng hàm roipoly trong Matlab.
Ghép hào quang vào ảnh núi phú sĩ.
Dùng hàm lọc băng thông thấp để làm trơn cạnh của ánh hào quang.
Ghép Ferguson vào ảnh núi phú sĩ và đè lên hào quang.
Ghép ngọn núi phú sĩ đã được cắt ra lúc đầu vào ảnh núi phú sĩ và đè lên Ferguson.
Những hình ảnh được cắt ra từ ảnh gốc bằng hàm roipoly:
Hào quang được trích ra từ ảnh gốc
Hình Ferguson trích ra từ ảnh gốc
Hình Núi Phú Sĩ trích ra từ ảnh gốc.
2.1 Thủ thuật ghép ảnh:
Sau khi chúng ta dùng hàm roipoly để cắt những đối tượng cần thiết trong từng ảnh, chúng ta sẽ thu được những ảnh nhị phân. Trong đó, đối tượng cần thiết là những pixel màu trắng. Những pixel màu trắng trong ảnh nhị phân được xem như là cờ hiệu để biết được đó có phải là pixel thuộc đối tượng cần thiết hay không.
Duyệt toàn bộ ảnh nhị phân, từ tọa độ của pixel trắng trong ảnh nhị phân chúng ta lần tới đúng tọa độ đó trong ảnh màu. Sau đó, chèn giá trị màu của của pixel đó vào vị trí chúng ta cần chèn trong ảnh nền. Tiến hành tương tự như vậy cho tất cả các pixel trắng trong ảnh nhị phân, chúng ta có được kết quả là đối tượng được chèn vào ảnh nền.
Một số điểm cần lưu ý trong kỹ thuật chèn hình ảnh là:
Xét tới thứ tự của ảnh được chèn: chúng ta phải xác định trước đối tượng nào đứng trước, đối tượng nào đứng sau để có thứ tự chèn hợp lý.
Chọn hình nền đủ lớn và có quang cảnh thoáng.
2.2 Hàm lọc băng thông thấp làm mịn ảnh:
Như trình bày trong phần lý thuyết thì hàm lọc băng thông thấp có tác dụng làm mờ các chi tiết góc cạnh. Có rất nhiều loại hàm lọc khác nhau, nhưng chúng tôi sử dụng hàm lọc Gauss vì tính hội tụ nhanh của nó.
Tuy nhiên, trong quá trình lọc cũng nên lưu ý đến việc đặt hàm lọc ở đâu và như thế nào để cho ra kết quả tối ưu nhất. Ví dụ, trong bài làm của nhóm chúng tôi, nếu tiến hành lọc ảnh xong mới ghép sẽ xảy ra hiện tượng ảnh bị nhòe toàn bộ như hình bên dưới:
Vì vậy sau khi ghép hào quang vào, chúng tôi tiến hành lọc ngay. Vì thực sự chỉ có một đối tượng là ánh hào quang cần được lọc mà thôi.
Kết luận:
Bài tập trên ứng dụng hai kỹ thuật cơ bản của xử lý ảnh là ghép hình và làm trơn cạnh của một vùng trong ảnh. Trong đó, kỹ thuật ghép hình dùng hàm ROIPOLY để cắt các đối tượng cần ghép ra, lưu vết đối tượng trong một ảnh nhị phân và ghép vào ảnh nền dựa vào vết đó. Kỹ thuật làm trơn cạnh trong ảnh dựa vào biến đổi dữ liệu ảnh thành dữ liệu trong miền tần số. Thông qua việc hiệu chỉnh các tín hiệu hình sin và khôi phục lại dữ liệu ảnh nhờ lý thuyết tích chập, chúng ta thu được ảnh cần hiệu chỉnh. Ngoài ra chúng ta cần có một chiến lược ghép ảnh hợp lý, xác định thứ tự ghép ảnh và thứ tự lọc để có thể có một bức ảnh theo ý muốn.
Tài liệu tham khảo:
[1] Gonzalez, Woods, Digital Image Processing Using Matlab (2nd).
[2] William K.Pratt, Digital Image Processing, Third Edition, Pixelsoft, Inc, LosAltos, California.