diff --git a/6_rendering/gtMakePoleFigure.m b/6_rendering/gtMakePoleFigure.m
index 2840f493df0273a20316efa8a1af6fabadfb6aa6..25c20d109d78e4a687eef6711cd707013b170faa 100644
--- a/6_rendering/gtMakePoleFigure.m
+++ b/6_rendering/gtMakePoleFigure.m
@@ -21,13 +21,14 @@ function [hfig, save_density, valid_poles_weight] = gtMakePoleFigure(varargin)
% 'poles' = poles from gtReadBoundariesStructure
% [poles, weigths] = gtReadBoundariesStructure(boundaries_structure, vThr, DoSym)
%
-% 'range' = {'sst'}/'phi90'/'phi360'
-% standard stereographic triangle, or quarter or whole hemisphere.
+% 'range' = 'sst'/'phi90'/'phi360'/{'hex'}
+% standard stereographic triangle, quarter or whole
+% hemisphere or hexagonal sst
%
% 'label_poles' = 'all'/'bcc'/'main'/{'none'}
% all hkl poles or just the corners of the sst.
%
-% 'mark_poles' = {'all'}/'bcc'/'main'/'none'
+% 'mark_poles' = 'all'/'bcc'/{'main'}/'none'
%
% 'weights' = can supply a vector the same length as poles to weight the
% pole data {[]}
@@ -38,7 +39,7 @@ function [hfig, save_density, valid_poles_weight] = gtMakePoleFigure(varargin)
%
% 'plot_figure' = can choose not to plot figure, and use the returned
% density values (non normalised) for something
-% {true}
+% {true}/false
%
% Version 002 26-06-2012 by LNervo
% Add phaseid as input argument for multiple phases materials
@@ -54,7 +55,7 @@ appdefaults.poles = [];
appdefaults.density = [];
appdefaults.valid_poles_weight = [];
appdefaults.weights = [];
-appdefaults.range = 'sst'; %or phi90, phi360
+appdefaults.range = 'hex'; %or phi90, phi360, sst
appdefaults.label_poles = 'none';
appdefaults.mark_poles = 'all';
appdefaults.increment = 0.01;
@@ -77,6 +78,7 @@ valid_poles_weight = app.valid_poles_weight;
parameters = [];
load('parameters.mat');
spacegroup = parameters.cryst(phaseid).spacegroup;
+crystal_system = parameters.cryst(phaseid).crystal_system;
if ~isempty(poles)
@@ -86,7 +88,7 @@ if ~isempty(poles)
end
% before building figure, take only those poles that fall in the range 0-phimax,
- % 0-psimax (allowing for seach radius)
+ % 0-psimax (allowing for search radius)
if strcmp(app.range, 'sst')
% for standard stereographic triangle
phimax=pi/4;
@@ -114,6 +116,15 @@ if ~isempty(poles)
weights=weights(dum);
valid_poles_weight=sum(weights(poles(:,3)>=0));
total_solid_angle=2*pi;
+ elseif strcmp(app.range, 'hex')
+ % for hexagonal sst
+ phimax=pi/6;
+ psimax=pi/2;
+ dum=find(poles(:,1)>-sin(searchR) & poles(:,2)>-sin(searchR) & poles(:,3)>cos(psimax)-sin(searchR));
+ poles=poles(dum, :);
+ weights=weights(dum);
+ valid_poles_weight=sum(weights(poles(:,1)>=0 & poles(:,2)>=0 & poles(:,3)>=0));
+ total_solid_angle=pi/12;
else
disp('unrecognised plot range :-(')
return
@@ -167,6 +178,11 @@ else
phimax=2*pi;
psimax=pi/2;
total_solid_angle=2*pi;
+ elseif strcmp(app.range, 'hex')
+ % for standard stereographic triangle
+ phimax=pi/6; % 30 degrees
+ psimax=pi/2;
+ total_solid_angle=pi/12;
else
disp('unrecognised plot range :-(')
return
@@ -183,29 +199,38 @@ if app.plot_figure
% have: density as a function of phi and psi.
% Can convert to x,y,z, and use surf.
% modify for steorographic projection
- phi=[0:inc:phimax phimax]; %around z-axis
- psi=[0:inc:psimax psimax]; %away from z axis
+ phi=[0:inc:phimax phimax]; % around z-axis
+ psi=[0:inc:psimax psimax]; % away from z axis
%radius distorted for stereographic proj
r=tan(psi/2);
[phi,r]=meshgrid(phi,r);
[x,y]=pol2cart(phi,r);
if app.new_figure
- hfig=figure();
+ hfig = figure();
end
+ % plot the poles
surf(x,y,density')
hold on
shading('interp')
+ axis auto
axis equal
+ if strcmp(app.range, 'phi360')
+ app.mark_poles = 'all';
+ elseif strcmp(app.range, 'hex')
+ app.mark_poles = 'main';
+ app.label_poles = 'none';
+ end
+
if ~strcmp(app.mark_poles, 'none') %draw poles
% Get symmetry operators
symm = gtCrystGetSymmetryOperators([], spacegroup);
% cubic materials
- if spacegroup==225 || spacegroup==229
+ if strcmp(crystal_system, 'cubic')
if strcmp(app.mark_poles, 'main')
uvw_poles=[0 0 1; 1 0 1; 1 1 1];
@@ -219,6 +244,7 @@ if app.plot_figure
else
disp('unrecognised "mark_poles" option - should be all/main/none')
end
+
uvw_poles2=[];
for i=1:length(symm)
uvw_poles2=[uvw_poles2; uvw_poles*symm(i).g];
@@ -226,39 +252,30 @@ if app.plot_figure
uvw_poles=[uvw_poles2; -uvw_poles2];
% hexagonal materials
- elseif spacegroup==194 || spacegroup==167 || spacegroup==663
+ elseif strcmp(crystal_system, 'hexagonal')
uvw_poles_hex=[...
- 1 0 -1 0;...
- 0 0 0 1;...
- 0 1 -1 2;...
- 0 1 -1 1;...
- 1 1 -2 3;...
- 1 1 -2 0];
-
- % permute
- uvw_poles_hex2=[];
- for ii = 1:size(uvw_poles_hex,1)
- uvw_poles_hex2 = [uvw_poles_hex2 ; gtGetReflections(uvw_poles_hex(ii,:), spacegroup)];
+ 0 0 0 1;...
+ 2 -1 -1 0;...
+ 1 0 -1 0];
+ if strcmp(app.mark_poles, 'all')
+
+ % permute
+ uvw_poles_hex2=[];
+ for ii = 1:size(uvw_poles_hex,1)
+ uvw_poles_hex2 = [uvw_poles_hex2 ; gtGetReflections(uvw_poles_hex(ii,:), spacegroup)];
+ end
+
+ % convert to cartesian
+ uvw_poles = gtCrystHKL2Cartesian(uvw_poles_hex2', gtCrystHKL2CartesianMatrix(parameters.cryst(phaseid).latticepar))';
end
- uvw_poles_hex=uvw_poles_hex2;
-
- % convert to cartesian
- uvw_poles(:,1)= uvw_poles_hex(:,1) + 0.5 * uvw_poles_hex(:,2);
- uvw_poles(:,2)= 3^0.5/2 *uvw_poles_hex(:,2);
- uvw_poles(:,3)= uvw_poles_hex(:,4);
- % account for cell parameters
- uvw_poles(:,1)=uvw_poles(:,1)*2/(sqrt(3)*parameters.cryst(phaseid).latticepar(1));
- uvw_poles(:,2)=uvw_poles(:,2)*2/(sqrt(3)*parameters.cryst(phaseid).latticepar(1));
- uvw_poles(:,3)=uvw_poles(:,3)/parameters.cryst(phaseid).latticepar(3);
- % add invert through origin
- uvw_poles=[uvw_poles; -uvw_poles];
+ z_uvw = gtCrystHKL2Cartesian(uvw_poles_hex', gtCrystHKL2CartesianMatrix(parameters.cryst(phaseid).latticepar))';
else
- disp('spacegroup not recognised...')
+ disp('crystal system not recognised...')
end
% get the relevent uvw_poles
if strcmp(app.range, 'sst')
- %for standard stereographic triangle
+ %for standard stereographic triangle for cubic
dum=(uvw_poles(:,2)>=0 & uvw_poles(:,2)<=uvw_poles(:,1) & uvw_poles(:,3)>=uvw_poles(:,1));
uvw_poles=uvw_poles(dum, :);
elseif strcmp(app.range, 'phi90')
@@ -268,24 +285,58 @@ if app.plot_figure
elseif strcmp(app.range, 'phi360')
%for 360 pole figure - reduce only using psi
uvw_poles=uvw_poles(uvw_poles(:,3)>=0, :);
+
+ [~,ind_uvw,~]=findDifferentRowIntoMatrix(uvw_poles, z_uvw, false);
+ z_hkl=uvw_poles(ind_uvw,:);
+ [~,ind_hkl,~]=findDifferentRowIntoMatrix(z_uvw, z_hkl, false);
+ elseif strcmp(app.range, 'hex')
+ %for the unitary triangle for hexagonal
+ dum=(uvw_poles(:,1)>=0 & uvw_poles(:,2)<=uvw_poles(:,1)*tand(30) & uvw_poles(:,3)>=0);
+ uvw_poles=uvw_poles(dum, :);
+
+ [~,ind_uvw,~]=findDifferentRowIntoMatrix(uvw_poles, z_uvw, false);
+ z_hkl=uvw_poles(ind_uvw,:);
+ [~,ind_hkl,~]=findDifferentRowIntoMatrix(z_uvw, z_hkl, false);
+ else
+ disp('range not recognised...')
end
% plot these uvw/hkl poles
phi_uvw=atan2(uvw_poles(:,2), uvw_poles(:,1));
psi_uvw=acos(uvw_poles(:,3)./(sqrt(sum(uvw_poles.*uvw_poles, 2))));
r_uvw=tan(psi_uvw/2);
+
dummy=find(r_uvw>1);
r_uvw(dummy)=[];
phi_uvw(dummy)=[];
+ psi_uvw(dummy)=[];
+ uvw_poles(dummy,:)=[];
+
[x_uvw,y_uvw]=pol2cart(phi_uvw, r_uvw);
dummy=find(x_uvw<-1 | x_uvw>1 | y_uvw<-1 | y_uvw>1);
x_uvw(dummy)=[];
y_uvw(dummy)=[];
-
+ uvw_poles(dummy,:)=[];
+
+ % drawing the poles
plot3(x_uvw, y_uvw, ones(size(x_uvw))*(max(density(:))+1), 'k*')
+ if strcmp(app.range, 'hex')
+ x_labels = x_uvw(ind_uvw(ind_hkl))+0.01;
+ y_labels = y_uvw(ind_uvw(ind_hkl))-0.02;
+ alignment_h = {'left','right','left'};
+ alignment_v = {'top','top','bottom'};
+
+ % drawing the labels
+ for ii=1:size(z_uvw,1)
+ text(x_labels(ii), y_labels(ii), ['[' num2str(uvw_poles_hex(ii,:)) ']'],...
+ 'HorizontalAlignment',alignment_h{ii},...
+ 'VerticalAlignment', alignment_v{ii})
+ end
+ end
end
+
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%% plot rings on phi360 figures
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
@@ -299,8 +350,8 @@ if app.plot_figure
ring_z=ones(size(ring_phi))*(max(density(:))+1);
plot3(ring_x, ring_y, ring_z, 'w-');
% label rings
- %h=text(ring_x(1)+0.01, 0, (max(density(:))+1), num2str(ring));
- %set(h, 'color', 'w');
+ h=text(ring_x(1)+0.01, 0, (max(density(:))+1), num2str(ring));
+ set(h, 'color', 'w');
end
end
@@ -311,10 +362,11 @@ if app.plot_figure
set(gca, 'XTickLabel','');
set(gca, 'YTickLabel','');
set(gca, 'GridLineStyle','none');
- set(gca, 'Ycolor', 'w');
set(gca, 'Xcolor', 'w');
- colorbar();
-
+ set(gca, 'Ycolor', 'w');
+ c=colorbar();
+ set(c, 'Location', 'EastOutside')
+
%change to x-y view
set(gca, 'CameraPosition', [(min(x(:))+max(x(:)))/2 (min(y(:))+max(y(:)))/2 5*(max(density(:))+1)]);
set(gca, 'CameraUpVector', [0 1 0]);
@@ -365,7 +417,7 @@ if app.plot_figure
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% label poles
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
- if ~strcmp(app.label_poles, 'none')
+ if ~strcmp(app.label_poles, 'none') && strcmp(crystal_system, 'cubic')
label_z=max(density(:))+1;
@@ -408,7 +460,7 @@ if app.plot_figure
end
else
- hfig=[];
+ hfig = [];
end
end % end of function
diff --git a/6_rendering/gtRandomAxisCmap.m b/6_rendering/gtRandomAxisCmap.m
new file mode 100644
index 0000000000000000000000000000000000000000..1124e58289adce9d90adea776f5581fe8d0056b3
--- /dev/null
+++ b/6_rendering/gtRandomAxisCmap.m
@@ -0,0 +1,103 @@
+function [axis, axismap] = gtRandomAxisCmap(phaseid, pl_norm, r_vectors)
+% GTRANDOMAXISCMAP Make a colourmap according to an axis orientation
+% [axis, axismap] = gtRandomAxisCmap(phaseid, pl_norm, r_vectors)
+% ---------------------------------------------------------------
+% Loads parameters.mat
+% Saves 6_rendering/axis_[pl_norm].mat and 6_rendering/axismap_[pl_norm].mat
+%
+% INPUT:
+% phaseid = phase number <int> {1}
+% pl_norm = hkil of the plane normal wrt calculate the poles <double 1x4>
+% {c_axis = [0 0 0 1], a_axis = [1 0 -1 0], ...}
+% r_vectors = r_vectors <double Nx4>
+%
+% OUTPUT:
+% axis = list of axis direction of all the grains
+% axismap = axis colourmap calculated from r_vectors
+%
+%
+% r_vectors -> axis -> phi/psi -> hsl -> rgb colormap
+% use psi to give brightness; phi to give colour
+%
+% Version 003 19-10-2012 by YGuilhem
+% Factorized using gtVectorOrientationColor
+%
+% Version 002 28-06-2012 by LNervo
+% Update to version 2 of parameters; cleaning formatting
+
+
+if ~exist('phaseid','var') || isempty(phaseid)
+ phaseid = 1;
+end
+if ~exist('pl_norm','var') || isempty(pl_norm)
+ pl_norm = [0 0 0 1];
+end
+
+pl_norm_str = num2str(pl_norm);
+pl_norm_str = strrep(pl_norm_str,' ','');
+
+parameters = [];
+load('parameters.mat');
+
+if size(r_vectors,2) == 3
+ r_vectors(:,2:4) = r_vectors;
+ r_vectors(:,1) = 1:size(r_vectors,1);
+end
+
+all_hkils=[];
+for i = 1:size(pl_norm,1)
+ all_hkils = [all_hkils ; gtGetReflections(pl_norm(i,:), parameters.cryst(phaseid).spacegroup)];
+end
+
+% going to cartesian reciprocal space + normalisation of the cartesian
+% space
+all_hkls(:, 1) = all_hkils(:, 1) + 0.5 * all_hkils(:, 2);
+all_hkls(:, 2) = sqrt(3)/2 * all_hkils(:, 2);
+all_hkls(:, 3) = all_hkils(:, 4);
+
+all_hkls(:, 1) = all_hkls(:, 1)*2/(sqrt(3)*parameters.cryst(phaseid).latticepar(1));
+all_hkls(:, 2) = all_hkls(:, 2)*2/(sqrt(3)*parameters.cryst(phaseid).latticepar(1));
+all_hkls(:, 3) = all_hkls(:, 3)/parameters.cryst(phaseid).latticepar(3);
+
+% normalise hkls vector
+tmp = sqrt(sum((all_hkls.*all_hkls),2));
+normalised_hkls = all_hkls./(repmat(tmp,1,3));
+%%%%% can be replaced with gtHexToCart
+
+
+axis = zeros(size(r_vectors, 1), 4);
+axis(:, 1) = r_vectors(:, 1);
+for i=1:size(r_vectors,1)
+ crys2sam = Rod2g(r_vectors(i,2:4));
+ all_normals = (crys2sam * normalised_hkls')';
+ if all_normals(2,3)>=0
+ axis(i,2)=all_normals(2,1);
+ axis(i,3)=all_normals(2,2);
+ axis(i,4)=all_normals(2,3);
+ else
+ axis(i,2)=all_normals(1,1);
+ axis(i,3)=all_normals(1,2);
+ axis(i,4)=all_normals(1,3);
+ end
+end
+
+
+fname = fullfile('6_rendering',['axis_' pl_norm_str '.mat']);
+save(fname,'axis'); % list of the axis
+disp(['Saved axis list in ' fname])
+
+% Get the orientation color defined by HSL color code, translated to RGB
+rgb = gtVectorOrientationColor(axis(:, 2:4));
+
+% Add black background and expand list length of max(grainid)+1
+axismap = zeros(max(axis(:, 1))+1, 3);
+for i=1:length(axis)
+ axismap(axis(i,1)+1, :) = rgb(i, :);
+end
+
+% Save the axis colourmap
+fname = fullfile('6_rendering',['axismap_' pl_norm_str '.mat']);
+save(fname,'axismap');
+disp(['Saved colour map in ' fname])
+
+end % end of function