diff --git a/Device/Histo/SimulationResult.cpp b/Device/Histo/SimulationResult.cpp
index 6fab60ed2a726085f2132d3243bea1780ddf5a38..5108a31b2ef35f8d2d706b778b7d1d844e49e23c 100644
--- a/Device/Histo/SimulationResult.cpp
+++ b/Device/Histo/SimulationResult.cpp
@@ -89,26 +89,20 @@ std::vector<AxisInfo> SimulationResult::axisInfo(Axes::Units units) const
 
 const IUnitConverter& SimulationResult::converter() const
 {
-    if (!m_unit_converter)
-        throw std::runtime_error(
-            "Error in SimulationResult::converter: unit converter was not initialized");
+    ASSERT(m_unit_converter);
     return *m_unit_converter;
 }
 
 double& SimulationResult::operator[](size_t i)
 {
-    if (m_data)
-        return (*m_data)[i];
-    throw std::runtime_error("Error in SimulationResult::operator[]: "
-                             "no data initialized");
+    ASSERT(m_data);
+    return (*m_data)[i];
 }
 
 const double& SimulationResult::operator[](size_t i) const
 {
-    if (m_data)
-        return (*m_data)[i];
-    throw std::runtime_error("Error in SimulationResult::operator[]: "
-                             "no data initialized");
+    ASSERT(m_data);
+    return (*m_data)[i];
 }
 
 size_t SimulationResult::size() const
@@ -116,6 +110,16 @@ size_t SimulationResult::size() const
     return m_data ? m_data->getAllocatedSize() : 0;
 }
 
+double SimulationResult::max() const
+{
+    ASSERT(m_data);
+    double result = 0;
+    for (size_t i=0; i<size(); ++i)
+        if ((*m_data)[i]>result)
+            result = (*m_data)[i];
+    return result;
+}
+
 #ifdef BORNAGAIN_PYTHON
 PyObject* SimulationResult::array(Axes::Units units) const
 {
diff --git a/Device/Histo/SimulationResult.h b/Device/Histo/SimulationResult.h
index ba5c6126b7380926ef41efd281eb317fbae70f37..d5011cb88841acf02016e524708f65fecf619bfe 100644
--- a/Device/Histo/SimulationResult.h
+++ b/Device/Histo/SimulationResult.h
@@ -65,6 +65,7 @@ public:
     double& operator[](size_t i);
     const double& operator[](size_t i) const;
     size_t size() const;
+    double max() const;
     bool empty() const { return size() == 0; }
 
     //! returns intensity data as Python numpy array
diff --git a/auto/Wrap/doxygenDevice.i b/auto/Wrap/doxygenDevice.i
index a5bdaf8d91d8103424946a7530a74b51cb9d97c8..0896fd73285617fc92988420c74bf8191db50fbc 100644
--- a/auto/Wrap/doxygenDevice.i
+++ b/auto/Wrap/doxygenDevice.i
@@ -2870,6 +2870,9 @@ Returns underlying unit converter.
 %feature("docstring")  SimulationResult::size "size_t SimulationResult::size() const
 ";
 
+%feature("docstring")  SimulationResult::max "double & SimulationResult::max() const
+";
+
 %feature("docstring")  SimulationResult::empty "bool SimulationResult::empty() const
 ";
 
diff --git a/auto/Wrap/libBornAgainDevice.py b/auto/Wrap/libBornAgainDevice.py
index a7e8a2f4c37c61a3486483fb93c59d34e04520a7..6a4529adb44be29552a1703d0a1bb676eaaf2a74 100644
--- a/auto/Wrap/libBornAgainDevice.py
+++ b/auto/Wrap/libBornAgainDevice.py
@@ -6124,6 +6124,14 @@ class SimulationResult(object):
         """
         return _libBornAgainDevice.SimulationResult_size(self)
 
+    def max(self):
+        r"""
+        max(SimulationResult self) -> double
+        double & SimulationResult::max() const
+
+        """
+        return _libBornAgainDevice.SimulationResult_max(self)
+
     def empty(self):
         r"""
         empty(SimulationResult self) -> bool
diff --git a/auto/Wrap/libBornAgainDevice_wrap.cpp b/auto/Wrap/libBornAgainDevice_wrap.cpp
index 6aaa54411afb1de906dc21f0eccd487fb5f9e684..cf6146414700ecf6c3a811ff41b7129b4732df12 100644
--- a/auto/Wrap/libBornAgainDevice_wrap.cpp
+++ b/auto/Wrap/libBornAgainDevice_wrap.cpp
@@ -45213,6 +45213,29 @@ fail:
 }
 
 
+SWIGINTERN PyObject *_wrap_SimulationResult_max(PyObject *SWIGUNUSEDPARM(self), PyObject *args) {
+  PyObject *resultobj = 0;
+  SimulationResult *arg1 = (SimulationResult *) 0 ;
+  void *argp1 = 0 ;
+  int res1 = 0 ;
+  PyObject *swig_obj[1] ;
+  double result;
+  
+  if (!args) SWIG_fail;
+  swig_obj[0] = args;
+  res1 = SWIG_ConvertPtr(swig_obj[0], &argp1,SWIGTYPE_p_SimulationResult, 0 |  0 );
+  if (!SWIG_IsOK(res1)) {
+    SWIG_exception_fail(SWIG_ArgError(res1), "in method '" "SimulationResult_max" "', argument " "1"" of type '" "SimulationResult const *""'"); 
+  }
+  arg1 = reinterpret_cast< SimulationResult * >(argp1);
+  result = (double)((SimulationResult const *)arg1)->max();
+  resultobj = SWIG_From_double(static_cast< double >(result));
+  return resultobj;
+fail:
+  return NULL;
+}
+
+
 SWIGINTERN PyObject *_wrap_SimulationResult_empty(PyObject *SWIGUNUSEDPARM(self), PyObject *args) {
   PyObject *resultobj = 0;
   SimulationResult *arg1 = (SimulationResult *) 0 ;
@@ -48688,6 +48711,11 @@ static PyMethodDef SwigMethods[] = {
 		"size_t SimulationResult::size() const\n"
 		"\n"
 		""},
+	 { "SimulationResult_max", _wrap_SimulationResult_max, METH_O, "\n"
+		"SimulationResult_max(SimulationResult self) -> double\n"
+		"double & SimulationResult::max() const\n"
+		"\n"
+		""},
 	 { "SimulationResult_empty", _wrap_SimulationResult_empty, METH_O, "\n"
 		"SimulationResult_empty(SimulationResult self) -> bool\n"
 		"bool SimulationResult::empty() const\n"