diff --git a/Core/Instrument/Beam.cpp b/Core/Instrument/Beam.cpp
index a13b7ca05a65683f8ed15e545f3e704b93dcc2a8..a906568f0cd8fd46e83f05e0a91c02e92df1cae7 100644
--- a/Core/Instrument/Beam.cpp
+++ b/Core/Instrument/Beam.cpp
@@ -82,6 +82,14 @@ void Beam::setPolarization(const kvector_t bloch_vector)
     m_polarization = calculatePolarization(bloch_vector);
 }
 
+kvector_t Beam::getBlochVector() const
+{
+    double x = 2.0*m_polarization(0, 1).real();
+    double y = 2.0*m_polarization(1, 0).imag();
+    double z = (m_polarization(0, 0) - m_polarization(1, 1)).real();
+    return { x, y, z };
+}
+
 Eigen::Matrix2cd Beam::calculatePolarization(const kvector_t bloch_vector) const
 {
     Eigen::Matrix2cd result;
diff --git a/Core/Instrument/Beam.h b/Core/Instrument/Beam.h
index da0f57a68156f3f4733c4bcc9e55b35b48607696..a52d5bbb36926f8679dea7835ce885023ef7fe2e 100644
--- a/Core/Instrument/Beam.h
+++ b/Core/Instrument/Beam.h
@@ -47,6 +47,8 @@ public:
     //! Sets the polarization density matrix according to the given Bloch vector
     void setPolarization(const kvector_t bloch_vector);
 
+    kvector_t getBlochVector() const;
+
 #ifndef SWIG
     //! Returns the polarization density matrix (in spin basis along z-axis)
     Eigen::Matrix2cd getPolarization() const  { return m_polarization; }
diff --git a/GUI/coregui/Models/GUIObjectBuilder.cpp b/GUI/coregui/Models/GUIObjectBuilder.cpp
index 3a6d226a564d6d9de25e1364793960519dcd3c75..d3b92c390036b5469ab61e82614e8425505daaf1 100644
--- a/GUI/coregui/Models/GUIObjectBuilder.cpp
+++ b/GUI/coregui/Models/GUIObjectBuilder.cpp
@@ -130,7 +130,8 @@ SessionItem* GUIObjectBuilder::populateDocumentModel(DocumentModel* p_document_m
     Q_ASSERT(p_options_item);
     if (simulation.getOptions().isIntegrate()) {
         p_options_item->setComputationMethod(Constants::SIMULATION_MONTECARLO);
-        p_options_item->setNumberOfMonteCarloPoints(static_cast<int>(simulation.getOptions().getMcPoints()));
+        p_options_item->setNumberOfMonteCarloPoints(
+                    static_cast<int>(simulation.getOptions().getMcPoints()));
     }
     if (simulation.getOptions().useAvgMaterials()) {
         p_options_item->setFresnelMaterialMethod(Constants::AVERAGE_LAYER_MATERIAL);
@@ -196,7 +197,10 @@ void GUIObjectBuilder::visit(const MultiLayer* p_sample)
     SessionItem* p_multilayer_item =
             m_sampleModel->insertNewItem(Constants::MultiLayerType);
     p_multilayer_item->setItemName(p_sample->getName().c_str());
-    p_multilayer_item->setItemValue(MultiLayerItem::P_CROSS_CORR_LENGTH, p_sample->crossCorrLength());
+    p_multilayer_item->setItemValue(MultiLayerItem::P_CROSS_CORR_LENGTH,
+                                    p_sample->crossCorrLength());
+    p_multilayer_item->setVectorItem(MultiLayerItem::P_EXTERNAL_FIELD,
+                                     p_sample->externalField());
     m_levelToParentItem[depth()] = p_multilayer_item;
     m_itemToSample[p_multilayer_item] = p_sample;
 }
@@ -571,7 +575,7 @@ void GUIObjectBuilder::visit(const RotationEuler* p_sample)
     Q_ASSERT(transformation_item);
     SessionItem* p_rotationItem = transformation_item->setGroupProperty(
                 TransformationItem::P_ROT, Constants::EulerRotationType);
-    p_rotationItem->setItemValue(EulerRotationItem::P_ALPHA,  Units::rad2deg(p_sample->getAlpha()) );
+    p_rotationItem->setItemValue(EulerRotationItem::P_ALPHA, Units::rad2deg(p_sample->getAlpha()) );
     p_rotationItem->setItemValue(EulerRotationItem::P_BETA, Units::rad2deg(p_sample->getBeta()) );
     p_rotationItem->setItemValue(EulerRotationItem::P_GAMMA, Units::rad2deg(p_sample->getGamma()) );
     m_levelToParentItem[depth()] = transformation_item;
diff --git a/GUI/coregui/Models/TransformFromDomain.cpp b/GUI/coregui/Models/TransformFromDomain.cpp
index d87453008a57327b0045b32b381b7e711a241697..98b89e2776768399e8630c9c6f48c452b0025019 100644
--- a/GUI/coregui/Models/TransformFromDomain.cpp
+++ b/GUI/coregui/Models/TransformFromDomain.cpp
@@ -235,6 +235,9 @@ void TransformFromDomain::setItemFromSample(BeamItem* beamItem, const GISASSimul
             setItemFromSample(azimuthalAngle, distributions[i]);
         }
     }
+
+    // polarization parameters
+    beamItem->setVectorItem(BeamItem::P_POLARIZATION, beam.getBlochVector());
 }
 
 void TransformFromDomain::setInstrumentDetectorFromSample(InstrumentItem* instrumentItem,
diff --git a/auto/Wrap/libBornAgainCore.py b/auto/Wrap/libBornAgainCore.py
index 4b3294a5195857b9e1529d98f1b301887ab5a271..00cef83fa46acd5c0943d115abc7cfb4e96eed5f 100644
--- a/auto/Wrap/libBornAgainCore.py
+++ b/auto/Wrap/libBornAgainCore.py
@@ -3774,6 +3774,11 @@ class Beam(INode):
         return _libBornAgainCore.Beam_setPolarization(self, bloch_vector)
 
 
+    def getBlochVector(self):
+        """getBlochVector(Beam self) -> kvector_t"""
+        return _libBornAgainCore.Beam_getBlochVector(self)
+
+
     def getWavelength(self):
         """
         getWavelength(Beam self) -> double
diff --git a/auto/Wrap/libBornAgainCore_wrap.cpp b/auto/Wrap/libBornAgainCore_wrap.cpp
index c3d8bf40da988fede9c741a7812f2ac9753be5c9..c1c9db018e811aad5d0912015aeafa07d7587522 100644
--- a/auto/Wrap/libBornAgainCore_wrap.cpp
+++ b/auto/Wrap/libBornAgainCore_wrap.cpp
@@ -35214,6 +35214,28 @@ fail:
 }
 
 
+SWIGINTERN PyObject *_wrap_Beam_getBlochVector(PyObject *SWIGUNUSEDPARM(self), PyObject *args) {
+  PyObject *resultobj = 0;
+  Beam *arg1 = (Beam *) 0 ;
+  void *argp1 = 0 ;
+  int res1 = 0 ;
+  PyObject * obj0 = 0 ;
+  kvector_t result;
+  
+  if (!PyArg_ParseTuple(args,(char *)"O:Beam_getBlochVector",&obj0)) SWIG_fail;
+  res1 = SWIG_ConvertPtr(obj0, &argp1,SWIGTYPE_p_Beam, 0 |  0 );
+  if (!SWIG_IsOK(res1)) {
+    SWIG_exception_fail(SWIG_ArgError(res1), "in method '" "Beam_getBlochVector" "', argument " "1"" of type '" "Beam const *""'"); 
+  }
+  arg1 = reinterpret_cast< Beam * >(argp1);
+  result = ((Beam const *)arg1)->getBlochVector();
+  resultobj = SWIG_NewPointerObj((new kvector_t(static_cast< const kvector_t& >(result))), SWIGTYPE_p_BasicVector3DT_double_t, SWIG_POINTER_OWN |  0 );
+  return resultobj;
+fail:
+  return NULL;
+}
+
+
 SWIGINTERN PyObject *_wrap_Beam_getWavelength(PyObject *SWIGUNUSEDPARM(self), PyObject *args) {
   PyObject *resultobj = 0;
   Beam *arg1 = (Beam *) 0 ;
@@ -111444,6 +111466,7 @@ static PyMethodDef SwigMethods[] = {
 		"Sets the polarization density matrix according to the given Bloch vector. \n"
 		"\n"
 		""},
+	 { (char *)"Beam_getBlochVector", _wrap_Beam_getBlochVector, METH_VARARGS, (char *)"Beam_getBlochVector(Beam self) -> kvector_t"},
 	 { (char *)"Beam_getWavelength", _wrap_Beam_getWavelength, METH_VARARGS, (char *)"\n"
 		"Beam_getWavelength(Beam self) -> double\n"
 		"\n"