Skip to content
Snippets Groups Projects
Commit f63ce796 authored by Peter Reischig's avatar Peter Reischig Committed by Nicola Vigano
Browse files

Deleted previous indexing functions in folder zUtil_Index.

git-svn-id: https://svn.code.sf.net/p/dct/code/trunk@507 4c865b51-4357-4376-afb4-474e03ccb993
parent c77025fe
No related branches found
No related tags found
No related merge requests found
Showing
with 0 additions and 2143 deletions
% function gtDrawGrainCubes(grains,sample,numbered,hlight)
%
% Draws a figure representing grain size and orientation by cubes
% (unit cells in a cubic lattice) based on the Rodrigues vectors.
%
% USAGE gtDrawGrainCubes(grains,sample)
% gtDrawGrainCubes(grains([1:10]),sample,1,[1,3,5])
%
% INPUT grains: grain structure from indexing
% sample: parameters and definition of the sample reference system
% (sample=gtSampleGeoInSampleSystem)
% numbered: if true, grain id-s are shown
% hlight: grain id-s to be highlighted
%
function gtDrawGrainCubes2(grains,sample,numbered,hlight)
if ~exist('hlight','var')
hlight=[];
end
if ~exist('numbered','var')
numbered=[];
end
nof_grains=length(grains);
% cube_vertices=[0 0 0; 1 0 0; 1 1 0; 0 1 0; 0 0 1; 1 0 1; 1 1 1; 0 1 1]-0.5*repmat([1 1 1],8,1);
% cube_faces=[1 2 6 5; 2 3 7 6; 3 4 8 7; 4 1 5 8; 1 2 3 4; 5 6 7 8];
strainscale=50; % 50
fs=14;
figure
colormap(jet(1000))
caxis([-0.01 0.01])
cb=colorbar('FontSize',fs,'ytick',-0.01:0.005:0.01);
%cb=colorbar('ytick',-0.01:0.002:0.01);
%cbpos=get(cb,'Position');
cbpos=[0.8 0.2 0.03 0.6];
set(cb,'Position',cbpos);
set(gca,'xtick',[],'xcolor','w')
set(gca,'ytick',[],'ycolor','w')
set(gca,'ztick',[],'zcolor','w')
view(3);
axhor=10*floor((sample.rad+30)/10);
axver=10*floor((max(sample.top,sample.bot)+50)/10);
axis([-axhor axhor -axhor axhor -axver axver])
axis equal;
%xlabel('X')
%ylabel('Y')
%zlabel('Z')
hold on
get(gca,'CameraPosition');
get(gca,'CameraViewAngle');
%camviewangle=
%campos=
t=1:360;
circx=cosd(t)*sample.rad;
circy=sind(t)*sample.rad;
circbotz(1:360)=sample.bot;
circtopz(1:360)=sample.top;
plot3(circx,circy,circbotz,'k-','Linewidth',1)
plot3(circx,circy,circtopz,'k-','Linewidth',1)
plot3([0 0],[0 0],[sample.bot sample.top],'k.-.','Linewidth',1)
plot3([-sample.rad -sample.rad],[0 0],[sample.bot sample.top],'k.-','Linewidth',1)
plot3([0 0],[-sample.rad -sample.rad],[sample.bot sample.top],'k.-','Linewidth',1)
plot3([sample.rad sample.rad],[0 0],[sample.bot sample.top],'k.-','Linewidth',1)
plot3([0 0],[sample.rad sample.rad],[sample.bot sample.top],'k.-','Linewidth',1)
% beam direction:
%plot3([-170 170],[0 0],[0 0],'r-','Linewidth',2)
cube_vertices=[0 0 0; 1 0 0; 1 1 0; 0 1 0; 0 0 1; 1 0 1; 1 1 1; 0 1 1]-0.5*repmat([1 1 1],8,1);
cube_faces=[2 3 7 6; 4 1 5 8; 1 2 6 5; 3 4 8 7; 5 6 7 8; 1 2 3 4];
% directions: +yz -yz -xz +xz +xy -xy
% +x -x -y +y +z -z
cube_planes=[1 0 0; -1 0 0; 0 -1 0; 0 1 0; 0 0 1; 0 0 -1];
for i=1:nof_grains
if ismember(i,hlight)
hl=1;
else
hl=false;
end
sfPlotGrainCube(grains{i},cube_vertices,cube_faces,cube_planes,hl,strainscale)
if numbered
text(grains{i}.center(1),-axhor,grains{i}.center(3),num2str(i),'Color','c')
text(axhor,grains{i}.center(2),grains{i}.center(3),num2str(i),'Color','c')
end
end
end
function sfPlotGrainCube(grain,cube_vertices,cube_faces,cube_planes,hl,strainscale)
ga=grain.stat.bbysmean/2;
gc=grain.center;
% if hl
% facecolors=cool(8);
% facecolors=facecolors(end-5:end,:);
% else
% facecolors=copper(8);
% facecolors=facecolors(end-5:end,:);
% end
% ELONGATION ALONG LOADING AXIS
%
ggrain=Rod2g(grain.R_vector);
grain_vertices=ga*cube_vertices*ggrain';
grain_vertices_st=repmat(gc,8,1)+((eye(3)+grain.strain.ST*strainscale)*grain_vertices')';
facecolor=grain.strain.Eps(3);
facecolor=repmat(facecolor,6,1);
%facecolor=facecolor.*repmat([1.0; 0.98; 0.96; 0.94; 0.92; 0.90],1,3);
% SHEAR STRAIN ON THE CUBE CELL FACES
% T=Rod2g(grain.R_vector); % T rotation tensor in lab system (align the axes with the crystal axes);
%
% ECRYST=T'*grain.strain.ST*T; % strain tensor in the crystal reference
%
% facecolor=zeros(6,1);
% facecolor([1,2])=sqrt(ECRYST(1,2)^2+ECRYST(1,3)^2); % xy & xz
% facecolor([3,4])=sqrt(ECRYST(2,1)^2+ECRYST(2,3)^2); % yx & yz
% facecolor([5,6])=sqrt(ECRYST(3,1)^2+ECRYST(3,2)^2); % zx & zy
%
% ST=grain.strain.ST;
%
% grain_vertices=ga*cube_vertices*T';
% grain_vertices_st=repmat(gc,8,1)+((eye(3)+ST*strainscale)*grain_vertices')';
%
patch('Vertices',grain_vertices_st,'Faces',cube_faces,...
'FaceVertexCData',facecolor,'FaceColor','flat');
% x0=[1 0 0];
% y0=[0 1 0];
% z0=[0 0 1];
%
% xcr=T*x0';
% ycr=T*y0';
% zcr=T*z0';
%
%
% quiver3(gc(1),gc(2),gc(3),xcr(1),xcr(2),xcr(3),200,'r')
% quiver3(gc(1),gc(2),gc(3),ycr(1),ycr(2),ycr(3),200,'b')
% quiver3(gc(1),gc(2),gc(3),zcr(1),zcr(2),zcr(3),200,'k')
%
% grain_vertices=ga*cube_vertices*T';
% grain_vertices_st=repmat(gc,8,1)+((eye(3))*grain_vertices')';
%
% patch('Vertices',grain_vertices_st,'Faces',cube_faces,...
% 'FaceVertexCData',facecolor,'FaceColor','flat');
%
%
% edgevec=[0 1 0];
% facepos=[1 0 0];
% edgevec_g=edgevec*T';
% facepos_g=gc+ga/2*facepos*T';
% quiver3(facepos_g(1),facepos_g(2),facepos(3),edgevec_g(1),edgevec_g(2),edgevec_g(3),100,'c');
end
%
% INPUT: Two data series of segmented and paired diffraction spots 180 degrees
% offset from Spotpair or Calibration tables.
%
%
%
% OUTPUT: Diffraction spots grouped into grains. Location, orientation
% and volume of grains are determined.
%
%
%
function [grain,allgrainstat,probl,singlesfit,errest]=gtINDEXTER(strategy,flag_update,grain)
load('parameters.mat');
if ~exist('strategy','var')
strategy=[];
end
if ~exist('grain','var')
grain=[];
end
% if ~exist('flag_cont','var')
% flag_cont=false;
% end
if ~exist('flag_update','var')
flag_update=false;
end
tic
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% SET PARAMETERS
% Sample coordinate system defined in pixels by:
sample=gtSampleGeoInSampleSystem(parameters);
spacegroup=parameters.cryst.spacegroup;
latticepar=parameters.cryst.latticepar;
% allowed angles and their multiplicity according of plane families
[ACM,ACMU]=gtAngConsMat(parameters.cryst); % (in degrees)
reflections=parameters.cryst.hklsp;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Load data from database
if isfield(parameters.acq, 'pair_tablename')
pairtable=parameters.acq.pair_tablename;
else % what should be standard convention
pairtable=sprintf('%sspotpairs',parameters.acq.name);
end
gtDBConnect();
% All the pair data are loaded in 'tot':
% if flag_cont % continue indexing
% mymcmd=sprintf(['SELECT pairID, difAID, difBID, theta, eta, omega, '...
% 'avintint, avbbXsize, avbbYsize, samcentXA, samcentYA, samcentZA, '...
% 'samcentXB, samcentYB, samcentZB, ldirX, ldirY, ldirZ, plX, plY, plZ, '...
% 'thetatype FROM %s WHERE isnull(grainid) ORDER BY avbbYsize desc'],pairtable);
% disp('Continuation does not work yet...')
% return
% else % (re)start indexing
mymcmd=sprintf(['SELECT pairID, difAID, difBID, theta, eta, omega, '...
'avintint, avbbXsize, avbbYsize, samcentXA, samcentYA, samcentZA, '...
'samcentXB, samcentYB, samcentZB, ldirX, ldirY, ldirZ, plX, plY, plZ, '...
'thetatype FROM %s ORDER BY avbbYsize desc'],pairtable);
%end
[tot.pairid,tot.aid,tot.bid,tot.theta,tot.eta,tot.omega,...
tot.int,tot.bbxs,tot.bbys, ...
tot.ca(:,1),tot.ca(:,2),tot.ca(:,3),tot.cb(:,1),tot.cb(:,2),tot.cb(:,3),...
tot.dir(:,1),tot.dir(:,2),tot.dir(:,3), ...
tot.pl(:,1),tot.pl(:,2),tot.pl(:,3),tot.thetatype]=mym(mymcmd);
nof_lines=length(tot.pairid) ;
% these line id-s will not correspond to pairID-s !!!
tot.id=(1:nof_lines)';
% Make sure directions are normalized
% (e.g. acosd function after rounding errors from database)
pllength=sqrt(tot.pl(:,1).^2 + tot.pl(:,2).^2 + tot.pl(:,3).^2) ;
tot.pl=[tot.pl(:,1)./pllength tot.pl(:,2)./pllength tot.pl(:,3)./pllength] ;
dirlength=sqrt(tot.dir(:,1).^2 + tot.dir(:,2).^2 + tot.dir(:,3).^2) ;
tot.dir=[tot.dir(:,1)./dirlength tot.dir(:,2)./dirlength tot.dir(:,3)./dirlength] ;
for i=1:nof_lines
tot.hkl(i,:)=reflections(tot.thetatype(i),:);
end
% Exclude Friedel pairs which are discarded
if isfield(parameters.index,'discard')
consider=true(nof_lines,1);
discard=[];
for i=1:length(parameters.index.discard)
tmp=find(parameters.index.discard(i)==tot.pairid);
if ~isempty(tmp)
discard=[discard, tmp]; % line ID belonging to pair ID
end
end
consider(discard)=false;
tot=gtIndexSelectLines(consider,tot);
nof_lines=length(tot.pairid) ;
tot.id=(1:nof_lines)';
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Execute strategy and create output
if isempty(strategy)
if isfield(parameters,'index')
if isfield(parameters.index,'strategy')
strategy=parameters.index.strategy;
else
strategy=gtIndexDefaultStrategy;
end
else
strategy=gtIndexDefaultStrategy;
end
end
[grain,singleids]=gtIndexExecuteStrategy(strategy,grain,tot,sample,ACM,spacegroup,latticepar);
if isempty(grain)
disp(' ')
disp(' No grain has been found.')
return
end
singles=gtIndexSelectLines(singleids,tot);
nof_grains=length(grain);
% Create final grain output
disp(' ')
disp('CREATING OUTPUT...')
tic
if flag_update % && ~flag_cont
mym(sprintf('UPDATE %s SET grainID=null', pairtable)) ;
end
grain=gtIndexCreateOutput(grain,(1:nof_grains),tot,strategy.r,spacegroup,ACM,latticepar,pairtable,flag_update);
% Intermediate saving before singles
fname=gtLastFileName('indexter_lastrun','new');
save(fname,'grain','singles','strategy');
disp(' ')
disp(['Results have been saved in ', fname])
toc
% Deal with singles
singlesfit=[];
if ~isempty(strategy.s)
if ~isempty(singles.id) %if no singles, then skip
disp(' ')
disp('TREATING SINGLES...')
tic
[grain,changedg,singlesfit]=gtIndexFitSinglesToGrains(singles,grain,strategy.s,spacegroup,latticepar);
toc
disp(' ')
disp('RECREATING OUTPUT...')
disp(' ')
tic
grain=gtIndexCreateOutput(grain,changedg,tot,strategy.r,spacegroup,ACM,latticepar,pairtable,flag_update);
toc
end
else
singlesfit.pairid=[];
end
% Collect lines which have been assigned
% grouped=false(nof_lines,1);
% for i=1:nof_grains
% grouped(grain{i}.lid)=true;
% end
% sum(grouped)
grouped=[];
for i=1:nof_grains
grouped=[grouped; [grain{i}.lid', grain{i}.pairid', repmat(i,length(grain{i}.lid),1)]];
end
%size(grouped,1)
%length(unique(grouped(:,1)))
grouped=sortrows(grouped,1);
group=gtIndexSelectLines(grouped(:,1),tot);
sfpairid=[];
for i=1:length(singlesfit.pairid)
pid=singlesfit.pairid(i);
sfpairid(i,:)=[pid, grouped((grouped(:,2)==pid),3)];
end
if ~isempty(sfpairid)
singlesfit.pairid=sortrows(sfpairid,1);
%singlesfit.pairid=[singlesfit.pairid', grouped(singlesfit.pairid',2)];
else
singlesfit.pairid=[];
end
% Determine unique planes
disp(' ')
grain=gtIndexUniquePlaneNormals(grain,ACM);
unipl=gtAllGrainValues(grain,'uni','nof_unipl',1);
% Grain statistics
[allgrainstat,probl]=gtIndexAllGrainStat(grain);
% Estimate errors:
errest=gtEstimateErrors(grain);
% Save results
disp(' ')
save(fname,'grain','group','singles','singlesfit','allgrainstat','probl','errest','strategy');
disp(['Results have been saved in ', fname])
disp(' ')
disp(' ')
disp('Number of Friedel pairs:'), disp(length(tot.id)) ;
disp('Number of pairs indexed:'), disp(length(group.id)) ;
disp('Number of pairs left unassigned:'), disp(length(tot.id)-length(group.id)) ;
disp('Number of single pairs:'), disp(length(singles.id)) ;
disp('Number of grains found:'), disp(length(grain)) ;
disp('Average number of pairs per grain:'), disp(length(group.id)/length(grain)) ;
disp('Average number of unique pairs per grain:'), disp(mean(unipl)) ;
disp(' ')
figure('name','Eta vs theta from grouped pairs')
plot(group.theta,group.eta,'b.')
end % of function
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% SUB-FUNCTIONS
function lname=sfLatestFileName(w)
fnames=dir('indexter_lastrun*.mat')
nmb=[];
for i=1:length(fnames)
tmp=str2num(fnames(i).name(17:end-4))
if isempty(tmp)
nmb(i)=0;
else
nmb(i)=tmp;
end
end
[maxval,maxloc]=max(nmb);
if isnan(maxval)
maxval=1;
end
switch w
case 'existing'
lname=fnames(maxloc).name;
case 'new'
lname=sprintf('indexter_lastrun%d.mat',maxval+1);
end
end
function [newgrain,remaining]=gtIndexAddToGrains(remaining,grain,tot,...
sample,ACM,tol_b)
disp(' ')
disp(' ADDING PAIRS TO GRAINS...')
disp(' ')
disp(tol_b)
tic
nof_grains=length(grain);
newgrain=grain;
for i=1:nof_grains
nof_lines=length(grain{i}.lid);
pgood=true(nof_lines+1,1);
pgood(nof_lines+1)=false;
ptry=nof_lines+1;
goodfit=false(length(remaining),1);
% Distance tolerance
if grain{i}.merged
tol_b.dist=max(max(grain{i}.stat.distlines),tol_b.distmax);
else
tol_b.dist=tol_b.distmax;
end
glines=gtIndexSelectLines(grain{i}.lid,tot);
for j=1:length(remaining)
%cand=gtIndexSelectLines([grain{i}.lid,remaining(j)],tot);
%goodfit(j)=gtIndexTryAddLineToGroup(ptry,pgood,cand,tol_b,sample,ACM);
tryline=gtIndexSelectLines(remaining(j),tot);
goodfit(j)=gtIndexCheckGroupCons(tryline,glines,grain{i}.center,tol_b,sample,ACM);
end
if any(goodfit)
lineids=[grain{i}.lid'; remaining(goodfit)];
newlines=gtIndexSelectLines(lineids,tot); % reload original lines in the new grain
newgrain{i}=gtIndexCreateGrainOutputAdd(newlines);
newgrain{i}.id=i;
newgrain{i}.merged=grain{i}.merged;
disp(['Pairs added to grain #' num2str(i) ':'])
disp(remaining(goodfit))
end
remaining(goodfit)=[];
end
disp(' ')
toc
end
% Grain statistics for output and check
function [allgrainstat,probl]=gtIndexAllGrainStat(grain)
nof_grains=length(grain);
allgrainstat.intrel=[];
allgrainstat.intrat=[];
allgrainstat.bbxsrel=[];
allgrainstat.bbysrel=[];
allgrainstat.bbxsrat=[];
allgrainstat.bbysrat=[];
allgrainstat.distcom=[];
allgrainstat.distlines=[];
allgrainstat.angplanes=[];
allgrainstat.Rdist=[];
%grouped=false(nof_lines,1);
% Detected problems:
probl.distcom=[];
probl.Rdist=[];
probl.int=[];
probl.bbxs=[];
probl.bbys=[];
for i=1:nof_grains
allgrainstat.distcom=[allgrainstat.distcom, grain{i}.stat.distcom];
allgrainstat.distlines=[allgrainstat.distlines, grain{i}.stat.distlines];
allgrainstat.angplanes=[allgrainstat.angplanes, grain{i}.stat.angplanes];
allgrainstat.Rdist=[allgrainstat.Rdist, grain{i}.stat.Rdist];
allgrainstat.intrel=[allgrainstat.intrel, grain{i}.stat.intrel];
allgrainstat.intrat=[allgrainstat.intrat, grain{i}.stat.intrat];
allgrainstat.bbxsrel=[allgrainstat.bbxsrel, grain{i}.stat.bbxsrel];
allgrainstat.bbysrel=[allgrainstat.bbysrel, grain{i}.stat.bbysrel];
allgrainstat.bbxsrat=[allgrainstat.bbxsrat, grain{i}.stat.bbxsrat];
allgrainstat.bbysrat=[allgrainstat.bbysrat, grain{i}.stat.bbysrat];
% Mark outliers
if ~isempty(grain{i}.stat.outl.int)
tmp=[repmat(i,length(grain{i}.stat.outl.int),1) grain{i}.stat.outl.int'];
probl.int=[probl.int; tmp]; % mark as problem
end
if ~isempty(grain{i}.stat.outl.bbxs)
tmp=[repmat(i,length(grain{i}.stat.outl.bbxs),1) grain{i}.stat.outl.bbxs'];
probl.bbxs=[probl.bbxs; tmp]; % mark as problem
end
if ~isempty(grain{i}.stat.outl.bbys)
tmp=[repmat(i,length(grain{i}.stat.outl.bbys),1) grain{i}.stat.outl.bbys'];
probl.bbys=[probl.bbys; tmp]; % mark as problem
end
if ~isempty(grain{i}.stat.outl.distcom)
tmp=[repmat(i,length(grain{i}.stat.outl.distcom),1) grain{i}.stat.outl.distcom'];
probl.distcom=[probl.distcom; tmp]; % mark as problem
end
if ~isempty(grain{i}.stat.outl.Rdist)
tmp=[repmat(i,length(grain{i}.stat.outl.Rdist),1) grain{i}.stat.outl.Rdist'];
probl.Rdist=[probl.Rdist; tmp]; % mark as problem
end
end
nof_bins=nof_grains;
figure('name','Relative intensities in the grains')
hist(allgrainstat.intrel,nof_bins)
figure('name','Relative bbox X size in the grains')
hist(allgrainstat.bbxsrel,nof_bins)
figure('name','Relative bbox Y size in the grains')
hist(allgrainstat.bbysrel,nof_bins)
figure('name','Max. intensity ratio in the grains')
hist(allgrainstat.intrat,max(1,round(nof_bins/10)))
figure('name','Max. bbox X size ratio in the grains')
hist(allgrainstat.bbxsrat,max(1,round(nof_bins/10)))
figure('name','Max. bbox Y size ratio in the grains')
hist(allgrainstat.bbysrat,max(1,round(nof_bins/10)))
figure('name','Angular deviation of planes')
hist(allgrainstat.angplanes,nof_bins)
figure('name','Difference between planes and grain orientation in Rodrigues space')
hist(allgrainstat.Rdist,nof_bins)
figure('name','Distance between diffraction paths')
hist(allgrainstat.distlines,nof_bins)
figure('name','Distance of diffraction paths to grain center of mass')
hist(allgrainstat.distcom,nof_bins)
end
function [grain,singles]=gtIndexBuildGrains(remaining_in,tot,sample,ACM,tol_b)
disp(' ')
disp(' BUILDING GRAINS...')
disp(' ')
disp(tol_b)
tic
remaining=remaining_in;
nof_grains=0;
singles=[];
grain=[];
while length(remaining)>=1
act=gtIndexSelectLines(remaining(1),tot);
cand_rem=gtIndexSelectLines(remaining(2:end),tot);
tol_b.dist=min(act.bbys/tol_b.distsf,tol_b.distmax);
[paircons,cand0,isecs]=gtIndexCheckPairCons(act,cand_rem,tol_b,sample,ACM);
if length(paircons)>3
cand_paircons=gtIndexSelectLines([act.id; cand0.id],tot);
%[pgroup1,grainoutput1]=gtIndexFindGroupInCand(cand_paircons,tol_b,sample,ACM);
[pgroup,grainoutput]=gtIndexFindGroupInCand(cand_paircons,isecs,tol_b,sample,ACM);
if sum(pgroup)>1
nof_grains=nof_grains+1;
grain{nof_grains}=grainoutput;
grain{nof_grains}.id=nof_grains;
disp(['Found new grain #' num2str(nof_grains)])
disp(grain{nof_grains}.pairid)
todel=cand_paircons.id(pgroup);
else
singles=[singles; act.id];
todel=act.id;
end
else
singles=[singles; act.id];
todel=act.id;
end
for i=1:length(todel)
remaining(remaining==todel(i))=[];
end
end
disp(' ')
toc
end
%
% Gives OK, if all the lines are pair consistent with a given line.
%
% act: the actual line (structure)
% cand: candidates (structure) to be tested for pair consistency with 'act'
% tol: tolerances (structure)
% sample: sample geometry parameters (structure)
% ACM: angular consistency matrix
function ok=gtIndexCheckAllPairCons(act,cand,tol,sample,ACM)
ok=false;
%paircons=(1:length(cand.id))'; % initial indices of consistent lines
%progress=length(cand.id);
% bbox y size close enough?
cons=(act.bbys/tol.bbys<cand.bbys) & (act.bbys*tol.bbys>cand.bbys);
if ~all(cons)
return
end
%cand=gtIndexSelectLines(cons,cand);
%paircons(~cons)=[];
%progress=[progress; length(cand.id)];
% bbox x size close enough?
cons=(act.bbxs/tol.bbxs<cand.bbxs) & (act.bbxs*tol.bbxs>cand.bbxs);
if ~all(cons)
return
end
%cand=gtIndexSelectLines(cons,cand);
%paircons(~cons)=[];
%progress=[progress; length(cand.id)];
% integrated intensities close enough?
cons=(act.int/tol.int<cand.int) & (act.int*tol.int>cand.int);
if ~all(cons)
return
end
%cand=gtIndexSelectLines(cons,cand);
%paircons(~cons)=[];
%progress=[progress; length(cand.id)];
% spatial distance small enough?
cons=gt2LinesDist([act.ca act.dir],[cand.ca cand.dir])<tol.dist;
if ~all(cons)
return
end
%cand=gtIndexSelectLines(cons,cand);
%paircons(~cons)=[];
%progress=[progress; length(cand.id)];
% is intersection in the sample volume?%
% NOTE: intersection of 2 lines might actually be far from the grain center
% if they are near being parallel. This could be taken into account here...
% On the other hand, there should be no parallel lines in a grain set.
pofinters=gt2LinesIntersection([act.ca act.dir],[cand.ca cand.dir]);
cons=gtIsPointInSample(pofinters,sample,tol);
if ~all(cons)
return
end
%cand=gtIndexSelectLines(cons,cand);
%pofinters(~cons,:)=[];
%paircons(~cons)=[];
%progress=[progress; length(cand.id)];
% plane angles consistent?
cons=false(length(cand.id),1);
angle=zeros(length(cand.id),1);
for i=1:length(cand.id)
ACV=ACM{act.thetatype,cand.thetatype(i)};
ang=gt2PlanesAngle(act.pl,cand.pl(i,:));
angle(i)=ang;
cons(i)=min(abs(ACV-ang))<tol.ang;
end
if ~all(cons)
return
end
%cand=gtIndexSelectLines(cons,cand);
%pofinters(~cons,:)=[];
%angle(~cons)=[];
%paircons(~cons)=[];
%progress=[progress; length(cand.id)];
ok=true;
end
function groupcons=gtIndexCheckGroupCons(tryline,goodl,goodpofinters,tol,sample,ACM)
groupcons=false;
paircons=gtIndexCheckAllPairCons(tryline,goodl,tol,sample,ACM);
if paircons
if pointtolinedist(goodpofinters,[tryline.ca tryline.dir])<tol.dist;
groupcons=true;
end
end
end
%
% Checks pair consistency in a set of lines.
%
% act: the actual line (structure)
% cand: candidates (structure) to be pair consistent with it
% tol: tolerances (structure)
% sample: sample geometry parameters (structure)
% ACM: angular consistency matrix
function [paircons,cand,pofinters]=gtIndexCheckPairCons(act,cand,tol,sample,ACM)
paircons=(1:length(cand.id))'; % initial indices of consistent lines
%progress=length(cand.id);
% integrated intensities close enough?
cons=(act.int/tol.int<cand.int) & (act.int*tol.int>cand.int);
cand=gtIndexSelectLines(cons,cand);
paircons(~cons)=[];
%progress=[progress; length(cand.id)];
% bbox x size close enough?
cons=(act.bbxs/tol.bbxs<cand.bbxs) & (act.bbxs*tol.bbxs>cand.bbxs);
cand=gtIndexSelectLines(cons,cand);
paircons(~cons)=[];
%progress=[progress; length(cand.id)];
% bbox y size close enough?
cons=(act.bbys/tol.bbys<cand.bbys) & (act.bbys*tol.bbys>cand.bbys);
cand=gtIndexSelectLines(cons,cand);
paircons(~cons)=[];
%progress=[progress; length(cand.id)];
% spatial distance small enough?
cons=gt2LinesDist([act.ca act.dir],[cand.ca cand.dir])<tol.dist;
cand=gtIndexSelectLines(cons,cand);
paircons(~cons)=[];
%progress=[progress; length(cand.id)];
% is intersection in the sample volume?%
% NOTE: intersection of 2 lines might actually be far from the grain center
% if they are near being parallel. This could be taken into account here...
% On the other hand, there should be no parallel lines in a grain set.
pofinters=gt2LinesIntersection([act.ca act.dir],[cand.ca cand.dir]);
cons=gtIsPointInSample(pofinters,sample,tol);
cand=gtIndexSelectLines(cons,cand);
pofinters(~cons,:)=[];
paircons(~cons)=[];
%progress=[progress; length(cand.id)];
% plane angles consistent?
cons=false(length(cand.id),1);
angle=zeros(length(cand.id),1);
for i=1:length(cand.id)
ACV=ACM{act.thetatype,cand.thetatype(i)};
ang=gt2PlanesAngle(act.pl,cand.pl(i,:));
angle(i)=ang;
cons(i)=min(abs(ACV-ang))<tol.ang;
end
cand=gtIndexSelectLines(cons,cand);
pofinters(~cons,:)=[];
%angle(~cons)=[];
paircons(~cons)=[];
%progress=[progress; length(cand.id)];
end
% Creates 'grain' output structure from indexed (grouped) pair-lines.
%
function grain=gtIndexCreateGrainOutputAdd(grainl) % ,spacegroup)
grain.id=[];
% Pair ID
grain.pairid=grainl.pairid';
% Line ID-s in the indexing (normally the pairs ordered by intensities)
% i.e. NOT the pairID-s!
grain.lid=grainl.id';
grain.center=pofintersectexp([grainl.ca grainl.dir])'; % center of mass
grain.stat.intmean=mean(grainl.int); % mean intensity value
grain.stat.bbxsmean=mean(grainl.bbxs); % mean
grain.stat.bbysmean=mean(grainl.bbys); % mean
grain.stat.distlines=[];
% grain.stat.angplanes=[];
for gi=1:length(grainl.id)
for gj=gi+1:length(grainl.id)
grain.stat.distlines=[grain.stat.distlines,...
gt2LinesDist([grainl.ca(gi,:) grainl.dir(gi,:)],...
[grainl.ca(gj,:) grainl.dir(gj,:)])];
%
% ACV=ACM{grainl.thetatype(gi),grainl.thetatype(gj)};
% ang=gt2PlanesAngle(grainl.pl(gi,:),grainl.pl(gj,:));
% ang=min(abs(ACV-ang));
% grain.stat.angplanes=[grain.stat.angplanes, ang];
end
end
%grain.R_vector=[0, 0, 0];
%grain.merged=false;
end
% Creates 'grain' output structure from indexed (grouped) pair-lines.
%
function grain=gtIndexCreateGrainOutputBuild(grainl) % ,spacegroup)
grain.id=[];
% Pair ID
grain.pairid=grainl.pairid';
% Line ID-s in the indexing (normally the pairs ordered by intensities)
% i.e. NOT the pairID-s!
grain.lid=grainl.id';
grain.center=pofintersectexp([grainl.ca grainl.dir])'; % center of mass
grain.stat.intmean=mean(grainl.int); % mean intensity value
grain.stat.bbxsmean=mean(grainl.bbxs); % mean
grain.stat.bbysmean=mean(grainl.bbys); % mean
grain.R_vector=[0, 0, 0];
grain.merged=false;
end
% Creates 'grain' output structure from indexed (grouped) pair-lines.
%
function grain=gtIndexCreateGrainOutputFinal(grainl,ACM,spacegroup,latticepar,tol_outl)
nof_lines=length(grainl.id);
% Tolerance factor for marking outliers (*std)
%tol_outl=3;
% Order lines by thetatype and plane normal
sortM(:,1)=grainl.thetatype';
sortM(:,2)=grainl.pl(:,3);
sortM(:,3)=(1:nof_lines)';
sortM=sortrows(sortM,[1 -2]);
grainl=gtIndexSelectLines(sortM(:,3),grainl);
% ID and Rodrigues vector
grain.id=[];
% Center of mass in sample system; for coordinate definition, look at
% gtSampleGeoInSampleSystem
grain.center=pofintersectexp([grainl.ca grainl.dir])'; % center of mass
grain.R_vector=gtRodriguesTest(grainl.pl,grainl.hkl,spacegroup,latticepar);
% Pair ID
grain.nof_pairs=length(grainl.pairid);
grain.pairid=grainl.pairid';
% Thetatype, hkl vectors and plane normals
grain.thetatype=grainl.thetatype';
grain.hkl=grainl.hkl;
grain.pl=grainl.pl;
grain.pla=NaN(1,nof_lines);
% Diffraction angles
grain.theta=grainl.theta';
grain.eta=grainl.eta';
grain.omega=grainl.omega';
% Line ID-s in the indexing (normally the pairs ordered by intensities)
% i.e. these are NOT the pairID-s!
grain.lid=grainl.id';
% Grain statistics for tolerances
% Difference of each plane and the average grain orientation in Rodrigues space:
for i=1:length(grainl.id)
Rv=gtRodriguesVectors2(grainl.pl(i,:),grainl.hkl(i,:),spacegroup,latticepar);
grain.stat.Rdist(i)=min(pointtolinedist(grain.R_vector,Rv));
end
% Mean and standard deviation of distance in Rodrigues space
grain.stat.Rdistmean=mean(grain.stat.Rdist);
grain.stat.Rdiststd=std(grain.stat.Rdist);
% Line distances to the center of mass
grain.stat.distcom=pointtolinedist(grain.center,[grainl.ca grainl.dir])';
grain.stat.distcommean=mean(grain.stat.distcom);
grain.stat.distcomstd=std(grain.stat.distcom);
% Line distances and plane angles pair-wise
grain.stat.distlines=[];
grain.stat.angplanes=[];
for gi=1:length(grainl.id)
for gj=gi+1:length(grainl.id)
grain.stat.distlines=[grain.stat.distlines,...
gt2LinesDist([grainl.ca(gi,:) grainl.dir(gi,:)],...
[grainl.ca(gj,:) grainl.dir(gj,:)])];
ACV=ACM{grainl.thetatype(gi),grainl.thetatype(gj)};
ang=gt2PlanesAngle(grainl.pl(gi,:),grainl.pl(gj,:));
ang=min(abs(ACV-ang));
grain.stat.angplanes=[grain.stat.angplanes, ang];
end
end
% Intensity
grain.stat.int=grainl.int';
grain.stat.intmean=mean(grainl.int); % mean intensity value
grain.stat.intstd=std(grainl.int); % std intensity
grain.stat.intrel=(grainl.int/grain.stat.intmean)'; % relative intensity
grain.stat.intrat=max(grainl.int)/min(grainl.int); % extreme intensity ratio
% Bounding box x size
grain.stat.bbxs=grainl.bbxs';
grain.stat.bbxsmean=mean(grainl.bbxs); % mean
grain.stat.bbxsstd=std(grainl.bbxs); % std
grain.stat.bbxsrel=(grainl.bbxs/grain.stat.bbxsmean)'; % relative bbox size
grain.stat.bbxsrat=max(grainl.bbxs)/min(grainl.bbxs); % extreme ratio
% Bounding box y size
grain.stat.bbys=grainl.bbys';
grain.stat.bbysmean=mean(grainl.bbys); % mean
grain.stat.bbysstd=std(grainl.bbys); % std
grain.stat.bbysrel=(grainl.bbys/grain.stat.bbysmean)'; % relative bbox size
grain.stat.bbysrat=max(grainl.bbys)/min(grainl.bbys); % extreme ratio
% Mark outliers
grain.stat.outl.int=gtOutliers(grainl.int', tol_outl, grain.stat.intmean, grain.stat.intstd);
grain.stat.outl.bbxs=gtOutliers(grainl.bbxs', tol_outl, grain.stat.bbxsmean, grain.stat.bbxsstd);
grain.stat.outl.bbys=gtOutliers(grainl.bbys', tol_outl, grain.stat.bbysmean, grain.stat.bbysstd);
grain.stat.outl.distcom=gtOutliers(grain.stat.distcom, tol_outl, grain.stat.distcommean, grain.stat.distcomstd);
grain.stat.outl.Rdist=gtOutliers(grain.stat.Rdist, tol_outl, grain.stat.Rdistmean, grain.stat.Rdiststd);
grain.stat.outl.int=grain.pairid(grain.stat.outl.int);
grain.stat.outl.bbxs=grain.pairid(grain.stat.outl.bbxs);
grain.stat.outl.bbys=grain.pairid(grain.stat.outl.bbys);
grain.stat.outl.distcom=grain.pairid(grain.stat.outl.distcom);
grain.stat.outl.Rdist=grain.pairid(grain.stat.outl.Rdist);
end % of function
% Creates 'grain' output structure from indexed (grouped) pair-lines.
%
function grain=gtIndexCreateGrainOutputMerge(grainl)
% nof_lines=length(grainl.id);
%
% % Tolerance factor for marking outliers (*std)
% tol_outl=3;
% % modif sabine
%
% %ACM
%
% % fin modif sabine
% % Order lines by thetatype and plane normal
% sortM(:,1)=grainl.thetatype';
% sortM(:,2)=grainl.pl(:,3);
% sortM(:,3)=(1:nof_lines)';
%
% sortM=sortrows(sortM,[1 -2]);
%
% grainl=gtIndexSelectLines(sortM(:,3),grainl);
% ID and Rodrigues vector
grain.id=[];
grain.merged=true;
% Center of mass in sample system; for coordinate definition, look at
% gtSampleGeoInSampleSystem
grain.center=pofintersectexp([grainl.ca grainl.dir])'; % center of mass
%grain.R_vector=NaN(1,3);
% Line ID-s in the indexing (normally the pairs ordered by intensities)
% i.e. NOT the pairID-s!
grain.lid=grainl.id';
% % Pair ID
% grain.pairid=grainl.pairid';
%
% % Thetatype, hkl vectors and plane normals
% grain.thetatype=grainl.thetatype';
%grain.hkl=grainl.hkl;
%grain.pl=grainl.pl;
% grain.pla=NaN(1,nof_lines);
% Unique plane normals (multiples are averaged)
%grain=gtINUniquePlaneNormals(grainl,ACM);
% Rodrigues vectors (multiple) of the plane normals of a line
%for i=1:length(grainl.id)
% Rv=gtRodriguesVectors(grainl.pl(i,:),grainl.hkl(i,:),spacegroup);
% grain.lRv(1:size(Rv,1),1:size(Rv,2),i)=Rv;
%end
%grain.lRv=grainl.lRv;
% grain.theta=grainl.theta';
% grain.eta=grainl.eta';
% grain.omega=grainl.omega';
%
% % Grain statistics for tolerances
% grain.stat.int=grainl.int;
grain.stat.intmean=mean(grainl.int); % mean intensity value
% grain.stat.intstd=std(grainl.int); % std intensity
% grain.stat.intrel=(grainl.int/grain.stat.intmean)'; % relative intensity
% grain.stat.intrat=max(grainl.int)/min(grainl.int); % extreme intensity ratio
% grain.stat.bbxs=grainl.bbxs;
grain.stat.bbxsmean=mean(grainl.bbxs); % mean
% grain.stat.bbxsstd=std(grainl.bbxs); % std
% grain.stat.bbxsrel=(grainl.bbxs/grain.stat.bbxsmean)'; % relative bbox size
% grain.stat.bbxsrat=max(grainl.bbxs)/min(grainl.bbxs); % extreme ratio
% grain.stat.bbys=grainl.bbys;
grain.stat.bbysmean=mean(grainl.bbys); % mean
% grain.stat.bbysstd=std(grainl.bbys); % std
% grain.stat.bbysrel=(grainl.bbys/grain.stat.bbysmean)'; % relative bbox size
% grain.stat.bbysrat=max(grainl.bbys)/min(grainl.bbys); % extreme ratio
% Line distances to the center of mass
% grain.stat.distcom=pointtolinedist(grain.center,[grainl.ca grainl.dir])';
% grain.stat.distcommean=mean(grain.stat.distcom);
% grain.stat.distcomstd=std(grain.stat.distcom);
%
% % Line distances and plane angles pair-wise
grain.stat.distlines=[];
% grain.stat.angplanes=[];
for gi=1:length(grainl.id)
for gj=gi+1:length(grainl.id)
grain.stat.distlines=[grain.stat.distlines,...
gt2LinesDist([grainl.ca(gi,:) grainl.dir(gi,:)],...
[grainl.ca(gj,:) grainl.dir(gj,:)])];
%
% ACV=ACM{grainl.thetatype(gi),grainl.thetatype(gj)};
% ang=gt2PlanesAngle(grainl.pl(gi,:),grainl.pl(gj,:));
% ang=min(abs(ACV-ang));
% grain.stat.angplanes=[grain.stat.angplanes, ang];
end
end
%
% % Outliers
% grain.stat.outl.int=gtOutliers(grainl.int', tol_outl, grain.stat.intmean, grain.stat.intstd);
% grain.stat.outl.bbxs=gtOutliers(grainl.bbxs', tol_outl, grain.stat.bbxsmean, grain.stat.bbxsstd);
% grain.stat.outl.bbys=gtOutliers(grainl.bbys', tol_outl, grain.stat.bbysmean, grain.stat.bbysstd);
% grain.stat.outl.distcom=gtOutliers(grain.stat.distcom, tol_outl, grain.stat.distcommean, grain.stat.distcomstd);
%
% grain.stat.outl.int=grain.pairid(grain.stat.outl.int);
% grain.stat.outl.bbxs=grain.pairid(grain.stat.outl.bbxs);
% grain.stat.outl.bbys=grain.pairid(grain.stat.outl.bbys);
% grain.stat.outl.distcom=grain.pairid(grain.stat.outl.distcom);
end % of function
function newgrain=gtIndexCreateOutput(grain,whichg,tot,strategy_r,spacegroup,ACM,latticepar,pairtable,flag_update)
for i=whichg
% Create grain structure
glines=gtIndexSelectLines(grain{i}.lid,tot);
grain{i}=gtIndexCreateGrainOutputFinal(glines,ACM,spacegroup,latticepar,strategy_r);
grain{i}.id=i;
% Get difspot ID-s
difAID=[];
difBID=[];
for j=1:length(grain{i}.pairid)
mysqlcmd=sprintf('SELECT difAID,difBID FROM %s WHERE pairID=%d',pairtable,grain{i}.pairid(j));
[difAID(j),difBID(j)]=mym(mysqlcmd);
end
grain{i}.difspots=[difAID, difBID];
end
% Order grains by approximate volume
bbys=gtAllGrainValues(grain,'stat','bbysmean',1);
bbxs=gtAllGrainValues(grain,'stat','bbxsmean',1);
M(:,1)=gtAllGrainValues(grain,'id',[],1);
M(:,2)=bbys.*bbxs;
M=sortrows(M,-2);
for i=1:length(grain)
newgrain{i}=grain{M(i,1)};
newgrain{i}.id=i;
end
% Update database table
if flag_update
for i=1:length(newgrain)
for j=1:length(newgrain{i}.pairid)
mysqlcmd=sprintf('UPDATE %s SET grainID=%d WHERE pairID=%d',pairtable,i,newgrain{i}.pairid(j));
mym(mysqlcmd);
end
end
end
end
function strategy=gtIndexDefaultStrategy()
% No. of iterations (linear extension of tolerances from 'beg' to 'end'):
strategy.iter=5;
% Tolerances for building grains in the first iteration:
strategy.b.beg.ang=0.2; % angular difference in degrees
strategy.b.beg.int=100; % max. intensity ratio
strategy.b.beg.bbxs=5; % max. bbox x size ratio
strategy.b.beg.bbys=1.4; % max. bbox y size ratio
strategy.b.beg.distsf=3; % size factor for maximum distance
strategy.b.beg.distmax=2; % max. absolut distance
strategy.b.beg.ming=4; % min. no. of Friedel pairs in a grain
strategy.b.beg.geo=10; % safety margin around sample volume
% Tolerances for building grains in the last iteration:
strategy.b.end.ang=1.2;
strategy.b.end.int=100;
strategy.b.end.bbxs=5;
strategy.b.end.bbys=1.4;
strategy.b.end.distsf=3;
strategy.b.end.distmax=20;
strategy.b.end.ming=4;
strategy.b.end.geo=10;
% Tolerances for merging grains in the first iteration:
strategy.m.beg.bbxs=strategy.b.beg.bbxs; % max. bbox x size ratio
strategy.m.beg.bbys=strategy.b.beg.bbys; % max. bbox y size ratio
strategy.m.beg.distsf=2; % size factor for maximum distance
strategy.m.beg.int=strategy.b.beg.int; % max. intensity ratio
% Tolerances for merging grains in the last iteration:
strategy.m.end.bbxs=strategy.b.end.bbxs;
strategy.m.end.bbys=strategy.b.end.bbys;
strategy.m.end.distsf=2;
strategy.m.end.int=strategy.b.end.int;
% Tolerances for assigning Friedel pairs to existing grains:
%strategy.a=strategy.b;
% Tolerance for marking outliers in final grain sets (in std):
strategy.r=3;
strategy.s=3;
end
function [grain,remaining]=gtIndexExecuteStrategy(strategy,grain,tot,sample,ACM,spacegroup,latticepar)
% Define indexing strategy:
% b = build
% m = merge
% l = loosen
% a = add
% f = fit singles
if ~isfield(strategy,'proc')
strategy.proc=[];
end
% Creat strategy.proc () if not given
if isempty(strategy.proc)
if strategy.iter<=1
strategy.iter=1;
strategy.proc{1,1}='b';
strategy.proc{2,1}='m';
strategy.proc{3,1}='a';
else
strfields=fieldnames(strategy.b.beg);
for i=1:length(strfields)
dif=strategy.b.end.(strfields{i})-strategy.b.beg.(strfields{i});
strategy.l.b.(strfields{i})=dif/(strategy.iter-1);
end
strfields=fieldnames(strategy.m.beg);
for i=1:length(strfields)
dif=strategy.m.end.(strfields{i})-strategy.m.beg.(strfields{i});
strategy.l.m.(strfields{i})=dif/(strategy.iter-1);
end
sp{1,1}='b';
sp{2,1}='m';
sp{3,1}='l';
sp{4,1}='a';
for i=1:strategy.iter-1
strategy.proc=[strategy.proc; sp];
end
strategy.proc{end+1}='b';
strategy.proc{end+1}='m';
strategy.proc{end+1}='a';
end
end
tol.b=strategy.b.beg;
tol.a=strategy.b.beg;
tol.m=strategy.m.beg;
%tol.r=strategy.r.beg;
%tol.f=strategy.f.beg;
remaining=1:length(tot.pairid); % list of all lines not yet grouped
for i=1:length(strategy.proc)
step=strategy.proc{i};
action=step(1);
switch action
case 'l' % Loosen tolerance values
strfields=fieldnames(strategy.b.beg);
for j=1:length(strfields)
tol.b.(strfields{j})=tol.b.(strfields{j})+strategy.(step).b.(strfields{j});
end
tol.a=tol.b;
strfields=fieldnames(strategy.m.beg);
for j=1:length(strfields)
tol.m.(strfields{j})=tol.m.(strfields{j})+strategy.(step).m.(strfields{j});
end
case 'b' % Build grains
[grainB,remaining]=gtIndexBuildGrains(remaining,tot,sample,ACM,tol.(step));
grain=[grain,grainB];
for j=1:length(grain)
grain{j}.id=j;
end
case 'm' % Merge grains
if ~isempty(grain)
grain=gtIndexMergeGrains(grain,tot,tol.(step),spacegroup,latticepar);
end
case 'a' % Add lines to grains
if ~isempty(grain)
[grain,remaining]=gtIndexAddToGrains(remaining,grain,tot,sample,ACM,tol.(step));
end
end
end
end % of function
function [goods,grainoutput]=gtIndexFindGroupInCand(cand,isecs,tol,sample,ACM)
% cand: all the lines that possibly belong to the same grain as the 1st line
% (they are pair-consistent)
% 1st line in "cand" is the actual line investigated
% isecs: intersections (length(cand)-1))
nof_cand=length(cand.id);
grainoutput=[];
goods=false(nof_cand,1);
goods(1)=true;
for c1=2:nof_cand % loop to find a 2nd line of the grain
basel=gtIndexSelectLines([1,c1],cand);
for c2=c1+1:nof_cand % loop to find a 3rd line of the grain
tryl=gtIndexSelectLines(c2,cand);
ok=gtIndexCheckGroupCons(tryl,basel,isecs(c1-1,:),tol,sample,ACM);
if ok
% loop to find a 4th,5th,6th,... line of the grain recursively
cs=sfAddLineToGroup([1,c1,c2],cand,nof_cand,tol,sample,ACM,1);
if length(cs)>=tol.ming % grain was found and cs output contains all lines
csf=sfCheckRest(cs,cand,nof_cand,tol,sample,ACM);
goods(csf)=true;
grainl=gtIndexSelectLines(csf,cand); % lines in the grain
grainoutput=gtIndexCreateGrainOutputBuild(grainl);
return
end
end
end
end
end
function csnew=sfAddLineToGroup(cs,cand,nof_cand,tol,sample,ACM,mylevel)
basel=gtIndexSelectLines(cs,cand);
baseisec=pofintersectexp([basel.ca basel.dir])';
csnew=cs;
for i=cs(end)+1:nof_cand
tryl=gtIndexSelectLines(i,cand);
ok=gtIndexCheckGroupCons(tryl,basel,baseisec,tol,sample,ACM);
if ok
if (length(cs)+1)>=tol.ming
csnew=[cs,i];
return
else
csext=sfAddLineToGroup([cs,i],cand,nof_cand,tol,sample,ACM,mylevel+1);
if length(csext)>=tol.ming
csnew=csext;
return
end
end
end
end
end
function csf=sfCheckRest(cs,cand,nof_cand,tol,sample,ACM)
basel=gtIndexSelectLines(cs,cand);
baseisec=pofintersectexp([basel.ca basel.dir])';
csf=cs;
for i=max(cs):nof_cand
tryl=gtIndexSelectLines(i,cand);
ok=gtIndexCheckGroupCons(tryl,basel,baseisec,tol,sample,ACM);
if ok
csf=[csf,i];
basel=gtIndexSelectLines(csf,cand);
baseisec=pofintersectexp([basel.ca basel.dir])';
end
end
end
function [match,d]=gtIndexMatchGrains(grain1,grain2,cent,tol_matchgrains)
if ~exist('cent','var')
cent=false;
end
% Relative position: vectorfrom pos. 1 to pos. 2
%dc=[0 0 0]
disp('Applied relative drift:')
dc=[9.2443, -13.7827, 16.2131]
% bug: uniqueness is not fulfilled
if ~exist('tol_matchgrains','var')
tol_matchgrains.dist=20; % 50
tol_matchgrains.Rdist=0.05; % 0.05
tol_matchgrains.bbxs=3; % 3
tol_matchgrains.bbys=2; % 3
tol_matchgrains.int=1e10; % 1e10
end
nof_grains1=length(grain1);
%nof_grains2=length(grain2);
if cent
cz1=gtAllGrainValues(grain1,'center',[],3);
ids1=(1:length(grain1))';
S=sortrows([ids1,cz1],2);
cgrain1=[];
for i=floor(nof_grains1/3):nof_grains1*2/3
cgrain1=[cgrain1,grain1(S(i,1))];
end
grain1=cgrain1;
nof_grains1=length(grain1);
end
match=NaN(nof_grains1,2);
match(:,1)=(1:nof_grains1)';
% Load relevant info into the gr1 structure (vectors and matrices):
gr1.id=gtAllGrainValues(grain1,'id',[],1);
gr1.center(:,1)=gtAllGrainValues(grain1,'center',[],1);
gr1.center(:,2)=gtAllGrainValues(grain1,'center',[],2);
gr1.center(:,3)=gtAllGrainValues(grain1,'center',[],3);
gr1.R_vector(:,1)=gtAllGrainValues(grain1,'R_vector',[],1);
gr1.R_vector(:,2)=gtAllGrainValues(grain1,'R_vector',[],2);
gr1.R_vector(:,3)=gtAllGrainValues(grain1,'R_vector',[],3);
gr1.bbxs(:,1)=gtAllGrainValues(grain1,'stat','bbxsmean',1);
gr1.bbys(:,1)=gtAllGrainValues(grain1,'stat','bbysmean',1);
gr1.int(:,1)=gtAllGrainValues(grain1,'stat','intmean',1);
% Load relevant info into the gr2 structure (vectors and matrices):
gr2.id=gtAllGrainValues(grain2,'id',[],1);
gr2.center(:,1)=gtAllGrainValues(grain2,'center',[],1);
gr2.center(:,2)=gtAllGrainValues(grain2,'center',[],2);
gr2.center(:,3)=gtAllGrainValues(grain2,'center',[],3);
gr2.R_vector(:,1)=gtAllGrainValues(grain2,'R_vector',[],1);
gr2.R_vector(:,2)=gtAllGrainValues(grain2,'R_vector',[],2);
gr2.R_vector(:,3)=gtAllGrainValues(grain2,'R_vector',[],3);
gr2.bbxs(:,1)=gtAllGrainValues(grain2,'stat','bbxsmean',1);
gr2.bbys(:,1)=gtAllGrainValues(grain2,'stat','bbysmean',1);
gr2.int(:,1)=gtAllGrainValues(grain2,'stat','intmean',1);
for i=1:nof_grains1
if i==85
%keyboard
end
act=sfSelectGrains(i,gr1); % actual grain
%tol_matchgrains.dist=min([act.bbxs, act.bbys])/2;
tomatch=sfCheckMatching(act,gr2,dc,tol_matchgrains); % grains to be matched
if ~isempty(tomatch)
ltm=length(tomatch);
if ltm>1
%actminsize=min([act.bbxs, act.bbys]);
%dcenter=sum((repmat(act.center,ltm,1)-gr2.center(tomatch,:)).^2,2)/(actminsize^2);
dR_vector=sum((repmat(act.center+dc,ltm,1)-gr2.R_vector(tomatch,:)).^2,2);
%dbbxs=sum((act.bbxs-gr2.bbxs(tomatch,:)).^2,2)/(act.bbxs^2);
%dbbys=sum((act.bbys-gr2.bbys(tomatch,:)).^2,2)/(act.bbys^2);
%dint=sum((act.int-gr2.int(tomatch,:)).^2,2)/(act.int^2);
%dd=dcenter+dR_vector+dbbxs+dbbys+dint;
dd=dR_vector;
sortM=[];
sortM(:,1)=dd;
sortM(:,2)=tomatch;
sortM(:,3)=gr2.id(tomatch);
sortM=sortrows(sortM,1);
tomatch=sortM(:,2);
disp(['Multiple matches found for grain #' num2str(i) ':'])
disp(tomatch)
end
match(i,2)=tomatch(1);
end
end
for i=1:nof_grains1
if isnan(match(i,2))
d.dist(i,:)=NaN(1,3);
d.Rdist(i,:)=NaN(1,3);
d.int(i,1)=NaN;
d.bbys(i,1)=NaN;
d.bbxs(i,1)=NaN;
else
d.dist(i,:)=gr2.center(match(i,2),:)-gr1.center(i,:)-dc;
d.Rdist(i,:)=gr2.R_vector(match(i,2),:)-gr1.R_vector(i,:);
d.bbxs(i,1)=gr2.bbxs(match(i,2))/gr1.bbxs(i);
d.bbys(i,1)=gr2.bbys(match(i,2))/gr1.bbys(i);
d.int(i,1)=gr2.int(match(i,2))/gr1.int(i);
end
end
mean(d.dist(~isnan(d.dist(:,1)),:))
disp('Matches found:')
disp(sum(~isnan(match(:,2))));
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% SUB_FUNCTIONS
function grk=sfSelectGrains(keeplist,gr)
grk.id=gr.id(keeplist);
grk.int=gr.int(keeplist);
grk.center=gr.center(keeplist,:);
grk.R_vector=gr.R_vector(keeplist,:);
grk.bbys=gr.bbys(keeplist);
grk.bbxs=gr.bbxs(keeplist);
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function tomatch=sfCheckMatching(act,cand,dc,tol_matchgrains)
if act.id==38
%keyboard
end
tomatch=(1:length(cand.id))'; % initial indices of grain candidates to be merged
% Delete those indices that don't meet the following constraints:
% Bounding box Y size close enough?
cons=(cand.bbys>act.bbys/tol_matchgrains.bbys) & ...
(cand.bbys<tol_matchgrains.bbys*act.bbys);
cand=sfSelectGrains(cons,cand);
tomatch(~cons)=[];
% Spatial distance of centers close enough?
dvec=repmat(act.center+dc,length(cand.id),1)-cand.center;
dist=sqrt(sum(dvec.*dvec,2));
cons=dist<tol_matchgrains.dist;
cand=sfSelectGrains(cons,cand);
tomatch(~cons)=[];
% Bounding box X size close enough?
cons=(cand.bbxs>act.bbxs/tol_matchgrains.bbxs) & ...
(cand.bbxs<tol_matchgrains.bbxs*act.bbxs);
cand=sfSelectGrains(cons,cand);
tomatch(~cons)=[];
% Intensity close enough?
cons=(cand.int>act.int/tol_matchgrains.int) & ...
(cand.int<tol_matchgrains.int*act.int);
cand=sfSelectGrains(cons,cand);
tomatch(~cons)=[];
% Rodriguez vector close enough?
dRvec=repmat(act.R_vector,length(cand.id),1)-cand.R_vector;
Rdist=sqrt(sum(dRvec.*dRvec,2));
cons=Rdist<tol_matchgrains.Rdist;
%cand=sfSelectGrains(cons,cand);
tomatch(~cons)=[];
% cons=false(length(cand.id),1);
% for i=1:length(cand.id)
%
% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% % [pof, good_plnorms, dsqr]=plot_rodrigues_consistancy_test([l1.pl;l2.pl],[l1.hkl;l2.hkl],0);
% [Rvec, good_plnorms]=gtRodriguesTest([act.pl;grain{cand.id(i)}.pl],[act.hkl;grain{cand.id(i)}.hkl],0,spacegroup,0);
%
% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% cons(i)=all(good_plnorms);
% end
% tomatch(~cons)=[];
end
function [match,d]=gtIndexMatchGrains_sab(grain1,grain2,cent,tol_matchgrains)
% Relative position: vectorfrom pos. 1 to pos. 2
dc=[10 10 10]
%dc =[0 0 0]
%dc=[9.2443, -13.7827, 16.2131]
% bug: uniqueness is not fulfilled
if ~exist('tol_matchgrains','var')
tol_matchgrains.dist=50; % 50
tol_matchgrains.Rdist=0.25; % 0.05
tol_matchgrains.bbxs=50; % 3
tol_matchgrains.bbys=50; % 3
tol_matchgrains.int=100e10; % 1e10
end
nof_grains1=length(grain1);
%nof_grains2=length(grain2);
if cent
cz1=gtAllGrainValues(grain1,'center',[],3);
ids1=(1:length(grain1))';
S=sortrows([ids1,cz1],2);
cgrain1=[];
for i=floor(nof_grains1/3):nof_grains1*2/3
cgrain1=[cgrain1,grain1(S(i,1))];
end
grain1=cgrain1;
nof_grains1=length(grain1);
end
match=NaN(nof_grains1,2);
match(:,1)=(1:nof_grains1)';
% Load relevant info into the gr1 structure (vectors and matrices):
gr1.id=gtAllGrainValues(grain1,'id',[],1);
gr1.center(:,1)=gtAllGrainValues(grain1,'center',[],1);
gr1.center(:,2)=gtAllGrainValues(grain1,'center',[],2);
gr1.center(:,3)=gtAllGrainValues(grain1,'center',[],3);
gr1.R_vector(:,1)=gtAllGrainValues(grain1,'R_vector',[],1);
gr1.R_vector(:,2)=gtAllGrainValues(grain1,'R_vector',[],2);
gr1.R_vector(:,3)=gtAllGrainValues(grain1,'R_vector',[],3);
gr1.bbxs(:,1)=gtAllGrainValues(grain1,'stat','bbxsmean',1);
gr1.bbys(:,1)=gtAllGrainValues(grain1,'stat','bbysmean',1);
gr1.int(:,1)=gtAllGrainValues(grain1,'stat','intmean',1);
% Load relevant info into the gr2 structure (vectors and matrices):
gr2.id=gtAllGrainValues(grain2,'id',[],1);
gr2.center(:,1)=gtAllGrainValues(grain2,'center',[],1);
gr2.center(:,2)=gtAllGrainValues(grain2,'center',[],2);
gr2.center(:,3)=gtAllGrainValues(grain2,'center',[],3);
gr2.R_vector(:,1)=gtAllGrainValues(grain2,'R_vector',[],1);
gr2.R_vector(:,2)=gtAllGrainValues(grain2,'R_vector',[],2);
gr2.R_vector(:,3)=gtAllGrainValues(grain2,'R_vector',[],3);
gr2.bbxs(:,1)=gtAllGrainValues(grain2,'stat','bbxsmean',1);
gr2.bbys(:,1)=gtAllGrainValues(grain2,'stat','bbysmean',1);
gr2.int(:,1)=gtAllGrainValues(grain2,'stat','intmean',1);
for i=1:nof_grains1
if i==85
%keyboard
end
act=sfSelectGrains(i,gr1); % actual grain
%tol_matchgrains.dist=min([act.bbxs, act.bbys])/2;
tomatch=sfCheckMatching(act,gr2,dc,tol_matchgrains); % grains to be matched
if ~isempty(tomatch)
ltm=length(tomatch);
if ltm>1
%actminsize=min([act.bbxs, act.bbys]);
%dcenter=sum((repmat(act.center,ltm,1)-gr2.center(tomatch,:)).^2,2)/(actminsize^2);
dR_vector=sum((repmat(act.center+dc,ltm,1)-gr2.R_vector(tomatch,:)).^2,2);
%dbbxs=sum((act.bbxs-gr2.bbxs(tomatch,:)).^2,2)/(act.bbxs^2);
%dbbys=sum((act.bbys-gr2.bbys(tomatch,:)).^2,2)/(act.bbys^2);
%dint=sum((act.int-gr2.int(tomatch,:)).^2,2)/(act.int^2);
%dd=dcenter+dR_vector+dbbxs+dbbys+dint;
dd=dR_vector;
sortM=[];
sortM(:,1)=dd;
sortM(:,2)=tomatch;
sortM(:,3)=gr2.id(tomatch);
sortM=sortrows(sortM,1);
tomatch=sortM(:,2);
disp(['Multiple matches found for grain #' num2str(i) ':'])
disp(tomatch)
end
match(i,2)=tomatch(1);
end
end
for i=1:nof_grains1
if isnan(match(i,2))
d.dist(i,:)=NaN(1,3);
d.Rdist(i,:)=NaN(1,3);
d.int(i,1)=NaN;
d.bbys(i,1)=NaN;
d.bbxs(i,1)=NaN;
else
d.dist(i,:)=gr2.center(match(i,2),:)-gr1.center(i,:)-dc;
d.Rdist(i,:)=gr2.R_vector(match(i,2),:)-gr1.R_vector(i,:);
d.bbxs(i,1)=gr2.bbxs(match(i,2))/gr1.bbxs(i);
d.bbys(i,1)=gr2.bbys(match(i,2))/gr1.bbys(i);
d.int(i,1)=gr2.int(match(i,2))/gr1.int(i);
end
end
mean(d.dist(~isnan(d.dist(:,1)),:))
disp('Matches found:')
disp(sum(~isnan(match(:,2))));
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% SUB_FUNCTIONS
function grk=sfSelectGrains(keeplist,gr)
grk.id=gr.id(keeplist);
grk.int=gr.int(keeplist);
grk.center=gr.center(keeplist,:);
grk.R_vector=gr.R_vector(keeplist,:);
grk.bbys=gr.bbys(keeplist);
grk.bbxs=gr.bbxs(keeplist);
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function tomatch=sfCheckMatching(act,cand,dc,tol_matchgrains)
if act.id==38
%keyboard
end
tomatch=(1:length(cand.id))'; % initial indices of grain candidates to be merged
% Delete those indices that don't meet the following constraints:
% Bounding box Y size close enough?
cons=(cand.bbys>act.bbys/tol_matchgrains.bbys) & ...
(cand.bbys<tol_matchgrains.bbys*act.bbys);
cand=sfSelectGrains(cons,cand);
tomatch(~cons)=[];
% Spatial distance of centers close enough?
dvec=repmat(act.center+dc,length(cand.id),1)-cand.center;
dist=sqrt(sum(dvec.*dvec,2));
cons=dist<tol_matchgrains.dist;
cand=sfSelectGrains(cons,cand);
tomatch(~cons)=[];
% Bounding box X size close enough?
cons=(cand.bbxs>act.bbxs/tol_matchgrains.bbxs) & ...
(cand.bbxs<tol_matchgrains.bbxs*act.bbxs);
cand=sfSelectGrains(cons,cand);
tomatch(~cons)=[];
% Intensity close enough?
cons=(cand.int>act.int/tol_matchgrains.int) & ...
(cand.int<tol_matchgrains.int*act.int);
cand=sfSelectGrains(cons,cand);
tomatch(~cons)=[];
% Rodriguez vector close enough?
dRvec=repmat(act.R_vector,length(cand.id),1)-cand.R_vector;
Rdist=sqrt(sum(dRvec.*dRvec,2));
cons=Rdist<tol_matchgrains.Rdist;
%cand=sfSelectGrains(cons,cand);
tomatch(~cons)=[];
% cons=false(length(cand.id),1);
% for i=1:length(cand.id)
%
% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% % [pof, good_plnorms, dsqr]=plot_rodrigues_consistancy_test([l1.pl;l2.pl],[l1.hkl;l2.hkl],0);
% [Rvec, good_plnorms]=gtRodriguesTest([act.pl;grain{cand.id(i)}.pl],[act.hkl;grain{cand.id(i)}.hkl],0,spacegroup,0);
%
% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% cons(i)=all(good_plnorms);
% end
% tomatch(~cons)=[];
end
function newgrain=gtIndexMergeGrains(grain,tot,tol_merge,spacegroup,latticepar)
disp(' ')
disp(' MERGING GRAINS...')
disp(' ')
disp(tol_merge)
tic
nof_newgrains=0;
nof_grains=length(grain);
remaining=1:nof_grains;
% Load relevant info into the gr structure (vectors and matrices):
gr.id=gtAllGrainValues(grain,'id',[],1);
gr.center(:,1)=gtAllGrainValues(grain,'center',[],1);
gr.center(:,2)=gtAllGrainValues(grain,'center',[],2);
gr.center(:,3)=gtAllGrainValues(grain,'center',[],3);
%gr.R_vector(:,1)=gtAllGrainValues(grain,'R_vector',[],1);
%gr.R_vector(:,2)=gtAllGrainValues(grain,'R_vector',[],2);
%gr.R_vector(:,3)=gtAllGrainValues(grain,'R_vector',[],3);
gr.bbxs(:,1)=gtAllGrainValues(grain,'stat','bbxsmean',1);
gr.bbys(:,1)=gtAllGrainValues(grain,'stat','bbysmean',1);
gr.int(:,1)=gtAllGrainValues(grain,'stat','intmean',1);
while length(remaining)>=1
act=sfSelectGrains(remaining(1),gr); % actual grain
cand=sfSelectGrains(remaining(2:end),gr); % remaining candidates for merge
tol_merge.dist=min([act.bbxs, act.bbys])/tol_merge.distsf;
tomerge=sfCheckMerge(act,cand,grain,tot,tol_merge,spacegroup,latticepar); % grains to be merged
nof_newgrains=nof_newgrains+1;
if isempty(tomerge)
newgrain{nof_newgrains}=grain{act.id}; % grain is unchanged
newgrain{nof_newgrains}.id=nof_newgrains;
else
lineids=grain{act.id}.lid;
for i=1:length(tomerge)
lineids=[lineids, grain{cand.id(tomerge(i))}.lid];
end
lineids=lineids';
newglines=gtIndexSelectLines(lineids,tot); % reload original lines in the new grain
newgrain{nof_newgrains}=gtIndexCreateGrainOutputMerge(newglines);
newgrain{nof_newgrains}.id=nof_newgrains;
%newgrain{nof_newgrains}.merged=true;
disp(['Grains merged as new grain #' num2str(nof_newgrains) ':'])
disp(act.id)
disp(cand.id(tomerge(:)))
end
remaining(1)=[];
remaining(tomerge)=[];
% for i=1:length(todel)
% remaining(remaining==todel(i))=[];
% end
%
% if nof_grains>=dongrains
% break
% end
end
disp(' ')
toc
end
function grk=sfSelectGrains(keeplist,gr)
grk.id=gr.id(keeplist);
grk.int=gr.int(keeplist);
grk.center=gr.center(keeplist,:);
%grk.R_vector=gr.R_vector(keeplist,:);
grk.bbys=gr.bbys(keeplist);
grk.bbxs=gr.bbxs(keeplist);
end
function tomerge=sfCheckMerge(act,cand,grain,tot,tol_merge,spacegroup,latticepar)
tomerge=(1:length(cand.id))'; % initial indices of grain candidates to be merged
% Delete those indices that don't meet the following constraints:
% Bounding box Y size close enough?
cons=(cand.bbys>act.bbys/tol_merge.bbys) & ...
(cand.bbys<tol_merge.bbys*act.bbys);
cand=sfSelectGrains(cons,cand);
tomerge(~cons)=[];
% Spatial distance of centres close enough?
dvec=repmat(act.center,length(cand.id),1)-cand.center;
dist=sqrt(sum(dvec.*dvec,2));
cons=dist<tol_merge.dist;
cand=sfSelectGrains(cons,cand);
tomerge(~cons)=[];
% Bounding box X size close enough?
cons=(cand.bbxs>act.bbxs/tol_merge.bbxs) & ...
(cand.bbxs<tol_merge.bbxs*act.bbxs);
cand=sfSelectGrains(cons,cand);
tomerge(~cons)=[];
% Intensity close enough?
cons=(cand.int>act.int/tol_merge.int) & ...
(cand.int<tol_merge.int*act.int);
cand=sfSelectGrains(cons,cand);
tomerge(~cons)=[];
% Rodriguez vector close enough?
% dRvec=repmat(act.R_vector,length(cand.id),1)-cand.R_vector;
% Rdist=sqrt(sum(dRvec.*dRvec,2));
% cons=Rdist<tol_merge.Rdist;
% cand=sfSelectGrains(cons,cand);
% tomerge(~cons)=[];
l1=gtIndexSelectLines(grain{act.id}.lid,tot);
cons=false(length(cand.id),1);
for i=1:length(cand.id)
l2=gtIndexSelectLines(grain{cand.id(i)}.lid,tot);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% [pof, good_plnorms, dsqr]=plot_rodrigues_consistancy_test([l1.pl;l2.pl],[l1.hkl;l2.hkl],0);
[Rvec, good_plnorms]=gtRodriguesTest([l1.pl;l2.pl],[l1.hkl;l2.hkl],spacegroup,latticepar);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
cons(i)=all(good_plnorms);
end
%cand=sfSelectGrains(cons,cand);
tomerge(~cons)=[];
end
% to be deleted
% new macro: gtIndexSelectLines
function linekept=gtIndexSelectLines(keeplist,line)
linekept.id=line.id(keeplist);
linekept.pairid=line.pairid(keeplist);
linekept.aid=line.aid(keeplist);
linekept.bid=line.bid(keeplist);
linekept.theta=line.theta(keeplist);
linekept.eta=line.eta(keeplist);
linekept.omega=line.omega(keeplist);
linekept.int=line.int(keeplist);
linekept.bbxs=line.bbxs(keeplist);
linekept.bbys=line.bbys(keeplist);
linekept.ca=line.ca(keeplist,:);
linekept.cb=line.cb(keeplist,:);
linekept.dir=line.dir(keeplist,:);
linekept.pl=line.pl(keeplist,:);
linekept.thetatype=line.thetatype(keeplist);
linekept.hkl=line.hkl(keeplist,:);
%linekept.lRv=line.lRv(:,:,keeplist);
end
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment