Skip to content
Snippets Groups Projects
Commit d6c1b0c7 authored by Laura Nervo's avatar Laura Nervo
Browse files

gtDrawGrainUnitCell : vectorize function to plot multiple datasets. 

gtPlotGrainUnitCell : better commented and optimized functionalities
gtComputeGrainUnitCell : added to compute the unit cell basing on grain center, size, strain, orientation

Signed-off-by: default avatarLaura Nervo <laura.nervo@esrf.fr>
parent 923bcef9
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
6_rendering/gtComputeGrainUnitCell.m 0 → 100644
+ 43
− 0
View file @ d6c1b0c7
function vertices = gtComputeGrainUnitCell(rvec, vertices, ga, gc, strainT, scale)
% GTCOMPUTEGRAINUNITCELL
% vertices = gtComputeGrainUnitCell(rvec, vertices, [ga], [gc], [strainT], [scale])
% ---------------------------------------------------------------------------------
% Trasforms reference unti cell in grain unit cell according to the grain center,
% grain orientation, grain size, strain tensor and strain scaling factor.
%
% INPUT:
% rvec = Rodriguez vector <double 1x3>
% vertices = original vertices of the unit cell <Nx3>
% ga = size of the grain <double>
% gc = grain center <double 1x3> (mm)
% strainT = strain tensor <double 3x3>
% scale = scale factor for strain <double>
%
% OUTPUT:
% vertices = modified vertices of the unit cell <Nx3>
%
% Version 002 06-06-2013 by LNervo
% compute orientation matrix
g_grain = gtMathsRod2OriMat(rvec');
% compute positions in the Lab ref system
vertices = gtVectorCryst2Lab(vertices, g_grain);
% apply grain size
if exist('ga','var') && ~isempty(ga) && ga ~= 0
vertices = ga*vertices;
end
% apply grain center
if exist('gc','var') && ~isempty(gc)
vertices = repmat(gc,size(vertices,1),1) + (eye(3)*vertices')';
end
% apply grain strain
if exist('strainT','var') && ~isempty(strainT) && ~any(isnan(strainT(:)))
if exist('scale','var') && ~isempty(scale) && scale ~= 0
strainT = strainT*scale;
end
vertices = vertices + (strainT*vertices')';
end
end % end of function
6_rendering/gtDrawGrainUnitCells.m 100755 → 100644
+ 170
− 185
View file @ d6c1b0c7
function [f_handle, cb_handle] = gtDrawGrainUnitCells(grains, varargin)
% [f_handle, cb_handle] = gtDrawGrainUnitCells(grains, varargin)
% --------------------------------------------------------------
function [f_handle, p_handle] = gtDrawGrainUnitCells(grains, varargin)
% [f_handle, p_handle] = gtDrawGrainUnitCells(grains, varargin)
% -------------------------------------------------------------
% Generalized cases to plot more than one dataset.
%
% Draws a figure representing grain size, location and orientation by
% Draws a unit cell representing grain size, location and orientation by
% cubes or hexagonal prismes based on the their Rodrigues vectors.
%
% INPUT:
% grains = grain structure from indexing <cell>
% grains = <cell array> grain cell-structure from indexing
%
% OPTIONAL INPUT (varargin):
% grainids = <double> list of grain IDs of interest
% hf = figure handle {}
% ids = <double> list of grain IDs of interest; by default all
% the grains are drawn
% cmap = <double> can be a component of the strain tensor
% {grain.strain.strainT(3,3)} or a colour map
% latticepar = <double 1x6> lattice parameters
% labgeo = <struct> laboratory geometry as in parameters
% phaseid = <int> phase number {1}
% type = <string> unit cell type {'cubic'}/'hexagonal'
% hlight = <double> grain id-s to be highlighted {}
% scale = <double> number to scale grain size; empty to use bbox info in grain {}
% view = <double 1x2> 3D angular view {[0 90]}
% showbar = <logical> flag to show the colourbar {false}
% colormap = <string> type of the colormap 'jet',{'cool'},'gray','autumn','hsv',...
% only if it is not with true colours (RGB)
% If needed, one can set the number of
% levels: i.e. 'colormap','jet(10)'
% hlight = <double> grain id-s to be highlighted {[]}
% linecolor = <string> color for edges
% strain = <double> scale factor for the strain values, if wanted {0}
% zero to use the unstrained unit cell
% strainT = <double> component of the strainT tensor in 'strain'
% {[3 3]}
% (zero to use the unstrained unit cell)
% alpha = <double> transparency value for patches {1}
% pxsize = <double> pixel size (mm/px); if empty, taken from
% parameters.labgeo
% ratio = <double> c/a ratio; if empty, taken from
% parameters.cryst
% phaseid = <int> phase number {1}
% type = <string> unit cell type; if empty, taken from
% parameters.cryst
% label = <string> name for entry in the legend <string>
%
% hf = figure handle if already existing {[]}
% zoom = <double> zoom in of this quantity {1.2}
% view = <double> 3D angular view {[0 90]} (1x2)
% legend = <logical> display or not the legend {false}
% legendpos = <double> location of the legend box {'NorthEast'}
%
% OUTPUT:
% f_handle = handle of the figure
% cb_handle = handle of the colorbar
% p_handle = handle of the patches
%
%
% Usage:
% gtDrawGrainUnitCells({grain}, 'ids', {1:10}, 'cmap', {gtRandCmap(10)}, 'pxsize', {0.0014})
%
%
% Version 006 06-06-2013 by LNervo
% Added function to compute unit cell vertices, modified
% gtPlotGrainUnitCell
%
% Version 005 11-03-2013 by LNervo
% Use cell arrays to plot multiple datasets
%
% Version 004 15-11-2012 by LNervo
% Use gtDrawSampleGeometry(labgeo) to create the graphics with axis and reference
% system
% Use gtDrawSampleGeometry(labgeo) to create the graphics with axis and reference
% system
%
% Version 003 23-10-2012 by LNervo
% Add varargin
% modify parameters/global settings
app.grainids = [];
app.hf = [];
app.cmap = [];
app.latticepar = [];
app.labgeo = [];
app.phaseid = 1;
app.type = [];
app.hlight = [];
app.scale = [];
app.view = [0 90];
app.showbar = false;
app.colormap = 'cool';
app.strain = 0;
app.strainT = [3 3];
app.zoom = 1.2;
% Add varargin
app.ids = arrayfun(@(num) 1:numel(grains{num}), 1:numel(grains), 'UniformOutput', false);
app.cmap{1} = [];
app.hlight{1} = [];
app.linecolor = {'k','r','y','g','b'};
app.strain = cellfun(@(num) 0, grains, 'UniformOutput', false);
app.alpha = cellfun(@(num) 1, grains, 'UniformOutput', false);
app.pxsize{1} = [];
app.ratio{1} = [];
app.phaseid = cellfun(@(num) 1, grains, 'UniformOutput', false);
app.type{1} = [];
app.label{1} = [];
app.hf = [];
app.zoom = 1.2;
app.view = [0 90];
app.legend = false;
app.legendpos = 'NorthEast';
app = parse_pv_pairs(app, varargin);
if ~exist('grains','var') || isempty(grains)
gtError('gtDrawGrainUnitCells:missingArgument','grain is not given as input argument...Quitting')
end
if isempty(app.grainids)
app.grainids = 1:length(grains);
end
if isempty(app.cmap)
if isfield(grains{1},'strain')
app.cmap = gtIndexAllGrainValues(grains,'strain','strainT',app.strainT(1),app.strainT(2));
disp(['Using strainT(' num2str(app.strainT(1)) ',' num2str(app.strainT(2)) ') values as colorscale'])
app.showbar = true;
else
app.cmap = zeros(length(grains),1);
app.showbar = false;
end
else
if size(app.cmap,2)==1
app.showbar = true;
else
app.showbar = false;
end
if size(app.cmap,1)==1 && size(app.cmap,2)>1
app.cmap = repmat(app.cmap,length(grains),1);
end
% remove bkg color
if all(app.cmap(1,:) == 0)
app.cmap(1,:) = [];
end
end
grains = grains(app.grainids);
app.cmap = app.cmap(app.grainids,:);
app.linecolor = app.linecolor(1:numel(grains));
parameters = [];
load('parameters.mat');
if isempty(app.latticepar)
app.latticepar = parameters.cryst(app.phaseid).latticepar;
disp('load latticepar from parameters.mat')
end
if isempty(app.labgeo)
app.labgeo = parameters.labgeo;
disp('load labgeo from parameters.mat')
end
if isempty(app.type)
app.type = parameters.cryst(app.phaseid).crystal_system;
disp('load crystal_system from parameters.mat')
if ~strcmpi(app.type,'cubic') && ~strcmpi(app.type, 'hexagonal')
app.type = 'cubic';
disp('Default shape for unit cells: cubic')
end
end
labgeo = parameters.labgeo;
if isempty(app.pxsize{1})
app.pxsize{1} = (labgeo.pixelsizeu + labgeo.pixelsizev)/2;
end
if isempty(app.hf)
hf = gtDrawSampleGeometry(app.labgeo.samenvrad, app.labgeo.samenvbot, app.labgeo.samenvtop);
hf = gtDrawSampleGeometry(labgeo.samenvrad, labgeo.samenvbot, labgeo.samenvtop);
else
hf = app.hf;
end
hold on
nof_grains = length(grains);
pixelsize = (app.labgeo.pixelsizeu + app.labgeo.pixelsizev)/2;
disp('Using pixelsize and samenv* from parameters.mat:labgeo')
if app.showbar
colormap(app.colormap)
disp(['Using colormap ' app.colormap '...'])
mmin = min(app.cmap);
mmax = max(app.cmap);
ticks = mmin:(mmax-mmin)/10:mmax;
cb = colorbar('FontSize',10,...
'YMinorTick','on');
end
set(gcf, 'Visible', 'off');
% vertices and faces calculation
if strcmpi(app.type, 'cubic')
data = gtCubicUnitCell('ratio',app.latticepar(3)/app.latticepar(1),'draw',false);
vertices = data.vertices(1:8,:);
faces_sides = data.faces_sides;
faces_end = data.faces_end;
hkl = data.hkl;
elseif strcmpi(app.type, 'hexagonal')
data = gtHexagonalUnitCell('ratio',app.latticepar(3)/app.latticepar(1),'draw',false);
vertices = data.vertices(1:12,:);
faces_sides = data.faces_sides;
faces_end = data.faces_end;
hkl = data.hkl;
end
hp = [];
[Bmat, Amat] = gtCrystHKL2CartesianMatrix(app.latticepar);
pl_norm = gtCrystHKL2Cartesian(hkl, Bmat); % hkl plane normals in cryst reference
% loop over datasets
for ii=1:numel(grains)
if app.strain~=0
disp('Plotting strained unit cells for grains with strainT tensor...')
disp(['Strain values are scaled by ' num2str(app.strain)])
end
if ii>1
if numel(app.type)==1, app.type{ii} = app.type{1}; end
if numel(app.ratio)==1, app.ratio{ii} = app.ratio{1}; end
if numel(app.pxsize)==1, app.pxsize{ii} = app.pxsize{1}; end
if numel(app.cmap)==1, app.cmap{ii} = app.cmap{1}; end
if numel(app.hlight)==1, app.hlight{ii} = app.hlight{1}; end
if numel(app.label)==1, app.label{ii} = app.label{1}; end
end
% draw unit cells
for ii = 1:nof_grains
if ismember(ii, app.hlight)
hl = true;
else
hl = false;
cryst = parameters.cryst(app.phaseid{ii});
if isempty(app.type{ii}), app.type{ii} = cryst.crystal_system; end
if isempty(app.ratio{ii}), app.ratio{ii} = cryst.latticepar(3)/cryst.latticepar(1); end
if isempty(app.ids{ii}), app.ids{ii} = 1:length(grains{ii}); end
if isempty(app.cmap{ii}), app.cmap{ii} = zeros(length(grains{ii}),1); end
ids = app.ids{ii};
grain = grains{ii};
cmap = app.cmap{ii};
% if single value (RGB), copy it for all the grains
if size(cmap,1)==1 && size(cmap,2)>1
cmap = repmat(cmap,length(grain),1);
end
sfPlotGrainCell(grains{ii}, vertices, faces_sides, hl, pixelsize, app.cmap(ii,:), app.scale, app.strain);
sfPlotGrainCell(grains{ii}, vertices, faces_end, hl, pixelsize, app.cmap(ii,:), app.scale, app.strain);
end
% remove bkg color
if all(cmap(1,:) == [0 0 0])
cmap = cmap(2:end,:);
end
% get selected grains for current datasets
grain = grain(ids);
cmap = cmap(ids,:);
% vertices and faces calculation : TO DO LIST : unify in one function
if strcmpi(app.type{ii}, 'cubic')
data{ii} = gtCubicUnitCell('ratio',app.ratio{ii},'draw',false);
elseif strcmpi(app.type{ii}, 'hexagonal')
data{ii} = gtHexagonalUnitCell('ratio',app.ratio{ii},'draw',false);
end
% draw unit cells
hold(gca,'on');
p = [];
for jj = 1:length(grain)
% grain size
if isfield(grain{jj}, 'stat')
if isfield(grain{jj}.stat, 'size_int')
grainsize = 0.5*app.pxsize{ii}*grain{jj}.stat.size_int;
elseif isfield(grain{jj}.stat, 'bbxsmean') && isfield(grain{jj}.stat, 'bbysmean')
grainsize = 0.5*app.pxsize{ii}*(grain{jj}.stat.bbxsmean + grain{jj}.stat.bbysmean)/2;
end
else
grainsize = 10;
end
% grain strain
if ~isfield(grain{jj}, 'strain')
strainT = NaN(3);
else
strainT = grain{jj}.strain.strainT;
end
% highlight grain?
if ismember(jj, app.hlight{ii}), hl = true; else hl = false; end
% computes vertices and caxis coordinates
vertices{jj} = gtComputeGrainUnitCell(grain{jj}.R_vector, data{ii}.vertices, ...
grainsize, grain{jj}.center, strainT, app.strain{ii});
caxis{jj} = gtComputeGrainUnitCell(grain{jj}.R_vector, data{ii}.caxis, ...
grainsize, grain{jj}.center, strainT, app.strain{ii});
% draws patches for prismatic planes
p{end+1} = gtPlotGrainUnitCell(...
vertices{jj}, caxis{jj}, data{ii}.faces_sides, ...
cmap(jj,:), hl, app.linecolor{ii}, app.alpha{ii});
set(p{end}, 'UserData', [ii ids(jj)]);
set(p{end}, 'Tag', sprintf('data_%02d_grain_%04d',ii,grain{jj}.id));
% draws patches for basal planes
p{end+1} = gtPlotGrainUnitCell(...
vertices{jj}, caxis{jj}, data{ii}.faces_end, ...
cmap(jj,:), hl, app.linecolor{ii}, app.alpha{ii});
set(p{end}, 'UserData', [ii ids(jj)]);
set(p{end}, 'Tag', sprintf('data_%02d_grain_%04d',ii,grain{jj}.id));
end % end for jj % loop over grains
app.cmap{ii} = cmap;
hp{ii} = p;
end % end for ii % loop over datasets
app.label = app.type;
axis equal;
zoom('out');
......@@ -177,57 +200,19 @@ zoom(app.zoom)
zoom('off')
view(app.view);
f_handle = gcf;
if nargout == 2
cb_handle = cb;
end
end % of main function
%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% SUB-FUNCTIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function sfPlotGrainCell(grain, vertices, faces, hl, pixelsize, cmap, scale, strainscale)
if isempty(scale)
ga = 0.5*pixelsize*(grain.stat.bbxsmean + grain.stat.bbysmean)/2;
else
ga = 0.5*pixelsize*(grain.stat.bbxsmean + grain.stat.bbysmean)/2*scale;
if app.legend
leg = legend('Parent',gcf,'String',app.label,'Location',app.legendpos);
end
if isfield(grain,'centre')
gc = grain.centre;
else
gc = grain.center;
end
% Get or recompute orientation matrix g
if isfield(grain, 'g')
g_grain = grain.g;
else
g_grain = gtMathsRod2OriMat(grain.R_vector);
end
% Express vertices in cartesian SAMPLE CS
grain_vertices = ga * gtVectorCryst2Lab(vertices, g_grain);
grain_vertices_st = repmat(gc,size(vertices,1),1) + (eye(3)*grain_vertices')';
set(gcf, 'Visible', 'on');
f_handle = gcf;
app.hf = f_handle;
if isfield(grain,'strain') && ~any(isnan(grain.strain.strainT(:))) && strainscale ~= 0
grain_vertices_st = grain_vertices_st + (grain.strain.strainT*strainscale*grain_vertices')';
if nargout == 2
p_handle = hp;
app.hp = hp;
end
facecolor = repmat(cmap,size(faces,1),1);
if hl
patch('Vertices', grain_vertices_st, 'Faces', faces , ...
'FaceVertexCData', facecolor, 'FaceColor', 'flat','EdgeColor','r','LineWidth',2);
else
patch('Vertices', grain_vertices_st, 'Faces', faces , ...
'FaceVertexCData', facecolor, 'FaceColor', 'flat');
end
print_structure(app, 'general settings:')
end % of sfPlotGrainCell
end % end of function
6_rendering/gtPlotGrainUnitCell.m
+ 36
− 28
View file @ d6c1b0c7
function hp = gtPlotGrainUnitCell(gc, rvec, strainT, vertices, faces, hl, ga, cmap, linecolor, strainscale, alpha)
g_grain = Rod2g(rvec);
grain_vertices = ga*vertices*g_grain';
grain_vertices_st = repmat(gc,size(vertices,1),1) + (eye(3)*grain_vertices')';
if ~any(isnan(strainT(:))) && strainscale ~= 0
grain_vertices_st = grain_vertices_st + (strainT*strainscale*grain_vertices')';
end
function [hp,hl] = gtPlotGrainUnitCell(vertices, caxis, faces, cmap, hlight, linecolor, alpha)
% GTPLOTGRAINUNITCELL Draws the grain unit cell given its properties
%
% [hp,hl] = gtPlotGrainUnitCell(vertices, caxis, faces, cmap, hlight, linecolor, alpha)
% -------------------------------------------------------------------------------------
% INPUT:
% vertices = vertices coordinates <double Nx3>
% caxis = coordinates of the c-axis <double 2x3>
% faces = faces <double>
% cmap = facecolor property <double 1x3>
% hlight = flag to highlight the unit cell (thicker edges) <logical>
% linecolor = linecolor property (edges)
% alpha = facealpha property <double>
%
% OUTPUT:
% hp = patch handle
% hl = c-axis handle
%
% Version 002 06-06-2013 by LNervo
facecolor = repmat(cmap,size(faces,1),1);
% draws the patches
hp = patch('vertices', vertices, 'faces', faces, ...
'FaceVertexCData', facecolor, 'FaceColor', 'flat', ...
'EdgeColor', linecolor, 'FaceAlpha', alpha);
% draw the caxis
if ~isempty(caxis)
hl = line(caxis(:,1),caxis(:,2),caxis(:,3),'Color','k');
else
hl = [];
end
if hl
hp = patch('vertices', grain_vertices_st, 'faces', faces, ...
'FaceVertexCData',facecolor,'FaceColor','flat','EdgeColor',linecolor,...
'LineWidth',2,'FaceAlpha',alpha);
if hlight
set(hp, 'LineWidth',2)
if hlight, set(findobj(gcf,'type','line'), 'LineWidth',2); end
else
hp = patch('vertices', grain_vertices_st, 'faces', faces, ...
'FaceVertexCData',facecolor,'FaceColor','flat','EdgeColor',linecolor,...
'LineWidth',1,'FaceAlpha',alpha);
set(hp, 'LineWidth',1)
if hlight, set(findobj(gcf,'type','line'), 'LineWidth',1); end
end
%
% if hl
% hp = drawMesh(grain_vertices_st, faces, ...
% 'FaceVertexCData',facecolor,'FaceColor','flat','EdgeColor',linecolor,...
% 'LineWidth',2,'FaceAlpha',alpha);
% else
% hp = drawMesh(grain_vertices_st, faces, ...
% 'FaceVertexCData',facecolor,'FaceColor','flat','EdgeColor',linecolor,...
% 'LineWidth',1,'FaceAlpha',alpha);
% end
end % of gtPlotGrainUnitCell
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