Merge pull request #128 from ralfHielscher/master

MTEX 4.2.0

EBSDAnalysis/@EBSD/EBSD.m

@@ -51,9 +51,9 @@
 
       if nargin == 0, return; end            
 
-      ebsd.rotations = rotation(rot(:));
+      ebsd.rotations = rotation(rot);
       ebsd = ebsd.init(phases,CSList);      
-      ebsd.id = (1:length(phases)).';
+      ebsd.id = (1:numel(phases)).';
 
       % extract additional properties
       ebsd.prop = get_option(varargin,'options',struct);

EBSDAnalysis/@EBSD/KAM.m

@@ -1,5 +1,5 @@
 function kam = KAM(ebsd,varargin)
-% intergranular average misorientation angle per orientation
+% intragranular average misorientation angle per orientation
 %
 % Syntax
 %

EBSDAnalysis/@EBSD/calcGrains.m

@@ -76,7 +76,7 @@
 % save 
 grains.prop.GOS = GOS;
 grains.prop.meanRotation = meanRotation;
-mis2mean = inv(rotation(q)) .* meanRotation(grainId(:));
+mis2mean = inv(rotation(q(:))) .* meanRotation(grainId(:));
 
 end
 

EBSDAnalysis/@EBSD/plot.m

@@ -18,9 +18,9 @@
 %   mtexdata forsterite
 %   plot(ebsd)
 %
-%   plot(ebsd('Forsterite'))
+%   plot(ebsd('Forsterite'),ebsd('Forsterite').orientations)
 %
-%   plot(ebsd,ebsd.mad)
+%   plot(ebsd,ebsd.mad,'micronbar','off')
 %
 % See also
 % EBSD/plot

EBSDAnalysis/@grain2d/checkInside.m

@@ -1,9 +1,10 @@
-function isInside = checkInside( grains, xy )
+function isInside = checkInside(grains, xy, varargin)
 % check for points or grains to be inside a big grain
 %
 % Syntax
 %   isInside = checkInside(hostGrains, [x,y])
 %   isInside = checkInside(hostGrains, inclusionGrains)
+%   isInside = checkInside(hostGrains, ebsd)
 %
 % Input
 %  hostGrains      - @grain2d
@@ -13,6 +14,8 @@
 % Output
 %  isInside - numInclusionGrains x numHostGrains matrix
 %
+% Options
+%
 % Example
 %  mtexdata small
 %  grains = calcGrains(ebsd('indexed'))
@@ -21,6 +24,7 @@
 % See also
 % EBSD/findByLocation grain2d/findByOrientation
 
+ignoreInklusions = check_option(varargin,'ignoreInklusions');
 
 if isa(xy,'grain2d') 
 
@@ -31,6 +35,9 @@
 
   % a grain should not contain itself
   xy(ind,:) = NaN;
+
+  ignoreInklusions = true;
+
 elseif isa(xy,'EBSD')
 
   % extract unit cell
@@ -52,16 +59,14 @@
 
 poly = grains.poly;
 V = grains.V;
+incl = grains.inclusionId;
 
 for i = 1:length(poly)
 
-  p = poly{i};
-
-  firstloop = find(p(2:end)==p(1),1);
-
-  p = p(1:firstloop);
-  Vx = V(p,1); Vy = V(p,2);
+  p = poly{i};  
+  if ignoreInklusions, p = p(1:end-incl(i)); end
 
+  Vx = V(p,1); Vy = V(p,2);  
   isInside(:,i) = inpolygon(xy(:,1),xy(:,2),Vx,Vy);
 
 end

EBSDAnalysis/@grain2d/grain2d.m

@@ -7,7 +7,11 @@
     id=[]      % id of each grain    
     grainSize = [] % number of measurements per grain
   end
-
+
+  properties (Hidden = true)
+    inclusionId = [];
+  end
+
   % general properties
   properties    
     boundary = grainBoundary % boundary of the grains
@@ -19,9 +23,10 @@
     V                % vertices with x,y coordinates
     scanUnit         % unit of the vertice coordinates
     id2ind           % 
-    GOS              % intergranular average misorientation angle    
+    GOS              % intragranular average misorientation angle    
     x                % x coordinates of the vertices of the grains
     y                % y coordinates of the vertices of the grains
+    tripelPoints     % tripel points
   end
 
   properties (Dependent = true, Access = protected)
@@ -66,11 +71,11 @@
       grains.id = (1:numel(grains.phaseId)).';
       grains.grainSize = full(sum(I_DG,1)).';
 
-      grains.boundary = grainBoundary(V,F,I_FDext,ebsd);
+      grains.boundary = grainBoundary(V,F,I_FDext,ebsd,grains.phaseId);
       grains.boundary.scanUnit = ebsd.scanUnit;
-      grains.innerBoundary = grainBoundary(V,F,I_FDint,ebsd);
+      grains.innerBoundary = grainBoundary(V,F,I_FDint,ebsd,grains.phaseId);
 
-      grains.poly = calcPolygons(I_FDext * I_DG,F,V);
+      [grains.poly, grains.inclusionId]  = calcPolygons(I_FDext * I_DG,F,V);
 
 
       function [I_FDext,I_FDint] = calcBoundary
@@ -116,6 +121,10 @@
 
     function grains = set.V(grains,V)
       grains.boundary.V = V;
+
+      % update V in tripel points
+      tP = grains.tripelPoints;
+      grains.tripelPoints.V = V(tP.id,:);
     end
 
     function idV = get.idV(grains)
@@ -148,6 +157,14 @@
       unit = grains.boundary.scanUnit;
     end
 
+    function tP = get.tripelPoints(grains)
+      tP = grains.boundary.tripelPoints;
+    end
+
+    function grains = set.tripelPoints(grains,tP)
+      grains.boundary.tripelPoints = tP;
+    end
+
   end
 
 end

EBSDAnalysis/@grain2d/hasHole.m

@@ -2,10 +2,10 @@
 % test if a grain has a hole or not
 %
 % Input
-% grains - @Grain2d
+% grains - @grain2d
 %
 % Output
 % h  - logical array, |true| if a grain has hole
 %
 
-h = cellfun(@(x) nnz(x == x(1))>2,grains.poly);
+h = grains.inclusionId>0;

EBSDAnalysis/@grain2d/perimeter.m

@@ -1,9 +1,11 @@
 function  peri = perimeter(grains,varargin)
-% calculates the perimeter of a grain without holes
+% calculates the perimeter of a grain without or without inclusions
 %
 % Syntax
 %
-% grains.perimeter
+%   grains.perimeter
+%   perimter(grains)
+%   perimter(grains,'withInclusions')
 %
 % Input
 %  grains - @grain2d
@@ -14,8 +16,15 @@
 % See also
 % grain2d/equivalentPerimeter
 
-% ignore holes
-poly = cellfun(@(x) x(1:(1+find(x(2:end) == x(1),1))),grains.poly,'uniformOutput',false);
+poly = grains.poly;
 
-peri =  cellfun(@(ind) sum(sqrt(sum(diff(grains.V(ind,:)).^2,2))),poly);
-
+% remove inclusions
+if ~check_option(varargin,'withInclusion')
+  incl = grains.inclusionId;
+  for i = find(incl>0).'
+    poly{i} = poly{i}(1:end-incl(i));
+  end
+end
+
+V = grains.V;
+peri =  cellfun(@(ind) sum(sqrt(sum(diff(V(ind,:)).^2,2))),poly);

EBSDAnalysis/@grain2d/private/calcPolygons.m

@@ -1,10 +1,11 @@
-function poly = calcPolygons(I_FG,F,V)
+function [poly,inclId] = calcPolygons(I_FG,F,V)
 %
 % Input:
 %  I_FG - incidence matrix faces to grains
 %  F    - list of faces
 %  V    - list of vertices
 
+inclId = zeros(size(I_FG,2),1);
 poly = cell(size(I_FG,2),1);
 
 % for all grains
@@ -19,12 +20,36 @@
       EC{c} = fliplr(EC{c});
     end
   end
-
+    
   % this is needed
   for c=2:numel(EC), EC{c} = [EC{c} EC{1}(1)]; end
 
   poly{k} = [EC{:}];
+  inclId(k) = length(poly{k}) - length(EC{1});
 
 end
 
 end
+
+function test
+
+[x,y] = meshgrid(1:3);
+
+cs = crystalSymmetry('1')
+
+ori = orientation('axis',xvector,'angle',[0,90*degree],cs);
+
+ori_id = [ 0 1 1 ; ...
+  1 0 1;...
+  1 1 1];
+
+ori = ori(ori_id + 1);
+
+ebsd = EBSD(ori,ori_id,{cs},'options',struct('x',x,'y',y))
+
+end
+
+
+%
+%
+

EBSDAnalysis/@grain2d/subSet.m

@@ -28,6 +28,7 @@
 grains = subSet@dynProp(grains,ind);
 
 grains.poly = grains.poly(ind);
+grains.inclusionId = grains.inclusionId(ind);
 grains.id = grains.id(ind);
 grains.phaseId = grains.phaseId(ind);
 grains.grainSize = grains.grainSize(ind);

EBSDAnalysis/@grain2d/subsasgn.m

@@ -30,6 +30,7 @@
       grains.phaseId = subsasgn(grains.phaseId,s(1),b.phaseId);
       grains.grainSize = subsasgn(grains.grainSize,s(1),b.grainSize);
       grains.poly = subsasgn(grains.poly,s(1),b.poly);
+      grains.inclusionId = subsasgn(grains.inclusionId,s(1),b.inclusionId);
       grains.CSList = b.CSList;
       grains.phaseMap = b.phaseMap;
 

EBSDAnalysis/@grainBoundary/grainBoundary.m

@@ -20,8 +20,9 @@
 
   % general properties
   properties
-    V = []          % vertices x,y coordinates        
+    V = []          % vertices x,y coordinates            
     scanUnit = 'um' % unit of the vertice coordinates
+    tripelPoints    % tripel points
   end
 
   properties (Dependent = true)
@@ -34,11 +35,11 @@
     segmentId      % connected component id
     segmentSize    % number of faces that form a segment
     x              % x coordinates of the vertices of the grains
-    y              % y coordinates of the vertices of the grains
+    y              % y coordinates of the vertices of the grains    
   end
 
   methods
-    function gB = grainBoundary(V,F,I_FD,ebsd)
+    function gB = grainBoundary(V,F,I_FD,ebsd,grainsPhaseId)
 
       if nargin == 0, return; end
 
@@ -52,17 +53,17 @@
 
       % scale fid down to 1:length(gB)
       d = diff([0;fId]);      
-      fId = cumsum(d>0) + (d==0)*length(gB.F);
+      fId = cumsum(d>0) + (d==0)*size(gB.F,1);
 
-      gB.ebsdId = zeros(length(gB.F),2);
+      gB.ebsdId = zeros(size(gB.F,1),2);
       gB.ebsdId(fId) = eId;      
 
       % compute grainId
-      gB.grainId = zeros(length(gB.F),2);
+      gB.grainId = zeros(size(gB.F,1),2);
       gB.grainId(fId) = ebsd.grainId(eId);
 
       % compute phaseId
-      gB.phaseId = zeros(length(gB.F),2);
+      gB.phaseId = zeros(size(gB.F,1),2);
       isNotBoundary = gB.ebsdId>0;
       gB.phaseId(isNotBoundary) = ebsd.phaseId(gB.ebsdId(isNotBoundary));
       gB.phaseMap = ebsd.phaseMap;
@@ -76,12 +77,34 @@
       gB.grainId(doSort,:) = fliplr(gB.grainId(doSort,:));
 
       % compute misrotations
-      gB.misrotation = rotation(idquaternion(length(gB.F),1));
+      gB.misrotation = rotation(idquaternion(size(gB.F,1),1));
       isNotBoundary = all(gB.ebsdId,2);
       gB.misrotation(isNotBoundary) = ...
         inv(ebsd.rotations(gB.ebsdId(isNotBoundary,2))) ...
         .* ebsd.rotations(gB.ebsdId(isNotBoundary,1));
-
+
+      % compute tripel points
+      I_VF = gB.I_VF;
+      I_VG = (I_VF * gB.I_FG)==2;
+      % tripel points are those with exactly 3 neigbouring grains and 3
+      % boundary segments
+      itP = full(sum(I_VG,2)==3 & sum(I_VF,2)==3);
+      [tpGrainId,~] = find(I_VG(itP,:).');
+      tpGrainId = reshape(tpGrainId,3,[]).';      
+      tpPhaseId = full(grainsPhaseId(tpGrainId));
+
+      % compute ebsdId
+      % first step: compute faces at the tripel point
+      % clean up incidence matrix
+      %I_FD(~any(I_FD,2),:) = [];
+      % incidence matrix between tripel points and voronoi cells
+      %I_TD = I_VF(itP,:) * I_FD;
+      [tPBoundaryId,~] = find(I_VF(itP,:).');
+      tPBoundaryId = reshape(tPBoundaryId,3,[]).';
+
+      gB.tripelPoints = tripelPointList(find(itP),gB.V(itP,:),...
+        tpGrainId,tPBoundaryId,tpPhaseId,gB.phaseMap,gB.CSList);
+
     end
 
     function gB = cat(dim,varargin)
@@ -146,10 +169,14 @@
     end   
 
     function A_F = get.A_F(gB)
-      I_VF = gB.I_VF;           
-      A_F = I_VF.' * I_VF;
+      I_VF = gB.I_VF; %#ok<PROP>
+      A_F = I_VF.' * I_VF; %#ok<PROP>
     end
 
+    %function tp = get.tripelPoints(gB)
+
+    %end
+
     function segmentId = get.segmentId(gB)
       segmentId = connectedComponents(gB.A_F).';
     end
@@ -187,11 +214,17 @@
 
     function out = hasPhaseId(gB,phaseId,phaseId2)
 
+      if isempty(phaseId), out = false(size(gB)); return; end
+
       if nargin == 2
         out = any(gB.phaseId == phaseId,2);
 
         % not indexed phase should include outer border as well
         if phaseId > 0 && ischar(gB.CSList{phaseId}), out = out | any(gB.phaseId == 0,2); end
+
+      elseif isempty(phaseId2)
+
+        out = false(size(gB));
 
       elseif phaseId == phaseId2