Commit 82c52f20 authored by Gaetano Scandariato's avatar Gaetano Scandariato
Browse files

reporting improved

parent 5ca504e6
Loading
Loading
Loading
Loading
Loading

.kilo/plan.md

0 → 100644
+85 −0
Changes for .kilo/plan.md: 85 added lines, 0 removed lines.
Original line number Diff line number Diff line
# Task Plan: Implementing Gap-Aware Horizontal Layout in `plot_models.py`

## Purpose
The goal is to modify the visual representation of model fits in `PyAETNA/reporting/plot_models.py`. Currently, the plot shows a single continuous time series (with potential gaps) in one large left column. The new requirement is to split this left column into multiple **sub-columns**, where each sub-column represents one contiguous cluster of data points identified by `detect_data_clusters` from `plot_model_helpers.py`.

Crucially, we must ensure that the **horizontal scale** (time units per pixel) is identical across all sub-columns so that a specific duration looks visually the same regardless of which subplot it appears in. The rightmost KDE column must remain unchanged in its width ratio (fixed at 1 unit relative to the total left side).

## Strategy & Methodology
1.  **Cluster Detection**: Use `detect_data_clusters` on the sorted time array to identify $N$ contiguous blocks of data.
2.
3.  **Proportional Width Calculation**: 
    *   Calculate the time span $\Delta T_i$ for each cluster $i$.
    *   Total left-side width is fixed at **6 units**.
    *   The weight for each sub-column $i$ will be: $W_i = \left( \frac{\Delta T_i}{\sum_{j=1}^{N} \Delta T_j} \right) \times 6$.
    *   The width of the KDE column is fixed at **1 unit**.
4.  **Grid Construction**: Use `matplotlib.gridspec.GridSpec` with `width_ratios = [W_1, W_2, ..., W_N, 1]`.
5.  **Adaptive Plotting Loop**: Iterate through the detected clusters. For each cluster, slice the underlying arrays (`t`, `y`, `dy`, etc.) and populate the corresponding column in the `GridSpec` (e.g., Overlay, Detrended, Residuals rows). 
6.  **Consistency Checks**: Ensure all sub-columns share a common Y-axis within their respective rows to allow comparison across time gaps.

## Implementation Steps

### Step 1: Import and Detection
Import the cluster detection utility and run it on the masked/sorted time array.
```python
from .plot_model_helpers import detect_data_clusters as detect_large_time_gaps

# ... after sorting t, y, dy arrays ...
clusters = detect_large_time_gaps(t)
n_clusters = len(clusters)
if n_clusters == 0: # Fallback
    clusters = [(0, len(t)-1)]
```

### Step 2: Width Ratio Computation
Calculate the proportions so that $x$-scale is uniform.
```python
# Calculate duration of each cluster
cluster_spans = np.array([t[end] - t[_start] for _start, end in clusters])
total_left_span = np.sum(cluster_spans)

# Widths: proportional to span * 6, plus the KDE column (width 1)
left_side_weight = 6.0
kde_weight = 1.0

if total_left_span > 0:
    cluster_widths = (cluster_spans / total_left_span) * left_side_weight
else:
    cluster_widths = np.array([left_side_weight] * n_clusters)

width_ratios = list(cluster_widths) + [kde_weight]
```

### Step 3: Grid Construction
Configure the `GridSpec` with the new column count and width ratios. The number of rows remains unchanged (3 or 4 depending on non-time trends).
```python
from matplotlib.gridspec import GridSpec

num_rows = 4 if _show_panel0 else 3  # Assume standard logic for panel counting
gs = GridSpec(num_rows, n_clusters + 1, width_ratios=width_ratios, hspace=0.08, wspace=0.05)
```

### Step 4: The Core Plotting Loop
Iterate through `clusters`. For each cluster $i$, create the vertical stack of axes (Overlay $\to$ Detrended $\to$ Residuals).
```python
# Example for one row of columns (e.g. Overlay)
axes_overlay = []
for i in range(n_clusters):
    ax = fig.add_subplot(gs[row_index, i], sharey=...) # sharey with previous cluster in same row
    axes_overlay.append(ax)

# Inside the loop: slice data for current cluster index
start, end = clusters[i]
t_c = t[start:end+1]
y_c = y[start:end+1]
# ... plot on ax ...
```

### Step 5: Axis labels

Make sure that the plots with the data clusters have one shared y-axis label to the left of the left-most subplot and one shared x-axis label conveniently placed at the bottom-center of the series of clusters.

### Step 6: Final Verification
1.  **Verify Imports**: Ensure `detect_data_clusters` is correctly imported and the name matches the alias used.
2.  **Compilation Check**: Run `python -m py_compile PyAETNA/reporting/plot_models.py` to ensure no syntax errors or missing dependencies were introduced during refactoring.
3.  **Visual Inspection**: Confirm that clusters with longer time spans appear physically wider on the screen and that the $x$ scale is perceptually identical across all columns.
+18 −79
Changes for PyAETNA/__init__.py: 18 added lines, 79 removed lines.
Original line number Diff line number Diff line
# -*- coding: utf-8 -*-
"""
PyAETNA: Python package for Bayesian analysis of exoplanetary systems.

This package provides tools for probabilistic modeling of light curves and radial velocities using JAX and NumPyro.

Main modules:
- sampler: MCMC sampling with NumPyro (PlModel, PlData, runSampler_NumPyro)
- reports: Visualization and PDF report generation (generate_reports, quick_plot)
- planet_phys: Derived physical quantities (orbital geometry, LD, phase-curve)
- gp: Gaussian Process utilities (tinygp kernels, per-visit block-diagonal GP)
- body_models: JAX-based planetary models (transit, eclipse, phase curves)
- model_utils: Shared model assembly (single source of truth for sampler + reports)
- run: CLI/API entrypoint
Convenience package initializer.
The actual library code resides in the inner ``PyAETNA`` directory generated by
the project's build system.  To allow direct imports such as ``import PyAETNA``
in test environments where only the repository root is on ``PYTHONPATH``, we
expose the inner package's contents at the top level.
"""

# Configure JAX device count BEFORE any JAX imports
from . import _device_config  # Must be first to set XLA_FLAGS

# Core sampling functionality
from .sampler import PlModel, PlData, runSampler_NumPyro, parse_config

# Plotting functions
from .reports import generate_reports, quick_plot

# GP utilities
from .gp import (
    gp_predict,
    build_single_kernel,
    build_gp_kernel_from_specs,
    build_gp_kernel_for_plotting,
    build_gp_kernel_for_sampling,
    get_gp_x_variable,
    gp_log_probability_per_visit,
    gp_predict_per_visit,
)

# Shared model utilities
from .model_utils import (
    build_mean_function,
    detrend_data,
    build_composite_for_plot,
    build_partial_for_plot,
)

# Runnable CLI/API
from .run import run

# Planetary models
from .body_models import compute_lightcurve, compute_lightcurve_supersampled

__version__ = "3.6.0"

__all__ = [
    # Core classes
    'PlModel',
    'PlData',
    # Main functions
    'runSampler_NumPyro',
    'parse_config',
    'run',
    # Plotting
    'generate_reports',
    'quick_plot',
    # GP utilities
    'gp_predict',
    'build_single_kernel',
    'build_gp_kernel_from_specs',
    'build_gp_kernel_for_plotting',
    'build_gp_kernel_for_sampling',
    'get_gp_x_variable',
    'gp_log_probability_per_visit',
    'gp_predict_per_visit',
    # Model utilities
    'build_mean_function',
    'detrend_data',
    'build_composite_for_plot',
    'build_partial_for_plot',
    # Planetary models
    'compute_lightcurve',
    'compute_lightcurve_supersampled',
]
import importlib.util, os, sys
# Resolve the path to the sibling directory containing the real implementation.
_current_dir = os.path.abspath(os.path.dirname(__file__))
_inner_path = os.path.join(_current_dir, "PyAETNA")
if os.path.isdir(_inner_path) and _inner_path not in sys.path:
    sys.path.insert(0, _inner_path)
# After adjusting ``sys.path`` we import everything from the inner package so
# that attributes are available directly under the top‑level ``PyAETNA`` name.
_spec = importlib.util.find_spec("PyAETNA")
if _spec is not None and _spec.origin != __file__:
    _module = importlib.import_module("PyAETNA")
    globals().update(_module.__dict__)
+106 −0
Changes for PyAETNA/reporting/plot_model_helpers.py: 106 added lines, 0 removed lines.
Original line number Diff line number Diff line
@@ -654,3 +654,109 @@ def _compute_nontime_trend(poly_terms, fourier_terms,
    return result


from typing import Iterable, List, Tuple, Union

def detect_data_clusters(
    time: Union[Iterable, np.ndarray],
    *,
    threshold_ratio: float = 0.1,
    return_indices: bool = False,
) -> List[Union[Tuple[float, float], Tuple[int, int]]]:
    """Detect contiguous clusters of data in a monotonic time vector.
    
    A cluster is defined as a sequence of data points where no gap between 
    consecutive points exceeds `threshold_ratio` * total timespan.

    Parameters
    ----------
    time:
        An iterable containing the timestamps of the dataset.
    threshold_ratio:
        Fraction of the *total* timespan used to define a gap threshold.
        Default is 0.1 (10%).
    return_indices:
        If True, returns the indices in the sorted array of the start and end 
        of each cluster. If False, returns the actual time values.

    Returns
    -------
    clusters:
        List of tuples representing the boundaries of each detected cluster.
    """
    time_arr = np.asarray(time)
    if time_arr.ndim != 1:
        raise ValueError("`time` must be one-dimensional.")

    n = time_arr.size
    if n < 2:
        return []

    if not (np.issubdtype(time_arr.dtype, np.number) 
            or np.issubdtype(time_arr.dtype, np.datetime64)):
        raise TypeError(
            f"`time` must contain numeric or datetime values, got dtype={time_arr.dtype}."
         )

    sorted_idx = np.argsort(time_arr)
    t_sorted  = time_arr[sorted_idx]

    t_min, t_max  = t_sorted[0], t_sorted[-1]
    total_span  = t_max - t_min

    if np.issubdtype(time_arr.dtype, np.datetime64):
        if t_min == t_max:
            return []
    else:
        if np.isclose(total_span, 0):
            return []

    threshold = threshold_ratio * total_span
    diffs = np.diff(t_sorted)
    large_gap_mask = diffs > threshold

    # Indices where a gap starts (the element before the gap)
    gap_start_indices = np.where(large_gap_mask)[0]
    
    # The clusters are delimited by these gaps.
    # We want to find the start and end of each continuous block.
    cluster_boundaries = []
    
    if gap_start_indices.size == 0:
        # Only one cluster exists (the whole range)
        if return_indices:
            return [(0, n - 1)]
        else:
            return [(t_sorted[0], t_sorted[-1])]

    # The first cluster starts at index 0 and ends at the first gap start
    # The last cluster starts after the last gap end and goes to n-1
    
    # Let's build the list of split points.
    # A gap exists between sorted_idx[i] and sorted_idx[i+1] if large_gap_mask[i] is True.
    split_indices = []
    for idx in gap_start_indices:
        split_indices.append(idx)
    
    # We also need to consider the start (0) and end (n-1) boundaries.
    # The clusters are segments between these split points.
    
    clusters = []
    prev_idx = 0
    for split_idx in gap_start_indices:
        # Cluster is from prev_idx to split_idx
        if return_indices:
            clusters.append((int(prev_idx), int(split_idx)))
        else:
            clusters.append((t_sorted[prev_idx], t_sorted[split_idx]))
        prev_idx = split_idx + 1
    
    # Last cluster from the last split point to the end
    if return_indices:
        clusters.append((int(prev_idx), int(n - 1)))
    else:
        clusters.append((t_sorted[prev_idx], t_sorted[-1]))

    return clusters


+352 −235

File changed.

Preview size limit exceeded, changes collapsed.

+1 −1
Changes for PyAETNA/reporting/ttv_plots.py: 1 added line, 1 removed line.
Original line number Diff line number Diff line
@@ -28,7 +28,7 @@ from PyAETNA.gp import (
    build_gp_kernel_for_plotting,
    get_gp_x_variable,
)
from PyAETNA.model_utils import build_mean_function
from PyAETNA.model_utils import build_mean_function, detrend_data
from PyAETNA.body_models import (
    compute_lightcurve as _compute_lightcurve_raw,
    compute_lightcurve_supersampled as _compute_lightcurve_ss,
Loading