在 Matlab 中过滤包含 NaN 的图像？
Filter image that contains NaNs in Matlab?
我有一个代表一些数据的二维数组 (doubles)，里面有一堆 NaNs。数据的等高线图如下所示：
所有空白都是NaNs，灰色菱形供参考，填充轮廓显示我的数据的形状。当我使用 imfilt 过滤数据时，NaNs 会显着影响数据，所以我们最终得到这样的结果：
您可以看到支持集显着收缩。我不能使用它，因为它已经咀嚼了边缘上的一些更有趣的变化（出于我实验的特定原因，这些边缘很重要）。
是否有在 NaNs 岛内进行过滤的功能，它处理类似于矩形过滤 windows 边缘的边缘，而不是仅仅消除边缘？有点像 nanmean 函数，除了卷积图像？
这是我的过滤器代码：
filtWidth = 7;
imageFilter=fspecial('gaussian',filtWidth,filtSigma);
%convolve them
dataFiltered = imfilter(rfVals,imageFilter,'symmetric','conv');
以及绘制等高线图的代码：
figure
contourf(dataFiltered); hold on
plot([-850 0 850 0 -850], [0 850 0 -850 0], 'Color', [.7 .7 .7],'LineWidth', 1); %the square (limits are data-specific)
axis equal
在 Mathworks 文件交换 (ndanfilter.m) 中有一些代码接近我想要的，但我相信它只插入 NaNs 散落在 图像的内部，而不是显示这种岛型效果的数据。
注意：我刚发现 nanconv.m，它完全符合我的要求，用法非常直观（对图像进行卷积，忽略 NaN，很像 nanmean 的作品）。我已将这部分作为我接受的答案，并包括与其他答案的性能比较。
相关问题
Gaussian filtering a image with Nan in Python
好吧，不用你的绘图功能，我仍然可以给你一个解决方案。您要做的是找到所有新的 NaN's 并将其替换为原始的未过滤数据（假设它是正确的）。虽然它没有被过滤，但它比减少轮廓图像的域要好。
% Toy Example Data
rfVals= rand(100,100);
rfVals(1:2,:) = nan;
rfVals(:,1:2) = nan;
% Create and Apply Filter
filtWidth = 3;
imageFilter=fspecial('gaussian',filtWidth,filtWidth);
dataFiltered = imfilter(rfVals,imageFilter,'symmetric','conv');
sum(sum(isnan( dataFiltered ) ) )
% Replace New NaN with Unfiltered Data
newnan = ~isnan( rfVals) & isnan( dataFiltered );
dataFiltered( newnan ) = rfVals( newnan );
sum(sum(isnan( rfVals) ) )
sum(sum(isnan( dataFiltered ) ) )
使用以下代码检测新的 NaN。您也可以使用 xor 函数。
newnan = ~isnan( rfVals) & isnan( dataFiltered );
然后这一行将 dataFiltered 中的索引设置为 rfVals
中的值
dataFiltered( newnan ) = rfVals( newnan );
结果
从控制台和我的代码中打印的行，您可以看到 dataFiltered 中 NaN 的数量从 688 减少到 396，NaN 中的数量也是如此 rfVals.
ans =
688
ans =
396
ans =
396
备选方案 1
您还可以通过指定较小的内核并在之后合并它来在边缘附近使用较小的过滤器，但是如果您只想要使用最少代码的有效数据，我的主要解决方案就可以了。
备选方案 2
另一种方法是 pad/replace NaN 值为零或您想要的某个常量，以便它可以工作，然后截断它。但是对于信号 processing/filtering，您可能需要我的主要解决方案。
一种方法是在执行过滤之前使用 scatteredInterpolant (or TriScatteredInterp 在旧版 MATLAB 中使用最近邻插值替换 NaN 值），然后再次用 NaN 值替换这些点。这类似于使用 'replicate' 参数过滤完整的二维数组，而不是 'symmetric' 参数作为 imfilter 的边界选项（即您正在复制而不是反射锯齿状 NaN 边界处的值）。
代码如下所示：
% Make your filter:
filtWidth = 7;
imageFilter = fspecial('gaussian', filtWidth, filtWidth);
% Interpolate new values for Nans:
nanMask = isnan(rfVals);
[r, c] = find(~nanMask);
[rNan, cNan] = find(nanMask);
F = scatteredInterpolant(c, r, rfVals(~nanMask), 'nearest');
interpVals = F(cNan, rNan);
data = rfVals;
data(nanMask) = interpVals;
% Filter the data, replacing Nans afterward:
dataFiltered = imfilter(data, imageFilter, 'replicate', 'conv');
dataFiltered(nanMask) = nan;
我最终使用的技术是 Matlab 文件交换中的函数 nanconv.m。它完全符合我的要求：它以忽略 NaNs 的方式运行过滤器，就像 Matlab 的内置函数 nanmean 所做的那样。这很难从函数的文档中破译，它有点神秘。
以下是我的使用方法：
filtWidth = 7;
filtSigma = 5;
imageFilter=fspecial('gaussian',filtWidth,filtSigma);
dataFiltered = nanconv(data,imageFilter, 'nanout');
我正在粘贴下面的 nanconv 函数（它被 the BSD license 覆盖）。当我有机会时，我会 post 图片等，只是想 post 我最终为那些对我所做的事情感到好奇的人所做的事情。
与其他答案的比较
使用 结果在直觉上看起来非常好，但边缘有一些数量上的波动，这是一个令人担忧的问题。实际上，超出边缘的图像外推导致我的数据边缘出现许多虚假的高值。
使用 用原始数据替换缺失的位，看起来也很不错（特别是对于非常小的过滤器），但是（按设计）你最终会在边缘得到未过滤的数据，这是我的申请有问题。
nanconv
function c = nanconv(a, k, varargin)
% NANCONV Convolution in 1D or 2D ignoring NaNs.
% C = NANCONV(A, K) convolves A and K, correcting for any NaN values
% in the input vector A. The result is the same size as A (as though you
% called 'conv' or 'conv2' with the 'same' shape).
%
% C = NANCONV(A, K, 'param1', 'param2', ...) specifies one or more of the following:
% 'edge' - Apply edge correction to the output.
% 'noedge' - Do not apply edge correction to the output (default).
% 'nanout' - The result C should have NaNs in the same places as A.
% 'nonanout' - The result C should have ignored NaNs removed (default).
% Even with this option, C will have NaN values where the
% number of consecutive NaNs is too large to ignore.
% '2d' - Treat the input vectors as 2D matrices (default).
% '1d' - Treat the input vectors as 1D vectors.
% This option only matters if 'a' or 'k' is a row vector,
% and the other is a column vector. Otherwise, this
% option has no effect.
%
% NANCONV works by running 'conv2' either two or three times. The first
% time is run on the original input signals A and K, except all the
% NaN values in A are replaced with zeros. The 'same' input argument is
% used so the output is the same size as A. The second convolution is
% done between a matrix the same size as A, except with zeros wherever
% there is a NaN value in A, and ones everywhere else. The output from
% the first convolution is normalized by the output from the second
% convolution. This corrects for missing (NaN) values in A, but it has
% the side effect of correcting for edge effects due to the assumption of
% zero padding during convolution. When the optional 'noedge' parameter
% is included, the convolution is run a third time, this time on a matrix
% of all ones the same size as A. The output from this third convolution
% is used to restore the edge effects. The 'noedge' parameter is enabled
% by default so that the output from 'nanconv' is identical to the output
% from 'conv2' when the input argument A has no NaN values.
%
% See also conv, conv2
%
% AUTHOR: Benjamin Kraus (bkraus@bu.edu, ben@benkraus.com)
% Copyright (c) 2013, Benjamin Kraus
% $Id: nanconv.m 4861 2013-05-27 03:16:22Z bkraus $
% Process input arguments
for arg = 1:nargin-2
switch lower(varargin{arg})
case 'edge'; edge = true; % Apply edge correction
case 'noedge'; edge = false; % Do not apply edge correction
case {'same','full','valid'}; shape = varargin{arg}; % Specify shape
case 'nanout'; nanout = true; % Include original NaNs in the output.
case 'nonanout'; nanout = false; % Do not include NaNs in the output.
case {'2d','is2d'}; is1D = false; % Treat the input as 2D
case {'1d','is1d'}; is1D = true; % Treat the input as 1D
end
end
% Apply default options when necessary.
if(exist('edge','var')~=1); edge = false; end
if(exist('nanout','var')~=1); nanout = false; end
if(exist('is1D','var')~=1); is1D = false; end
if(exist('shape','var')~=1); shape = 'same';
elseif(~strcmp(shape,'same'))
error([mfilename ':NotImplemented'],'Shape ''%s'' not implemented',shape);
end
% Get the size of 'a' for use later.
sza = size(a);
% If 1D, then convert them both to columns.
% This modification only matters if 'a' or 'k' is a row vector, and the
% other is a column vector. Otherwise, this argument has no effect.
if(is1D);
if(~isvector(a) || ~isvector(k))
error('MATLAB:conv:AorBNotVector','A and B must be vectors.');
end
a = a(:); k = k(:);
end
% Flat function for comparison.
o = ones(size(a));
% Flat function with NaNs for comparison.
on = ones(size(a));
% Find all the NaNs in the input.
n = isnan(a);
% Replace NaNs with zero, both in 'a' and 'on'.
a(n) = 0;
on(n) = 0;
% Check that the filter does not have NaNs.
if(any(isnan(k)));
error([mfilename ':NaNinFilter'],'Filter (k) contains NaN values.');
end
% Calculate what a 'flat' function looks like after convolution.
if(any(n(:)) || edge)
flat = conv2(on,k,shape);
else flat = o;
end
% The line above will automatically include a correction for edge effects,
% so remove that correction if the user does not want it.
if(any(n(:)) && ~edge); flat = flat./conv2(o,k,shape); end
% Do the actual convolution
c = conv2(a,k,shape)./flat;
% If requested, replace output values with NaNs corresponding to input.
if(nanout); c(n) = NaN; end
% If 1D, convert back to the original shape.
if(is1D && sza(1) == 1); c = c.'; end
end
nanfilter在过滤的时候做的事情和nanconv完全一样，只要filter是一样的。如果您在使用 nanfilter 之前获取 nan 值，然后将其添加回过滤后的矩阵，您将获得与使用选项 'nanout' 从 nanconv 获得的结果相同的结果，只要您使用相同的筛选。
我有一个代表一些数据的二维数组 (doubles)，里面有一堆 NaNs。数据的等高线图如下所示：
所有空白都是NaNs，灰色菱形供参考，填充轮廓显示我的数据的形状。当我使用 imfilt 过滤数据时，NaNs 会显着影响数据，所以我们最终得到这样的结果：
您可以看到支持集显着收缩。我不能使用它，因为它已经咀嚼了边缘上的一些更有趣的变化（出于我实验的特定原因，这些边缘很重要）。
是否有在 NaNs 岛内进行过滤的功能，它处理类似于矩形过滤 windows 边缘的边缘，而不是仅仅消除边缘？有点像 nanmean 函数，除了卷积图像？
这是我的过滤器代码：
filtWidth = 7;
imageFilter=fspecial('gaussian',filtWidth,filtSigma);
%convolve them
dataFiltered = imfilter(rfVals,imageFilter,'symmetric','conv');
以及绘制等高线图的代码：
figure
contourf(dataFiltered); hold on
plot([-850 0 850 0 -850], [0 850 0 -850 0], 'Color', [.7 .7 .7],'LineWidth', 1); %the square (limits are data-specific)
axis equal
在 Mathworks 文件交换 (ndanfilter.m) 中有一些代码接近我想要的，但我相信它只插入 NaNs 散落在 图像的内部，而不是显示这种岛型效果的数据。
注意：我刚发现 nanconv.m，它完全符合我的要求，用法非常直观（对图像进行卷积，忽略 NaN，很像 nanmean 的作品）。我已将这部分作为我接受的答案，并包括与其他答案的性能比较。
相关问题
Gaussian filtering a image with Nan in Python
好吧，不用你的绘图功能，我仍然可以给你一个解决方案。您要做的是找到所有新的 NaN's 并将其替换为原始的未过滤数据（假设它是正确的）。虽然它没有被过滤，但它比减少轮廓图像的域要好。
% Toy Example Data
rfVals= rand(100,100);
rfVals(1:2,:) = nan;
rfVals(:,1:2) = nan;
% Create and Apply Filter
filtWidth = 3;
imageFilter=fspecial('gaussian',filtWidth,filtWidth);
dataFiltered = imfilter(rfVals,imageFilter,'symmetric','conv');
sum(sum(isnan( dataFiltered ) ) )
% Replace New NaN with Unfiltered Data
newnan = ~isnan( rfVals) & isnan( dataFiltered );
dataFiltered( newnan ) = rfVals( newnan );
sum(sum(isnan( rfVals) ) )
sum(sum(isnan( dataFiltered ) ) )
使用以下代码检测新的 NaN。您也可以使用 xor 函数。
newnan = ~isnan( rfVals) & isnan( dataFiltered );
然后这一行将 dataFiltered 中的索引设置为 rfVals
dataFiltered( newnan ) = rfVals( newnan );
结果
从控制台和我的代码中打印的行，您可以看到 dataFiltered 中 NaN 的数量从 688 减少到 396，NaN 中的数量也是如此 rfVals.
ans =
688
ans =
396
ans =
396
备选方案 1
您还可以通过指定较小的内核并在之后合并它来在边缘附近使用较小的过滤器，但是如果您只想要使用最少代码的有效数据，我的主要解决方案就可以了。
备选方案 2
另一种方法是 pad/replace NaN 值为零或您想要的某个常量，以便它可以工作，然后截断它。但是对于信号 processing/filtering，您可能需要我的主要解决方案。
一种方法是在执行过滤之前使用 scatteredInterpolant (or TriScatteredInterp 在旧版 MATLAB 中使用最近邻插值替换 NaN 值），然后再次用 NaN 值替换这些点。这类似于使用 'replicate' 参数过滤完整的二维数组，而不是 'symmetric' 参数作为 imfilter 的边界选项（即您正在复制而不是反射锯齿状 NaN 边界处的值）。
代码如下所示：
% Make your filter:
filtWidth = 7;
imageFilter = fspecial('gaussian', filtWidth, filtWidth);
% Interpolate new values for Nans:
nanMask = isnan(rfVals);
[r, c] = find(~nanMask);
[rNan, cNan] = find(nanMask);
F = scatteredInterpolant(c, r, rfVals(~nanMask), 'nearest');
interpVals = F(cNan, rNan);
data = rfVals;
data(nanMask) = interpVals;
% Filter the data, replacing Nans afterward:
dataFiltered = imfilter(data, imageFilter, 'replicate', 'conv');
dataFiltered(nanMask) = nan;
我最终使用的技术是 Matlab 文件交换中的函数 nanconv.m。它完全符合我的要求：它以忽略 NaNs 的方式运行过滤器，就像 Matlab 的内置函数 nanmean 所做的那样。这很难从函数的文档中破译，它有点神秘。
以下是我的使用方法：
filtWidth = 7;
filtSigma = 5;
imageFilter=fspecial('gaussian',filtWidth,filtSigma);
dataFiltered = nanconv(data,imageFilter, 'nanout');
我正在粘贴下面的 nanconv 函数（它被 the BSD license 覆盖）。当我有机会时，我会 post 图片等，只是想 post 我最终为那些对我所做的事情感到好奇的人所做的事情。
与其他答案的比较
使用
使用
nanconv
function c = nanconv(a, k, varargin)
% NANCONV Convolution in 1D or 2D ignoring NaNs.
% C = NANCONV(A, K) convolves A and K, correcting for any NaN values
% in the input vector A. The result is the same size as A (as though you
% called 'conv' or 'conv2' with the 'same' shape).
%
% C = NANCONV(A, K, 'param1', 'param2', ...) specifies one or more of the following:
% 'edge' - Apply edge correction to the output.
% 'noedge' - Do not apply edge correction to the output (default).
% 'nanout' - The result C should have NaNs in the same places as A.
% 'nonanout' - The result C should have ignored NaNs removed (default).
% Even with this option, C will have NaN values where the
% number of consecutive NaNs is too large to ignore.
% '2d' - Treat the input vectors as 2D matrices (default).
% '1d' - Treat the input vectors as 1D vectors.
% This option only matters if 'a' or 'k' is a row vector,
% and the other is a column vector. Otherwise, this
% option has no effect.
%
% NANCONV works by running 'conv2' either two or three times. The first
% time is run on the original input signals A and K, except all the
% NaN values in A are replaced with zeros. The 'same' input argument is
% used so the output is the same size as A. The second convolution is
% done between a matrix the same size as A, except with zeros wherever
% there is a NaN value in A, and ones everywhere else. The output from
% the first convolution is normalized by the output from the second
% convolution. This corrects for missing (NaN) values in A, but it has
% the side effect of correcting for edge effects due to the assumption of
% zero padding during convolution. When the optional 'noedge' parameter
% is included, the convolution is run a third time, this time on a matrix
% of all ones the same size as A. The output from this third convolution
% is used to restore the edge effects. The 'noedge' parameter is enabled
% by default so that the output from 'nanconv' is identical to the output
% from 'conv2' when the input argument A has no NaN values.
%
% See also conv, conv2
%
% AUTHOR: Benjamin Kraus (bkraus@bu.edu, ben@benkraus.com)
% Copyright (c) 2013, Benjamin Kraus
% $Id: nanconv.m 4861 2013-05-27 03:16:22Z bkraus $
% Process input arguments
for arg = 1:nargin-2
switch lower(varargin{arg})
case 'edge'; edge = true; % Apply edge correction
case 'noedge'; edge = false; % Do not apply edge correction
case {'same','full','valid'}; shape = varargin{arg}; % Specify shape
case 'nanout'; nanout = true; % Include original NaNs in the output.
case 'nonanout'; nanout = false; % Do not include NaNs in the output.
case {'2d','is2d'}; is1D = false; % Treat the input as 2D
case {'1d','is1d'}; is1D = true; % Treat the input as 1D
end
end
% Apply default options when necessary.
if(exist('edge','var')~=1); edge = false; end
if(exist('nanout','var')~=1); nanout = false; end
if(exist('is1D','var')~=1); is1D = false; end
if(exist('shape','var')~=1); shape = 'same';
elseif(~strcmp(shape,'same'))
error([mfilename ':NotImplemented'],'Shape ''%s'' not implemented',shape);
end
% Get the size of 'a' for use later.
sza = size(a);
% If 1D, then convert them both to columns.
% This modification only matters if 'a' or 'k' is a row vector, and the
% other is a column vector. Otherwise, this argument has no effect.
if(is1D);
if(~isvector(a) || ~isvector(k))
error('MATLAB:conv:AorBNotVector','A and B must be vectors.');
end
a = a(:); k = k(:);
end
% Flat function for comparison.
o = ones(size(a));
% Flat function with NaNs for comparison.
on = ones(size(a));
% Find all the NaNs in the input.
n = isnan(a);
% Replace NaNs with zero, both in 'a' and 'on'.
a(n) = 0;
on(n) = 0;
% Check that the filter does not have NaNs.
if(any(isnan(k)));
error([mfilename ':NaNinFilter'],'Filter (k) contains NaN values.');
end
% Calculate what a 'flat' function looks like after convolution.
if(any(n(:)) || edge)
flat = conv2(on,k,shape);
else flat = o;
end
% The line above will automatically include a correction for edge effects,
% so remove that correction if the user does not want it.
if(any(n(:)) && ~edge); flat = flat./conv2(o,k,shape); end
% Do the actual convolution
c = conv2(a,k,shape)./flat;
% If requested, replace output values with NaNs corresponding to input.
if(nanout); c(n) = NaN; end
% If 1D, convert back to the original shape.
if(is1D && sza(1) == 1); c = c.'; end
end
nanfilter在过滤的时候做的事情和nanconv完全一样，只要filter是一样的。如果您在使用 nanfilter 之前获取 nan 值，然后将其添加回过滤后的矩阵，您将获得与使用选项 'nanout' 从 nanconv 获得的结果相同的结果，只要您使用相同的筛选。