Skip to content
Snippets Groups Projects
Commit aa8d52cf authored by Laura Nervo's avatar Laura Nervo Committed by Nicola Vigano
Browse files

Colour Maps : unified functions gtAaxisCmap and gtCaxisCmap into... 

Colour Maps : unified functions gtAaxisCmap and gtCaxisCmap into gtRandomAxisCmap, giving phaseid, plane normal and r_vectors list.
Pole Figure : created the case 'hex' to draw only the SST for hexagonal crystals....Moving towards the creation of the Inverse Pole Figure for any crystal system....

Signed-off-by: default avatarLaura Nervo <laura.nervo@esrf.fr>

git-svn-id: https://svn.code.sf.net/p/dct/code/trunk@967 4c865b51-4357-4376-afb4-474e03ccb993
parent 03dd9f16
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
6_rendering/gtMakePoleFigure.m
+ 96
− 44
View file @ aa8d52cf
......@@ -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
6_rendering/gtRandomAxisCmap.m 0 → 100644
+ 103
− 0
View file @ aa8d52cf
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
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment