clc;
clear all;
ns=2;
alpha=1.3	;
number_sample=500;
 R=[1 , 0.2+0.1i;
    0.2-0.1i, 1];
X = wgn(ns,number_sample,0,'complex');
A=chol(R);
D=inv(A)';
xx=inv(D)*X;
%%%%%%%%%%%%%%%%%%產生列為隨機變量%%%%%%%%%%%
if alpha==2
     LV=ones(1,number_sample);
else LV=stblrnd(alpha/2,-1,cos(alpha*pi/4)^(2/alpha),0,1,number_sample) ;
end 
%%%%%%%%%%%%%%%%%%產生次高斯隨機矢量%%%%%%%%%%%
for t=1:ns
    y(t,:)=sqrt(LV(1,:)).*xx(t,:);
end
gam=(R(1,1)/2)^alpha
x=-5:0.1:5;
p_z = stblpdf(x,alpha,0,gam,0);
plot(x,p_z);hold on
[a,b]=stblest(y(2,:))