From f999bbf40fcc2850ce9afc6429fcd4d4b95109e9 Mon Sep 17 00:00:00 2001 From: Thomas Robitaille Date: Fri, 24 Jul 2026 06:50:50 +0000 Subject: [PATCH 1/3] Read each spot spectrum from its own HDF5 group instead of the parent sphere group --- src/sources/source_type.f90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/sources/source_type.f90 b/src/sources/source_type.f90 index 2b0a710d..ed41ea47 100644 --- a/src/sources/source_type.f90 +++ b/src/sources/source_type.f90 @@ -177,7 +177,7 @@ subroutine source_read(group, s) s%spot(i)%cost = cos(spot_size * deg2rad) g_spot = mp_open_group(group, spot_names(i)) - call set_spectrum(group, s%spot(i)%freq_type, s%spot(i)%spectrum, s%spot(i)%temperature) + call set_spectrum(g_spot, s%spot(i)%freq_type, s%spot(i)%spectrum, s%spot(i)%temperature) call mp_close_group(g_spot) if(s%freq_type == 3) call error("source_read", "Spot cannot have LTE spectrum") From 2acb119c2ff2dd07ad36ed33c1be18e05feefcc9 Mon Sep 17 00:00:00 2001 From: Thomas Robitaille Date: Fri, 24 Jul 2026 06:52:49 +0000 Subject: [PATCH 2/3] Only sample from the sphere spectrum in source_emit when the photon was not emitted from a spot --- src/sources/source_type.f90 | 54 ++++++++++++++++++++----------------- 1 file changed, 30 insertions(+), 24 deletions(-) diff --git a/src/sources/source_type.f90 b/src/sources/source_type.f90 index ed41ea47..abbdf839 100644 --- a/src/sources/source_type.f90 +++ b/src/sources/source_type.f90 @@ -456,19 +456,22 @@ subroutine source_emit(src,p,nu) p%emiss_var_frac = jnu_var_frac(p%icell%ic, p%dust_id) call dust_sample_emit_probability(d(p%dust_id),p%emiss_var_id,p%emiss_var_frac,nu,p%energy) end select - end if - select case(src%freq_type) - case(1) - p%energy = interpolate_pdf(src%spectrum, nu, bounds_error=.false., fill_value=0._dp) - case(2) - p%energy = normalized_B_nu(nu, src%temperature) - case(3) - p%dust_id = select_dust_specific_energy_rho(p%icell) - p%emiss_var_id = jnu_var_id(p%icell%ic, p%dust_id) - p%emiss_var_frac = jnu_var_frac(p%icell%ic, p%dust_id) - call dust_sample_emit_probability(d(p%dust_id),p%emiss_var_id,p%emiss_var_frac,nu,p%energy) - end select + else + + select case(src%freq_type) + case(1) + p%energy = interpolate_pdf(src%spectrum, nu, bounds_error=.false., fill_value=0._dp) + case(2) + p%energy = normalized_B_nu(nu, src%temperature) + case(3) + p%dust_id = select_dust_specific_energy_rho(p%icell) + p%emiss_var_id = jnu_var_id(p%icell%ic, p%dust_id) + p%emiss_var_frac = jnu_var_frac(p%icell%ic, p%dust_id) + call dust_sample_emit_probability(d(p%dust_id),p%emiss_var_id,p%emiss_var_frac,nu,p%energy) + end select + + end if else @@ -486,19 +489,22 @@ subroutine source_emit(src,p,nu) p%emiss_var_frac = jnu_var_frac(p%icell%ic, p%dust_id) call dust_sample_j_nu(d(p%dust_id),p%emiss_var_id,p%emiss_var_frac,p%nu) end select - end if - select case(src%freq_type) - case(1) - p%nu = sample_pdf(src%spectrum) - case(2) - call random_planck_frequency(p%nu, src%temperature) - case(3) - p%dust_id = select_dust_specific_energy_rho(p%icell) - p%emiss_var_id = jnu_var_id(p%icell%ic, p%dust_id) - p%emiss_var_frac = jnu_var_frac(p%icell%ic, p%dust_id) - call dust_sample_j_nu(d(p%dust_id),p%emiss_var_id,p%emiss_var_frac,p%nu) - end select + else + + select case(src%freq_type) + case(1) + p%nu = sample_pdf(src%spectrum) + case(2) + call random_planck_frequency(p%nu, src%temperature) + case(3) + p%dust_id = select_dust_specific_energy_rho(p%icell) + p%emiss_var_id = jnu_var_id(p%icell%ic, p%dust_id) + p%emiss_var_frac = jnu_var_frac(p%icell%ic, p%dust_id) + call dust_sample_j_nu(d(p%dust_id),p%emiss_var_id,p%emiss_var_frac,p%nu) + end select + + end if end if From 41879c3fac9bea8f41f92b61670275d474136ec0 Mon Sep 17 00:00:00 2001 From: Thomas Robitaille Date: Mon, 27 Jul 2026 12:33:20 +0000 Subject: [PATCH 3/3] Add regression test that a spot emits with its own spectrum rather than the sphere spectrum --- hyperion/model/tests/test_spot_source.py | 70 ++++++++++++++++++++++++ 1 file changed, 70 insertions(+) create mode 100644 hyperion/model/tests/test_spot_source.py diff --git a/hyperion/model/tests/test_spot_source.py b/hyperion/model/tests/test_spot_source.py new file mode 100644 index 00000000..bf0c5a04 --- /dev/null +++ b/hyperion/model/tests/test_spot_source.py @@ -0,0 +1,70 @@ +from __future__ import print_function, division + +import os +import tempfile + +import numpy as np +import pytest + +from .. import Model +from ...util.functions import random_id + + +@pytest.mark.requires_hyperion_binaries +def test_spot_uses_its_own_spectrum(): + # Regression test: photons emitted from a spot must use the spot's own + # spectrum, not the parent sphere's. We give the sphere and the spot + # disjoint emission bands (no dust, so photons escape directly), so any + # flux in the spot's band can only come from the spot having used its own + # spectrum. Before the fix, source_emit overwrote the spot-sampled + # frequency with the sphere spectrum, so the spot band got zero flux. + + m = Model() + m.set_cartesian_grid([-1e12, 1e12], [-1e12, 1e12], [-1e12, 1e12]) + + nu = np.logspace(np.log10(3e12), np.log10(1e15), 300) + fnu_sphere = np.where((nu > 1e13) & (nu < 2e13), 1., 0.) # ~15-30 micron + fnu_spot = np.where((nu > 3e14) & (nu < 6e14), 1., 0.) # ~0.5-1 micron (disjoint) + + s = m.add_spherical_source() + s.radius = 1e11 + s.position = (0., 0., 0.) + s.luminosity = 1. + s.spectrum = (nu, fnu_sphere) + + spot = s.add_spot() + spot.longitude = 0. + spot.latitude = 0. + spot.radius = s.radius + spot.luminosity = 1. + spot.spectrum = (nu, fnu_spot) + + sed = m.add_peeled_images(sed=True, image=False) + sed.set_viewing_angles([45.], [45.]) + sed.set_wavelength_range(60, 0.1, 100.) + sed.set_aperture_radii(1, 1e12, 1e12) + + m.set_n_initial_iterations(0) + m.set_n_photons(imaging=20000) + + tmpdir = tempfile.mkdtemp() + m.write(os.path.join(tmpdir, random_id())) + mo = m.run() + + wav, nufnu = mo.get_sed() + wav = np.array(wav) + nufnu = np.squeeze(np.array(nufnu)) + + sphere_band = (wav > 10.) & (wav < 40.) # ~15-30 micron + spot_band = (wav > 0.4) & (wav < 1.2) # ~0.5-1 micron + + sphere_flux = np.nansum(nufnu[sphere_band]) + spot_flux = np.nansum(nufnu[spot_band]) + + # the model actually ran and the sphere emitted in its band + assert sphere_flux > 0. + # the spot emitted in its OWN band (exactly zero before the fix) + assert spot_flux > 0. + assert spot_flux > 0.2 * sphere_flux + # and essentially nothing leaked between the two disjoint bands + assert spot_flux + sphere_flux > 0.9 * np.nansum(nufnu)