Skip to content

Add two genie datat sets to the test data #26

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
kjmeagher merged 1 commit into from
Dec 19, 2023
Merged

Add two genie datat sets to the test data #26

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
135 changes: 76 additions & 59 deletions contrib/book_simweights_testdata.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -7,12 +7,12 @@
"""Script to generate the test data used by simweights testing."""

import os.path
import shutil
import sys
import tarfile
import tempfile
from pathlib import Path

from I3Tray import I3Tray
from icecube import dataclasses, hdfwriter, icetray, rootwriter, simclasses, tableio # noqa: F401


Expand All @@ -28,66 +28,74 @@ def fake_event_header(frame: dict) -> None:
fake_event_header.event_id = 0


filelist = [
"/data/sim/IceCube/2015/filtered/level2/CORSIKA-in-ice/12602/0000000-0000999/Level2_IC86.2015_corsika.012602.000000.i3.bz2",
"/data/sim/IceCube/2015/filtered/level2/CORSIKA-in-ice/20014/0000000-0000999/Level2_IC86.2015_corsika.020014.000000.i3.bz2",
"/data/sim/IceCube/2015/filtered/level2/CORSIKA-in-ice/20021/0000000-0000999/Level2_IC86.2015_corsika.020021.000000.i3.bz2",
"/data/sim/IceCube/2016/filtered/level2/CORSIKA-in-ice/20208/0000000-0000999/Level2_IC86.2016_corsika.020208.000001.i3.bz2",
"/data/sim/IceCube/2016/filtered/level2/CORSIKA-in-ice/20243/0000000-0000999/Level2_IC86.2016_corsika.020243.000001.i3.bz2",
"/data/sim/IceCube/2016/filtered/level2/CORSIKA-in-ice/20263/0000000-0000999/Level2_IC86.2016_corsika.020263.000000.i3.zst",
"/data/sim/IceCube/2016/filtered/level2/CORSIKA-in-ice/20777/0000000-0000999/Level2_IC86.2016_corsika.020777.000000.i3.zst",
"/data/sim/IceCube/2016/filtered/level2/CORSIKA-in-ice/20778/0000000-0000999/Level2_IC86.2016_corsika.020778.000000.i3.zst",
"/data/sim/IceCube/2016/filtered/level2/CORSIKA-in-ice/20780/0000000-0000999/Level2_IC86.2016_corsika.020780.000000.i3.zst",
"/data/sim/IceCube/2016/filtered/level2/CORSIKA-in-ice/21889/0000000-0000999/Level2_IC86.2016_corsika.021889.000000.i3.zst",
"/data/sim/IceCube/2011/filtered/level2/neutrino-generator/10634/00000-00999/Level2_IC86.2011_nugen_NuMu.010634.000000.i3.bz2",
"/data/sim/IceCube/2011/filtered/level2/neutrino-generator/10692/00000-00999/Level2_IC86.2011_nugen_NuE.010692.000000.i3.bz2",
"/data/sim/IceCube/2012/filtered/level2/neutrino-generator/11029/00000-00999/Level2_nugen_numu_IC86.2012.011029.000000.i3.bz2",
"/data/sim/IceCube/2012/filtered/level2/neutrino-generator/11065/00000-00999/Level2_IC86.2012_nugen_NuTau.011065.000001.i3.bz2",
"/data/sim/IceCube/2012/filtered/level2/neutrino-generator/11069/00000-00999/Level2_nugen_numu_IC86.2012.011069.000000.i3.bz2",
"/data/sim/IceCube/2012/filtered/level2/neutrino-generator/11070/00000-00999/Level2_nugen_numu_IC86.2012.011070.000000.i3.bz2",
"/data/sim/IceCube/2012/filtered/level2/neutrino-generator/11297/00000-00999/Level2_nugen_nutau_IC86.2012.011297.000000.i3.bz2",
"/data/sim/IceCube/2012/filtered/level2/neutrino-generator/11374/00000-00999/clsim-base-4.0.3.0.99_eff/Level2_IC86.2012_nugen_numu.011374.000050.clsim-base-4.0.3.0.99_eff.i3.bz2",
"/data/sim/IceCube/2012/filtered/level2/neutrino-generator/11477/00000-00999/clsim-base-4.0.3.0.99_eff/Level2_IC86.2012_nugen_nutau.011477.000000.clsim-base-4.0.3.0.99_eff.i3.bz2",
"/data/sim/IceCube/2012/filtered/level2/neutrino-generator/11836/00000-00999/clsim-base-4.0.3.0.99_eff/Level2_IC86.2012_nugen_nutau.011836.000000.clsim-base-4.0.3.0.99_eff.i3.bz2",
"/data/sim/IceCube/2012/filtered/level2/neutrino-generator/12646/0000000-0000999/clsim-base-4.0.5.0.99_eff/Level2_IC86.2012_nugen_nue.012646.000000.clsim-base-4.0.5.0.99_eff.i3.bz2",
"/data/sim/IceCube/2016/filtered/level2/neutrino-generator/20878/0000000-0000999/Level2_IC86.2016_NuMu.020878.000000.i3.zst",
"/data/sim/IceCube/2016/filtered/level2/neutrino-generator/20885/0000000-0000999/Level2_IC86.2016_NuE.020885.000000.i3.zst",
"/data/sim/IceCube/2016/filtered/level2/neutrino-generator/20895/0000000-0000999/Level2_IC86.2016_NuTau.020895.000000.i3.zst",
"/data/ana/CosmicRay/IceTop_level3/sim/IC86.2012/SIBYLL2.1/p/12360_v1s/Level3_IC86.2012_SIBYLL2.1_p_12360_E6.0_0.i3.bz2",
"/data/sim/IceCubeUpgrade/genie/step3/140021/upgrade_genie_step3_140021_000000.i3.zst",
]
filelist = {
"corsika": [
"/data/sim/IceCube/2015/filtered/level2/CORSIKA-in-ice/12602/0000000-0000999/Level2_IC86.2015_corsika.012602.000000.i3.bz2",
"/data/sim/IceCube/2015/filtered/level2/CORSIKA-in-ice/20014/0000000-0000999/Level2_IC86.2015_corsika.020014.000000.i3.bz2",
"/data/sim/IceCube/2015/filtered/level2/CORSIKA-in-ice/20021/0000000-0000999/Level2_IC86.2015_corsika.020021.000000.i3.bz2",
"/data/sim/IceCube/2016/filtered/level2/CORSIKA-in-ice/20208/0000000-0000999/Level2_IC86.2016_corsika.020208.000001.i3.bz2",
"/data/sim/IceCube/2016/filtered/level2/CORSIKA-in-ice/20243/0000000-0000999/Level2_IC86.2016_corsika.020243.000001.i3.bz2",
"/data/sim/IceCube/2016/filtered/level2/CORSIKA-in-ice/20263/0000000-0000999/Level2_IC86.2016_corsika.020263.000000.i3.zst",
"/data/sim/IceCube/2016/filtered/level2/CORSIKA-in-ice/20777/0000000-0000999/Level2_IC86.2016_corsika.020777.000000.i3.zst",
"/data/sim/IceCube/2016/filtered/level2/CORSIKA-in-ice/20778/0000000-0000999/Level2_IC86.2016_corsika.020778.000000.i3.zst",
"/data/sim/IceCube/2016/filtered/level2/CORSIKA-in-ice/20780/0000000-0000999/Level2_IC86.2016_corsika.020780.000000.i3.zst",
"/data/sim/IceCube/2016/filtered/level2/CORSIKA-in-ice/21889/0000000-0000999/Level2_IC86.2016_corsika.021889.000000.i3.zst",
],
"nugen": [
"/data/sim/IceCube/2011/filtered/level2/neutrino-generator/10634/00000-00999/Level2_IC86.2011_nugen_NuMu.010634.000000.i3.bz2",
"/data/sim/IceCube/2011/filtered/level2/neutrino-generator/10692/00000-00999/Level2_IC86.2011_nugen_NuE.010692.000000.i3.bz2",
"/data/sim/IceCube/2012/filtered/level2/neutrino-generator/11029/00000-00999/Level2_nugen_numu_IC86.2012.011029.000000.i3.bz2",
"/data/sim/IceCube/2012/filtered/level2/neutrino-generator/11065/00000-00999/Level2_IC86.2012_nugen_NuTau.011065.000001.i3.bz2",
"/data/sim/IceCube/2012/filtered/level2/neutrino-generator/11069/00000-00999/Level2_nugen_numu_IC86.2012.011069.000000.i3.bz2",
"/data/sim/IceCube/2012/filtered/level2/neutrino-generator/11070/00000-00999/Level2_nugen_numu_IC86.2012.011070.000000.i3.bz2",
"/data/sim/IceCube/2012/filtered/level2/neutrino-generator/11297/00000-00999/Level2_nugen_nutau_IC86.2012.011297.000000.i3.bz2",
"/data/sim/IceCube/2012/filtered/level2/neutrino-generator/11374/00000-00999/clsim-base-4.0.3.0.99_eff/Level2_IC86.2012_nugen_numu.011374.000050.clsim-base-4.0.3.0.99_eff.i3.bz2",
"/data/sim/IceCube/2012/filtered/level2/neutrino-generator/11477/00000-00999/clsim-base-4.0.3.0.99_eff/Level2_IC86.2012_nugen_nutau.011477.000000.clsim-base-4.0.3.0.99_eff.i3.bz2",
"/data/sim/IceCube/2012/filtered/level2/neutrino-generator/11836/00000-00999/clsim-base-4.0.3.0.99_eff/Level2_IC86.2012_nugen_nutau.011836.000000.clsim-base-4.0.3.0.99_eff.i3.bz2",
"/data/sim/IceCube/2012/filtered/level2/neutrino-generator/12646/0000000-0000999/clsim-base-4.0.5.0.99_eff/Level2_IC86.2012_nugen_nue.012646.000000.clsim-base-4.0.5.0.99_eff.i3.bz2",
"/data/sim/IceCube/2016/filtered/level2/neutrino-generator/20878/0000000-0000999/Level2_IC86.2016_NuMu.020878.000000.i3.zst",
"/data/sim/IceCube/2016/filtered/level2/neutrino-generator/20885/0000000-0000999/Level2_IC86.2016_NuE.020885.000000.i3.zst",
"/data/sim/IceCube/2016/filtered/level2/neutrino-generator/20895/0000000-0000999/Level2_IC86.2016_NuTau.020895.000000.i3.zst",
],
"icetop": [
"/data/ana/CosmicRay/IceTop_level3/sim/IC86.2012/SIBYLL2.1/p/12360_v1s/Level3_IC86.2012_SIBYLL2.1_p_12360_E6.0_0.i3.bz2",
],
"genie": [
"/data/sim/IceCubeUpgrade/genie/step3/141828/upgrade_genie_step3_141828_000000.i3.zst",
"/data/sim/IceCube/2023/generated/GENIE/22590/0000000-0000999/GENIE_NuMu_IceCubeUpgrade_v58.22590.000000.i3.zst",
],
}
keys = {
"corsika": [
"PolyplopiaPrimary",
"I3PrimaryInjectorInfo",
"I3CorsikaInfo",
"CorsikaWeightMap",
"I3CorsikaWeight",
],
"nugen": ["I3MCWeightDict"],
"genie": ["I3GenieInfo", "I3GenieResult", "I3MCWeightDict"],
"icetop": ["I3TopInjectorInfo", "MCPrimary"],
}
streams = {
"corsika": ["InIceSplit"],
"nugen": ["InIceSplit", "in_ice"],
"genie": ["NullSplit"],
"icetop": ["IceTopSplit"],
}


if "notemp" in sys.argv:
outdir = "/scratch/kmeagher/simweights/"
outdir = Path("/scratch/kmeagher/simweights/")
else:
tempdir = tempfile.TemporaryDirectory(prefix="simweights_testdata_")
outdir = Path(tempdir.name)

for filename in filelist:
for simtype, filename in ((i, x) for i in filelist for x in filelist[i]):
basename = Path(filename).name.replace(".i3.zst", "").replace(".i3.bz2", "").replace(".i3.gz", "")
assert basename != Path(filename).name

split = False
if "corsika" in basename:
keys = [
"PolyplopiaPrimary",
"I3PrimaryInjectorInfo",
"I3CorsikaInfo",
"CorsikaWeightMap",
"I3CorsikaWeight",
]
streams = ["InIceSplit"]
elif "genie" in basename:
keys = ["I3GenieInfo", "I3GenieResult", "I3MCWeightDict"]
streams = ["NullSplit"]
split = True
elif "nugen" in basename or "Nu" in basename:
keys = ["I3MCWeightDict"]
streams = ["InIceSplit", "in_ice"]
else:
keys = ["I3TopInjectorInfo", "MCPrimary"]
streams = ["IceTopSplit"]

split = simtype == "genie"
outfile = outdir / basename

if Path(outfile.name + ".hdf5").exists():
Expand All @@ -96,10 +104,10 @@ def fake_event_header(frame: dict) -> None:

print(f"Booking : {filename}")
print(f" outfile: {outfile}")
print(f" keys : {keys}")
print(f" streams: {streams}")
print(f" keys : {keys[simtype]}")
print(f" streams: {streams[simtype]}")

tray = I3Tray()
tray = icetray.I3Tray()
tray.Add("I3Reader", FileNameList=[filename])

if split:
Expand All @@ -115,20 +123,29 @@ def fake_event_header(frame: dict) -> None:
hdfwriter.I3HDFTableService(str(outfile) + ".hdf5"),
rootwriter.I3ROOTTableService(str(outfile) + ".root"),
],
SubEventStreams=streams,
keys=keys,
SubEventStreams=streams[simtype],
keys=keys[simtype],
)
tray.Add("Keep", keys=keys[simtype])
tray.Add(
"I3Writer",
Filename=str(outfile) + ".i3.zst",
Streams=[icetray.I3Frame.Simulation, icetray.I3Frame.DAQ],
DropOrphanStreams=[icetray.I3Frame.Physics],
)

tray.Execute()
del tray


tarfilename = "/data/user/kmeagher/simweights_testdata_test.tar.gz"
print(f"Writing tarfile {tarfilename}")

shutil.copy("upgrade_genie_step3_140021_000000.root", outdir)
shutil.copy("upgrade_genie_step3_140021_000000.hdf5", outdir)

with tarfile.open(tarfilename, "w:gz") as tar:
for f in os.listdir(outdir):
print(f"Adding {f} to tarball")
tar.add(outdir / f, arcname=f)

print("Done!")
print("Finished writing tarfile {tarfilename}")
12 changes: 9 additions & 3 deletions tests/test_genie_datasets.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -53,12 +53,12 @@ def cmp_dataset(self, fname):
with self.subTest(lib=str(fobj)):
w = GenieWeighter(fobj)

pdf0 = w.surface.spectra[14][0].dists[0]
pdf0 = next(iter(w.surface.spectra.values()))[0].dists[0]
np.testing.assert_allclose(1 / pdf0.v, global_probability_scale * solid_angle * injection_area, 1e-5)

np.testing.assert_allclose(w.get_weight_column("wght"), genie_weight)

power_law = w.surface.spectra[14][0].dists[2]
power_law = next(iter(w.surface.spectra.values()))[0].dists[2]
energy_term = 1 / power_law.pdf(w.get_weight_column("energy"))
np.testing.assert_allclose(energy_term, energy_factor)

Expand All @@ -70,9 +70,15 @@ def cmp_dataset(self, fname):
for fobj in fobjs:
fobj.close()

def test_NuE(self):
def test_140021(self):
self.cmp_dataset("upgrade_genie_step3_140021_000000")

def test_141828(self):
self.cmp_dataset("upgrade_genie_step3_141828_000000")

def test_22590(self):
self.cmp_dataset("GENIE_NuMu_IceCubeUpgrade_v58.22590.000000")


if __name__ == "__main__":
unittest.main()