BornAgain  1.19.0
Simulate and fit neutron and x-ray scattering at grazing incidence
DepthProbeSimulation Class Reference
Inheritance diagram for DepthProbeSimulation:
[legend]
Collaboration diagram for DepthProbeSimulation:
[legend]

Public Member Functions

 DepthProbeSimulation ()
 
 ~DepthProbeSimulation () override
 
void accept (INodeVisitor *visitor) const final
 Calls the INodeVisitor's visit method. More...
 
void addParameterDistribution (const ParameterDistribution &par_distr)
 
void addParameterDistribution (const std::string &param_name, const IDistribution1D &distribution, size_t nbr_samples, double sigma_factor=0.0, const RealLimits &limits=RealLimits())
 
const IBackgroundbackground () const
 
Beambeam ()
 
const Beambeam () const
 
DepthProbeSimulationclone () const override
 
SimulationResult convertData (const OutputData< double > &data, bool put_masked_areas_to_zero=true)
 Convert user data to SimulationResult object for later drawing in various axes units. More...
 
int copyNumber (const INode *node) const
 Returns copyNumber of child, which takes into account existence of children with same name. More...
 
ParameterPoolcreateParameterTree () const
 Creates new parameter pool, with all local parameters and those of its children. More...
 
std::unique_ptr< IUnitConvertercreateUnitConverter () const
 
IDetectordetector ()
 
const IDetectordetector () const
 
std::string displayName () const
 Returns display name, composed from the name of node and it's copy number. More...
 
const IAxisgetAlphaAxis () const
 Returns a pointer to incident angle axis. More...
 
std::vector< const INode * > getChildren () const
 Returns a vector of children. More...
 
IDetectorgetDetector ()
 
const IDetectorgetDetector () const
 
const DistributionHandlergetDistributionHandler () const
 
const std::string & getName () const
 
SimulationOptionsgetOptions ()
 
const SimulationOptionsgetOptions () const
 
const IAxisgetZAxis () const
 Returns a pointer to z-position axis. More...
 
Instrumentinstrument ()
 
const Instrumentinstrument () const
 
size_t intensityMapSize () const override
 Returns the total number of the intensity values in the simulation result. More...
 
virtual void onChange ()
 Action to be taken in inherited class when a parameter has changed. More...
 
RealParameterparameter (const std::string &name) const
 Returns parameter with given 'name'. More...
 
ParameterPoolparameterPool () const
 Returns pointer to the parameter pool. More...
 
std::string parametersToString () const
 Returns multiline string representing available parameters. More...
 
INodeparent ()
 
const INodeparent () const
 
virtual void prepareSimulation ()
 Put into a clean state for running a simulation. More...
 
std::vector< const INode * > progeny () const
 Returns a vector of all descendants. More...
 
void registerChild (INode *node)
 
RealParameterregisterParameter (const std::string &name, double *parpointer)
 
void registerVector (const std::string &base_name, kvector_t *p_vec, const std::string &units="nm")
 
void removeParameter (const std::string &name)
 
void removeVector (const std::string &base_name)
 
SimulationResult result () const override
 Returns the results of the simulation in a format that supports unit conversion and export to numpy arrays. More...
 
void runMPISimulation ()
 Run a simulation in a MPI environment. More...
 
void runSimulation ()
 Run a simulation, possibly averaged over parameter distributions. More...
 
const MultiLayersample () const
 
void setBackground (const IBackground &bg)
 
void setBeamParameters (double lambda, int nbins, double alpha_i_min, double alpha_i_max, const IFootprintFactor *beam_shape=nullptr)
 Sets beam parameters with alpha_i of the beam defined in the range. More...
 
void setDetectorResolutionFunction (const IResolutionFunction2D &resolution_function)
 
void setInstrument (const Instrument &instrument_)
 
void setName (const std::string &name)
 
void setOptions (const SimulationOptions &options)
 
void setParameterValue (const std::string &name, double value)
 
virtual void setParent (const INode *newParent)
 
void setSample (const MultiLayer &sample)
 The MultiLayer object will not be owned by the ISimulation object. More...
 
void setSampleBuilder (const std::shared_ptr< ISampleBuilder > &sample_builder)
 
void setTerminalProgressMonitor ()
 Initializes a progress monitor that prints to stdout. More...
 
void setVectorValue (const std::string &base_name, kvector_t value)
 
void setZSpan (size_t n_bins, double z_min, double z_max)
 Set z positions for intensity calculations. More...
 
void subscribe (ProgressHandler::Callback_t inform)
 
virtual void transferToCPP ()
 Used for Python overriding of clone (see swig/tweaks.py) More...
 
virtual std::string treeToString () const
 Returns multiline string representing tree structure below the node. More...
 

Static Public Member Functions

static std::string XComponentName (const std::string &base_name)
 
static std::string YComponentName (const std::string &base_name)
 
static std::string ZComponentName (const std::string &base_name)
 

Protected Member Functions

const SimulationOptionsoptions () const
 
ProgressHandlerprogress ()
 
virtual void transferResultsToIntensityMap ()
 Creates the appropriate data structure (e.g. More...
 
virtual void updateIntensityMap ()
 

Protected Attributes

const size_t m_NP
 
std::vector< double > m_P
 

Private Member Functions

 DepthProbeSimulation (const DepthProbeSimulation &other)
 
void addBackgroundIntensity (size_t start_ind, size_t n_elements) override
 
void addDataToCache (double weight) override
 
void checkCache () const
 
std::unique_ptr< OutputData< double > > createIntensityData () const
 Creates intensity data from simulation elements. More...
 
std::vector< DepthProbeElementgenerateSimulationElements (const Beam &beam)
 Generate simulation elements for given beam. More...
 
std::unique_ptr< IComputationgenerateSingleThreadedComputation (size_t start, size_t n_elements) override
 Generate a single threaded computation for a given range of simulation elements. More...
 
double incidentAngle (size_t index) const
 
void initialize ()
 Initializes simulation. More...
 
void initSimulationElementVector () override
 Initializes the vector of ISimulation elements. More...
 
void moveDataFromCache () override
 
void normalize (size_t start_ind, size_t n_elements) override
 Normalize the detector counts to beam intensity, to solid angle, and to exposure angle. More...
 
size_t numberOfSimulationElements () const override
 Gets the number of elements this simulation needs to calculate. More...
 
std::vector< double > rawResults () const override
 
void runSingleSimulation (size_t batch_start, size_t batch_size, double weight=1.0)
 Runs a single simulation with fixed parameter values. More...
 
void setBeamParameters (double lambda, const IAxis &alpha_axis, const IFootprintFactor *beam_shape)
 Sets beam parameters with alpha_i of the beam defined in the range. More...
 
void setRawResults (const std::vector< double > &raw_data) override
 
void validateParametrization (const ParameterDistribution &par_distr) const override
 Checks the distribution validity for simulation. More...
 
void validityCheck () const
 Checks if simulation data is ready for retrieval. More...
 

Private Attributes

std::unique_ptr< IAxism_alpha_axis
 
std::unique_ptr< IBackgroundm_background
 
std::vector< std::valarray< double > > m_cache
 
DistributionHandler m_distribution_handler
 
Instrument m_instrument
 
std::string m_name
 
SimulationOptions m_options
 
const INodem_parent {nullptr}
 
std::unique_ptr< ParameterPoolm_pool
 parameter pool (kind of pointer-to-implementation) More...
 
ProgressHandler m_progress
 
SampleProvider m_sample_provider
 
std::vector< DepthProbeElementm_sim_elements
 
std::unique_ptr< IAxism_z_axis
 

Detailed Description

Definition at line 33 of file DepthProbeSimulation.h.

Constructor & Destructor Documentation

◆ DepthProbeSimulation() [1/2]

DepthProbeSimulation::DepthProbeSimulation ( )

Definition at line 32 of file DepthProbeSimulation.cpp.

32  : ISimulation()
33 {
34  initialize();
35 }
void initialize()
Initializes simulation.

References initialize().

Referenced by clone().

Here is the call graph for this function:

◆ ~DepthProbeSimulation()

DepthProbeSimulation::~DepthProbeSimulation ( )
overridedefault

◆ DepthProbeSimulation() [2/2]

DepthProbeSimulation::DepthProbeSimulation ( const DepthProbeSimulation other)
private

Definition at line 100 of file DepthProbeSimulation.cpp.

101  : ISimulation(other), m_sim_elements(other.m_sim_elements), m_cache(other.m_cache)
102 {
103  if (other.m_alpha_axis)
104  m_alpha_axis.reset(other.m_alpha_axis->clone());
105  if (other.m_z_axis)
106  m_z_axis.reset(other.m_z_axis->clone());
107  for (auto iter = m_sim_elements.begin(); iter != m_sim_elements.end(); ++iter)
108  iter->setZPositions(m_alpha_axis.get());
109  initialize();
110 }
std::unique_ptr< IAxis > m_z_axis
std::vector< std::valarray< double > > m_cache
std::unique_ptr< IAxis > m_alpha_axis
std::vector< DepthProbeElement > m_sim_elements

References initialize(), m_alpha_axis, m_sim_elements, and m_z_axis.

Here is the call graph for this function:

Member Function Documentation

◆ accept()

void DepthProbeSimulation::accept ( INodeVisitor visitor) const
inlinefinalvirtual

Calls the INodeVisitor's visit method.

Implements INode.

Definition at line 40 of file DepthProbeSimulation.h.

40 { visitor->visit(this); }
virtual void visit(const BasicLattice2D *)
Definition: INodeVisitor.h:151

◆ addBackgroundIntensity()

void DepthProbeSimulation::addBackgroundIntensity ( size_t  start_ind,
size_t  n_elements 
)
overrideprivatevirtual

Implements ISimulation.

Definition at line 239 of file DepthProbeSimulation.cpp.

240 {
241  if (background())
242  throw std::runtime_error(
243  "Error: nonzero background is not supported by DepthProbeSimulation");
244 }
const IBackground * background() const
Definition: ISimulation.h:74

References ISimulation::background().

Here is the call graph for this function:

◆ addDataToCache()

void DepthProbeSimulation::addDataToCache ( double  weight)
overrideprivatevirtual

Implements ISimulation.

Definition at line 246 of file DepthProbeSimulation.cpp.

247 {
248  checkCache();
249  for (size_t i = 0, size = m_sim_elements.size(); i < size; ++i)
250  m_cache[i] += m_sim_elements[i].getIntensities() * weight;
251 }

References checkCache(), m_cache, and m_sim_elements.

Here is the call graph for this function:

◆ addParameterDistribution() [1/2]

void ISimulation::addParameterDistribution ( const ParameterDistribution par_distr)
inherited

Definition at line 272 of file ISimulation.cpp.

273 {
274  validateParametrization(par_distr);
276 }
void addParameterDistribution(const std::string &param_name, const IDistribution1D &distribution, size_t nbr_samples, double sigma_factor=0.0, const RealLimits &limits=RealLimits())
add a sampled parameter distribution
virtual void validateParametrization(const ParameterDistribution &) const
Checks the distribution validity for simulation.
Definition: ISimulation.h:139
DistributionHandler m_distribution_handler
Definition: ISimulation.h:159

References DistributionHandler::addParameterDistribution(), ISimulation::m_distribution_handler, and ISimulation::validateParametrization().

Here is the call graph for this function:

◆ addParameterDistribution() [2/2]

void ISimulation::addParameterDistribution ( const std::string &  param_name,
const IDistribution1D distribution,
size_t  nbr_samples,
double  sigma_factor = 0.0,
const RealLimits limits = RealLimits() 
)
inherited

Definition at line 264 of file ISimulation.cpp.

267 {
268  ParameterDistribution par_distr(param_name, distribution, nbr_samples, sigma_factor, limits);
269  addParameterDistribution(par_distr);
270 }
void addParameterDistribution(const std::string &param_name, const IDistribution1D &distribution, size_t nbr_samples, double sigma_factor=0.0, const RealLimits &limits=RealLimits())
A parametric distribution function, for use with any model parameter.

Referenced by StandardSimulations::MiniGISASBeamDivergence(), and TransformToDomain::setBeamDistribution().

◆ background()

const IBackground* ISimulation::background ( ) const
inlineinherited

◆ beam() [1/2]

◆ beam() [2/2]

const Beam& ISimulation::beam ( ) const
inlineinherited

Definition at line 59 of file ISimulation.h.

59 { return m_instrument.beam(); }

References Instrument::beam(), and ISimulation::m_instrument.

Here is the call graph for this function:

◆ checkCache()

void DepthProbeSimulation::checkCache ( ) const
private

Definition at line 193 of file DepthProbeSimulation.cpp.

194 {
195  if (m_sim_elements.size() != m_cache.size())
196  throw std::runtime_error("Error in DepthProbeSimulation: the sizes of simulation element "
197  "vector and of its cache are different");
198 }

References m_cache, and m_sim_elements.

Referenced by addDataToCache(), and moveDataFromCache().

◆ clone()

DepthProbeSimulation * DepthProbeSimulation::clone ( ) const
overridevirtual

Implements ISimulation.

Definition at line 39 of file DepthProbeSimulation.cpp.

40 {
41  return new DepthProbeSimulation(*this);
42 }

References DepthProbeSimulation().

Here is the call graph for this function:

◆ convertData()

SimulationResult ISimulation::convertData ( const OutputData< double > &  data,
bool  put_masked_areas_to_zero = true 
)
inherited

Convert user data to SimulationResult object for later drawing in various axes units.

User data will be cropped to the ROI defined in the simulation, amplitudes in areas corresponding to the masked areas of the detector will be set to zero.

Definition at line 308 of file ISimulation.cpp.

310 {
311  auto converter = UnitConverterUtils::createConverter(*this);
312  auto roi_data = UnitConverterUtils::createOutputData(*converter, converter->defaultUnits());
313 
314  if (roi_data->hasSameDimensions(data)) {
315  // data is already cropped to ROI
316  if (put_masked_areas_to_zero) {
317  detector().iterate(
318  [&](IDetector::const_iterator it) {
319  (*roi_data)[it.roiIndex()] = data[it.roiIndex()];
320  },
321  /*visit_masked*/ false);
322  } else {
323  roi_data->setRawDataVector(data.getRawDataVector());
324  }
325 
326  } else if (detHasSameDimensions(detector(), data)) {
327  // exp data has same shape as the detector, we have to put orig data to smaller roi map
328  detector().iterate(
329  [&](IDetector::const_iterator it) {
330  (*roi_data)[it.roiIndex()] = data[it.detectorIndex()];
331  },
332  /*visit_masked*/ !put_masked_areas_to_zero);
333 
334  } else {
335  throw std::runtime_error("FitObject::init_dataset() -> Error. Detector and exp data have "
336  "different shape.");
337  }
338 
339  return SimulationResult(*roi_data, *converter);
340 }
void iterate(std::function< void(const_iterator)> func, bool visit_masks=false) const
Definition: IDetector.cpp:197
IDetector & detector()
Definition: ISimulation.h:63
std::vector< T > getRawDataVector() const
Returns copy of raw data vector.
Definition: OutputData.h:334
An iterator for SimulationArea.
Wrapper around OutputData<double> that also provides unit conversions.
std::unique_ptr< IUnitConverter > createConverter(const ISimulation &simulation)
std::unique_ptr< OutputData< double > > createOutputData(const IUnitConverter &converter, Axes::Units units)
Returns zero-valued output data array in specified units.

References UnitConverterUtils::createConverter(), UnitConverterUtils::createOutputData(), ISimulation::detector(), SimulationAreaIterator::detectorIndex(), OutputData< T >::getRawDataVector(), IDetector::iterate(), and SimulationAreaIterator::roiIndex().

Here is the call graph for this function:

◆ copyNumber()

int INode::copyNumber ( const INode node) const
inherited

Returns copyNumber of child, which takes into account existence of children with same name.

Definition at line 94 of file INode.cpp.

95 {
96  if (node->parent() != this)
97  return -1;
98 
99  int result(-1), count(0);
100  for (auto child : getChildren()) {
101 
102  if (child == nullptr)
103  throw std::runtime_error("INode::copyNumber() -> Error. Nullptr as child.");
104 
105  if (child == node)
106  result = count;
107 
108  if (child->getName() == node->getName())
109  ++count;
110  }
111 
112  return count > 1 ? result : -1;
113 }
const INode * parent() const
Definition: INode.cpp:84
virtual std::vector< const INode * > getChildren() const
Returns a vector of children.
Definition: INode.cpp:63
const std::string & getName() const

References INode::getChildren(), IParametricComponent::getName(), and INode::parent().

Referenced by INode::displayName().

Here is the call graph for this function:

◆ createIntensityData()

std::unique_ptr< OutputData< double > > DepthProbeSimulation::createIntensityData ( ) const
private

Creates intensity data from simulation elements.

Definition at line 267 of file DepthProbeSimulation.cpp.

268 {
269  std::unique_ptr<OutputData<double>> result = std::make_unique<OutputData<double>>();
270  result->addAxis(*getAlphaAxis());
271  result->addAxis(*getZAxis());
272 
273  std::vector<double> rawData;
274  rawData.reserve(getAlphaAxis()->size() * getZAxis()->size());
275  for (size_t i = 0, size = m_sim_elements.size(); i < size; ++i) {
276  const std::valarray<double>& fixed_angle_result = m_sim_elements[i].getIntensities();
277  rawData.insert(rawData.end(), std::begin(fixed_angle_result), std::end(fixed_angle_result));
278  }
279  result->setRawDataVector(rawData);
280 
281  return result;
282 }
const IAxis * getAlphaAxis() const
Returns a pointer to incident angle axis.
SimulationResult result() const override
Returns the results of the simulation in a format that supports unit conversion and export to numpy a...
const IAxis * getZAxis() const
Returns a pointer to z-position axis.
virtual size_t size() const =0
retrieve the number of bins

References getAlphaAxis(), getZAxis(), m_sim_elements, result(), and IAxis::size().

Referenced by result().

Here is the call graph for this function:

◆ createParameterTree()

ParameterPool * INode::createParameterTree ( ) const
virtualinherited

Creates new parameter pool, with all local parameters and those of its children.

Reimplemented from IParametricComponent.

Definition at line 126 of file INode.cpp.

127 {
128  std::unique_ptr<ParameterPool> result(new ParameterPool);
129 
130  for (const INode* child : progeny()) {
131  const std::string path = NodeUtils::nodePath(child, parent()) + "/";
132  child->parameterPool()->copyToExternalPool(path, result.get());
133  }
134 
135  return result.release();
136 }
Base class for tree-like structures containing parameterized objects.
Definition: INode.h:49
std::vector< const INode * > progeny() const
Returns a vector of all descendants.
Definition: INode.cpp:68
Container with parameters for IParametricComponent object.
Definition: ParameterPool.h:29
std::string nodePath(const INode *node, const INode *root=nullptr)
Returns path composed of node's displayName, with respect to root node.
Definition: NodeUtils.cpp:89

References NodeUtils::nodePath(), INode::parent(), and INode::progeny().

Referenced by ISimulation::runSimulation(), validateParametrization(), OffSpecularSimulation::validateParametrization(), and SpecularSimulation::validateParametrization().

Here is the call graph for this function:

◆ createUnitConverter()

std::unique_ptr< IUnitConverter > DepthProbeSimulation::createUnitConverter ( ) const

Definition at line 95 of file DepthProbeSimulation.cpp.

96 {
97  return std::make_unique<DepthProbeConverter>(beam(), *m_alpha_axis, *m_z_axis);
98 }
Beam & beam()
Definition: ISimulation.h:58

References ISimulation::beam(), m_alpha_axis, and m_z_axis.

Referenced by result().

Here is the call graph for this function:

◆ detector() [1/2]

◆ detector() [2/2]

const IDetector& ISimulation::detector ( ) const
inlineinherited

Definition at line 64 of file ISimulation.h.

64 { return m_instrument.detector(); }

References Instrument::detector(), and ISimulation::m_instrument.

Here is the call graph for this function:

◆ displayName()

std::string INode::displayName ( ) const
inherited

Returns display name, composed from the name of node and it's copy number.

Definition at line 115 of file INode.cpp.

116 {
117  std::string result = getName();
118  if (m_parent) {
119  int index = m_parent->copyNumber(this);
120  if (index >= 0)
121  result = result + std::to_string(index);
122  }
123  return result;
124 }
const INode * m_parent
Definition: INode.h:83
int copyNumber(const INode *node) const
Returns copyNumber of child, which takes into account existence of children with same name.
Definition: INode.cpp:94

References INode::copyNumber(), IParametricComponent::getName(), and INode::m_parent.

Referenced by NodeUtils::nodePath().

Here is the call graph for this function:

◆ generateSimulationElements()

std::vector< DepthProbeElement > DepthProbeSimulation::generateSimulationElements ( const Beam beam)
private

Generate simulation elements for given beam.

Definition at line 152 of file DepthProbeSimulation.cpp.

153 {
154  std::vector<DepthProbeElement> result;
155 
156  const double wavelength = beam.wavelength();
157  const double angle_shift = beam.direction().alpha();
158 
159  const size_t axis_size = getAlphaAxis()->size();
160  result.reserve(axis_size);
161  for (size_t i = 0; i < axis_size; ++i) {
162  double result_angle = incidentAngle(i) + angle_shift;
163  result.emplace_back(wavelength, -result_angle, getZAxis());
164  if (!alpha_limits.isInRange(result_angle))
165  result.back().setCalculationFlag(false); // false = exclude from calculations
166  }
167  return result;
168 }
Direction direction() const
Definition: Beam.h:45
double wavelength() const
Definition: Beam.h:43
double incidentAngle(size_t index) const
double alpha() const
Definition: Direction.h:29

References Direction::alpha(), ISimulation::beam(), Beam::direction(), getAlphaAxis(), getZAxis(), incidentAngle(), result(), IAxis::size(), and Beam::wavelength().

Referenced by initSimulationElementVector().

Here is the call graph for this function:

◆ generateSingleThreadedComputation()

std::unique_ptr< IComputation > DepthProbeSimulation::generateSingleThreadedComputation ( size_t  start,
size_t  n_elements 
)
overrideprivatevirtual

Generate a single threaded computation for a given range of simulation elements.

Parameters
startIndex of the first element to include into computation
n_elementsNumber of elements to process

Implements ISimulation.

Definition at line 171 of file DepthProbeSimulation.cpp.

172 {
173  ASSERT(start < m_sim_elements.size() && start + n_elements <= m_sim_elements.size());
174  const auto& begin = m_sim_elements.begin() + static_cast<long>(start);
175  return std::make_unique<DepthProbeComputation>(*sample(), options(), progress(), begin,
176  begin + static_cast<long>(n_elements));
177 }
#define ASSERT(condition)
Definition: Assert.h:31
const MultiLayer * sample() const
ProgressHandler & progress()
Definition: ISimulation.h:125
const SimulationOptions & options() const
Definition: ISimulation.h:124

References ASSERT, m_sim_elements, ISimulation::options(), ISimulation::progress(), and ISimulation::sample().

Here is the call graph for this function:

◆ getAlphaAxis()

const IAxis * DepthProbeSimulation::getAlphaAxis ( ) const

Returns a pointer to incident angle axis.

Definition at line 71 of file DepthProbeSimulation.cpp.

72 {
73  if (!m_alpha_axis)
74  throw std::runtime_error("Error in DepthProbeSimulation::getAlphaAxis: incident angle axis "
75  "was not initialized.");
76  return m_alpha_axis.get();
77 }

References m_alpha_axis.

Referenced by createIntensityData(), generateSimulationElements(), numberOfSimulationElements(), rawResults(), setRawResults(), and validityCheck().

◆ getChildren()

std::vector< const INode * > ISimulation::getChildren ( ) const
virtualinherited

Returns a vector of children.

Reimplemented from INode.

Definition at line 254 of file ISimulation.cpp.

255 {
256  std::vector<const INode*> result;
257  result.push_back(&instrument());
259  if (m_background)
260  result.push_back(m_background.get());
261  return result;
262 }
SampleProvider m_sample_provider
Definition: ISimulation.h:158
const Instrument & instrument() const
Definition: ISimulation.h:55
virtual SimulationResult result() const =0
Returns the results of the simulation in a format that supports unit conversion and export to numpy a...
std::vector< const INode * > getChildren() const override
Returns a vector of children.

References SampleProvider::getChildren(), ISimulation::instrument(), ISimulation::m_background, ISimulation::m_sample_provider, and ISimulation::result().

Here is the call graph for this function:

◆ getDetector() [1/2]

IDetector* ISimulation::getDetector ( )
inlineinherited

Definition at line 61 of file ISimulation.h.

61 { return m_instrument.getDetector(); }
IDetector * getDetector()
Definition: Instrument.cpp:96

References Instrument::getDetector(), and ISimulation::m_instrument.

Referenced by ISimulation2D::detector2D(), TransformFromDomain::setDetector(), and TransformFromDomain::setDetectorMasks().

Here is the call graph for this function:

◆ getDetector() [2/2]

const IDetector* ISimulation::getDetector ( ) const
inlineinherited

Definition at line 62 of file ISimulation.h.

62 { return m_instrument.getDetector(); }

References Instrument::getDetector(), and ISimulation::m_instrument.

Here is the call graph for this function:

◆ getDistributionHandler()

const DistributionHandler& ISimulation::getDistributionHandler ( ) const
inlineinherited

Definition at line 88 of file ISimulation.h.

88 { return m_distribution_handler; }

References ISimulation::m_distribution_handler.

Referenced by TransformFromDomain::setGISASBeamItem().

◆ getName()

◆ getOptions() [1/2]

SimulationOptions& ISimulation::getOptions ( )
inlineinherited

Definition at line 92 of file ISimulation.h.

92 { return m_options; }
SimulationOptions m_options
Definition: ISimulation.h:156

References ISimulation::m_options.

◆ getOptions() [2/2]

◆ getZAxis()

const IAxis * DepthProbeSimulation::getZAxis ( ) const

Returns a pointer to z-position axis.

Definition at line 79 of file DepthProbeSimulation.cpp.

80 {
81  if (!m_z_axis)
82  throw std::runtime_error("Error in DepthProbeSimulation::getZAxis: position axis "
83  "was not initialized.");
84  return m_z_axis.get();
85 }

References m_z_axis.

Referenced by createIntensityData(), generateSimulationElements(), initSimulationElementVector(), rawResults(), and setRawResults().

◆ incidentAngle()

double DepthProbeSimulation::incidentAngle ( size_t  index) const
private

Definition at line 262 of file DepthProbeSimulation.cpp.

263 {
264  return m_alpha_axis->bin(index).center();
265 }

References m_alpha_axis.

Referenced by generateSimulationElements().

◆ initialize()

void DepthProbeSimulation::initialize ( )
private

Initializes simulation.

Definition at line 215 of file DepthProbeSimulation.cpp.

216 {
217  setName("DepthProbeSimulation");
218 
219  // allow for negative inclinations in the beam of specular simulation
220  // it is required for proper averaging in the case of divergent beam
221  auto inclination = beam().parameter("InclinationAngle");
222  inclination->setLimits(RealLimits::limited(-M_PI_2, M_PI_2));
223 }
#define M_PI_2
Definition: Constants.h:45
RealParameter * parameter(const std::string &name) const
Returns parameter with given 'name'.
void setName(const std::string &name)
static RealLimits limited(double left_bound_value, double right_bound_value)
Creates an object bounded from the left and right.
Definition: RealLimits.cpp:125
RealParameter & setLimits(const RealLimits &limits)

References ISimulation::beam(), RealLimits::limited(), M_PI_2, IParametricComponent::parameter(), RealParameter::setLimits(), and IParametricComponent::setName().

Referenced by DepthProbeSimulation().

Here is the call graph for this function:

◆ initSimulationElementVector()

void DepthProbeSimulation::initSimulationElementVector ( )
overrideprivatevirtual

Initializes the vector of ISimulation elements.

Implements ISimulation.

Definition at line 143 of file DepthProbeSimulation.cpp.

144 {
146 
147  if (!m_cache.empty())
148  return;
149  m_cache.resize(m_sim_elements.size(), std::valarray<double>(0.0, getZAxis()->size()));
150 }
std::vector< DepthProbeElement > generateSimulationElements(const Beam &beam)
Generate simulation elements for given beam.

References ISimulation::beam(), generateSimulationElements(), getZAxis(), m_cache, and m_sim_elements.

Here is the call graph for this function:

◆ instrument() [1/2]

Instrument& ISimulation::instrument ( )
inlineinherited

Definition at line 56 of file ISimulation.h.

56 { return m_instrument; }

References ISimulation::m_instrument.

◆ instrument() [2/2]

◆ intensityMapSize()

size_t DepthProbeSimulation::intensityMapSize ( ) const
overridevirtual

Returns the total number of the intensity values in the simulation result.

Implements ISimulation.

Definition at line 87 of file DepthProbeSimulation.cpp.

88 {
89  if (!m_z_axis || !m_alpha_axis)
90  throw std::runtime_error("Error in DepthProbeSimulation::intensityMapSize: attempt to "
91  "access non-initialized data.");
92  return m_z_axis->size() * m_alpha_axis->size();
93 }

References m_alpha_axis, and m_z_axis.

◆ moveDataFromCache()

void DepthProbeSimulation::moveDataFromCache ( )
overrideprivatevirtual

Implements ISimulation.

Definition at line 253 of file DepthProbeSimulation.cpp.

254 {
255  checkCache();
256  for (size_t i = 0, size = m_sim_elements.size(); i < size; ++i)
257  m_sim_elements[i].setIntensities(std::move(m_cache[i]));
258  m_cache.clear();
259  m_cache.shrink_to_fit();
260 }

References checkCache(), m_cache, and m_sim_elements.

Here is the call graph for this function:

◆ normalize()

void DepthProbeSimulation::normalize ( size_t  start_ind,
size_t  n_elements 
)
overrideprivatevirtual

Normalize the detector counts to beam intensity, to solid angle, and to exposure angle.

Parameters
start_indIndex of the first element to operate on
n_elementsNumber of elements to process

Implements ISimulation.

Definition at line 225 of file DepthProbeSimulation.cpp.

226 {
227  const double beam_intensity = beam().intensity();
228  for (size_t i = start_ind, stop_point = start_ind + n_elements; i < stop_point; ++i) {
229  auto& element = m_sim_elements[i];
230  const double alpha_i = -element.getAlphaI();
231  const auto footprint = beam().footprintFactor();
232  double intensity_factor = beam_intensity;
233  if (footprint != nullptr)
234  intensity_factor = intensity_factor * footprint->calculate(alpha_i);
235  element.setIntensities(element.getIntensities() * intensity_factor);
236  }
237 }
double intensity() const
Returns the beam intensity in neutrons/sec.
Definition: Beam.h:42
const IFootprintFactor * footprintFactor() const
Returns footprint factor.
Definition: Beam.cpp:101
virtual double calculate(double alpha) const =0
Calculate footprint correction coefficient from the beam incident angle alpha.

References ISimulation::beam(), IFootprintFactor::calculate(), Beam::footprintFactor(), Beam::intensity(), and m_sim_elements.

Here is the call graph for this function:

◆ numberOfSimulationElements()

size_t DepthProbeSimulation::numberOfSimulationElements ( ) const
overrideprivatevirtual

Gets the number of elements this simulation needs to calculate.

Implements ISimulation.

Definition at line 44 of file DepthProbeSimulation.cpp.

45 {
46  return getAlphaAxis()->size();
47 }

References getAlphaAxis(), and IAxis::size().

Here is the call graph for this function:

◆ onChange()

◆ options()

◆ parameter()

RealParameter * IParametricComponent::parameter ( const std::string &  name) const
inherited

Returns parameter with given 'name'.

Definition at line 86 of file IParametricComponent.cpp.

87 {
88  return m_pool->parameter(name);
89 }
std::unique_ptr< ParameterPool > m_pool
parameter pool (kind of pointer-to-implementation)
QString const & name(EShape k)
Definition: particles.cpp:21

References IParametricComponent::m_pool, and RealSpace::Particles::name().

Referenced by initialize(), SpecularSimulation::initialize(), Lattice3D::initialize(), IParticle::registerAbundance(), ParticleLayout::registerParticleDensity(), IParticle::registerPosition(), Layer::registerThickness(), Lattice2D::setRotationEnabled(), and DistributionLogNormal::setUnits().

Here is the call graph for this function:

◆ parameterPool()

ParameterPool* IParametricComponent::parameterPool ( ) const
inlineinherited

Returns pointer to the parameter pool.

Definition at line 39 of file IParametricComponent.h.

39 { return m_pool.get(); } // has non-const usages!

References IParametricComponent::m_pool.

Referenced by INode::INode(), IParametricComponent::IParametricComponent(), pyfmt2::argumentList(), SampleBuilderNode::borrow_builder_parameters(), SampleBuilderNode::reset(), and IDistribution1D::setUnits().

◆ parametersToString()

std::string IParametricComponent::parametersToString ( ) const
inherited

Returns multiline string representing available parameters.

Definition at line 43 of file IParametricComponent.cpp.

44 {
45  std::ostringstream result;
46  std::unique_ptr<ParameterPool> P_pool(createParameterTree());
47  result << *P_pool << "\n";
48  return result.str();
49 }
virtual ParameterPool * createParameterTree() const
Creates new parameter pool, with all local parameters and those of its children.

References IParametricComponent::createParameterTree().

Here is the call graph for this function:

◆ parent() [1/2]

INode * INode::parent ( )
inherited

Definition at line 89 of file INode.cpp.

90 {
91  return const_cast<INode*>(m_parent);
92 }

References INode::m_parent.

◆ parent() [2/2]

◆ prepareSimulation()

void ISimulation::prepareSimulation ( )
virtualinherited

Put into a clean state for running a simulation.

Reimplemented in SpecularSimulation, OffSpecularSimulation, ISimulation2D, and GISASSimulation.

Definition at line 182 of file ISimulation.cpp.

183 {
186  throw std::runtime_error(
187  "Error in ISimulation::prepareSimulation(): non-default materials of"
188  " several different types are used in the sample provided");
189  gsl_set_error_handler_off();
190 }
void updateSample()
Generates new sample if sample builder defined.
const MultiLayer * sample() const
Returns current sample.
bool ContainsCompatibleMaterials(const MultiLayer &multilayer)
Returns true if the multilayer contains non-default materials of one type only.

References MultiLayerUtils::ContainsCompatibleMaterials(), ISimulation::m_sample_provider, SampleProvider::sample(), and SampleProvider::updateSample().

Referenced by ISimulation2D::prepareSimulation(), SpecularSimulation::prepareSimulation(), and ISimulation::runSimulation().

Here is the call graph for this function:

◆ progeny()

std::vector< const INode * > INode::progeny ( ) const
inherited

Returns a vector of all descendants.

Definition at line 68 of file INode.cpp.

69 {
70  std::vector<const INode*> result;
71  result.push_back(this);
72  for (const auto* child : getChildren()) {
73  for (const auto* p : child->progeny())
74  result.push_back(p);
75  }
76  return result;
77 }

References INode::getChildren().

Referenced by INode::createParameterTree(), and ParticleDistribution::generateParticles().

Here is the call graph for this function:

◆ progress()

ProgressHandler& ISimulation::progress ( )
inlineprotectedinherited

◆ rawResults()

std::vector< double > DepthProbeSimulation::rawResults ( ) const
overrideprivatevirtual

Implements ISimulation.

Definition at line 284 of file DepthProbeSimulation.cpp.

285 {
286  validityCheck();
287  const size_t z_size = getZAxis()->size();
288  const size_t alpha_size = getAlphaAxis()->size();
289 
290  std::vector<double> result;
291  result.reserve(alpha_size * z_size);
292  for (size_t i = 0; i < alpha_size; ++i) {
293  if (m_sim_elements[i].size() != z_size)
294  throw std::runtime_error("Error in DepthProbeSimulation::rawResults: simulation "
295  "element size is not equal to the size of the position axis");
296  const auto& intensities = m_sim_elements[i].getIntensities();
297  result.insert(result.end(), std::begin(intensities), std::end(intensities));
298  }
299 
300  return result;
301 }
void validityCheck() const
Checks if simulation data is ready for retrieval.

References getAlphaAxis(), getZAxis(), m_sim_elements, result(), IAxis::size(), and validityCheck().

Here is the call graph for this function:

◆ registerChild()

void INode::registerChild ( INode node)
inherited

Definition at line 57 of file INode.cpp.

58 {
59  ASSERT(node);
60  node->setParent(this);
61 }
virtual void setParent(const INode *newParent)
Definition: INode.cpp:79

References ASSERT, and INode::setParent().

Referenced by Beam::Beam(), Crystal::Crystal(), IDetector::IDetector(), Instrument::Instrument(), InterferenceFunction2DLattice::InterferenceFunction2DLattice(), InterferenceFunction2DParaCrystal::InterferenceFunction2DParaCrystal(), InterferenceFunction2DSuperLattice::InterferenceFunction2DSuperLattice(), InterferenceFunctionFinite2DLattice::InterferenceFunctionFinite2DLattice(), Particle::Particle(), ParticleDistribution::ParticleDistribution(), ParticleLayout::addAndRegisterAbstractParticle(), ParticleCoreShell::addAndRegisterCore(), MultiLayer::addAndRegisterInterface(), MultiLayer::addAndRegisterLayer(), ParticleCoreShell::addAndRegisterShell(), Layer::addLayout(), ParticleComposition::addParticlePointer(), ISimulation::initialize(), MesoCrystal::initialize(), Beam::operator=(), Instrument::operator=(), IParticle::rotate(), ParticleLayout::setAndRegisterInterferenceFunction(), ISimulation::setBackground(), InterferenceFunction1DLattice::setDecayFunction(), InterferenceFunction2DLattice::setDecayFunction(), Instrument::setDetector(), IDetector::setDetectorResolution(), Beam::setFootprintFactor(), Particle::setFormFactor(), InterferenceFunctionFinite3DLattice::setLattice(), InterferenceFunctionRadialParaCrystal::setProbabilityDistribution(), InterferenceFunction2DParaCrystal::setProbabilityDistributions(), ConvolutionDetectorResolution::setResolutionFunction(), IParticle::setRotation(), LayerInterface::setRoughness(), and InterferenceFunction2DSuperLattice::setSubstructureIFF().

Here is the call graph for this function:

◆ registerParameter()

RealParameter & IParametricComponent::registerParameter ( const std::string &  name,
double *  parpointer 
)
inherited

Definition at line 51 of file IParametricComponent.cpp.

52 {
53  return m_pool->addParameter(
54  new RealParameter(name, data, getName(), [&]() -> void { onChange(); }));
55 }
virtual void onChange()
Action to be taken in inherited class when a parameter has changed.
Wraps a parameter of type double.
Definition: RealParameter.h:31

References IParametricComponent::getName(), IParametricComponent::m_pool, RealSpace::Particles::name(), and IParametricComponent::onChange().

Referenced by BasicLattice2D::BasicLattice2D(), Beam::Beam(), CylindersInBABuilder::CylindersInBABuilder(), DetectionProperties::DetectionProperties(), HexagonalLattice2D::HexagonalLattice2D(), IInterferenceFunction::IInterferenceFunction(), INode::INode(), InterferenceFunction1DLattice::InterferenceFunction1DLattice(), InterferenceFunction2DParaCrystal::InterferenceFunction2DParaCrystal(), InterferenceFunctionHardDisk::InterferenceFunctionHardDisk(), InterferenceFunctionRadialParaCrystal::InterferenceFunctionRadialParaCrystal(), InterferenceFunctionTwin::InterferenceFunctionTwin(), Lattice2D::Lattice2D(), LayerRoughness::LayerRoughness(), MultiLayer::MultiLayer(), ParticleDistribution::ParticleDistribution(), PlainMultiLayerBySLDBuilder::PlainMultiLayerBySLDBuilder(), ResolutionFunction2DGaussian::ResolutionFunction2DGaussian(), ResonatorBuilder::ResonatorBuilder(), SquareLattice2D::SquareLattice2D(), TriangularRippleBuilder::TriangularRippleBuilder(), IParticle::registerAbundance(), ParticleLayout::registerParticleDensity(), Layer::registerThickness(), IParametricComponent::registerVector(), ParticleLayout::registerWeight(), and Lattice2D::setRotationEnabled().

Here is the call graph for this function:

◆ registerVector()

void IParametricComponent::registerVector ( const std::string &  base_name,
kvector_t p_vec,
const std::string &  units = "nm" 
)
inherited

Definition at line 57 of file IParametricComponent.cpp.

59 {
60  registerParameter(XComponentName(base_name), &((*p_vec)[0])).setUnit(units);
61  registerParameter(YComponentName(base_name), &((*p_vec)[1])).setUnit(units);
62  registerParameter(ZComponentName(base_name), &((*p_vec)[2])).setUnit(units);
63 }
static std::string XComponentName(const std::string &base_name)
static std::string ZComponentName(const std::string &base_name)
RealParameter & registerParameter(const std::string &name, double *parpointer)
static std::string YComponentName(const std::string &base_name)
RealParameter & setUnit(const std::string &name)
MVVM_MODEL_EXPORT std::string base_name(const std::string &path)
Provide the filename of a file path.
Definition: fileutils.cpp:78

References ModelView::Utils::base_name(), IParametricComponent::registerParameter(), RealParameter::setUnit(), IParametricComponent::XComponentName(), IParametricComponent::YComponentName(), and IParametricComponent::ZComponentName().

Referenced by Beam::Beam(), DetectionProperties::DetectionProperties(), InterferenceFunctionTwin::InterferenceFunctionTwin(), MultiLayer::MultiLayer(), Lattice3D::initialize(), and IParticle::registerPosition().

Here is the call graph for this function:

◆ removeParameter()

void IParametricComponent::removeParameter ( const std::string &  name)
inherited

Definition at line 91 of file IParametricComponent.cpp.

92 {
93  m_pool->removeParameter(name);
94 }

References IParametricComponent::m_pool, and RealSpace::Particles::name().

Referenced by IParticle::registerAbundance(), ParticleLayout::registerParticleDensity(), Layer::registerThickness(), IParametricComponent::removeVector(), and Lattice2D::setRotationEnabled().

Here is the call graph for this function:

◆ removeVector()

void IParametricComponent::removeVector ( const std::string &  base_name)
inherited

Definition at line 96 of file IParametricComponent.cpp.

References ModelView::Utils::base_name(), IParametricComponent::removeParameter(), IParametricComponent::XComponentName(), IParametricComponent::YComponentName(), and IParametricComponent::ZComponentName().

Referenced by IParticle::registerPosition().

Here is the call graph for this function:

◆ result()

SimulationResult DepthProbeSimulation::result ( ) const
overridevirtual

Returns the results of the simulation in a format that supports unit conversion and export to numpy arrays.

Implements ISimulation.

Definition at line 49 of file DepthProbeSimulation.cpp.

50 {
51  validityCheck();
52  auto data = createIntensityData();
53  return SimulationResult(*data, *createUnitConverter());
54 }
std::unique_ptr< IUnitConverter > createUnitConverter() const
std::unique_ptr< OutputData< double > > createIntensityData() const
Creates intensity data from simulation elements.

References createIntensityData(), createUnitConverter(), and validityCheck().

Referenced by createIntensityData(), generateSimulationElements(), and rawResults().

Here is the call graph for this function:

◆ runMPISimulation()

void ISimulation::runMPISimulation ( )
inherited

Run a simulation in a MPI environment.

Definition at line 226 of file ISimulation.cpp.

227 {
228  MPISimulation ompi;
229  ompi.runSimulation(this);
230 }
void runSimulation(ISimulation *simulation)

References MPISimulation::runSimulation().

Here is the call graph for this function:

◆ runSimulation()

void ISimulation::runSimulation ( )
inherited

Run a simulation, possibly averaged over parameter distributions.

Run simulation with possible averaging over parameter distributions.

Definition at line 193 of file ISimulation.cpp.

194 {
196 
197  const size_t total_size = numberOfSimulationElements();
198  size_t param_combinations = m_distribution_handler.getTotalNumberOfSamples();
199 
200  m_progress.reset();
201  m_progress.setExpectedNTicks(param_combinations * total_size);
202 
203  // restrict calculation to current batch
204  const size_t n_batches = m_options.getNumberOfBatches();
205  const size_t current_batch = m_options.getCurrentBatch();
206 
207  const size_t batch_start = getStartIndex(n_batches, current_batch, total_size);
208  const size_t batch_size = getNumberOfElements(n_batches, current_batch, total_size);
209  if (batch_size == 0)
210  return;
211 
212  if(param_combinations == 1)
213  runSingleSimulation(batch_start, batch_size, 1.);
214  else{
215  std::unique_ptr<ParameterPool> param_pool(createParameterTree());
216  for (size_t index = 0; index < param_combinations; ++index) {
217  double weight = m_distribution_handler.setParameterValues(param_pool.get(), index);
218  runSingleSimulation(batch_start, batch_size, weight);
219  }
221  }
224 }
double setParameterValues(ParameterPool *p_parameter_pool, size_t index)
set the parameter values of the simulation object to a specific combination of values,...
void setParameterToMeans(ParameterPool *p_parameter_pool) const
Sets mean distribution values to the parameter pool.
size_t getTotalNumberOfSamples() const
get the total number of parameter value combinations (product of the individual sizes of each paramet...
ParameterPool * createParameterTree() const
Creates new parameter pool, with all local parameters and those of its children.
Definition: INode.cpp:126
void runSingleSimulation(size_t batch_start, size_t batch_size, double weight=1.0)
Runs a single simulation with fixed parameter values.
virtual size_t numberOfSimulationElements() const =0
Gets the number of elements this simulation needs to calculate.
virtual void moveDataFromCache()=0
virtual void prepareSimulation()
Put into a clean state for running a simulation.
virtual void transferResultsToIntensityMap()
Creates the appropriate data structure (e.g.
Definition: ISimulation.h:114
void setExpectedNTicks(size_t n)
unsigned getNumberOfBatches() const
unsigned getCurrentBatch() const

References INode::createParameterTree(), SimulationOptions::getCurrentBatch(), SimulationOptions::getNumberOfBatches(), DistributionHandler::getTotalNumberOfSamples(), ISimulation::m_distribution_handler, ISimulation::m_options, ISimulation::m_progress, ISimulation::moveDataFromCache(), ISimulation::numberOfSimulationElements(), ISimulation::prepareSimulation(), ProgressHandler::reset(), ISimulation::runSingleSimulation(), ProgressHandler::setExpectedNTicks(), DistributionHandler::setParameterToMeans(), DistributionHandler::setParameterValues(), and ISimulation::transferResultsToIntensityMap().

Referenced by JobWorker::start().

Here is the call graph for this function:

◆ runSingleSimulation()

void ISimulation::runSingleSimulation ( size_t  batch_start,
size_t  batch_size,
double  weight = 1.0 
)
privateinherited

Runs a single simulation with fixed parameter values.

If desired, the simulation is run in several threads.

Definition at line 280 of file ISimulation.cpp.

281 {
283 
284  const size_t n_threads = m_options.getNumberOfThreads();
285  ASSERT(n_threads > 0);
286 
287  std::vector<std::unique_ptr<IComputation>> computations;
288 
289  for (size_t i_thread = 0; i_thread < n_threads;
290  ++i_thread) { // Distribute computations by threads
291  const size_t thread_start = batch_start + getStartIndex(n_threads, i_thread, batch_size);
292  const size_t thread_size = getNumberOfElements(n_threads, i_thread, batch_size);
293  if (thread_size == 0)
294  break;
295  computations.push_back(generateSingleThreadedComputation(thread_start, thread_size));
296  }
297  runComputations(computations);
298 
299  normalize(batch_start, batch_size);
300  addBackgroundIntensity(batch_start, batch_size);
301  addDataToCache(weight);
302 }
virtual void addDataToCache(double weight)=0
virtual std::unique_ptr< IComputation > generateSingleThreadedComputation(size_t start, size_t n_elements)=0
Generate a single threaded computation for a given range of simulation elements.
virtual void addBackgroundIntensity(size_t start_ind, size_t n_elements)=0
virtual void initSimulationElementVector()=0
Initializes the vector of ISimulation elements.
virtual void normalize(size_t start_ind, size_t n_elements)=0
Normalize the detector counts to beam intensity, to solid angle, and to exposure angle.
unsigned getNumberOfThreads() const

References ISimulation::addBackgroundIntensity(), ISimulation::addDataToCache(), ASSERT, ISimulation::generateSingleThreadedComputation(), SimulationOptions::getNumberOfThreads(), ISimulation::initSimulationElementVector(), ISimulation::m_options, and ISimulation::normalize().

Referenced by ISimulation::runSimulation().

Here is the call graph for this function:

◆ sample()

const MultiLayer * ISimulation::sample ( ) const
inherited

Definition at line 238 of file ISimulation.cpp.

239 {
240  return m_sample_provider.sample();
241 }

References ISimulation::m_sample_provider, and SampleProvider::sample().

Referenced by ISimulation::ISimulation(), generateSingleThreadedComputation(), ISimulation2D::generateSingleThreadedComputation(), SpecularSimulation::generateSingleThreadedComputation(), ISimulation::setSample(), and validityCheck().

Here is the call graph for this function:

◆ setBackground()

void ISimulation::setBackground ( const IBackground bg)
inherited

Definition at line 248 of file ISimulation.cpp.

249 {
250  m_background.reset(bg.clone());
252 }
virtual IBackground * clone() const =0
void registerChild(INode *node)
Definition: INode.cpp:57

References IBackground::clone(), ISimulation::m_background, and INode::registerChild().

Referenced by ISimulation::ISimulation(), and StandardSimulations::ConstantBackgroundGISAS().

Here is the call graph for this function:

◆ setBeamParameters() [1/2]

void DepthProbeSimulation::setBeamParameters ( double  lambda,
const IAxis alpha_axis,
const IFootprintFactor beam_shape 
)
private

Sets beam parameters with alpha_i of the beam defined in the range.

Definition at line 112 of file DepthProbeSimulation.cpp.

114 {
115  if (lambda <= 0.0)
116  throw std::runtime_error(
117  "Error in DepthProbeSimulation::setBeamParameters: wavelength must be positive.");
118  if (alpha_axis.lowerBound() < 0.0)
119  throw std::runtime_error(
120  "Error in DepthProbeSimulation::setBeamParameters: minimum value on "
121  "angle axis is negative.");
122  if (alpha_axis.lowerBound() >= alpha_axis.upperBound())
123  throw std::runtime_error(
124  "Error in DepthProbeSimulation::setBeamParameters: maximal value on "
125  "angle axis is less or equal to the minimal one.");
126  if (alpha_axis.size() == 0)
127  throw std::runtime_error(
128  "Error in DepthProbeSimulation::setBeamParameters: angle axis is empty");
129 
130  SpecularDetector1D detector(alpha_axis);
132  m_alpha_axis.reset(alpha_axis.clone());
133 
134  // beam is initialized with zero-valued angles
135  // Zero-valued incident alpha is required for proper
136  // taking into account beam resolution effects
137  instrument().setBeamParameters(lambda, zero_alpha_i, zero_phi_i);
138 
139  if (beam_shape)
140  beam().setFootprintFactor(*beam_shape);
141 }
void setFootprintFactor(const IFootprintFactor &shape_factor)
Sets footprint factor to the beam.
Definition: Beam.cpp:106
virtual IAxis * clone() const =0
clone function
virtual double upperBound() const =0
Returns value of last point of axis.
virtual double lowerBound() const =0
Returns value of first point of axis.
void setBeamParameters(double wavelength, double alpha_i, double phi_i)
Sets the beam wavelength and incoming angles.
Definition: Instrument.cpp:76
void setDetector(const IDetector &detector)
Sets the detector (axes can be overwritten later)
Definition: Instrument.cpp:52
1D detector for specular simulations.

References ISimulation::beam(), IAxis::clone(), ISimulation::detector(), ISimulation::instrument(), IAxis::lowerBound(), m_alpha_axis, Instrument::setBeamParameters(), Instrument::setDetector(), Beam::setFootprintFactor(), IAxis::size(), and IAxis::upperBound().

Here is the call graph for this function:

◆ setBeamParameters() [2/2]

void DepthProbeSimulation::setBeamParameters ( double  lambda,
int  nbins,
double  alpha_i_min,
double  alpha_i_max,
const IFootprintFactor beam_shape = nullptr 
)

Sets beam parameters with alpha_i of the beam defined in the range.

Definition at line 56 of file DepthProbeSimulation.cpp.

58 {
59  FixedBinAxis axis("alpha_i", static_cast<size_t>(nbins), alpha_i_min, alpha_i_max);
60  setBeamParameters(lambda, axis, beam_shape);
61 }
void setBeamParameters(double lambda, int nbins, double alpha_i_min, double alpha_i_max, const IFootprintFactor *beam_shape=nullptr)
Sets beam parameters with alpha_i of the beam defined in the range.
Axis with fixed bin size.
Definition: FixedBinAxis.h:23

◆ setDetectorResolutionFunction()

void ISimulation::setDetectorResolutionFunction ( const IResolutionFunction2D resolution_function)
inherited

Definition at line 177 of file ISimulation.cpp.

178 {
179  detector().setResolutionFunction(resolution_function);
180 }
void setResolutionFunction(const IResolutionFunction2D &resFunc)
Definition: IDetector.cpp:112

References ISimulation::detector(), and IDetector::setResolutionFunction().

Referenced by StandardSimulations::MiniGISASDetectorResolution().

Here is the call graph for this function:

◆ setInstrument()

void ISimulation::setInstrument ( const Instrument instrument_)
inherited

Definition at line 158 of file ISimulation.cpp.

159 {
160  m_instrument = instrument_;
162 }
virtual void updateIntensityMap()
Definition: ISimulation.h:119

References ISimulation::m_instrument, and ISimulation::updateIntensityMap().

Here is the call graph for this function:

◆ setName()

void IParametricComponent::setName ( const std::string &  name)
inlineinherited

Definition at line 69 of file IParametricComponent.h.

69 { m_name = name; }

References IParametricComponent::m_name, and RealSpace::Particles::name().

Referenced by BasicLattice2D::BasicLattice2D(), Beam::Beam(), ConvolutionDetectorResolution::ConvolutionDetectorResolution(), Crystal::Crystal(), DetectionProperties::DetectionProperties(), DistributionHandler::DistributionHandler(), FormFactorCoreShell::FormFactorCoreShell(), FormFactorCrystal::FormFactorCrystal(), FormFactorDecoratorMaterial::FormFactorDecoratorMaterial(), FormFactorDecoratorPositionFactor::FormFactorDecoratorPositionFactor(), FormFactorDecoratorRotation::FormFactorDecoratorRotation(), FormFactorWeighted::FormFactorWeighted(), HexagonalLattice2D::HexagonalLattice2D(), IDetector::IDetector(), INode::INode(), Instrument::Instrument(), InterferenceFunction1DLattice::InterferenceFunction1DLattice(), InterferenceFunction2DLattice::InterferenceFunction2DLattice(), InterferenceFunction2DParaCrystal::InterferenceFunction2DParaCrystal(), InterferenceFunction2DSuperLattice::InterferenceFunction2DSuperLattice(), InterferenceFunction3DLattice::InterferenceFunction3DLattice(), InterferenceFunctionFinite2DLattice::InterferenceFunctionFinite2DLattice(), InterferenceFunctionFinite3DLattice::InterferenceFunctionFinite3DLattice(), InterferenceFunctionHardDisk::InterferenceFunctionHardDisk(), InterferenceFunctionNone::InterferenceFunctionNone(), InterferenceFunctionRadialParaCrystal::InterferenceFunctionRadialParaCrystal(), InterferenceFunctionTwin::InterferenceFunctionTwin(), ISampleBuilder::ISampleBuilder(), IsGISAXSDetector::IsGISAXSDetector(), Lattice3D::Lattice3D(), Layer::Layer(), LayerInterface::LayerInterface(), LayerRoughness::LayerRoughness(), MultiLayer::MultiLayer(), ParticleCoreShell::ParticleCoreShell(), ParticleDistribution::ParticleDistribution(), ParticleLayout::ParticleLayout(), RectangularDetector::RectangularDetector(), ResolutionFunction2DGaussian::ResolutionFunction2DGaussian(), SampleBuilderNode::SampleBuilderNode(), SphericalDetector::SphericalDetector(), SquareLattice2D::SquareLattice2D(), Layer::clone(), LayersWithAbsorptionBuilder::createSampleByIndex(), Basic2DParaCrystalBuilder::createSampleByIndex(), ParticleInVacuumBuilder::createSampleByIndex(), SimpleMagneticRotationBuilder::createSampleByIndex(), initialize(), GISASSimulation::initialize(), OffSpecularSimulation::initialize(), SpecularSimulation::initialize(), SpecularDetector1D::initialize(), MesoCrystal::initialize(), Particle::initialize(), ParticleComposition::initialize(), Beam::operator=(), SampleBuilderNode::operator=(), SampleBuilderNode::reset(), and SampleBuilderNode::setSBN().

Here is the call graph for this function:

◆ setOptions()

void ISimulation::setOptions ( const SimulationOptions options)
inlineinherited

Definition at line 90 of file ISimulation.h.

90 { m_options = options; }

References ISimulation::m_options, and ISimulation::options().

Here is the call graph for this function:

◆ setParameterValue()

void IParametricComponent::setParameterValue ( const std::string &  name,
double  value 
)
inherited

Definition at line 65 of file IParametricComponent.cpp.

66 {
67  if (name.find('*') == std::string::npos && name.find('/') == std::string::npos) {
68  m_pool->setParameterValue(name, value);
69  } else {
70  std::unique_ptr<ParameterPool> P_pool{createParameterTree()};
71  if (name.find('*') != std::string::npos)
72  P_pool->setMatchedParametersValue(name, value);
73  else
74  P_pool->setParameterValue(name, value);
75  }
76 }
int setMatchedParametersValue(const std::string &wildcards, double value)
Sets value of the nonzero parameters that match pattern ('*' allowed), or throws.

References IParametricComponent::createParameterTree(), IParametricComponent::m_pool, RealSpace::Particles::name(), and ParameterPool::setMatchedParametersValue().

Referenced by AsymRippleBuilder::buildSample(), and IParametricComponent::setVectorValue().

Here is the call graph for this function:

◆ setParent()

void INode::setParent ( const INode newParent)
virtualinherited

Reimplemented in SampleProvider.

Definition at line 79 of file INode.cpp.

80 {
81  m_parent = newParent;
82 }

References INode::m_parent.

Referenced by INode::registerChild(), SampleProvider::setBuilder(), and SampleProvider::setParent().

◆ setRawResults()

void DepthProbeSimulation::setRawResults ( const std::vector< double > &  raw_data)
overrideprivatevirtual

Implements ISimulation.

Definition at line 303 of file DepthProbeSimulation.cpp.

304 {
305  validityCheck();
306  const size_t z_size = getZAxis()->size();
307  const size_t alpha_size = getAlphaAxis()->size();
308 
309  if (raw_results.size() != z_size * alpha_size)
310  throw std::runtime_error(
311  "Error in DepthProbeSimulation::setRawResults: the vector to set is of invalid size");
312 
313  const double* raw_array = raw_results.data();
314  for (size_t i = 0; i < alpha_size; ++i) {
315  std::valarray<double> fixed_angle_result(raw_array, z_size);
316  m_sim_elements[i].setIntensities(std::move(fixed_angle_result));
317  raw_array += z_size;
318  }
319 }

References getAlphaAxis(), getZAxis(), m_sim_elements, IAxis::size(), and validityCheck().

Here is the call graph for this function:

◆ setSample()

void ISimulation::setSample ( const MultiLayer sample)
inherited

The MultiLayer object will not be owned by the ISimulation object.

Definition at line 233 of file ISimulation.cpp.

234 {
236 }
void setSample(const MultiLayer &multilayer)

References ISimulation::m_sample_provider, ISimulation::sample(), and SampleProvider::setSample().

Referenced by ISimulation::ISimulation().

Here is the call graph for this function:

◆ setSampleBuilder()

void ISimulation::setSampleBuilder ( const std::shared_ptr< ISampleBuilder > &  sample_builder)
inherited

Definition at line 243 of file ISimulation.cpp.

244 {
245  m_sample_provider.setBuilder(sample_builder);
246 }
void setBuilder(const std::shared_ptr< ISampleBuilder > &sample_builder)

References ISimulation::m_sample_provider, and SampleProvider::setBuilder().

Here is the call graph for this function:

◆ setTerminalProgressMonitor()

void ISimulation::setTerminalProgressMonitor ( )
inherited

Initializes a progress monitor that prints to stdout.

Definition at line 166 of file ISimulation.cpp.

167 {
168  m_progress.subscribe([](size_t percentage_done) -> bool {
169  if (percentage_done < 100)
170  std::cout << std::setprecision(2) << "\r... " << percentage_done << "%" << std::flush;
171  else // wipe out
172  std::cout << "\r... 100%\n";
173  return true;
174  });
175 }
void subscribe(ProgressHandler::Callback_t callback)

References ISimulation::m_progress, and ProgressHandler::subscribe().

Here is the call graph for this function:

◆ setVectorValue()

void IParametricComponent::setVectorValue ( const std::string &  base_name,
kvector_t  value 
)
inherited

Definition at line 78 of file IParametricComponent.cpp.

79 {
83 }
T z() const
Returns z-component in cartesian coordinate system.
Definition: BasicVector3D.h:67
T y() const
Returns y-component in cartesian coordinate system.
Definition: BasicVector3D.h:65
T x() const
Returns x-component in cartesian coordinate system.
Definition: BasicVector3D.h:63
void setParameterValue(const std::string &name, double value)

References ModelView::Utils::base_name(), IParametricComponent::setParameterValue(), BasicVector3D< T >::x(), IParametricComponent::XComponentName(), BasicVector3D< T >::y(), IParametricComponent::YComponentName(), BasicVector3D< T >::z(), and IParametricComponent::ZComponentName().

Here is the call graph for this function:

◆ setZSpan()

void DepthProbeSimulation::setZSpan ( size_t  n_bins,
double  z_min,
double  z_max 
)

Set z positions for intensity calculations.

Negative z's correspond to the area under sample surface. The more negative z is, the deeper layer corresponds to it.

Definition at line 63 of file DepthProbeSimulation.cpp.

64 {
65  if (z_max <= z_min)
66  throw std::runtime_error("Error in DepthProbeSimulation::setZSpan: maximum on-axis value "
67  "is less or equal to the minimum one");
68  m_z_axis = std::make_unique<FixedBinAxis>("z", n_bins, z_min, z_max);
69 }

References m_z_axis.

◆ subscribe()

void ISimulation::subscribe ( ProgressHandler::Callback_t  inform)
inlineinherited

Definition at line 94 of file ISimulation.h.

94 { m_progress.subscribe(inform); }

References ISimulation::m_progress, and ProgressHandler::subscribe().

Referenced by JobWorker::start().

Here is the call graph for this function:

◆ transferResultsToIntensityMap()

virtual void ISimulation::transferResultsToIntensityMap ( )
inlineprotectedvirtualinherited

Creates the appropriate data structure (e.g.

2D intensity map) from the calculated SimulationElement objects

Reimplemented in OffSpecularSimulation.

Definition at line 114 of file ISimulation.h.

114 {}

Referenced by ISimulation::runSimulation(), ISimulation2D::setRawResults(), and SpecularSimulation::setRawResults().

◆ transferToCPP()

virtual void ICloneable::transferToCPP ( )
inlinevirtualinherited

Used for Python overriding of clone (see swig/tweaks.py)

Definition at line 34 of file ICloneable.h.

◆ treeToString()

std::string INode::treeToString ( ) const
virtualinherited

Returns multiline string representing tree structure below the node.

Definition at line 52 of file INode.cpp.

53 {
54  return NodeUtils::nodeToString(this);
55 }
std::string nodeToString(const INode *node)
Returns multiline string representing tree structure starting from given node.
Definition: NodeUtils.cpp:81

References NodeUtils::nodeToString().

Here is the call graph for this function:

◆ updateIntensityMap()

virtual void ISimulation::updateIntensityMap ( )
inlineprotectedvirtualinherited

Reimplemented in OffSpecularSimulation.

Definition at line 119 of file ISimulation.h.

119 {}

Referenced by ISimulation2D::setDetectorParameters(), and ISimulation::setInstrument().

◆ validateParametrization()

void DepthProbeSimulation::validateParametrization ( const ParameterDistribution par_distr) const
overrideprivatevirtual

Checks the distribution validity for simulation.

Reimplemented from ISimulation.

Definition at line 200 of file DepthProbeSimulation.cpp.

201 {
202  const bool zero_mean = par_distr.getDistribution()->getMean() == 0.0;
203  if (zero_mean)
204  return;
205 
206  std::unique_ptr<ParameterPool> parameter_pool(createParameterTree());
207  const std::vector<RealParameter*> names =
208  parameter_pool->getMatchedParameters(par_distr.getMainParameterName());
209  for (const auto par : names)
210  if (par->getName().find("InclinationAngle") != std::string::npos && !zero_mean)
211  throw std::runtime_error("Error in DepthProbeSimulation: parameter distribution of "
212  "beam inclination angle should have zero mean.");
213 }
virtual double getMean() const =0
Returns the distribution-specific mean.
const IDistribution1D * getDistribution() const
std::string getMainParameterName() const
get the main parameter's name

References INode::createParameterTree(), ParameterDistribution::getDistribution(), ParameterDistribution::getMainParameterName(), and IDistribution1D::getMean().

Here is the call graph for this function:

◆ validityCheck()

void DepthProbeSimulation::validityCheck ( ) const
private

Checks if simulation data is ready for retrieval.

Definition at line 179 of file DepthProbeSimulation.cpp.

180 {
181  const MultiLayer* current_sample = sample();
182  if (!current_sample)
183  throw std::runtime_error(
184  "Error in DepthProbeSimulation::validityCheck: no sample found in the simulation.");
185 
186  const size_t data_size = m_sim_elements.size();
187  if (data_size != getAlphaAxis()->size())
188  throw std::runtime_error(
189  "Error in DepthProbeSimulation::validityCheck: length of simulation "
190  "element vector is not equal to the number of inclination angles");
191 }
Our sample model: a stack of layers one below the other.
Definition: MultiLayer.h:41

References getAlphaAxis(), m_sim_elements, and ISimulation::sample().

Referenced by rawResults(), result(), and setRawResults().

Here is the call graph for this function:

◆ XComponentName()

std::string IParametricComponent::XComponentName ( const std::string &  base_name)
staticinherited

Definition at line 103 of file IParametricComponent.cpp.

104 {
105  return base_name + "X";
106 }

References ModelView::Utils::base_name().

Referenced by Lattice3D::initialize(), IParticle::registerPosition(), IParametricComponent::registerVector(), IParametricComponent::removeVector(), IParametricComponent::setVectorValue(), and VectorParameterTranslator::translate().

Here is the call graph for this function:

◆ YComponentName()

std::string IParametricComponent::YComponentName ( const std::string &  base_name)
staticinherited

Definition at line 108 of file IParametricComponent.cpp.

109 {
110  return base_name + "Y";
111 }

References ModelView::Utils::base_name().

Referenced by IParametricComponent::registerVector(), IParametricComponent::removeVector(), IParametricComponent::setVectorValue(), and VectorParameterTranslator::translate().

Here is the call graph for this function:

◆ ZComponentName()

std::string IParametricComponent::ZComponentName ( const std::string &  base_name)
staticinherited

Definition at line 113 of file IParametricComponent.cpp.

114 {
115  return base_name + "Z";
116 }

References ModelView::Utils::base_name().

Referenced by IParametricComponent::registerVector(), IParametricComponent::removeVector(), IParametricComponent::setVectorValue(), and VectorParameterTranslator::translate().

Here is the call graph for this function:

Member Data Documentation

◆ m_alpha_axis

std::unique_ptr<IAxis> DepthProbeSimulation::m_alpha_axis
private

◆ m_background

std::unique_ptr<IBackground> ISimulation::m_background
privateinherited

◆ m_cache

std::vector<std::valarray<double> > DepthProbeSimulation::m_cache
private

◆ m_distribution_handler

DistributionHandler ISimulation::m_distribution_handler
privateinherited

◆ m_instrument

◆ m_name

std::string IParametricComponent::m_name
privateinherited

◆ m_NP

const size_t INode::m_NP
protectedinherited

Definition at line 88 of file INode.h.

Referenced by INode::INode().

◆ m_options

◆ m_P

std::vector<double> INode::m_P
protectedinherited

Definition at line 89 of file INode.h.

Referenced by INode::INode(), and IFootprintFactor::setWidthRatio().

◆ m_parent

const INode* INode::m_parent {nullptr}
privateinherited

Definition at line 83 of file INode.h.

Referenced by INode::displayName(), INode::parent(), and INode::setParent().

◆ m_pool

◆ m_progress

ProgressHandler ISimulation::m_progress
privateinherited

◆ m_sample_provider

◆ m_sim_elements

◆ m_z_axis

std::unique_ptr<IAxis> DepthProbeSimulation::m_z_axis
private

The documentation for this class was generated from the following files: