Merge pull request #106 from boschresearch/105-meanstress-transform-r…

…-goal

Preserve the mean stress intervals index transforming a histogram

src/pylife/strength/meanstress.py

@@ -422,15 +422,19 @@ def _validate(self):
 
     def fkm_goodman(self, haigh, R_goal):
         ranges = HaighDiagram.fkm_goodman(haigh).transform(self._obj, R_goal)['range']
-        return self._rebin_results(ranges).load_collective
+        return self._rebin_results(ranges, R_goal).load_collective
 
-    def _rebin_results(self, ranges):
+    def _rebin_results(self, ranges, R_goal):
 
         def resulting_intervals():
             ranges_max = ranges.max()
             binsize = np.hypot(self._binsize_x, self._binsize_y) / np.sqrt(2.)
             bincount = int(np.ceil(ranges_max / binsize))
-            return pd.IntervalIndex.from_breaks(np.linspace(0, ranges_max, bincount+1), name="range")
+            range_bins = np.linspace(0, ranges_max, bincount+1)
+            means_bins = range_bins * (1. + R_goal) / (2.0 * (1. - R_goal))
+            range_itv = pd.IntervalIndex.from_breaks(range_bins, name="range")
+            means_itv = pd.IntervalIndex.from_breaks(means_bins, name="mean")
+            return range_itv, means_itv
 
         def aggregate_on_projection(projection, itvs, ranges, obj):
             level_values = tuple(projection)
@@ -452,18 +456,20 @@ def sum_intervals(iv, ranges, obj):
             new_idx = pd.IntervalIndex(pd.interval_range(0.,  0., 0), name='range')
             return pd.Series([], index=new_idx, name='cycles', dtype=np.float64)
 
-        intv_idx = resulting_intervals()
+        range_itv_idx, means_itg_idx = resulting_intervals()
 
         if len(self._remaining_names) > 0:
             remaining_idx = self._obj.index.to_frame(index=False).groupby(self._remaining_names).first().index
             result = (
                 remaining_idx.to_frame(index=False)
-                .apply(lambda p: aggregate_on_projection(p, intv_idx, ranges, self._obj), axis=1)
+                .apply(lambda p: aggregate_on_projection(p, range_itv_idx, ranges, self._obj), axis=1)
             )
             result.index = remaining_idx
-            result = result.stack().reorder_levels(['range'] + self._remaining_names)
+            result.columns = pd.MultiIndex.from_arrays([result.columns, means_itg_idx])
+            result = result.stack(['range', 'mean']).reorder_levels(['range', 'mean'] + self._remaining_names)
         else:
-            result = intv_idx.to_series().apply(lambda iv: sum_intervals(iv, ranges, self._obj))
+            result = range_itv_idx.to_series().apply(lambda iv: sum_intervals(iv, ranges, self._obj))
+            result.index = pd.MultiIndex.from_arrays([result.index, means_itg_idx])
 
         return result
 

tests/strength/test_meanstress.py

@@ -197,6 +197,26 @@ def test_five_segment_multiple_M_sm():
     np.testing.assert_array_almost_equal(res.amplitude, np.ones_like(res))
 
 
+#GH-105
+@pytest.mark.parametrize("R_goal", [-1., 0., -1./3., 1./3.])
+def test_fkm_goodman_hist_R_goal(R_goal):
+    rg = pd.IntervalIndex.from_breaks(np.linspace(0., 2., 25), closed='left')
+    mn = pd.IntervalIndex.from_breaks(np.linspace(-1./12., 23./12., 25), closed='left')
+
+    mat = pd.Series(np.zeros(24*24), name='cycles',
+                    index=pd.MultiIndex.from_product([rg, mn], names=['range', 'mean']))
+    mat.loc[(7./6., 7./6.)] = 1.
+    mat.loc[(4./3., 2./3.)] = 3.
+    mat.loc[(2. - 1e-9, 0.)] = 5.
+
+    haigh = pd.Series({'M': 0.5, 'M2': 0.5/3.})
+    result = mat.meanstress_transform.fkm_goodman(haigh, R_goal).R
+
+    expected = pd.Series(R_goal, index=result.index, name=result.name)
+
+    pd.testing.assert_series_equal(result, expected)
+
+
 @pytest.mark.parametrize("R_goal, expected", [  # all calculated by pencil on paper
     (-1., 2.0),
     (0., 4./3.),
@@ -216,8 +236,10 @@ def test_fkm_goodman_hist_range_mean(R_goal, expected):
     haigh = pd.Series({'M': 0.5, 'M2': 0.5/3.})
     res = mat.meanstress_transform.fkm_goodman(haigh, R_goal).to_pandas()
     test_interval = pd.Interval(expected-1./96., expected+1./96.)
-    assert res.loc[res.index.overlaps(test_interval)].sum() == 9
-    assert res.loc[np.logical_not(res.index.overlaps(test_interval))].sum() == 0
+
+    mask = res.index.get_level_values('range').overlaps(test_interval)
+    assert res.loc[mask].sum() == 9
+    assert res.loc[~mask].sum() == 0
 
 
 @pytest.mark.parametrize("R_goal, expected", [  # all calculated by pencil on paper
@@ -240,8 +262,10 @@ def test_fkm_goodman_hist_from_to(R_goal, expected):
     res = mat.meanstress_transform.fkm_goodman(haigh, R_goal).to_pandas()
 
     test_interval = pd.Interval(expected-1./96., expected+1./96.)
-    assert res.loc[res.index.overlaps(test_interval)].sum() == 9
-    assert res.loc[np.logical_not(res.index.overlaps(test_interval))].sum() == 0
+
+    mask = res.index.get_level_values('range').overlaps(test_interval)
+    assert res.loc[mask].sum() == 9
+    assert res.loc[~mask].sum() == 0
 
 
 def test_meanstress_transform_additional_index():
@@ -270,15 +294,18 @@ def test_meanstress_transform_additional_index():
     haigh = pd.Series({'M': 0.5, 'M2': 0.5/3.})
     res = mat.meanstress_transform.fkm_goodman(haigh, -1.0).to_pandas()
 
-    assert set(res.index.names) == {'range', 'node_id'}
+    assert set(res.index.names) == {'range', 'mean', 'node_id'}
+
+    null_itv = pd.Interval(0.0, 0.0)
 
-    assert res.loc[(2.0, 1)] == 9
-    assert res.loc[(2.0, 2)] == 18
-    assert res.loc[(2.0, 3)] == 36
+    assert res.loc[(2.0, null_itv, 1)] == 9
+    assert res.loc[(2.0, null_itv, 2)] == 18
+    assert res.loc[(2.0, null_itv, 3)] == 36
 
     assert res.sum() == 9 + 18 + 36
     assert res.min() == 0
 
+
 def test_meanstress_transform_two_additional_indices():
     fr = pd.IntervalIndex.from_breaks(np.linspace(-25./24., 1., 49), closed='left', name='from')
     to = pd.IntervalIndex.from_breaks(np.linspace(-1./24., 2., 49), closed='left', name='to')
@@ -310,13 +337,15 @@ def test_meanstress_transform_two_additional_indices():
     haigh = pd.Series({'M': 0.5, 'M2': 0.5/3.})
     res = mat.meanstress_transform.fkm_goodman(haigh, -1.0).to_pandas()
 
-    assert set(res.index.names) == {'range', 'node_id', 'element_id'}
+    assert set(res.index.names) == {'range', 'mean', 'node_id', 'element_id'}
+
+    null_itv = pd.Interval(0.0, 0.0)
 
-    assert res.loc[(2.0, 1, 1)] == 9
-    assert res.loc[(2.0, 2, 1)] == 18
+    assert res.loc[(2.0, null_itv, 1, 1)] == 9
+    assert res.loc[(2.0, null_itv, 2, 1)] == 18
 
-    assert res.loc[(2.0, 1, 2)] == 36
-    assert res.loc[(2.0, 2, 2)] == 72
+    assert res.loc[(2.0, null_itv, 1, 2)] == 36
+    assert res.loc[(2.0, null_itv, 2, 2)] == 72
 
     assert res.sum() == 9 + 18 + 36 + 72
     assert res.min() == 0
@@ -343,8 +372,9 @@ def test_fkm_goodman_hist_range_mean_nonzero(R_goal, expected):
 
     test_interval = pd.Interval(expected-1./96., expected+1./96.)
 
-    assert res.loc[res.index.overlaps(test_interval)].sum() == 9
-    assert res.loc[np.logical_not(res.index.overlaps(test_interval))].sum() == 0
+    mask = res.index.get_level_values('range').overlaps(test_interval)
+    assert res.loc[mask].sum() == 9
+    assert res.loc[~mask].sum() == 0
 
     binsize = res.index.get_level_values('range').length.min()
     np.testing.assert_approx_equal(binsize, 2./24., significant=1)

tests/strength/test_meanstress_matrix.py

@@ -69,8 +69,9 @@ def test_meanstress_collective_fkm_goodman_single_ms_sens(ms_sens):
     expected_interval = pd.Interval(expected - 1./96., expected + 1./96.)
     res = rf.meanstress_transform.fkm_goodman(haigh, R_goal).to_pandas()
 
-    assert res.loc[res.index.overlaps(expected_interval)].sum() == 9
-    assert res.loc[~res.index.overlaps(expected_interval)].sum() == 0
+    mask = res.index.get_level_values('range').overlaps(expected_interval)
+    assert res.loc[mask].sum() == 9
+    assert res.loc[~mask].sum() == 0
 
 
 def test_meanstress_collective_fkm_goodman_multiple_ms_sens():