function [z,h]=akdest2d(x,y,A,k,h)

% AKDEST2D    Adaptive Kernel Density Estimator for bi-variate data.
%
%             AKDEST2D(X,Y,A,K,H) evaluates the bi-variate Adaptive Kernel 
%             Density Estimate with kernel given by the function in K.
%
%             'KNORM' is the multivariate normal distribution.
%             'KEPAN' is the Epanechnikov kernal.
%
%             (X,Y): The point for which we want the density estimate.
%             A: The bivariate data is stored in the rows of A.
%             H: The optional parameter H represents the window width.
%                If omitted an 'optimal' value is used.
%
%             Christian C. Beardah 1994-95

if nargin==4,
  h=0;
end;

[n d]=size(A); % Establish dimensions of A

alpha=1/d;

xd=A(:,1);
yd=A(:,2);

if h==0,
  h=[hns(xd,k), hns(yd,k)];
end;

[f,h]=kdest2d(xd,yd,A,k,h);

l=log(f);
t=sum(l);

g=exp(t/n);

lambda=(f/g).^(-alpha);
h1=h(1)*lambda;
h2=h(2)*lambda;

my_sum = 0;
for i=1:n;
   if strcmp(k,'knorm')==1;
      temp = knorm((x-A(i,1)*ones(size(x)))/h1(i),(y-A(i,2)*ones(size(y)))/h2(i));
   elseif strcmp(k,'kepan')==1;
      temp = kepan((x-A(i,1)*ones(size(x)))/h1(i),(y-A(i,2)*ones(size(y)))/h2(i));
   elseif strcmp(k,'krect')==1;
      temp = krect((x-A(i,1)*ones(size(x)))/h1(i),(y-A(i,2)*ones(size(y)))/h2(i));
   elseif strcmp(k,'ktria')==1;
      temp = ktria((x-A(i,1)*ones(size(x)))/h1(i),(y-A(i,2)*ones(size(y)))/h2(i));
   elseif strcmp(k,'kbiwe')==1;
      temp = kbiwe((x-A(i,1)*ones(size(x)))/h1(i),(y-A(i,2)*ones(size(y)))/h2(i));
   elseif strcmp(k,'ktriw')==1;
      temp = ktriw((x-A(i,1)*ones(size(x)))/h1(i),(y-A(i,2)*ones(size(y)))/h2(i));
   end;
   my_sum = my_sum+((h1(i)*h2(i))^(-1)*temp);
end;


%my_sum=(h1(1)*h2(1))^(-1)*...
%    feval(k,(x-A(1,1)*ones(size(x)))/h1(1),(y-A(1,2)*ones(size(y)))/h2(1));

%for i=2:n,
%  my_sum=my_sum+(h1(i)*h2(i))^(-1)*...
%      feval(k,(x-A(i,1)*ones(size(x)))/h1(i),(y-A(i,2)*ones(size(y)))/h2(i));
%end;    

z=my_sum/n;
 