function y=kernel_density(x,xdat,b,method)

%畫核密度函數圖，其中b為窗寬
%xmin=min(xdat);
%xmax=max(xdat);
%band=(xmax-xmin)/1000;
%x=xmin:band:xmax;
n=length(x);
m=length(xdat);
for i=1:n
    xxb=(xdat-x(i))/b;
    xxb=xxb(abs(xxb)<=1);
    if strcmp(method,'epan')
        kx=(1-xxb.^2)*3/4;
    elseif strcmp(method,'gauss')
        kx=exp(-xxb.^2/2)/sqrt(2*pi);
    elseif strcmp(method,'tri3')
        kx=(1-xxb.^3).^3;
    elseif strcmp(method,'tri2')
        kx=(35/32)*(1-xxb.^2).^3;
    elseif strcmp(method,'qua')
        kx=(15/16)*(1-xxb.^2).^2;
    elseif strcmp(method,'cosi')
        kx=pi*cos(pi*xxb/2)/4;
    end
    yy(i)=sum(kx)/(m*b);
end

if nargout==0
    plot(x,yy,'r')
else
    y=yy;
end