In order to make full use of the Matlab interface to S-GeMS some knowledge of the use of S-GeMS is essential. The book Applied Geostatistics with SGeMS (Remy, Boucher and Wu, Cambridge University Press, 2009), written by the developers of S-GeMS is highly recommended.
The Matlab interface to S-GeMS relies on a feature of S-GeMS that allow S-GeMS to read and execute a series of Python commands from the command line, without the need to load the graphical user interface, as for example:
sgems -s sgems_python_script.py
The Matlab interface consists of methods and functions to automatically create such a Python script, execute the script using S-GeMS and load the simulated/estimated results into Matlab
One function (Section 4.7, “sgems_grid”) handles these actions allowing simulation on grids in the following manner:
Define a parameter file (Section 4.6, “sgems_get_par”, Section 4.12, “sgems_read_xml”)
Write a python script that (Section 4.8, “sgems_grid_py”)
sets up a grid or pointset where simulation or estimation is performed
performs the simulation/estimation
export the simulated/estimated data
Load the data into Matlab (Section 4.14, “sgems_write”)
For a complete list of S-GeMS related commands on mGstat see Section 4, “S-GeMS functions”
This section contains a rather detailed explanation of using S-GeMS to perform simulation. Much more compact example can be found in the following chapters.
Unconditional and conditional sequential simulation can be performed using Section 4.7, “sgems_grid” :
S = sgems_grid(S);
Where S is a Matlab data structure containing all the information needed to setup and run S-GeMS
A number of different simulation algorithms are available in S-GeMS The behavior of each algorithm is controlled through an XML file. Such an XML file can for example be exported from S-GeMS by choosing to save a parameter file for a specific algorithm.
Such an XML formatted parameter is needed to perform any kind of simulation. A number of 'default' parameter files available using the Section 4.6, “sgems_get_par” function. For example to obtain a default parameter file for sequential Gaussian simulation use
S = sgems_get_par('sgsim')
S =
xml_file: 'sgsim.par'
XML: [1x1 struct]
As can be seen the adds the name of the XML file (S.xml_file) as well as a XML data structure in the S-GeMS matlab structure S.XML.
All supported simulation/estimation types can be found calling sgems_get_par without arguments:
>> sgems_get_par sgems_get_par : available SGeMS type dssim sgems_get_par : available SGeMS type filtersim_cate sgems_get_par : available SGeMS type filtersim_cont sgems_get_par : available SGeMS type lusim sgems_get_par : available SGeMS type sgsim sgems_get_par : available SGeMS type snesim_std
Now all parameters for 'sgsim' simulation can be set directly from the Matlab command line. To see the number of fields in the XML file (refer to the S-GeMS book described above for the meaning of all parameters):
>> S.XML.parameters
ans =
algorithm: [1x1 struct]
Grid_Name: [1x1 struct]
Property_Name: [1x1 struct]
Nb_Realizations: [1x1 struct]
Seed: [1x1 struct]
Kriging_Type: [1x1 struct]
Trend: [1x1 struct]
Local_Mean_Property: [1x1 struct]
Assign_Hard_Data: [1x1 struct]
Hard_Data: [1x1 struct]
Max_Conditioning_Data: [1x1 struct]
Search_Ellipsoid: [1x1 struct]
Use_Target_Histogram: [1x1 struct]
nonParamCdf: [1x1 struct]
Variogram: [1x1 struct]
To see the number of realization :
>> S.XML.parameters.Nb_Realizations
ans =
value: 10
To set the number of realization to 20 do:
>> S.XML.parameters.Nb_Realizations.value=20;
One also need to define the grid used for simulation. This is done through the S.dim data structure:
%grid size S.dim.nx=70; S.dim.ny=60; S.dim.nz=1; % grid cell size S.dim.dx=1; S.dim.dy=1; S.dim.dz=1; % grid origin S.dim.x0=0; S.dim.y0=0; S.dim.z0=0;
All the values listed above for the S.dim data structure are default, thus if they are not set, they are assumed as listed.
Unconditional simulation is now performed using:
>> S=sgems_grid(S);
sgems_grid : Trying to run SGeMS using sgsim.py, output to SGSIM.out
'import site' failed; use -v for traceback
Executing script...
working on realization 1
|# | 5%working on realization 2
|## | 10%working on realization 3
|### | 15%working on realization 4
|#### | 20%working on realization 5
|##### | 25%working on realization 6
|###### | 30%working on realization 7
|####### | 35%working on realization 8
|######## | 40%working on realization 9
|######### | 45%working on realization 10
|########## | 50%working on realization 11
|########### | 55%working on realization 12
|############ | 60%working on realization 13
|############# | 65%working on realization 14
|############## | 70%working on realization 15
|############### | 75%working on realization 16
|################ | 80%working on realization 17
|################# | 85%working on realization 18
|################## | 90%working on realization 19
|################### | 95%working on realization 20
|####################| 100%
sgems_read : Reading GRID data from SGSIM.sgems
sgems_grid : SGeMS ran successfully
S =
xml_file: 'sgsim.par'
XML: [1x1 struct]
dim: [1x1 struct]
data: [4200x20 double]
O: [1x1 struct]
x: [1x70 double]
y: [1x60 double]
z: 1
D: [4-D double]
As seen above the following field have been added to the S-GeMS matlab structure:
S.x,
S.y,
S.z,
S.data and
S.D.
S.x,
S.y,
S.z are 3 arrays defining the grid.
S.data, is the simulated data as exported from S-GeMS. Note the each realization is returned as a list of size nx*ny*nz.
S.D, is but a rearrangement of S.data into a 4D dimensional data structure, of size (nx,ny,nz,nsim). To visualize for example the 3rd realization use for example:
imagesc(S.x,S.y,S.D(:,:,1,3));
Conditional simulation can be performed by setting the S.d_obs parameter. For example:
S.d_obs=[18 13 0 0; 5 5 0 1; 2 28 0 1]; S=sgems_grid(S); imagesc(S.x,S.y,S.D(:,:,1,3));
A simple example of unconditional sequential simulation.
S=sgems_get_par('sgsim');
S.XML.parameters.Nb_Realizations.value=12;
S=sgems_grid(S);
for i=1:S.XML.parameters.Nb_Realizations.value;
subplot(4,3,i);
imagesc(S.x,S.y,S.D(:,:,1,i));
end
A simple example of conditional sequential simulation.
S=sgems_get_par('sgsim');
% Define observed data=
S.d_obs=[18 13 0 0; 5 5 0 1; 2 28 0 1];
S.XML.parameters.Nb_Realizations.value=12;
S=sgems_grid(S);
for i=1:S.XML.parameters.Nb_Realizations.value;
subplot(4,3,i);
imagesc(S.x,S.y,S.D(:,:,1,i));axis image;
end
A simple example of unconditional SNESIM AND FILTERSIM simulation.
S1=sgems_get_par('snesim_std'); %
% Note that S1.ti_file is automatically set.
% simply change this to point to another training to use.
S1.XML.parameters.Nb_Realizations.value=4;
S2=sgems_get_par('filtersim_cont');
S2.XML.parameters.Nb_Realizations.value=4;
S1=sgems_grid(S1);
S2=sgems_grid(S2);
for i=1:S1.XML.parameters.Nb_Realizations.value;
subplot(S1.XML.parameters.Nb_Realizations.value,2,i);
imagesc(S1.x,S1.y,S1.D(:,:,1,i));axis image;
subplot(S1.XML.parameters.Nb_Realizations.value,2,i+S1.XML.parameters.Nb_Realizations.value);
imagesc(S2.x,S2.y,S2.D(:,:,1,i));axis image;
end
Using a JPG file from FLICKR as training image:
% LOAD IMAGE AND CONVERT TO SGeMS binary TRAINING IMAGE
file_img='1609350318_7300f07360_m.jpg'; % pattern
f_out=sgems_image2ti(file_img);
TI=sgems_read(f_out);
% SETUP FILTERSIM
S=sgems_get_par('filtersim_cont');
S.ti_file=f_out;
S.XML.parameters.PropertySelector_Training.grid=TI.grid_name;
S.XML.parameters.PropertySelector_Training.property=TI.property{1};
S.XML.parameters.Nb_Realizations.value=1;
S.dim.x=[1:1:200];
S.dim.y=[1:1:200];
S.dim.z=[0];
% RUN SIMULATION
S=sgems_grid(S);
% VISUALIZE RELIZATION
subplot(1,2,1);
imagesc(TI.x,TI.y,TI.D(:,:,:,1)');axis image;
subplot(1,2,2);
imagesc(S.x,S.y,S.D(:,:,:,1)');axis image;
 |
Example of converting an image and using it for continuous filtersim simulation
Demonstration simulation of mGstat supported simulation algorithms can be performed using Section 4.5, “sgems_demo”. To see a list of supported simulation algorithms use:
>> sgems_get_par sgems_get_par : available SGeMS type dssim sgems_get_par : available SGeMS type filtersim_cate sgems_get_par : available SGeMS type filtersim_cont sgems_get_par : available SGeMS type lusim sgems_get_par : available SGeMS type sgsim sgems_get_par : available SGeMS type snesim_std
To run a demonstration of continuous filtersim simulation using the 'filtersim_cont' algorithm do
>> sgems_demo('filtersim_cont');
This will perform both unconditional and conditional simulation, and visualize the results as for example here below.
 |
Unconditional simulation
 |
Conditional simulation
 |
E-type on conditional simulations
Running Section 4.5, “sgems_demo” without arguments will run the demonstration using all supported simulation algorithms.
