defmodule Cepheid.Normalization.Magnitude do
  @moduledoc """
  Functions for normalizing the magnitude of light curves.
  """

  import Nx.Defn

  @doc """
  Computes range- and phase-normalized visual magnitudes using a lambertian sphere phase function.

  ## Parameters
  - `observed_magnitudes`: Nx.Tensor of observed apparent magnitudes
  - `observed_phases`: Nx.Tensor of phase angles in radians
  - `observed_ranges`: Nx.Tensor of ranges
  - `std_range`: Standard range for normalization (same units as observed_ranges)
  - `std_phase`: Standard phase for normalization in radians

  ## Returns
  - Nx.Tensor of normalized magnitudes
  """
  defn normalize_to_range_and_phase(
         observed_magnitudes,
         observed_phases,
         observed_ranges,
         std_range,
         std_phase
       ) do
    # lambertian phase function: (1 + cos(α)) / 2
    phase_factor_obs = (1.0 + Nx.cos(observed_phases)) / 2.0
    phase_factor_std = (1.0 + Nx.cos(std_phase)) / 2.0

    # phase correction: -2.5 * log10(phase_factor_obs / phase_factor_std)
    phase_correction =
      phase_factor_obs
      |> Nx.divide(phase_factor_std)
      |> Nx.log10()
      |> Nx.multiply(-2.5)

    # range correction: 5 * log10(range / standard_range)
    range_correction =
      observed_ranges
      |> Nx.divide(std_range)
      |> Nx.log10()
      |> Nx.multiply(5.0)

    # normalized magnitude = observed - phase_correction - range_correction
    observed_magnitudes
    |> Nx.subtract(phase_correction)
    |> Nx.subtract(range_correction)
  end
end