diff --git a/docs/conf.py b/docs/conf.py
index 197fc66..0461c65 100644
--- a/docs/conf.py
+++ b/docs/conf.py
@@ -144,6 +144,7 @@
     # Other articles
     "ai_2023": ("https://arxiv.org/abs/2303.10171%s", "Ai et al., 2023%s"),
     "borsanyi_2016": ("https://arxiv.org/abs/1606.07494%s", "Borsanyi et al., 2016%s"),
+    "smith_2019": ("https://arxiv.org/abs/1908.00546%s", "Smith & Caldwell, 2019%s"),
     "caprini_2020": ("https://arxiv.org/abs/1910.13125%s", "Caprini et al., 2020%s"),
     "giese_2020": ("https://arxiv.org/abs/2004.06995%s", "Giese et al., 2020%s"),
     "giese_2021": ("https://arxiv.org/abs/2010.09744%s", "Giese et al., 2021%s"),
diff --git a/examples/const_cs_gw.py b/examples/const_cs_gw.py
index 037c3fe..8d0c365 100644
--- a/examples/const_cs_gw.py
+++ b/examples/const_cs_gw.py
@@ -54,6 +54,9 @@ def main():
             spectrum_kwargs={
                 "r_star": r_star
                 # "z": z
+                # "Tn": 100,
+                # "g_star": 100,
+                # "gs_star": 100
             }
             # bubble_kwargs={"allow_invalid": False}, allow_bubble_failure=True
         )
@@ -150,9 +153,9 @@ def main():
     fig.tight_layout()
     fig2.tight_layout()
     fig3.tight_layout()
-    utils.save_and_show(fig, "const_cs_gw.png")
-    utils.save_and_show(fig2, "const_cs_gw_v.png")
-    utils.save_and_show(fig3, "const_cs_gw_omgw0.png")
+    utils.save_and_show(fig, "const_cs_gw")
+    utils.save_and_show(fig2, "const_cs_gw_v")
+    utils.save_and_show(fig3, "const_cs_gw_omgw0")
 
 
 if __name__ == "__main__":
diff --git a/examples/noise.py b/examples/noise.py
new file mode 100644
index 0000000..f8f117e
--- /dev/null
+++ b/examples/noise.py
@@ -0,0 +1,44 @@
+from matplotlib import pyplot as plt
+import numpy as np
+
+from examples.utils import save_and_show
+from pttools.omgw0 import noise
+
+
+def main():
+    fig: plt.Figure = plt.figure()
+    axs = fig.subplots(2, 2)
+
+    f = np.logspace(-4, -1, 50)
+    ax1 = axs[0, 0]
+    P_oms = noise.P_oms()
+    ax1.plot([f.min(), f.max()], [P_oms, P_oms], label=r"$P_\text{oms}$")
+    ax1.plot(f, noise.P_acc(f), label=r"$P_\text{acc}$")
+    ax1.set_ylabel("$P(f)$")
+
+    ax2 = axs[0, 1]
+    ax2.plot(f, noise.S_AE(f), label="$S_{A,E}$")
+    ax2.plot(f, noise.S_AE_approx(f), label=r"$S_{A,E,\text{approx}}$")
+    ax2.plot(f, noise.S_gb(f), label=r"$S_\text{gb}$")
+    ax2.set_ylabel("$S(f)$")
+
+    ax3 = axs[1, 0]
+    ax3.plot(f, noise.omega_ins(f), label=r"$\Omega_\text{ins}$")
+    ax3.plot(f, noise.omega_eb(f), label=r"$\Omega_\text{eb}$")
+    ax3.plot(f, noise.omega_gb(f), label=r"$\Omega_\text{gb}$")
+    ax3.plot(f, noise.omega_noise(f), label=r"$\Omega_\text{noise}$")
+    ax3.set_ylabel(r"$\Omega(f)$")
+
+    for ax in axs.flat:
+        ax.set_xlabel(r"$f(\text{Hz})$")
+        ax.set_xscale("log")
+        ax.set_yscale("log")
+        ax.legend()
+    fig.tight_layout()
+
+    return fig
+
+
+if __name__ == '__main__':
+    fig = main()
+    save_and_show(fig, "noise")
diff --git a/pttools/models/bag.py b/pttools/models/bag.py
index 2e00535..7ebf7df 100644
--- a/pttools/models/bag.py
+++ b/pttools/models/bag.py
@@ -35,6 +35,7 @@ class BagModel(AnalyticModel):
     DEFAULT_LABEL_LATEX = "Bag model"
     DEFAULT_LABEL_UNICODE = DEFAULT_LABEL_LATEX
     DEFAULT_NAME = "bag"
+    TEMPERATURE_IS_PHYSICAL = False
 
     def __init__(
             self,
diff --git a/pttools/models/base.py b/pttools/models/base.py
index 37a1a08..ce8e597 100644
--- a/pttools/models/base.py
+++ b/pttools/models/base.py
@@ -13,7 +13,10 @@
 
 
 class BaseModel(abc.ABC):
-    """The base for both Model and ThermoModel"""
+    """The base for both Model and ThermoModel
+
+    All temperatures must be in units of GeV for the frequency conversion in Spectrum to work.
+    """
     DEFAULT_LABEL_LATEX: str = None
     DEFAULT_LABEL_UNICODE: str = None
     DEFAULT_NAME: str = None
@@ -21,6 +24,9 @@ class BaseModel(abc.ABC):
     DEFAULT_T_MIN: float = 1e-3
     DEFAULT_T_MAX: float = np.inf
 
+    #: Whether the temperature is in proper physics units
+    TEMPERATURE_IS_PHYSICAL: bool = None
+
     def __init__(
             self,
             name: str = None,
@@ -29,13 +35,15 @@ def __init__(
             label_latex: str = None,
             label_unicode: str = None,
             gen_cs2: bool = True,
-            gen_cs2_neg: bool = True):
+            gen_cs2_neg: bool = True,
+            temperature_is_physical: bool = None):
         self.name = self.DEFAULT_NAME if name is None else name
         self.label_latex = self.DEFAULT_LABEL_LATEX if label_latex is None else label_latex
         self.label_unicode = self.DEFAULT_LABEL_UNICODE if label_unicode is None else label_unicode
         self.T_min = self.DEFAULT_T_MIN if T_min is None else T_min
         self.T_max = self.DEFAULT_T_MAX if T_max is None else T_max
         self.restrict_to_valid = restrict_to_valid
+        self.temperature_is_physical = self.TEMPERATURE_IS_PHYSICAL if temperature_is_physical is None else temperature_is_physical
 
         if self.name is None:
             raise ValueError("The model must have a name.")
diff --git a/pttools/models/const_cs.py b/pttools/models/const_cs.py
index 7d2d9cd..f4763e6 100644
--- a/pttools/models/const_cs.py
+++ b/pttools/models/const_cs.py
@@ -29,6 +29,7 @@ class ConstCSModel(AnalyticModel):
     DEFAULT_LABEL_LATEX = "Constant $c_s$ model"
     DEFAULT_LABEL_UNICODE = "Constant cₛ model"
     DEFAULT_NAME = "const_cs"
+    TEMPERATURE_IS_PHYSICAL = False
 
     def __init__(
             self,
diff --git a/pttools/models/const_cs_thermo.py b/pttools/models/const_cs_thermo.py
index 79fba0d..8ed695e 100644
--- a/pttools/models/const_cs_thermo.py
+++ b/pttools/models/const_cs_thermo.py
@@ -13,6 +13,7 @@ class ConstCSThermoModel(ThermoModel):
     DEFAULT_LABEL_LATEX = "Constant $c_s$ thermo-model"
     DEFAULT_LABEL_UNICODE = "Constant cₛ thermo-model"
     DEFAULT_NAME = "const_cs_thermo"
+    TEMPERATURE_IS_PHYSICAL = False
 
     GEFF_DATA_LOG_TEMP = np.linspace(-1, 3, 1000)
     GEFF_DATA_TEMP = 10**GEFF_DATA_LOG_TEMP
diff --git a/pttools/models/full.py b/pttools/models/full.py
index b49cfaf..b238b3f 100644
--- a/pttools/models/full.py
+++ b/pttools/models/full.py
@@ -49,7 +49,9 @@ def __init__(
             V_s=V_s, V_b=V_b,
             T_min=thermo.t_min, T_max=thermo.t_max,
             name=name, label_latex=label_latex, label_unicode=label_unicode,
-            gen_critical=False, gen_cs2=False, gen_cs2_neg=False, implicit_V=True)
+            gen_critical=False, gen_cs2=False, gen_cs2_neg=False, implicit_V=True,
+            temperature_is_physical=thermo.TEMPERATURE_IS_PHYSICAL
+        )
 
         self.temp_spline_s = splrep(
             np.log10(self.w(self.thermo.GEFF_DATA_TEMP, Phase.SYMMETRIC)), self.thermo.GEFF_DATA_LOG_TEMP
diff --git a/pttools/models/model.py b/pttools/models/model.py
index 8389a4d..ff4bdb5 100644
--- a/pttools/models/model.py
+++ b/pttools/models/model.py
@@ -43,6 +43,7 @@ def __init__(
             gen_cs2: bool = True,
             gen_cs2_neg: bool = True,
             implicit_V: bool = False,
+            temperature_is_physical: bool = None,
             allow_invalid: bool = False):
 
         if implicit_V:
@@ -62,6 +63,13 @@ def __init__(
             # if V_s == V_b:
             #     logger.warning("The bubble will not expand, when V_s <= V_b. Got: V_b = V_s = %s.", V_s)
 
+        self.temperature_is_physical = self.TEMPERATURE_IS_PHYSICAL if temperature_is_physical is None else temperature_is_physical
+        if self.temperature_is_physical is None:
+            raise ValueError(
+                "It has not been specified whether the temperature scale for the model is physical. "
+                "Please specify it in the model definition."
+            )
+
         self.T_ref: float = T_ref
         self.V_s: float = V_s
         self.V_b: float = V_b
diff --git a/pttools/models/sm.py b/pttools/models/sm.py
index 7a1d658..e0d9982 100644
--- a/pttools/models/sm.py
+++ b/pttools/models/sm.py
@@ -21,6 +21,8 @@ class StandardModel(ThermoModel):
     DEFAULT_LABEL_LATEX = "Standard Model"
     DEFAULT_LABEL_UNICODE = DEFAULT_LABEL_LATEX
     DEFAULT_NAME = "standard_model"
+    TEMPERATURE_IS_PHYSICAL = True
+
     # Copied from the ArXiv file som_eos.tex
     GEFF_DATA = np.array([
         [0.00, 10.71, 1.00228],
diff --git a/pttools/omgw0/noise.py b/pttools/omgw0/noise.py
index f206ed7..b2ebfb5 100644
--- a/pttools/omgw0/noise.py
+++ b/pttools/omgw0/noise.py
@@ -13,7 +13,8 @@ def signal_to_noise_ratio(
     $$\rho = \sqrt{T_{\text{obs}} \int_{f_\text{min}}^{f_\text{max}} df \frac{
     h^2 \Omega_{\text{gw},0}^2}{
     h^2 \Omega_{\text{n}}^2}}
-    :gowling_2021:`\ ` eq. 3.12
+    :gowling_2021:`\ ` eq. 3.12,
+    :smith_2019:`\ ` p. 1
 
     :param f: frequencies
     :param signal: signal array
@@ -27,22 +28,33 @@ def signal_to_noise_ratio(
 def ft(L: th.FloatOrArr = const.LISA_ARM_LENGTH) -> th.FloatOrArr:
     r"""Transfer frequency
     $$f_t = \frac{c}{2\pi L}$$
+    :gowling_2021:`\ ` p. 12
     """
     return const.c / (2*np.pi*L)
 
+FT_LISA = ft()
+_ft_func = ft
 
-# def N_AE(f: th.FloatOrArr, L: th.FloatOrArr = const.LISA_ARM_LENGTH) -> th.FloatOrArrNumba:
-#     r"""A and E channels of LISA instrument noise
-#     $$N_A = N_E = ...$$
-#     :gowling_2021:`\ ` eq. 3.4
-#     """
-#     f_frac = f/ft
-#     cos_f_frac = np.cos(f_frac)
-#     W = 1 - np.exp(-2j*f_frac)
-#     return ((4 + 2*cos_f_frac)*P_oms(L) + 8*(1 + cos_f_frac + cos_f_frac**2 * P_acc(f, L))) * np.abs(W)**2
+
+def N_AE(f: th.FloatOrArr, ft: th.FloatOrArr = FT_LISA, L: th.FloatOrArr = const.LISA_ARM_LENGTH, W_abs2: th.FloatOrArr = None) -> th.FloatOrArrNumba:
+    r"""A and E channels of LISA instrument noise
+    $$N_A = N_E = ...$$
+    :gowling_2021:`\ ` eq. 3.4
+    """
+    cos_f_frac = np.cos(f/ft)
+    if W_abs2 is None:
+        W_abs2 = np.abs(W(f, ft))**2
+    return ((4 + 2*cos_f_frac)*P_oms(L) + 8*(1 + cos_f_frac + cos_f_frac**2) * P_acc(f, L)) * W_abs2
 
 
 def omega(f: th.FloatOrArr, S: th.FloatOrArr) -> th.FloatOrArr:
+    r"""Convert an effective noise power spectral density (aka. sensitivity) $S$
+    to a fractional GW energy density power spectrum $\Omega$
+    $$\Omega = \frac{4 \pi^2}{3 H_0^2} f^3 S(f)$$
+    :gowling_2021:`\ ` eq. 3.8,
+    :gowling_2023:`\ ` eq. 3.8,
+    :smith_2019:`\ ` p. 1
+    """
     return 4*np.pi**2 / (3*const.H0_HZ**2) * f**3 * S
 
 
@@ -58,10 +70,11 @@ def omega_ins(f: th.FloatOrArr) -> th.FloatOrArr:
     r"""LISA instrument noise
     $$\Omega_\text{ins} = \left( \frac{4 \pi^2}{3 H_0^2} f^3 S_A(f)$$
     """
-    return omega(f, S_AE(f))
+    return omega(f=f, S=S_AE(f))
 
 
 def omega_noise(f: th.FloatOrArr) -> th.FloatOrArr:
+    r"""$$\Omega_\text{noise} = \Omega_\text{ins} + \Omega_\text{eb} + \Omega_\text{gb}$$"""
     return omega_ins(f) + omega_eb(f) + omega_gb(f)
 
 
@@ -79,13 +92,46 @@ def P_oms(L: th.FloatOrArr = const.LISA_ARM_LENGTH) -> th.FloatOrArr:
     $$P_\text{oms}(f) = \left( \frac{1.5 \cdot 10^{-11} \text{m}}{L} \right)^2 \text{Hz}^{-1}$$
     :gowling_2021:`\ ` eq. 3.2
     This is white noise and therefore independent of the frequency.
+    Note that there is a typo on :gowling_2021:`\ ` p. 12:
+    the correct $L = 2.5 \cdot 10^9 \text{m}$.
     """
     return (1.5e-11 / L)**2
 
 
-def S_AE(f: th.FloatOrArr, L: th.FloatOrArr = const.LISA_ARM_LENGTH, both_channels: bool = True) -> th.FloatOrArr:
+def R_AE(f: th.FloatOrArr, ft: th.FloatOrArr = FT_LISA, W_abs2: th.FloatOrArr = None) -> th.FloatOrArr:
+    r"""Gravitational wave response function for the A and E channels
+    :gowling_2021:`\ ` eq. 3.6
+    """
+    if W_abs2 is None:
+        W_abs2 = np.abs(W(f, ft))**2
+    return 9/20 * W_abs2 / (1 + (3*f/(4*ft))**2)
+
+
+def S(N: th.FloatOrArr, R: th.FloatOrArr) -> th.FloatOrArr:
+    r"""Noise power spectral density
+    $$S = \frac{N}{\mathcal{R}}$$
+    :gowling_2021:`\ ` eq. 3.1
+    """
+    return N / R
+
+
+def S_AE(f: th.FloatOrArr, ft: th.FloatOrArr = FT_LISA, L: th.FloatOrArr = const.LISA_ARM_LENGTH, both_channels: bool = True) -> th.FloatOrArr:
+    r"""Noise power spectral density for the LISA A and E channels
+    $$S_A = S_E = \frac{N_A}{\mathcal{R}_A}$$
+    :gowling_2021:`\ ` eq. 3.7
+    """
+    # The W_abs2 cancels and can therefore be set to unity
+    ret = S(N=N_AE(f=f, ft=ft, L=L, W_abs2=1), R=R_AE(f=f, ft=ft, W_abs2=1))
+    if both_channels:
+        return 1/np.sqrt(2) * ret
+    return ret
+
+
+def S_AE_approx(f: th.FloatOrArr, L: th.FloatOrArr = const.LISA_ARM_LENGTH, both_channels: bool = True) -> th.FloatOrArr:
     r"""Approximate noise power spectral density for the LISA A and E channels
-    $$S_A = S_E$
+    $$S_A = S_E = \frac{N_A}{\mathcal{R}_A}
+    \approx \frac{40}{3} (P_\text{oms} + 4P_\text{acc}) \left( 1 + \frac{3f}{4f_t} \right^2$$
+    :gowling_2021:`\ ` eq. 3.7
     """
     ret = 40/3 * (P_oms(L) + 4*P_acc(f, L)) * (1 + (3*f/(4*ft(L)))**2)
     if both_channels:
@@ -104,3 +150,11 @@ def S_gb(
         d: float = 1680) -> th.FloatOrArr:
     r"""Noise power spectral density for galactic binaries"""
     return A * (1e-3 / f)**(-7/3) * np.exp(-(f/f_ref_gb)**a - b*f*np.sin(c*f)) * (1 + np.tanh(d*(fk - f)))
+
+
+def W(f: th.FloatOrArr, ft: th.FloatOrArr) -> th.FloatOrArr:
+    r"""Round trip modulation
+    $$W(f,f_t) = 1 - e^{-2i \frac{f}{f_t}}$$
+    :gowling_2021:`\ ` p. 12
+    """
+    return 1 - np.exp(-2j * f / ft)
diff --git a/pttools/omgw0/omgw0_ssm.py b/pttools/omgw0/omgw0_ssm.py
index 10186d4..2493c29 100644
--- a/pttools/omgw0/omgw0_ssm.py
+++ b/pttools/omgw0/omgw0_ssm.py
@@ -8,20 +8,70 @@
 """
 
 import functools
+import logging
 
 import numpy as np
 
 from pttools.bubble.boundary import Phase
+from pttools.bubble.bubble import Bubble
 import pttools.bubble.ke_frac_approx as K
 import pttools.omgw0.suppression as sup
 from pttools.ssmtools.const import NPTDEFAULT, NptType
-import pttools.ssmtools.spectrum as ssm
+import pttools.ssmtools as ssm
 import pttools.type_hints as th
 from . import const
 from . import noise
 
+logger = logging.getLogger(__name__)
+
 
 class Spectrum(ssm.SSMSpectrum):
+    r"""A spectrum object that includes $\Omega_{\text{gw},0}$"""
+    def __init__(
+            self,
+            bubble: Bubble,
+            y: np.ndarray = None,
+            z_st_thresh: float = ssm.Z_ST_THRESH,
+            nuc_type: ssm.NucType = ssm.DEFAULT_NUC_TYPE,
+            nt: int = 10000,
+            r_star: float = None,
+            lifetime_multiplier: float = 1,
+            compute: bool = True,
+            Tn: float = None,
+            g_star: float = None,
+            gs_star: float = None
+            ):
+        """
+        :param bubble: the Bubble object
+        :param y: $z = kR*$ array
+        :param z_st_thresh: for $z$ values above z_sh_tresh, use approximation rather than doing the sine transform integral.
+        :param nuc_type: nucleation type
+        :param nt: number of points in the t array
+        :param r_star: $r_*$
+        :param lifetime_multiplier: used for computing the source lifetime factor
+        :param compute: whether to compute the spectrum immediately
+        :param Tn: $T_n$, nucleation temperature override
+        :param g_star: $g_*$, degrees of freedom override at the time of GW production
+        :param gs_star: $g_{s,*}$ degrees of freedom override for entropy at the time of GW production
+        """
+        super().__init__(
+            bubble=bubble,
+            y=y,
+            z_st_thresh=z_st_thresh,
+            nuc_type=nuc_type,
+            nt=nt,
+            r_star=r_star,
+            lifetime_multiplier=lifetime_multiplier,
+            compute=compute
+        )
+        self.override_necessary = not self.bubble.model.temperature_is_physical
+        self.Tn_manual_override = Tn is not None
+        self.g_star_manual_override = g_star is not None
+        self.gs_star_manual_override = gs_star is not None
+        self.Tn_override = 100 if Tn is None else Tn
+        self.g_star_override = 100 if g_star is None else gs_star
+        self.gs_star_override = self.g_star_override if gs_star is None else gs_star
+
     def f(self, z: np.ndarray = None) -> th.FloatOrArr:
         if z is None:
             z = self.y
@@ -30,7 +80,7 @@ def f(self, z: np.ndarray = None) -> th.FloatOrArr:
     @functools.cached_property
     def f_star0(self) -> float:
         return f_star0(
-            T_n=self.bubble.Tn,
+            Tn=self.Tn,
             g_star=self.g_star
         )
 
@@ -39,13 +89,29 @@ def F_gw0(self, g0: float = const.G0, gs0: float = const.GS0) -> float:
             g_star=self.g_star,
             g0=g0,
             gs0=gs0,
-            gs_star=self.bubble.model.gs(w=self.bubble.va_enthalpy_density, phase=Phase.BROKEN)
+            gs_star=self.gs_star
         )
 
     @functools.cached_property
     def g_star(self) -> float:
+        if self.override_necessary or self.gs_star_manual_override:
+            return self.g_star_override
+        return self.g_star_computed
+
+    @functools.cached_property
+    def g_star_computed(self):
         return self.bubble.model.gp(w=self.bubble.va_enthalpy_density, phase=Phase.BROKEN)
 
+    @functools.cached_property
+    def gs_star(self) -> float:
+        if self.override_necessary or self.gs_star_manual_override:
+            return self.gs_star_override
+        return self.gs_star_computed
+
+    @functools.cached_property
+    def gs_star_computed(self) -> float:
+        return self.bubble.model.gs(w=self.bubble.va_enthalpy_density, phase=Phase.BROKEN)
+
     def noise(self) -> np.ndarray:
         return noise.omega_noise(self.f())
 
@@ -66,6 +132,12 @@ def signal_to_noise_ratio(self) -> float:
     def suppression_factor(self, method: sup.SuppressionMethod = sup.SuppressionMethod.DEFAULT) -> float:
         return sup.get_suppression_factor(vw=self.bubble.v_wall, alpha=self.bubble.alpha_n, method=method)
 
+    @functools.cached_property
+    def Tn(self) -> float:
+        if self.override_necessary or self.Tn_manual_override:
+            return self.Tn_override
+        return self.bubble.Tn
+
 
 def f(z: th.FloatOrArr, r_star: th.FloatOrArr, f_star0: th.FloatOrArr) -> th.FloatOrArr:
     r"""Convert the dimensionless wavenumber $z$ to frequency today by taking into account the redshift.
@@ -82,16 +154,16 @@ def f0(rs: th.FloatOrArr, T_n: th.FloatOrArr = const.T_default, g_star: float =
     return f_star0(T_n, g_star) / rs
 
 
-def f_star0(T_n: th.FloatOrArr, g_star: th.FloatOrArr = 100) -> th.FloatOrArr:
+def f_star0(Tn: th.FloatOrArr, g_star: th.FloatOrArr = 100) -> th.FloatOrArr:
     r"""
     Conversion factor between the frequencies at the time of the nucleation and frequencies today.
     $$f_{*,0} = 2.6 \cdot 10^{-6} \text{Hz} \left( \frac{T_n}{100 \text{GeV}} \right) \left( \frac{g_*}{100} \right)^{\frac{1}{6}}$$,
     :gowling_2021:`\ ` eq. 2.13
-    :param T_n: Nucleation temperature
+    :param Tn: Nucleation temperature
     :param g_star: Degrees of freedom at the time the GWs were produced. The default value is from the article.
     :return:
     """
-    return const.fs0_ref * (T_n / 100) * (g_star / 100)**(1/6)
+    return const.fs0_ref * (Tn / 100) * (g_star / 100)**(1 / 6)
 
 
 def F_gw0(
@@ -165,7 +237,7 @@ def omgw0_bag(
 
 
 def r_star(H_n: th.FloatOrArr, R_star: th.FloatOrArr) -> th.FloatOrArr:
-    """
+    r"""
     $$r_* = H_n R_*$$
     :gowling_2021:`\ ` eq. 2.2
     """
diff --git a/pttools/ssmtools/spectrum.py b/pttools/ssmtools/spectrum.py
index 354ee23..582b673 100644
--- a/pttools/ssmtools/spectrum.py
+++ b/pttools/ssmtools/spectrum.py
@@ -41,20 +41,19 @@ def __init__(
             y: np.ndarray = None,
             z_st_thresh: float = const.Z_ST_THRESH,
             nuc_type: NucType = DEFAULT_NUC_TYPE,
-            # method: ssm.Method = ssm.Method.E_CONSERVING,
-            # de_method: ssm.DE_Method = ssm.DE_Method.STANDARD,
-            # nxi: int = 5000,
             nt: int = 10000,
-            # nq: int = 1000,
             r_star: float = None,
             lifetime_multiplier: float = 1,
             compute: bool = True):
         r"""
-        :param bubble: Fluid simulation object
+        :param bubble: the Bubble object
         :param y: $z = kR*$ array
         :param z_st_thresh: for $z$ values above z_sh_tresh, use approximation rather than doing the sine transform integral.
         :param nuc_type: nucleation type
-        :param source_duration: $\frac{\Delta \eta_\text{v}}{\eta_*}$
+        :param nt: number of points in the t array
+        :param r_star: $r_*$
+        :param lifetime_multiplier: used for computing the source lifetime factor
+        :param compute: whether to compute the spectrum immediately
         """
         # Parameters
         self.bubble = bubble