function se = BOGaussKernel(L,theta,sigma)
% BOGaussKernel - calculate 2D Gaussian kernel
% 
% Boguslaw Obara, http://boguslawobara.net/
%
% Help: 
%       se = BOGaussKernel(L,theta,sigma)
%       'L'         - input length of kernel
%       'theta'     - input theta value
%       'sigma'     - input sigma value
%       'se'        - output kernel matrix
% Example:
%       L=9.0; theta=pi/4; sigma=2.0;
%       SE=BOGaussKernel(L,theta,sigma);
%
% Version:
%   0.1 - 24/09/2005 First implementation

%% Settings
x_size=0; sum=0.0;
start=-(L-1); stop=(L-1);
size=abs(start)+abs(stop)+1;
K=zeros(size,size);
se=zeros(size,size);
%%
for x=start:stop
    for y=start:stop
        u=double(x)*cos(theta)-double(y)*sin(theta);
        v=double(x)*sin(theta)+double(y)*cos(theta);	
        if (abs(u)<=3.0*sigma && abs(v)<=L/2.0) 
            K(x+abs(start)+1,y+abs(start)+1)=-exp(-(u*u)/(2.0*sigma*sigma));
            x_size=x_size+1;
        else
            K(x+abs(start)+1,y+abs(start)+1)=0.0;
        end
    end
end
%%
for x=start:stop
    for y=start:stop
        sum=sum+K(x+abs(start)+1,y+abs(start)+1);
    end
end
me=sum/double((x_size));
%%
for x=start:stop
    for y=start:stop
        u=double(x)*cos(theta)-double(y)*sin(theta);
        v=double(x)*sin(theta)+double(y)*cos(theta);	
        if abs(u)<=3*sigma && abs(v)<=L/2
            se(x+abs(start)+1,y+abs(start)+1)=(K(x+abs(start)+1,y+abs(start)+1)-me);
        else
            se(x+abs(start)+1,y+abs(start)+1)=0.0;
        end
    end
end
%%
end