Skip to content

Commit

Permalink
Simplify the StationParameters class
Browse files Browse the repository at this point in the history
Internal dictionaries `params_dict`, `params_err_dict`, and
`is_outlier_dict` are removed.

Instead, the class now has a `get_spectral_parameters` method that
returns a dictionary of spectral parameters. This method is used in
the `ssp_inversion` and `ssp_sqlite_output` modules.
  • Loading branch information
claudiodsf committed Mar 29, 2024
1 parent 435e9ea commit aa87add
Show file tree
Hide file tree
Showing 6 changed files with 81 additions and 82 deletions.
98 changes: 42 additions & 56 deletions sourcespec/ssp_data_types.py
Original file line number Diff line number Diff line change
Expand Up @@ -190,13 +190,13 @@ def __init__(self, param_id, name=None, units=None, value=None,
self.confidence_level = confidence_level
self.outlier = False

def value_uncertainty(self):
"""Return value and uncertainty as 3-element tuple."""
def compact_uncertainty(self):
"""Return uncertainty in a compact form."""
if self.lower_uncertainty is not None:
uncertainty = (self.lower_uncertainty, self.upper_uncertainty)
else:
uncertainty = (self.uncertainty, self.uncertainty)
return (self.value, *uncertainty)
return (self.lower_uncertainty, self.upper_uncertainty)
if self.uncertainty is not None:
return (self.uncertainty, self.uncertainty)
return (np.nan, np.nan)


class StationParameters(OrderedAttribDict):
Expand All @@ -208,55 +208,34 @@ class StationParameters(OrderedAttribDict):
objects.
"""

def __init__(self, param_id, instrument_type=None,
def __init__(self, station_id, instrument_type=None,
latitude=None, longitude=None,
hypo_dist_in_km=None, epi_dist_in_km=None, azimuth=None):
self.param_id = param_id
self.station_id = station_id
self.instrument_type = instrument_type
self.latitude = latitude
self.longitude = longitude
self.hypo_dist_in_km = hypo_dist_in_km
self.epi_dist_in_km = epi_dist_in_km
self.azimuth = azimuth
self.params_dict = {}
self.params_err_dict = {}
self.is_outlier_dict = {}

def __setattr__(self, attr, value):
if isinstance(value, SpectralParameter):
parname = attr
par = value
self.params_dict[parname] = par.value
if par.uncertainty is not None:
self.params_err_dict[parname] = (
par.uncertainty, par.uncertainty
)
else:
self.params_err_dict[parname] = (
par.lower_uncertainty, par.upper_uncertainty
)
self.is_outlier_dict[parname] = par.outlier
self[attr] = value

def rebuild_dictionaries(self):
"""Rebuild spectral parameters dictionaries."""
for key, value in self.items():
if not isinstance(value, SpectralParameter):
continue
parname = key
par = value
self.params_dict[parname] = par.value
if par.uncertainty is not None:
self.params_err_dict[parname] = (
par.uncertainty, par.uncertainty
)
elif par.lower_uncertainty is not None:
self.params_err_dict[parname] = (
par.lower_uncertainty, par.upper_uncertainty
)
else:
self.params_err_dict[parname] = (np.nan, np.nan)
self.is_outlier_dict[parname] = par.outlier
self.Mw = None
self.fc = None
self.t_star = None
self.Mo = None
self.radius = None
self.ssd = None
self.Qo = None
self.Er = None
self.sigma_a = None
self.Ml = None

def get_spectral_parameters(self):
"""Return a dictionary of spectral parameters."""
return {
key: value
for key, value in self.items()
if isinstance(value, SpectralParameter)
}


class SummaryStatistics(OrderedAttribDict):
Expand Down Expand Up @@ -294,7 +273,9 @@ def compact_uncertainty(self):
"""Return uncertainty in a compact form."""
if self.lower_uncertainty is not None:
return (self.lower_uncertainty, self.upper_uncertainty)
return (self.uncertainty, self.uncertainty)
if self.uncertainty is not None:
return (self.uncertainty, self.uncertainty)
return (np.nan, np.nan)


class SummarySpectralParameter(OrderedAttribDict):
Expand Down Expand Up @@ -342,19 +323,23 @@ def __init__(self):
def value_array(self, key, filter_outliers=False):
"""Return an array of values for the given key."""
vals = np.array([
x.params_dict.get(key, np.nan)
for x in self.station_parameters.values()
stat_par[key].value
if key in stat_par and stat_par[key] is not None
else np.nan
for stat_par in self.station_parameters.values()
])
if filter_outliers:
outliers = self.outlier_array(key)
vals = vals[~outliers]
return vals

def error_array(self, key, filter_outliers=False):
"""Return an array of errors for the given key."""
"""Return an array of errors (two columns) for the given key."""
errs = np.array([
x.params_err_dict.get(key, np.nan)
for x in self.station_parameters.values()
stat_par[key].compact_uncertainty()
if key in stat_par and stat_par[key] is not None
else (np.nan, np.nan)
for stat_par in self.station_parameters.values()
])
if filter_outliers:
outliers = self.outlier_array(key)
Expand All @@ -365,9 +350,11 @@ def outlier_array(self, key):
"""Return an array of outliers for the given key."""
return np.array(
[
stat_par[key].outlier
if key in stat_par
# if we cannot find the given key, we assume outlier=True
x.is_outlier_dict.get(key, True)
for x in self.station_parameters.values()
else True
for stat_par in self.station_parameters.values()
]
)

Expand Down Expand Up @@ -408,7 +395,6 @@ def find_outliers(self, key, n):
for stat_id, outl in zip(station_ids, outliers):
stat_par = self.station_parameters[stat_id]
stat_par[key].outlier = outl
stat_par.rebuild_dictionaries()

def mean_values(self):
"""Return a dictionary of mean values."""
Expand Down
7 changes: 4 additions & 3 deletions sourcespec/ssp_html_report.py
Original file line number Diff line number Diff line change
Expand Up @@ -223,9 +223,10 @@ def _summary_value_and_err_text(value, error, fmt):

def _station_value_and_err_text(par, key, fmt):
"""Format station value and error text."""
try:
_par = par[key]
except KeyError:
# par[key] can be None even if key is in par (e.g., if key is 'Ml' and
# local magnitude is not computed)
_par = par[key] if key in par else None
if _par is None:
return '', ''
if isinstance(_par, SpectralParameter):
value = _par.value
Expand Down
16 changes: 11 additions & 5 deletions sourcespec/ssp_inversion.py
Original file line number Diff line number Diff line change
Expand Up @@ -295,7 +295,7 @@ def _spec_inversion(config, spec, spec_weight):
)

station_pars = StationParameters(
param_id=spec.id, instrument_type=spec.stats.instrtype,
station_id=spec.id, instrument_type=spec.stats.instrtype,
latitude=stla, longitude=stlo,
hypo_dist_in_km=spec.stats.hypo_dist,
epi_dist_in_km=spec.stats.epi_dist,
Expand All @@ -319,7 +319,7 @@ def _spec_inversion(config, spec, spec_weight):
travel_time = spec.stats.travel_times[config.wave_type[0]]
# seismic moment
station_pars.Mo = SpectralParameter(
param_id='Mw', value=mag_to_moment(Mw), format_spec='{:.3e}')
param_id='Mo', value=mag_to_moment(Mw), format_spec='{:.3e}')
# source radius in meters
station_pars.radius = SpectralParameter(
param_id='radius', value=source_radius(fc, beta, k_coeff),
Expand Down Expand Up @@ -389,8 +389,14 @@ def _spec_inversion(config, spec, spec_weight):

def _synth_spec(config, spec, station_pars):
"""Return a stream with one or more synthetic spectra."""
par = station_pars.params_dict
par_err = station_pars.params_err_dict
par = {
x.param_id: x.value
for x in station_pars.get_spectral_parameters().values()
}
par_err = {
x.param_id: x.compact_uncertainty()
for x in station_pars.get_spectral_parameters().values()
}
spec_st = Stream()
params_opt = [par[key] for key in ('Mw', 'fc', 't_star')]

Expand Down Expand Up @@ -505,7 +511,7 @@ def spectral_inversion(config, spec_st, weight_st):
logger.warning(msg)
continue
spec_st += _synth_spec(config, spec, station_pars)
sspec_output.station_parameters[station_pars.param_id] = station_pars
sspec_output.station_parameters[station_pars.station_id] = station_pars

logger.info('Inverting spectra: done')
logger.info('---------------------------------------------------')
Expand Down
5 changes: 2 additions & 3 deletions sourcespec/ssp_local_magnitude.py
Original file line number Diff line number Diff line change
Expand Up @@ -171,9 +171,8 @@ def local_magnitude(config, st, proc_st, sspec_output):

# Make sure that the param_Ml object is always defined, even when Ml
# is not computed (i.e., "continue" below)
try:
param_Ml = station_pars.Ml
except KeyError:
param_Ml = station_pars.Ml
if param_Ml is None:
param_Ml = SpectralParameter(
param_id='Ml', value=np.nan, format_spec='{:.2f}')
station_pars.Ml = param_Ml
Expand Down
10 changes: 4 additions & 6 deletions sourcespec/ssp_radiated_energy.py
Original file line number Diff line number Diff line change
Expand Up @@ -193,15 +193,13 @@ def radiated_energy_and_apparent_stress(
# Make sure that the param_Er and param_sigma_a objects are always
# defined, even when Er and/or sigma_a are not computed
# (i.e., "continue" below)
try:
param_Er = station_pars.Er
except KeyError:
param_Er = station_pars.Er
if param_Er is None:
param_Er = SpectralParameter(
param_id='Er', value=np.nan, format_spec='{:.3e}')
station_pars.Er = param_Er
try:
param_sigma_a = station_pars.sigma_a
except KeyError:
sigma_a = station_pars.sigma_a
if sigma_a is None:
param_sigma_a = SpectralParameter(
param_id='sigma_a', value=np.nan, format_spec='{:.3e}')
station_pars.sigma_a = param_sigma_a
Expand Down
27 changes: 18 additions & 9 deletions sourcespec/ssp_sqlite_output.py
Original file line number Diff line number Diff line change
Expand Up @@ -155,23 +155,32 @@ def _write_stations_table(cursor, db_file, sspec_output, config):
# Insert new line
t = (
statId, evid, runid,
*par.Mo.value_uncertainty(),
par.Mo.value,
*par.Mo.compact_uncertainty(),
int(par.Mo.outlier),
*par.Mw.value_uncertainty(),
par.Mw.value,
*par.Mw.compact_uncertainty(),
int(par.Mw.outlier),
*par.fc.value_uncertainty(),
par.fc.value,
*par.fc.compact_uncertainty(),
int(par.fc.outlier),
*par.t_star.value_uncertainty(),
par.t_star.value,
*par.t_star.compact_uncertainty(),
int(par.t_star.outlier),
*par.Qo.value_uncertainty(),
par.Qo.value,
*par.Qo.compact_uncertainty(),
int(par.Qo.outlier),
*par.ssd.value_uncertainty(),
par.ssd.value,
*par.ssd.compact_uncertainty(),
int(par.ssd.outlier),
*par.radius.value_uncertainty(),
par.radius.value,
*par.radius.compact_uncertainty(),
int(par.radius.outlier),
*par.Er.value_uncertainty(),
par.Er.value,
*par.Er.compact_uncertainty(),
int(par.Er.outlier),
*par.sigma_a.value_uncertainty(),
par.sigma_a.value,
*par.sigma_a.compact_uncertainty(),
int(par.sigma_a.outlier),
par.hypo_dist_in_km,
par.azimuth
Expand Down

0 comments on commit aa87add

Please sign in to comment.