diff --git a/Core/Algorithms/inc/DWBASimulation.h b/Core/Algorithms/inc/DWBASimulation.h
index 76329ccdd749dc25737ca1f132569bb49dea4e96..764f187864d90524e62520cda36a0f07ce2872e9 100644
--- a/Core/Algorithms/inc/DWBASimulation.h
+++ b/Core/Algorithms/inc/DWBASimulation.h
@@ -46,12 +46,24 @@ class DWBASimulation : public ISimulation
//! Returns output data containing calculated intensity.
const OutputData<double>& getDWBAIntensity() const
- { return m_dwba_intensity; }
+ {
+ if (mp_polarization_output) return getPolarizationData();
+ return m_dwba_intensity;
+ }
+
+ //! Returns output data containing calculated polarized intensity.
+ const OutputData<Eigen::Matrix2cd>& getPolarizedDWBAIntensity() const
+ { return *mp_polarization_output; }
//! Adds intensity to current dwba intensity
void addDWBAIntensity(const OutputData<double>& data_to_add)
{ m_dwba_intensity += data_to_add; }
+ //! Adds polarized intensity to current polarized dwba intensity
+ void addPolarizedDWBAIntensity(const OutputData<Eigen::Matrix2cd>
+ &data_to_add)
+ { (*mp_polarization_output) += data_to_add; }
+
virtual DWBASimulation *clone() const;
// ---------------------------------
@@ -85,7 +97,10 @@ protected:
//! Returns the wavelength of the incoming beam
double getWaveLength() const;
- OutputData<double> m_dwba_intensity;
+ //! apply beam polarization to get specific polarized intensity map
+ const OutputData<double>& getPolarizationData() const;
+
+ mutable OutputData<double> m_dwba_intensity;
OutputData<Eigen::Matrix2cd> *mp_polarization_output;
cvector_t m_ki;
double m_alpha_i;
diff --git a/Core/Algorithms/inc/LayerDWBASimulation.h b/Core/Algorithms/inc/LayerDWBASimulation.h
index a398baaf9464eccbe8b2773faec603480cf019b8..cacb65702fe34ab3f72f40982c182851a67e605a 100644
--- a/Core/Algorithms/inc/LayerDWBASimulation.h
+++ b/Core/Algorithms/inc/LayerDWBASimulation.h
@@ -19,6 +19,8 @@
#include "DWBASimulation.h"
#include "IDoubleToComplexFunction.h"
+class MagneticCoefficientsMap;
+
//! Base class for LayerDecoratorDWBASimulation, DiffuseDWBASimulation.
class LayerDWBASimulation : public DWBASimulation
@@ -38,6 +40,7 @@ class LayerDWBASimulation : public DWBASimulation
const IDoubleToPairOfComplexMap& rt_map);
void setKzAndRTFunctions(const IDoubleToComplexMap& kz_function,
const IDoubleToPairOfComplexMap& rt_map);
+ void setMagneticCoefficientsMap(const MagneticCoefficientsMap& coeff_map);
protected:
Bin1DCVector getKfBin(double wavelength,
@@ -45,6 +48,7 @@ class LayerDWBASimulation : public DWBASimulation
const Bin1D& phi_bin) const;
IDoubleToComplexMap *mp_kz_function;
IDoubleToPairOfComplexMap *mp_RT_function;
+ MagneticCoefficientsMap *mp_coeff_map;
};
#endif /* LAYERDWBASIMULATION_H_ */
diff --git a/Core/Algorithms/inc/MultiLayerDWBASimulation.h b/Core/Algorithms/inc/MultiLayerDWBASimulation.h
index 54242889b9b03d63f4e7a35c8c7a0b9dff10d500..af2082b37de7e0fd0d9c14359459de4ff0fc16c4 100644
--- a/Core/Algorithms/inc/MultiLayerDWBASimulation.h
+++ b/Core/Algorithms/inc/MultiLayerDWBASimulation.h
@@ -47,6 +47,9 @@ public:
virtual void run();
protected:
+ //! calculates intensity map for samples with magnetization
+ void runMagnetic();
+
std::set<double> getAlphaList() const;
std::map<size_t, LayerDWBASimulation*> m_layer_dwba_simulation_map;
MultiLayer *mp_multi_layer;
diff --git a/Core/Algorithms/src/DWBASimulation.cpp b/Core/Algorithms/src/DWBASimulation.cpp
index 99cdc45a832ff0273803dd7f2e119a04f5400353..33ae8c194b27bad22c9d3feb1348af042595c300 100644
--- a/Core/Algorithms/src/DWBASimulation.cpp
+++ b/Core/Algorithms/src/DWBASimulation.cpp
@@ -41,6 +41,7 @@ void DWBASimulation::init(const Simulation& simulation)
kvector_t ki_real(m_ki.x().real(), m_ki.y().real(), m_ki.z().real());
m_alpha_i = std::asin(ki_real.z()/ki_real.mag());
m_sim_params = simulation.getSimulationParameters();
+
// initialize polarization output if needed
if (checkPolarizationPresent()) {
mp_polarization_output = new OutputData<Eigen::Matrix2cd>();
@@ -94,4 +95,17 @@ double DWBASimulation::getWaveLength() const
return 2*M_PI/real_ki.mag();
}
-
+const OutputData<double>& DWBASimulation::getPolarizationData() const
+{
+ Eigen::Matrix2cd pol_density = mp_simulation->getInstrument()
+ .getBeam().getPolarization();
+ OutputData<double>::iterator it = m_dwba_intensity.begin();
+ OutputData<Eigen::Matrix2cd>::const_iterator mat_it =
+ mp_polarization_output->begin();
+ while (it != m_dwba_intensity.end()) {
+ Eigen::Matrix2cd mat = pol_density * (*mat_it);
+ *it = std::abs(mat.trace());
+ ++it, ++mat_it;
+ }
+ return m_dwba_intensity;
+}
diff --git a/Core/Algorithms/src/LayerDWBASimulation.cpp b/Core/Algorithms/src/LayerDWBASimulation.cpp
index e3f87fd7b553edd636cabaf4358fe56b0b2e9859..0e65ac3f131a4d47deb3f75a4154d69bcfa9bed1 100644
--- a/Core/Algorithms/src/LayerDWBASimulation.cpp
+++ b/Core/Algorithms/src/LayerDWBASimulation.cpp
@@ -14,11 +14,15 @@
// ************************************************************************** //
#include "LayerDWBASimulation.h"
+
+#include "MagneticCoefficientsMap.h"
+
#include <cassert>
LayerDWBASimulation::LayerDWBASimulation()
- : mp_kz_function(0)
- , mp_RT_function(0)
+: mp_kz_function(0)
+, mp_RT_function(0)
+, mp_coeff_map(0)
{
}
@@ -26,6 +30,7 @@ LayerDWBASimulation::~LayerDWBASimulation()
{
delete mp_kz_function;
delete mp_RT_function;
+ delete mp_coeff_map;
}
void LayerDWBASimulation::setKzFunction(const IDoubleToComplexMap& kz_function)
@@ -46,6 +51,12 @@ void LayerDWBASimulation::setKzAndRTFunctions(const IDoubleToComplexMap& kz_func
setReflectionTransmissionFunction(rt_map);
}
+void LayerDWBASimulation::setMagneticCoefficientsMap(
+ const MagneticCoefficientsMap& coeff_map)
+{
+ mp_coeff_map = coeff_map.clone();
+}
+
Bin1DCVector LayerDWBASimulation::getKfBin(double wavelength, const Bin1D& alpha_bin, const Bin1D& phi_bin) const
{
assert(mp_kz_function != NULL);
diff --git a/Core/Algorithms/src/MultiLayerDWBASimulation.cpp b/Core/Algorithms/src/MultiLayerDWBASimulation.cpp
index 5d8185b7f12ef1b4c9ad0cba05a366491d7cecff..0e046b5db9b7de44333e1abd94410fb4fe75f040 100644
--- a/Core/Algorithms/src/MultiLayerDWBASimulation.cpp
+++ b/Core/Algorithms/src/MultiLayerDWBASimulation.cpp
@@ -16,11 +16,13 @@
#include "MultiLayerDWBASimulation.h"
#include "SpecularMatrix.h"
+#include "SpecularMagnetic.h"
#include "MultiLayer.h"
#include "DoubleToComplexInterpolatingFunction.h"
#include "MultiLayerRoughnessDWBASimulation.h"
#include "DoubleToComplexMap.h"
#include "MessageService.h"
+#include "MagneticCoefficientsMap.h"
MultiLayerDWBASimulation::MultiLayerDWBASimulation(
@@ -148,6 +150,62 @@ void MultiLayerDWBASimulation::run()
}
}
+void MultiLayerDWBASimulation::runMagnetic()
+{
+ SpecularMagnetic specularCalculator;
+
+ kvector_t m_ki_real(m_ki.x().real(), m_ki.y().real(), m_ki.z().real());
+
+ m_dwba_intensity.setAllTo(0.0);
+ mp_polarization_output->setAllTo(Eigen::Matrix2cd::Zero());
+ double lambda = 2*M_PI/m_ki_real.mag();
+
+ // collect all alpha angles and calculate Fresnel coefficients
+ typedef std::pair<double, SpecularMagnetic::MultiLayerCoeff_t>
+ doubleFresnelPair_t;
+ std::vector<doubleFresnelPair_t> doubleFresnel_buffer;
+ std::set<double> alpha_set = getAlphaList();
+ doubleFresnel_buffer.reserve(alpha_set.size());
+
+ double angle;
+ kvector_t kvec;
+ SpecularMagnetic::MultiLayerCoeff_t coeffs;
+ for (std::set<double>::const_iterator it =
+ alpha_set.begin(); it != alpha_set.end(); ++it) {
+ angle = *it;
+ kvec.setLambdaAlphaPhi(lambda, -angle, 0.0);
+ specularCalculator.execute(*mp_multi_layer, kvec, coeffs);
+ doubleFresnel_buffer.push_back( doubleFresnelPair_t(angle,coeffs) );
+ }
+
+ // run through layers and add DWBA calculated from these layers
+ for(size_t i_layer=0;
+ i_layer<mp_multi_layer->getNumberOfLayers(); ++i_layer) {
+ msglog(MSG::DEBUG) << "MultiLayerDWBASimulation::runMagnetic()"
+ "-> Layer " << i_layer;
+ MagneticCoefficientsMap coeff_map;
+
+ for(std::vector<doubleFresnelPair_t >::const_iterator it=
+ doubleFresnel_buffer.begin();
+ it!=doubleFresnel_buffer.end(); ++it) {
+ double angle = (*it).first;
+ const SpecularMagnetic::LayerMatrixCoeff& coeff = (*it).second[i_layer];
+ coeff_map[angle] = coeff;
+ }
+
+ // layer DWBA simulation
+ std::map<size_t, LayerDWBASimulation*>::const_iterator pos =
+ m_layer_dwba_simulation_map.find(i_layer);
+ if(pos != m_layer_dwba_simulation_map.end() ) {
+ LayerDWBASimulation *p_layer_dwba_sim = pos->second;
+ p_layer_dwba_sim->setMagneticCoefficientsMap(coeff_map);
+ p_layer_dwba_sim->run();
+ addPolarizedDWBAIntensity(
+ p_layer_dwba_sim->getPolarizedDWBAIntensity() );
+ }
+ } // i_layer
+}
+
std::set<double> MultiLayerDWBASimulation::getAlphaList() const
{
std::set<double> result;
diff --git a/Core/Tools/inc/DoubleToComplexMap.h b/Core/Tools/inc/DoubleToComplexMap.h
index 0c6847fb029b8404162cbaf94a289fef165dd950..f145bf4152ab900604cf80c2251ee0142b920e48 100644
--- a/Core/Tools/inc/DoubleToComplexMap.h
+++ b/Core/Tools/inc/DoubleToComplexMap.h
@@ -16,8 +16,6 @@
#ifndef DOUBLETOCOMPLEXUNORDEREDMAP_H
#define DOUBLETOCOMPLEXUNORDEREDMAP_H
-//#include "Exceptions.h"
-//#include "Types.h"
#include "Utils.h"
#include "IDoubleToComplexFunction.h"
diff --git a/Core/Tools/inc/MagneticCoefficientsMap.h b/Core/Tools/inc/MagneticCoefficientsMap.h
new file mode 100644
index 0000000000000000000000000000000000000000..6734fba9cf04dad8239e48d1578caab9cd98fb4d
--- /dev/null
+++ b/Core/Tools/inc/MagneticCoefficientsMap.h
@@ -0,0 +1,52 @@
+// ************************************************************************** //
+//
+// BornAgain: simulate and fit scattering at grazing incidence
+//
+//! @file Tools/inc/MagneticCoefficientsMap.h
+//! @brief Defines class MagneticCoefficientsMap.
+//!
+//! @homepage http://apps.jcns.fz-juelich.de/BornAgain
+//! @license GNU General Public License v3 or higher (see COPYING)
+//! @copyright Forschungszentrum Jülich GmbH 2013
+//! @authors Scientific Computing Group at MLZ Garching
+//! @authors C. Durniak, G. Pospelov, W. Van Herck, J. Wuttke
+//
+// ************************************************************************** //
+
+#ifndef MAGNETICCOEFFICIENTSMAP_H_
+#define MAGNETICCOEFFICIENTSMAP_H_
+
+#include "Utils.h"
+#include "SpecularMagnetic.h"
+
+//! Map from angles (double) to matrix coefficients for magnetic
+//! DWBA calculation.
+
+class MagneticCoefficientsMap
+{
+ public:
+ typedef Utils::UnorderedMap<double, SpecularMagnetic::LayerMatrixCoeff>
+ container_t;
+
+ MagneticCoefficientsMap(){}
+ MagneticCoefficientsMap(const container_t& value_map) : m_value_map(value_map) {}
+
+ SpecularMagnetic::LayerMatrixCoeff& operator[] (double key) {
+ return m_value_map[key];
+ }
+ MagneticCoefficientsMap *clone() const {
+ return new MagneticCoefficientsMap(m_value_map);
+ }
+ const SpecularMagnetic::LayerMatrixCoeff& evaluate(double value) const {
+ return m_value_map.find(value);
+ }
+ private:
+ container_t m_value_map;
+};
+
+//! Double to pair of complex unordered map.
+
+
+
+
+#endif /* MAGNETICCOEFFICIENTSMAP_H_ */