feat: Add new alarm forecast and temporal evaluation functions and tutorials with updated data and plotting utilities.

Author: KennyGraham1

csep/__init__.py

@@ -54,6 +54,8 @@
     'load_gridded_forecast',
     'load_catalog_forecast',
     'load_alarm_forecast',
+    'load_temporal_forecast',
+    'compute_temporal_observations',
     'utc_now_datetime',
     'strptime_to_utc_datetime',
     'datetime_to_utc_epoch',
@@ -640,3 +642,111 @@ def load_alarm_forecast(fname, name=None, score_field='alarm_score',
 
     return forecast
 
+
+def load_temporal_forecast(fname, time_col='time', probability_col='probability',
+                           delimiter=',', **kwargs):
+    """ Load temporal probability forecast from CSV file.
+
+    This function loads time-series probability forecasts where each row represents
+    a time window (e.g., day) with an associated probability of earthquake occurrence.
+    Observations should be computed separately from a catalog using
+    compute_temporal_observations().
+
+    Args:
+        fname (str): Path to temporal forecast CSV file
+        time_col (str): Name of time/index column (default: 'time')
+        probability_col (str): Name of probability column (default: 'probability')
+        delimiter (str): CSV delimiter character. Default: ','
+
+    Returns:
+        dict: Dictionary containing:
+            - 'times': Time indices (numpy array or DatetimeIndex)
+            - 'probabilities': Forecast probabilities (numpy array)
+            - 'metadata': Additional information (dict)
+
+    Raises:
+        FileNotFoundError: If the CSV file does not exist
+        ValueError: If required columns are missing
+
+    Example:
+        >>> # Load temporal forecast
+        >>> data = csep.load_temporal_forecast('daily_forecast.csv')
+        >>> catalog = csep.load_catalog('events.csv')
+        >>>
+        >>> # Compute observations from catalog
+        >>> observations = csep.compute_temporal_observations(
+        ...     catalog, data['times'], magnitude_min=4.0
+        ... )
+        >>>
+        >>> # Evaluate with temporal Molchan diagram
+        >>> from csep.utils import plots
+        >>> ax, ass, sigma = plots.plot_temporal_Molchan_diagram(
+        ...     data['probabilities'], observations, name='DailyM4+'
+        ... )
+    """
+    # Sanity checks
+    if not os.path.exists(fname):
+        raise FileNotFoundError(
+            f"Could not locate file {fname}. Unable to load temporal forecast.")
+
+    times, probabilities, metadata = readers.temporal_forecast_csv(
+        filename=fname,
+        time_col=time_col,
+        probability_col=probability_col,
+        delimiter=delimiter
+    )
+
+    return {
+        'times': times,
+        'probabilities': probabilities,
+        'metadata': metadata
+    }
+
+
+def compute_temporal_observations(catalog, times, magnitude_min=None,
+                                   start_time=None, time_delta=None):
+    """ Compute binary observations from a catalog for temporal forecast evaluation.
+
+    This function counts whether one or more earthquakes occurred in each time
+    window, creating a binary observation vector suitable for temporal ROC
+    and Molchan evaluation.
+
+    Args:
+        catalog: A CSEPCatalog object containing observed earthquakes
+        times (array-like): Array representing forecast time windows.
+                           Can be datetime objects OR integer indices (1, 2, 3, ...).
+                           If integers, use start_time and time_delta.
+        magnitude_min (float, optional): Minimum magnitude threshold.
+        start_time (str or datetime, optional): Start time of the forecast.
+                                               Required if times are integer indices.
+                                               Example: '2024-01-01'
+        time_delta (str, optional): Duration of each time window.
+                                   Required if times are integer indices.
+                                   Examples: '1D' (1 day), '1H' (1 hour), '7D' (1 week)
+
+    Returns:
+        numpy.ndarray: Binary array (1 if event occurred in window, 0 otherwise)
+
+    Example:
+        >>> # With datetime times in CSV
+        >>> data = csep.load_temporal_forecast('forecast.csv')
+        >>> catalog = csep.load_catalog('events.csv')
+        >>> observations = csep.compute_temporal_observations(
+        ...     catalog, data['times'], magnitude_min=4.0
+        ... )
+        >>>
+        >>> # With integer indices (1, 2, 3, ...)
+        >>> observations = csep.compute_temporal_observations(
+        ...     catalog, data['times'],
+        ...     start_time='2024-01-01',
+        ...     time_delta='1D',
+        ...     magnitude_min=4.0
+        ... )
+    """
+    return readers.compute_temporal_observations(
+        catalog=catalog,
+        times=times,
+        magnitude_min=magnitude_min,
+        start_time=start_time,
+        time_delta=time_delta
+    )