%% computes and displays bivariate normal densities
% Also relies on the function MVN in the separate file, MVN.m
% 

corr = 0.4; % correlation coefficient between X1 and X2. 
G = ones(2,2);   % set Var(X1) = Var(X2) = 1
G(1,2) = corr;
G(2,1) = corr;

Ginv = inv(G);
detG = det(G);
D = size(G);
N = D(1);
factor = (2*pi)^(N/2) *sqrt(det(G))  % miltiplier in front of the exp(...) in the MVN density


K = 40;
D = K+1;
fGrid = ones(D,D);
for xi=0:K,
 for yj = 0:K,  
  x = [xi*2/K, yj*2/K]';    
  [f] = exp(-(xi+yj)/2);
  fGrid(xi+1,yj+1) = f;  
 end
end
xvec = (0:K)*2/K;
surf(xvec,xvec,fGrid)