Skip to content

Commit

Permalink
Formatting
Browse files Browse the repository at this point in the history
  • Loading branch information
JBorrow committed Sep 11, 2024
1 parent 725b600 commit aa3c9f5
Show file tree
Hide file tree
Showing 7 changed files with 76 additions and 54 deletions.
2 changes: 0 additions & 2 deletions swiftsimio/optional_packages.py
Original file line number Diff line number Diff line change
Expand Up @@ -92,5 +92,3 @@ def x(func):

# For additional CUDA API access
cuda = None


1 change: 0 additions & 1 deletion swiftsimio/visualisation/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -6,4 +6,3 @@
from .slice import slice_scatter as slice
from .slice import slice_gas, slice_gas_pixel_grid
from .smoothing_length import generate_smoothing_lengths

18 changes: 9 additions & 9 deletions swiftsimio/visualisation/power_spectrum.py
Original file line number Diff line number Diff line change
Expand Up @@ -199,7 +199,7 @@ def render_to_deposit(
"""

# Get the positions and masses
folding = 2.0**folding
folding = 2.0 ** folding
positions = data.coordinates
quantity = getattr(data, project)

Expand Down Expand Up @@ -244,10 +244,10 @@ def render_to_deposit(
units = 1.0 / (
data.metadata.boxsize[0] * data.metadata.boxsize[1] * data.metadata.boxsize[2]
)
units.convert_to_units(1.0 / data.metadata.boxsize.units**3)
units.convert_to_units(1.0 / data.metadata.boxsize.units ** 3)

units *= quantity.units
new_cosmo_factor = quantity.cosmo_factor / (coord_cosmo_factor**3)
new_cosmo_factor = quantity.cosmo_factor / (coord_cosmo_factor ** 3)

return cosmo_array(
deposition, comoving=comoving, cosmo_factor=new_cosmo_factor, units=units
Expand Down Expand Up @@ -399,7 +399,7 @@ def folded_depositions_to_power_spectrum(

if folding != final_folding:
cutoff_wavenumber = (
2.0**folding * np.min(depositions[folding].shape) / np.min(box_size)
2.0 ** folding * np.min(depositions[folding].shape) / np.min(box_size)
)

if cutoff_above_wavenumber_fraction is not None:
Expand All @@ -424,7 +424,7 @@ def folded_depositions_to_power_spectrum(
corrected_wavenumber_centers[prefer_bins] = folded_wavenumber_centers[
prefer_bins
].to(corrected_wavenumber_centers.units)
folding_tracker[prefer_bins] = 2.0**folding
folding_tracker[prefer_bins] = 2.0 ** folding

contributed_counts[prefer_bins] = folded_counts[prefer_bins]
elif transition == "average":
Expand Down Expand Up @@ -457,7 +457,7 @@ def folded_depositions_to_power_spectrum(

# For debugging, we calculate an effective fold number.
folding_tracker[use_bins] = (
(folding_tracker * existing_weight + (2.0**folding) * new_weight)
(folding_tracker * existing_weight + (2.0 ** folding) * new_weight)
/ transition_norm
)[use_bins]

Expand Down Expand Up @@ -538,7 +538,7 @@ def deposition_to_power_spectrum(
deposition.shape == cross_deposition.shape
), "Depositions must have the same shape"

folding = 2.0**folding
folding = 2.0 ** folding

box_size_folded = box_size[0] / folding
npix = deposition.shape[0]
Expand All @@ -560,7 +560,7 @@ def deposition_to_power_spectrum(
else:
conj_fft = fft.conj()

fourier_amplitudes = (fft * conj_fft).real * box_size_folded**3
fourier_amplitudes = (fft * conj_fft).real * box_size_folded ** 3

# Calculate k-value spacing (centered FFT)
dk = 2 * np.pi / (box_size_folded)
Expand Down Expand Up @@ -592,7 +592,7 @@ def deposition_to_power_spectrum(
divisor[zero_mask] = 1

# Correct for folding
binned_amplitudes *= folding**3
binned_amplitudes *= folding ** 3

# Correct units and names
wavenumbers = unyt.unyt_array(
Expand Down
50 changes: 35 additions & 15 deletions swiftsimio/visualisation/ray_trace.py
Original file line number Diff line number Diff line change
Expand Up @@ -14,11 +14,15 @@

from swiftsimio.objects import cosmo_array
from swiftsimio.reader import __SWIFTParticleDataset, SWIFTDataset
from swiftsimio.visualisation.projection_backends.kernels import kernel_gamma, kernel_double_precision as kernel
from swiftsimio.visualisation.projection_backends.kernels import (
kernel_gamma,
kernel_double_precision as kernel,
)

from swiftsimio.accelerated import jit, prange, NUM_THREADS
import unyt


@jit(nopython=True, fastmath=True)
def core_panels(
x: np.float64,
Expand Down Expand Up @@ -114,7 +118,9 @@ def core_panels(
and particle_cell_y >= 0
and particle_cell_y <= maximal_array_index
):
output[particle_cell_x, particle_cell_y, panel] += (m[i] * inverse_cell_area)
output[particle_cell_x, particle_cell_y, panel] += (
m[i] * inverse_cell_area
)
continue

normalisation = 0.0
Expand Down Expand Up @@ -157,12 +163,11 @@ def core_panels(

r = math.sqrt(distance_x_2 + distance_y_2)

output[cell_x, cell_y, panel] += (
kernel(r, kernel_width) * normalisation
)
output[cell_x, cell_y, panel] += kernel(r, kernel_width) * normalisation

return output


@jit(nopython=True, fastmath=True)
def core_panels_parallel(
x: np.float64,
Expand Down Expand Up @@ -250,9 +255,9 @@ def panel_pixel_grid(
raise AttributeError(
f'Comoving quantity "{project}" is not compatible with physical coordinates!'
)
m = m.value
m = m.value

# This provides a default 'slice it all' mask.
# This provides a default 'slice it all' mask.
if mask is None:
mask = np.s_[:]

Expand All @@ -275,15 +280,14 @@ def panel_pixel_grid(
z_min = unyt.unyt_quantity(0.0, units=box_z.units)
z_max = box_z


x_range = x_max - x_min
y_range = y_max - y_min

# Deal with non-cubic boxes:
# we always use the maximum of x_range and y_range to normalise the coordinates
# empty pixels in the resulting square image are trimmed afterwards
max_range = max(x_range, y_range)

try:
hsml = data.smoothing_lengths
except AttributeError:
Expand Down Expand Up @@ -322,6 +326,7 @@ def panel_pixel_grid(
max_z=z_max,
)


def panel_gas(
data: SWIFTDataset,
resolution: int,
Expand Down Expand Up @@ -373,18 +378,24 @@ def panel_gas(
)



# --- Functions that actually perform the 'ray tracing'.


def transfer_function(value, width, center):
"""
A simple gaussian transfer function centered around a specific value.
"""
return 1 / (width * np.sqrt(2.0 * np.pi)) * np.exp(-0.5 * ((value - center) / width) ** 2)
return (
1
/ (width * np.sqrt(2.0 * np.pi))
* np.exp(-0.5 * ((value - center) / width) ** 2)
)


@jit(fastmath=True, nopython=True)
def integrate_ray_numba_specific(input: np.array, red: float, green: float, blue: float, center: float, width: float):
def integrate_ray_numba_specific(
input: np.array, red: float, green: float, blue: float, center: float, width: float
):
"""
Given a ray, integrate the transfer function along it
"""
Expand All @@ -393,10 +404,16 @@ def integrate_ray_numba_specific(input: np.array, red: float, green: float, blue
color = np.array([red, green, blue], dtype=np.float32)

for i in input:
value += color * 1 / (width * np.sqrt(2.0 * np.pi)) * np.exp(-0.5 * ((i - center) / width) ** 2)
value += (
color
* 1
/ (width * np.sqrt(2.0 * np.pi))
* np.exp(-0.5 * ((i - center) / width) ** 2)
)

return value / len(input)


@jit(fastmath=True, nopython=True)
def integrate_ray_numba_nocolor(input: np.array, center: float, width: float):
"""
Expand All @@ -407,7 +424,11 @@ def integrate_ray_numba_nocolor(input: np.array, center: float, width: float):

for i in input:
value *= 0.99
value += 1 / (width * np.sqrt(2.0 * np.pi)) * np.exp(-0.5 * ((i - center) / width) ** 2)
value += (
1
/ (width * np.sqrt(2.0 * np.pi))
* np.exp(-0.5 * ((i - center) / width) ** 2)
)

return np.float32(value / len(input))

Expand Down Expand Up @@ -445,7 +466,6 @@ def integrate_ray_numba_nocolor(input: np.array, center: float, width: float):
# return output



# # %%
# import matplotlib.pyplot as plt
# import swiftascmaps
Expand Down
8 changes: 2 additions & 6 deletions swiftsimio/visualisation/smoothing_length/__init__.py
Original file line number Diff line number Diff line change
@@ -1,10 +1,6 @@
from .sph import get_hsml as sph_get_hsml
from .nearest_neighbours import (
get_hsml as nearest_neighbours_get_hsml,
)
from .generate import (
generate_smoothing_lengths,
)
from .nearest_neighbours import get_hsml as nearest_neighbours_get_hsml
from .generate import generate_smoothing_lengths

backends_get_hsml = {
"sph": sph_get_hsml,
Expand Down
18 changes: 10 additions & 8 deletions swiftsimio/visualisation/volume_render.py
Original file line number Diff line number Diff line change
Expand Up @@ -845,14 +845,14 @@ def render_gas(
* data.metadata.boxsize[1]
* data.metadata.boxsize[2]
)
units.convert_to_units(1.0 / data.metadata.boxsize.units**3)
units.convert_to_units(1.0 / data.metadata.boxsize.units ** 3)

comoving = data.gas.coordinates.comoving
coord_cosmo_factor = data.gas.coordinates.cosmo_factor
if project is not None:
units *= getattr(data.gas, project).units
project_cosmo_factor = getattr(data.gas, project).cosmo_factor
new_cosmo_factor = project_cosmo_factor / coord_cosmo_factor**3
new_cosmo_factor = project_cosmo_factor / coord_cosmo_factor ** 3
else:
new_cosmo_factor = coord_cosmo_factor ** (-3)

Expand Down Expand Up @@ -946,7 +946,8 @@ def visualise_render(
]

images = [
array([color[0] * x, color[1] * x, color[2] * x]).T for color, x in zip(colors, images)
array([color[0] * x, color[1] * x, color[2] * x]).T
for color, x in zip(colors, images)
]

if return_type == "all":
Expand All @@ -960,9 +961,7 @@ def visualise_render(


def visualise_render_options(
centers: list[float],
widths: Union[list[float], float],
cmap: str = "viridis",
centers: list[float], widths: Union[list[float], float], cmap: str = "viridis"
) -> tuple["plt.Figure", "plt.Axes"]:
"""
Creates a figure of your rendering options. The y-axis is the output value
Expand Down Expand Up @@ -998,10 +997,13 @@ def visualise_render_options(

for center, width, color in zip(centers, widths, colors):
xs = linspace(center - 5.0 * width, center + 5.0 * width, 100)
ys = [exp(-0.5 * ((center - x) / width)**2) / (width * sqrt(2.0 * pi)) for x in xs]
ys = [
exp(-0.5 * ((center - x) / width) ** 2) / (width * sqrt(2.0 * pi))
for x in xs
]

ax.axvline(center, color=color, linestyle="--")

ax.plot(xs, ys, color=color)

return fig, ax
return fig, ax
33 changes: 20 additions & 13 deletions tests/test_visualisation.py
Original file line number Diff line number Diff line change
Expand Up @@ -416,25 +416,24 @@ def test_panel_rendering(filename):
res = 1024

# Test the panel rendering
panel = panel_gas(
data,
resolution=res,
panels=N_depth,
project="masses",
)
panel = panel_gas(data, resolution=res, panels=N_depth, project="masses")

assert panel.shape[-1] == N_depth

import matplotlib.pyplot as plt
from matplotlib.colors import LogNorm

plt.imsave("panels_added.png", plt.get_cmap()(LogNorm(vmin=10**6, vmax=10**6.5)(np.sum(panel, axis=-1))))

projected = project_gas(
data, res, "masses", backend="renormalised"
plt.imsave(
"panels_added.png",
plt.get_cmap()(LogNorm(vmin=10 ** 6, vmax=10 ** 6.5)(np.sum(panel, axis=-1))),
)

plt.imsave("projected.png", plt.get_cmap()(LogNorm(vmin=10**6, vmax=10**6.5)(projected)),)
projected = project_gas(data, res, "masses", backend="renormalised")

plt.imsave(
"projected.png",
plt.get_cmap()(LogNorm(vmin=10 ** 6, vmax=10 ** 6.5)(projected)),
)

fullstack = np.zeros((res, res))

Expand All @@ -443,11 +442,19 @@ def test_panel_rendering(filename):

offset = 32

plt.imsave("stacked.png", plt.get_cmap()(LogNorm()(fullstack[offset:-offset, offset:-offset])))
plt.imsave(
"stacked.png",
plt.get_cmap()(LogNorm()(fullstack[offset:-offset, offset:-offset])),
)

assert np.isclose(panel.sum(axis=-1)[offset:-offset, offset:-offset], projected[offset:-offset, offset:-offset], rtol=0.1).all()
assert np.isclose(
panel.sum(axis=-1)[offset:-offset, offset:-offset],
projected[offset:-offset, offset:-offset],
rtol=0.1,
).all()
return


def test_periodic_boundary_wrapping():
"""
Test that periodic boundary wrapping works.
Expand Down

0 comments on commit aa3c9f5

Please sign in to comment.