q916_qfunction.m

using Monte Carlo integeration calculate Q function

y=input('Enter y, argument of Q function >');
M=input('Enter M, argument of experiments >');
N=input('Enter N, argument of trials >');
u=(pi/2)*rand(N,M);    %Generate random numbers
uu=(1/pi).*exp(-(y^2)./(2*(sin(u)).^2)); %Define functon
data=zeros(N,M);%Initialize array
data1=zeros(N,1);
%The following four lines of code determine
%M estimates as a function of j, 0<j<=N
data(1,:)=uu(1,:);
for j=2:N
    data(j,:)=(pi/2)*sum(uu(1:j,:))/j;
end
est=data(N,:);   %N estimates of Q(y)
est1=sum(est)/M;
for k=1:N
    data1(k,1)=sum(data(k,:))/M;
end;
real=qfunc(y)*ones(1,N);
n=100:N;
plot(n,data1(100:N),'k',n,real(100:N),'k--')  %plot resulte
xlabel('Number of Trials')
ylabel('Estimate of Q function')
%End