diff --git a/Core/Binning/SimulationElement.cpp b/Core/Binning/SimulationElement.cpp
index 307c4f012f74f335ebaf624633686852b59aa4d7..8aef522613003862b6d46313304e3aed235b9af4 100644
--- a/Core/Binning/SimulationElement.cpp
+++ b/Core/Binning/SimulationElement.cpp
@@ -25,25 +25,28 @@ SimulationElement::SimulationElement(double wavelength, double alpha_i, double p
, m_phi_i(phi_i)
, m_intensity(0.0)
, mP_pixel(std::move(pixel))
- , m_contains_specular(false)
{
initPolarization();
}
SimulationElement::SimulationElement(const SimulationElement& other)
: m_wavelength(other.m_wavelength), m_alpha_i(other.m_alpha_i), m_phi_i(other.m_phi_i)
- , m_intensity(other.m_intensity), m_contains_specular(other.m_contains_specular)
+ , m_intensity(other.m_intensity)
{
mP_pixel.reset(other.mP_pixel->clone());
+ if (other.m_specular_data)
+ m_specular_data.reset(new SpecularData(*other.m_specular_data));
m_polarization = other.m_polarization;
m_analyzer_operator = other.m_analyzer_operator;
}
SimulationElement::SimulationElement(const SimulationElement& other, double x, double y)
: m_wavelength(other.m_wavelength), m_alpha_i(other.m_alpha_i), m_phi_i(other.m_phi_i)
- , m_intensity(other.m_intensity), m_contains_specular(other.m_contains_specular)
+ , m_intensity(other.m_intensity)
{
mP_pixel.reset(other.mP_pixel->createZeroSizePixel(x, y));
+ if (other.m_specular_data)
+ m_specular_data.reset(new SpecularData(*other.m_specular_data));
m_polarization = other.m_polarization;
m_analyzer_operator = other.m_analyzer_operator;
}
@@ -56,7 +59,7 @@ SimulationElement::SimulationElement(SimulationElement&& other) noexcept
, m_polarization(std::move(other.m_polarization))
, m_analyzer_operator(std::move(other.m_analyzer_operator))
, mP_pixel(std::move(other.mP_pixel))
- , m_contains_specular(other.m_contains_specular)
+ , m_specular_data(std::move(other.m_specular_data))
{
}
@@ -108,7 +111,7 @@ void SimulationElement::swapContent(SimulationElement &other)
std::swap(m_polarization, other.m_polarization);
std::swap(m_analyzer_operator, other.m_analyzer_operator);
std::swap(mP_pixel, other.mP_pixel);
- std::swap(m_contains_specular, other.m_contains_specular);
+ std::swap(m_specular_data, other.m_specular_data);
}
void SimulationElement::initPolarization()
@@ -127,14 +130,14 @@ double SimulationElement::getPhi(double x, double y) const
return getKf(x,y).phi();
}
-bool SimulationElement::containsSpecularWavevector() const
+void SimulationElement::setSpecular()
{
- return m_contains_specular;
+ m_specular_data.reset(new SpecularData);
}
-void SimulationElement::setSpecular(bool contains_specular)
+void SimulationElement::setSpecular(std::unique_ptr<SpecularData> specular_data)
{
- m_contains_specular = contains_specular;
+ m_specular_data = std::move(specular_data);
}
double SimulationElement::getIntegrationFactor(double x, double y) const {
@@ -152,3 +155,26 @@ void addElementsWithWeight(std::vector<SimulationElement>::const_iterator first,
for (std::vector<SimulationElement>::const_iterator it = first; it != last; ++it, ++result)
result->addIntensity(it->getIntensity() * weight);
}
+
+SpecularData::SpecularData() : data_type_used(DATA_TYPE::Invalid) {}
+
+SpecularData::SpecularData(MatrixVector coefficients)
+ : data(std::move(coefficients))
+ , data_type_used(DATA_TYPE::Matrix)
+{}
+
+SpecularData::SpecularData(ScalarVector coefficients)
+ : data(std::move(coefficients))
+ , data_type_used(DATA_TYPE::Scalar)
+{}
+
+const ILayerRTCoefficients& SpecularData::operator[](size_t index) const
+{
+ if (data_type_used == DATA_TYPE::Invalid)
+ throw std::runtime_error(
+ "Error in SpecularData::operator[]: attempt to access uninitialized data");
+ if (data_type_used == DATA_TYPE::Scalar)
+ return (*boost::get<ScalarVector>(&data))[index];
+ else
+ return (*boost::get<MatrixVector>(&data))[index];
+}
diff --git a/Core/Binning/SimulationElement.h b/Core/Binning/SimulationElement.h
index 09017df9ad2540da692fa095b9aed4eb660e87a9..a2c965dc4cea2be6707ee3145870affb3ba3a6ae 100644
--- a/Core/Binning/SimulationElement.h
+++ b/Core/Binning/SimulationElement.h
@@ -20,10 +20,14 @@
#include "EigenCore.h"
#include "Vectors3D.h"
#include "IPixel.h"
+#include "MatrixRTCoefficients.h"
+#include "ScalarRTCoefficients.h"
+#include <boost/variant.hpp>
#include <memory>
#include <vector>
class IPixel;
+class SpecularData;
//! Data stucture containing both input and output of a single detector cell.
//! @ingroup simulation
@@ -80,11 +84,12 @@ public:
double getAlpha(double x, double y) const;
double getPhi(double x, double y) const;
- //! check if element contains given wavevector
- bool containsSpecularWavevector() const;
+ //! check if element corresponds to specular peak
+ SpecularData* specularData() const {return m_specular_data.get();}
- //! indicate that this element contains the specular wavevector
- void setSpecular(bool contains_specular);
+ //! Turn on specular data
+ void setSpecular();
+ void setSpecular(std::unique_ptr<SpecularData> specular_data);
private:
void swapContent(SimulationElement &other);
@@ -99,7 +104,7 @@ private:
Eigen::Matrix2cd m_analyzer_operator; //!< polarization analyzer operator
#endif
std::unique_ptr<IPixel> mP_pixel;
- bool m_contains_specular;
+ std::unique_ptr<SpecularData> m_specular_data;
};
@@ -109,4 +114,29 @@ void addElementsWithWeight(std::vector<SimulationElement>::const_iterator first,
std::vector<SimulationElement>::iterator result,
double weight);
+//! Helper class for SimulationElement to carry specular information
+//! @ingroup simulation
+
+class SpecularData
+{
+ // FIXME: find a better way to carry the specular data in SimulationElement
+ using ScalarVector = std::vector<ScalarRTCoefficients>;
+ using MatrixVector = std::vector<MatrixRTCoefficients>;
+public:
+ SpecularData();
+
+ SpecularData(MatrixVector coefficients);
+
+ SpecularData(ScalarVector coefficients);
+
+ const ILayerRTCoefficients& operator[](size_t index) const;
+
+ bool isInited() const {return data_type_used != DATA_TYPE::Invalid;}
+
+private:
+ enum class DATA_TYPE { Invalid = -1, Scalar, Matrix };
+ boost::variant<ScalarVector, MatrixVector> data;
+ DATA_TYPE data_type_used;
+};
+
#endif // SIMULATIONELEMENT_H
diff --git a/Core/Computation/SpecularComputationTerm.cpp b/Core/Computation/SpecularComputationTerm.cpp
index b7dc91548a7aaa5f781a0d5e22701a8b1185c6ea..954e51680477c30c66fd18e89a00290db2a720d1 100644
--- a/Core/Computation/SpecularComputationTerm.cpp
+++ b/Core/Computation/SpecularComputationTerm.cpp
@@ -17,12 +17,13 @@ void SpecularComputationTerm::eval(ProgressHandler*,
return;
for (auto it = begin_it; it != end_it; ++it)
- if (it->containsSpecularWavevector())
+ if (it->specularData())
evalSingle(it);
}
void SpecularComputationTerm::evalSingle(const std::vector<SimulationElement>::iterator& iter) const
{
- complex_t R = mp_fresnel_map->getInCoefficients(*iter, 0)->getScalarR();
- iter->setIntensity(std::norm(R));
+ mp_fresnel_map->fillSpecularData(*iter);
+ const ILayerRTCoefficients& layer_data = (*iter->specularData())[0];
+ iter->setIntensity(std::norm(layer_data.getScalarR()));
}
diff --git a/Core/Instrument/IDetector2D.cpp b/Core/Instrument/IDetector2D.cpp
index d4e807e8b1f53fddc96152f8ce2f6e57e3697793..7f933e97896b54bca6aa310d5fcfddb349fa47e8 100644
--- a/Core/Instrument/IDetector2D.cpp
+++ b/Core/Instrument/IDetector2D.cpp
@@ -111,7 +111,7 @@ std::vector<SimulationElement> IDetector2D::createSimulationElements(const Beam
sim_element.setPolarization(beam_polarization);
sim_element.setAnalyzerOperator(analyzer_operator);
if (it.index()==spec_index) {
- sim_element.setSpecular(true);
+ sim_element.setSpecular();
}
result.push_back(std::move(sim_element));
}
diff --git a/Core/Instrument/SpecularDetector1D.cpp b/Core/Instrument/SpecularDetector1D.cpp
index bb6013eaa643ec12a7b048dbf1577e530eb1821e..5aa6d7bd47c051a8a9322a4e28d626a190a00efc 100644
--- a/Core/Instrument/SpecularDetector1D.cpp
+++ b/Core/Instrument/SpecularDetector1D.cpp
@@ -80,7 +80,7 @@ std::vector<SimulationElement> SpecularDetector1D::createSimulationElements(cons
std::make_unique<SpecularPixel>(alpha_i));
sim_element.setPolarization(beam_polarization);
sim_element.setAnalyzerOperator(analyzer_operator);
- sim_element.setSpecular(true);
+ sim_element.setSpecular();
result.push_back(std::move(sim_element));
}
diff --git a/Core/Multilayer/IFresnelMap.h b/Core/Multilayer/IFresnelMap.h
index ab72fc3cf60fca0955c3367ab96c0fd1494a9a7c..5f274abb9fba872ffe75ec1063f843fd00fa6ba2 100644
--- a/Core/Multilayer/IFresnelMap.h
+++ b/Core/Multilayer/IFresnelMap.h
@@ -23,6 +23,7 @@
class ILayerRTCoefficients;
class MultiLayer;
class SimulationElement;
+struct SpecularData;
//! Holds the necessary information to calculate the radiation wavefunction in every layer
//! for different incoming (outgoing) angles of the beam in the top layer
@@ -43,6 +44,9 @@ public:
virtual const ILayerRTCoefficients* getInCoefficients(
const SimulationElement& sim_element, size_t layer_index) const =0;
+ //! Fills simulation element specular data
+ virtual void fillSpecularData(SimulationElement& sim_element) const = 0;
+
//! Sets the multilayer to be used for the Fresnel calculations.
virtual void setMultilayer(const MultiLayer& multilayer);
diff --git a/Core/Multilayer/MatrixFresnelMap.cpp b/Core/Multilayer/MatrixFresnelMap.cpp
index f0459a02f0885466b4b466558870cd40d987117f..546bebf2d4cbc2e61d86e202d4a648247f58258c 100644
--- a/Core/Multilayer/MatrixFresnelMap.cpp
+++ b/Core/Multilayer/MatrixFresnelMap.cpp
@@ -46,6 +46,17 @@ MatrixFresnelMap::getInCoefficients(const SimulationElement& sim_element, size_t
return getCoefficients(sim_element.getKi(), layer_index, *mP_multilayer, m_hash_table_in);
}
+void MatrixFresnelMap::fillSpecularData(SimulationElement& sim_element) const
+{
+ const auto& kvec = sim_element.getKi();
+ std::vector<MatrixRTCoefficients> coef_vector;
+ if (m_use_cache)
+ coef_vector = getCoefficientsFromCache(kvec, *mP_multilayer, m_hash_table_in);
+ else
+ coef_vector = calculateCoefficients(*mP_multilayer, kvec);
+ sim_element.setSpecular(std::make_unique<SpecularData>(std::move(coef_vector)));
+}
+
const ILayerRTCoefficients* MatrixFresnelMap::getCoefficients(kvector_t kvec, size_t layer_index,
const MultiLayer& multilayer,
CoefficientHash& hash_table) const
diff --git a/Core/Multilayer/MatrixFresnelMap.h b/Core/Multilayer/MatrixFresnelMap.h
index 2076f7e046962c2ffa05eb0453090ece3420b43e..bbbb5cd6520f4e09c02d5abd9a91e6a5fe7deb6f 100644
--- a/Core/Multilayer/MatrixFresnelMap.h
+++ b/Core/Multilayer/MatrixFresnelMap.h
@@ -44,6 +44,9 @@ public:
void setMultilayer(const MultiLayer& multilayer) final override;
+ //! Fills simulation element specular data
+ virtual void fillSpecularData(SimulationElement& sim_element) const override;
+
typedef std::unordered_map<kvector_t, std::vector<MatrixRTCoefficients>, HashKVector>
CoefficientHash;
diff --git a/Core/Multilayer/ScalarFresnelMap.cpp b/Core/Multilayer/ScalarFresnelMap.cpp
index f64371cf963a6be7afd4e97282dfb512b59d0225..1276f51cf72a1ae22bd5bcf6238159e7657f22c8 100644
--- a/Core/Multilayer/ScalarFresnelMap.cpp
+++ b/Core/Multilayer/ScalarFresnelMap.cpp
@@ -42,6 +42,17 @@ const ILayerRTCoefficients* ScalarFresnelMap::getInCoefficients(
return getCoefficients(sim_element.getKi(), layer_index);
}
+void ScalarFresnelMap::fillSpecularData(SimulationElement& sim_element) const
+{
+ const auto& kvec = sim_element.getKi();
+ std::vector<ScalarRTCoefficients> coef_vector;
+ if (m_use_cache)
+ coef_vector = getCoefficientsFromCache(kvec);
+ else
+ coef_vector = calculateCoefficients(*mP_multilayer, kvec);
+ sim_element.setSpecular(std::make_unique<SpecularData>(std::move(coef_vector)));
+}
+
const ScalarRTCoefficients* ScalarFresnelMap::getCoefficients(
kvector_t kvec, size_t layer_index) const
{
diff --git a/Core/Multilayer/ScalarFresnelMap.h b/Core/Multilayer/ScalarFresnelMap.h
index 60fb90fa80ff206b73a18671c5a2d116b1b2812b..9f258660425cca80f530b5df5f5b6d3d4ee2e0e8 100644
--- a/Core/Multilayer/ScalarFresnelMap.h
+++ b/Core/Multilayer/ScalarFresnelMap.h
@@ -44,6 +44,9 @@ public:
virtual const ILayerRTCoefficients* getInCoefficients(
const SimulationElement& sim_element, size_t layer_index) const final override;
+ //! Fills simulation element specular data
+ virtual void fillSpecularData(SimulationElement& sim_element) const override;
+
private:
const ScalarRTCoefficients* getCoefficients(kvector_t kvec, size_t layer_index) const;
const std::vector<ScalarRTCoefficients>& getCoefficientsFromCache(kvector_t kvec) const;