diff --git a/modules/performRegionalEventSimulation/regionalGroundMotion/ComputeIntensityMeasure.py b/modules/performRegionalEventSimulation/regionalGroundMotion/ComputeIntensityMeasure.py index 4d0e9c881..3f18e309b 100644 --- a/modules/performRegionalEventSimulation/regionalGroundMotion/ComputeIntensityMeasure.py +++ b/modules/performRegionalEventSimulation/regionalGroundMotion/ComputeIntensityMeasure.py @@ -1106,21 +1106,25 @@ def export_im( # noqa: C901, D103, PLR0912 json.dump(res, f, indent=2) # export the event grid and station csv files if csv_flag: - # output EventGrid.csv + # output EventGrid.csv. Site characterization values (Vs30, z1pt0, + # z2pt5) come from the SimCenterSiteModel.csv that CreateStation + # wrote earlier in the workflow station_name = [ 'site' + str(stations[j]['ID']) + '.csv' for j in range(len(stations)) ] lat = [stations[j]['lat'] for j in range(len(stations))] lon = [stations[j]['lon'] for j in range(len(stations))] - # vs30 = [stations[j]['vs30'] for j in range(len(stations))] - # zTR = [stations[j]['DepthToRock'] for j in range(len(stations))] + vs30 = [stations[j].get('vs30') for j in range(len(stations))] + z1pt0 = [stations[j].get('z1pt0') for j in range(len(stations))] + z2pt5 = [stations[j].get('z2pt5') for j in range(len(stations))] df = pd.DataFrame( # noqa: PD901 { 'GP_file': station_name, 'Longitude': lon, 'Latitude': lat, - # 'Vs30': vs30, - # 'DepthToRock': zTR + 'Vs30': vs30, + 'z1pt0': z1pt0, + 'z2pt5': z2pt5, } ) # if cur_eq[2]: diff --git a/modules/performRegionalEventSimulation/regionalGroundMotion/CreateStation.py b/modules/performRegionalEventSimulation/regionalGroundMotion/CreateStation.py index c505a25c3..28ce25457 100644 --- a/modules/performRegionalEventSimulation/regionalGroundMotion/CreateStation.py +++ b/modules/performRegionalEventSimulation/regionalGroundMotion/CreateStation.py @@ -82,6 +82,39 @@ def get_label(options, labels, label_name): # noqa: D103 print(f'WARNING: Could not identify the label for the {label_name}') # noqa: T201, RET503 +# Canonical OpenSHA basin-depth provider NAMEs (matches `getSourceName()`). +# Used to recognize Z1.0 / Z2.5 Type strings that select a specific OpenSHA +# basin-depth model rather than the default sequenced fallback. +_OPENSHA_BASIN_DEPTH_NAMES = frozenset([ + 'USGS SF Bay Area Velocity Model Release 21.1', + 'USGS Bay Area Velocity Model Release 8.3.0', + 'SCEC Community Velocity Model Version 4, Iteration 26, Basin Depth', + 'SCEC Community Velocity Model Version 4, Iteration 26, M01 w/ Taper, Basin Depth', + 'SCEC CCA, Iteration 6, Basin Depth', + 'SCEC Community Velocity Model Version 4 Basin Depth', + 'SCEC/Harvard Community Velocity Model Version 11.9.x Basin Depth', + 'SCEC CyberShake Study 18.8 Stitched Basin Depth', + 'SCEC CyberShake Study 24.8 Stitched Basin Depth', + 'SCEC Community Velocity Model Version 2 Basin Depth', +]) + + +def _z_opensha_dispatch(z_config, tag_key): + """Decide whether to fetch Z1.0/Z2.5 from OpenSHA and which model. + + Returns a (model_name, should_fetch) tuple: + - ('OpenSHA default model' Type with tag==2) → (None, True): sequenced + - (Type is a basin-depth NAME) → (NAME, True): specific + - (anything else) → (None, False): no fetch + """ + cfg_type = z_config.get('Type', '') + if cfg_type in _OPENSHA_BASIN_DEPTH_NAMES: + return cfg_type, True + if cfg_type == 'OpenSHA default model' and z_config.get(tag_key) == 2: # noqa: PLR2004 + return None, True + return None, False + + class Station: """A class for stations in an earthquake scenario""" # noqa: D400 @@ -253,46 +286,33 @@ def create_stations( # noqa: C901, PLR0912, PLR0915 nan_loc = [] else: nan_loc = list(range(len(selected_stn.index))) - if 'Global Vs30' in vs30Config['Type']: - vs30_tag = 1 - elif 'Thompson' in vs30Config['Type']: - vs30_tag = 2 - elif 'NCM' in vs30Config['Type']: - vs30_tag = 3 - else: - vs30_tag = 0 - if len(nan_loc) and vs30_tag == 1: - print('CreateStation: Interpolating global Vs30 map for defined stations.') # noqa: T201 - selected_stn.loc[nan_loc, vs30_label] = get_vs30_global( - selected_stn.iloc[ # noqa: PD011 - nan_loc, list(selected_stn.keys()).index(lat_label) - ].values.tolist(), - selected_stn.iloc[ # noqa: PD011 - nan_loc, list(selected_stn.keys()).index(lon_label) - ].values.tolist(), - ) - if len(nan_loc) and vs30_tag == 2: # noqa: PLR2004 - print('CreateStation: Interpolating Thompson Vs30 map for defined stations.') # noqa: T201 - selected_stn.loc[nan_loc, vs30_label] = get_vs30_thompson( - selected_stn.iloc[ # noqa: PD011 - nan_loc, list(selected_stn.keys()).index(lat_label) - ].values.tolist(), - selected_stn.iloc[ # noqa: PD011 - nan_loc, list(selected_stn.keys()).index(lon_label) - ].values.tolist(), - ) - if len(nan_loc) and vs30_tag == 3: # noqa: PLR2004 - print('CreateStation: Fetch National Crustal Model Vs for defined stations.') # noqa: T201 - selected_stn.loc[nan_loc, vs30_label] = get_vs30_ncm( - selected_stn.iloc[ # noqa: PD011 - nan_loc, list(selected_stn.keys()).index(lat_label) - ].values.tolist(), - selected_stn.iloc[ # noqa: PD011 - nan_loc, list(selected_stn.keys()).index(lon_label) - ].values.tolist(), + # Vs30 dispatch. + # + # All non-User-specified Vs30 paths route through OpenSHA. The + # vs30Config['Type'] string is either a canonical OpenSHA provider + # NAME (e.g. 'CGS/Wills VS30 Map (2015)') or a legacy substring + # ('Global Vs30 ...', 'Thompson ...') that get_site_vs30_from_opensha + # maps to a canonical NAME for backward compatibility. + # + # The NCM branch below is not exposed in the R2D UI; the underlying + # USGS web service has not been verified since 2021 (see get_vs30_ncm). + if 'NCM' in vs30Config['Type']: + # Not reachable from the R2D UI + if len(nan_loc): + print('CreateStation: Fetch National Crustal Model Vs for defined stations.') # noqa: T201 + selected_stn.loc[nan_loc, vs30_label] = get_vs30_ncm( + selected_stn.iloc[ # noqa: PD011 + nan_loc, list(selected_stn.keys()).index(lat_label) + ].values.tolist(), + selected_stn.iloc[ # noqa: PD011 + nan_loc, list(selected_stn.keys()).index(lon_label) + ].values.tolist(), + ) + elif len(nan_loc): + print( # noqa: T201 + f'CreateStation: Fetch Vs30 from OpenSHA ' + f'({vs30Config["Type"]}) for defined stations.' ) - if len(nan_loc) and vs30_tag == 0: - print('CreateStation: Fetch OpenSHA Vs30 map for defined stations.') # noqa: T201 selected_stn.loc[nan_loc, vs30_label] = get_site_vs30_from_opensha( selected_stn.iloc[ # noqa: PD011 nan_loc, list(selected_stn.keys()).index(lat_label) @@ -300,6 +320,7 @@ def create_stations( # noqa: C901, PLR0912, PLR0915 selected_stn.iloc[ # noqa: PD011 nan_loc, list(selected_stn.keys()).index(lon_label) ].values.tolist(), + vs30model=vs30Config['Type'], ) # Get zTR @@ -456,55 +477,65 @@ def create_stations( # noqa: C901, PLR0912, PLR0915 'chi', ]: user_param_list.pop(user_param_list.index(cur_param)) - # If z1pt0 is OpenSHA default model, use parallel processing to get z1pt0 - if z1Config['Type'] == 'OpenSHA default model': - z1_tag = z1Config['z1_tag'] - if z1_tag == 2: # noqa: PLR2004 - # num_cores = z1Config.get('num_cores', multiprocessing.cpu_count()) - num_cores = z1Config.get('num_cores', 1) - if num_cores == 1: - z1pt0_results = [ - get_site_z1pt0_from_opensha(lat, lon) + # Z1.0 / Z2.5 fetch from OpenSHA. + # + # z*Config['Type'] is one of: + # - 'User-specified' → values come from the site file or a constant + # - 'OpenSHA default model' → legacy: uses z*_tag for empirical / sequenced / null + # - → fetch that specific provider's value + # + # For both 'OpenSHA default model' (with tag=2) and a named OpenSHA + # model, we fetch all values up front (optionally in parallel) so the + # per-station loop below can just index into the results array. + _z1_model, _z1_fetch = _z_opensha_dispatch(z1Config, 'z1_tag') + if _z1_fetch: + num_cores = z1Config.get('num_cores', 1) + if num_cores == 1: + z1pt0_results = [ + get_site_z1pt0_from_opensha(lat, lon, model_name=_z1_model) + for lat, lon in zip( + selected_stn['Latitude'].tolist(), + selected_stn['Longitude'].tolist(), + ) + ] + else: + with tqdm_joblib( + tqdm(desc='Get z1pt0 from openSHA', total=selected_stn.shape[0]) + ) as progress_bar: # noqa: F841 + z1pt0_results = Parallel(n_jobs=num_cores)( + delayed(get_site_z1pt0_from_opensha)( + lat, lon, model_name=_z1_model + ) for lat, lon in zip( selected_stn['Latitude'].tolist(), selected_stn['Longitude'].tolist(), ) - ] - else: - with tqdm_joblib( - tqdm(desc='Get z1pt0 from openSHA', total=selected_stn.shape[0]) - ) as progress_bar: - z1pt0_results = Parallel(n_jobs=num_cores)( - delayed(get_site_z1pt0_from_opensha)(lat, lon) - for lat, lon in zip( - selected_stn['Latitude'].tolist(), - selected_stn['Longitude'].tolist(), - ) + ) + + _z25_model, _z25_fetch = _z_opensha_dispatch(z25Config, 'z25_tag') + if _z25_fetch: + num_cores = z25Config.get('num_cores', 1) + if num_cores == 1: + z2pt5_results = [ + get_site_z2pt5_from_opensha(lat, lon, model_name=_z25_model) + for lat, lon in zip( + selected_stn['Latitude'].tolist(), + selected_stn['Longitude'].tolist(), + ) + ] + else: + with tqdm_joblib( + tqdm(desc='Get z2pt5 from openSHA', total=selected_stn.shape[0]) + ) as progress_bar: # noqa: F841 + z2pt5_results = Parallel(n_jobs=num_cores)( + delayed(get_site_z2pt5_from_opensha)( + lat, lon, model_name=_z25_model ) - if z25Config['Type'] == 'OpenSHA default model': - z25_tag = z25Config['z25_tag'] - if z25_tag == 2: # noqa: PLR2004 - # num_cores = z25Config.get('num_cores', multiprocessing.cpu_count()) - num_cores = z25Config.get('num_cores', 1) - if num_cores == 1: - z2pt5_results = [ - get_site_z2pt5_from_opensha(lat, lon) for lat, lon in zip( selected_stn['Latitude'].tolist(), selected_stn['Longitude'].tolist(), ) - ] - else: - with tqdm_joblib( - tqdm(desc='Get z2pt5 from openSHA', total=selected_stn.shape[0]) - ) as progress_bar: # noqa: F841 - z2pt5_results = Parallel(n_jobs=num_cores)( - delayed(get_site_z2pt5_from_opensha)(lat, lon) - for lat, lon in zip( - selected_stn['Latitude'].tolist(), - selected_stn['Longitude'].tolist(), - ) - ) + ) ground_failure_input_keys = set() for ind in tqdm(range(selected_stn.shape[0]), desc='Stations'): @@ -559,6 +590,14 @@ def create_stations( # noqa: C901, PLR0912, PLR0915 elif z1_tag == 0: z1pt0 = get_z1(tmp.get('Vs30')) tmp.update({'z1pt0': z1pt0}) + elif z1Config['Type'] in _OPENSHA_BASIN_DEPTH_NAMES: + # Specific OpenSHA basin-depth model; fall back to the + # empirical Chiou & Youngs (2013) equation if the chosen + # model has no data at this site. + z1pt0 = z1pt0_results[ind] + if np.isnan(z1pt0): + z1pt0 = get_z1(tmp.get('Vs30')) + tmp.update({'z1pt0': z1pt0}) if (z25Config['Type'] == 'User-specified') and stn.get('z2pt5'): tmp.update({'z2pt5': stn.get('z2pt5')}) @@ -578,6 +617,14 @@ def create_stations( # noqa: C901, PLR0912, PLR0915 elif z25_tag == 0: z2pt5 = get_z25(tmp['z1pt0']) tmp.update({'z2pt5': z2pt5}) + elif z25Config['Type'] in _OPENSHA_BASIN_DEPTH_NAMES: + # Specific OpenSHA basin-depth model; fall back to the + # empirical Campbell & Bozorgnia (2013) equation if the + # chosen model has no data at this site. + z2pt5 = z2pt5_results[ind] + if np.isnan(z2pt5): + z2pt5 = get_z25(tmp['z1pt0']) + tmp.update({'z2pt5': z2pt5}) if 'DepthToRock' in stn.index: tmp.update({'DepthToRock': stn.get('DepthToRock')}) @@ -713,32 +760,6 @@ def create_gridded_stations( ) -def get_vs30_global(lat, lon): - """Interpolate global Vs30 at given latitude and longitude - Input: - lat: list of latitude - lon: list of longitude - Output: - vs30: list of vs30 - """ # noqa: D205, D400 - import os - import pickle - - from scipy import interpolate - - # Loading global Vs30 data - cwd = os.path.dirname(os.path.realpath(__file__)) # noqa: PTH120 - with open(cwd + '/database/site/global_vs30_4km.pkl', 'rb') as f: # noqa: PTH123 - vs30_global = pickle.load(f) # noqa: S301 - # Interpolation function (linear) - interpFunc = interpolate.interp2d( # noqa: N806 - vs30_global['Longitude'], vs30_global['Latitude'], vs30_global['Vs30'] - ) - vs30 = [float(interpFunc(x, y)) for x, y in zip(lon, lat)] - # return - return vs30 # noqa: DOC201, RET504, RUF100 - - def parallel_interpolation(func, lat, lon): """Interpolate data in parallel Input: @@ -751,40 +772,6 @@ def parallel_interpolation(func, lat, lon): return func(lat, lon) -def get_vs30_thompson(lat, lon): - """Interpolate global Vs30 at given latitude and longitude - Input: - lat: list of latitude - lon: list of longitude - Output: - vs30: list of vs30 - """ # noqa: D205, D400 - import os - import pickle - - from scipy import interpolate - - # Loading Thompson Vs30 data - cwd = os.path.dirname(os.path.realpath(__file__)) # noqa: PTH120 - with open(cwd + '/database/site/thompson_vs30_4km.pkl', 'rb') as f: # noqa: PTH123 - vs30_thompson = pickle.load(f) # noqa: S301 - # Interpolation function (linear) - vs30_thompson['Vs30'][vs30_thompson['Vs30'] < 0.1] = 760 # noqa: PLR2004 - interpFunc = interpolate.interp2d( # noqa: N806 - vs30_thompson['Longitude'], vs30_thompson['Latitude'], vs30_thompson['Vs30'] - ) - vs30 = [float(interpFunc(x, y)) for x, y in zip(lon, lat)] - - num_zeros = len([x for x in vs30 if x == 0]) - if num_zeros > 0: - # Thompson's map gives zero values for water-covered region and outside CA -> use 760 for default - print( # noqa: T201 - f'CreateStation: Warning - approximate 760 m/s for {num_zeros} sites not supported by Thompson Vs30 map (water/outside CA).' - ) - # return - return vs30 # noqa: DOC201, RET504, RUF100 - - def get_z1(vs30): """Compute z1 based on the prediction equation by Chiou and Youngs (2013) (unit of vs30 is meter/second and z1 is meter).""" return np.exp( @@ -866,6 +853,13 @@ def export_site_prop(stn_file, output_dir, filename): def get_zTR_ncm(lat, lon): # noqa: N802 """Call USGS National Crustal Model services for zTR https://earthquake.usgs.gov/nshmp/ncm + + NOTE: not currently reachable from the R2D UI (no UI option produces + a BedrockDepth Type containing 'NCM' / 'National Crustal Model'). The + request URL below appears to repeat `/ws/nshmp/ncm/` — looks like a + pre-existing typo from 2021; the endpoint may have moved or been + retired. Verify before re-exposing. + Input: lat: list of latitude lon: list of longitude @@ -878,6 +872,7 @@ def get_zTR_ncm(lat, lon): # noqa: N802 # Looping over sites for cur_lat, cur_lon in zip(lat, lon): + # `/ws/nshmp/ncm/ws/nshmp/ncm/` looks duplicated — likely a typo. url_geology = f'https://earthquake.usgs.gov/ws/nshmp/ncm/ws/nshmp/ncm/geologic-framework?location={cur_lat}%2C{cur_lon}' # geological data (depth to bedrock) r1 = requests.get(url_geology) # noqa: S113 @@ -900,6 +895,12 @@ def get_zTR_ncm(lat, lon): # noqa: N802 def get_vsp_ncm(lat, lon, depth): """Call USGS National Crustal Model services for Vs30 profile https://earthquake.usgs.gov/nshmp/ncm + + NOTE: only reachable through `get_vs30_ncm`, which is itself not + exposed in the R2D UI (no UI option produces a Vs30 Type containing + 'NCM'). The request URL has the same `/ws/nshmp/ncm/` repetition as + `get_zTR_ncm`; endpoint status unverified. + Input: lat: list of latitude lon: list of longitude @@ -914,6 +915,7 @@ def get_vsp_ncm(lat, lon, depth): # Looping over sites for cur_lat, cur_lon in zip(lat, lon): + # `/ws/nshmp/ncm/ws/nshmp/ncm/` looks duplicated — likely a typo. url_geophys = f'https://earthquake.usgs.gov/ws/nshmp/ncm/ws/nshmp/ncm/geophysical?location={cur_lat}%2C{cur_lon}&depths={depthMin}%2C{depthInc}%2C{depthMax}' r1 = requests.get(url_geophys) # noqa: S113 cur_res = r1.json() @@ -955,7 +957,13 @@ def compute_vs30_from_vsp(depthp, vsp): def get_vs30_ncm(lat, lon): - """Fetch Vs30 at given latitude and longitude from NCM + """Fetch Vs30 at given latitude and longitude from NCM. + + NOTE: not currently reachable from the R2D UI — the dispatch in the + Vs30 section triggers on `'NCM' in vs30Config['Type']`, but no R2D + Vs30 option produces such a Type string. The underlying USGS + /ws/nshmp/ncm endpoint status is unverified. + Input: lat: list of latitude lon: list of longitude diff --git a/modules/performRegionalEventSimulation/regionalGroundMotion/FetchOpenSHA.py b/modules/performRegionalEventSimulation/regionalGroundMotion/FetchOpenSHA.py index dcffbe103..a68667726 100644 --- a/modules/performRegionalEventSimulation/regionalGroundMotion/FetchOpenSHA.py +++ b/modules/performRegionalEventSimulation/regionalGroundMotion/FetchOpenSHA.py @@ -1006,7 +1006,54 @@ def get_IM( # noqa: C901, N802, D103 return res, station_info +# Canonical OpenSHA Vs30 provider NAME strings (matches `getSourceName()`). +# Used to validate / normalize user-supplied Vs30 Type strings. +_OPENSHA_VS30_NAMES = ( + 'Thompson VS30 Map (2022)', + 'Thompson VS30 Map (2020)', + 'Thompson VS30 Map (2018)', + 'CGS/Wills VS30 Map (2015)', + 'CGS/Wills Site Classification Map (2006)', + 'Global Vs30 from Topographic Slope (Wald & Allen 2008)', + 'Vs30 from various Community Velocity Models', +) + + +def _resolve_vs30_model_name(vs30model): + """Map a Vs30 Type string to a canonical OpenSHA provider NAME. + + Accepts a canonical NAME, a canonical NAME with a UI scope suffix + (e.g. " [California only]"), or a legacy substring from older R2D + configurations. Returns the canonical NAME so the rest of the code + can do exact name matching against `getSourceName()`. + """ + if vs30model in _OPENSHA_VS30_NAMES: + return vs30model + # Match by prefix to accept UI labels like " [California only]". + for canonical in _OPENSHA_VS30_NAMES: + if vs30model.startswith(canonical): + return canonical + # Legacy substrings (pre-2026 R2D configurations). + if 'Global Vs30' in vs30model: + return 'Global Vs30 from Topographic Slope (Wald & Allen 2008)' + if 'Wills' in vs30model and '2015' in vs30model: + return 'CGS/Wills VS30 Map (2015)' + if 'Wills' in vs30model and '2006' in vs30model: + return 'CGS/Wills Site Classification Map (2006)' + if 'Thompson' in vs30model and '2018' in vs30model: + return 'Thompson VS30 Map (2018)' + if 'Thompson' in vs30model and '2020' in vs30model: + return 'Thompson VS30 Map (2020)' + if 'Thompson' in vs30model: + # Unqualified or 2022-era "Thompson" string → newest map. + return 'Thompson VS30 Map (2022)' + # Fall through with the original string; OpenSHA will simply find no + # matching provider and the caller will see NaN values. + return vs30model + + def get_site_vs30_from_opensha(lat, lon, vs30model='CGS/Wills VS30 Map (2015)'): # noqa: D103 + vs30model = _resolve_vs30_model_name(vs30model) sites = ArrayList() # noqa: F405 num_sites = len(lat) for i in range(num_sites): @@ -1026,10 +1073,11 @@ def get_site_vs30_from_opensha(lat, lon, vs30model='CGS/Wills VS30 Map (2015)'): else: # noqa: RET508 continue - # The Wills 2015 map returns NaN for offshore sites; fall back to the - # global topographic-slope Vs30 model. Look up the fallback provider by - # name so an upstream reorder of the provider list does not silently - # route this code to the wrong dataset. + # The selected Vs30 map may return NaN for cells outside its coverage + # (e.g. Wills 2015 / Thompson outside California, offshore cells). Fall + # back to the global topographic-slope model for those sites. Look up + # the fallback provider by name so an upstream reorder of the provider + # list does not silently route this code to the wrong dataset. if any([np.isnan(x) for x in vs30]): # noqa: C419 fallback_name = 'Global Vs30 from Topographic Slope (Wald & Allen 2008)' fallback_provider_data = None @@ -1050,27 +1098,45 @@ def get_site_vs30_from_opensha(lat, lon, vs30model='CGS/Wills VS30 Map (2015)'): return vs30 -def get_site_z1pt0_from_opensha(lat, lon): # noqa: D103 +def _fetch_z_from_opensha(lat, lon, data_type, model_name): + """Fetch Z1.0 or Z2.5 (in km) at a single site from OpenSHA. + + When `model_name` is None or 'OpenSHA default model', returns the value + from the first provider in OpenSHA's default order that has non-NaN + data for the site (sequenced fallback). When `model_name` is a specific + OpenSHA basin-depth provider NAME (matching `getSourceName()`), returns + the value from that provider only; NaN if the provider has no data + for the site or is not in the list. + """ sites = ArrayList() # noqa: F405 sites.add(Site(Location(lat, lon))) # noqa: F405 siteDataProviders = OrderedSiteDataProviderList.createSiteDataProviderDefaults() # noqa: N806, F405 siteData = siteDataProviders.getAllAvailableData(sites) # noqa: N806 + + use_default_sequence = ( + model_name is None or model_name == 'OpenSHA default model' + ) + z = float('nan') for data in siteData: - if data.getValue(0).getDataType() == 'Depth to Vs = 1.0 km/sec': - z1pt0 = float(data.getValue(0).getValue()) - if not np.isnan(z1pt0): - break - return z1pt0 * 1000.0 + if str(data.getValue(0).getDataType()) != data_type: + continue + if not use_default_sequence and str(data.getSourceName()) != model_name: + continue + candidate = float(data.getValue(0).getValue()) + if not np.isnan(candidate): + z = candidate + break + if not use_default_sequence: + # Asked for a specific model and it returned NaN here. + break + return z -def get_site_z2pt5_from_opensha(lat, lon): # noqa: D103 - sites = ArrayList() # noqa: F405 - sites.add(Site(Location(lat, lon))) # noqa: F405 - siteDataProviders = OrderedSiteDataProviderList.createSiteDataProviderDefaults() # noqa: N806, F405 - siteData = siteDataProviders.getAllAvailableData(sites) # noqa: N806 - for data in siteData: - if data.getValue(0).getDataType() == 'Depth to Vs = 2.5 km/sec': - z2pt5 = float(data.getValue(0).getValue()) - if not np.isnan(z2pt5): - break - return z2pt5 * 1000.0 +def get_site_z1pt0_from_opensha(lat, lon, model_name=None): # noqa: D103 + z = _fetch_z_from_opensha(lat, lon, 'Depth to Vs = 1.0 km/sec', model_name) + return z * 1000.0 + + +def get_site_z2pt5_from_opensha(lat, lon, model_name=None): # noqa: D103 + z = _fetch_z_from_opensha(lat, lon, 'Depth to Vs = 2.5 km/sec', model_name) + return z * 1000.0 diff --git a/modules/performRegionalEventSimulation/regionalGroundMotion/ScenarioForecast.py b/modules/performRegionalEventSimulation/regionalGroundMotion/ScenarioForecast.py index a80670221..bcac4af9d 100644 --- a/modules/performRegionalEventSimulation/regionalGroundMotion/ScenarioForecast.py +++ b/modules/performRegionalEventSimulation/regionalGroundMotion/ScenarioForecast.py @@ -122,11 +122,11 @@ # shutil.rmtree(os.environ.get('OQ_DATADIR')) # os.makedirs(f"{os.environ.get('OQ_DATADIR')}") - # untar site databases + # untar site databases. global_zTR_4km is the depth-to-bedrock map used + # by ground-failure code (liquefaction/landslide); all Vs30 paths + # go through OpenSHA, so no Vs30 pickles need to be staged here. site_database = [ - 'global_vs30_4km.tar.gz', 'global_zTR_4km.tar.gz', - 'thompson_vs30_4km.tar.gz', ] print('HazardSimulation: Extracting site databases.') # noqa: T201 cwd = os.path.dirname(os.path.realpath(__file__)) # noqa: PTH120 diff --git a/modules/performRegionalEventSimulation/regionalGroundMotion/database/CMakeLists.txt b/modules/performRegionalEventSimulation/regionalGroundMotion/database/CMakeLists.txt index 0276b2ded..79a1c1cef 100644 --- a/modules/performRegionalEventSimulation/regionalGroundMotion/database/CMakeLists.txt +++ b/modules/performRegionalEventSimulation/regionalGroundMotion/database/CMakeLists.txt @@ -1,2 +1,3 @@ add_subdirectory(gmdb) add_subdirectory(groundfailure) +add_subdirectory(site) diff --git a/modules/performRegionalEventSimulation/regionalGroundMotion/database/site/CMakeLists.txt b/modules/performRegionalEventSimulation/regionalGroundMotion/database/site/CMakeLists.txt new file mode 100644 index 000000000..22e255992 --- /dev/null +++ b/modules/performRegionalEventSimulation/regionalGroundMotion/database/site/CMakeLists.txt @@ -0,0 +1,2 @@ +simcenter_add_file(NAME global_vs30_4km.tar.gz) +simcenter_add_file(NAME global_zTR_4km.tar.gz) diff --git a/modules/performRegionalEventSimulation/regionalGroundMotion/database/site/thompson_vs30_4km.tar.gz b/modules/performRegionalEventSimulation/regionalGroundMotion/database/site/thompson_vs30_4km.tar.gz deleted file mode 100644 index 439c6c914..000000000 Binary files a/modules/performRegionalEventSimulation/regionalGroundMotion/database/site/thompson_vs30_4km.tar.gz and /dev/null differ