From 090b557a6fc8739886e5259ee9329c39baa3ec43 Mon Sep 17 00:00:00 2001 From: Duncan Date: Mon, 29 Apr 2024 15:54:56 +0200 Subject: [PATCH] Attempts to fix #27, currently fails due to keyword arguments in interp_kind. --- zodipy/_interpolate_source.py | 14 +++++++------- zodipy/zodipy.py | 14 ++++---------- 2 files changed, 11 insertions(+), 17 deletions(-) diff --git a/zodipy/_interpolate_source.py b/zodipy/_interpolate_source.py index 6b981fa..03a14af 100644 --- a/zodipy/_interpolate_source.py +++ b/zodipy/_interpolate_source.py @@ -27,13 +27,13 @@ def get_source_parameters_kelsall_comp( if not bandpass.frequencies.unit.is_equivalent(model.spectrum.unit): bandpass.switch_convention() - spectrum = ( + center_freqs = ( model.spectrum.to_value(u.Hz) if model.spectrum.unit.is_equivalent(u.Hz) else model.spectrum.to_value(u.micron) ) - interpolator = partial(interpolator, x=spectrum) + interpolator = partial(interpolator, x=center_freqs) source_parameters: dict[ComponentLabel | str, dict[str, Any]] = {} for comp_label in model.comps: @@ -52,11 +52,11 @@ def get_source_parameters_kelsall_comp( source_parameters[comp_label]["albedo"] = albedo if model.phase_coefficients is not None: - phase_coefficients = interpolator(y=np.asarray(model.phase_coefficients))( - bandpass.frequencies.value - ) - phase_coefficients = interpolator(y=np.asarray(model.phase_coefficients))( - bandpass.frequencies.value + phase_coefficients = np.array( + [ + interpolator(y=np.asarray(model.phase_coefficients[i]))(bandpass.frequencies.value) + for i in range(len(model.phase_coefficients)) + ] ) else: phase_coefficients = np.repeat(np.zeros((3, 1)), repeats=bandpass.frequencies.size, axis=-1) diff --git a/zodipy/zodipy.py b/zodipy/zodipy.py index f7ed6b5..4e4b548 100644 --- a/zodipy/zodipy.py +++ b/zodipy/zodipy.py @@ -9,7 +9,7 @@ import healpy as hp import numpy as np from astropy.coordinates import solar_system_ephemeris -from scipy.interpolate import interp1d +from scipy.interpolate import CubicSpline from zodipy._bandpass import get_bandpass_interpolation_table, validate_and_get_bandpass from zodipy._constants import SPECIFIC_INTENSITY_UNITS @@ -45,10 +45,6 @@ class Zodipy: Defaults to DIRBE. gauss_quad_degree (int): Order of the Gaussian-Legendre quadrature used to evaluate the line-of-sight integral in the simulations. Default is 50 points. - interp_kind (str): Interpolation kind used to interpolate relevant model parameters. - Defaults to 'linear'. For more information on available interpolation methods, - please visit the [Scipy documentation]( - https://docs.scipy.org/doc/scipy/reference/generated/scipy.interpolate.interp1d.html). extrapolate (bool): If `True` all spectral quantities in the selected model are extrapolated to the requested frequencies or wavelengths. If `False`, an exception is raised on requested frequencies/wavelengths outside of the @@ -72,7 +68,6 @@ def __init__( model: str = "dirbe", gauss_quad_degree: int = 50, extrapolate: bool = False, - interp_kind: str = "linear", ephemeris: str = "de432s", solar_cut: u.Quantity[u.deg] | None = None, solar_cut_fill_value: float = np.nan, @@ -81,16 +76,15 @@ def __init__( self.model = model self.gauss_quad_degree = gauss_quad_degree self.extrapolate = extrapolate - self.interp_kind = interp_kind self.ephemeris = ephemeris self.solar_cut = solar_cut.to(u.rad) if solar_cut is not None else solar_cut self.solar_cut_fill_value = solar_cut_fill_value self.n_proc = n_proc self._interpolator = partial( - interp1d, - kind=self.interp_kind, - fill_value="extrapolate" if self.extrapolate else np.nan, + CubicSpline, + extrapolate=True, + bc_type="natural", ) self._ipd_model = model_registry.get_model(model) self._gauss_points_and_weights = np.polynomial.legendre.leggauss(gauss_quad_degree)