Signed-off-by: Wolfgang Ludwig <wolfgang.ludwig@esrf.fr>

zUtil_Topotomo/calc_pole_tilts_id11_index2.m

+function [phx, phy] = calc_pole_tilts_id11_index(parameters, gid, thetatypes)
+% assumes initial tilts of the goniometer to be at 0, 0  - uses ID11
+% goniometer geometry with maximum tilt of 20 deg for samrx and 10 deg for samry
+
+load(sprintf('%s/4_grains/phase_01/index.mat',parameters.acq.dir));
+maxvalue = tand(20);
+
+%% Sample reference system
+LabX = [1 0 0]'; samgeo.dirx = [1 0 0];
+LabY = [0 1 0]'; samgeo.diry = [0 1 0];
+LabZ = [0 0 1]'; samgeo.dirz = [0 0 1];
+
+%% Directions of instrument axis at sample rotation omega = 0 (...due to sample omega = 0 corresponds to instrument rotation = -90 deg)
+instrgeo.dirx = [  0 -1  0];
+instrgeo.diry = [  1  0  0];
+instrgeo.dirz = [  0  0  1]; 
+
+%% Rotation Matrix  (based on angle & axis  -> Rodrigues Rotation Formula)  v_rotated = R(angle, axis) * v
+
+R = @(angle,axis)gtMathsRotationTensor(angle, gtMathsRotationMatrixComp(axis, 'col'));
+ 
+%% Loop through grains and their reflections...
+
+for i = 1:length(gid)
+    
+    if exist(sprintf('%s/4_grains/phase_01/grain_%04d.mat',parameters.acq.dir, gid(i)), 'file')
+        g=load(sprintf('%s/4_grains/phase_01/grain_%04d.mat',parameters.acq.dir, gid(i)));
+    else
+        g=gtCalculateGrain(grain{gid(i)},parameters);
+    end
+    
+    ind = find(g.allblobs.pl(:,3) > maxvalue & ismember(g.allblobs.thetatype, thetatypes));
+    
+    if ~isempty(ind)
+        [pl,ia,ic] = unique(g.allblobs.pl(ind,:), 'rows', 'stable');
+        ind = ind(ia);  
+      
+        for j = 1:length(ind)
+
+            hkl = g.allblobs.hkl(ind(j), :); 
+            theta = g.allblobs.theta(ind(j),:);
+            G_sam  = g.allblobs.pl(ind(j),:);
+            diffrz = g.allblobs.omega(ind(j))-90;
+            uv = g.allblobs.detector.uvw(ind(j),1:2);
+
+            %% GInstr =  Plane normal in Instrument reference system:
+            % samrx is upper tilt and rotates around LabX at diffrz = 0 (omega = 90)
+            % ramry is lower tile and rotates around LabY at diffrz = 0 (omega = 90)
+
+            G_instr = gtGeoSam2Sam(G_sam, samgeo, instrgeo, 1, 1);
+
+            % first solve for upper tilt:  RX * G_sam = [r 0 s]   ->  G_sam(2)*cos(phx) - G_sam(3)*sin(phx) = 0
+            phx = atand(G_instr(2)/G_instr(3));
+            tmp = R(phx, LabX) * G_instr';
+
+            % now solve for lower tilt: RY * tmp = [0 0 1], RY = [cos(phy) 0  sin(phy); 0 1 0;  -sin(phy) 0 cos(phy)];  -> 
+            % phy = atand(-tmp(1)/tmp(3)); 
+
+            phy = atand(-tmp(1)/tmp(3)); 
+     
+            % testing
+            % G_aligned = R(phy, LabY) * tmp;
+            if (abs(phy) < 10 & abs(phx) <20)
+                disp(sprintf('Grain %d:  Found %d %d %d reflection:  diffry = %f samrx = %f   samry = %f\n', gid(i), hkl(1), hkl(2), hkl(3), theta, phx, phy));
+                disp(sprintf('Positions in DCT scan: diffrz = %f, u = %f v = %f\n', diffrz, uv(1), uv(2)));
+            end
+        end
+    end
+end
+
+        
+        
+