diff --git a/6_rendering/gtCrystAxisPoles.m b/6_rendering/gtCrystAxisPoles.m
index f7ee0b6b1c4884da6c721984a3cb9f2f7f45e20f..6c674aaf76eec3196267f6ba59f5d5e658d66c2a 100644
--- a/6_rendering/gtCrystAxisPoles.m
+++ b/6_rendering/gtCrystAxisPoles.m
@@ -1,26 +1,33 @@
-function [poles, cmap, valids] = gtCrystAxisPoles(phaseid, pl_norm, r_vectors, axesRot, save_flag)
+function [poles, cmap, valids] = gtCrystAxisPoles(phaseid, pl_norm, varargin)
% GTCRYSTAXISPOLES Makes a colourmap according to a crystallographic axis
-% orientation and gives the colormap.
-% [poles, cmap, valids] = gtCrystAxisPoles([phaseid], [pl_norm], r_vectors, [axesRot], [save_flag])
-% -------------------------------------------------------------------------------------------------
+% orientation and gives the colormap
+% [poles, cmap, valids] = gtCrystAxisPoles([phaseid], [pl_norm], varargin)
+% ------------------------------------------------------------------------
% Calculates poles wrt a crystallographic direction using r_vectors.
-% Loads parameters.mat:cryst(phaseid):spacegroup, latticepar, crystal_system
+% Loads parameters.mat:cryst(phaseid):latticepar, crystal_system
%
% INPUT:
-% phaseid = <int> phase number {1}
-% pl_norm = <double> hkil of the plane normal wrt calculate the poles (1,3(4))
-% {[0 0 0 2]},[2 -1 -1 0] for hexagonal, [0 0 1],[0 1 0] for cubic...
-% r_vectors = <double> Rodrigues vectors (N,3)
-% axesRot = <double> vector to trasform sample reference. It rotates
-% poles basing on component order in this vector {[1 2 3]}
-% save_flag = <logical> flag to save file {false}
+% phaseid = <int> phase number {1}
+% pl_norm = <double> hkil of the plane normal wrt calculate the poles (1,3(4))
+% {[0 0 0 2]},[2 -1 -1 0] for hexagonal, [0 0 1],[0 1 0] for cubic...
+%
+% OPTIONAL INPUT (varargin as a list of pairs, see parse_pv_pairs.m)
+% 'crystal_system' = <string> Legacy to set (force) the crystal system
+% for old datasets
+% 'latticepar' = <double> Lattice parameters to convert from hex to
+% cartesian base
+% 'r_vectors' = <double> Rodrigues vectors as N x 3 array
+% 'axesRot' = <double> vector to trasform sample reference. It rotates
+% poles basing on component order in this vector {[1 2 3]}
+% 'background' = <logical> Keep or remove the background color in first row {true}
+% 'save' = <logical> flag to save file {false}
%
% OUTPUT:
-% poles = <double> list of poles direction of all the grains (M,4)
-% We keep the index of the reference r_vector (first column)
-% NB. M is different from N because of the crystal symmetry.
-% cmap = <double> poles colourmap (N,4)
-% valids = <double> poles list (only one per orientation) (N,3)
+% poles = <double> list of poles direction of all the grains (M,4)
+% We keep the index of the reference r_vector (first column)
+% NB. M is different from N because of the crystal symmetry.
+% cmap = <double> poles colourmap (N,4)
+% valids = <double> poles list (only one per orientation) (N,3)
%
%
% Version 005 18-09-2013 by LNervo
@@ -36,78 +43,103 @@ function [poles, cmap, valids] = gtCrystAxisPoles(phaseid, pl_norm, r_vectors, a
% Version 002 28-06-2012 by LNervo
% Update to version 2 of parameters; cleaning formatting
+% Default parameters
+par.crystal_system = '';
+par.latticepar = [];
+par.r_vectors = [];
+par.axesRot = [1 2 3];
+par.background = true;
+par.save = false;
-if ~exist('r_vectors','var') || isempty(r_vectors)
- gtError('gtCrystAxisPoles:missingArgument','The argument r_vectors is empty...')
-end
+par = parse_pv_pairs(par, varargin);
+
+% Set default phaseid
if ~exist('phaseid','var') || isempty(phaseid)
phaseid = 1;
end
-parameters = [];
-load('parameters.mat');
-crystal_system = parameters.cryst(phaseid).crystal_system;
-latticepar = parameters.cryst(phaseid).latticepar;
-
+% Set default crystal direction if pl_norm empty
if ~exist('pl_norm','var') || isempty(pl_norm)
- if strcmp(crystal_system, 'hexagonal')
- pl_norm = [0 0 0 2];
- else
- pl_norm = [0 0 1];
- end
-end
-if ~exist('axesRot','var') || isempty(axesRot)
- axesRot = [1 2 3];
-end
-if ~exist('save_flag','var') || isempty(save_flag)
- save_flag = false;
+ pl_norm = [0 0 1];
end
-pl_norm_str = num2str(pl_norm);
-pl_norm_str = strrep(pl_norm_str,' ','');
+r_vectors = [];
+phaseDir = fullfile('4_grains', sprintf('phase_%02d', phaseid));
+if ~isempty(par.r_vectors)
+ r_vectors = par.r_vectors;
+elseif exist(fullfile(phaseDir, 'r_vectors.mat'), 'file')
+ load(fullfile(phaseDir, 'r_vectors.mat'));
+elseif exist(fullfile(phaseDir, 'index.mat'), 'file')
+ grain = [];
+ load(fullfile(phaseDir, 'index.mat'), 'grain');
+ r_vectors = gtIndexAllGrainValues(grain, 'R_vector', [], 1, 1:3);
+else
+ gtError('gtCrystAxisPoles:missing_file', 'Could not find r_vectors.mat or index.mat file!');
+end
+% Resize r_vectors if needed
if size(r_vectors,2) == 4
r_vectors = r_vectors(:,2:4);
end
-% take all the symmetric equivalents
-all_hkils = gtCrystSignedHKLs(pl_norm, gtCrystGetSymmetryOperators(crystal_system));
-% finds the original pl_norm among the list of its symmetries
-[~,ind,~] = findDifferentRowIntoMatrix(all_hkils, pl_norm);
-% going to cartesian reciprocal space + normalisation of the cartesian space
-normalised_hkls = gtHex2Cart(all_hkils, latticepar);
+% Get the crystal system
+if isempty(par.crystal_system)
+ parameters = load('parameters.mat');
+ par.crystal_system = parameters.parameters.cryst(phaseid).crystal_system;
+end
+% Get the lattice parameters
+if isempty(par.latticepar)
+ parameters = load('parameters.mat');
+ par.latticepar = parameters.parameters.cryst(phaseid).latticepar;
+end
-poles = zeros(0,4);
-valid_poles = zeros(0,4);
+% convert to Miller-Bravais notation if the case
+if strcmp(par.crystal_system, 'hexagonal') && size(pl_norm,2) == 3
+ pl_norm(4) = pl_norm(3);
+ pl_norm(3) = -pl_norm(1)-pl_norm(2);
+end
% Compute all orientation matrices g
+% Gcry = g Gsam --> Gsam = Gcry inv(g)
+% g = U = gtMathsRod2OriMat(r_vectors.')
all_g = gtMathsRod2OriMat(r_vectors.');
-% Express normalised hkls in cartesian SAMPLE CS
-all_normals = gtVectorCryst2Lab(normalised_hkls, all_g);
+% take all the symmetric equivalents
+pl_norm_symm = gtCrystSignedHKLs(pl_norm, gtCrystGetSymmetryOperators(par.crystal_system));
+% finds the original pl_norm among the list of its symmetries
+[~,ind,~] = findDifferentRowIntoMatrix(pl_norm_symm, pl_norm);
+% going to cartesian reciprocal space + normalisation of the cartesian space
+% PLnorm = pl_norm normalised and symmetrised
+PLnorm = gtHex2Cart(pl_norm_symm, par.latticepar);
+% Express normalised hkls (PLs) in cartesian SAMPLE CS
+PLs = gtVectorCryst2Lab(PLnorm, all_g);
+
+
+poles = zeros(0,4);
+valid_poles = zeros(0,4);
for ii = 1:length(all_g)
- normals = all_normals(:, :, ii);
+ normals = PLs(:, :, ii);
n_bkp = normals;
% rotate
- if axesRot(1) < 0
- normals(:,1) = -n_bkp(:,axesRot(1));
+ if par.axesRot(1) < 0
+ normals(:,1) = -n_bkp(:,par.axesRot(1));
else
- normals(:,1) = n_bkp(:,axesRot(1));
+ normals(:,1) = n_bkp(:,par.axesRot(1));
end
- if axesRot(2) < 0
- normals(:,2) = -n_bkp(:,axesRot(2));
+ if par.axesRot(2) < 0
+ normals(:,2) = -n_bkp(:,par.axesRot(2));
else
- normals(:,2) = n_bkp(:,axesRot(2));
+ normals(:,2) = n_bkp(:,par.axesRot(2));
end
- if axesRot(3) < 0
- normals(:,3) = -n_bkp(:,axesRot(3));
+ if par.axesRot(3) < 0
+ normals(:,3) = -n_bkp(:,par.axesRot(3));
else
- normals(:,3) = n_bkp(:,axesRot(3));
+ normals(:,3) = n_bkp(:,par.axesRot(3));
end
- for jj = 1:size(all_normals, 1)
+ for jj = 1:size(PLs, 1)
% takes the positive poles and saves the grainid which it refers to
if (normals(jj, 3) >= 0)
poles(end+1, 1) = ii;
@@ -130,13 +162,22 @@ clear tmp
if nargout > 1
cmap = gtCrystAxisCmap(valid_poles);
+
+ if ~par.background
+ cmap(1,:) = [];
+ end
+
if nargout > 2
valids = valid_poles;
end
end
-if save_flag
+if par.save
+
+ pl_norm_str = num2str(pl_norm);
+ pl_norm_str = strrep(pl_norm_str,' ','');
+
% saves the list of poles to be used in the pole figure
fname = fullfile('6_rendering',['poles_' pl_norm_str '.mat']);
save(fname,'poles'); % list of the poles