lap_fus.m

matlab code for image fusion used by canga.

function [fusion] = lap_fus(im1,im2,ns,consistency,a)
%
%lap_fus.m, v 1.2 2002/03/14 11:07:40
%===========================================================================
%               Eduardo Fernandez Canga - University of Bath
%
%                        Copyright (c) 2002
%===========================================================================
%
%                   [fusion] = lap_fus(im1,im2,ns,consistency)  
%
%                   im1: input image
%                   im2: input image
%                    ns: number of scales (default 4)
%           consistency: apply consistency
%                           1=yes (default)
%                           0=no 
%
%===========================================================================

if nargin < 5, a=6/16;end
if nargin < 4, consistency=1;end
if nargin < 3, ns=4;end

if a>1 | a<0
    error('Error: Invalid filter')
end

if any (size(im1)~=size(im2))
    error('Error: Different Size Images')
end

r=2^ns;
sx=size(im1,1);     % check image row size to know if it is  
x=mod(sx,r);        % divisible by 2^number of escales
if x                % if not, add rows
    im1=rowpad(im1,r-x);
    im2=rowpad(im2,r-x);
end
sy=size(im1,2);     % check image col size to know if it is  
y=mod(sy,r);        % divisible by 2^number of escales
if y                % if not, add cols
    im1=colpad(im1,r-y);
    im2=colpad(im2,r-y);    
end

% Laplacian Filter
r=[(1/4-a/2) 1/4 a 1/4 (1/4-a/2)];
w=conv2(r,r');

g1_1=reduce(im1,w);     % apply filter and downsample in each image
g1_2=reduce(im2,w);     % to obtain the next level of gaussian pyramid

g0_1=expand(g1_1,w);    % upsample and apply filter in each image
g0_2=expand(g1_2,w);    % to estimate the previous level of gaussian pyramid

L1=im1-g0_1;            % substract estimation from previous image
L2=im2-g0_2;            % to obtain current level of the laplacian pyramid

decision=abs(L1)>abs(L2);       % obtain decision mask for current level fusion
if consistency==1               % apply consistency check if it is active
    mask=ones(3,3)/9;           % use 3x3 window mask
    decision=impad(double(decision),1);
    decision=round(conv2(decision,mask,'valid'));
end
L=L1.*decision+L2.*~decision;   % perform fusion using decision mask

ns=ns-1;                % decrease scale number
if ns > 0               % apply recursively till top level is reached
    fus=lap_fus(g1_1,g1_2,ns);
else
    fus=(g1_1+g1_2)/2;  % top level fused by averaging
end
fus=expand(fus,w)+L;    % inverse transform of the laplacian pyramid
fusion=fus(1:sx,1:sy);  % select valid area of the final fusion


function [y] = reduce(x,w)
% performs filtering and downsampling by two
y=downsamp(conv2(impad(x,2),w,'valid'),2);

function [y] = expand(x,w)
% performs upsampling by two and filtering
y=conv2(impad(upsamp(x,2),2),4*w,'valid');

function [y] = downsamp(x,r)
[a b]=size(x);
y = x(1:r:a,1:r:b);

function [y] = upsamp(x,r)
[a b]=size(x);
y=zeros(a*r,b*r);
y(1:r:a*r,1:r:b*r) = x;