Update src/seismic_hazard_forecasting.py

fix orientation issue attempt 1
This commit is contained in:
2026-06-22 19:06:19 +02:00
parent b7ee8d52ab
commit 0ef41d72f6
+8 -7
View File
@@ -699,18 +699,22 @@ verbose: {verbose}")
fxy = xy_kde[0] fxy = xy_kde[0]
logger.debug(f"Normalization check; sum of all f(x,y) values = {np.sum(fxy)}") logger.debug(f"Normalization check; sum of all f(x,y) values = {np.sum(fxy)}")
xx, yy = np.meshgrid(x_range, y_range, indexing='ij') # grid points xx, yy = np.meshgrid(x_range, y_range) # grid points
# set every grid point to be a receiver # set every grid point to be a receiver
grid_shape = xx.shape
x_rx = xx.flatten() x_rx = xx.flatten()
y_rx = yy.flatten() y_rx = yy.flatten()
num_points = x_rx.size
distances = np.zeros(shape=(num_points, grid_shape[0], grid_shape[1]))
# compute distance matrix for each receiver # compute distance matrix for each receiver
distances = np.zeros(shape=(nx * ny, nx, ny)) #distances = np.zeros(shape=(nx * ny, nx, ny))
rx_lat = np.zeros(nx * ny) rx_lat = np.zeros(nx * ny)
rx_lon = np.zeros(nx * ny) rx_lon = np.zeros(nx * ny)
for i in range(nx * ny): for i in range(num_points):
# Compute the squared distances directly using NumPy's vectorized operations # Compute the squared distances directly using NumPy's vectorized operations
squared_distances = (xx - x_rx[i]) ** 2 + (yy - y_rx[i]) ** 2 squared_distances = (xx - x_rx[i]) ** 2 + (yy - y_rx[i]) ** 2
distances[i] = np.sqrt(squared_distances) distances[i] = np.sqrt(squared_distances)
@@ -726,7 +730,7 @@ verbose: {verbose}")
if exclude_low_fxy: if exclude_low_fxy:
indices = list(np.where(fxy.flatten() > thresh_fxy)[0]) indices = list(np.where(fxy.flatten() > thresh_fxy)[0])
else: else:
indices = range(0, len(distances)) indices = np.arange(num_points)
if use_AOI: if use_AOI:
# Filter out receivers outside the AOI; Find indices where values are OUTSIDE the AOI # Filter out receivers outside the AOI; Find indices where values are OUTSIDE the AOI
@@ -836,9 +840,6 @@ verbose: {verbose}")
iml_grid_hd[iml_grid_hd == 0.0] = np.nan # change zeroes back to nan iml_grid_hd[iml_grid_hd == 0.0] = np.nan # change zeroes back to nan
# correct orientation is AOI is used
if use_AOI:
iml_grid_hd = np.rot90(iml_grid_hd, k=-1)
#vmin_hd = min(x for x in iml_grid_hd.flatten() if not isnan(x)) #vmin_hd = min(x for x in iml_grid_hd.flatten() if not isnan(x))
vmax_hd = max(x for x in iml_grid_hd.flatten() if not isnan(x)) vmax_hd = max(x for x in iml_grid_hd.flatten() if not isnan(x))