Repository: StingraySoftware/stingray
Pull request #798 (2024-02-08 11:38:42)

Fix count rate estimation error for particularly bad GTIs
In `StingrayTimeseries.fill_bad_time_intervals`, when the buffer size is small and count rates are low, there was the possibility that the count rate estimation failed with a numpy error. E.g., the test `test_no_counts_in_buffer` I added failed. Here, I catch the condition and warn that this is happening. Also, I made the estimation more robust when only one of the GTIs on the side of the bad time intervals has good data.

pep8speaks: Hello @matteobachetti! Thanks for opening this PR. We checked the lines you've touched for [PEP 8](https://www.python.org/dev/peps/pep-0008) issues, and found:

* In the file [`stingray/base.py`](https://github.com/StingraySoftware/stingray/blob/db5674d1456055bfc6533335d0765601dc01d1c8/stingray/base.py):

> [Line 2294:101](https://github.com/StingraySoftware/stingray/blob/db5674d1456055bfc6533335d0765601dc01d1c8/stingray/base.py#L2294): [E501](https://duckduckgo.com/?q=pep8%20E501) line too long (116 > 100 characters)


codecov[bot]: ## [Codecov](https://app.codecov.io/gh/StingraySoftware/stingray/pull/798?src=pr&el=h1&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=StingraySoftware) Report
Attention: `2 lines` in your changes are missing coverage. Please review.
> Comparison is base [(`5285c23`)](https://app.codecov.io/gh/StingraySoftware/stingray/commit/5285c23221ed81c946007e8ae116c84fdc9947d4?el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=StingraySoftware) 96.35% compared to head [(`db5674d`)](https://app.codecov.io/gh/StingraySoftware/stingray/pull/798?src=pr&el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=StingraySoftware) 81.90%.

| [Files](https://app.codecov.io/gh/StingraySoftware/stingray/pull/798?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=StingraySoftware) | Patch % | Lines |
|---|---|---|
| [stingray/base.py](https://app.codecov.io/gh/StingraySoftware/stingray/pull/798?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=StingraySoftware#diff-c3RpbmdyYXkvYmFzZS5weQ==) | 77.77% | [2 Missing :warning: ](https://app.codecov.io/gh/StingraySoftware/stingray/pull/798?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=StingraySoftware) |

<details><summary>Additional details and impacted files</summary>


```diff
@@             Coverage Diff             @@
##             main     #798       +/-   ##
===========================================
- Coverage   96.35%   81.90%   -14.46%     
===========================================
  Files          43       43               
  Lines        8784     8792        +8     
===========================================
- Hits         8464     7201     -1263     
- Misses        320     1591     +1271     
```



</details>

[:umbrella: View full report in Codecov by Sentry](https://app.codecov.io/gh/StingraySoftware/stingray/pull/798?src=pr&el=continue&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=StingraySoftware).   
:loudspeaker: Have feedback on the report? [Share it here](https://about.codecov.io/codecov-pr-comment-feedback/?utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=StingraySoftware).


Patch:
diff --git a/docs/notebooks b/docs/notebooks
index ccc49119..d3443305 160000
--- a/docs/notebooks
+++ b/docs/notebooks
@@ -1,1 +1,1 @@
-Subproject commit ccc491197b11b018304ca4397af7968a6111a2bb
+Subproject commit d344330557a96105e83304fe4b931287395201cd
diff --git a/stingray/base.py b/stingray/base.py
index f9734725..a3d4a4f7 100644
--- a/stingray/base.py
+++ b/stingray/base.py
@@ -1,4 +1,5 @@
 """Base classes"""
+
 from __future__ import annotations
 
 from collections.abc import Iterable
@@ -2285,12 +2286,27 @@ class StingrayTimeseries(StingrayObject):
                 nevents = local_new_times.size
             else:
                 low_time_arr = filtered_times[max(filt_low_idx - buffer_size, 0) : filt_low_idx]
+                low_time_arr = low_time_arr[low_time_arr > bti[0] - buffer_size]
                 high_time_arr = filtered_times[filt_hig_idx : buffer_size + filt_hig_idx]
+                high_time_arr = high_time_arr[high_time_arr < bti[1] + buffer_size]
+
+                if len(low_time_arr) > 0 and (filt_low_t - low_time_arr[0]) > 0:
+                    ctrate_low = np.count_nonzero(low_time_arr) / (filt_low_t - low_time_arr[0])
+                else:
+                    ctrate_low = np.nan
+                if len(high_time_arr) > 0 and (high_time_arr[-1] - filt_hig_t) > 0:
+                    ctrate_high = np.count_nonzero(high_time_arr) / (high_time_arr[-1] - filt_hig_t)
+                else:
+                    ctrate_high = np.nan
 
-                ctrate = (
-                    np.count_nonzero(low_time_arr) / (filt_low_t - low_time_arr[0])
-                    + np.count_nonzero(high_time_arr) / (high_time_arr[-1] - filt_hig_t)
-                ) / 2
+                if not np.isfinite(ctrate_low) and not np.isfinite(ctrate_high):
+                    warnings.warn(
+                        f"No valid data around to simulate the time series in interval "
+                        f"{bti[0]:g}-{bti[1]:g}. Skipping. Please check that the buffer size is "
+                        f"adequate."
+                    )
+                    continue
+                ctrate = np.nanmean([ctrate_low, ctrate_high])
                 nevents = rs.poisson(ctrate * (bti[1] - bti[0]))
                 local_new_times = rs.uniform(bti[0], bti[1], nevents)
             new_times.append(local_new_times)


Test patch:
diff --git a/stingray/tests/test_base.py b/stingray/tests/test_base.py
index ecdc8a18..1766bec7 100644
--- a/stingray/tests/test_base.py
+++ b/stingray/tests/test_base.py
@@ -1327,6 +1327,18 @@ class TestFillBTI(object):
         for attr in ["time", "energy", "blablas"]:
             assert np.allclose(getattr(new_masked, attr), getattr(filt_masked, attr))
 
+    def test_no_counts_in_buffer(self):
+        ev_like_filt = copy.deepcopy(self.ev_like)
+        # I introduce a small gap in the GTIs
+        ev_like_filt.gti = np.asarray([[0, 490], [491, 498], [500, 505], [510, 520], [522, 700]])
+
+        # I empty out two GTIs
+        bad = (ev_like_filt.time > 490) & (ev_like_filt.time < 510)
+        ev_like_filt = ev_like_filt.apply_mask(~bad)
+
+        with pytest.warns(UserWarning, match="simulate the time series in interval 498-500"):
+            ev_like_filt.fill_bad_time_intervals(max_length=3, buffer_size=2)
+
     def test_lc_like(self):
         lc_like_filt = copy.deepcopy(self.lc_like)
         # I introduce a small gap in the GTIs