Contents

% ----------------------------------------------------------%%   create_17K_2p0A%% ----------------------------------------------------------
horace_on()

clear, close all

!==================================================================!
!                      HORACE                                      !
!------------------------------------------------------------------!
!  Visualisation of multi-dimensional neutron spectroscopy data    !
!                                                                  !
!  T.G.Perring, J van Duijn, R.A.Ewings         November 2008      !
!------------------------------------------------------------------!
! Matlab  code: $Revision:: 558  $ ($Date:: 2011-07-19 11:06:05 $) !
! Mex code:    Disabled  or not supported on this platform         !
!------------------------------------------------------------------!

Define working directories

% ----------------------------------------------------------% Define working directories% ----------------------------------------------------------% Windows% dirroot  = '\\serdon-nv\illdata\133\in5\exp_4-01-1230\';% dirfinal = [dirroot,'processed\spe\'];% dirraw  =  [dirroot,'rawdata\'];% linux
dirroot  = '/usr/illdata/133/in5/exp_4-01-1230/';
dirfinal = [dirroot,'processed/spe/'];
dirraw  =  [dirroot,'rawdata/'];

indir      = dirfinal;                             % directory of spe files
par_file   = [dirfinal,'in5psd_318_384_240.par'];  % detector parameter file
sqw_file   = [dirfinal,'tmo_2p0A_17K.sqw'];        % output sqw file

Initial parameters

Lambda = 2.0;                          % Wavelenght [A]
efix   = 81.805/Lambda^2 ;             % Ei [meV]
emode  = 1;                            % Direct geometry (2: inverse)
alatt  = [5.297, 5.831, 7.403];        % Crystal parameters (from literature)
angdeg = [90.0,90.0,90.0];
u      = [0,1,0];                      % 1st reference axis: u=// to incident beam
v      = [0,0,1];                      % 2nd reference axis: v=perp to incident beam
omega=0;dpsi=0;                        % Unuseful (?)
gl=2.37;   % Estimated from the Braggs at level
gs=1.003;  %   --"--% Or left to zero and orientation corrected on .SQW files later% (see horace_rotation_corrections).

Sample rotation offset

Omega = Omega_0 when Ki is parallel to the first reference axis (u)

% ----------------------------------------------------------% Sample rotation offset (determined from Braggs positions)% ----------------------------------------------------------
Offset = 250.32 - 180.0;

Proceed

Give all sample rotation angles (as measured) and corresponding .spe files (generated in Lamp at present). Sort them (not compulsory), check Omega with the angle recorded in the raw datafiles and proceed to the projection in crystal reciprocal space (generate .sqw file).

% Sample rotation angles as given from the measurements
Omega  = [240:-2:234,232:-2:120,121:2:139,141:2:233]';
% Run numbers corresponding to the sample rotation
irun   = [110017:110020,110023:110079,110080:110089,110091:110137]';
[Omega, ind] = sort(Omega);

% Sample rotation without offset
psi = Offset - Omega;

irun = irun(ind);
nfiles = numel(irun);
spe_file = cell(1,nfiles);

% Cell matrix of spe files to usefor k=1:nfiles
    spe_file{k}=[indir,'sqw_',num2str(irun(k,1)),'_',num2str(irun(k,1)),'_318_384_240.spe'];
end% Check consistency (compare the sample rotation angle with the value in each file)for k = 1:nfiles
    rawFile = [dirraw,num2str(irun(k,1)),'.nxs'];
    SROT=0;
    SRot = hdf5read(rawFile, '/entry0/instrument/SRot/value');
    fprintf('%d)\tPsi = %5.2f ( Om = %5.2f / SRot = %5.2f) -> %s \n',k, psi(k),Omega(k),SRot,spe_file{k});
end

1)  Psi = -49.68 ( Om = 120.00 / SRot = 120.01) -> /usr/illdata/133/in5/exp_4-01-1230/processed/spe/sqw_110079_110079_318_384_240.spe 
2)  Psi = -50.68 ( Om = 121.00 / SRot = 120.99) -> /usr/illdata/133/in5/exp_4-01-1230/processed/spe/sqw_110080_110080_318_384_240.spe 
3)  Psi = -51.68 ( Om = 122.00 / SRot = 122.01) -> /usr/illdata/133/in5/exp_4-01-1230/processed/spe/sqw_110078_110078_318_384_240.spe 
4)  Psi = -52.68 ( Om = 123.00 / SRot = 122.99) -> /usr/illdata/133/in5/exp_4-01-1230/processed/spe/sqw_110081_110081_318_384_240.spe 

...

116)  Psi = -165.68 ( Om = 236.00 / SRot = 236.01) -> /usr/illdata/133/in5/exp_4-01-1230/processed/spe/sqw_110019_110019_318_384_240.spe 
117)  Psi = -167.68 ( Om = 238.00 / SRot = 238.01) -> /usr/illdata/133/in5/exp_4-01-1230/processed/spe/sqw_110018_110018_318_384_240.spe 
118)  Psi = -169.68 ( Om = 240.00 / SRot = 239.99) -> /usr/illdata/133/in5/exp_4-01-1230/processed/spe/sqw_110017_110017_318_384_240.spe 

% Generate the .SQW file
gen_sqw (spe_file, par_file, sqw_file, efix, emode, alatt, angdeg, u, v, ...
         psi , omega, dpsi, gl, gs);

Published with MATLAB® R2012b

Contents

% CORRECTION_SAMPLE_MISALIGNMENT Correcting for sample misalignments in Horace%% Corrections of the sample mis-alignments and lattice parameter% refinement. (from the script on the Horace website )% To be done after the creation of the .sqw file%% 1) Check which Bragg peak to deal with (in 3 directions)% 2) Once the notional Braggs are defined, refine the Braggs% 3) Compute the corrections% 4) Apply corrections to the .sqw file%
clear, close all

Parameters to enter

sqw_file   = 'spe/Alcur_155K.sqw';  % .sqw file to correct% projection axes
proj.u = [1,1,0];
proj.v = [0,0,1];
proj.type = 'rrr';
proj.uoffset = [0,0,0,0];
% Lattice parameters
alatt  = [8.22,8.22, 11.02];
angdeg = [90,90,120];

Look for Braggs

%Make a series slices in order to work out what Bragg positions we have.
dQ = 0.025;
H = 0; dH = dQ;
K = 0; dK = dQ;
L = 0; dL = dQ;
E = 0; dE = 0.125;
alignment_slice1=cut_sqw(sqw_file,proj, dH, dK,[L-dL, L+dL],[E-dE, E+dE],'-nopix');
alignment_slice2=cut_sqw(sqw_file,proj, [H-dH, H+dH], dK,dL,[E-dE, E+dE],'-nopix');
alignment_slice3=cut_sqw(sqw_file,proj, dH, [K-dK, K+dK],dL,[E-dE, E+dE],'-nopix');


plot(compact(alignment_slice1)); keep_figure; lc 01
plot(compact(alignment_slice2)); keep_figure; lc 01
plot(compact(alignment_slice3)); keep_figure; lc 01

Define notional Bragg peaks

Notional Bragg peaks = rlu

bp=[1,1,0; 2,2,0; 0,0,6; 1,1,2; 1,1,-2];

Refine the Bragg peak positions

%Get true Bragg peak positions

% longitudinal extend, bin width, thickness
% transversal  extend, bin width, thickness
% dE
% 
% longitudinal cuts: 
%  wl = cut_sqw(sqw_file, proj, [-longit. extend/2, longit. bin, lomngit. extend/2], ...
%                               [-transverse thick, transverse bin, +transverse thick],  ...
%                               [-transverse thick, transverse bin, +transverse thick], [-dE, +dE], '-nopix'); 
%
% transverse cuts:
%  wt1= cut_sqw(sqw_file, proj, [-longit. thick/2,+longit. thick/2], ...
%                               [-trans. extend/2, trans. bin,+trans. extend/2] , ...
%                               [-trans. thick/2,+trans. thick/2], [-dE, +dE], '-nopix');
%
%  wt2= cut_sqw(sqw_file, proj, [-longit. thick/2,+longit. thick/2], ...
%                               [-trans. thick/2,+trans. thick/2], [-dE, +dE], ...
%                               [-trans. extend/2, trans. bin,+trans. extend/2], '-nopix');
%
%

[rlu0,~,wcut,wpeak]=bragg_positions(sqw_file, bp, ...
                                        0.5, 0.0125, 0.25,...
                                        0.5, 0.0125, 0.25, 0.1, 'gauss','bin_ab');

%Check how well the function did:
bragg_positions_view(wcut,wpeak)

--------------------------------------------------------------------------------
Peak 1:  [1  1  0]    scan: 1 (radial scan)
--------------------------------------------------------------------------------
Peak 1:  [1  1  0]    scan: 2 (transverse scan)
--------------------------------------------------------------------------------
Peak 1:  [1  1  0]    scan: 3 (transverse scan)
--------------------------------------------------------------------------------
--------------------------------------------------------------------------------
Peak 2:  [2  2  0]    scan: 1 (radial scan)
--------------------------------------------------------------------------------
Peak 2:  [2  2  0]    scan: 2 (transverse scan)
-------------------------------------------------------------------------------- 
Peak 2:  [2  2  0]    scan: 3 (transverse scan)
--------------------------------------------------------------------------------
--------------------------------------------------------------------------------
Peak 3:  [0  0  6]    scan: 1 (radial scan)
--------------------------------------------------------------------------------
Peak 3:  [0  0  6]    scan: 2 (transverse scan)
--------------------------------------------------------------------------------
Peak 3:  [0  0  6]    scan: 3 (transverse scan)
--------------------------------------------------------------------------------
--------------------------------------------------------------------------------
Peak 4:  [1  1  2]    scan: 1 (radial scan)
--------------------------------------------------------------------------------
Peak 4:  [1  1  2]    scan: 2 (transverse scan)
--------------------------------------------------------------------------------
Peak 4:  [1  1  2]    scan: 3 (transverse scan)
--------------------------------------------------------------------------------
--------------------------------------------------------------------------------
Peak 5:  [1  1 -2]    scan: 1 (radial scan)
--------------------------------------------------------------------------------
Peak 5:  [1  1 -2]    scan: 2 (transverse scan)
--------------------------------------------------------------------------------
Peak 5:  [1  1 -2]    scan: 3 (transverse scan)
--------------------------------------------------------------------------------

Compute the corrections

%Determine corrections to lattice and orientation:
[rlu_corr,alatt,angdeg] = refine_crystal(rlu0, alatt, angdeg, bp,'fix_angdeg','fix_alatt_ratio');

Apply changes to sqw file

choice = questdlg('Apply changes to sqw file?', 'Confirm dialog', 'Ok', 'Cancel','Cancel');
if strcmp(choice,'Ok')
    change_crystal_horace(sqw_file, rlu_corr);
end

Published with MATLAB® R2014b

Contents

% ---------------------------------------------------------------------%  cuts% ---------------------------------------------------------------------
clear, close all% Windows% dirroot  = '\\serdon-nv\illdata\133\in5\exp_4-01-1230\';% dirfinal = [dirroot,'processed\spe\'];% linux
dirroot  = '/usr/illdata/133/in5/exp_4-01-1230/';
dirfinal = [dirroot,'processed/spe/'];

% Name of the .sqw file
data_source  = [dirfinal,'tmo_2p0A_17K.sqw'];    % output sqw file

Projection axes

Define a projection axes for the cut. Could be any orientation in the spanned reciprocal space, in particular, different from the one used to generate the .sqw file.

proj.u = [0,1,0];         % 2 perpendicular axes, the 3rd is calculated
proj.v = [0,0,1];
proj.type = 'rrr';        % units: 'r' means r.l.u, 'a' mean A-1
proj.uoffset = [0,0,0,0]; % Origin offset of the cut

3D cut in the [u,v,w] plane at fixed energy

With u=[010] and v = [001] this is volume [HKL,hw=0]

%                              u    v    w       hw
w1 = cut_sqw(data_source,proj,0.05,0.05,0.05,[-0.1,0.1]);
plot(w1)    % 3D plots with sliceomatic
view(0,0)

% Optionally smooth  before plotting
w2=d3d(w1);
w3=smooth(w2,[2,2,2]);
plot(w3);
view(0,90)

lc 01e-5% Color axis in 3D plots

3D cut in the [u,v,hw]: equatorial plane

With u=[010] and v = [001] this is volume [0KL,hw]

%                              u    v         w      hw
w1 = cut_sqw(data_source,proj,0.05,0.05,[-0.05,0.05],0);
plot(w1)    % 3D plots with sliceomatic
view(0,0)

% Optionally smooth  before plotting
w2=d3d(w1);
w3=smooth(w2,[2,2,2]);
plot(w3);
view(0,90)

3D cut in the [u,w,hw]: a vertical plane

With u=[010] and v = [001] this is volume [H0L,hw]

%                              u         v       w   hw
w1 = cut_sqw(data_source,proj,0.05,[-0.05,0.05],0.05,0);
plot(w1)    % 3D plots with sliceomatic
view(0,0)

% Optionally smooth  before plotting
w2=d3d(w1);
w3=smooth(w2,[2,2,2]);
plot(w3);
view(0,90)

2D cut in the [u,0,hw] direction

With u=[010] and v = [001] this is plane [0K0]

%                              u         v            w      hw
w1 = cut_sqw(data_source,proj,0.05,[-0.05,0.05],[-0.05,0.05],0);
plot(w1)   

% Optionally smooth  before plotting
w2=d2d(w1);
w3=smooth(w2,[2,2]);
plot(w3);

lz 01e-5% Color axis in 2D plots

2D cut in the [u,v,fixed hw]: iso-energy planes

With u=[010] and v = [001] this is plane [0KL] at hw = 1meV

%                             u     v         w          hw
w1 = cut_sqw(data_source,proj,0.05,0.05,[-0.05,0.05],[0.95,1.05]);
plot(w1)    

% Optionally smooth  before plotting
w2=d2d(w1);
w3=smooth(w2,[2,2]);
plot(w3);

1D cut in the [fixed u,fixed v, hw]: iso-energy planes

With u=[010] and v = [001] this is at position [001] u v w hw

w1 = cut_sqw(data_source,proj,[-0.05,0.05],[0.95,1.05],[-0.05,0.05],0);
plot(w1)    

Save the plot

Plot type can be 'pdf', 'jpeg', 'png', 'tif', 'eps',...

printplot2('figure.pdf',gcf,'pdf');

Published with MATLAB� R2012b