diff --git a/Device/Beam/Beam.cpp b/Device/Beam/Beam.cpp
index 7a5b11dc79431f9c081d9d639e0fc5be093d82cf..40847f5a8204672d3b83206f0f880ab5ef0b6bbc 100644
--- a/Device/Beam/Beam.cpp
+++ b/Device/Beam/Beam.cpp
@@ -14,6 +14,7 @@
 
 #include "Device/Beam/Beam.h"
 #include "Base/Math/Constants.h"
+#include "Base/Spin/SpinMatrix.h"
 #include "Base/Util/Assert.h"
 #include "Device/Beam/FootprintGauss.h"
 #include <heinz/Complex.h>
diff --git a/Device/Beam/Beam.h b/Device/Beam/Beam.h
index ef1458db054b13fa2fe562d2e0d0234c871ef7a3..5d7031a4e3040f5ab835fb4e6278e358507036da 100644
--- a/Device/Beam/Beam.h
+++ b/Device/Beam/Beam.h
@@ -19,11 +19,8 @@
 #include "Fit/Param/RealLimits.h"
 #include "Param/Node/INode.h"
 
-#ifndef SWIG
-#include "Base/Spin/SpinMatrix.h"
-#endif
-
 class IFootprintFactor;
+class SpinMatrix;
 
 //! An incident neutron or x-ray beam.
 //! @ingroup beam
@@ -64,10 +61,8 @@ public:
     //! Throws if limits are violated.
     void setWavelengthGuarded(double value);
 
-#ifndef SWIG
     //! Returns the polarization density matrix (in spin basis along z-axis)
     SpinMatrix getPolarization() const;
-#endif
 
     void setWavelength(double wavelength);
     void setDirection(const Direction& direction);
diff --git a/Device/Detector/DetectorContext.h b/Device/Detector/DetectorContext.h
index 678d916bcd34edc1ccdf71becd4cee0c9f93355d..54b05188992b1a1fff4851963a87088db9d2d285 100644
--- a/Device/Detector/DetectorContext.h
+++ b/Device/Detector/DetectorContext.h
@@ -17,13 +17,10 @@
 #define BORNAGAIN_DEVICE_DETECTOR_DETECTORCONTEXT_H
 
 #include "Base/Pixel/IPixel.h"
+#include "Base/Spin/SpinMatrix.h"
 #include <memory>
 #include <vector>
 
-#ifndef SWIG
-#include "Base/Spin/SpinMatrix.h"
-#endif
-
 class IDetector2D;
 
 //! Holds precalculated information for faster DiffuseElement generation.
@@ -45,9 +42,7 @@ public:
 private:
     void setup_context(const IDetector2D* detector);
 
-#ifndef SWIG
     SpinMatrix m_analyzer_operator;
-#endif
     std::vector<std::unique_ptr<IPixel>> m_pixels; //! All unmasked pixels inside ROI.
     std::vector<size_t> m_active_indices; //! The sequence of bin indices (unmasked, in ROI)
 };
diff --git a/Device/Pol/PolFilter.cpp b/Device/Pol/PolFilter.cpp
index 2f5662549ec951e7278ec9f705972f1f5cf56ba1..635d583b4459c71d8aa2bbdd4fcc3d8ed4f08aaf 100644
--- a/Device/Pol/PolFilter.cpp
+++ b/Device/Pol/PolFilter.cpp
@@ -13,6 +13,7 @@
 //  ************************************************************************************************
 
 #include "Device/Pol/PolFilter.h"
+#include "Base/Spin/SpinMatrix.h"
 #include "Fit/Param/RealLimits.h"
 #include <heinz/Complex.h>
 
diff --git a/Device/Pol/PolFilter.h b/Device/Pol/PolFilter.h
index d30463b3f3872b55fba9e631d7e00f95e3c3f87c..c6f5e91f282e95e379b7f2e86581fd5925fc92fe 100644
--- a/Device/Pol/PolFilter.h
+++ b/Device/Pol/PolFilter.h
@@ -19,9 +19,7 @@
 #include "Param/Node/INode.h"
 #include <heinz/Vectors3D.h>
 
-#ifndef SWIG
-#include "Base/Spin/SpinMatrix.h"
-#endif
+class SpinMatrix;
 
 //! Detector properties (efficiency, transmission).
 //! @ingroup detector
@@ -43,10 +41,8 @@ public:
     //! Sets the polarization analyzer characteristics of the detector
     void setDirEffTra(R3 direction, double efficiency, double total_transmission);
 
-#ifndef SWIG
     //! Return the polarization density matrix (in spin basis along z-axis)
     SpinMatrix matrix() const;
-#endif
 
     //! Retrieve the analyzer characteristics
     R3 polDirection() const;
diff --git a/Resample/Particle/IReParticle.cpp b/Resample/Particle/IReParticle.cpp
index d378f9bfd749d682257b91fb17eb375d19eda0f3..380407ba6441e36530a10d93833b839738a7a48b 100644
--- a/Resample/Particle/IReParticle.cpp
+++ b/Resample/Particle/IReParticle.cpp
@@ -13,6 +13,7 @@
 //  ************************************************************************************************
 
 #include "Resample/Particle/IReParticle.h"
+#include "Base/Spin/SpinMatrix.h"
 #include "Base/Vector/WavevectorInfo.h"
 #include <memory>
 
diff --git a/Resample/Particle/IReParticle.h b/Resample/Particle/IReParticle.h
index 225a43bb066bf31054e8d0bb90f7c9bdfd788d55..f47d00e9a9e516dd373c18353001497cf85471a1 100644
--- a/Resample/Particle/IReParticle.h
+++ b/Resample/Particle/IReParticle.h
@@ -20,13 +20,13 @@
 #ifndef BORNAGAIN_RESAMPLE_PARTICLE_IREPARTICLE_H
 #define BORNAGAIN_RESAMPLE_PARTICLE_IREPARTICLE_H
 
-#include "Base/Spin/SpinMatrix.h"
 #include "Base/Types/ICloneable.h"
 #include <heinz/Complex.h>
 #include <heinz/Vectors3D.h>
 
 class IRotation;
 class Material;
+class SpinMatrix;
 class WavevectorInfo;
 
 //! Abstract base class for reprocessed particles.
@@ -59,14 +59,11 @@ public:
     //! Returns the z-coordinate of the lowest point in this shape after a given rotation
     virtual double topZ(const IRotation* rotation) const = 0;
 
-#ifndef SWIG
     //! Returns scattering amplitude for matrix interactions
     virtual SpinMatrix thePolFF(const WavevectorInfo& wavevectors) const;
 
     //! Returns the total volume of the particle of this form factor's shape
     virtual double volume() const;
-
-#endif // SWIG
 };
 
 #endif // BORNAGAIN_RESAMPLE_PARTICLE_IREPARTICLE_H
diff --git a/Resample/Particle/ReCoreShell.cpp b/Resample/Particle/ReCoreShell.cpp
index b49d8a3d7c28313e3f0a9d5183e58f5910d72b6d..f1f5d2d8d84715174e87c37472bda255532173ee 100644
--- a/Resample/Particle/ReCoreShell.cpp
+++ b/Resample/Particle/ReCoreShell.cpp
@@ -13,6 +13,7 @@
 //  ************************************************************************************************
 
 #include "Resample/Particle/ReCoreShell.h"
+#include "Base/Spin/SpinMatrix.h"
 #include "Base/Vector/WavevectorInfo.h" // debug
 
 ReCoreShell::ReCoreShell(IReParticle* core, IReParticle* shell)
diff --git a/Resample/Particle/ReCoreShell.h b/Resample/Particle/ReCoreShell.h
index 71935dc166063f8862208e0649c1fd8ce18cd0a2..b953dbf9906b5eff020dd7f42fda20291608a4ca 100644
--- a/Resample/Particle/ReCoreShell.h
+++ b/Resample/Particle/ReCoreShell.h
@@ -42,11 +42,7 @@ public:
     void setAmbientMaterial(const Material& material) override;
 
     complex_t theFF(const WavevectorInfo& wavevectors) const override;
-
-#ifndef SWIG
-    //! Calculates and returns a polarized form factor calculation in DWBA
     SpinMatrix thePolFF(const WavevectorInfo& wavevectors) const override;
-#endif
 
 protected:
     std::unique_ptr<IReParticle> m_core;
diff --git a/Resample/Particle/ReMesocrystal.cpp b/Resample/Particle/ReMesocrystal.cpp
index f4cd2f2880c451373e0c9daab317e1a9d5759147..a276572623a79d1ec17cbb966d9b3e4f001a9e53 100644
--- a/Resample/Particle/ReMesocrystal.cpp
+++ b/Resample/Particle/ReMesocrystal.cpp
@@ -14,6 +14,7 @@
 
 #include "Resample/Particle/ReMesocrystal.h"
 #include "Base/Math/Constants.h"
+#include "Base/Spin/SpinMatrix.h"
 #include "Base/Vector/WavevectorInfo.h"
 #include "Resample/Particle/ReParticle.h"
 
diff --git a/Resample/Particle/ReMesocrystal.h b/Resample/Particle/ReMesocrystal.h
index 3f974285b04bebf0ab911db3a5154195309fceea..8cf2ed003bd4d4342f6fd55a12f88f1408bea572 100644
--- a/Resample/Particle/ReMesocrystal.h
+++ b/Resample/Particle/ReMesocrystal.h
@@ -50,9 +50,7 @@ public:
     double topZ(const IRotation* rotation) const override;
 
     complex_t theFF(const WavevectorInfo& wavevectors) const override;
-#ifndef SWIG
     SpinMatrix thePolFF(const WavevectorInfo& wavevectors) const override;
-#endif
 
 private:
     void calculateLargestReciprocalDistance();
diff --git a/Resample/Particle/ReParticle.h b/Resample/Particle/ReParticle.h
index 3cd8edb8704bd37bd6b7bc15d031374053a2e175..868e3e5d88c2cefa5a0e20630b0586597849e980 100644
--- a/Resample/Particle/ReParticle.h
+++ b/Resample/Particle/ReParticle.h
@@ -47,9 +47,7 @@ public:
     void setRotMatrix(const RotMatrix& rotMatrix);
 
     complex_t theFF(const WavevectorInfo& wavevectors) const override;
-#ifndef SWIG
     SpinMatrix thePolFF(const WavevectorInfo& wavevectors) const override;
-#endif
 
     double volume() const override;
 
diff --git a/Sample/Material/Material.cpp b/Sample/Material/Material.cpp
index b4764e9041c00ed437c2517796ca66002b68aa63..89ddb6c1ed84c442beea059d140f8b61b82bd3f0 100644
--- a/Sample/Material/Material.cpp
+++ b/Sample/Material/Material.cpp
@@ -13,6 +13,7 @@
 //  ************************************************************************************************
 
 #include "Sample/Material/Material.h"
+#include "Base/Spin/SpinMatrix.h"
 #include "Base/Util/Assert.h"
 #include "Base/Vector/RotMatrix.h"
 #include "Base/Vector/WavevectorInfo.h"
diff --git a/Sample/Material/Material.h b/Sample/Material/Material.h
index 40a6c2ca3d7ba7c813b055277307c599ec32aaad..8ddf81e3767a01cfc8d28173bbfc400adc0979e7 100644
--- a/Sample/Material/Material.h
+++ b/Sample/Material/Material.h
@@ -25,6 +25,7 @@
 #endif // SWIG
 
 class RotMatrix;
+class SpinMatrix;
 class WavevectorInfo;
 
 //! A wrapper for underlying material implementation
@@ -86,10 +87,8 @@ public:
     //! sld (in \f$nm^{-2}\f$) being the scattering length density
     complex_t scalarSubtrSLD(const WavevectorInfo& wavevectors) const;
 
-#ifndef SWIG
     //! Returns (\f$ \pi/\lambda^2 \f$ - sld) matrix with magnetization corrections
     SpinMatrix polarizedSubtrSLD(const WavevectorInfo& wavevectors) const;
-#endif
 
     Material rotatedMaterial(const RotMatrix& transform) const;
 
diff --git a/Sample/Particle/IFormFactor.h b/Sample/Particle/IFormFactor.h
index 98c9e7269c32610b13256d64deec700ed610d908..e3b7a9472a2e7901239a621067be343a010a5465 100644
--- a/Sample/Particle/IFormFactor.h
+++ b/Sample/Particle/IFormFactor.h
@@ -20,10 +20,8 @@
 #include <heinz/Complex.h>
 #include <heinz/Vectors3D.h>
 
-#ifndef SWIG
 #include "Base/Spin/SpinMatrix.h"
 #include "Base/Spin/Spinor.h"
-#endif
 
 class IRotation;
 class IShape3D;
@@ -44,12 +42,10 @@ public:
     IFormFactor* clone() const override = 0;
 
     virtual complex_t theFF(const WavevectorInfo& wavevectors) const;
+    virtual SpinMatrix thePolFF(const WavevectorInfo& wavevectors) const;
 
     std::string shapeName() const;
 
-#ifndef SWIG
-    virtual SpinMatrix thePolFF(const WavevectorInfo& wavevectors) const;
-#endif
     virtual double bottomZ(const IRotation* rotation) const;
     virtual double topZ(const IRotation* rotation) const;
 
@@ -67,12 +63,10 @@ public:
     //! Default implementation only allows rotations along z-axis
     virtual bool canSliceAnalytically(const IRotation* rot) const;
 
-#ifndef SWIG
     //! Returns scattering amplitude for complex scattering wavevector q=k_i-k_f in case
     //! of matrix interactions. Default implementation calls formfactor_at_bottom(q) and
     //! multiplies with the unit matrix.
     virtual SpinMatrix formfactor_pol(C3 q) const;
-#endif
 
 protected:
     //! IShape3D object, used to retrieve vertices (which may be approximate in the case
diff --git a/auto/Wrap/doxygenResample.i b/auto/Wrap/doxygenResample.i
index 152c2b290dc595a89b62210c870ff466c7e02afc..14c8601bbfe14d5dd5ee266bdacaa9cbd2b51ec8 100644
--- a/auto/Wrap/doxygenResample.i
+++ b/auto/Wrap/doxygenResample.i
@@ -463,7 +463,7 @@ Returns scattering amplitude for complex wavevectors ki, kf.
 
 %feature("docstring")  ReCoreShell::thePolFF "SpinMatrix ReCoreShell::thePolFF(const WavevectorInfo &wavevectors) const override
 
-Calculates and returns a polarized form factor calculation in DWBA. 
+Returns scattering amplitude for matrix interactions. 
 ";
 
 
diff --git a/auto/Wrap/doxygenSample.i b/auto/Wrap/doxygenSample.i
index a8cc3f1e3cf972c194d8b301d6f2ccc993ffabdc..8c2f6c8ba3162fc259c79a38acb0b55b32bec8c5 100644
--- a/auto/Wrap/doxygenSample.i
+++ b/auto/Wrap/doxygenSample.i
@@ -1189,10 +1189,10 @@ Returns a clone of this  ISampleNode object.
 %feature("docstring")  IFormFactor::theFF "complex_t IFormFactor::theFF(const WavevectorInfo &wavevectors) const
 ";
 
-%feature("docstring")  IFormFactor::shapeName "std::string IFormFactor::shapeName() const
+%feature("docstring")  IFormFactor::thePolFF "SpinMatrix IFormFactor::thePolFF(const WavevectorInfo &wavevectors) const
 ";
 
-%feature("docstring")  IFormFactor::thePolFF "SpinMatrix IFormFactor::thePolFF(const WavevectorInfo &wavevectors) const
+%feature("docstring")  IFormFactor::shapeName "std::string IFormFactor::shapeName() const
 ";
 
 %feature("docstring")  IFormFactor::bottomZ "double IFormFactor::bottomZ(const IRotation *rotation) const
diff --git a/auto/Wrap/libBornAgainDevice.py b/auto/Wrap/libBornAgainDevice.py
index c1e1a1084c1a893dc0ab8eaa2af787f3597b5713..20deb845406f0ccfd41a6e2378d8c3b8e0b4a72f 100644
--- a/auto/Wrap/libBornAgainDevice.py
+++ b/auto/Wrap/libBornAgainDevice.py
@@ -3015,6 +3015,16 @@ class Beam(libBornAgainParam.INode):
         """
         return _libBornAgainDevice.Beam_setWavelengthGuarded(self, value)
 
+    def getPolarization(self):
+        r"""
+        getPolarization(Beam self) -> SpinMatrix
+        SpinMatrix Beam::getPolarization() const
+
+        Returns the polarization density matrix (in spin basis along z-axis) 
+
+        """
+        return _libBornAgainDevice.Beam_getPolarization(self)
+
     def setWavelength(self, wavelength):
         r"""
         setWavelength(Beam self, double wavelength)
diff --git a/auto/Wrap/libBornAgainDevice_wrap.cpp b/auto/Wrap/libBornAgainDevice_wrap.cpp
index 07d97e3d00afe9ae171c2364631b3689dbb7ceed..891f31997fe12ea36272f061bc4864eb1a5f9e08 100644
--- a/auto/Wrap/libBornAgainDevice_wrap.cpp
+++ b/auto/Wrap/libBornAgainDevice_wrap.cpp
@@ -3141,73 +3141,74 @@ namespace Swig {
 #define SWIGTYPE_p_ResolutionFunction2DGaussian swig_types[41]
 #define SWIGTYPE_p_SimulationResult swig_types[42]
 #define SWIGTYPE_p_SphericalDetector swig_types[43]
-#define SWIGTYPE_p_Vec3T_double_t swig_types[44]
-#define SWIGTYPE_p_Vec3T_int_t swig_types[45]
-#define SWIGTYPE_p_Vec3T_std__complexT_double_t_t swig_types[46]
-#define SWIGTYPE_p_VerticalLine swig_types[47]
-#define SWIGTYPE_p_allocator_type swig_types[48]
-#define SWIGTYPE_p_bool swig_types[49]
-#define SWIGTYPE_p_char swig_types[50]
-#define SWIGTYPE_p_const_iterator swig_types[51]
-#define SWIGTYPE_p_corr_matrix_t swig_types[52]
-#define SWIGTYPE_p_difference_type swig_types[53]
-#define SWIGTYPE_p_double swig_types[54]
-#define SWIGTYPE_p_first_type swig_types[55]
-#define SWIGTYPE_p_int swig_types[56]
-#define SWIGTYPE_p_iterator swig_types[57]
-#define SWIGTYPE_p_key_type swig_types[58]
-#define SWIGTYPE_p_long_long swig_types[59]
-#define SWIGTYPE_p_mapped_type swig_types[60]
-#define SWIGTYPE_p_p_ICoordSystem swig_types[61]
-#define SWIGTYPE_p_p_PyObject swig_types[62]
-#define SWIGTYPE_p_parameters_t swig_types[63]
-#define SWIGTYPE_p_second_type swig_types[64]
-#define SWIGTYPE_p_short swig_types[65]
-#define SWIGTYPE_p_signed_char swig_types[66]
-#define SWIGTYPE_p_size_type swig_types[67]
-#define SWIGTYPE_p_std__allocatorT_Vec3T_double_t_t swig_types[68]
-#define SWIGTYPE_p_std__allocatorT_Vec3T_std__complexT_double_t_t_t swig_types[69]
-#define SWIGTYPE_p_std__allocatorT_double_t swig_types[70]
-#define SWIGTYPE_p_std__allocatorT_int_t swig_types[71]
-#define SWIGTYPE_p_std__allocatorT_std__complexT_double_t_t swig_types[72]
-#define SWIGTYPE_p_std__allocatorT_std__pairT_double_double_t_t swig_types[73]
-#define SWIGTYPE_p_std__allocatorT_std__pairT_std__string_const_double_t_t swig_types[74]
-#define SWIGTYPE_p_std__allocatorT_std__string_t swig_types[75]
-#define SWIGTYPE_p_std__allocatorT_std__vectorT_double_std__allocatorT_double_t_t_t swig_types[76]
-#define SWIGTYPE_p_std__allocatorT_std__vectorT_int_std__allocatorT_int_t_t_t swig_types[77]
-#define SWIGTYPE_p_std__allocatorT_unsigned_long_t swig_types[78]
-#define SWIGTYPE_p_std__arrayT_double_3_t swig_types[79]
-#define SWIGTYPE_p_std__arrayT_std__complexT_double_t_3_t swig_types[80]
-#define SWIGTYPE_p_std__complexT_double_t swig_types[81]
-#define SWIGTYPE_p_std__functionT_void_fSimulationAreaIterator_const_RF_t swig_types[82]
-#define SWIGTYPE_p_std__invalid_argument swig_types[83]
-#define SWIGTYPE_p_std__lessT_std__string_t swig_types[84]
-#define SWIGTYPE_p_std__mapT_std__string_double_std__lessT_std__string_t_std__allocatorT_std__pairT_std__string_const_double_t_t_t swig_types[85]
-#define SWIGTYPE_p_std__pairT_double_double_t swig_types[86]
-#define SWIGTYPE_p_std__vectorT_AxisInfo_std__allocatorT_AxisInfo_t_t swig_types[87]
-#define SWIGTYPE_p_std__vectorT_INode_const_p_std__allocatorT_INode_const_p_t_t swig_types[88]
-#define SWIGTYPE_p_std__vectorT_ParaMeta_std__allocatorT_ParaMeta_t_t swig_types[89]
-#define SWIGTYPE_p_std__vectorT_Vec3T_double_t_std__allocatorT_Vec3T_double_t_t_t swig_types[90]
-#define SWIGTYPE_p_std__vectorT_Vec3T_std__complexT_double_t_t_std__allocatorT_Vec3T_std__complexT_double_t_t_t_t swig_types[91]
-#define SWIGTYPE_p_std__vectorT_double_std__allocatorT_double_t_t swig_types[92]
-#define SWIGTYPE_p_std__vectorT_int_std__allocatorT_int_t_t swig_types[93]
-#define SWIGTYPE_p_std__vectorT_size_t_std__allocatorT_size_t_t_t swig_types[94]
-#define SWIGTYPE_p_std__vectorT_std__complexT_double_t_std__allocatorT_std__complexT_double_t_t_t swig_types[95]
-#define SWIGTYPE_p_std__vectorT_std__pairT_double_double_t_std__allocatorT_std__pairT_double_double_t_t_t swig_types[96]
-#define SWIGTYPE_p_std__vectorT_std__string_std__allocatorT_std__string_t_t swig_types[97]
-#define SWIGTYPE_p_std__vectorT_std__unique_ptrT_DiffuseElement_t_std__allocatorT_std__unique_ptrT_DiffuseElement_t_t_t swig_types[98]
-#define SWIGTYPE_p_std__vectorT_std__vectorT_double_std__allocatorT_double_t_t_std__allocatorT_std__vectorT_double_std__allocatorT_double_t_t_t_t swig_types[99]
-#define SWIGTYPE_p_std__vectorT_std__vectorT_int_std__allocatorT_int_t_t_std__allocatorT_std__vectorT_int_std__allocatorT_int_t_t_t_t swig_types[100]
-#define SWIGTYPE_p_std__vectorT_unsigned_int_std__allocatorT_unsigned_int_t_t swig_types[101]
-#define SWIGTYPE_p_std__vectorT_unsigned_long_std__allocatorT_unsigned_long_t_t swig_types[102]
-#define SWIGTYPE_p_swig__SwigPyIterator swig_types[103]
-#define SWIGTYPE_p_unsigned_char swig_types[104]
-#define SWIGTYPE_p_unsigned_int swig_types[105]
-#define SWIGTYPE_p_unsigned_long_long swig_types[106]
-#define SWIGTYPE_p_unsigned_short swig_types[107]
-#define SWIGTYPE_p_value_type swig_types[108]
-static swig_type_info *swig_types[110];
-static swig_module_info swig_module = {swig_types, 109, 0, 0, 0, 0};
+#define SWIGTYPE_p_SpinMatrix swig_types[44]
+#define SWIGTYPE_p_Vec3T_double_t swig_types[45]
+#define SWIGTYPE_p_Vec3T_int_t swig_types[46]
+#define SWIGTYPE_p_Vec3T_std__complexT_double_t_t swig_types[47]
+#define SWIGTYPE_p_VerticalLine swig_types[48]
+#define SWIGTYPE_p_allocator_type swig_types[49]
+#define SWIGTYPE_p_bool swig_types[50]
+#define SWIGTYPE_p_char swig_types[51]
+#define SWIGTYPE_p_const_iterator swig_types[52]
+#define SWIGTYPE_p_corr_matrix_t swig_types[53]
+#define SWIGTYPE_p_difference_type swig_types[54]
+#define SWIGTYPE_p_double swig_types[55]
+#define SWIGTYPE_p_first_type swig_types[56]
+#define SWIGTYPE_p_int swig_types[57]
+#define SWIGTYPE_p_iterator swig_types[58]
+#define SWIGTYPE_p_key_type swig_types[59]
+#define SWIGTYPE_p_long_long swig_types[60]
+#define SWIGTYPE_p_mapped_type swig_types[61]
+#define SWIGTYPE_p_p_ICoordSystem swig_types[62]
+#define SWIGTYPE_p_p_PyObject swig_types[63]
+#define SWIGTYPE_p_parameters_t swig_types[64]
+#define SWIGTYPE_p_second_type swig_types[65]
+#define SWIGTYPE_p_short swig_types[66]
+#define SWIGTYPE_p_signed_char swig_types[67]
+#define SWIGTYPE_p_size_type swig_types[68]
+#define SWIGTYPE_p_std__allocatorT_Vec3T_double_t_t swig_types[69]
+#define SWIGTYPE_p_std__allocatorT_Vec3T_std__complexT_double_t_t_t swig_types[70]
+#define SWIGTYPE_p_std__allocatorT_double_t swig_types[71]
+#define SWIGTYPE_p_std__allocatorT_int_t swig_types[72]
+#define SWIGTYPE_p_std__allocatorT_std__complexT_double_t_t swig_types[73]
+#define SWIGTYPE_p_std__allocatorT_std__pairT_double_double_t_t swig_types[74]
+#define SWIGTYPE_p_std__allocatorT_std__pairT_std__string_const_double_t_t swig_types[75]
+#define SWIGTYPE_p_std__allocatorT_std__string_t swig_types[76]
+#define SWIGTYPE_p_std__allocatorT_std__vectorT_double_std__allocatorT_double_t_t_t swig_types[77]
+#define SWIGTYPE_p_std__allocatorT_std__vectorT_int_std__allocatorT_int_t_t_t swig_types[78]
+#define SWIGTYPE_p_std__allocatorT_unsigned_long_t swig_types[79]
+#define SWIGTYPE_p_std__arrayT_double_3_t swig_types[80]
+#define SWIGTYPE_p_std__arrayT_std__complexT_double_t_3_t swig_types[81]
+#define SWIGTYPE_p_std__complexT_double_t swig_types[82]
+#define SWIGTYPE_p_std__functionT_void_fSimulationAreaIterator_const_RF_t swig_types[83]
+#define SWIGTYPE_p_std__invalid_argument swig_types[84]
+#define SWIGTYPE_p_std__lessT_std__string_t swig_types[85]
+#define SWIGTYPE_p_std__mapT_std__string_double_std__lessT_std__string_t_std__allocatorT_std__pairT_std__string_const_double_t_t_t swig_types[86]
+#define SWIGTYPE_p_std__pairT_double_double_t swig_types[87]
+#define SWIGTYPE_p_std__vectorT_AxisInfo_std__allocatorT_AxisInfo_t_t swig_types[88]
+#define SWIGTYPE_p_std__vectorT_INode_const_p_std__allocatorT_INode_const_p_t_t swig_types[89]
+#define SWIGTYPE_p_std__vectorT_ParaMeta_std__allocatorT_ParaMeta_t_t swig_types[90]
+#define SWIGTYPE_p_std__vectorT_Vec3T_double_t_std__allocatorT_Vec3T_double_t_t_t swig_types[91]
+#define SWIGTYPE_p_std__vectorT_Vec3T_std__complexT_double_t_t_std__allocatorT_Vec3T_std__complexT_double_t_t_t_t swig_types[92]
+#define SWIGTYPE_p_std__vectorT_double_std__allocatorT_double_t_t swig_types[93]
+#define SWIGTYPE_p_std__vectorT_int_std__allocatorT_int_t_t swig_types[94]
+#define SWIGTYPE_p_std__vectorT_size_t_std__allocatorT_size_t_t_t swig_types[95]
+#define SWIGTYPE_p_std__vectorT_std__complexT_double_t_std__allocatorT_std__complexT_double_t_t_t swig_types[96]
+#define SWIGTYPE_p_std__vectorT_std__pairT_double_double_t_std__allocatorT_std__pairT_double_double_t_t_t swig_types[97]
+#define SWIGTYPE_p_std__vectorT_std__string_std__allocatorT_std__string_t_t swig_types[98]
+#define SWIGTYPE_p_std__vectorT_std__unique_ptrT_DiffuseElement_t_std__allocatorT_std__unique_ptrT_DiffuseElement_t_t_t swig_types[99]
+#define SWIGTYPE_p_std__vectorT_std__vectorT_double_std__allocatorT_double_t_t_std__allocatorT_std__vectorT_double_std__allocatorT_double_t_t_t_t swig_types[100]
+#define SWIGTYPE_p_std__vectorT_std__vectorT_int_std__allocatorT_int_t_t_std__allocatorT_std__vectorT_int_std__allocatorT_int_t_t_t_t swig_types[101]
+#define SWIGTYPE_p_std__vectorT_unsigned_int_std__allocatorT_unsigned_int_t_t swig_types[102]
+#define SWIGTYPE_p_std__vectorT_unsigned_long_std__allocatorT_unsigned_long_t_t swig_types[103]
+#define SWIGTYPE_p_swig__SwigPyIterator swig_types[104]
+#define SWIGTYPE_p_unsigned_char swig_types[105]
+#define SWIGTYPE_p_unsigned_int swig_types[106]
+#define SWIGTYPE_p_unsigned_long_long swig_types[107]
+#define SWIGTYPE_p_unsigned_short swig_types[108]
+#define SWIGTYPE_p_value_type swig_types[109]
+static swig_type_info *swig_types[111];
+static swig_module_info swig_module = {swig_types, 110, 0, 0, 0, 0};
 #define SWIG_TypeQuery(name) SWIG_TypeQueryModule(&swig_module, &swig_module, name)
 #define SWIG_MangledTypeQuery(name) SWIG_MangledTypeQueryModule(&swig_module, &swig_module, name)
 
@@ -33794,6 +33795,29 @@ fail:
 }
 
 
+SWIGINTERN PyObject *_wrap_Beam_getPolarization(PyObject *SWIGUNUSEDPARM(self), PyObject *args) {
+  PyObject *resultobj = 0;
+  Beam *arg1 = (Beam *) 0 ;
+  void *argp1 = 0 ;
+  int res1 = 0 ;
+  PyObject *swig_obj[1] ;
+  SpinMatrix result;
+  
+  if (!args) SWIG_fail;
+  swig_obj[0] = args;
+  res1 = SWIG_ConvertPtr(swig_obj[0], &argp1,SWIGTYPE_p_Beam, 0 |  0 );
+  if (!SWIG_IsOK(res1)) {
+    SWIG_exception_fail(SWIG_ArgError(res1), "in method '" "Beam_getPolarization" "', argument " "1"" of type '" "Beam const *""'"); 
+  }
+  arg1 = reinterpret_cast< Beam * >(argp1);
+  result = ((Beam const *)arg1)->getPolarization();
+  resultobj = SWIG_NewPointerObj((new SpinMatrix(static_cast< const SpinMatrix& >(result))), SWIGTYPE_p_SpinMatrix, SWIG_POINTER_OWN |  0 );
+  return resultobj;
+fail:
+  return NULL;
+}
+
+
 SWIGINTERN PyObject *_wrap_Beam_setWavelength(PyObject *SWIGUNUSEDPARM(self), PyObject *args) {
   PyObject *resultobj = 0;
   Beam *arg1 = (Beam *) 0 ;
@@ -47686,6 +47710,13 @@ static PyMethodDef SwigMethods[] = {
 		"Check for limits, set the value if within limits. Throws if limits are violated. \n"
 		"\n"
 		""},
+	 { "Beam_getPolarization", _wrap_Beam_getPolarization, METH_O, "\n"
+		"Beam_getPolarization(Beam self) -> SpinMatrix\n"
+		"SpinMatrix Beam::getPolarization() const\n"
+		"\n"
+		"Returns the polarization density matrix (in spin basis along z-axis) \n"
+		"\n"
+		""},
 	 { "Beam_setWavelength", _wrap_Beam_setWavelength, METH_VARARGS, "\n"
 		"Beam_setWavelength(Beam self, double wavelength)\n"
 		"void Beam::setWavelength(double wavelength)\n"
@@ -49563,6 +49594,7 @@ static swig_type_info _swigt__p_RectangularPixel = {"_p_RectangularPixel", "Rect
 static swig_type_info _swigt__p_ResolutionFunction2DGaussian = {"_p_ResolutionFunction2DGaussian", "ResolutionFunction2DGaussian *", 0, 0, (void*)0, 0};
 static swig_type_info _swigt__p_SimulationResult = {"_p_SimulationResult", "SimulationResult *", 0, 0, (void*)0, 0};
 static swig_type_info _swigt__p_SphericalDetector = {"_p_SphericalDetector", "SphericalDetector *", 0, 0, (void*)0, 0};
+static swig_type_info _swigt__p_SpinMatrix = {"_p_SpinMatrix", "SpinMatrix *", 0, 0, (void*)0, 0};
 static swig_type_info _swigt__p_Vec3T_double_t = {"_p_Vec3T_double_t", "std::vector< Vec3< double > >::value_type *|R3 *|Vec3< double > *", 0, 0, (void*)0, 0};
 static swig_type_info _swigt__p_Vec3T_int_t = {"_p_Vec3T_int_t", "I3 *|Vec3< int > *", 0, 0, (void*)0, 0};
 static swig_type_info _swigt__p_Vec3T_std__complexT_double_t_t = {"_p_Vec3T_std__complexT_double_t_t", "Vec3< std::complex< double > > *|C3 *|std::vector< Vec3< std::complex< double > > >::value_type *", 0, 0, (void*)0, 0};
@@ -49674,6 +49706,7 @@ static swig_type_info *swig_type_initial[] = {
   &_swigt__p_ResolutionFunction2DGaussian,
   &_swigt__p_SimulationResult,
   &_swigt__p_SphericalDetector,
+  &_swigt__p_SpinMatrix,
   &_swigt__p_Vec3T_double_t,
   &_swigt__p_Vec3T_int_t,
   &_swigt__p_Vec3T_std__complexT_double_t_t,
@@ -49785,6 +49818,7 @@ static swig_cast_info _swigc__p_RectangularPixel[] = {  {&_swigt__p_RectangularP
 static swig_cast_info _swigc__p_ResolutionFunction2DGaussian[] = {  {&_swigt__p_ResolutionFunction2DGaussian, 0, 0, 0},{0, 0, 0, 0}};
 static swig_cast_info _swigc__p_SimulationResult[] = {  {&_swigt__p_SimulationResult, 0, 0, 0},{0, 0, 0, 0}};
 static swig_cast_info _swigc__p_SphericalDetector[] = {  {&_swigt__p_SphericalDetector, 0, 0, 0},{0, 0, 0, 0}};
+static swig_cast_info _swigc__p_SpinMatrix[] = {  {&_swigt__p_SpinMatrix, 0, 0, 0},{0, 0, 0, 0}};
 static swig_cast_info _swigc__p_Vec3T_double_t[] = {  {&_swigt__p_Vec3T_double_t, 0, 0, 0},{0, 0, 0, 0}};
 static swig_cast_info _swigc__p_Vec3T_int_t[] = {  {&_swigt__p_Vec3T_int_t, 0, 0, 0},{0, 0, 0, 0}};
 static swig_cast_info _swigc__p_Vec3T_std__complexT_double_t_t[] = {  {&_swigt__p_Vec3T_std__complexT_double_t_t, 0, 0, 0},{0, 0, 0, 0}};
@@ -49896,6 +49930,7 @@ static swig_cast_info *swig_cast_initial[] = {
   _swigc__p_ResolutionFunction2DGaussian,
   _swigc__p_SimulationResult,
   _swigc__p_SphericalDetector,
+  _swigc__p_SpinMatrix,
   _swigc__p_Vec3T_double_t,
   _swigc__p_Vec3T_int_t,
   _swigc__p_Vec3T_std__complexT_double_t_t,
diff --git a/auto/Wrap/libBornAgainSample.py b/auto/Wrap/libBornAgainSample.py
index ce3144af323f1937cb1ed7ab5fd6e1a346012262..57488771ace307432f3092504f5fe6b48eeeb0fc 100644
--- a/auto/Wrap/libBornAgainSample.py
+++ b/auto/Wrap/libBornAgainSample.py
@@ -2891,6 +2891,16 @@ class Material(object):
         """
         return _libBornAgainSample.Material_scalarSubtrSLD(self, wavevectors)
 
+    def polarizedSubtrSLD(self, wavevectors):
+        r"""
+        polarizedSubtrSLD(Material self, WavevectorInfo const & wavevectors) -> SpinMatrix
+        SpinMatrix Material::polarizedSubtrSLD(const WavevectorInfo &wavevectors) const
+
+        Returns (  $ \\pi/\\lambda^2 $ - sld) matrix with magnetization corrections. 
+
+        """
+        return _libBornAgainSample.Material_polarizedSubtrSLD(self, wavevectors)
+
     def rotatedMaterial(self, transform):
         r"""
         rotatedMaterial(Material self, RotMatrix transform) -> Material
@@ -3062,6 +3072,14 @@ class IFormFactor(ISampleNode):
         """
         return _libBornAgainSample.IFormFactor_theFF(self, wavevectors)
 
+    def thePolFF(self, wavevectors):
+        r"""
+        thePolFF(IFormFactor self, WavevectorInfo const & wavevectors) -> SpinMatrix
+        SpinMatrix IFormFactor::thePolFF(const WavevectorInfo &wavevectors) const
+
+        """
+        return _libBornAgainSample.IFormFactor_thePolFF(self, wavevectors)
+
     def shapeName(self):
         r"""
         shapeName(IFormFactor self) -> std::string
@@ -3131,6 +3149,16 @@ class IFormFactor(ISampleNode):
 
         """
         return _libBornAgainSample.IFormFactor_canSliceAnalytically(self, rot)
+
+    def formfactor_pol(self, q):
+        r"""
+        formfactor_pol(IFormFactor self, C3 q) -> SpinMatrix
+        SpinMatrix IFormFactor::formfactor_pol(C3 q) const
+
+        Returns scattering amplitude for complex scattering wavevector q=k_i-k_f in case of matrix interactions. Default implementation calls formfactor_at_bottom(q) and multiplies with the unit matrix. 
+
+        """
+        return _libBornAgainSample.IFormFactor_formfactor_pol(self, q)
     def __disown__(self):
         self.this.disown()
         _libBornAgainSample.disown_IFormFactor(self)
diff --git a/auto/Wrap/libBornAgainSample_wrap.cpp b/auto/Wrap/libBornAgainSample_wrap.cpp
index ba0e2ba77c3023f6d588f0f921fba1ecd96fbd5a..10d27264012ff74de7a0bdcdc32f4411e2aa6017 100644
--- a/auto/Wrap/libBornAgainSample_wrap.cpp
+++ b/auto/Wrap/libBornAgainSample_wrap.cpp
@@ -3200,71 +3200,72 @@ namespace Swig {
 #define SWIGTYPE_p_SawtoothRippleGauss swig_types[100]
 #define SWIGTYPE_p_SawtoothRippleLorentz swig_types[101]
 #define SWIGTYPE_p_SimpleSelectionRule swig_types[102]
-#define SWIGTYPE_p_SquareLattice2D swig_types[103]
-#define SWIGTYPE_p_TruncatedCube swig_types[104]
-#define SWIGTYPE_p_TruncatedSphere swig_types[105]
-#define SWIGTYPE_p_TruncatedSpheroid swig_types[106]
-#define SWIGTYPE_p_Vec3T_double_t swig_types[107]
-#define SWIGTYPE_p_Vec3T_int_t swig_types[108]
-#define SWIGTYPE_p_Vec3T_std__complexT_double_t_t swig_types[109]
-#define SWIGTYPE_p_WavevectorInfo swig_types[110]
-#define SWIGTYPE_p_allocator_type swig_types[111]
-#define SWIGTYPE_p_char swig_types[112]
-#define SWIGTYPE_p_difference_type swig_types[113]
-#define SWIGTYPE_p_first_type swig_types[114]
-#define SWIGTYPE_p_int swig_types[115]
-#define SWIGTYPE_p_key_type swig_types[116]
-#define SWIGTYPE_p_long_long swig_types[117]
-#define SWIGTYPE_p_mapped_type swig_types[118]
-#define SWIGTYPE_p_p_PyObject swig_types[119]
-#define SWIGTYPE_p_second_type swig_types[120]
-#define SWIGTYPE_p_short swig_types[121]
-#define SWIGTYPE_p_signed_char swig_types[122]
-#define SWIGTYPE_p_size_type swig_types[123]
-#define SWIGTYPE_p_std__allocatorT_INode_const_p_t swig_types[124]
-#define SWIGTYPE_p_std__allocatorT_INode_p_t swig_types[125]
-#define SWIGTYPE_p_std__allocatorT_Vec3T_double_t_t swig_types[126]
-#define SWIGTYPE_p_std__allocatorT_Vec3T_std__complexT_double_t_t_t swig_types[127]
-#define SWIGTYPE_p_std__allocatorT_double_t swig_types[128]
-#define SWIGTYPE_p_std__allocatorT_int_t swig_types[129]
-#define SWIGTYPE_p_std__allocatorT_std__complexT_double_t_t swig_types[130]
-#define SWIGTYPE_p_std__allocatorT_std__pairT_double_double_t_t swig_types[131]
-#define SWIGTYPE_p_std__allocatorT_std__pairT_std__string_const_double_t_t swig_types[132]
-#define SWIGTYPE_p_std__allocatorT_std__string_t swig_types[133]
-#define SWIGTYPE_p_std__allocatorT_std__vectorT_double_std__allocatorT_double_t_t_t swig_types[134]
-#define SWIGTYPE_p_std__allocatorT_std__vectorT_int_std__allocatorT_int_t_t_t swig_types[135]
-#define SWIGTYPE_p_std__allocatorT_unsigned_long_t swig_types[136]
-#define SWIGTYPE_p_std__arrayT_double_3_t swig_types[137]
-#define SWIGTYPE_p_std__arrayT_std__complexT_double_t_3_t swig_types[138]
-#define SWIGTYPE_p_std__complexT_double_t swig_types[139]
-#define SWIGTYPE_p_std__invalid_argument swig_types[140]
-#define SWIGTYPE_p_std__lessT_std__string_t swig_types[141]
-#define SWIGTYPE_p_std__mapT_std__string_double_std__lessT_std__string_t_std__allocatorT_std__pairT_std__string_const_double_t_t_t swig_types[142]
-#define SWIGTYPE_p_std__pairT_double_double_t swig_types[143]
-#define SWIGTYPE_p_std__vectorT_INode_const_p_std__allocatorT_INode_const_p_t_t swig_types[144]
-#define SWIGTYPE_p_std__vectorT_INode_p_std__allocatorT_INode_p_t_t swig_types[145]
-#define SWIGTYPE_p_std__vectorT_IParticle_const_p_std__allocatorT_IParticle_const_p_t_t swig_types[146]
-#define SWIGTYPE_p_std__vectorT_Material_const_p_std__allocatorT_Material_const_p_t_t swig_types[147]
-#define SWIGTYPE_p_std__vectorT_ParaMeta_std__allocatorT_ParaMeta_t_t swig_types[148]
-#define SWIGTYPE_p_std__vectorT_ParticleLayout_const_p_std__allocatorT_ParticleLayout_const_p_t_t swig_types[149]
-#define SWIGTYPE_p_std__vectorT_Vec3T_double_t_std__allocatorT_Vec3T_double_t_t_t swig_types[150]
-#define SWIGTYPE_p_std__vectorT_Vec3T_std__complexT_double_t_t_std__allocatorT_Vec3T_std__complexT_double_t_t_t_t swig_types[151]
-#define SWIGTYPE_p_std__vectorT_double_std__allocatorT_double_t_t swig_types[152]
-#define SWIGTYPE_p_std__vectorT_int_std__allocatorT_int_t_t swig_types[153]
-#define SWIGTYPE_p_std__vectorT_std__complexT_double_t_std__allocatorT_std__complexT_double_t_t_t swig_types[154]
-#define SWIGTYPE_p_std__vectorT_std__pairT_double_double_t_std__allocatorT_std__pairT_double_double_t_t_t swig_types[155]
-#define SWIGTYPE_p_std__vectorT_std__string_std__allocatorT_std__string_t_t swig_types[156]
-#define SWIGTYPE_p_std__vectorT_std__vectorT_double_std__allocatorT_double_t_t_std__allocatorT_std__vectorT_double_std__allocatorT_double_t_t_t_t swig_types[157]
-#define SWIGTYPE_p_std__vectorT_std__vectorT_int_std__allocatorT_int_t_t_std__allocatorT_std__vectorT_int_std__allocatorT_int_t_t_t_t swig_types[158]
-#define SWIGTYPE_p_std__vectorT_unsigned_long_std__allocatorT_unsigned_long_t_t swig_types[159]
-#define SWIGTYPE_p_swig__SwigPyIterator swig_types[160]
-#define SWIGTYPE_p_unsigned_char swig_types[161]
-#define SWIGTYPE_p_unsigned_int swig_types[162]
-#define SWIGTYPE_p_unsigned_long_long swig_types[163]
-#define SWIGTYPE_p_unsigned_short swig_types[164]
-#define SWIGTYPE_p_value_type swig_types[165]
-static swig_type_info *swig_types[167];
-static swig_module_info swig_module = {swig_types, 166, 0, 0, 0, 0};
+#define SWIGTYPE_p_SpinMatrix swig_types[103]
+#define SWIGTYPE_p_SquareLattice2D swig_types[104]
+#define SWIGTYPE_p_TruncatedCube swig_types[105]
+#define SWIGTYPE_p_TruncatedSphere swig_types[106]
+#define SWIGTYPE_p_TruncatedSpheroid swig_types[107]
+#define SWIGTYPE_p_Vec3T_double_t swig_types[108]
+#define SWIGTYPE_p_Vec3T_int_t swig_types[109]
+#define SWIGTYPE_p_Vec3T_std__complexT_double_t_t swig_types[110]
+#define SWIGTYPE_p_WavevectorInfo swig_types[111]
+#define SWIGTYPE_p_allocator_type swig_types[112]
+#define SWIGTYPE_p_char swig_types[113]
+#define SWIGTYPE_p_difference_type swig_types[114]
+#define SWIGTYPE_p_first_type swig_types[115]
+#define SWIGTYPE_p_int swig_types[116]
+#define SWIGTYPE_p_key_type swig_types[117]
+#define SWIGTYPE_p_long_long swig_types[118]
+#define SWIGTYPE_p_mapped_type swig_types[119]
+#define SWIGTYPE_p_p_PyObject swig_types[120]
+#define SWIGTYPE_p_second_type swig_types[121]
+#define SWIGTYPE_p_short swig_types[122]
+#define SWIGTYPE_p_signed_char swig_types[123]
+#define SWIGTYPE_p_size_type swig_types[124]
+#define SWIGTYPE_p_std__allocatorT_INode_const_p_t swig_types[125]
+#define SWIGTYPE_p_std__allocatorT_INode_p_t swig_types[126]
+#define SWIGTYPE_p_std__allocatorT_Vec3T_double_t_t swig_types[127]
+#define SWIGTYPE_p_std__allocatorT_Vec3T_std__complexT_double_t_t_t swig_types[128]
+#define SWIGTYPE_p_std__allocatorT_double_t swig_types[129]
+#define SWIGTYPE_p_std__allocatorT_int_t swig_types[130]
+#define SWIGTYPE_p_std__allocatorT_std__complexT_double_t_t swig_types[131]
+#define SWIGTYPE_p_std__allocatorT_std__pairT_double_double_t_t swig_types[132]
+#define SWIGTYPE_p_std__allocatorT_std__pairT_std__string_const_double_t_t swig_types[133]
+#define SWIGTYPE_p_std__allocatorT_std__string_t swig_types[134]
+#define SWIGTYPE_p_std__allocatorT_std__vectorT_double_std__allocatorT_double_t_t_t swig_types[135]
+#define SWIGTYPE_p_std__allocatorT_std__vectorT_int_std__allocatorT_int_t_t_t swig_types[136]
+#define SWIGTYPE_p_std__allocatorT_unsigned_long_t swig_types[137]
+#define SWIGTYPE_p_std__arrayT_double_3_t swig_types[138]
+#define SWIGTYPE_p_std__arrayT_std__complexT_double_t_3_t swig_types[139]
+#define SWIGTYPE_p_std__complexT_double_t swig_types[140]
+#define SWIGTYPE_p_std__invalid_argument swig_types[141]
+#define SWIGTYPE_p_std__lessT_std__string_t swig_types[142]
+#define SWIGTYPE_p_std__mapT_std__string_double_std__lessT_std__string_t_std__allocatorT_std__pairT_std__string_const_double_t_t_t swig_types[143]
+#define SWIGTYPE_p_std__pairT_double_double_t swig_types[144]
+#define SWIGTYPE_p_std__vectorT_INode_const_p_std__allocatorT_INode_const_p_t_t swig_types[145]
+#define SWIGTYPE_p_std__vectorT_INode_p_std__allocatorT_INode_p_t_t swig_types[146]
+#define SWIGTYPE_p_std__vectorT_IParticle_const_p_std__allocatorT_IParticle_const_p_t_t swig_types[147]
+#define SWIGTYPE_p_std__vectorT_Material_const_p_std__allocatorT_Material_const_p_t_t swig_types[148]
+#define SWIGTYPE_p_std__vectorT_ParaMeta_std__allocatorT_ParaMeta_t_t swig_types[149]
+#define SWIGTYPE_p_std__vectorT_ParticleLayout_const_p_std__allocatorT_ParticleLayout_const_p_t_t swig_types[150]
+#define SWIGTYPE_p_std__vectorT_Vec3T_double_t_std__allocatorT_Vec3T_double_t_t_t swig_types[151]
+#define SWIGTYPE_p_std__vectorT_Vec3T_std__complexT_double_t_t_std__allocatorT_Vec3T_std__complexT_double_t_t_t_t swig_types[152]
+#define SWIGTYPE_p_std__vectorT_double_std__allocatorT_double_t_t swig_types[153]
+#define SWIGTYPE_p_std__vectorT_int_std__allocatorT_int_t_t swig_types[154]
+#define SWIGTYPE_p_std__vectorT_std__complexT_double_t_std__allocatorT_std__complexT_double_t_t_t swig_types[155]
+#define SWIGTYPE_p_std__vectorT_std__pairT_double_double_t_std__allocatorT_std__pairT_double_double_t_t_t swig_types[156]
+#define SWIGTYPE_p_std__vectorT_std__string_std__allocatorT_std__string_t_t swig_types[157]
+#define SWIGTYPE_p_std__vectorT_std__vectorT_double_std__allocatorT_double_t_t_std__allocatorT_std__vectorT_double_std__allocatorT_double_t_t_t_t swig_types[158]
+#define SWIGTYPE_p_std__vectorT_std__vectorT_int_std__allocatorT_int_t_t_std__allocatorT_std__vectorT_int_std__allocatorT_int_t_t_t_t swig_types[159]
+#define SWIGTYPE_p_std__vectorT_unsigned_long_std__allocatorT_unsigned_long_t_t swig_types[160]
+#define SWIGTYPE_p_swig__SwigPyIterator swig_types[161]
+#define SWIGTYPE_p_unsigned_char swig_types[162]
+#define SWIGTYPE_p_unsigned_int swig_types[163]
+#define SWIGTYPE_p_unsigned_long_long swig_types[164]
+#define SWIGTYPE_p_unsigned_short swig_types[165]
+#define SWIGTYPE_p_value_type swig_types[166]
+static swig_type_info *swig_types[168];
+static swig_module_info swig_module = {swig_types, 167, 0, 0, 0, 0};
 #define SWIG_TypeQuery(name) SWIG_TypeQueryModule(&swig_module, &swig_module, name)
 #define SWIG_MangledTypeQuery(name) SWIG_MangledTypeQueryModule(&swig_module, &swig_module, name)
 
@@ -7974,6 +7975,41 @@ complex_t SwigDirector_IFormFactor::theFF(WavevectorInfo const &wavevectors) con
 }
 
 
+SpinMatrix SwigDirector_IFormFactor::thePolFF(WavevectorInfo const &wavevectors) const {
+  void *swig_argp ;
+  int swig_res = 0 ;
+  
+  SpinMatrix c_result;
+  swig::SwigVar_PyObject obj0;
+  obj0 = SWIG_NewPointerObj(SWIG_as_voidptr(&wavevectors), SWIGTYPE_p_WavevectorInfo,  0 );
+  if (!swig_get_self()) {
+    Swig::DirectorException::raise("'self' uninitialized, maybe you forgot to call IFormFactor.__init__.");
+  }
+#if defined(SWIG_PYTHON_DIRECTOR_VTABLE)
+  const size_t swig_method_index = 7;
+  const char *const swig_method_name = "thePolFF";
+  PyObject *method = swig_get_method(swig_method_index, swig_method_name);
+  swig::SwigVar_PyObject result = PyObject_CallFunctionObjArgs(method ,(PyObject *)obj0, NULL);
+#else
+  swig::SwigVar_PyObject swig_method_name = SWIG_Python_str_FromChar("thePolFF");
+  swig::SwigVar_PyObject result = PyObject_CallMethodObjArgs(swig_get_self(), (PyObject *) swig_method_name ,(PyObject *)obj0, NULL);
+#endif
+  if (!result) {
+    PyObject *error = PyErr_Occurred();
+    if (error) {
+      Swig::DirectorMethodException::raise("Error detected when calling 'IFormFactor.thePolFF'");
+    }
+  }
+  swig_res = SWIG_ConvertPtr(result,&swig_argp,SWIGTYPE_p_SpinMatrix,  0  | 0);
+  if (!SWIG_IsOK(swig_res)) {
+    Swig::DirectorTypeMismatchException::raise(SWIG_ErrorType(SWIG_ArgError(swig_res)), "in output value of type '""SpinMatrix""'");
+  }
+  c_result = *(reinterpret_cast< SpinMatrix * >(swig_argp));
+  if (SWIG_IsNewObj(swig_res)) delete reinterpret_cast< SpinMatrix * >(swig_argp);
+  return (SpinMatrix) c_result;
+}
+
+
 double SwigDirector_IFormFactor::bottomZ(IRotation const *rotation) const {
   double c_result = SwigValueInit< double >() ;
   
@@ -7983,7 +8019,7 @@ double SwigDirector_IFormFactor::bottomZ(IRotation const *rotation) const {
     Swig::DirectorException::raise("'self' uninitialized, maybe you forgot to call IFormFactor.__init__.");
   }
 #if defined(SWIG_PYTHON_DIRECTOR_VTABLE)
-  const size_t swig_method_index = 7;
+  const size_t swig_method_index = 8;
   const char *const swig_method_name = "bottomZ";
   PyObject *method = swig_get_method(swig_method_index, swig_method_name);
   swig::SwigVar_PyObject result = PyObject_CallFunctionObjArgs(method ,(PyObject *)obj0, NULL);
@@ -8016,7 +8052,7 @@ double SwigDirector_IFormFactor::topZ(IRotation const *rotation) const {
     Swig::DirectorException::raise("'self' uninitialized, maybe you forgot to call IFormFactor.__init__.");
   }
 #if defined(SWIG_PYTHON_DIRECTOR_VTABLE)
-  const size_t swig_method_index = 8;
+  const size_t swig_method_index = 9;
   const char *const swig_method_name = "topZ";
   PyObject *method = swig_get_method(swig_method_index, swig_method_name);
   swig::SwigVar_PyObject result = PyObject_CallFunctionObjArgs(method ,(PyObject *)obj0, NULL);
@@ -8047,7 +8083,7 @@ double SwigDirector_IFormFactor::volume() const {
     Swig::DirectorException::raise("'self' uninitialized, maybe you forgot to call IFormFactor.__init__.");
   }
 #if defined(SWIG_PYTHON_DIRECTOR_VTABLE)
-  const size_t swig_method_index = 9;
+  const size_t swig_method_index = 10;
   const char *const swig_method_name = "volume";
   PyObject *method = swig_get_method(swig_method_index, swig_method_name);
   swig::SwigVar_PyObject args = PyTuple_New(0);
@@ -8079,7 +8115,7 @@ double SwigDirector_IFormFactor::radialExtension() const {
     Swig::DirectorException::raise("'self' uninitialized, maybe you forgot to call IFormFactor.__init__.");
   }
 #if defined(SWIG_PYTHON_DIRECTOR_VTABLE)
-  const size_t swig_method_index = 10;
+  const size_t swig_method_index = 11;
   const char *const swig_method_name = "radialExtension";
   PyObject *method = swig_get_method(swig_method_index, swig_method_name);
   swig::SwigVar_PyObject args = PyTuple_New(0);
@@ -8112,7 +8148,7 @@ complex_t SwigDirector_IFormFactor::formfactor_at_bottom(C3 q) const {
     Swig::DirectorException::raise("'self' uninitialized, maybe you forgot to call IFormFactor.__init__.");
   }
 #if defined(SWIG_PYTHON_DIRECTOR_VTABLE)
-  const size_t swig_method_index = 11;
+  const size_t swig_method_index = 12;
   const char *const swig_method_name = "formfactor_at_bottom";
   PyObject *method = swig_get_method(swig_method_index, swig_method_name);
   swig::SwigVar_PyObject result = PyObject_CallFunctionObjArgs(method ,(PyObject *)obj0, NULL);
@@ -8142,7 +8178,7 @@ std::string SwigDirector_IFormFactor::pythonConstructor() const {
     Swig::DirectorException::raise("'self' uninitialized, maybe you forgot to call IFormFactor.__init__.");
   }
 #if defined(SWIG_PYTHON_DIRECTOR_VTABLE)
-  const size_t swig_method_index = 12;
+  const size_t swig_method_index = 13;
   const char *const swig_method_name = "pythonConstructor";
   PyObject *method = swig_get_method(swig_method_index, swig_method_name);
   swig::SwigVar_PyObject args = PyTuple_New(0);
@@ -8177,7 +8213,7 @@ bool SwigDirector_IFormFactor::canSliceAnalytically(IRotation const *rot) const
     Swig::DirectorException::raise("'self' uninitialized, maybe you forgot to call IFormFactor.__init__.");
   }
 #if defined(SWIG_PYTHON_DIRECTOR_VTABLE)
-  const size_t swig_method_index = 13;
+  const size_t swig_method_index = 14;
   const char *const swig_method_name = "canSliceAnalytically";
   PyObject *method = swig_get_method(swig_method_index, swig_method_name);
   swig::SwigVar_PyObject result = PyObject_CallFunctionObjArgs(method ,(PyObject *)obj0, NULL);
@@ -8201,6 +8237,41 @@ bool SwigDirector_IFormFactor::canSliceAnalytically(IRotation const *rot) const
 }
 
 
+SpinMatrix SwigDirector_IFormFactor::formfactor_pol(C3 q) const {
+  void *swig_argp ;
+  int swig_res = 0 ;
+  
+  SpinMatrix c_result;
+  swig::SwigVar_PyObject obj0;
+  obj0 = SWIG_NewPointerObj(SWIG_as_voidptr(new C3((const C3 &)q)), SWIGTYPE_p_Vec3T_std__complexT_double_t_t, SWIG_POINTER_OWN |  0 );
+  if (!swig_get_self()) {
+    Swig::DirectorException::raise("'self' uninitialized, maybe you forgot to call IFormFactor.__init__.");
+  }
+#if defined(SWIG_PYTHON_DIRECTOR_VTABLE)
+  const size_t swig_method_index = 15;
+  const char *const swig_method_name = "formfactor_pol";
+  PyObject *method = swig_get_method(swig_method_index, swig_method_name);
+  swig::SwigVar_PyObject result = PyObject_CallFunctionObjArgs(method ,(PyObject *)obj0, NULL);
+#else
+  swig::SwigVar_PyObject swig_method_name = SWIG_Python_str_FromChar("formfactor_pol");
+  swig::SwigVar_PyObject result = PyObject_CallMethodObjArgs(swig_get_self(), (PyObject *) swig_method_name ,(PyObject *)obj0, NULL);
+#endif
+  if (!result) {
+    PyObject *error = PyErr_Occurred();
+    if (error) {
+      Swig::DirectorMethodException::raise("Error detected when calling 'IFormFactor.formfactor_pol'");
+    }
+  }
+  swig_res = SWIG_ConvertPtr(result,&swig_argp,SWIGTYPE_p_SpinMatrix,  0  | 0);
+  if (!SWIG_IsOK(swig_res)) {
+    Swig::DirectorTypeMismatchException::raise(SWIG_ErrorType(SWIG_ArgError(swig_res)), "in output value of type '""SpinMatrix""'");
+  }
+  c_result = *(reinterpret_cast< SpinMatrix * >(swig_argp));
+  if (SWIG_IsNewObj(swig_res)) delete reinterpret_cast< SpinMatrix * >(swig_argp);
+  return (SpinMatrix) c_result;
+}
+
+
 #ifdef __cplusplus
 extern "C" {
 #endif
@@ -36278,6 +36349,39 @@ fail:
 }
 
 
+SWIGINTERN PyObject *_wrap_Material_polarizedSubtrSLD(PyObject *SWIGUNUSEDPARM(self), PyObject *args) {
+  PyObject *resultobj = 0;
+  Material *arg1 = (Material *) 0 ;
+  WavevectorInfo *arg2 = 0 ;
+  void *argp1 = 0 ;
+  int res1 = 0 ;
+  void *argp2 = 0 ;
+  int res2 = 0 ;
+  PyObject *swig_obj[2] ;
+  SpinMatrix result;
+  
+  if (!SWIG_Python_UnpackTuple(args, "Material_polarizedSubtrSLD", 2, 2, swig_obj)) SWIG_fail;
+  res1 = SWIG_ConvertPtr(swig_obj[0], &argp1,SWIGTYPE_p_Material, 0 |  0 );
+  if (!SWIG_IsOK(res1)) {
+    SWIG_exception_fail(SWIG_ArgError(res1), "in method '" "Material_polarizedSubtrSLD" "', argument " "1"" of type '" "Material const *""'"); 
+  }
+  arg1 = reinterpret_cast< Material * >(argp1);
+  res2 = SWIG_ConvertPtr(swig_obj[1], &argp2, SWIGTYPE_p_WavevectorInfo,  0  | 0);
+  if (!SWIG_IsOK(res2)) {
+    SWIG_exception_fail(SWIG_ArgError(res2), "in method '" "Material_polarizedSubtrSLD" "', argument " "2"" of type '" "WavevectorInfo const &""'"); 
+  }
+  if (!argp2) {
+    SWIG_exception_fail(SWIG_ValueError, "invalid null reference " "in method '" "Material_polarizedSubtrSLD" "', argument " "2"" of type '" "WavevectorInfo const &""'"); 
+  }
+  arg2 = reinterpret_cast< WavevectorInfo * >(argp2);
+  result = ((Material const *)arg1)->polarizedSubtrSLD((WavevectorInfo const &)*arg2);
+  resultobj = SWIG_NewPointerObj((new SpinMatrix(static_cast< const SpinMatrix& >(result))), SWIGTYPE_p_SpinMatrix, SWIG_POINTER_OWN |  0 );
+  return resultobj;
+fail:
+  return NULL;
+}
+
+
 SWIGINTERN PyObject *_wrap_Material_rotatedMaterial(PyObject *SWIGUNUSEDPARM(self), PyObject *args) {
   PyObject *resultobj = 0;
   Material *arg1 = (Material *) 0 ;
@@ -37320,6 +37424,51 @@ fail:
 }
 
 
+SWIGINTERN PyObject *_wrap_IFormFactor_thePolFF(PyObject *SWIGUNUSEDPARM(self), PyObject *args) {
+  PyObject *resultobj = 0;
+  IFormFactor *arg1 = (IFormFactor *) 0 ;
+  WavevectorInfo *arg2 = 0 ;
+  void *argp1 = 0 ;
+  int res1 = 0 ;
+  void *argp2 = 0 ;
+  int res2 = 0 ;
+  PyObject *swig_obj[2] ;
+  Swig::Director *director = 0;
+  bool upcall = false;
+  SpinMatrix result;
+  
+  if (!SWIG_Python_UnpackTuple(args, "IFormFactor_thePolFF", 2, 2, swig_obj)) SWIG_fail;
+  res1 = SWIG_ConvertPtr(swig_obj[0], &argp1,SWIGTYPE_p_IFormFactor, 0 |  0 );
+  if (!SWIG_IsOK(res1)) {
+    SWIG_exception_fail(SWIG_ArgError(res1), "in method '" "IFormFactor_thePolFF" "', argument " "1"" of type '" "IFormFactor const *""'"); 
+  }
+  arg1 = reinterpret_cast< IFormFactor * >(argp1);
+  res2 = SWIG_ConvertPtr(swig_obj[1], &argp2, SWIGTYPE_p_WavevectorInfo,  0  | 0);
+  if (!SWIG_IsOK(res2)) {
+    SWIG_exception_fail(SWIG_ArgError(res2), "in method '" "IFormFactor_thePolFF" "', argument " "2"" of type '" "WavevectorInfo const &""'"); 
+  }
+  if (!argp2) {
+    SWIG_exception_fail(SWIG_ValueError, "invalid null reference " "in method '" "IFormFactor_thePolFF" "', argument " "2"" of type '" "WavevectorInfo const &""'"); 
+  }
+  arg2 = reinterpret_cast< WavevectorInfo * >(argp2);
+  director = SWIG_DIRECTOR_CAST(arg1);
+  upcall = (director && (director->swig_get_self()==swig_obj[0]));
+  try {
+    if (upcall) {
+      result = ((IFormFactor const *)arg1)->IFormFactor::thePolFF((WavevectorInfo const &)*arg2);
+    } else {
+      result = ((IFormFactor const *)arg1)->thePolFF((WavevectorInfo const &)*arg2);
+    }
+  } catch (Swig::DirectorException&) {
+    SWIG_fail;
+  }
+  resultobj = SWIG_NewPointerObj((new SpinMatrix(static_cast< const SpinMatrix& >(result))), SWIGTYPE_p_SpinMatrix, SWIG_POINTER_OWN |  0 );
+  return resultobj;
+fail:
+  return NULL;
+}
+
+
 SWIGINTERN PyObject *_wrap_IFormFactor_shapeName(PyObject *SWIGUNUSEDPARM(self), PyObject *args) {
   PyObject *resultobj = 0;
   IFormFactor *arg1 = (IFormFactor *) 0 ;
@@ -37624,6 +37773,56 @@ fail:
 }
 
 
+SWIGINTERN PyObject *_wrap_IFormFactor_formfactor_pol(PyObject *SWIGUNUSEDPARM(self), PyObject *args) {
+  PyObject *resultobj = 0;
+  IFormFactor *arg1 = (IFormFactor *) 0 ;
+  C3 arg2 ;
+  void *argp1 = 0 ;
+  int res1 = 0 ;
+  void *argp2 ;
+  int res2 = 0 ;
+  PyObject *swig_obj[2] ;
+  Swig::Director *director = 0;
+  bool upcall = false;
+  SpinMatrix result;
+  
+  if (!SWIG_Python_UnpackTuple(args, "IFormFactor_formfactor_pol", 2, 2, swig_obj)) SWIG_fail;
+  res1 = SWIG_ConvertPtr(swig_obj[0], &argp1,SWIGTYPE_p_IFormFactor, 0 |  0 );
+  if (!SWIG_IsOK(res1)) {
+    SWIG_exception_fail(SWIG_ArgError(res1), "in method '" "IFormFactor_formfactor_pol" "', argument " "1"" of type '" "IFormFactor const *""'"); 
+  }
+  arg1 = reinterpret_cast< IFormFactor * >(argp1);
+  {
+    res2 = SWIG_ConvertPtr(swig_obj[1], &argp2, SWIGTYPE_p_Vec3T_std__complexT_double_t_t,  0  | 0);
+    if (!SWIG_IsOK(res2)) {
+      SWIG_exception_fail(SWIG_ArgError(res2), "in method '" "IFormFactor_formfactor_pol" "', argument " "2"" of type '" "C3""'"); 
+    }  
+    if (!argp2) {
+      SWIG_exception_fail(SWIG_ValueError, "invalid null reference " "in method '" "IFormFactor_formfactor_pol" "', argument " "2"" of type '" "C3""'");
+    } else {
+      C3 * temp = reinterpret_cast< C3 * >(argp2);
+      arg2 = *temp;
+      if (SWIG_IsNewObj(res2)) delete temp;
+    }
+  }
+  director = SWIG_DIRECTOR_CAST(arg1);
+  upcall = (director && (director->swig_get_self()==swig_obj[0]));
+  try {
+    if (upcall) {
+      result = ((IFormFactor const *)arg1)->IFormFactor::formfactor_pol(arg2);
+    } else {
+      result = ((IFormFactor const *)arg1)->formfactor_pol(arg2);
+    }
+  } catch (Swig::DirectorException&) {
+    SWIG_fail;
+  }
+  resultobj = SWIG_NewPointerObj((new SpinMatrix(static_cast< const SpinMatrix& >(result))), SWIGTYPE_p_SpinMatrix, SWIG_POINTER_OWN |  0 );
+  return resultobj;
+fail:
+  return NULL;
+}
+
+
 SWIGINTERN PyObject *_wrap_disown_IFormFactor(PyObject *SWIGUNUSEDPARM(self), PyObject *args) {
   PyObject *resultobj = 0;
   IFormFactor *arg1 = (IFormFactor *) 0 ;
@@ -65282,6 +65481,13 @@ static PyMethodDef SwigMethods[] = {
 		"Returns (  $ \\\\pi/\\\\lambda^2 $ - sld), sld (in  $nm^{-2}$) being the scattering length density \n"
 		"\n"
 		""},
+	 { "Material_polarizedSubtrSLD", _wrap_Material_polarizedSubtrSLD, METH_VARARGS, "\n"
+		"Material_polarizedSubtrSLD(Material self, WavevectorInfo const & wavevectors) -> SpinMatrix\n"
+		"SpinMatrix Material::polarizedSubtrSLD(const WavevectorInfo &wavevectors) const\n"
+		"\n"
+		"Returns (  $ \\\\pi/\\\\lambda^2 $ - sld) matrix with magnetization corrections. \n"
+		"\n"
+		""},
 	 { "Material_rotatedMaterial", _wrap_Material_rotatedMaterial, METH_VARARGS, "\n"
 		"Material_rotatedMaterial(Material self, RotMatrix transform) -> Material\n"
 		"Material Material::rotatedMaterial(const RotMatrix &transform) const\n"
@@ -65383,6 +65589,11 @@ static PyMethodDef SwigMethods[] = {
 		"complex_t IFormFactor::theFF(const WavevectorInfo &wavevectors) const\n"
 		"\n"
 		""},
+	 { "IFormFactor_thePolFF", _wrap_IFormFactor_thePolFF, METH_VARARGS, "\n"
+		"IFormFactor_thePolFF(IFormFactor self, WavevectorInfo const & wavevectors) -> SpinMatrix\n"
+		"SpinMatrix IFormFactor::thePolFF(const WavevectorInfo &wavevectors) const\n"
+		"\n"
+		""},
 	 { "IFormFactor_shapeName", _wrap_IFormFactor_shapeName, METH_O, "\n"
 		"IFormFactor_shapeName(IFormFactor self) -> std::string\n"
 		"std::string IFormFactor::shapeName() const\n"
@@ -65429,6 +65640,13 @@ static PyMethodDef SwigMethods[] = {
 		"Default implementation only allows rotations along z-axis. \n"
 		"\n"
 		""},
+	 { "IFormFactor_formfactor_pol", _wrap_IFormFactor_formfactor_pol, METH_VARARGS, "\n"
+		"IFormFactor_formfactor_pol(IFormFactor self, C3 q) -> SpinMatrix\n"
+		"SpinMatrix IFormFactor::formfactor_pol(C3 q) const\n"
+		"\n"
+		"Returns scattering amplitude for complex scattering wavevector q=k_i-k_f in case of matrix interactions. Default implementation calls formfactor_at_bottom(q) and multiplies with the unit matrix. \n"
+		"\n"
+		""},
 	 { "disown_IFormFactor", _wrap_disown_IFormFactor, METH_O, NULL},
 	 { "IFormFactor_swigregister", IFormFactor_swigregister, METH_O, NULL},
 	 { "IFormFactor_swiginit", IFormFactor_swiginit, METH_VARARGS, NULL},
@@ -71119,6 +71337,7 @@ static swig_type_info _swigt__p_SawtoothRippleBox = {"_p_SawtoothRippleBox", "Sa
 static swig_type_info _swigt__p_SawtoothRippleGauss = {"_p_SawtoothRippleGauss", "SawtoothRippleGauss *", 0, 0, (void*)0, 0};
 static swig_type_info _swigt__p_SawtoothRippleLorentz = {"_p_SawtoothRippleLorentz", "SawtoothRippleLorentz *", 0, 0, (void*)0, 0};
 static swig_type_info _swigt__p_SimpleSelectionRule = {"_p_SimpleSelectionRule", "SimpleSelectionRule *", 0, 0, (void*)0, 0};
+static swig_type_info _swigt__p_SpinMatrix = {"_p_SpinMatrix", "SpinMatrix *", 0, 0, (void*)0, 0};
 static swig_type_info _swigt__p_SquareLattice2D = {"_p_SquareLattice2D", "SquareLattice2D *", 0, 0, (void*)0, 0};
 static swig_type_info _swigt__p_TruncatedCube = {"_p_TruncatedCube", "TruncatedCube *", 0, 0, (void*)0, 0};
 static swig_type_info _swigt__p_TruncatedSphere = {"_p_TruncatedSphere", "TruncatedSphere *", 0, 0, (void*)0, 0};
@@ -71287,6 +71506,7 @@ static swig_type_info *swig_type_initial[] = {
   &_swigt__p_SawtoothRippleGauss,
   &_swigt__p_SawtoothRippleLorentz,
   &_swigt__p_SimpleSelectionRule,
+  &_swigt__p_SpinMatrix,
   &_swigt__p_SquareLattice2D,
   &_swigt__p_TruncatedCube,
   &_swigt__p_TruncatedSphere,
@@ -71455,6 +71675,7 @@ static swig_cast_info _swigc__p_SawtoothRippleBox[] = {  {&_swigt__p_SawtoothRip
 static swig_cast_info _swigc__p_SawtoothRippleGauss[] = {  {&_swigt__p_SawtoothRippleGauss, 0, 0, 0},{0, 0, 0, 0}};
 static swig_cast_info _swigc__p_SawtoothRippleLorentz[] = {  {&_swigt__p_SawtoothRippleLorentz, 0, 0, 0},{0, 0, 0, 0}};
 static swig_cast_info _swigc__p_SimpleSelectionRule[] = {  {&_swigt__p_SimpleSelectionRule, 0, 0, 0},{0, 0, 0, 0}};
+static swig_cast_info _swigc__p_SpinMatrix[] = {  {&_swigt__p_SpinMatrix, 0, 0, 0},{0, 0, 0, 0}};
 static swig_cast_info _swigc__p_SquareLattice2D[] = {  {&_swigt__p_SquareLattice2D, 0, 0, 0},{0, 0, 0, 0}};
 static swig_cast_info _swigc__p_TruncatedCube[] = {  {&_swigt__p_TruncatedCube, 0, 0, 0},{0, 0, 0, 0}};
 static swig_cast_info _swigc__p_TruncatedSphere[] = {  {&_swigt__p_TruncatedSphere, 0, 0, 0},{0, 0, 0, 0}};
@@ -71623,6 +71844,7 @@ static swig_cast_info *swig_cast_initial[] = {
   _swigc__p_SawtoothRippleGauss,
   _swigc__p_SawtoothRippleLorentz,
   _swigc__p_SimpleSelectionRule,
+  _swigc__p_SpinMatrix,
   _swigc__p_SquareLattice2D,
   _swigc__p_TruncatedCube,
   _swigc__p_TruncatedSphere,
diff --git a/auto/Wrap/libBornAgainSample_wrap.h b/auto/Wrap/libBornAgainSample_wrap.h
index 59f45853b19353fad632762b0b75b7b0f94586f2..d0cb8f32afe3bacb1a91475b8a600207cb91c579 100644
--- a/auto/Wrap/libBornAgainSample_wrap.h
+++ b/auto/Wrap/libBornAgainSample_wrap.h
@@ -76,6 +76,7 @@ public:
     virtual std::vector< INode const *, std::allocator< INode const * > > nodeChildren() const;
     virtual Material const *material() const;
     virtual complex_t theFF(WavevectorInfo const &wavevectors) const;
+    virtual SpinMatrix thePolFF(WavevectorInfo const &wavevectors) const;
     virtual double bottomZ(IRotation const *rotation) const;
     virtual double topZ(IRotation const *rotation) const;
     virtual double volume() const;
@@ -83,6 +84,7 @@ public:
     virtual complex_t formfactor_at_bottom(C3 q) const;
     virtual std::string pythonConstructor() const;
     virtual bool canSliceAnalytically(IRotation const *rot) const;
+    virtual SpinMatrix formfactor_pol(C3 q) const;
 
 /* Internal director utilities */
 public:
@@ -113,7 +115,7 @@ private:
       return method;
     }
 private:
-    mutable swig::SwigVar_PyObject vtable[14];
+    mutable swig::SwigVar_PyObject vtable[16];
 #endif
 
 };