enh: adding a criteria to have only candidates with minimum two points

fink_filters/rubin/livestream/filter_extragalactic_new_candidate/filter.py

@@ -1,5 +1,5 @@
 # Copyright 2019-2026 AstroLab Software
-# Author: Julien Peloton
+# Author: Julien Peloton, Camille Douzet
 #
 # Licensed under the Apache License, Version 2.0 (the "License");
 # you may not use this file except in compliance with the License.
@@ -12,19 +12,72 @@
 # WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 # See the License for the specific language governing permissions and
 # limitations under the License.
-"""Select LSST alerts new (< 48h first apparition) and potentially extragalactic"""
+"""Select LSST alerts new (< 5days first apparition), bright (mag < 24), potentially extragalactic with a fading or rising rate passing the cuts"""
 
 import pandas as pd
+
 import fink_filters.rubin.blocks as fb
+import fink_filters.rubin.utils as fu
+
+DESCRIPTION = "Select LSST alerts new (< 5days first apparition), bright (mag < 24), potentially extragalactic with a fading or rising rate passing the cuts"
+
+
+def has_two_points_same_band(
+    diaSource: pd.DataFrame, diaObject: pd.DataFrame
+) -> pd.Series:
+    """Check if an object has at least 2 detections in the same band as the current alert.
+
+    Parameters
+    ----------
+    diaSource: pd.DataFrame
+        Full diaSource section of an alert (dictionary exploded)
+    diaObject: pd.DataFrame
+        Full diaObject section of an alert (dictionary exploded)
+
+    Returns
+    -------
+    out: pd.Series of booleans
+        True if at least 2 detections exist in the current alert's band, False otherwise.
+    """
+    result = pd.Series(False, index=diaSource.index)
 
+    for b in ["u", "g", "r", "i", "z", "y"]:
+        mask = diaSource.band == b
+        col = f"{b}_psfFluxNdata"
+        try:
+            result[mask] = diaObject.loc[mask, col] >= 2
+        except KeyError:
+            continue
+    return result
 
-DESCRIPTION = (
-    "Select LSST alerts new (< 48h first apparition) and potentially extragalactic"
-)
+
+def get_latest_source_same_band(row) -> dict | None:
+    """Return the most recent previous detection in the same band as the current alert.
+
+    Parameters
+    ----------
+    row: pd.Series
+        A single row from a DataFrame combining diaSource and prvDiaSources columns.
+
+    Returns
+    -------
+    out: dict or None
+        The most recent previous diaSource dict in the same band, or None if none exists.
+    """
+    sources = row["prvDiaSources"]
+    if sources is None or len(sources) == 0:
+        return None
+    current_band = row["band"]
+    same_band = [s for s in sources if s["band"] == current_band]
+    if not same_band:
+        return None
+    return max(same_band, key=lambda s: s["midpointMjdTai"])
 
 
 def extragalactic_new_candidate(
     diaSource: pd.DataFrame,
+    diaObject: pd.DataFrame,
+    prvDiaSources: pd.Series,
     simbad_otype: pd.Series,
     mangrove_lum_dist: pd.Series,
     is_sso: pd.Series,
@@ -35,20 +88,22 @@ def extragalactic_new_candidate(
     legacydr8_zphot: pd.Series,
     firstDiaSourceMjdTaiFink: pd.Series,
 ) -> pd.Series:
-    """Select LSST alerts new (< 48h first apparition) and potentially extragalactic
+    """Select LSST alerts new (< 5days first apparition), bright (mag < 24), potentially extragalactic with a fading or rising rate passing the cuts
 
     Notes
     -----
-    Based on an extragalactic block and time cut
-
-
+    Based on an extragalactic block, time cut, sampling cut, and rate cut.
+    Rising alerts must have rate < -0.2 mag/day and last less than 3 days.
+    Fading alerts must have rate > 0.2 mag/day in r/i bands, or > 0.5 mag/day in g/u bands.
 
     Parameters
     ----------
     diaSource: pd.DataFrame
         Full diaSource section of an alert (dictionary exploded)
     diaObject: pd.DataFrame
         Full diaObject section of an alert (dictionary exploded)
+    prvDiaSources: pd.Series
+        Series of lists of previous diaSource dicts for each alert
     simbad_otype: pd.Series
         Series containing labels from `xm.simbad_otype`
     mangrove_lum_dist: pd.Series
@@ -95,10 +150,42 @@ def extragalactic_new_candidate(
         legacydr8_zphot,
     )
 
-    # 48h maximum
-    f_new = (diaSource.midpointMjdTai - firstDiaSourceMjdTaiFink) < 2.0
+    # 5 days maximum
+    f_new = (diaSource.midpointMjdTai - firstDiaSourceMjdTaiFink) < 5.0
+
+    f_bright = fu.flux_to_apparent_mag(diaSource.psfFlux) < 24
+
+    # Minimum 2 points
+
+    f_sampling = has_two_points_same_band(diaSource, diaObject)
+
+    # Check rising and fading rate
+    df = pd.concat([diaSource, prvDiaSources.rename("prvDiaSources")], axis=1)
+    prev = df.apply(get_latest_source_same_band, axis=1)
 
-    f_extragalactic_new = f_extragalactic_near_galaxy & f_new
+    prev_flux = prev.apply(lambda s: s["psfFlux"] if s is not None else float("nan"))
+    prev_time = prev.apply(
+        lambda s: s["midpointMjdTai"] if s is not None else float("nan")
+    )
+
+    delta_mag = fu.flux_to_apparent_mag(diaSource.psfFlux) - fu.flux_to_apparent_mag(
+        prev_flux
+    )
+    delta_time = diaSource.midpointMjdTai - prev_time
+    delta_time_rising = diaSource.midpointMjdTai - firstDiaSourceMjdTaiFink
+
+    rate = delta_mag / delta_time
+
+    f_rising = (rate < -0.2) & (delta_time_rising < 3)
+
+    f_fading_ri = ((diaSource.band == "r") | (diaSource.band == "i")) & (rate > 0.2)
+    f_fading_gu = ((diaSource.band == "g") | (diaSource.band == "u")) & (rate > 0.5)
+
+    f_rate = f_rising | f_fading_ri | f_fading_gu
+
+    f_extragalactic_new = (
+        f_extragalactic_near_galaxy & f_new & f_sampling & f_bright & f_rate
+    )
 
     return f_extragalactic_new