diff --git a/Device/Beam/Beam.cpp b/Device/Beam/Beam.cpp
index 6ef5c9dc0b0f2f8c52440000ded4c7d3250cdd49..491a446f0b193b33f31d742d13d1a0e12a90b4f9 100644
--- a/Device/Beam/Beam.cpp
+++ b/Device/Beam/Beam.cpp
@@ -147,12 +147,12 @@ void Beam::setPolarization(const R3 bloch_vector)
     m_beamPolarization = bloch_vector;
 }
 
-R3 Beam::getBlochVector() const
+R3 Beam::polVector() const
 {
     return m_beamPolarization;
 }
 
-SpinMatrix Beam::getPolarization() const
+SpinMatrix Beam::polMatrix() const
 {
     return SpinMatrix::FromBlochVector(m_beamPolarization);
 }
diff --git a/Device/Beam/Beam.h b/Device/Beam/Beam.h
index 9b4e72d197df2f5bd44a4cd75d20981a4d86fa03..99c8b747e74c3f116e9ef424b49bb394804da720 100644
--- a/Device/Beam/Beam.h
+++ b/Device/Beam/Beam.h
@@ -45,7 +45,11 @@ public:
     // Direction& direction() { return m_direction; }
     Direction direction() const { return {m_alpha, m_phi}; } // TODO -> const .. &
 
-    R3 getBlochVector() const;
+    //! Returns polarization density as Bloch vector
+    R3 polVector() const;
+    //! Returns the polarization density matrix (in spin basis along z-axis)
+    SpinMatrix polMatrix() const;
+
     //! Returns footprint factor.
     const IFootprintFactor* footprintFactor() const;
 
@@ -61,8 +65,6 @@ public:
     //! Throws if limits are violated.
     void setWavelengthGuarded(double value);
 
-    //! Returns the polarization density matrix (in spin basis along z-axis)
-    SpinMatrix getPolarization() const;
 
     void setWavelength(double wavelength);
     void setDirection(const Direction& direction);
diff --git a/Device/Pol/PolFilter.cpp b/Device/Pol/PolFilter.cpp
index 04e7ae4216e24fe8d94b0ff60e785da2144ae214..fe9b65b4a6ec8748231ab4b212c5907da5c529bb 100644
--- a/Device/Pol/PolFilter.cpp
+++ b/Device/Pol/PolFilter.cpp
@@ -49,14 +49,8 @@ SpinMatrix PolFilter::matrix() const
 {
     if (m_direction.mag() == 0.0 || m_efficiency == 0.0)
         return m_total_transmission * SpinMatrix::One();
-    double x = m_direction.x() / m_direction.mag();
-    double y = m_direction.y() / m_direction.mag();
-    double z = m_direction.z() / m_direction.mag();
-    double sum = m_total_transmission * 2.0;
-    double diff = m_total_transmission * m_efficiency * 2.0;
-    complex_t im(0.0, 1.0);
-    return {(sum + diff * z) / 2.0, diff * (x - im * y) / 2.0, diff * (x + im * y) / 2.0,
-            (sum - diff * z) / 2.0};
+    R3 v = m_direction.unit() * m_efficiency;
+    return 2* m_total_transmission * SpinMatrix::FromBlochVector(v);
 }
 
 R3 PolFilter::polDirection() const
diff --git a/GUI/View/FromDomain/FromDomain.cpp b/GUI/View/FromDomain/FromDomain.cpp
index 2c8c90b8abb6c287067cd5fd1b2136b7f9d9442a..31ecaa73ddd66adb4b6728e38332281dc15f311a 100644
--- a/GUI/View/FromDomain/FromDomain.cpp
+++ b/GUI/View/FromDomain/FromDomain.cpp
@@ -547,7 +547,7 @@ void GUI::Transform::FromDomain::setDetectorResolution(DetectorItem* detector_it
 void GUI::Transform::FromDomain::setPolarizerAnalyzer(Instrument2DItem* instrument_item,
                                                       const ISimulation2D& simulation)
 {
-    instrument_item->setPolarization(simulation.beam().getBlochVector());
+    instrument_item->setPolarization(simulation.beam().polVector());
 
     const IDetector& detector = simulation.detector();
     double total_transmission = detector.analyzer().totalTransmission();
diff --git a/Sim/Export/SimulationToPython.cpp b/Sim/Export/SimulationToPython.cpp
index d47a577da3b8c3b5cc2e64b0cee301a0545e65a3..e7e6c4172bc5d626de47a4bf8b8be9ea785e7202 100644
--- a/Sim/Export/SimulationToPython.cpp
+++ b/Sim/Export/SimulationToPython.cpp
@@ -260,7 +260,7 @@ std::string defineDetectorPolarizationAnalysis(const ISimulation* simulation)
 std::string defineBeamPolarization(const Beam& beam)
 {
     std::ostringstream result;
-    auto bloch_vector = beam.getBlochVector();
+    auto bloch_vector = beam.polVector();
     if (bloch_vector.mag() > 0.0) {
         std::string beam_polpair = "beam_polpair";
         result << indent() << beam_polpair << " = R3(" << Py::Fmt::printDouble(bloch_vector.x())
diff --git a/Sim/Simulation/ISimulation2D.cpp b/Sim/Simulation/ISimulation2D.cpp
index 08129ad3bc9d0950b1562493b11e6f97743758a3..ee3e65be207e597bf15395dbb49f61868cde7cbc 100644
--- a/Sim/Simulation/ISimulation2D.cpp
+++ b/Sim/Simulation/ISimulation2D.cpp
@@ -92,7 +92,7 @@ std::vector<std::unique_ptr<DiffuseElement>> ISimulation2D::generateElements(con
     const double wavelength = beam.wavelength();
     const double alpha_i = -beam.direction().alpha(); // Defined to be always positive in Beam
     const double phi_i = beam.direction().phi();
-    const SpinMatrix beam_polpair = beam.getPolarization();
+    const SpinMatrix beam_polpair = beam.polMatrix();
 
     const IDetector2D& detector = detector2D();
     const SpinMatrix analyzer_operator = detector.analyzer().matrix();
diff --git a/Sim/Simulation/SpecularSimulation.cpp b/Sim/Simulation/SpecularSimulation.cpp
index b3fb4c14c664903314a0f6313edd25efcfc44b71..73dd20361bf3ecf18fbdae2921abe91691757f02 100644
--- a/Sim/Simulation/SpecularSimulation.cpp
+++ b/Sim/Simulation/SpecularSimulation.cpp
@@ -213,5 +213,5 @@ void SpecularSimulation::moveDataFromCache()
 
 PolMatrices SpecularSimulation::polarizerPair() const
 {
-    return {beam().getPolarization(), detector().analyzer().matrix()};
+    return {beam().polMatrix(), detector().analyzer().matrix()};
 }
diff --git a/Tests/Unit/Sim/SpecularSimulationTest.cpp b/Tests/Unit/Sim/SpecularSimulationTest.cpp
index 370f643de33cac451f427669e47d6270a64ac34a..6693a27aee11048ff608d04a56a2bf8fc5294029 100644
--- a/Tests/Unit/Sim/SpecularSimulationTest.cpp
+++ b/Tests/Unit/Sim/SpecularSimulationTest.cpp
@@ -86,7 +86,7 @@ TEST_F(SpecularSimulationTest, SetAngularScan)
     EXPECT_EQ(0.0, sim4.beam().direction().alpha());
     EXPECT_EQ(0.0, sim4.beam().direction().phi());
 
-    EXPECT_EQ(sim4.beam().getBlochVector(), polarizer_dir);
+    EXPECT_EQ(sim4.beam().polVector(), polarizer_dir);
     EXPECT_EQ(sim4.detector().analyzer().polDirection(), analyzer_dir);
     EXPECT_EQ(sim4.detector().analyzer().polEfficiency(), 0.33);
     EXPECT_EQ(sim4.detector().analyzer().totalTransmission(), 0.22);
@@ -134,7 +134,7 @@ TEST_F(SpecularSimulationTest, SetQScan)
     EXPECT_EQ(0.0, sim3.beam().direction().alpha());
     EXPECT_EQ(0.0, sim3.beam().direction().phi());
 
-    EXPECT_EQ(sim3.beam().getBlochVector(), polarizer_dir);
+    EXPECT_EQ(sim3.beam().polVector(), polarizer_dir);
     EXPECT_EQ(sim3.detector().analyzer().polDirection(), analyzer_dir);
     EXPECT_EQ(sim3.detector().analyzer().polEfficiency(), 0.33);
     EXPECT_EQ(sim3.detector().analyzer().totalTransmission(), 0.22);
diff --git a/auto/Wrap/doxygenDevice.i b/auto/Wrap/doxygenDevice.i
index 2f4eb2126598287f51c277f56f54ca402119fe95..0a6f788bb683fc3af695248dd31962c253f33261 100644
--- a/auto/Wrap/doxygenDevice.i
+++ b/auto/Wrap/doxygenDevice.i
@@ -77,7 +77,14 @@ Returns the beam intensity in neutrons/sec.
 %feature("docstring")  Beam::direction "Direction Beam::direction() const
 ";
 
-%feature("docstring")  Beam::getBlochVector "R3 Beam::getBlochVector() const
+%feature("docstring")  Beam::polVector "R3 Beam::polVector() const
+
+Returns polarization density as Bloch vector. 
+";
+
+%feature("docstring")  Beam::polMatrix "SpinMatrix Beam::polMatrix() const
+
+Returns the polarization density matrix (in spin basis along z-axis) 
 ";
 
 %feature("docstring")  Beam::footprintFactor "const IFootprintFactor * Beam::footprintFactor() const
@@ -100,11 +107,6 @@ Check for limits, set the value if within limits. Throws if limits are violated.
 Check for limits, set the value if within limits. Throws if limits are violated. 
 ";
 
-%feature("docstring")  Beam::getPolarization "SpinMatrix Beam::getPolarization() const
-
-Returns the polarization density matrix (in spin basis along z-axis) 
-";
-
 %feature("docstring")  Beam::setWavelength "void Beam::setWavelength(double wavelength)
 ";
 
diff --git a/auto/Wrap/libBornAgainDevice.py b/auto/Wrap/libBornAgainDevice.py
index 28465d4ea6df848027250eb878e8b43dd558cd45..81e57e112ff4a31ade6b2d86f355ba98ee146007 100644
--- a/auto/Wrap/libBornAgainDevice.py
+++ b/auto/Wrap/libBornAgainDevice.py
@@ -2406,13 +2406,25 @@ class Beam(libBornAgainParam.INode):
         """
         return _libBornAgainDevice.Beam_direction(self)
 
-    def getBlochVector(self):
+    def polVector(self):
         r"""
-        getBlochVector(Beam self) -> R3
-        R3 Beam::getBlochVector() const
+        polVector(Beam self) -> R3
+        R3 Beam::polVector() const
+
+        Returns polarization density as Bloch vector. 
+
+        """
+        return _libBornAgainDevice.Beam_polVector(self)
+
+    def polMatrix(self):
+        r"""
+        polMatrix(Beam self) -> SpinMatrix
+        SpinMatrix Beam::polMatrix() const
+
+        Returns the polarization density matrix (in spin basis along z-axis) 
 
         """
-        return _libBornAgainDevice.Beam_getBlochVector(self)
+        return _libBornAgainDevice.Beam_polMatrix(self)
 
     def footprintFactor(self):
         r"""
@@ -2454,16 +2466,6 @@ 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 d8e94d7b0ffc1580eeb1b951a88aaba3e145b17b..ed7af056d4eb8e9ab6b35970c186acc3bbdea396 100644
--- a/auto/Wrap/libBornAgainDevice_wrap.cpp
+++ b/auto/Wrap/libBornAgainDevice_wrap.cpp
@@ -28501,7 +28501,7 @@ fail:
 }
 
 
-SWIGINTERN PyObject *_wrap_Beam_getBlochVector(PyObject *SWIGUNUSEDPARM(self), PyObject *args) {
+SWIGINTERN PyObject *_wrap_Beam_polVector(PyObject *SWIGUNUSEDPARM(self), PyObject *args) {
   PyObject *resultobj = 0;
   Beam *arg1 = (Beam *) 0 ;
   void *argp1 = 0 ;
@@ -28513,10 +28513,10 @@ SWIGINTERN PyObject *_wrap_Beam_getBlochVector(PyObject *SWIGUNUSEDPARM(self), P
   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_getBlochVector" "', argument " "1"" of type '" "Beam const *""'"); 
+    SWIG_exception_fail(SWIG_ArgError(res1), "in method '" "Beam_polVector" "', argument " "1"" of type '" "Beam const *""'"); 
   }
   arg1 = reinterpret_cast< Beam * >(argp1);
-  result = ((Beam const *)arg1)->getBlochVector();
+  result = ((Beam const *)arg1)->polVector();
   resultobj = SWIG_NewPointerObj((new R3(static_cast< const R3& >(result))), SWIGTYPE_p_Vec3T_double_t, SWIG_POINTER_OWN |  0 );
   return resultobj;
 fail:
@@ -28524,6 +28524,29 @@ fail:
 }
 
 
+SWIGINTERN PyObject *_wrap_Beam_polMatrix(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_polMatrix" "', argument " "1"" of type '" "Beam const *""'"); 
+  }
+  arg1 = reinterpret_cast< Beam * >(argp1);
+  result = ((Beam const *)arg1)->polMatrix();
+  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_footprintFactor(PyObject *SWIGUNUSEDPARM(self), PyObject *args) {
   PyObject *resultobj = 0;
   Beam *arg1 = (Beam *) 0 ;
@@ -28634,29 +28657,6 @@ 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 ;
@@ -41732,9 +41732,18 @@ static PyMethodDef SwigMethods[] = {
 		"Direction Beam::direction() const\n"
 		"\n"
 		""},
-	 { "Beam_getBlochVector", _wrap_Beam_getBlochVector, METH_O, "\n"
-		"Beam_getBlochVector(Beam self) -> R3\n"
-		"R3 Beam::getBlochVector() const\n"
+	 { "Beam_polVector", _wrap_Beam_polVector, METH_O, "\n"
+		"Beam_polVector(Beam self) -> R3\n"
+		"R3 Beam::polVector() const\n"
+		"\n"
+		"Returns polarization density as Bloch vector. \n"
+		"\n"
+		""},
+	 { "Beam_polMatrix", _wrap_Beam_polMatrix, METH_O, "\n"
+		"Beam_polMatrix(Beam self) -> SpinMatrix\n"
+		"SpinMatrix Beam::polMatrix() const\n"
+		"\n"
+		"Returns the polarization density matrix (in spin basis along z-axis) \n"
 		"\n"
 		""},
 	 { "Beam_footprintFactor", _wrap_Beam_footprintFactor, METH_O, "\n"
@@ -41765,13 +41774,6 @@ 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"