Commands we use in this page: mcstas
In this document, we demonstrate how
iFit/
iOptim and
iData can
all be used transparently to simulate and optimize
McStas instrument models.
To import a McStas simulation results, just use
>> results = iData('single_detector_file')
>> results = iData('directory') % import everything in the directory, including sim files (files may imported more than once)
>> results = iData('scan_directory/mcstas.dat') % for scans
>> results = iData('scan_directory/mcstas.sim') % for all files from the directory, uniquely.
which returns a single or an array of
iData
object(s), that you can plot with e.g.
subplot.
The resulting McStas data has an additional
Parameters field which holds the instrument parameters used for the simulation.
>> results(1).Parameters
Running a single
simulation
As iFit contains powerful import routines, that support McStas data
files, it is straight forward to launch a McStas instrument simulation
from Matlab and get back the simulation results. To do so, you will
obviously need to have
iFit
installed (refer to the
Install page), as
well as
McStas (refer to the
McStas web site). An easy solution is
to install our
ILL/CS
Live
DVD ready-to-run system with
McStas
pre-installed, and then put on top
Matlab
(requires a valid license).
The syntax for running a McStas instrument model is:
>> results = mcstas('instrument', parameters, options);
where
parameters is a
structure that holds parameter names and values, e.g.
parameters.RV=1
parameters.lambda=2.36
parameters.powder='Na2Ca3Al2F14.laz'
and both numerical values and strings are supported. Numerical values
should be single scalars, whereas strings are used as-is, that is they
may assign DEFINITION instrument parameters such as file names, arrays
and other types that are sent to the commands line. The unspecified
instrument parameters are used with their default values (when defined).
The parameters can also be given as a single string with members separated with the ';' character ans assigned with '=':
>> results = mcstas('instrument', 'RV=1; lambda=2.36; powder=Na2Ca3Al2F14.laz');
When the instrument is not found in the current directory, it is
searched in the Mcstas library, and eventually copied locally, which
triggers its compilation.
The last argument
options
controls simulation configuration.
- options.dir:
directory
where
to
store results. When not defined, a temporary
location is created, and removed after the simulation (string).
- options.ncount: number of
neutron events per iteration, e.g. 1e5
(double)
- options.mpi:
number
of
processors to use with MPI (integer. The simulation should be
re-compiled accordingly with options.compile='yes'
- options.seed:
random number seed to use for each iteration (double)
- options.gravitation: 0 or
1 to set gravitation handling in neutron propagation (boolean)
- options.mode:
should be 'simulate'. This is automatically set when no
'optimization' configuration parameter is set.
- options.compile: 0 or 1
to force re-compilation of the instrument (boolean)

In
the
following
examples, we define required structures on the fly, run
the
templateDIFF instrument
(the file must be accessible from Matlab current path) and get back
monitors:
>> monitors=mcstas('templateDIFF', struct('RV',1), struct('ncount',1e6))
monitors = array [1 5] iData object:
Index [Tag] [Dimension] [Title] [Last command] [Label]
1 id445199 [40 40] 'File Diff_Mono_XY_1299773321...' id445199=load(iData,... Diff_Mono_XY
2 id445212 [40 1] 'File Diff_Mono_Lambda_129977...' id445212=load(iData,... Diff_Mono_...
3 id445224 [40 1] 'File Diff_Sample_Lambda_1299...' id445224=load(iData,... Diff_Sampl...
4 id445236 [340 1] 'File Diff_BananaTheta_129977...' id445236=load(iData,... Diff_Banan...
5 id445248 [170 25] 'File Diff_BananaPSD_12997733...' id445248=load(iData,... Diff_Banan...
>> subplot(monitors)
The last line displays the simulation results. Usual data analysis
(display, mathematics, fit, export) can be performed (refer to the
iFit and
iData
documentation).
The resulting McStas data has an additional
Parameters field which holds the instrument parameters used for the simulation.
>> monitors(1).Parameters
Displaying the instrument geometry
To display the instrument geometry, use the syntax:
>> results = mcstas('instrument', 'RV=1','mode=display');
which will launch 'mcdisplay' in the background, and open the resulting
figure. No neutron trajectory is shown though, only the instrument
geometry.
Running series of
simulations (scans)
It is very simple to perform parameter scans, even in multi-dimensional
spaces. For this, the parameter values just need to be vectors of
numeric or vector cells, such as
- parameters.RV={0.5 1 1.5};
- parameters.L2= 2:0.5:4;
- parameters.L3= 1.3;
Scans with string parameters are also possible, such as with
- parameters.Powder={'Na2Ca3Al2F14.laz','Al.laz','Cd.laz'};
and then executes for instance:
>> [integral,monitors]=mcstas('templateDIFF', parameters, struct('ncount',1e4));
A live display of the on-going scan is shown when specifying
OutputFcn=
'fminplot' in the last
options argument (struct).
The
integral argument
returned is an
iData object containing the
integral of monitors as a function of the scanned parameters. It can be
directly plotted with e.g.
>> plot(integral); % simple plot, or median surface for multidimensional scans
>> scatter3(integral); % coloured points (possibly in 2D,3D)
>> plot3(integral); % coloured lines in 2D, or volume for multidimensional scans
>> slice(integral); % only for 3D, use slice-o-matic for volume inspection

The
monitors
argument returned is an array containing individual scan monitors.
You may extract and catenate a selection (slices) or it by e.g. the
iData/cat method.
In the following example, we launch a quick 2D scan:
>> [p,m]=mcstas('templateDIFF', struct('RV',[0.7 0.9 1.2],'L2',[1 1.3 1.5]),struct('ncount',1e4));
Then we may plot, as a 2D surface, the integral value for the last
(5th) monitor:
>> plot(p(:,:,5));

You
may
look
at the instrument parameters used for the scans, e.g.
>> get(m,'RV')
which are also available in the
Parameters field of the objects.
Let's finally look at the evolution of the 4-th monitor, which is a 1D
diffractogram, along an 'RV' parameter scan:
>> a = squeeze(m(2,:,:)); % extract an 'RV' line, with L2=1.3
>> b = a(:,4); % get 4th monitor on this line
>> setaxis(b,2,'RV'); % define a new axis along the line RV
>> c = cat(2, b); % assemble into a single iData object
>> plot(c,'tight interp'); % plot
This methodology also works with 2D monitors assembled along, e.g. a
scan line, which creates a 3D volume to explore with
slice, or
plot or
plot3.
Warning: Beware
of
the
number
of iterations required to scan large dimensionality spaces.
In practice you should avoid spaces above 2 or 3 scanned parameters.
Exceeding this will undoubtedly require very long computation times,
and large memory storage.
Optimizing instrument
parameters
The optimization of instrument parameters is performed as simply as a
single simulation, just by e.g. setting
options.mode to
'optimize'.
>> [solution, monitors] = mcstas('instrument', parameters, options);
The parameters, as a named structure, are initiated to their starting
value, from which the optimization proceeds.
Only numerical parameters can be
optimized. Also, it is
highly encouraged to use bounded parameters (that is define constraints, see
below), else the optimizer can get lost, produce unmeaningful results, or return
NaN parameter values.
More specifically, the
options
structure may contain, in addition to the optimization mode, the
following members:
- options.type:
'minimize' or 'maximize', which is the
default (string)
- options.monitors:
cell string of monitor names, or empty for all (cellstr)
- options.optimizer:
function name of the optimizer to use, default is fminpso (string or function
handle). See iOptim.
- options.mode:
should be 'optimize'.
This is automatically set whenever an 'optimization' configuration
field is set, such as monitors,
optimizer, TolFun, ...
as well as any other optimizer configuration fields such as
- options.TolFun ='1%'; % stop when criteria changes are smaller
that 1%, this is the default
- options.TolX ='0.1%'; % stop when parameter changes are smaller
that 0.1%
- options.Display='final';
- options.OutputFcn='fminplot';
% plots the criteria and parameters evolution during the
optimization. This is the default. Set to '' to override this choice. A view of the monitors used as criteria is also displayed.
The default optimizer configuration is used when not specified (see the