Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Feature/sort segments #2247

Open
wants to merge 4 commits into
base: develop
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
131 changes: 88 additions & 43 deletions EBSDAnalysis/@EBSD/calcGrains.m
Original file line number Diff line number Diff line change
Expand Up @@ -62,6 +62,7 @@
% See also
% GrainReconstruction GrainReconstructionAdvanced

% (1)
% subdivide the domain into cells according to the measurement locations,
% i.e. by Voronoi tessellation or unit cell
if isa(ebsd,'EBSDsquare') || isa(ebsd,'EBSDhex')
Expand All @@ -75,6 +76,7 @@
% D - cell array of cells
% I_FD - incidence matrix faces to vertices

% (2)
% determine which cells to connect
[A_Db,I_DG] = doSegmentation(I_FD,ebsd,varargin{:});
% A_db - neighboring cells with (inner) grain boundary
Expand All @@ -94,15 +96,58 @@
spdiags(ebsd.phaseId,0,length(ebsd),length(ebsd)),[],2));
phaseId(phaseId==0) = 1; % why this is needed?

% (3)
% compute boundary this gives
% I_FDext - faces x cells external grain boundaries
% I_FDint - faces x cells internal grain boundaries
[I_FDext, I_FDint, Fext, Fint] = calcBoundary;
% [I_FDext, I_FDint, Fext, Fint] = calcBoundary;

% TODO:
% should return
[isExt, isInt] = calcBoundary(A_Db, I_DG, I_FD);

isIntOrExt = isInt | isExt;

% reduce F
F = F(isIntOrExt,:);
I_FD = I_FD(isIntOrExt,:);
isInt = isInt(isIntOrExt);
isExt = isExt(isIntOrExt);

% reduce V
[inUse,~,F] = unique(F(:));
V = V(inUse,:);
F = reshape(F,[],2);

%(4)
% detect and sort segments
% sP - list of segment points
% s - starting point of segements in sP
% sfF - indices to F of the new face order
[s, sP, sf, sfF, sfD] = TriplepointSegments(F,V);


% * reorder F
% * split F into Fext and Fint
% * reorder and reduce I_FD and then split into I_FDext and I_FD_int

F = F(sfF,:);
isExt = isExt(sfF);
isInt = isInt(sfF);
F(sfD == 0,[1,2]) = F(sfD == 0,[2,1]);% flip direction
I_FD = I_FD(sfF,:);

Fext = F(isExt,:); % external boundary segments
Fint = F(isInt,:); % internal boundary segments
I_FDext = I_FD(isExt,:);
I_FDint = I_FD(isInt,:);

% does this still work? TODO
if check_option(varargin,'removeQuadruplePoints')
qAdded = removeQuadruplePoints;
end


% setup grains
try
poly = calcPolygonsC(I_FDext * I_DG,Fext,V);
Expand Down Expand Up @@ -258,48 +303,6 @@

end

function [I_FDext, I_FDint, Fext, Fint] = calcBoundary
% distinguish between interior and exterior grain boundaries

% cells that have a subgrain boundary, i.e. a boundary with a cell
% belonging to the same grain
sub = ((A_Db * I_DG) & I_DG)'; % grains x cell
[i,j] = find( diag(any(sub,1))*double(A_Db) ); % all adjacent to those
sub = any(sub(:,i) & sub(:,j),1); % pairs in a grain

% split grain boundaries A_Db into interior and exterior
A_Db_int = sparse(i(sub),j(sub),1,size(I_DG,1),size(I_DG,1));
A_Db_ext = A_Db - A_Db_int; % adjacent over grain boundary

% create incidence graphs
I_FDbg = diag( sum(I_FD,2)==1 ) * I_FD;
D_Fbg = diag(any(I_FDbg,2));

[ix,iy] = find(A_Db_ext);
D_Fext = diag(sum(abs(I_FD(:,ix)) & abs(I_FD(:,iy)),2)>0);

I_FDext = (D_Fext| D_Fbg)*I_FD;

[ix,iy] = find(A_Db_int);
D_Fsub = diag(sum(abs(I_FD(:,ix)) & abs(I_FD(:,iy)),2)>0);
I_FDint = D_Fsub*I_FD;

% remove empty lines from I_FD, F, and V
isExt = full(any(I_FDext,2));
I_FDext = I_FDext.'; I_FDext = I_FDext(:,isExt).';

isInt = full(any(I_FDint,2));
I_FDint = I_FDint.'; I_FDint = I_FDint(:,isInt).';

% remove vertices that are not needed anymore
[inUse,~,F] = unique(F(isExt | isInt,:));
V = V(inUse,:);
F = reshape(F,[],2);
Fext = F(isExt(isExt | isInt),:); % external boundary segments
Fint = F(isInt(isExt | isInt),:); % internal boundary segments

end

function gB = makeBoundary(F,I_FD)

% compute ebsdInd
Expand Down Expand Up @@ -449,6 +452,48 @@

end

function [isExt, isInt] = calcBoundary(A_Db, I_DG, I_FD)
% distinguish between interior and exterior grain boundaries

% cells that have a subgrain boundary, i.e. a boundary with a cell
% belonging to the same grain
sub = ((A_Db * I_DG) & I_DG)'; % grains x cell
[i,j] = find( diag(any(sub,1))*double(A_Db) ); % all adjacent to those
sub = any(sub(:,i) & sub(:,j),1); % pairs in a grain

% split grain boundaries A_Db into interior and exterior
A_Db_int = sparse(i(sub),j(sub),1,size(I_DG,1),size(I_DG,1));
A_Db_ext = A_Db - A_Db_int; % adjacent over grain boundary

% create incidence graphs
I_FDbg = diag( sum(I_FD,2)==1 ) * I_FD;
D_Fbg = diag(any(I_FDbg,2));

[ix,iy] = find(A_Db_ext);
D_Fext = diag(sum(abs(I_FD(:,ix)) & abs(I_FD(:,iy)),2)>0);

I_FDext = (D_Fext| D_Fbg)*I_FD;

[ix,iy] = find(A_Db_int);
D_Fsub = diag(sum(abs(I_FD(:,ix)) & abs(I_FD(:,iy)),2)>0);
I_FDint = D_Fsub*I_FD;

% remove empty lines from I_FD, F, and V
isExt = full(any(I_FDext,2));
%I_FDext = I_FDext.'; I_FDext = I_FDext(:,isExt).';

isInt = full(any(I_FDint,2));
%I_FDint = I_FDint.'; I_FDint = I_FDint(:,isInt).';

% remove vertices that are not needed anymore
%[inUse,~,F] = unique(F(isExt | isInt,:));
%V = V(inUse,:);
%F = reshape(F,[],2);
%Fext = F(isExt(isExt | isInt),:); % external boundary segments
%Fint = F(isInt(isExt | isInt),:); % internal boundary segments

end

function poly = calcPolygons(I_FG,F,V)
%
% Input:
Expand Down
Binary file added mex/TriplepointSegments.mexw64
Binary file not shown.
4 changes: 2 additions & 2 deletions mex/mex_install.m
Original file line number Diff line number Diff line change
Expand Up @@ -13,8 +13,8 @@ function mex_install(varargin)
'extern/jcvoronoi/',...
'extern/libDirectional/numerical',...
'SO3Fun/@SO3FunHarmonic/private/wignerTrafo',...
'tools/graph_tools/EulerCyclesC'
};
'tools/graph_tools/EulerCyclesC',...
'tools/graph_tools/TriplepointSegments'};

% TODO: Check for mex-Compiler

Expand Down
Loading