Update src/seismic_hazard_forecasting.py

temporary global variables since GUI not ready
This commit is contained in:
2026-06-19 15:38:35 +02:00
parent fb00131997
commit 037863e975
+12 -12
View File
@@ -5,6 +5,14 @@ from scipy.stats import bootstrap
from matplotlib.dates import DateFormatter, AutoDateLocator from matplotlib.dates import DateFormatter, AutoDateLocator
import numpy as np import numpy as np
global ncp_choice, tcp_max, torder_min, torder_max, AOI_lat, AOI_lon
ncp_choice = 'default'
tcp_max = 5
torder_min = 0
torder_max = 1
AOI_lat = np.array([51.48, 51.54])
AOI_lon = np.array([16.15, 16.24])
def plot_results(act_rate, bin_edges, bin_edges_dt, rt, boundaries, def plot_results(act_rate, bin_edges, bin_edges_dt, rt, boundaries,
bin_dur, unit, multiplicator, bin_dur, unit, multiplicator,
rate_forecast, rate_unc_high, rate_unc_low, rate_forecast, rate_unc_high, rate_unc_low,
@@ -130,7 +138,7 @@ def calc_rates(act_rate, cps):
for i in range(len(boundaries)-1)] for i in range(len(boundaries)-1)]
return rt, boundaries return rt, boundaries
def apply_beast(act_rate, *args): def apply_beast(act_rate):
""" """
Applies BEAST to the smmothed rate data using different smoothing windows. Applies BEAST to the smmothed rate data using different smoothing windows.
Input Input
@@ -183,7 +191,7 @@ def apply_beast(act_rate, *args):
return beast_result, np.sort(cps) return beast_result, np.sort(cps)
def bins_and_beast(dates, unit, bin_dur, multiplicator, *args): def bins_and_beast(dates, unit, bin_dur, multiplicator):
start_date = dates.min() start_date = dates.min()
end_date = dates.max() end_date = dates.max()
@@ -224,7 +232,7 @@ def bins_and_beast(dates, unit, bin_dur, multiplicator, *args):
act_rate = [count / ((bin_edges[i + 1] - bin_edges[i]) * multiplicator / bin_dur) act_rate = [count / ((bin_edges[i + 1] - bin_edges[i]) * multiplicator / bin_dur)
for i, count in enumerate(bin_counts)] for i, count in enumerate(bin_counts)]
out, cps = apply_beast(act_rate, *args) out, cps = apply_beast(act_rate)
if len(cps) > 0: if len(cps) > 0:
rt, boundaries = calc_rates(act_rate, cps) rt, boundaries = calc_rates(act_rate, cps)
@@ -334,9 +342,6 @@ def main(catalog_file, mc_file, pdf_file, m_file, m_select, mag_label, mc, m_max
# AOI_lat = np.array(AOI_extent[:2]) # AOI_lat = np.array(AOI_extent[:2])
# AOI_lon = np.array(AOI_extent[2:]) # AOI_lon = np.array(AOI_extent[2:])
AOI_lat = np.array([51.48, 51.54])
AOI_lon = np.array([16.15, 16.24])
# log user selections # log user selections
logger.debug(f"User input files\n Catalog: {catalog_file}\n Mc: {mc_file}\n Mag_PDF: {pdf_file}\n Mag: {m_file}") logger.debug(f"User input files\n Catalog: {catalog_file}\n Mc: {mc_file}\n Mag_PDF: {pdf_file}\n Mag: {m_file}")
@@ -587,11 +592,6 @@ verbose: {verbose}")
elif rate_select: elif rate_select:
logger.info(f"Activity rate modeling selected") logger.info(f"Activity rate modeling selected")
ncp_choice = 'default'
tcp_max = 5
torder_min = 0
torder_max = 1
time, mag_data, lat_dummy, lon_dummy, depth_dummy = read_mat_cat(catalog_file, mag_label=mag_label, output_datenum=True) time, mag_data, lat_dummy, lon_dummy, depth_dummy = read_mat_cat(catalog_file, mag_label=mag_label, output_datenum=True)
datenum_data = time # REMEMBER THE DECIMAL DENOTES DAYS datenum_data = time # REMEMBER THE DECIMAL DENOTES DAYS
@@ -617,7 +617,7 @@ verbose: {verbose}")
#-------split the data into bins and apply BEAST for changepoint detection-------------------- #-------split the data into bins and apply BEAST for changepoint detection--------------------
act_rate, bin_counts, bin_edges, bin_edges_dt, out, cps, rt, boundaries, bin_dur, time_unit = bins_and_beast( act_rate, bin_counts, bin_edges, bin_edges_dt, out, cps, rt, boundaries, bin_dur, time_unit = bins_and_beast(
np.array(datenum_data), time_unit, time_win_duration, multiplicator, tcp_max, torder_min, torder_max) np.array(datenum_data), time_unit, time_win_duration, multiplicator)
#------Forecasted rate is taken from BEAST or is equal to last value if no changepoints detected----- #------Forecasted rate is taken from BEAST or is equal to last value if no changepoints detected-----
if len(cps) > 0: if len(cps) > 0: