Lets Learn together... Happy Reading

" Two roads diverged in a wood, and I,
I took the one less traveled by,
And that has made all the difference "-Robert Frost

Recent Posts

UPSAMPLING IN MATLAB

Upsampling is the process of inserting zeros in between the signal value in order to increase the size of the matrix.  We will discuss about upsampling in both spatial and time domain.
1.1  Upsampling using MATLAB built-in function
1.2  Upsampling in 1D
1.3  Upsampling in 2D or image matrix
2.1  Upsampling a 1D signal
2.2  Upsampling a image matrix

UPSAMPLING IN SPATIAL DOMAIN:

Given: 1-D array of size 1xN
Upsample by factor of 2:  

Output: 1-D array of size 1x(2*N)
Where N represents length of the array, n represents the index starting from 0,1,2…,N
%UPSAMPLING USING MATLAB BUILT-IN FUNCTION 'UPSAMPLE'

A = 1:25;
B = upsample(A,2);

EXPLANATION:

The above MATLAB function will insert zeros in between the samples.

To upsample an array by ratio 2, update the output array as follows:
1.      If index(n) of the output array is divisible by 2(ratio), then update the output array with the value of the input array with index n/2
2.      Otherwise, insert zero

STEPS TO PERFORM:
1.      Consider an array A of size 1xM
2.      Obtain the upsample Ratio N
3.      Pre-allocate an array B of size 1x(M*N)
4.      If the index is divisible by N then update the array B with value of A else zero


MATLAB CODE:

%UPSAMPLING IN SPATIAL DOMAIN

A = 1:10;

%UPSAMPLING WITH RATIO 2
B = zeros([1 size(A,2)*2]);
B(1:2:end)=A;



%UPSAMPLING WITH RATIO N
N=3;
B=zeros([1 size(A,2)*N]);
B(1:N:end)=A;

EXPLANATION:
Let A = [1 2 3 4 5]
Let us upsample try to upsample A with ratio 2
Pre-allocate output matrix B = [0 0 0 0 0 0 0 0 0 0]
Substitute the value of the matrix A in indices divisible by the ratio i.e. 2 of matrix B
B(1:2:end) = A
Now B(1,3,5,7,9) =[1 2 3 4 5]
Thus B = [1 0 2 0 3 0 4 0 5 0]

NOTE:
Definition of upsampling is usually given assuming the index starts from zero. But in case of MATLAB, the index starts with one that is why the odd positions are considered instead of even.

STEPS TO PERFORM TO UPSAMPLE A 2D MATRIX:


MATLAB CODE:



%IMAGE UPSAMPLING

A = imread('cameraman.tif');

M = 2;
N = 3;

B = zeros([size(A,1)*M size(A,2)*N]);
B(1:M:end,1:N:end) = A;

figure,imagesc(B);colormap(gray);


sz = M*N;
H = fspecial('average',[sz sz]);
C = conv2(B,H,'same');

figure,imagesc(C);colormap(gray);


EXPLANATION:


The above image is the pixel representation of the zero inserted image. In each row, two zeros are inserted between the pixels and in the each column; single zero is inserted between the pixels. In spatial domain, inserting zeros will be quite visible so it’s advisable to perform any low pass filtering or approximation like spatial averaging or applying Gaussian.
In the above example, spatial averaging is done. In order to understand how spatial averaging is done using convolution check the following link:








like button Like "IMAGE PROCESSING" page

UPSAMPLING IN FREQUENCY DOMAIN



1.1  Upsampling using MATLAB built-in function
1.2  Upsampling in 1D
1.3  Upsampling in 2D or image matrix
2.1  Upsampling a 1D signal
2.2  Upsampling a image matrix

In Frequency domain, upsampling means nothing but the padding of zeros at the end of high frequency components on both sides of the signal.

STEPS TO PERFORM:

1.      Read an image
2.      Obtain the ratio to upsample
3.      Perform Fast Fourier Transform
4.      Shift the Low frequency components to the centre and High frequency components outside.
5.      Add zeros on both the sides of the image
6.      Shift the High frequency components to the centre and Low frequency components to the exterior (Inverse of fftshift)
7.      Perform Inverse Fast Fourier Transform
8.      Display the Upsampled Image

MATLAB CODE:

%UPSAMPLING IN FREQUENCY DOMAIN
%1D UPSAMPLING

FS = 100;
t  = 0:(1/FS):1;

A = 10*sin(2*pi*5*t);
figure,plot(t,A);
 



%FOURIER DOMAIN

FT = fft(A);
fq =linspace(-1/FS,1/FS,101);
figure,plot(fq,abs(FT));title('Before FFTSHIFT');


FT_s = fftshift(FT);
figure,plot(fq,abs(FT_s));title('After FFTSHIFT');


pad_zero = padarray(FT_s,[0 50]);

fq =linspace(-1/FS,1/FS,201);
figure,plot(fq,abs(pad_zero));title('After PADDING WITH ZEROS');


%INVERSE FOURIER TRANSFORM
IFT = ifft(ifftshift(pad_zero));

%UPSAMPLED SIGNAL
t1 = linspace(0,1,201);
figure,plot(t1,(IFT*2),'r',t,A,'g');




EXPLANATION:

Amplitude of the input original signal is 10 and the frequency is 5.
Similarly, the amplitude of the upsampled signal is 10 and the frequency is 5. The number of samples used to plot the signal is increased in the later case.



IMAGE UPSAMPLING IN FOURIER DOMAIN

MATLAB CODE:

%READ AN INPUT IMAGE
A = imread('cameraman.tif');

%RATIO
RatioM = 3;
RatioN = 2;

%UPSAMPLING OVER EACH ROW

mnrow = round(size(A,2)*(RatioM-1)/2);
% 1D FFT ON EACH ROW
row_fft = fft(A,[],2);

%PAD WITH ZEROS ON BOTH SIDES OF EACH ROW
pad_row = padarray(fftshift(row_fft,2),[0 mnrow]);

 
Logarthmic scale
%UPSAMPLING OVER EACH COLUMN
mncol = round(size(A,1)*(RatioN-1)/2);

% 1D FFT ON EACH COLUMN
col_fft = fft(pad_row,[],1);

%PAD WITH ZEROS ON BOTH SIDES OF EACH COLUMN
pad_col = padarray(fftshift(col_fft,1),[mncol 0]);
 
Logarthmic scale
%PERFORM 1D IFFT ON EACH COLUMN
ifft_col = ifft(ifftshift(pad_col,1),[],1);

%PERFORM 1D IFFT ON EACH ROW
ifft_col_row = ifft(ifftshift(ifft_col,2),[],2);

%DISPLAY THE IMAGE AFTER UPSAMPLING
res = abs(ifft_col_row);
res = uint8(res*(numel(res)/numel(A)));


figure,imagesc(res);




SIMPLE VERSION :

A = imread('cameraman.tif');

%RATIO
RatioM = 3;
RatioN = 2;



mnrow = round(size(A,2)*(RatioM-1)/2);
mncol = round(size(A,1)*(RatioN-1)/2);

%FFT ON 2D MATRIX
FT = fftshift(fft2(A));


%PADDING WITH ZEROS
pad_rc = padarray(FT,[mncol mnrow]);


%INVERSE FOURIER TRANSFORM
IFT = ifft2(ifftshift(pad_rc));

Img = uint8(abs(IFT)*(numel(IFT)/numel(A)));
figure,imagesc(Img);




like button Like "IMAGE PROCESSING" page

CONVOLUTION IN MATLAB

                 Let us try to understand convolution by performing spatial averaging on a matrix without using MATLAB built in function ‘conv2()’.

%CONVOLUTION IN MATLAB

%INPUT  MATRIX
A = zeros(5);
A(:) = 1:25;

%KERNEL
avg3 = ones(3)/9;

%CONVOLUTION
Result = conv2(A,avg3,'same');
display(Result);

Steps to be performed:





NOTE :
To define a kernel for spatial averaging, fill the kernel with ones and divide it by the number of elements in it.
For instance, consider kernel of size 4x4 , fill the  matrix with ones and divide it by 16. i.e the total number of elements in the matrix.

MATLAB CODE:
%INPUT  MATRIX
A = zeros(5);
A(:) = 1:25;

%KERNEL
avg3 = ones(3)/9;

%PAD THE MATRIX WITH ZEROS
B = padarray(A,[1 1]);
% PRE-ALLOCATE THE MATRIX
Output = zeros([size(A,1) size(A,2)]);

%PERFORM COONVOLUTION
for i = 1:size(B,1)-2
    for j = 1:size(B,2)-2
        Temp = B(i:i+2,j:j+2).*avg3;
        Output(i,j) = sum(Temp(:));
    end
 end

display(Output);

like button Like "IMAGE PROCESSING" page

Multilook Technique for speckle reduction




Consider a stack of images  affected by speckle noise of a same scene. 

In order to reduce the speckle, the availability of the stack of images can be utilized to obtain speckle reduced image.
In order to understand the multilook technique, let’s first generate noisy images and then apply speckle reduction.




Generate stack of noisy images:
1.      Read a noise free image
2.      Generate exponential noise with mean 1 and multiply it with the noise free image

3.      Repeat the process of generating exponential noise and multiplying  it with the noise free image to create stack of ‘n’ images

MATLAB code:
%Multilook - Speckle reduction

%Read a noisy free Image
I = imread('cameraman.tif');
I = double(I);
figure,imagesc(I);colormap(gray);title('Original Image');

[m, n] = size(I);
numofimg = 9;
Stack = zeros([m n numofimg]);

for i = 1:numofimg
    %Generate exponential noise
    noise=random('exp',1,m,n);
   
    %Multiply Noise free Image with noise
    Stack(:,:,i)=I.*noise;
   
end

Explanation:

Random noise with exponential distribution is multiplied with the noise free image to generate image affected by speckle.





Fig: Probability Density Function (PDF) of the noise generated for ‘n’ images that follows exponential distribution.

Multilook(ML) Technique:
Multilook image is the average of the stack of images. This technique is used widely in the field of Synthetic Aperture Radar(SAR) to reduce speckle on the same scene but different acquisition periods.
The stack of data is created for the same scene but with different time of acquisition and ML is applied to reduce the speckle considerably since the noise on the single image will be quite high.

MATLAB code:
ML_image = mean(Stack,3);
figure,imagesc(Stack(:,:,1));colormap(gray);title('Single Noisy Image');
figure,imagesc(ML_image);colormap(gray);title('Multilook Image');



Explanation:


For instance, in the above shown image, consider the red dot as the pixel value of the first pixel position of each image in the stack. Add all the pixel values represented in red dot and divide by the number of images in the stack. Similarly, perform the average for the second pixel position (green dot) and so on and so forth for the whole image.

Matlab code ‘mean(Stack,3) ‘ finds the mean of the image in 3rd dimension and the final result is the Multilooked Image with reduced speckle.


Moving Average Filter:
In order to reduce speckle further, a moving average filter can be used.
Moving average filter of window size 3x3 and 5x5 is applied on the multilook image and window of size 3x3 on a single speckled image in order to appreciate the multilook technique performance.

MATLAB code:

box = ones(3)/9;
ML_avg3 = conv2(ML_image,box,'same');
figure,imagesc(ML_avg3);colormap(gray);title('Moving Average - Window 3 x 3');
box = ones(5)/25;
ML_avg5 = conv2(ML_image,box,'same');
figure,imagesc(ML_avg5);colormap(gray);title('Moving Average - Window 5 x 5');
box = ones(3)/9;
avg_flt = conv2(Stack(:,:,1),box,'same');
figure,imagesc(avg_flt);colormap(gray);title('Moving Average - Single Image');







like button Like "IMAGE PROCESSING" page
Next Post Home

Link

Related Posts Plugin for WordPress, Blogger...
Google ping Hypersmash.com