calculate sample goniometer tilts for a grain with known orientation

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

zUtil_Topotomo/calc_pole_tilts_id06_index2.m

+function [samry, samrz] = calc_pole_tilts_id06_index2(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
+name = sprintf('%s/4_grains/phase_01/index.mat',parameters.acq.dir);
+load(name);
+maxvalue = cosd(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 (...depends on ftomo settings...!!)
+instrgeo.dirx = [  1  0  0];
+instrgeo.diry = [  0  1  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),:);
+            diffrx = g.allblobs.omega(ind(j));
+            image  = g.allblobs.omega(ind(j)) / parameters.labgeo.omstep;
+            uv = g.allblobs.detector.uvw(ind(j),1:2);
+
+            %% GInstr =  Plane normal in Instrument reference system:
+            % samrz is upper tilt and rotates around  -LabY at diffrx = 0 
+            % samry is lower tilt and rotates around   LabX at diffrx = 0 
+
+            G_instr = gtGeoSam2Sam(G_sam, samgeo, instrgeo, 1, 1);
+
+            % first solve for upper tilt: RUT * G_sam = [0 s r]; (RUT = [0 -1 0])
+            % ->  G_sam(1)*cos(ut) - G_sam(3)*sin(ut) = 0
+            ut = atand(G_instr(1)/G_instr(3));    
+            tmp = R(ut, -LabY) * G_instr';
+
+            % now solve for lower tilt: RLT * tmp = [0 0 1];  (RLT = [1 0 0]) 
+            % cos(lt)*G_sam(2) - sin(lt)*G_sam(3)=0
+            
+            lt = atand(tmp(2)/tmp(3)); 
+     
+            % now translate this to the ID06 instrument axes...
+            samrz = ut;
+            samry = lt;
+            % testing
+            % G_aligned = R(lt, LabX) * tmp
+            
+            
+            
+            if (abs(ut) < 20 & abs(lt) <20)
+                disp(sprintf('Grain %d:  Found %d %d %d %d reflection:  diffrz = %f mv samry %f   samrz %f\n', gid(i), hkl(1), hkl(2), hkl(3), hkl(4), theta, lt, ut));
+                disp(sprintf('Positions in DCT scan: image number = %d, u = %f v = %f\n', round(image), uv(1), uv(2)));
+            end
+        end
+    end
+end
+
+        
+        
+