mirror of
https://github.com/slsdetectorgroup/aare.git
synced 2026-04-20 15:44:37 +02:00
lunin_l
8f8173feb6
To improve codebase quality and reduce human error, this PR introduces the pre-commit framework. This ensures that all code adheres to project standards before it is even committed, maintaining a consistent style and catching common mistakes early. Key Changes: - Code Formatting: Automated C++ formatting using clang-format (based on the project's .clang-format file). - Syntax Validation: Basic checks for file integrity and syntax. - Spell Check: Automated scanning for typos in source code and comments. - CMake Formatting: Standardization of CMakeLists.txt and .cmake configuration files. - GitHub Workflow: Added a CI action that validates every Pull Request against the pre-commit configuration to ensure compliance. The configuration includes a [ci] block to handle automated fixes within the PR. Currently, this is disabled. If we want the CI to automatically commit formatting fixes back to the PR branch, this can be toggled to true in .pre-commit-config.yaml. ```yaml ci: autofix_commit_msg: [pre-commit] auto fixes from pre-commit hooks autofix_prs: false autoupdate_schedule: monthly ``` The last large commit with the fit functions, for example, was not formatted according to the clang-format rules. This PR would allow to avoid similar mistakes in the future. Python fomat with `ruff` for tests and sanitiser for `.ipynb` notebooks can be added as well.
923 KiB
923 KiB
In [1]:
import numpy as np from scipy.stats import multivariate_normal import matplotlib.pyplot as plt from aare import ClusterFinder, calculate_eta2, Interpolator import boost_histogram as bh from aare._aare import Cluster2x2d, ClusterVector_Cluster2x2d import pickle import os from test_Interpolation import create_photon_hit_with_gaussian_distribution, create_2x2cluster_from_frame, create_3x3cluster_from_frame, calculate_eta_distribution, photon_hit_in_euclidean_space
In [2]:
### plotting functions def plot_eta_distribution(eta_distribution): plt.imshow( eta_distribution[:, :, 0].view().transpose(), origin='lower', extent=[eta_distribution.axes[0].edges[0], eta_distribution.axes[0].edges[-1], eta_distribution.axes[1].edges[0], eta_distribution.axes[1].edges[-1]], aspect='auto') plt.colorbar() plt.xlabel('eta_x') plt.ylabel('eta_y') def plot_uniform_eta_distribution(uniform_etax, uniform_etay, eta_extent): fig, axes = plt.subplots(1, 2, figsize=(10, 5)) # (rows, cols) im1 = axes[0].imshow( uniform_etax[:, :, 0].view().transpose(), origin='lower', extent=[eta_extent[0], eta_extent[1], eta_extent[0], eta_extent[1]], aspect='auto' ) axes[0].set_title("uniform distribution of etax") fig.colorbar(im1, ax=axes[0], fraction=0.046, pad=0.04) axes[0].set_xlabel("eta_x") axes[0].set_ylabel("eta_y") im2 = axes[1].imshow( uniform_etay[:, :, 0].view().transpose(), origin='lower', extent=[eta_extent[0], eta_extent[1], eta_extent[0], eta_extent[1]], aspect='auto' ) axes[1].set_title("uniform distribution of etay") fig.colorbar(im2, ax=axes[1], fraction=0.046, pad=0.04) axes[1].set_xlabel("eta_x") axes[1].set_ylabel("eta_y") plt.tight_layout() plt.show() def plot_photon_hit_distribution(eta_distribution, marg_CDF_etax, cond_CDF_etay): normalized_eta_distribution = eta_distribution.values()[:,:, 0] / eta_distribution.values()[: , :, 0].sum() #normnalize n_samples = 1000 rng = np.random.default_rng(42) flat_indices = rng.choice( normalized_eta_distribution.size, size=n_samples, p=normalized_eta_distribution.ravel()) x_indices, y_indices = np.unravel_index(flat_indices, normalized_eta_distribution.shape) photon_position_x = marg_CDF_etax[x_indices, y_indices] photon_position_y = cond_CDF_etay[x_indices, y_indices] fig, axes = plt.subplots(1, 2, figsize=(10, 5)) # (rows, cols) axes[0].scatter(eta_distribution.axes[0].edges[:-1][x_indices], eta_distribution.axes[1].edges[:-1][y_indices], s=1, alpha=0.5) axes[0].set_title("Eta distribution (eta_x, eta_y)") #plt.gca().set_aspect('equal') axes[0].set_xlabel("eta_x") axes[0].set_ylabel("eta_y") #axes[0].set_xlim(0,1) #axes[0].set_ylim(0,1) axes[1].scatter(photon_position_x, photon_position_y, s=1, alpha=0.5) axes[1].set_title("uniform Photon positions") #plt.gca().set_aspect('equal') axes[1].set_xlabel("x") axes[1].set_ylabel("y") #axes[1].set_xlim(0,1) #axes[1].set_ylim(0,1)
In [3]:
pixel_width = 1e-4 values = np.arange(0.5*pixel_width, 0.1, pixel_width) num_pixels = values.size print(f"num pixels: ({values.size} x {values.size})") X, Y = np.meshgrid(values, values) data_points = np.stack([X.ravel(), Y.ravel()], axis=1)
num pixels: (1000 x 1000)
In [4]:
# mean is in bottom right quadrant variance = 10*pixel_width covariance_matrix = np.array([[variance, 0], [0, variance]]) mean = 650*pixel_width mean = np.array([mean, mean]) base_frame = create_photon_hit_with_gaussian_distribution(mean, variance, data_points) pixels_per_superpixel = int(num_pixels*0.5) plt.imshow(base_frame) plt.colorbar() plt.axvline(x=pixels_per_superpixel, color='black', linestyle='--', linewidth=0.5) plt.axhline(y=pixels_per_superpixel, color='black', linestyle='--', linewidth=0.5) plt.xlabel('x_0') plt.ylabel('y_0')
Out[4]:
Text(0, 0.5, 'y_0')
In [5]:
### create eta distribution num_frames = 5000 random_number_generator = np.random.default_rng(42) eta_distribution = calculate_eta_distribution(num_frames, pixels_per_superpixel, random_number_generator, bh.axis.Regular(100, -0.1, 0.6), bh.axis.Regular(100, -0.1, 0.6)) test_data_path = os.getenv("AARE_TEST_DATA") + "/eta_distributions" filename = test_data_path + "/eta_distribution_2x2cluster_gaussian.pkl" with open(filename, "wb") as f: pickle.dump(eta_distribution, f)
In [7]:
filename = os.getenv("AARE_TEST_DATA") + "/eta_distributions/eta_distribution_2x2cluster_gaussian.pkl" with open(filename, "rb") as f: eta_distribution = pickle.load(f)
In [8]:
plot_eta_distribution(eta_distribution)
In [9]:
interpolator = Interpolator(eta_distribution, eta_distribution.axes[0].edges, eta_distribution.axes[1].edges, eta_distribution.axes[2].edges[:-1]) #interpolator = Interpolator(eta_distribution.axes[0].edges, eta_distribution.axes[1].edges, eta_distribution.axes[2].edges[:-1]) #interpolator.rosenblatttransform(eta_distribution.values())
In [18]:
marg_CDF_etax = interpolator.get_ietax() cond_CDF_etay = interpolator.get_ietay() plot_photon_hit_distribution(eta_distribution, marg_CDF_etax, cond_CDF_etay)
In [10]:
marg_CDF_etax = interpolator.get_ietax() cond_CDF_etay = interpolator.get_ietay() #plot_uniform_eta_distribution(uniform_etax, uniform_etay, [eta_distribution.axes[0].edges[0], eta_distribution.axes[0].edges[-1]])
In [11]:
cluster = create_2x2cluster_from_frame(base_frame, pixels_per_superpixel) print("cluster.x:", cluster.x) print("cluster.y:", cluster.y) clustervec = ClusterVector_Cluster2x2d() clustervec.push_back(cluster) eta = calculate_eta2(cluster) print("eta: ", eta) bin_size = (eta_distribution.axes[0].edges[-1] - eta_distribution.axes[0].edges[0])/eta_distribution.axes[0].edges.shape[0] bin_index_x = int((eta[0] - eta_distribution.axes[0].edges[0])/bin_size) bin_index_y = int((eta[1] - eta_distribution.axes[1].edges[0])/bin_size) print("distance x:", marg_CDF_etax[bin_index_x, bin_index_y, 0]) print("distance y:", cond_CDF_etay[bin_index_x, bin_index_y, 0]) photon_hit = interpolator.interpolate(clustervec) print(photon_hit)
cluster.x: 1 cluster.y: 1 eta: (0.3519491547025806, 0.3519491547025806, 71553511.0) distance x: 0.7083333333333334 distance y: 0.7564102564102564 [(0.78166653, 0.81284991, 71553511.)]
In [12]:
# scale to cluster pixel width #cluster_center = num_pixels*0.75*pixel_width cluster_center = (pixels_per_superpixel*pixel_width*(cluster.x + 0.5), pixels_per_superpixel*pixel_width*(cluster.y + 0.5)) #scaled_photon_hit = photon_hit_in_euclidean_space(cluster_center, pixels_per_superpixel, photon_hit) print(photon_hit[0]) scaled_photon_hit = ((photon_hit[0][0])*pixels_per_superpixel*pixel_width, (photon_hit[0][1])*pixels_per_superpixel*pixel_width) print(f"previous center: ({cluster_center[0]}, {cluster_center[1]})") print(f"interpolated center: ({photon_hit[0][0]},{photon_hit[0][1]})") print(f"scaled interpolated center: ({scaled_photon_hit[0]},{scaled_photon_hit[1]})") print(f"actual center: ({mean},{mean})")
(0.7816665310999981, 0.8128499138725078, 71553511.0) previous center: (0.07500000000000001, 0.07500000000000001) interpolated center: (0.7816665310999981,0.8128499138725078) scaled interpolated center: (0.03908332655499991,0.040642495693625394) actual center: ([0.065 0.065],[0.065 0.065])
In [13]:
plt.imshow(base_frame) plt.colorbar() plt.axvline(x=pixels_per_superpixel, color='black', linestyle='--', linewidth=0.5) plt.axhline(y=pixels_per_superpixel, color='black', linestyle='--', linewidth=0.5) plt.scatter(mean/pixel_width, mean/pixel_width, color='red', s=40, marker='x', label='Actual photon hit') plt.scatter(scaled_photon_hit[0]/pixel_width, scaled_photon_hit[1]/pixel_width, color = 'blue', s=40, marker='x', label='interpolated photon hit') plt.xlabel('x_0') plt.ylabel('y_0')
Out[13]:
Text(0, 0.5, 'y_0')
3x3 Cluster¶
In [14]:
from aare._aare import Cluster3x3d, ClusterVector_Cluster3x3d
In [15]:
# mean is in center quadrant variance = 10*pixel_width covariance_matrix = np.array([[variance, 0], [0, variance]]) mean_x = (1 + 0.8)*(num_pixels/3)*pixel_width mean_y = (1 + 0.2)*(num_pixels/3)*pixel_width mean = np.array([mean_x, mean_y]) base_frame = create_photon_hit_with_gaussian_distribution(mean, variance, data_points) pixels_per_superpixel = int(num_pixels/3) plt.imshow(base_frame) plt.colorbar() plt.axvline(x=pixels_per_superpixel, color='black', linestyle='--', linewidth=0.5) plt.axhline(y=pixels_per_superpixel, color='black', linestyle='--', linewidth=0.5) plt.axvline(x=2*pixels_per_superpixel, color='black', linestyle='--', linewidth=0.5) plt.axhline(y=2*pixels_per_superpixel, color='black', linestyle='--', linewidth=0.5) plt.xlabel('x_0') plt.ylabel('y_0')
Out[15]:
Text(0, 0.5, 'y_0')
In [ ]:
### calculate eta distribution num_frames = 5000 random_number_generator = np.random.default_rng(42) eta_distribution = calculate_eta_distribution(num_frames, pixels_per_superpixel, random_number_generator, bh.axis.Regular(100, -0.1, 1.0), bh.axis.Regular(100, -0.1, 1.0), False) filename = os.getenv("AARE_TEST_DATA") + "/eta_distributions/eta_distribution_3x3cluster_gaussian.pkl" with open(filename, "wb") as f: pickle.dump(eta_distribution, f)
In [ ]:
filename = os.getenv("AARE_TEST_DATA") + "/eta_distributions/eta_distribution_3x3cluster_gaussian.pkl" with open(filename, "rb") as f: eta_distribution = pickle.load(f)
In [7]:
plot_eta_distribution(eta_distribution)
In [ ]:
In [14]:
#interpolator = Interpolator(eta_distribution, eta_distribution.axes[0].edges, eta_distribution.axes[1].edges, eta_distribution.axes[2].edges[:-1]) interpolator = Interpolator(eta_distribution.axes[0].edges, eta_distribution.axes[1].edges, eta_distribution.axes[2].edges[:-1]) interpolator.rosenblatttransform(eta_distribution.values()) marg_CDF_etax = interpolator.get_ietax() cond_CDF_etay = interpolator.get_ietay()
In [15]:
plot_photon_hit_distribution(eta_distribution, marg_CDF_etax, cond_CDF_etay)
In [ ]:
plot_uniform_eta_distribution(uniform_etax, uniform_etay, [eta_distribution.axes[0].edges[0], eta_distribution.axes[0].edges[-1]])
In [16]:
cluster = create_3x3cluster_from_frame(base_frame, pixels_per_superpixel) eta = calculate_eta2(cluster) print(eta) print(eta_distribution.axes[0].edges) clustervec = ClusterVector_Cluster3x3d() #print(marg_CDF_etax) print(marg_CDF_etax[93, 0,0]) clustervec.push_back(cluster) photon_hit = interpolator.interpolate(clustervec)
(0.5499154381017892, 0.44933336606613455, 48591988.0) [-0.1 -0.093 -0.086 -0.079 -0.072 -0.065 -0.058 -0.051 -0.044 -0.037 -0.03 -0.023 -0.016 -0.009 -0.002 0.005 0.012 0.019 0.026 0.033 0.04 0.047 0.054 0.061 0.068 0.075 0.082 0.089 0.096 0.103 0.11 0.117 0.124 0.131 0.138 0.145 0.152 0.159 0.166 0.173 0.18 0.187 0.194 0.201 0.208 0.215 0.222 0.229 0.236 0.243 0.25 0.257 0.264 0.271 0.278 0.285 0.292 0.299 0.306 0.313 0.32 0.327 0.334 0.341 0.348 0.355 0.362 0.369 0.376 0.383 0.39 0.397 0.404 0.411 0.418 0.425 0.432 0.439 0.446 0.453 0.46 0.467 0.474 0.481 0.488 0.495 0.502 0.509 0.516 0.523 0.53 0.537 0.544 0.551 0.558 0.565 0.572 0.579 0.586 0.593 0.6 ] 0.9999999999999998 dX: 1 dY: 0 u: 1 v: 0
In [17]:
## scale to cluster pixel width #cluster_center = (1+0.5)*pixels_per_superpixel*pixel_width cluster_center = (pixels_per_superpixel*pixel_width*(cluster.x + 0.5), pixels_per_superpixel*pixel_width*(cluster.y + 0.5)) print(photon_hit[0]) scaled_photon_hit = (photon_hit[0][0]*pixels_per_superpixel*pixel_width, photon_hit[0][1]*pixels_per_superpixel*pixel_width) print(f"previous center: ({cluster_center[0]}, {cluster_center[1]})") print(f"interpolated center: ({photon_hit[0][0]},{photon_hit[0][1]})") print(f"scaled interpolated center: ({scaled_photon_hit[0]},{scaled_photon_hit[1]})") print(f"actual center: ({mean},{mean})")
(1.5000000000000002, 1.5, 48591988.0) previous center: (0.04995000000000001, 0.04995000000000001) interpolated center: (1.5000000000000002,1.5) scaled interpolated center: (0.04995000000000001,0.04995) actual center: ([0.06 0.04],[0.06 0.04])
In [18]:
plt.imshow(base_frame) plt.colorbar() plt.axvline(x=pixels_per_superpixel, color='black', linestyle='--', linewidth=0.5) plt.axhline(y=pixels_per_superpixel, color='black', linestyle='--', linewidth=0.5) plt.axvline(x=2*pixels_per_superpixel, color='black', linestyle='--', linewidth=0.5) plt.axhline(y=2*pixels_per_superpixel, color='black', linestyle='--', linewidth=0.5) plt.scatter(mean_x/pixel_width, mean_y/pixel_width, color='red', s=40, marker='x', label='Actual photon hit') plt.scatter(scaled_photon_hit[0]/pixel_width, scaled_photon_hit[1]/pixel_width, color = 'blue', s=40, marker='x', label='interpolated photon hit') plt.xlabel('x_0') plt.ylabel('y_0')
Out[18]:
Text(0, 0.5, 'y_0')
Interpolation with Eta3¶
In [6]:
from aare._aare import Cluster3x3d, ClusterVector_Cluster3x3d, calculate_eta3, calculate_cross_eta3
In [4]:
# mean is in center quadrant variance = 10*pixel_width covariance_matrix = np.array([[variance, 0], [0, variance]]) mean = (1 + 0.3)*(num_pixels/3)*pixel_width mean = np.array([mean, mean]) base_frame = create_photon_hit_with_gaussian_distribution(mean, variance, data_points) pixels_per_superpixel = int(num_pixels/3) plt.imshow(base_frame) plt.colorbar() plt.axvline(x=pixels_per_superpixel, color='black', linestyle='--', linewidth=0.5) plt.axhline(y=pixels_per_superpixel, color='black', linestyle='--', linewidth=0.5) plt.axvline(x=2*pixels_per_superpixel, color='black', linestyle='--', linewidth=0.5) plt.axhline(y=2*pixels_per_superpixel, color='black', linestyle='--', linewidth=0.5) plt.xlabel('x_0') plt.ylabel('y_0')
Out[4]:
Text(0, 0.5, 'y_0')
In [7]:
num_frames = 2000 random_number_generator = np.random.default_rng(42) eta_distribution = calculate_eta_distribution(num_frames, pixels_per_superpixel, random_number_generator, calculate_eta3, False)
In [8]:
plot_eta_distribution(eta_distribution)
In [9]:
num_frames = 2000 random_number_generator = np.random.default_rng(42) eta_distribution = calculate_eta_distribution(num_frames, pixels_per_superpixel, random_number_generator, calculate_cross_eta3, False)
In [10]:
plot_eta_distribution(eta_distribution)
In [ ]: