MisfitTarget, TargetGroup, MisfitResult)
logger.debug('Excluding potential target %s: %s' % ( target.string_id(), reason))
''' Generate targets to compare ratios or log ratios of two seismic phases.
misfit = | a_obs / (a_obs + b_obs) - a_syn / (a_syn + b_syn) |
or
misfit = | log(a_obs / (a_obs + b_obs) + waterlevel) - log(a_syn / (a_syn + b_syn) + waterlevel) | '''
optional=True, help='select GF store based on event-station geometry.')
default=True, help='If true, compare synthetic and observed log ratios')
default=0.01, help='Waterlevel added to both ratios when comparing on logarithmic ' 'scale, to avoid log(0)')
logger.debug('Selecting phase ratio targets...') origin = event targets = []
for st in ds.get_stations(): blacklisted = False for measure in [self.measure_a, self.measure_b]: for cha in measure.channels: if ds.is_blacklisted((st.nsl() + (cha,))): blacklisted = True
if self.store_id_selector: store_id = self.store_id_selector.get_store_id( event, st, cha) else: store_id = self.store_id
target = PhaseRatioTarget( codes=st.nsl(), lat=st.lat, lon=st.lon, north_shift=st.north_shift, east_shift=st.east_shift, depth=st.depth, interpolation=self.interpolation, store_id=store_id, measure_a=self.measure_a, measure_b=self.measure_b, manual_weight=self.weight, normalisation_family=self.normalisation_family, path=self.path or default_path, backazimuth=0.0, fit_log_ratio=self.fit_log_ratio, fit_log_ratio_waterlevel=self.fit_log_ratio_waterlevel)
if blacklisted: log_exclude(target, 'blacklisted') continue
if self.distance_min is not None and \ target.distance_to(origin) < self.distance_min: log_exclude(target, 'distance < distance_min') continue
if self.distance_max is not None and \ target.distance_to(origin) > self.distance_max: log_exclude(target, 'distance > distance_max') continue
if self.distance_3d_min is not None and \ target.distance_3d_to(origin) < self.distance_3d_min: log_exclude(target, 'distance_3d < distance_3d_min') continue
if self.distance_3d_max is not None and \ target.distance_3d_to(origin) > self.distance_3d_max: log_exclude(target, 'distance_3d > distance_3d_max') continue
if self.depth_min is not None and \ target.depth < self.depth_min: log_exclude(target, 'depth < depth_min') continue
if self.depth_max is not None and \ target.depth > self.depth_max: log_exclude(target, 'depth > depth_max') continue
bazi, _ = target.azibazi_to(origin) target.backazimuth = bazi target.set_dataset(ds) targets.append(target)
return targets
''' Target to compare ratios or log ratios of two seismic phases.
misfit = | a_obs / (a_obs + b_obs) - a_syn / (a_syn + b_syn) |
'''
3, String.T(), help='network, station, location codes.')
help='ID of Green\'s function store to use for the computation.')
default=True, help='if true, compare synthetic and observed log ratios')
default=0.01, help='Waterlevel added to both ratios when comparing on logarithmic ' 'scale, to avoid log(0)')
gf.Location.__init__(self, **kwargs) MisfitTarget.__init__(self, **kwargs)
def get_plot_classes(cls):
return '.'.join(x for x in (self.path,) + self.codes)
return self.prepare_modelling(engine, source, None)
modelling_targets = [] for measure in [self.measure_a, self.measure_b]: modelling_targets.extend(measure.get_modelling_targets( self.codes, self.lat, self.lon, self.depth, self.store_id, self.backazimuth))
return modelling_targets
self, engine, source, modelling_targets, modelling_results):
ds = self.get_dataset()
try: imt = 0 amps = [] for measure in [self.measure_a, self.measure_b]: nmt_this = measure.get_nmodelling_targets() amp_obs, _ = measure.evaluate( engine, source, modelling_targets[imt:imt+nmt_this], dataset=ds)
amp_syn, _ = measure.evaluate( engine, source, modelling_targets[imt:imt+nmt_this], trs=[r.trace.pyrocko_trace() for r in modelling_results[imt:imt+nmt_this]])
amps.append((amp_obs, amp_syn))
imt += nmt_this
(a_obs, a_syn), (b_obs, b_syn) = amps
eps = self.fit_log_ratio_waterlevel if self.fit_log_ratio: res_a = num.log(a_obs / (a_obs + b_obs) + eps) \ - num.log(a_syn / (a_syn + b_syn) + eps) else: res_a = a_obs / (a_obs + b_obs) - a_syn / (a_syn + b_syn)
misfit = num.abs(res_a) norm = 1.0
result = PhaseRatioResult( misfits=num.array([[misfit, norm]], dtype=num.float), a_obs=a_obs, b_obs=b_obs, a_syn=a_syn, b_syn=b_syn)
return result
except dataset.NotFound as e: logger.debug(str(e)) return gf.SeismosizerError('No waveform data: %s' % str(e))
PhaseRatioTargetGroup PhaseRatioTarget PhaseRatioResult '''.split() |