mirror/halo_comparison
1
0
Fork 0
mirror of https://github.com/Findus23/halo_comparison.git synced 2024-09-19 16:03:50 +02:00
Code Releases Activity

Merge branch 'main' of github.com:Findus23/halo_comparison

Browse source
This commit is contained in:
glatterf42 2022-08-19 11:46:01 +02:00
commit d83da6d742
2 changed files with 201 additions and 203 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

374
PlotDaubechies.py
View file

@ -8,7 +8,6 @@ from typing import List
import matplotlib.pyplot as plt
import numpy as np
# two-fold upsampling -- https://cnx.org/contents/xsppCgXj@8.18:H_wA16rf@16/Upsampling
from matplotlib.axes import Axes
from matplotlib.figure import Figure
@ -101,38 +100,38 @@ g_Haar = np.array([0.7071067811865476, -0.7071067811865476])
xhaar, phihaar, psihaar = cascade_algorithm(h_Haar, g_Haar, maxit)
# DB2 -- http://wavelets.pybytes.com/wavelet/db2/
h_DB2 = np.array([0.4829629131, 0.8365163037, 0.2241438680, -0.1294095226])
g_DB2 = np.array([-0.1294095226, -0.2241438680, 0.8365163037, -0.4829629131])
# DB4 -- http://wavelets.pybytes.com/wavelet/db2/
h_DB4 = np.array([0.4829629131, 0.8365163037, 0.2241438680, -0.1294095226])
g_DB4 = np.array([-0.1294095226, -0.2241438680, 0.8365163037, -0.4829629131])
xdb2, phidb2, psidb2 = cascade_algorithm(h_DB2, g_DB2, maxit)
xdb4, phidb4, psidb4 = cascade_algorithm(h_DB4, g_DB4, maxit)
# DB3 -- http://wavelets.pybytes.com/wavelet/db3/
h_DB3 = np.array(
[
0.3326705529509569,
0.8068915093133388,
0.4598775021193313,
-0.13501102001039084,
-0.08544127388224149,
0.035226291882100656,
]
)
g_DB3 = np.array(
[
0.035226291882100656,
0.08544127388224149,
-0.13501102001039084,
-0.4598775021193313,
0.8068915093133388,
-0.3326705529509569,
]
)
# h_DB3 = np.array(
# [
# 0.3326705529509569,
# 0.8068915093133388,
# 0.4598775021193313,
# -0.13501102001039084,
# -0.08544127388224149,
# 0.035226291882100656,
# ]
# )
# g_DB3 = np.array(
# [
# 0.035226291882100656,
# 0.08544127388224149,
# -0.13501102001039084,
# -0.4598775021193313,
# 0.8068915093133388,
# -0.3326705529509569,
# ]
# )
#
# xdb3, phidb3, psidb3 = cascade_algorithm(h_DB3, g_DB3, maxit)
xdb3, phidb3, psidb3 = cascade_algorithm(h_DB3, g_DB3, maxit)
# DB4 -- http://wavelets.pybytes.com/wavelet/db4/
h_DB4 = np.array(
# DB8 -- http://wavelets.pybytes.com/wavelet/db4/
h_DB8 = np.array(
[
0.23037781330885523,
0.7148465705525415,
@ -144,7 +143,7 @@ h_DB4 = np.array(
-0.010597401784997278,
]
)
g_DB4 = np.array(
g_DB8 = np.array(
[
-0.010597401784997278,
-0.032883011666982945,
@ -157,138 +156,144 @@ g_DB4 = np.array(
]
)
xdb4, phidb4, psidb4 = cascade_algorithm(h_DB4, g_DB4, maxit)
# DB8 -- http://wavelets.pybytes.com/wavelet/db8/
h_DB8 = np.array(
[
0.05441584224308161,
0.3128715909144659,
0.6756307362980128,
0.5853546836548691,
-0.015829105256023893,
-0.2840155429624281,
0.00047248457399797254,
0.128747426620186,
-0.01736930100202211,
-0.04408825393106472,
0.013981027917015516,
0.008746094047015655,
-0.00487035299301066,
-0.0003917403729959771,
0.0006754494059985568,
-0.00011747678400228192,
]
)
g_DB8 = np.array(
[
-0.00011747678400228192,
-0.0006754494059985568,
-0.0003917403729959771,
0.00487035299301066,
0.008746094047015655,
-0.013981027917015516,
-0.04408825393106472,
0.01736930100202211,
0.128747426620186,
-0.00047248457399797254,
-0.2840155429624281,
0.015829105256023893,
0.5853546836548691,
-0.6756307362980128,
0.3128715909144659,
-0.05441584224308161,
]
)
xdb8, phidb8, psidb8 = cascade_algorithm(h_DB8, g_DB8, maxit)
# DB16 --
h_DB16 = np.array(
[
0.0031892209253436892,
0.03490771432362905,
0.1650642834886438,
0.43031272284545874,
0.6373563320829833,
0.44029025688580486,
-0.08975108940236352,
-0.3270633105274758,
-0.02791820813292813,
0.21119069394696974,
0.027340263752899923,
-0.13238830556335474,
-0.006239722752156254,
0.07592423604445779,
-0.007588974368642594,
-0.036888397691556774,
0.010297659641009963,
0.013993768859843242,
-0.006990014563390751,
-0.0036442796214883506,
0.00312802338120381,
0.00040789698084934395,
-0.0009410217493585433,
0.00011424152003843815,
0.00017478724522506327,
-6.103596621404321e-05,
-1.394566898819319e-05,
1.133660866126152e-05,
-1.0435713423102517e-06,
-7.363656785441815e-07,
2.3087840868545578e-07,
-2.1093396300980412e-08,
]
)
g_DB16 = np.array(
[
-2.1093396300980412e-08,
-2.3087840868545578e-07,
-7.363656785441815e-07,
1.0435713423102517e-06,
1.133660866126152e-05,
1.394566898819319e-05,
-6.103596621404321e-05,
-0.00017478724522506327,
0.00011424152003843815,
0.0009410217493585433,
0.00040789698084934395,
-0.00312802338120381,
-0.0036442796214883506,
0.006990014563390751,
0.013993768859843242,
-0.010297659641009963,
-0.036888397691556774,
0.007588974368642594,
0.07592423604445779,
0.006239722752156254,
-0.13238830556335474,
-0.027340263752899923,
0.21119069394696974,
0.02791820813292813,
-0.3270633105274758,
0.08975108940236352,
0.44029025688580486,
-0.6373563320829833,
0.43031272284545874,
-0.1650642834886438,
0.03490771432362905,
-0.0031892209253436892,
]
)
# # DB8 -- http://wavelets.pybytes.com/wavelet/db8/
# h_DB8 = np.array(
# [
# 0.05441584224308161,
# 0.3128715909144659,
# 0.6756307362980128,
# 0.5853546836548691,
# -0.015829105256023893,
# -0.2840155429624281,
# 0.00047248457399797254,
# 0.128747426620186,
# -0.01736930100202211,
# -0.04408825393106472,
# 0.013981027917015516,
# 0.008746094047015655,
# -0.00487035299301066,
# -0.0003917403729959771,
# 0.0006754494059985568,
# -0.00011747678400228192,
# ]
# )
# g_DB8 = np.array(
# [
# -0.00011747678400228192,
# -0.0006754494059985568,
# -0.0003917403729959771,
# 0.00487035299301066,
# 0.008746094047015655,
# -0.013981027917015516,
# -0.04408825393106472,
# 0.01736930100202211,
# 0.128747426620186,
# -0.00047248457399797254,
# -0.2840155429624281,
# 0.015829105256023893,
# 0.5853546836548691,
# -0.6756307362980128,
# 0.3128715909144659,
# -0.05441584224308161,
# ]
# )
#
# xdb8, phidb8, psidb8 = cascade_algorithm(h_DB8, g_DB8, maxit)
xdb16, phidb16, psidb16 = cascade_algorithm(h_DB16, g_DB16, maxit)
# # DB16 --
# h_DB16 = np.array(
# [
# 0.0031892209253436892,
# 0.03490771432362905,
# 0.1650642834886438,
# 0.43031272284545874,
# 0.6373563320829833,
# 0.44029025688580486,
# -0.08975108940236352,
# -0.3270633105274758,
# -0.02791820813292813,
# 0.21119069394696974,
# 0.027340263752899923,
# -0.13238830556335474,
# -0.006239722752156254,
# 0.07592423604445779,
# -0.007588974368642594,
# -0.036888397691556774,
# 0.010297659641009963,
# 0.013993768859843242,
# -0.006990014563390751,
# -0.0036442796214883506,
# 0.00312802338120381,
# 0.00040789698084934395,
# -0.0009410217493585433,
# 0.00011424152003843815,
# 0.00017478724522506327,
# -6.103596621404321e-05,
# -1.394566898819319e-05,
# 1.133660866126152e-05,
# -1.0435713423102517e-06,
# -7.363656785441815e-07,
# 2.3087840868545578e-07,
# -2.1093396300980412e-08,
# ]
# )
# g_DB16 = np.array(
# [
# -2.1093396300980412e-08,
# -2.3087840868545578e-07,
# -7.363656785441815e-07,
# 1.0435713423102517e-06,
# 1.133660866126152e-05,
# 1.394566898819319e-05,
# -6.103596621404321e-05,
# -0.00017478724522506327,
# 0.00011424152003843815,
# 0.0009410217493585433,
# 0.00040789698084934395,
# -0.00312802338120381,
# -0.0036442796214883506,
# 0.006990014563390751,
# 0.013993768859843242,
# -0.010297659641009963,
# -0.036888397691556774,
# 0.007588974368642594,
# 0.07592423604445779,
# 0.006239722752156254,
# -0.13238830556335474,
# -0.027340263752899923,
# 0.21119069394696974,
# 0.02791820813292813,
# -0.3270633105274758,
# 0.08975108940236352,
# 0.44029025688580486,
# -0.6373563320829833,
# 0.43031272284545874,
# -0.1650642834886438,
# 0.03490771432362905,
# -0.0031892209253436892,
# ]
# )
#
# xdb16, phidb16, psidb16 = cascade_algorithm(h_DB16, g_DB16, maxit)
###################################
fig: Figure
fig, ax = plt.subplots(
4,
3,
2,
figsize=figsize_from_page_fraction(height_to_width=12 / 8),
figsize=figsize_from_page_fraction(),
# sharex="all", sharey="all"
)
labels = ["Haar", "DB2", "DB4", "DB8", "DB16"]
fig.subplots_adjust(hspace=0, wspace=0)
labels = ["Haar (DB2)", "DB4", "DB8", "DB16"]
prange = [[0., 1.], [0., 3.], [0., 7.]]
for p in prange:
p[0] -= 0.25
p[1] += 0.25
ax[0, 0].set_title("scaling functions $\\varphi$")
ax[0, 1].set_title("wavelets $\\psi$")
@ -296,28 +301,28 @@ ax[0, 1].set_title("wavelets $\\psi$")
ax[0, 0].plot(xhaar, phihaar, lw=1)
ax[0, 1].plot(xhaar, psihaar, lw=1)
ax[1, 0].plot(xdb2, phidb2, lw=1)
ax[1, 1].plot(xdb2, psidb2, lw=1)
ax[1, 0].plot(xdb4, phidb4, lw=1)
ax[1, 1].plot(xdb4, psidb4, lw=1)
ax[2, 0].plot(xdb4, phidb4, lw=1)
ax[2, 1].plot(xdb4, psidb4, lw=1)
ax[2, 0].plot(xdb8, phidb8, lw=1)
ax[2, 1].plot(xdb8, psidb8, lw=1)
ax[3, 0].plot(xdb8, phidb8, lw=1)
ax[3, 1].plot(xdb8, psidb8, lw=1)
# ax[4, 0].plot(xdb16, phidb16, lw=1)
# ax[4, 1].plot(xdb16, psidb16, lw=1)
for a in ax.flatten():
a.set_xlabel("t")
for i, a in enumerate(ax.flatten()):
if i in [5, 6]:
a.set_xlabel('t')
a.plot(prange[i // 2], [0.0, 0.0], 'k:', lw=1)
a.set_xlim(prange[i // 2])
a.set_xticks(np.arange(prange[i // 2][1] - 0.25 + 1))
a.set_yticks([-1.0, 0.0, 1.0])
def inset_label(ax: Axes, text: str):
def inset_label(ax: Axes, text: str, top=False):
ax.text(
0.75,
0.05,
0.98,
0.90 if top else 0.05,
text,
horizontalalignment="left",
verticalalignment="bottom",
horizontalalignment="right",
verticalalignment="top" if top else "bottom",
transform=ax.transAxes,
)
@ -325,13 +330,13 @@ def inset_label(ax: Axes, text: str):
for a, label in zip(ax[:, 0], labels):
text = r"$\varphi_{\textrm{LABEL}}$".replace("LABEL", label)
inset_label(a, text)
a.set_ylim([-1.0, 1.5])
a.set_ylim([-1.0, 1.499])
for a, label in zip(ax[:, 1], labels):
text = r"$\psi_{\textrm{LABEL}}$".replace("LABEL", label)
inset_label(a, text)
inset_label(a, text, top=True)
# a.set_ylabel('$\\psi_{\\rm ' + i + '}[t]$')
a.set_ylim([-2, 2])
a.set_ylim([-1.499, 1.999])
fig.tight_layout()
fig.savefig(Path(f"~/tmp/wavelets.pdf").expanduser())
@ -358,42 +363,32 @@ def fourier_wavelet(h, g, n):
# ax.plot([0, np.pi], [1., 1.], 'k:')
kh, fphih, fpsih = fourier_wavelet(h_Haar, g_Haar, 256)
ax.plot(kh, np.abs(fphih) ** 2, label=r"$\hat\varphi_\textrm{Haar}$", c="C0")
ax.plot(kh, np.abs(fphih) ** 2, label=r"$\hat\varphi_\textrm{Haar (DB2)}$", c="C0")
ax.plot(
kh,
np.abs(fpsih) ** 2,
label=r"$\hat\psi_\textrm{Haar}$",
label=r"$\hat\psi_\textrm{Haar (DB2)}$",
c="C0",
linestyle="dashed",
)
kdb2, fphidb2, fpsidb2 = fourier_wavelet(h_DB2, g_DB2, 256)
ax.plot(kdb2, np.abs(fphidb2) ** 2, label=r"$\hat\varphi_\textrm{DB2}$", c="C1")
kdb4, fphidb4, fpsidb4 = fourier_wavelet(h_DB4, g_DB4, 256)
ax.plot(kdb4, np.abs(fphidb4) ** 2, label=r"$\hat\varphi_\textrm{DB4}$", c="C1")
ax.plot(
kdb2,
np.abs(fpsidb2) ** 2,
kdb4,
np.abs(fpsidb4) ** 2,
label=r"$\hat\psi_\textrm{DB2}$",
c="C1",
linestyle="dashed",
)
kdb4, fphidb4, fpsidb4 = fourier_wavelet(h_DB4, g_DB4, 256)
ax.plot(kdb4, np.abs(fphidb4) ** 2, label=r"$\hat\varphi_\textrm{DB4}$", c="C2")
ax.plot(
kdb4,
np.abs(fpsidb4) ** 2,
label=r"$\hat\psi_\textrm{DB4}$",
c="C2",
linestyle="dashed",
)
kdb8, fphidb8, fpsidb8 = fourier_wavelet(h_DB8, g_DB8, 256)
ax.plot(kdb8, np.abs(fphidb8) ** 2, label=r"$\hat\varphi_\textrm{DB8}$", c="C3")
ax.plot(kdb8, np.abs(fphidb8) ** 2, label=r"$\hat\varphi_\textrm{DB8}$", c="C2")
ax.plot(
kdb8,
np.abs(fpsidb8) ** 2,
label=r"$\hat\psi_\textrm{DB8}$",
c="C3",
label=r"$\hat\psi_\textrm{DB4}$",
c="C2",
linestyle="dashed",
)
@ -445,3 +440,4 @@ ax.set_xticklabels(
fig2.tight_layout()
fig2.savefig(Path(f"~/tmp/wavelet_power.pdf").expanduser())
plt.show()
f"C{i}"

30
spectra_plot.py
View file

@ -218,20 +218,22 @@ def create_plot(mode):
ax_end.set_xlim(right=k0 * resolutions[-1])
ax_end.set_ylim(0.8, 1.02)
if bottom_row:
if mode == "power":
ax_ics.legend(loc="lower left")
else:
lines: List[Line2D] = ax_ics.get_lines()
half_lines1 = []
half_lines2 = []
for line in lines:
if line.get_label().startswith("128"):
half_lines1.append(line)
else:
half_lines2.append(line)
lines: List[Line2D] = ax_ics.get_lines()
half_lines1 = []
half_lines2 = []
for line in lines:
lab = line.get_label()
if (
(mode == "cross" and lab.startswith("128"))
or
(mode == "power" and ("128" in lab or "256" in lab))
):
half_lines1.append(line)
else:
half_lines2.append(line)
ax_ics.legend(handles=half_lines1, loc="lower left")
ax_z1.legend(handles=half_lines2, loc="lower left")
ax_ics.legend(handles=half_lines1, loc="lower left")
ax_z1.legend(handles=half_lines2, loc="lower left")
if not bottom_row:
last_xtick: XTick = ax_ics.yaxis.get_major_ticks()[0]
@ -239,7 +241,7 @@ def create_plot(mode):
# fig.suptitle(f"Cross Spectra {time}") #Not needed for paper
# fig.tight_layout()
if mode=="cross":
if mode == "cross":
print(crossings)
crossings_df = pd.DataFrame(crossings, columns=combination_list, index=waveforms)
# print(crossings_df.to_markdown())