How Can I remove the Periodic Noise lines from Color image attached?
11 views (last 30 days)
Show older comments
Zain ul abideen ahsan
on 30 Nov 2020
Commented: Zain ul abideen ahsan
on 1 Dec 2020
Hy, I am trying to filter out the periodic lines in Colour image attached using FFT and IFFT and median filter. Can anybody help me in this please. The end result should be same colour picture but those grey lines gone or atlease faded away. Help will be appreciated
I am attaching the Picture and Code modified with the help of @Image Analyst.
% Demo to filter an image, with periodic ripple, in the Fourier domain.
clc; % Clear the command window.
close all; % Close all figures (except those of imtool.)
clear; % Erase all existing variables. Or clearvars if you want.
workspace; % Make sure the workspace panel is showing.
format long g;
format compact;
fontSize = 14;
%Read Image
grayImage = imread('Cereal_noise.png');
[rows, columns, numberOfColorBands] = size(grayImage);
if numberOfColorBands > 1
grayImage = rgb2gray(grayImage);
end
subplot(2, 3, 1);
imshow(grayImage, [0 255]);
title('Original Image', 'FontSize', fontSize);
set(gcf, 'units','normalized','outerposition',[0 0 1 1]); % Maximize figure.
% Compute the 2D fft.
frequencyImage = fftshift(fft2(grayImage));
% Take log magnitude so we can see it better in the display.
amplitudeImage = log(abs(frequencyImage));
minValue = min(min(amplitudeImage))
maxValue = max(max(amplitudeImage))
subplot(2, 3, 4);
imshow(amplitudeImage, []);
caption = sprintf('Notice the two spikes\nperpendicular to the periodic frequency');
title(caption, 'FontSize', fontSize);
axis on;
% zoom(10)
% Find the location of the big spikes.
amplitudeThreshold = 16;
brightSpikes = amplitudeImage > amplitudeThreshold; % Binary image.
subplot(2, 3, 5);
imshow(brightSpikes);
axis on;
title('Bright Spikes', 'FontSize', fontSize);
% Let user see the image.
promptMessage = sprintf('The image below shows the bright spikes.\nNow we will exclude the central spike.');
titleBarCaption = 'Continue?';
button = questdlg(promptMessage, titleBarCaption, 'Continue', 'Cancel', 'Continue');
if strcmpi(button, 'Cancel')
return;
end
% Exclude the central DC spike. Everything from row 115 to 143.
brightSpikes(115:143, :) = 0;
imshow(brightSpikes);
title('Bright spikes other than central spike', 'FontSize', fontSize);
promptMessage = sprintf('Now we will use these bright spikes to filter (mask) the spectrum.');
button = questdlg(promptMessage, titleBarCaption, 'Continue', 'Cancel', 'Continue');
if strcmpi(button, 'Cancel')
return;
end
% Filter/mask the spectrum.
frequencyImage(brightSpikes) = 0;
% Take log magnitude so we can see it better in the display.
amplitudeImage2 = log(abs(frequencyImage));
minValue = min(min(amplitudeImage2))
maxValue = max(max(amplitudeImage2))
subplot(2, 3, 5);
imshow(amplitudeImage2, [minValue maxValue]);
axis on;
title('Spikes zeroed out', 'FontSize', fontSize);
% zoom(10)
filteredImage = ifft2(fftshift(frequencyImage));
amplitudeImage3 = abs(filteredImage);
minValue = min(min(amplitudeImage3))
maxValue = max(max(amplitudeImage3))
subplot(2, 3, 6);
imshow(amplitudeImage3, [minValue maxValue]);
title('Filtered Image', 'FontSize', fontSize);
0 Comments
Accepted Answer
Image Analyst
on 30 Nov 2020
I would not do this in the Fourier domain since the lines are so perfect, sharp, and periodic in the spatial domain. How do you think you would do it? What are your ideas? How about you sum the image vertically to get a horizontal profile? By looking at each color channel in turn, you can determine if the striping is in just one color channel, or all of them. Then examine the sum in the brighter and darker areas and determine if the brighter ones are a multiplication factor brighter, or an additive factor brighter. Then adjust the stripes by adding or multiplying by the correct factor. So that should be a lot more precise than anything you do in the frequency domain. Give it a try and if you really can't figure it out, come back here to get rescued.
3 Comments
Image Analyst
on 30 Nov 2020
You can try this. But the problem is that your image is saturated. Clipped! So you don't know what the values were originally supposed to be so subtracting a constant won't be able to fix it in those areas where it was clipped. To reconstruct an image perfectly you need to not have data loss due to things like clipping/saturation.
% Demo to find stripes and correct them in the image. By Image Analyst, November 30, 2020.
clc; % Clear the command window.
close all; % Close all figures (except those of imtool.)
clearvars;
workspace; % Make sure the workspace panel is showing.
format long g;
format compact;
fontSize = 16;
fprintf('Beginning to run %s.m ...\n', mfilename);
%-----------------------------------------------------------------------------------------------------------------------------------
% Read in image.
folder = [];
baseFileName = 'cereal.png';
fullFileName = fullfile(folder, baseFileName);
% Check if file exists.
if ~exist(fullFileName, 'file')
% The file doesn't exist -- didn't find it there in that folder.
% Check the entire search path (other folders) for the file by stripping off the folder.
fullFileNameOnSearchPath = baseFileName; % No path this time.
if ~exist(fullFileNameOnSearchPath, 'file')
% Still didn't find it. Alert user.
errorMessage = sprintf('Error: %s does not exist in the search path folders.', fullFileName);
uiwait(warndlg(errorMessage));
return;
end
end
rgbImage = imread(fullFileName);
[rows, columns, numberOfColorChannels] = size(rgbImage)
% Display the test image full size.
subplot(2, 2, 1);
imshow(rgbImage, []);
axis('on', 'image');
caption = sprintf('Reference Image : "%s"', baseFileName);
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
drawnow;
hp = impixelinfo(); % Set up status line to see values when you mouse over the image.
% Set up figure properties:
% Enlarge figure to full screen.
hFig1 = gcf;
hFig1.Units = 'Normalized';
hFig1.WindowState = 'maximized';
% Get rid of tool bar and pulldown menus that are along top of figure.
% set(gcf, 'Toolbar', 'none', 'Menu', 'none');
% Give a name to the title bar.
hFig1.Name = 'Demo by Image Analyst';
[r, g, b] = imsplit(rgbImage);
% How about you sum the image vertically to get a horizontal profile?
horizontalProfile = mean(b(894:end, :));
subplot(2, 2, [2, 4]);
plot(horizontalProfile, 'b-');
grid on;
title('Horizontal Profile between rows 894 and 903', 'FontSize', fontSize, 'Interpreter', 'None');
drawnow;
% By looking at each color channel in turn, you can determine if the striping
% is in just one color channel, or all of them.
% Then examine the sum in the brighter and darker areas and determine
% if the brighter ones are a multiplication factor brighter, or an additive factor brighter.
% From the plot, it looks like it shows up in all 3 color channels and
% there is simple an additive 50 added to all color channels.
% Then adjust the stripes by adding or multiplying by the correct factor.
% but first we need to find out what columns are the brighter ones.
% first fine the min
bottomEnvelope = movmin(horizontalProfile, 20);
upperEnvelope = movmax(horizontalProfile, 20);
% Find the mean delta gray level between columns 200 and 300.
% This is where there is no saturation (over 255).
deltaGL = mean(upperEnvelope(200:300) - bottomEnvelope(200:300))
deltaGL = uint8(deltaGL); % Cast to uint8 since image was uint8.
hold on;
plot(bottomEnvelope, 'r-', 'LineWidth', 2);
plot(upperEnvelope, 'r-', 'LineWidth', 2);
% Compute midline
midline = (bottomEnvelope + upperEnvelope) / 2;
plot(midline, 'm-', 'LineWidth', 2);
% Determine what rows are above the midline.
% Those are the columns that are 50 gray levels brighter.
columnsToDim = horizontalProfile > midline;
% Darken those columns in each color channel.
r(:, columnsToDim) = r(:, columnsToDim) - deltaGL;
g(:, columnsToDim) = g(:, columnsToDim) - deltaGL;
b(:, columnsToDim) = b(:, columnsToDim) - deltaGL;
% Reconstruct output image.
rgbOut = cat(3, r, g, b);
subplot(2, 2, 3);
imshow(rgbOut, []);
caption = sprintf('Corrected Image where %d was subtracted from the stripes.', deltaGL);
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
drawnow;
hp = impixelinfo(); % Set up status line to see values when you mouse over the image.
More Answers (1)
sushanth govinahallisathyanarayana
on 30 Nov 2020
You can filter it in the frequency domain with a suitable low pass filter. Your code seems to indicate that you know where the fourier transform is spiking.
See Also
Categories
Find more on Spectral Measurements in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!