Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Fix shakemap, edits to allow updated HAZUS Bridges #348

Merged
merged 9 commits into from
Nov 13, 2024
4 changes: 2 additions & 2 deletions modules/Workflow/WorkflowApplications.json
Original file line number Diff line number Diff line change
Expand Up @@ -285,12 +285,12 @@
"ApplicationSpecificInputs": [
{
"id": "Directory",
"type": "path",
"type": "string",
"description": "Path to file containing folder of shake maps"
},
{
"id": "EventPath",
"type": "string",
"type": "path",
"description": "Path to the shake map event"
},
{
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -259,7 +259,7 @@ def read_pressure_data(file_names):
# index += 1

except: # noqa: E722
# sys.exit('Fatal Error!: the pressure filese have time gap')
# sys.exit('Fatal Error!: the pressure files have time gap')
index = 0 # Joint them even if they have a time gap

connected_time = np.concatenate((connected_time, time2[index:]))
Expand Down
41 changes: 30 additions & 11 deletions modules/createEVENT/shakeMapEvent/shakeMapEvent.py
Original file line number Diff line number Diff line change
Expand Up @@ -44,12 +44,17 @@
from pathlib import Path


class IntensityMeasureTypeError(Exception):
def __init__(self, im_type):
super().__init__(f'Intensity measure type {im_type} not found in grid data')


def create_shakemap_event(eventDirectory, eventPath, IMTypes): # noqa: D103, N803
IMTypesList = eval(IMTypes) # noqa: S307, N806

print('Creating shakemap event') # noqa: T201

xml_file_path = Path(eventDirectory) / eventPath / 'grid.xml'
xml_file_path = Path(eventPath) / 'grid.xml'

# Parse the XML file
tree = ET.parse(xml_file_path) # noqa: S314
Expand All @@ -62,6 +67,19 @@ def create_shakemap_event(eventDirectory, eventPath, IMTypes): # noqa: D103, N8
points = []
attributes = []

# Get the attribute_mapping
namespace = {'ns': 'http://earthquake.usgs.gov/eqcenter/shakemap'}
grid_fields = {}
for grid_field in root.findall('ns:grid_field', namespace):
index = grid_field.get('index')
name = grid_field.get('name')
units = grid_field.get('units')
grid_fields[name] = {'index': index, 'units': units}
attribute_mapping = {}
for im_type in ['PGA', 'PGV', 'MMI', 'PSA03', 'PSA10', 'PSA30']:
if im_type not in grid_fields:
raise IntensityMeasureTypeError(im_type)
attribute_mapping[im_type] = int(grid_fields[im_type]['index']) - 1
# Parse the grid data
for line in grid_data.text.strip().split('\n'):
values = line.split()
Expand All @@ -71,15 +89,6 @@ def create_shakemap_event(eventDirectory, eventPath, IMTypes): # noqa: D103, N8

# Store only the specified attributes
attr = {}
attribute_mapping = {
'PGA': 2,
'PGV': 3,
'MMI': 4,
'PSA03': 5,
'PSA10': 6,
'PSA30': 7,
}

for im_type in IMTypesList:
if im_type in attribute_mapping:
attr[im_type] = float(values[attribute_mapping[im_type]])
Expand All @@ -89,11 +98,21 @@ def create_shakemap_event(eventDirectory, eventPath, IMTypes): # noqa: D103, N8
# Create GeoDataFrame
gdf = gpd.GeoDataFrame(attributes, geometry=points, crs='EPSG:4326')

if 'PGA' in gdf.columns:
gdf['PGA'] = gdf['PGA'] / 100 # convert from pct g to g

if 'PSA03' in gdf.columns:
gdf['PSA03'] = gdf['PSA03'] / 100 # convert from pct g to g
gdf = gdf.rename(columns={'PSA03': 'SA_0.3'})
if 'PSA10' in gdf.columns:
gdf['PSA10'] = gdf['PSA10'] / 100 # convert from pct g to g
gdf = gdf.rename(columns={'PSA10': 'SA_1.0'})

# Display the first few rows
print('Saving shakemap to gpkg') # noqa: T201

# Save as a GeoPackage file
gdf_path = Path(eventDirectory) / 'EventGrid.gpkg'
gdf_path = Path(eventPath) / 'EventGrid.gpkg'
gdf.to_file(gdf_path, driver='GPKG')

return # noqa: PLR1711
Expand Down
37 changes: 20 additions & 17 deletions modules/performRegionalMapping/NearestNeighborEvents/NNE.py
Original file line number Diff line number Diff line change
Expand Up @@ -321,7 +321,8 @@ def find_neighbors( # noqa: C901, D103
for col in grid_df.columns
if col not in ['geometry', 'Longitude', 'Latitude']
]
event_count = len(im_columns)
# event_count = len(im_columns)
event_count = 1

# for each neighbor
for sample_j, nbr in enumerate(nbr_samples):
Expand All @@ -332,26 +333,28 @@ def find_neighbors( # noqa: C901, D103
nbr_index = ind_list[nbr]

# For GIS files, create a new CSV file
csv_filename = f'Site_{sample_j}.csv'
csv_filename = f'Site_{nbr_index}.csv'

csv_path = event_dir / csv_filename

# Create a CSV file with data from the GIS file
# Use actual data from the GIS file if available, otherwise use dummy data
im_columns = [
col
for col in grid_df.columns
if col not in ['geometry', 'Longitude', 'Latitude']
]

im_data = pd.DataFrame(
{
col: [grid_df.iloc[nbr_index][col]] * event_count
for col in im_columns
}
)
if not csv_path.exists():

# Create a CSV file with data from the GIS file
# Use actual data from the GIS file if available, otherwise use dummy data
im_columns = [
col
for col in grid_df.columns
if col not in ['geometry', 'Longitude', 'Latitude']
]

im_data = pd.DataFrame(
{
col: [grid_df.iloc[nbr_index][col]] * event_count
for col in im_columns
}
)

im_data.to_csv(csv_path, index=False)
im_data.to_csv(csv_path, index=False)
# save the collection file name and the IM row id
event_list.append(csv_filename + f'x{event_j}')

Expand Down
3 changes: 2 additions & 1 deletion modules/performSIMULATION/capacitySpectrum/CapacityModels.py
Original file line number Diff line number Diff line change
Expand Up @@ -88,6 +88,7 @@ def convert_story_rise(structureType, stories): # noqa: N803
rise = None

else:
rise = None
# First, check if we have valid story information
try:
stories = int(stories)
Expand Down Expand Up @@ -340,7 +341,7 @@ def __init__(self, general_info, dD=0.001): # noqa: N803
self.Dy = self.capacity_data[self.design_level][self.HAZUS_type]['Dy']
self.Ay = self.capacity_data[self.design_level][self.HAZUS_type]['Ay']
except KeyError:
msg = f'No capacity data for {self.HAZUS_type} and {self.design_level}'
msg = f'No capacity data for build class {self.HAZUS_type} and design level {self.design_level}'
raise KeyError(msg) # noqa: B904
self.cao_peterson_2006 = cao_peterson_2006(
self.Dy, self.Ay, self.Du, self.Au, dD
Expand Down
2 changes: 2 additions & 0 deletions modules/performSIMULATION/capacitySpectrum/DampingModels.py
Original file line number Diff line number Diff line change
Expand Up @@ -256,6 +256,8 @@ def get_beta(self, Dp, Ap): # noqa: N803
f'The base model {self.base_model} does not have a useful'
'get_kappa method.'
)
if Dp <= 0 or Ap <= 0:
return beta_elastic
Du = self.capacity.Du # noqa: N806
Ax = self.capacity.Ax # noqa: N806
B = self.capacity.B # noqa: N806
Expand Down
6 changes: 2 additions & 4 deletions modules/performSIMULATION/capacitySpectrum/runCMS.py
Original file line number Diff line number Diff line change
Expand Up @@ -97,10 +97,8 @@ def find_performance_point(cap_x, cap_y, dem_x, dem_y, dd=0.001):
elif dem_y_interp[0] < cap_y_interp[0]:
perf_x = 0.001 # x_interp[0]
perf_y = 0.001 # cap_y_interp[0]
# except IndexError as err:
# print('No performance point found; curves do not intersect.')
# print('IndexError: ')
# print(err)
else:
print('No performance point found; curves do not intersect.')

return perf_x, perf_y

Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -868,6 +868,8 @@ def run_residual_demand( # noqa: C901
edges_gdf['capacity'] = edges_gdf['lanes'] * 1800
edges_gdf['normal_capacity'] = edges_gdf['capacity']
edges_gdf['normal_maxspeed'] = edges_gdf['maxspeed']
edges_gdf['start_nid'] = edges_gdf['start_nid'].astype(int)
edges_gdf['end_nid'] = edges_gdf['end_nid'].astype(int)
# edges_gdf['fft'] = edges_gdf['length']/edges_gdf['maxspeed'] * 2.23694
edges_gdf.to_csv('edges.csv', index=False)
nodes_gdf = gpd.read_file(node_geojson)
Expand Down
6 changes: 5 additions & 1 deletion modules/tools/BRAILS/runBrailsTransp.py
Original file line number Diff line number Diff line change
Expand Up @@ -72,6 +72,7 @@ def runBrails( # noqa: N802, D103
minimumHAZUS, # noqa: N803
maxRoadLength, # noqa: N803
lengthUnit, # noqa: N803
saveTrafficSimulationAttr, # noqa: N803
):
# Initialize TranspInventoryGenerator:
invGenerator = TranspInventoryGenerator( # noqa: N806
Expand All @@ -83,7 +84,8 @@ def runBrails( # noqa: N802, D103

# Combine and format the generated inventory to SimCenter transportation network inventory json format
invGenerator.combineAndFormat_HWY(
minimumHAZUS=minimumHAZUS, maxRoadLength=maxRoadLength, lengthUnit=lengthUnit
minimumHAZUS=minimumHAZUS, maxRoadLength=maxRoadLength, lengthUnit=lengthUnit,
connectivity=saveTrafficSimulationAttr
)


Expand All @@ -99,6 +101,7 @@ def main(args): # noqa: D103
)
parser.add_argument('--maxRoadLength', default=100, type=float)
parser.add_argument('--lengthUnit', default='m', type=str)
parser.add_argument('--saveTrafficSimulationAttr', default=False, type=str2bool, nargs='?', const=True)

args = parser.parse_args(args)

Expand All @@ -115,6 +118,7 @@ def main(args): # noqa: D103
args.minimumHAZUS,
args.maxRoadLength,
args.lengthUnit,
args.saveTrafficSimulationAttr,
)

log_msg('BRAILS successfully generated the requested transportation inventory')
Expand Down
Loading