Skip to content
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

45 numpy vs np #46

Closed
wants to merge 2 commits into from
Closed

45 numpy vs np #46

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
143ec2b
replaced all occurences of numpy with np
shallmann Jun 18, 2022
190c2f9
remove unwanted file
shallmann Jun 18, 2022
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
6 changes: 3 additions & 3 deletions conda-recipes/nusigma/run_test.sh
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -4,12 +4,12 @@
check=$(
python <<-EOF | tail -n1
import nusigma
import numpy
import numpy as np

numpy.testing.assert_allclose(nusigma.nucross(1e5,1,"p","CC",1), 1.80246030548e-34)
np.testing.assert_allclose(nusigma.nucross(1e5,1,"p","CC",1), 1.80246030548e-34)
print("checkychecky")
EOF
)
if [ "$check" != "checkychecky" ]; then
exit 1
fi
fi
15 changes: 8 additions & 7 deletions notebooks/tutorials/001-Framework_Introduction.ipynb
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -41,9 +41,9 @@
}
],
"source": [
"import numpy\n",
"import numpy as np\n",
"from toise import factory\n",
"factory.set_kwargs(psi_bins={k: [0, numpy.pi] for k in ('tracks', 'cascades', 'radio')})"
"factory.set_kwargs(psi_bins={k: [0, np.pi] for k in ('tracks', 'cascades', 'radio')})"
]
},
{
Expand Down Expand Up @@ -216,7 +216,7 @@
"def make_components(aeffs):\n",
" aeff, muon_aeff = aeffs\n",
" \n",
" energy_threshold = numpy.inf\n",
" energy_threshold = np.inf\n",
" atmo = diffuse.AtmosphericNu.conventional(aeff, 1, hard_veto_threshold=energy_threshold)\n",
" atmo.prior = lambda v, **kwargs: -(v-1)**2/(2*0.1**2)\n",
" prompt = diffuse.AtmosphericNu.prompt(aeff, 1, hard_veto_threshold=energy_threshold)\n",
Expand Down Expand Up @@ -356,9 +356,9 @@
],
"source": [
"# also plot a profile likelihood\n",
"x = numpy.linspace(0., 3e-2, 11)\n",
"x = np.linspace(0., 3e-2, 11)\n",
"prof = llh.profile1d('astro', x)\n",
"ts = -2*(prof['LLH'] - numpy.nanmax(prof['LLH']))\n",
"ts = -2*(prof['LLH'] - np.nanmax(prof['LLH']))\n",
"\n",
"plt.plot(x, ts)\n",
"plt.axvline(limit, ls='--', color='grey')\n",
Expand Down Expand Up @@ -667,7 +667,8 @@
"hash": "58c0895f5e0218fe42ea0d65c34afa611258fa581a151421405042a8d733e1b8"
},
"kernelspec": {
"display_name": "Python 3.9.7 64-bit ('gen2-analysis': conda)",
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
Expand All @@ -680,7 +681,7 @@
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.9.9"
"version": "3.8.5"
}
},
"nbformat": 4,
Expand Down
20 changes: 10 additions & 10 deletions resources/scripts/2021_midscale/count_events.ipynb
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,7 @@
"metadata": {},
"outputs": [],
"source": [
"import numpy\n",
"import numpy as np\n",
"from toise import factory, diffuse, plotting, effective_areas\n",
"%matplotlib inline\n",
"import matplotlib.pyplot as plt"
Expand Down Expand Up @@ -59,12 +59,12 @@
],
"source": [
"# change the factory defaults\n",
"factory.set_kwargs(psi_bins={k: [0, numpy.pi] for k in ('tracks', 'cascades', 'radio')})\n",
"factory.set_kwargs(psi_bins={k: [0, np.pi] for k in ('tracks', 'cascades', 'radio')})\n",
"factory.set_kwargs(cos_theta=[-1,1])\n",
"\n",
"def make_components(aeffs):\n",
" nu, mu = aeffs\n",
" energy_threshold = numpy.inf\n",
" energy_threshold = np.inf\n",
" astro = diffuse.DiffuseAstro(nu, livetime=1)\n",
" return dict(astro=astro)\n",
"\n",
Expand All @@ -77,12 +77,12 @@
"bundle = factory.component_bundle({'IceCube': 0, 'Gen2-InIce': 0, 'Gen2-HCR': 1}, make_components)\n",
"components = bundle.get_components()\n",
"expectations = components['astro'].expectations(gamma=-2.)\n",
"total = numpy.zeros([len(energy_centers)]) # empty array for total numbers\n",
"total = np.zeros([len(energy_centers)]) # empty array for total numbers\n",
"for stream, counts in expectations.iteritems():\n",
" if 'shadow' in stream: # count both the shadowed and unshadowed tracks\n",
" total += counts\n",
"\n",
"total_num = numpy.sum(total[mask])\n",
"total_num = np.sum(total[mask])\n",
"print(\"total num above 100 TeV {}\".format(total_num))"
]
},
Expand All @@ -108,21 +108,21 @@
],
"metadata": {
"kernelspec": {
"display_name": "Python 2",
"display_name": "Python 3",
"language": "python",
"name": "python2"
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 2
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython2",
"version": "2.7.15"
"pygments_lexer": "ipython3",
"version": "3.8.5"
}
},
"nbformat": 4,
Expand Down
54 changes: 27 additions & 27 deletions resources/scripts/calculate_fom.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -49,7 +49,7 @@ def get_label(opts):
opts.spacing = 125.0

import sys, os
import numpy
import numpy as np


from toise import (
Expand Down Expand Up @@ -85,9 +85,9 @@ def survey_distance(phi, L0=1e45):
"""
phi = phi * (intflux(1e5) - intflux(1e-1)) # integrate from 100 GeV to 100 PeV
phi *= 1e3 / grb.units.GeV # TeV / cm^2 s -> erg / cm^2 s
dl = numpy.sqrt(L0 / (4 * numpy.pi * phi))
dl = np.sqrt(L0 / (4 * np.pi * phi))
dl *= grb.units.cm / 1e3 # cm -> Gpc
return numpy.where(numpy.isfinite(phi), dl, 0)
return np.where(np.isfinite(phi), dl, 0)


def survey_volume(sindec, phi, L0=1e45):
Expand All @@ -98,7 +98,7 @@ def survey_volume(sindec, phi, L0=1e45):
:returns: the volume in Gpc^3 in which standard candles of luminosity L0 would be detectable
"""
dl = survey_distance(phi, L0)
return ((sindec.max() - sindec.min()) * 2 * numpy.pi / 3.0) * ((dl ** 3).mean())
return ((sindec.max() - sindec.min()) * 2 * np.pi / 3.0) * ((dl ** 3).mean())


def print_result(value, **kwargs):
Expand Down Expand Up @@ -173,7 +173,7 @@ def make_components(aeff, zi):

dp = []

sindec = numpy.linspace(-1, 1, 21)
sindec = np.linspace(-1, 1, 21)
for zi in range(20):
bundle = factory.component_bundle(
{"IceCube": 0.0, "Gen2": opts.livetime}, partial(make_components, zi=zi)
Expand All @@ -197,12 +197,12 @@ def make_components(aeff, zi):
# print dp[-1]
# s = ps.expectations(**fixed)['tracks'].sum(axis=0)
# b = bkg.expectations['tracks'].sum(axis=0)
# if (~numpy.isnan(s/b)).any():
# if (~np.isnan(s/b)).any():
# pass
# # print s.cumsum()/s.sum()
# # print s/b
# # bkgs.append(bkg.expectations['tracks'][:,:20].sum())
dp = numpy.array(dp)[::-1]
dp = np.array(dp)[::-1]
print(dp)

volume = survey_volume(sindec, dp)
Expand All @@ -212,7 +212,7 @@ def make_components(aeff, zi):
elif opts.figure_of_merit == "ps_time_evolution":

livetime = opts.livetimes[0]
cos_theta = numpy.linspace(-1, 1, 21)
cos_theta = np.linspace(-1, 1, 21)
aeff = factory.create_aeff(opts, cos_theta=cos_theta)
energy_threshold = effective_areas.StepFunction(opts.veto_threshold, 90)
atmo = diffuse.AtmosphericNu.conventional(
Expand Down Expand Up @@ -241,8 +241,8 @@ def scale_livetime(component, livetime):
else:
return component

livetimes = numpy.linspace(*opts.livetimes)
dps = numpy.zeros(livetimes.size)
livetimes = np.linspace(*opts.livetimes)
dps = np.zeros(livetimes.size)
for i in tqdm(list(range(livetimes.size))):
lt = livetimes[i]
dps[i] = 1e-12 * pointsource.discovery_potential(
Expand All @@ -252,7 +252,7 @@ def scale_livetime(component, livetime):
)

if opts.outfile is not None:
numpy.savez(
np.savez(
opts.outfile,
livetime=livetimes,
discovery_potential=dps,
Expand All @@ -269,7 +269,7 @@ def scale_livetime(component, livetime):
if opts.outfile is None:
parser.error("You must supply an output file name")

aeff = factory.create_aeff(opts, cos_theta=numpy.linspace(-1, 1, 20))
aeff = factory.create_aeff(opts, cos_theta=np.linspace(-1, 1, 20))
energy_threshold = effective_areas.StepFunction(opts.veto_threshold, 90)
atmo = diffuse.AtmosphericNu.conventional(
aeff, opts.livetime, hard_veto_threshold=energy_threshold
Expand All @@ -283,7 +283,7 @@ def scale_livetime(component, livetime):

values = dict()

sindec = center(numpy.linspace(-1, 1, 20)[::-1])
sindec = center(np.linspace(-1, 1, 20)[::-1])
for zi in range(0, 19, 1):
ps = pointsource.SteadyPointSource(aeff, opts.livetime, zenith_bin=zi)
atmo_bkg = atmo.point_source_background(zenith_index=zi)
Expand Down Expand Up @@ -313,7 +313,7 @@ def scale_livetime(component, livetime):
pickle.dump(values, f, 2)

elif opts.figure_of_merit == "grb":
aeff = factory.create_aeff(opts, cos_theta=numpy.linspace(-1, 1, 21))
aeff = factory.create_aeff(opts, cos_theta=np.linspace(-1, 1, 21))
energy_threshold = effective_areas.StepFunction(opts.veto_threshold, 90)
atmo = diffuse.AtmosphericNu.conventional(
aeff, opts.livetime, hard_veto_threshold=energy_threshold
Expand All @@ -325,9 +325,9 @@ def scale_livetime(component, livetime):
astro.seed = 2
gamma = multillh.NuisanceParam(-2.3, 0.5, min=-2.7, max=-1.7)

z = 2 * numpy.ones(opts.livetime * 170 * 2)
t90 = numpy.ones(z.size) * 45.1
Eiso = 10 ** (53.5) * numpy.ones(z.size)
z = 2 * np.ones(opts.livetime * 170 * 2)
t90 = np.ones(z.size) * 45.1
Eiso = 10 ** (53.5) * np.ones(z.size)

pop = grb.GRBPopulation(aeff, z, Eiso)
atmo_bkg = atmo.point_source_background(
Expand Down Expand Up @@ -369,20 +369,20 @@ def components(aeff):
gzk = diffuse.AhlersGZK(aeff, 1.0)
return dict(atmo=atmo, prompt=prompt, astro=astro, gzk=gzk)

bundle = aeff_bundle(components, cos_theta=numpy.linspace(-1, 1, 21))
bundle = aeff_bundle(components, cos_theta=np.linspace(-1, 1, 21))
components = bundle.get_components()
components["gamma"] = multillh.NuisanceParam(-2.3, 0.5, min=-2.7, max=-1.7)
gzk = components.pop("gzk")
aeff = list(bundle.aeffs.values())[0]

pev = numpy.where(aeff.bin_edges[2][1:] > 5e7)[0][0]
pev = np.where(aeff.bin_edges[2][1:] > 5e7)[0][0]

def pev_events(observables):
return sum((v.sum(axis=0)[pev:].sum() for v in observables.values()))

ns = pev_events(gzk.expectations())
nb = pev_events(components["astro"].expectations(gamma=-2.3))
baseline = 5 * numpy.sqrt(nb) / ns
baseline = 5 * np.sqrt(nb) / ns

scale = pointsource.discovery_potential(
gzk, components, baseline=baseline, tolerance=1e-4, gamma=-2.3
Expand All @@ -401,7 +401,7 @@ def pev_events(observables):
if opts.outfile is None:
parser.error("You must supply an output file name")

aeff = factory.create_aeff(opts, cos_theta=numpy.linspace(-1, 1, 21))
aeff = factory.create_aeff(opts, cos_theta=np.linspace(-1, 1, 21))
energy_threshold = effective_areas.StepFunction(opts.veto_threshold, 90)
atmo = diffuse.AtmosphericNu.conventional(
aeff, opts.livetime, hard_veto_threshold=energy_threshold
Expand Down Expand Up @@ -433,9 +433,9 @@ def pev_events(observables):
)
)

numpy.savetxt(
np.savetxt(
opts.outfile,
numpy.vstack((energies, dps)).T,
np.vstack((energies, dps)).T,
header="# energy\tdiscovery flux [GeV cm-2 sr^-1 s^-1]",
)

Expand All @@ -462,7 +462,7 @@ def components(aeff):
)
astro_lo.seed = 2.0
astro = diffuse.DiffuseAstro(
aeff.restrict_energy_range(opts.energy_threshold, numpy.inf), 1
aeff.restrict_energy_range(opts.energy_threshold, np.inf), 1
)
astro.seed = 2.0
return dict(atmo=atmo, prompt=prompt, astro_lo=astro_lo, astro=astro)
Expand All @@ -471,7 +471,7 @@ def components(aeff):
astro.seed = 2.0
return dict(atmo=atmo, prompt=prompt, astro=astro)

bundle = aeff_bundle(components, cos_theta=numpy.linspace(-1, 1, 20))
bundle = aeff_bundle(components, cos_theta=np.linspace(-1, 1, 20))
components = bundle.get_components()
components["gamma"] = multillh.NuisanceParam(-2.3)
if two_component:
Expand Down Expand Up @@ -507,7 +507,7 @@ def ts_diff(gamma):
if plotit:
import pylab

g = numpy.linspace(g0 - 0.5, g0 + 0.5, 21)
g = np.linspace(g0 - 0.5, g0 + 0.5, 21)
pylab.plot(g, [ts_diff(g_) for g_ in g])
color = pylab.gca().lines[-1].get_color()

Expand Down Expand Up @@ -568,7 +568,7 @@ def ts(flux_norm, **fixed):
)
ns = exes["galactic"]["tracks"].sum()

print_result(numpy.sqrt(fit_ts), nb=nb, ns=ns)
print_result(np.sqrt(fit_ts), nb=nb, ns=ns)

else:
parser.error("Unknown figure of merit '%s'" % (opts.figure_of_merit))
15 changes: 8 additions & 7 deletions resources/scripts/data prep/fit_psf.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -5,7 +5,8 @@
from icecube.photospline.utils import pad_knots
from icecube.photospline import splinefitstable
from collections import defaultdict
import numpy, os
import numpy as np
import os


def create_psf(
Expand All @@ -31,9 +32,9 @@ def create_psf(
"alpha",
),
bins=(
numpy.linspace(3, 8, 51),
numpy.linspace(-1, 1, 21),
numpy.linspace(0, 25, int(25e2)),
np.linspace(3, 8, 51),
np.linspace(-1, 1, 21),
np.linspace(0, 25, int(25e2)),
),
)
del q
Expand Down Expand Up @@ -66,13 +67,13 @@ def fit_psf(h, smooth=1e-6):
power = [1.0, 1.0, 2.0]
nknots = [25, 25, 25]
knots = [
pad_knots(numpy.linspace(e[0] ** (1.0 / p), e[-1] ** (1.0 / p), n) ** p, o)
pad_knots(np.linspace(e[0] ** (1.0 / p), e[-1] ** (1.0 / p), n) ** p, o)
for n, o, p, e in zip(nknots, order, power, h.binedges)
]

z = h.bincontent.cumsum(axis=2)
w = 1.0 / h.squaredweights.cumsum(axis=2)
w[~numpy.isfinite(w)] = 0.0
w[~np.isfinite(w)] = 0.0

del h

Expand All @@ -94,7 +95,7 @@ def plot_psf(h, spline, cos_theta=0):
dashi.visual()

ax = plt.gca()
for logE in numpy.arange(4, 8):
for logE in np.arange(4, 8):

ei = h._h_binedges[0].searchsorted(logE) - 1
zi = h._h_binedges[1].searchsorted(cos_theta) - 1
Expand Down
Loading