EDA365电子论坛网
标题:
求助：关于非局部均值滤波程序的问题
[打印本页]
作者:
qUzalq
时间:
2020-8-28 14:36
标题:
求助：关于非局部均值滤波程序的问题
% 沿行插值
( f. B Z0 |; D3 X0 g5 u# l
Sd = cumsum(diff,1); % 行叠加
) X( P! v8 S% u2 c
% 沿列插值
' z. ~1 a- F7 O# I
Sd = cumsum(Sd,2); % 列叠加
9 K6 p7 t9 ]0 g9 j; W% U' p% A
为什么要进行叠加？
0 Q. c, E* U+ l9 q4 {6 `7 u9 v {( o/ f
还有这一步也不是很明白，怎么就得到每个像素点的权值了
]" e6 \! E( ~! [
SqDist = Sd(PatchSizeHalf+1:end-PatchSizeHalf,PatchSizeHalf+1:end-PatchSizeHalf)+Sd(1:end-2*PatchSizeHalf,1:end-2*PatchSizeHalf)-Sd(1:end-2*PatchSizeHalf,PatchSizeHalf+1:end-PatchSizeHalf)-Sd(PatchSizeHalf+1:end-PatchSizeHalf,1:end-2*PatchSizeHalf);
+ m# q/ P4 G$ x
作者:
NingW
时间:
2020-8-28 15:59
路过
作者:
mutougeda
时间:
2020-8-28 16:00
分享这个源代码，应该有用
8 ~* g6 a3 q, A! r& c+ g U8 z9 v) Y' K
' E- G( \6 J5 b5 Y0 g
function DeNimg = Non_Local_Means(Nimg,PSH,WSH,Sigma)
. D, l1 c. p& b8 K9 R
% Non_Local_Means滤波器
2 G/ U" u E' e2 Y# n
%函数输入:
6 |2 Y! C, x( J, }
% Nimg: 输入的图像矩阵 + 带噪声的
0 N# {/ t8 ~6 `% F/ s! E# g
% PSH: 扩展窗尺寸大小
& G9 r0 s0 |+ V# y% N9 B3 s6 C
% WSH: 窗尺寸大小
: x/ @# i# y7 V& ~2 [
% Sigma：方差
1 X/ R+ s. C+ u+ q! P
%函数输出:
; U6 u& y$ n& I* K @) `' R
% DeNimg: 重构滤波图像
( `0 U' l+ Z: y* p1 e- p
( Q) q/ ]3 ?7 Y2 ?* P+ U
if ~isa(Nimg,'double')
. r6 a! o, [- @1 W
Nimg = double(Nimg)/255;
+ K0 Y0 m, m: E
end
, X8 ?8 x/ w* q4 a. Z' K8 Z7 T1 m
2 {- z3 d& d3 V
% 图像维数
, p+ {' e" \6 E' z3 M
[Height,Width] = size(Nimg);
! |& l, x. Z% F/ n' j! S
u = zeros(Height,Width); % 初始化去噪图像矩阵
) d0 ^5 l% _$ U" w
M = u; % 初始化权值矩阵
: `7 O* {# \( h4 a; f8 B
Z = M; % 初始化叠加权值 accumlated weights
4 e& W: k6 J6 v& _
% 避免边界效应
) C* ]/ v0 `' J( C- O8 X6 ^7 `
PP = padarray(Nimg,[PSH,PSH],'symmetric','both');
0 U) Z, j# ]- [; F) |
PW = padarray(Nimg,[WSH,WSH],'symmetric','both');
7 h4 K: [+ O* m0 ~. W+ o
% padarray使用
( m8 d) X' N, P
% A =
5 W- ^- k# c+ c) R/ v7 }
% 1 3 4
/ P1 M8 U. D* M4 D5 V) ?! U" o8 x
% 2 3 4
0 p- s7 k+ m3 i# e' u4 m2 a/ b5 z
% 3 4 5
3 s* |5 _% q; S; A$ ]
% B = padarray(A, 2 * [1 1], 0, 'both')
" M8 @3 W5 `0 X8 u2 d; L1 A$ q
% 0 0 0 0 0 0 0
0 e: j0 C6 |( Y, y. \; L9 b+ Q
% 0 0 0 0 0 0 0
0 ?4 f- b+ @& p5 B( K- Y- i
% 0 0 1 3 4 0 0
! s& d1 c0 l$ C
% 0 0 2 3 4 0 0
( ~8 ` a6 S8 X. e' t N* a
% 0 0 3 4 5 0 0
# T! }8 \2 H- s
% 0 0 0 0 0 0 0
{$ X8 D! ^/ ?3 E+ Z9 s( P. g
% 0 0 0 0 0 0 0
. r/ _9 y) T6 F) F% h9 N* |
% 主循环
( @. V. Y) Y: Y* B5 n! x! b& P
for dx = -WSH:WSH
& Y, d! Q9 \$ }
for dy = -WSH:WSH
& }1 [: z* b9 H- R4 k S8 x" z; `
if dx ~= 0 || dy ~= 0
. R" A. F% }- N9 ~& J% {, M
Sd = integral_img(PP,dx,dy); % 插值图像
' |$ ~- k3 |% N
% 获取对应像素点的平方差矩阵
3 c6 q; ~3 G3 n0 z
SDist = Sd(PSH+1:end-PSH,PSH+1:end-PSH)+Sd(1:end-2*PSH,1:end-2*PSH)-Sd(1:end-2*PSH,PSH+1:end-PSH)-Sd(PSH+1:end-PSH,1:end-2*PSH);
4 C2 P9 c' n& e- }# S
% 计算每一个像素点的权值
+ G) B8 D4 Y1 r
w = exp(-SDist/(2*Sigma^2));
7 v! r8 F2 c" m/ y& n
% 得到相应的噪声点
3 } f( K T5 v# z2 X% _
v = PW((WSH+1+dx)
WSH+dx+Height),(WSH+1+dy)
WSH+dy+Width));
, l/ ~+ `* a7 n9 U; V5 G
% 更新去噪图像矩阵
3 [2 L/ Y' a( j& v* s: H
u = u+w.*v;
9 s; S/ B) Y7 Z# K; O! W
% 更新权值去噪图像矩阵
/ a$ d6 E1 I, x) X# C9 V
M = max(M,w);
5 R) C0 L8 Y" B* X( k1 z6 Y6 @
% 更新叠加权值 accumlated weights
E0 C2 K# T# `2 s$ V
Z = Z+w;
; n, U7 q4 v% d3 u
end
# i: a- K9 E- M8 e+ ^
end
2 Y1 g2 B/ e3 W
end
0 d) N, I1 Q$ |8 D( q1 {: U+ w, w
% 重构图像
! E0 M ~9 f$ ~/ {$ b
f = 1;
' z% F: L0 s" K0 [) o
u = u+f*M.*Nimg;
( `4 ~( y: u, [6 Y/ c0 L
u = u./(Z+f*M);
& X, G' y8 R0 [! b
DeNimg = u; % 重构去噪图像
6 S i1 P; _' @0 W
; v3 o' @5 M( L: V; ^
function Sd = integral_img(v,dx,dy)
+ ]3 C4 z: W/ W" B3 y: d2 ^
% 根据平方差，插值图像
8 @7 `4 c, t( }* s3 R
% 变换计算：tx = vx+dx; ty = vy+dy
- ?: H( q8 U+ }) p( E! ?
t = img_Shift(v,dx,dy);
3 U. o7 C v9 h2 N* i0 G
% 平方差图像
+ N" ?2 r2 R. R2 c1 @$ ?7 _
diff = (v-t).^2;
# w) @, g- ]$ X9 s2 x
% 沿行插值
' C1 U- l. Y3 N$ W" k
Sd = cumsum(diff,1); % 行叠加
( s. k n# D9 t! F
% 沿列插值
! W0 |% |$ t# t7 n
Sd = cumsum(Sd,2); % 列叠加
* A& N- a4 T& q W H5 S6 V
0 V6 H" R% Y. c/ i
function t = img_Shift(v,dx,dy)
\0 ?0 x8 w3 j) D0 ~! O
% 在xy坐标系下，进行图像变换操作
- C+ {7 M( ?: O
t = zeros(size(v));
- }& N, F, B* P8 Y- ~1 l6 r# S! Y0 _
type = (dx>0)*2+(dy>0);
* u1 p7 U/ S7 }& C) I5 W
switch type
" B& T+ V' v1 I1 h! r" G
case 0 % dx<0,dy<0: 向右下方移动
d. k0 R; j# @9 u
t(-dx+1:end,-dy+1:end) = v(1:end+dx,1:end+dy);
. l; D: @" w& {+ e
case 1 % dx<0,dy>0: 向左下方移动
6 e, q% g( q0 U8 G4 P
t(-dx+1:end,1:end-dy) = v(1:end+dx,dy+1:end);
8 A( l. `9 O# K
case 2 % dx>0,dy<0: 向右上方移动
3 ~1 K* [6 r: J7 X9 | _% \6 W
t(1:end-dx,-dy+1:end) = v(dx+1:end,1:end+dy);
/ s8 J* e/ p9 k4 i: F
case 3 % dx>0,dy>0: 向左上方移动
# W" U8 S1 e9 _8 m" G) w
t(1:end-dx,1:end-dy) = v(dx+1:end,dy+1:end);
x3 S, X! n+ p1 j
end
) S- U8 p1 }" k, @5 z: Z
欢迎光临 EDA365电子论坛网 (https://bbs.eda365.com/)
Powered by Discuz! X3.2
WSH+dx+Height),(WSH+1+dy)WSH+dy+Width));, l/ ~+ `* a7 n9 U; V5 G