gtCrystCalculateReflections : modified output (now as structure) with all the...

gtCrystCalculateReflections : modified output (now as structure) with all the info inside (to be imported in cryst), remove 'typelist' from input arguments (the calculation is made for both the unique and full lists of reflections)

Signed-off-by: Laura Nervo <laura.nervo@esrf.fr>

git-svn-id: https://svn.code.sf.net/p/dct/code/trunk@960 4c865b51-4357-4376-afb4-474e03ccb993

zUtil_Cryst/gtCrystCalculateReflections.m

-function [hkl, theta, hkil, sintheta] = gtCrystCalculateReflections(latticepar, mintheta, maxtheta, sg, typelist, filename)
+function list = gtCrystCalculateReflections(latticepar, mintheta, maxtheta, sg, filename)
 % GTCRYSTCALCULATEREFLECTIONS  Calculate the reflections list using python 
 %                              from xfab library in fable
-%     [hkl, theta, hkil, sintheta] = gtCrystCalculateReflections(latticepar, mintheta, maxtheta, sg, typelist, filename)
-%     ------------------------------------------------------------------------------------------------------------------
+%     list = gtCrystCalculateReflections(latticepar, mintheta, maxtheta, sg, filename)
+%     --------------------------------------------------------------------------------
 %
 %     INPUT: 
 %       latticepar      = lattice parameters <double 1x6>
@@ -10,15 +10,21 @@ function [hkl, theta, hkil, sintheta] = gtCrystCalculateReflections(latticepar,
 %       maxtheta        = maximum sin(theta)/lambda allowed <double>
 %       sg              = spacegroup number <int>
 %       typelist        = type of list <string> [{'all'}/'unique']
-%       filename        = filename for output <string>
 %
 %     OUTPUT:
-%       hkl      = reflection list <double>
-%       theta    = theta value for each reflection <double>
-%       hkil     = hkil indexes if needed <double>
-%       sintheta = sind(theta)/lambda <double>
+%       list = structure with fields (the same with 'sp' for the full lists)
+%           .hkl       = reflection list <double>
+%           .hkil      = hkil indexes if needed <double>
+%           .theta     = theta value for each reflection <double>
+%           .thetatype = theta type <double>
+%           .sintheta  = sind(theta)/lambda <double>
+%           .dspacing  = d-spacing <double>
+%           .mult      = multiplicity of each reflections <double>
 %
 %
+%     Version 002 12-12-2012 by LNervo
+%       Compact output as a structure
+%
 %     Version 001 03-12-2012 by LNervo
 
 
@@ -39,49 +45,78 @@ if ~exist('sg','var') || isempty(sg)
     load parameters;
     sg = parameters.cryst(1).spacegroup;
 end
-if ~exist('typelist','var') || isempty(typelist)
-    typelist = 'all';
-end
 if ~exist('filename','var') || isempty(filename)
     filename = 'reflections_';
 end
 
-latticepar_str = sprintf('%f ', latticepar);
+parameters = [];
+load('parameters.mat');
+energy = parameters.acq.energy;
+
+latticepar_str = sprintf('%f ',latticepar);
+latticepar_str = strtrim(latticepar_str);
+% python script to be used
 script_file = fullfile('/users/',getenv('USER'),'matlabDCT','zUtil_Python','reflections_list.py');
+
+% run the command for all the reflections
 cmd = sprintf('python26 %s %s %f %f %d %s %s',...
-              script_file,latticepar_str,mintheta,maxtheta,sg,typelist,filename);
-[~,msg]=unix(cmd);
-disp(msg)
+              script_file,latticepar_str,mintheta,maxtheta,sg,'all',filename);
+[~,msg]=unix(cmd); disp(msg)
 
-[C,Cmat] = gtReadTextFile([filename typelist '.txt'],'%f %f %f %f','\n','#',[1 3]);
-hkl = Cmat(:,1:3);
-sinth = Cmat(:,4);
+% run the command for only the reflections unique
+cmd = sprintf('python26 %s %s %f %f %d %s %s',...
+              script_file,latticepar_str,mintheta,maxtheta,sg,'unique',filename);
+[~,msg]=unix(cmd); disp(msg)
 
-if nargout > 1
-    load parameters;
-    energy = parameters.acq.energy;
-    theta = asind(sinth)*gtConvEnergyToWavelength(energy);
-    if nargout > 2
-      % for planes with Miller Indexes hkl -> hkil
-      hkil = hkl;
-      hkil(:,4) = hkl(:,3);
-      hkil(:,3) = -hkil(:,1)-hkil(:,2);
+% read the full list of reflections
+[~,Cmat] = gtReadTextFile([filename 'all.txt'],'%f %f %f %f','\n','#',[1 3]);
 
-      disp('  h   k   i   l   theta')
-      disp('-------------------------')
-      for ii=1:size(hkl,1)
-          fprintf('%3d %3d %3d %3d %7.3f\n',hkil(ii,:),theta(ii))
-      end
-      if nargout > 3
-          sintheta = sinth;
-      end
-    end
-else
-    disp('  h   k   l   sin(theta)/lambda')
-    disp('---------------------------------')
-    for ii=1:size(hkl,1)
-        fprintf('%3d %3d %3d %7.3f\n',hkl(ii,:),sinth(ii))
-    end
+% extracting information
+hklsp = Cmat(:,1:3);
+sinthetasp = Cmat(:,4);
+% getting 4-indexes notation
+hkilsp = hklsp;
+hkilsp(:,4) = hkilsp(:,3);
+hkilsp(:,3) = -hkilsp(:,1)-hkilsp(:,2);
+
+% extracting theta
+thetasp = asind(sinthetasp)*gtConvEnergyToWavelength(energy);
+% calculation of d-spacing
+Bmat = gtCrystHKL2CartesianMatrix(latticepar);
+dspacingsp = gtCrystDSpacing(hklsp',Bmat)';
+% calculation of the multiplicity
+[sintheta,ind,thetatypesp] = unique(sinthetasp,'stable');
+mult = [ind(2:end)',length(thetatypesp)+1]-ind';
+
+hkl       = hklsp(ind,:);
+theta     = thetasp(ind);
+hkil      = hkilsp(ind,:);
+thetatype = thetatypesp(ind);
+dspacing  = dspacingsp(ind);
+
+disp('  h   k   i   l     theta   thetatype  mult    d-spacing')
+disp('----------------------------------------------------------')
+for ii=1:size(ind,1)
+    fprintf( '%3d %3d %3d %3d %9.3f %5d %10d %10.3f\n',hkilsp(ind(ii),:),thetasp(ind(ii)),thetatypesp(ind(ii)),mult(ii),dspacingsp(ind(ii)) )
 end
 
+% creating the output
+list.latticepar  = latticepar;
+list.spacegroup  = sg;
+list.energy      = energy;
+list.hkl         = hkl';
+list.hklsp       = hklsp';
+list.hkil        = hkil';
+list.hkilsp      = hkilsp';
+list.theta       = theta';
+list.thetasp     = thetasp';
+list.thetatype   = thetatype';
+list.thetatypesp = thetatypesp';
+list.sintheta    = sintheta';
+list.sinthetasp  = sinthetasp';
+list.dspacing    = dspacing';
+list.dspacingsp  = dspacingsp';
+list.mult        = mult;
+
+
 end % end of function