Browse Source

pick residual fitting (Grond can now be used as a locator)

dev-scenario
Sebastian Heimann 2 years ago
committed by Sebastian Heimann
parent
commit
de65687338
  1. 1
      setup.py
  2. 9
      src/apps/grond.py
  3. 30
      src/scenario.py
  4. 1
      src/targets/__init__.py
  5. 1
      src/targets/phase_pick/__init__.py
  6. 3
      src/targets/phase_pick/plot.py
  7. 180
      src/targets/phase_pick/target.py
  8. 83
      test/test_locate.py

1
setup.py

@ -134,6 +134,7 @@ setup(
'grond.targets.waveform_oac',
'grond.targets.satellite',
'grond.targets.gnss_campaign',
'grond.targets.phase_pick',
'grond.problems',
'grond.problems.cmt',
'grond.problems.double_dc',

9
src/apps/grond.py

@ -418,6 +418,12 @@ def command_scenario(args):
'--gf-store-superdirs',
dest='gf_store_superdirs',
help='Comma-separated list of directories containing GF stores')
parser.add_option(
'--no-map',
dest='make_map',
default=True,
action='store_false',
help='suppress generation of map')
parser, options, args = cl_parse('scenario', args, setup)
@ -468,7 +474,8 @@ def command_scenario(args):
scenario.build(
force=options.force,
interactive=True,
gf_store_superdirs=gf_store_superdirs)
gf_store_superdirs=gf_store_superdirs,
make_map=options.make_map)
logger.info(CLIHints('scenario',
config=scenario.get_grond_config_path(),

30
src/scenario.py

@ -4,6 +4,7 @@ import shutil
import os.path as op
from pyrocko import gf, scenario, util
from pyrocko.gui import marker as pmarker
import math
import grond
@ -95,7 +96,12 @@ class GrondScenario(object):
for obs in self.observations],
source_generator=self.problem.get_scenario_source_generator())
def create_scenario(self, interactive=True, gf_store_superdirs=None):
def create_scenario(
self,
interactive=True,
gf_store_superdirs=None,
make_map=True):
logger.info('Creating scenario...')
scenario = self.get_scenario()
@ -128,11 +134,18 @@ class GrondScenario(object):
scenario.dump(filename=op.join(data_dir, 'scenario.yml'))
scenario.dump_data(path=data_dir)
scenario.make_map(op.join(self.project_dir, 'scenario_map.pdf'))
if make_map:
scenario.make_map(op.join(self.project_dir, 'scenario_map.pdf'))
shutil.move(op.join(data_dir, 'sources.yml'),
op.join(data_dir, 'scenario_sources.yml'))
markers = scenario.get_onsets()
marker_path = op.join(data_dir, 'picks', 'picks.markers')
if markers:
util.ensuredirs(marker_path)
pmarker.save_markers(markers, marker_path)
def get_grond_config(self):
engine_config = grond.EngineConfig(
gf_stores_from_pyrocko_config=False,
@ -163,13 +176,21 @@ class GrondScenario(object):
util.ensuredirs(config_path)
grond.write_config(self.get_grond_config(), config_path)
def build(self, force=False, interactive=False, gf_store_superdirs=None):
def build(
self,
force=False,
interactive=False,
gf_store_superdirs=None,
make_map=True):
logger.info('Building scenario...')
self.create_project_dir(force)
self.create_scenario(
interactive=interactive, gf_store_superdirs=gf_store_superdirs)
interactive=interactive,
gf_store_superdirs=gf_store_superdirs,
make_map=make_map)
self.create_grond_files()
@ -212,6 +233,7 @@ class WaveformObservation(Observation):
station_generator=scenario.targets.RandomStationGenerator(
nstations=self.nstations),
store_id=self.store_id,
tabulated_phases_from_store=True,
seismogram_quantity='displacement')
def get_grond_target_group(self):

1
src/targets/__init__.py

@ -4,3 +4,4 @@ from .waveform_phase_ratio import * # noqa
from .waveform_oac import * # noqa
from .satellite import * # noqa
from .gnss_campaign import * # noqa
from .phase_pick import * # noqa

1
src/targets/phase_pick/__init__.py

@ -0,0 +1 @@
from .target import * # noqa

3
src/targets/phase_pick/plot.py

@ -0,0 +1,3 @@
def get_plot_classes():
return []

180
src/targets/phase_pick/target.py

@ -0,0 +1,180 @@
import logging
import numpy as num
from pyrocko.guts import Float, Tuple, String
from pyrocko import gf
from ..base import (
MisfitTarget, TargetGroup, MisfitResult)
from grond.meta import has_get_plot_classes
guts_prefix = 'grond'
logger = logging.getLogger('grond.targets.phase_pick.target')
def log_exclude(target, reason):
logger.debug('Excluding potential target %s: %s' % (
target.string_id(), reason))
class PhasePickTargetGroup(TargetGroup):
'''
Generate targets to fit phase arrival times.
'''
distance_min = Float.T(optional=True)
distance_max = Float.T(optional=True)
distance_3d_min = Float.T(optional=True)
distance_3d_max = Float.T(optional=True)
depth_min = Float.T(optional=True)
depth_max = Float.T(optional=True)
store_id = gf.StringID.T(optional=True)
pick_synthetic_traveltime = gf.Timing.T(
help='Synthetic phase arrival definition.')
pick_phasename = String.T(
help='Name of phase in pick file.')
def get_targets(self, ds, event, default_path='none'):
logger.debug('Selecting phase pick targets...')
origin = event
targets = []
for st in ds.get_stations():
target = PhasePickTarget(
codes=st.nsl(),
lat=st.lat,
lon=st.lon,
north_shift=st.north_shift,
east_shift=st.east_shift,
depth=st.depth,
store_id=self.store_id,
manual_weight=self.weight,
normalisation_family=self.normalisation_family,
path=self.path or default_path,
pick_synthetic_traveltime=self.pick_synthetic_traveltime,
pick_phasename=self.pick_phasename)
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
target.set_dataset(ds)
targets.append(target)
return targets
class PhasePickResult(MisfitResult):
pass
@has_get_plot_classes
class PhasePickTarget(gf.Location, MisfitTarget):
'''
Target to fit phase arrival times.
'''
codes = Tuple.T(
3, String.T(),
help='network, station, location codes.')
store_id = gf.StringID.T(
help='ID of Green\'s function store (only used for earth model).')
pick_synthetic_traveltime = gf.Timing.T(
help='Synthetic phase arrival definition.')
pick_phasename = String.T(
help='Name of phase in pick file.')
can_bootstrap_weights = True
def __init__(self, **kwargs):
gf.Location.__init__(self, **kwargs)
MisfitTarget.__init__(self, **kwargs)
@classmethod
def get_plot_classes(cls):
from . import plot
plots = super(PhasePickTarget, cls).get_plot_classes()
plots.extend(plot.get_plot_classes())
return plots
def string_id(self):
return '.'.join(x for x in (self.path,) + self.codes)
def get_plain_targets(self, engine, source):
return self.prepare_modelling(engine, source, None)
def prepare_modelling(self, engine, source, targets):
return []
def get_times(self, engine, source):
tobs = None
tsyn = None
ds = self.get_dataset()
store = engine.get_store(self.store_id)
tsyn = source.time + store.t(
self.pick_synthetic_traveltime, source, self)
marker = ds.get_pick(
source.name,
self.codes[:3],
self.pick_phasename)
if marker:
tobs = marker.tmin
return tobs, tsyn
def finalize_modelling(
self, engine, source, modelling_targets, modelling_results):
ds = self.get_dataset() # noqa
tobs, tsyn = self.get_times(engine, source)
misfit = abs(tobs - tsyn)
norm = 1.0
result = PhasePickResult(
misfits=num.array([[misfit, norm]], dtype=num.float))
return result
__all__ = '''
PhasePickTargetGroup
PhasePickTarget
PhasePickResult
'''.split()

83
test/test_locate.py

@ -0,0 +1,83 @@
from __future__ import absolute_import
import os
import shutil
from pyrocko import guts
from grond import config
from grond import Environment
from . import common
from .common import grond, chdir
_multiprocess_can_split = True
def test_locate():
playground_dir = common.get_playground_dir()
common.get_test_data('gf_stores/crust2_ib/')
gf_stores_path = common.test_data_path('gf_stores')
with chdir(playground_dir):
scenario_dir = 'scenario_locate'
if os.path.exists(scenario_dir):
shutil.rmtree(scenario_dir)
grond('scenario', '--targets=waveforms', '--nevents=1',
'--nstations=6', '--gf-store-superdirs=%s' % gf_stores_path,
'--no-map',
scenario_dir)
with chdir(scenario_dir):
config_path = 'config/scenario.gronf'
quick_config_path = 'config/scenario_quick.gronf'
event_names = grond('events', config_path).strip().split('\n')
env = Environment([config_path] + event_names)
conf = env.get_config()
mod_conf = conf.clone()
target_groups = guts.load_all(string='''
--- !grond.PhasePickTargetGroup
normalisation_family: picks
path: pick.p
weight: 1.0
store_id: crust2_ib
distance_min: 10000.0
distance_max: 1000000.0
pick_synthetic_traveltime: '{stored:any_P}'
pick_phasename: 'any_P'
--- !grond.PhasePickTargetGroup
normalisation_family: picks
path: pick.s
weight: 1.0
store_id: crust2_ib
distance_min: 10000.0
distance_max: 1000000.0
pick_synthetic_traveltime: '{stored:any_S}'
pick_phasename: 'any_S'
''')
mod_conf.dataset_config.picks_paths = [
'data/scenario/picks/picks.markers']
mod_conf.target_groups = target_groups
# mod_conf.set_elements(
# 'analyser_configs[:].niterations', 100)
# mod_conf.set_elements(
# 'optimiser_config.sampler_phases[:].niterations', 100)
# mod_conf.set_elements(
# 'optimiser_config.nbootstrap', 5)
mod_conf.optimiser_config.sampler_phases[-1].niterations = 15000
mod_conf.set_basepath(conf.get_basepath())
config.write_config(mod_conf, quick_config_path)
grond('diff', config_path, quick_config_path)
grond('check', quick_config_path, *event_names)
grond('go', quick_config_path, *event_names)
rundir_paths = common.get_rundir_paths(config_path, event_names)
grond('report', *rundir_paths)
Loading…
Cancel
Save