Added gtIndexDifSpots: test version for fast reindexing of difspots.

Signed-off-by: Peter Reischig <peter.reischig@esrf.fr>

zUtil_Indexter/gtIndexDifSpots.m

+function gid = gtIndexDifSpots(grain, sp, parameters)
+% GTINDEXDIFSPOTS Finds grain ID-s for all diffraction spots. Function 
+% under development!!
+%
+% gid = gtIndexDifSpots(grain, sp, parameters)
+%
+% ---------------------------------------------------------
+% Under development! NOT properly tested yet!
+%
+% It tries to assign a grain ID to each diffraction spot. Those for which 
+% no acceptable grain is found or compatible with more than one grain 
+% (conflict) will have zero assigned.
+%
+% Run these functions beforehand:
+%   index = load('4_grains/phase_01/index.mat')
+%   grain = gtIndexAddInfo(index.grain, phaseid, parameters);
+%   grain = gtIndexFwdSimGrains(grain, parameters);
+%
+% INPUT
+%   grain      - all grain data (1xn or nx1 cell array) as in index.mat
+%
+%   Optional:
+%   sp         - diffraction spot data; if undefined or empty, 
+%                will be loaded from database
+%   parameters - as in parameters file; if undefined or empty, will be
+%                loaded from parameters file
+%
+% OUTPUT
+%   gid        - list of grain ID-s
+%
+ 
+
+%%%%%%%%%%%%%%%%%%
+%% Prepare input
+%%%%%%%%%%%%%%%%%%
+
+% Load parameters
+if ~exist('parameters','var') || isempty(parameters)
+    parameters = load('parameters.mat');
+    parameters = parameters.parameters;
+end
+
+% Get spot and pair data
+if ~exist('sp','var') || isempty(sp)
+    if isfield(parameters.acq, 'difA_name')
+        difspottable = [parameters.acq.difA_name 'difspot'];
+    else
+        difspottable = [parameters.acq.name 'difspot'];
+    end
+    
+    fprintf('Loading diff. spot data from difspottable %s ...\n', difspottable)
+    gtDBConnect;
+    mysqlcmd = sprintf('SELECT difspotID,CentroidX,CentroidY,CentroidImage,BoundingBoxXsize,BoundingBoxYsize,Integral FROM %s',difspottable);
+    [sp.id, sp.cu, sp.cv, sp.cw, sp.bbx, sp.bby, sp.int] = mym(mysqlcmd);
+end
+
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%% Vectorise fwd simulated grain spots 
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+gr.cu  = [];
+gr.cv  = [];
+gr.cw  = [];
+gr.id  = [];
+
+for ii = 1:length(grain)
+   gr.cu  = [gr.cu; grain{ii}.fsim.cu'];
+   gr.cv  = [gr.cv; grain{ii}.fsim.cv'];
+   gr.cw  = [gr.cw; grain{ii}.fsim.cw'];
+   gr.id  = [gr.id; grain{ii}.id(ones(1,length(grain{ii}.fsim.cu)))'];
+end
+    
+%%%%%%%%%%%%%%%%%%%%%%%
+%% Set tolerance limits
+%%%%%%%%%%%%%%%%%%%%%%%
+
+% Gather all grain data
+bbxsrat = gtIndexAllGrainValues(grain,'stat','bbxsrat',[],[]);
+bbysrat = gtIndexAllGrainValues(grain,'stat','bbysrat',[],[]);
+intrat  = gtIndexAllGrainValues(grain,'stat','intrat',[],[]);
+uvwA    = gtIndexAllGrainValues(grain,'fsimA','uvw',[],[]);
+uvwB    = gtIndexAllGrainValues(grain,'fsimB','uvw',[],[]);
+centA   = gtIndexAllGrainValues(grain,'centA',[],[],[]);
+centB   = gtIndexAllGrainValues(grain,'centB',[],[],[]);
+centimA = gtIndexAllGrainValues(grain,'centimA',[],[],[]);
+centimB = gtIndexAllGrainValues(grain,'centimB',[],[],[]);
+
+% Deviations between fwd simulation and measured spot positions
+du = [uvwA(1,:) - centA(1,:), uvwB(1,:) - centB(1,:)];
+dv = [uvwA(2,:) - centA(2,:), uvwB(2,:) - centB(2,:)];
+dw = [uvwA(3,:) - centimA,    uvwB(3,:) - centimB];
+
+% Set limits within 3std
+lim.u = mean(du) + 3*std(du);
+lim.v = mean(dv) + 3*std(dv);
+lim.w = mean(dw) + 3*std(dw);
+lim.bbxhi = mean(bbxsrat) + 3*std(bbxsrat);
+lim.bbyhi = mean(bbysrat) + 3*std(bbysrat);
+lim.inthi = mean(intrat)  + 3*std(intrat);
+lim.bbxlo = 1/lim.bbxhi; 
+lim.bbylo = 1/lim.bbyhi; 
+lim.intlo = 1/lim.inthi; 
+
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%
+%% Loop to find grain ID-s
+%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+gid  = zeros(1,length(sp.cu));
+conf = false(1,length(sp.cu));
+ok   = false(length(gr.cu),3);
+
+fprintf('Processing %d spots: 0000000',length(sp.cu))
+cnt = 0;
+
+tic
+
+% Loop through spots - could be a parallel parfor loop
+for ii = 1:length(sp.cu)
+    
+    % Display progress
+    cnt = cnt + 1;
+    if cnt == 100
+       fprintf('\b\b\b\b\b\b\b%07d',ii)
+       cnt = 0;
+    end
+    
+    % Check if position is within limits
+    ok(:,1) = abs(gr.cu - sp.cu(ii)) < lim.u;
+    ok(:,2) = abs(gr.cv - sp.cv(ii)) < lim.v;
+    ok(:,3) = abs(gr.cw - sp.cw(ii)) < lim.w;
+        
+    oks = all(ok, 2);
+    sumoks = sum(oks);
+    
+    % Check metadata, if exactly one candidate grain found
+    if sum(oks) == 1
+       
+        idg = gr.id(find(oks,1,'first'));
+        
+        bbxrat = sp.bbx(ii)/grain{idg}.stat.bbxsmean;
+        bbyrat = sp.bby(ii)/grain{idg}.stat.bbysmean;
+        intrat = sp.int(ii)/grain{idg}.stat.intmean;
+        
+        ok2(1) = bbxrat <= lim.bbxhi;
+        ok2(2) = bbxrat >= lim.bbxlo;
+        ok2(3) = bbyrat <= lim.bbyhi;
+        ok2(4) = bbyrat >= lim.bbylo;
+        ok2(5) = intrat <= lim.inthi;
+        ok2(6) = intrat >= lim.intlo;
+        
+        if all(ok2)
+            gid(ii) = idg;
+        end
+    elseif sumoks > 1
+        % spot in conflict 
+        conf(ii) = true; 
+    end
+end
+
+fprintf('\b\b\b\b\b\b\bdone.\n')
+
+toc
+
+fprintf('No. of spots analysed   : %d\n', length(gid))
+fprintf('No. of spots indexed    : %d\n', sum(gid>0))
+fprintf('No. of spots in conflict: %d\n', sum(conf))
+
+
+end % of function
\ No newline at end of file