Small MATLAB script

Unless otherwise stated, this page is written by DTU Nanolab internal

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