Skip to content

Commit

Permalink
possible fix for the failing tests and removal of unnecessary object
Browse files Browse the repository at this point in the history
  • Loading branch information
matthiasplum committed Oct 25, 2023
1 parent 0b5b481 commit a27dd7e
Show file tree
Hide file tree
Showing 5 changed files with 5 additions and 34 deletions.
8 changes: 1 addition & 7 deletions src/simweights/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -18,7 +18,6 @@
"CircleInjector",
"CorsikaWeighter",
"NaturalRateCylinder",
"IceTopNaturalRateCylinder",
"UniformSolidAngleCylinder",
"TIG1996",
"FixedFractionFlux",
Expand Down Expand Up @@ -64,11 +63,6 @@
from ._nugen_weighter import NuGenWeighter
from ._pdgcode import PDGCode
from ._powerlaw import PowerLaw
from ._spatial import (
CircleInjector,
IceTopNaturalRateCylinder,
NaturalRateCylinder,
UniformSolidAngleCylinder,
)
from ._spatial import CircleInjector, NaturalRateCylinder, UniformSolidAngleCylinder
from ._utils import corsika_to_pdg
from ._weighter import Weighter
4 changes: 2 additions & 2 deletions src/simweights/_icetop_weighter.py
Original file line number Diff line number Diff line change
Expand Up @@ -9,7 +9,7 @@

from ._generation_surface import GenerationSurface, generation_surface
from ._powerlaw import PowerLaw
from ._spatial import IceTopNaturalRateCylinder
from ._spatial import NaturalRateCylinder
from ._utils import get_column, get_table
from ._weighter import Weighter

Expand All @@ -25,7 +25,7 @@ def sframe_icetop_surface(table: Any) -> GenerationSurface:
get_column(table, "min_energy")[i],
get_column(table, "max_energy")[i],
)
spatial = IceTopNaturalRateCylinder(
spatial = NaturalRateCylinder(
0, # set cylinder height to 0 to get simple surface plane
get_column(table, "sampling_radius")[i],
np.cos(get_column(table, "max_zenith")[i]),
Expand Down
23 changes: 0 additions & 23 deletions src/simweights/_spatial.py
Original file line number Diff line number Diff line change
Expand Up @@ -126,29 +126,6 @@ def pdf(self, cos_zen: ArrayLike) -> NDArray[np.float64]:
)


class IceTopNaturalRateCylinder(CylinderBase):
r"""Angular distribution when zenith distribution matched the natural rate of an isotropic source.
For a given zenith angle the intensity of particles thrown was proportional to the cross-sectional area
perpendicular to the direction of the particle.
"""

def __init__(self, length: float, radius: float, cos_zen_min: float, cos_zen_max: float) -> None:
super().__init__(length, radius, cos_zen_min, cos_zen_max)
self.etendue = (
2 * np.pi * (self.cos_zen_max - self.cos_zen_min) * self._cap
) # This is needed to make the tests work
self._normalization = 1 / self.etendue

def pdf(self, cos_zen: ArrayLike) -> NDArray[np.float64]:
cosz = np.asfarray(cos_zen)
return np.piecewise(
cosz,
[(cosz >= self.cos_zen_min) & (cosz <= self.cos_zen_max)],
[self._normalization],
)


class CircleInjector:
"""Particles are generated on a circle perpendicular to the direction of the particle.
Expand Down
2 changes: 1 addition & 1 deletion tests/test_icetop_datasets.py
Original file line number Diff line number Diff line change
Expand Up @@ -31,7 +31,7 @@ def cmp_dataset(self, fname):
assert len(reffile["I3TopInjectorInfo"]) == 1
si = reffile["I3TopInjectorInfo"][0]
pri = reffile["MCPrimary"]
solid_angle = 2 * np.pi * (np.cos(si["min_zenith"]) - np.cos(si["max_zenith"]))
solid_angle = np.pi * (np.cos(si["min_zenith"]) ** 2 - np.cos(si["max_zenith"]) ** 2)
injection_area = np.pi * (si["sampling_radius"] * 1e2) ** 2
energy_integral = np.log(si["max_energy"] / si["min_energy"]) # assuming E^-1
energy_factor = energy_integral * pri["energy"]
Expand Down
2 changes: 1 addition & 1 deletion tests/test_icetop_weighter.py
Original file line number Diff line number Diff line change
Expand Up @@ -27,7 +27,7 @@ class TestIceTopWeighter(unittest.TestCase):
def test_icetop_corsika(self):
nevents = 10000
pdgid = 12
c1 = simweights.IceTopNaturalRateCylinder(0, 300, 0, 1)
c1 = simweights.NaturalRateCylinder(0, 300, 0, 1)
p1 = simweights.PowerLaw(0, 1e3, 1e4)

weight = np.zeros(nevents, dtype=particle_dtype)
Expand Down

0 comments on commit a27dd7e

Please sign in to comment.