function PoleFigureAngles(varargin)
hs=varargin{1};
qs=varargin{2};
rotation=0;
if length(varargin)>2
filename=varargin{3};
if length(varargin)==4
rotation=varargin{4};
end
end
if size(hs,1)==size(qs,1)
for j=1:size(hs,1)
% for each plane and substrate pair calculate the alpha and beta angles
h=hs(j,:);
% find all permutations of measurement plane
qperm=perms(qs(j,:));
% calculate the all permutations with both plus and minus sign
signperm=[1,1,1;unique(perms([-1,1,1]),'rows');unique(perms([-1,-1,1]),'rows');-1,-1,-1];
allperm=zeros(size(qperm,1)*size(signperm,1),3);
k=size(qperm,1);
for i =1: size(signperm,1)
allperm((i-1)*k+1:i*k,:)= qperm.*signperm(i,:);
end
% Find all unique permutations
uniqperm=unique(allperm,'rows');
alpha=zeros(length(uniqperm),1);
onplane=zeros(size(uniqperm));
%Calculate alpha angles and their projection on the plane
for i=1:length(uniqperm)
alpha(i)=90-acosd(h*uniqperm(i,:)'/(norm(h)*norm(uniqperm(i,:))));
onplane(i,:) = uniqperm(i,:)-h*uniqperm(i,:)'/norm(h)^2*h;
end
% Calculate beta angles from the first vector in onplane
beta=zeros(size(onplane,1),1);
betat=zeros(size(onplane,1),3);
for i=2:size(onplane,1)
betat(i,:)=cross(onplane(1,:),onplane(i,:));
if betat(i,:)*h'<0
beta(i)=360-atan2d(norm(cross(onplane(1,:),onplane(i,:))),...
dot(onplane(1,:),onplane(i,:)));
else
beta(i)=atan2d(norm(cross(onplane(1,:),onplane(i,:))),...
dot(onplane(1,:),onplane(i,:)));
end
end
% sort angles acorting to alpha and print to console
[A,I]=sort(alpha);
index=I(A>0);
angles=[alpha(index),beta(index)];
str=strcat(sprintf('alpha\t beta'),sprintf('\n%3.2f\t %3.2f',angles'));
disp(str);
% Plot angles on polar figure
if j==1
scat=figure;
else
figure(scat)
hold on
end
polarscatter(beta(index)*pi/180,90-alpha(index))
pax=gca;
pax.ThetaZeroLocation = 'top';
rlim([0 90])
rticks([0 30 60 90])
rticklabels({'90','60','30','0'});
% If a data set from 3D explore is given, overlay the angles on the
% data
if exist('filename','var')
if j==1
cont=figure;
hold on
axis square
axis off
polefig=importdata(filename,'\t',14);
sizes=textscan(polefig.textdata{14},'%[# number of data :] %f %f');
x=reshape((90-polefig.data(:,1)).*sind(polefig.data(:,2)-rotation)...
,[sizes{2},sizes{3}]);
y=reshape((90-polefig.data(:,1)).*cosd(polefig.data(:,2)-rotation)...
,[sizes{2},sizes{3}]);
z=reshape(log10(polefig.data(:,3)),[sizes{2},sizes{3}]);
z(z==-inf)=-1;
contourf(x,y,z,'linestyle','none')
colorbar
caxis([-1,max(max(z))])
colormap('jet')
else
figure(cont)
hold on
end
xa=(90-alpha(index)).*sind(beta(index));
ya=(90-alpha(index)).*cosd(beta(index));
scatter(xa,ya,'LineWidth',1,'MarkerEdgeColor','k')
end
end
end