Commit 8a618494 authored by Marius Isken's avatar Marius Isken
Browse files

StaticTarget levels // wip

parent be74892d
......@@ -344,6 +344,10 @@ def command_go(args):
parser.add_option(
'--parallel', dest='nparallel', type='int', default=1,
help='set number of events to process in parallel')
parser.add_option(
'--baraddur', dest='baraddur', default=False,
action='store_true',
help='start the Barad-dur plotting server')
parser, options, args = cl_parse('go', args, setup)
if len(args) < 2:
......@@ -358,6 +362,11 @@ def command_go(args):
else:
status = tuple(options.status.split(','))
if options.baraddur:
from grond.baraddur import BaraddurProcess
baraddur = BaraddurProcess(project_dir=op.abspath(op.curdir))
baraddur.start()
grond.go(
config,
event_names=event_names,
......@@ -365,6 +374,11 @@ def command_go(args):
status=status,
nparallel=options.nparallel)
if options.baraddur:
logger.info('Grond finished processing %d events. '
'Ctrl+C to kill Barad-dur')
baraddur.join()
def command_forward(args):
def setup(parser):
......
......@@ -19,5 +19,5 @@ class BaraddurProcess(_mp.Process):
if __name__ == '__main__':
p = BaraddurProcess('/home/marius/Development/testing/grond/rundir')
p = BaraddurProcess(project_dir='/home/marius/Development/testing/grond')
p.start()
......@@ -40,24 +40,42 @@ def makeColorGradient(misfits, fr=1., fg=.5, fb=1.,
class BaraddurRequestHandler(RequestHandler):
def initialize(self, config):
self.config = config
self.model = config.model
class BaraddurBokehHandler(BokehHandler):
def __init__(self, config, *args, **kwargs):
BokehHandler.__init__(self, *args, **kwargs)
self.config = config
self.model = config.model
self.nmodels = 0
self.source = ColumnDataSource()
def update_source(self, func):
new_nmodels, new_data = func(
skip_nmodels=self.nmodels,
keys=self.source.data.keys())
self.source.stream(new_data)
self.nmodels += new_nmodels
@gen.coroutine
def update_models(self):
return self.update_source(self.model.get_models)
@gen.coroutine
def update_misfits(self):
return self.update_source(self.model.get_misfits)
class BaraddurModel(object):
def __init__(self, rundir):
self.set_rundir(rundir)
def set_rundir(self, rundir):
logger.debug('Loading problem from %s' % rundir)
self.rundir = op.abspath(rundir)
self.problem = grond.core.load_problem_info(self.rundir)
self.parameters = self.problem.parameters
self.nparameters = self.problem.nparameters
self.targets = self.problem.targets
self.ntargets = self.problem.ntargets
def _jp(self, file):
......@@ -79,8 +97,8 @@ class BaraddurModel(object):
@property
def best_misfit(self):
_, data = self.get_misfits(keys='target_mean')
return data['target_mean'].min()
_, data = self.get_misfits(keys='misfit_mean')
return data['misfit_mean'].min()
@property
def niterations(self):
......@@ -151,9 +169,9 @@ class BaraddurModel(object):
misfits = combi.reshape((nmodels, self.ntargets, 2))
mf_dict = {}
for it in xrange(self.ntargets):
mf_dict['target_%03d' % it] = misfits[:, it, 0]
mf_dict['target_mean'] = num.mean(misfits, axis=1)[:, 0]
for it, target in enumerate(self.targets):
mf_dict[target.id] = misfits[:, it, 0]
mf_dict['misfit_mean'] = num.mean(misfits, axis=1)[:, 0]
mf_dict['niter'] = num.arange(nmodels, dtype=num.int)\
+ skip_nmodels + 1
......@@ -190,7 +208,7 @@ class Summary(BaraddurRequestHandler):
def get(self):
self.render('summary.html',
pages=pages,
problem=self.config.problem)
problem=self.model.problem)
class Misfit(BaraddurRequestHandler):
......@@ -198,10 +216,8 @@ class Misfit(BaraddurRequestHandler):
class MisfitsPlot(BaraddurBokehHandler):
def modify_document(self, doc):
self.nmodels = 0
self.source = ColumnDataSource()
self.source.add(num.ndarray(0, dtype=num.float32),
'target_mean')
'misfit_mean')
self.source.add(num.ndarray(0, dtype=num.float32),
'niter')
self.update_misfits()
......@@ -211,21 +227,13 @@ class Misfit(BaraddurRequestHandler):
plot_height=300,
x_axis_label='Iteration #',
y_axis_label='Misfit')
plot.scatter('niter', 'target_mean',
plot.scatter('niter', 'misfit_mean',
source=self.source, alpha=.4)
doc.add_root(plot)
doc.add_periodic_callback(self.update_misfits,
self.config.plot_update_small)
@gen.coroutine
def update_misfits(self):
new_nmodels, new_data = self.config.model.get_misfits(
skip_nmodels=self.nmodels,
keys=self.source.data.keys())
self.source.stream(new_data)
self.nmodels += new_nmodels
bokeh_handlers = {'misfit': MisfitsPlot}
@gen.coroutine
......@@ -236,7 +244,7 @@ class Misfit(BaraddurRequestHandler):
None,
app_path='/misfit',
url='plots'),
problem=self.config.problem)
problem=self.model.problem)
class Parameters(BaraddurRequestHandler):
......@@ -246,14 +254,14 @@ class Parameters(BaraddurRequestHandler):
ncols = 4
def modify_document(self, doc):
self.nmodels = 0
problem = self.config.problem
problem = self.model.problem
self.source = ColumnDataSource()
for p in ['niter'] + [p.name for p in problem.parameters]:
self.source.add(num.ndarray(0, dtype=num.float32), p)
for param in ['niter'] + [p.name for p in problem.parameters]:
self.source.add(
num.ndarray(0, dtype=num.float32),
param)
self.model = BaraddurModel(self.config.rundir)
self.update_parameters()
self.update_models()
plots = []
for par in problem.parameters:
......@@ -264,7 +272,8 @@ class Parameters(BaraddurRequestHandler):
fig.scatter('niter', par.name,
source=self.source, alpha=.4)
plots.append(fig)
plots += [None] * (self.ncols - (len(plots) % self.ncols))
if (len(plots) % self.ncols) != 0:
plots += [None] * (self.ncols - (len(plots) % self.ncols))
grid = layouts.gridplot(
plots,
......@@ -272,17 +281,9 @@ class Parameters(BaraddurRequestHandler):
ncols=self.ncols)
doc.add_root(grid)
doc.add_periodic_callback(self.update_parameters,
doc.add_periodic_callback(self.update_models,
self.config.plot_update_small)
@gen.coroutine
def update_parameters(self):
new_nmodels, new_data = self.config.model.get_models(
skip_nmodels=self.nmodels,
keys=self.source.data.keys())
self.source.stream(new_data)
self.nmodels += new_nmodels
bokeh_handlers = {'parameters': ParameterPlots}
@gen.coroutine
......@@ -294,7 +295,7 @@ class Parameters(BaraddurRequestHandler):
None,
app_path='/parameters',
url='plots'),
problem=self.config.problem)
problem=self.model.problem)
class Targets(BaraddurRequestHandler):
......@@ -302,45 +303,52 @@ class Targets(BaraddurRequestHandler):
class TargetContributionPlot(BaraddurBokehHandler):
def modify_document(self, doc):
self.nmodels = 0
self.source = ColumnDataSource()
self.update_contributions()
plot = figure(webgl=True,
responsive=True,
plot_height=300,
x_axis_label='Iteration #',
y_axis_label='Misfit')
doc.add_root(plot)
doc.add_periodic_callback(self.update_contributions, 1e3)
target_ids = [t.id for t in self.model.targets]
for param in ['niter'] + target_ids:
self.source.add(
num.ndarray(0, dtype=num.float32),
param)
self.update_misfits()
colors = ['blue', 'red']
for it, target in enumerate(self.model.targets):
plot.scatter('niter', target.id,
source=self.source,
color=colors[it])
@gen.coroutine
def update_contributions(self):
mx, misfits = grond.core.load_problem_data(
self.config.rundir, self.config.problem,
skip_models=self.nmodels)
# plot.patches(
# ['niter'] * self.model.ntargets,
# target_ids,
# source=self.source)
new_nmodels = mx.shape[0]
self.nmodels += new_nmodels
doc.add_root(plot)
doc.add_periodic_callback(self.update_misfits,
self.config.plot_update_small)
bokeh_handlers = {'contributions': TargetContributionPlot}
bokeh_handlers = {'targets': TargetContributionPlot}
@gen.coroutine
def get(self):
self.render('targets.html',
contribution_plot=autoload_server(
target_misfit_plot=autoload_server(
None,
app_path='/contributions',
app_path='/targets',
url='plots'),
pages=pages,
problem=self.config.problem)
problem=self.model.problem)
pages = OrderedDict([
('Rundirs', Rundirs),
('Summary', Misfit),
('Summary', Summary),
('Misfit', Misfit),
('Parameters', Parameters),
# ('Targets', Targets),
('Targets', Targets),
])
......@@ -375,10 +383,6 @@ class BaraddurConfig(Object):
self._rundir = None
self._model = None
@property
def problem(self):
return grond.core.load_problem_info(self.rundir)
@property
def rundir(self):
return self._rundir
......
......@@ -40,7 +40,10 @@
{% end %}
{% if len(models) == 0 %}
<tr>
<td colspan=4 style="text-align: center;"><i>No rundirs found</i></td>
<td colspan="6" style="text-align: center;">
<p style="font-style: italic;">No rundirs found</p>
<p style="font-style: italic; font-size: 90%; color: #999;">Refresh this page to reload available rundirs.</p>
</td>
</tr>
{% end %}
</tbody>
......
......@@ -15,7 +15,7 @@
<h2 class="content-subhead">Contributors</h3>
<p>Contributions of each target to the solution</p>
<div class="pure-u-1">
{% raw contribution_plot %}
{% raw target_misfit_plot %}
</div>
</div>
<div class="pure-u-1 pure-u-md-1 pure-u-lg-1-1">
......
......@@ -117,7 +117,7 @@ class Config(HasPaths):
synt = ds.synthetic_test
if synt:
synt.set_problem(problem)
problem.base_source = problem.unpack(synt.get_x())
problem.base_source = problem.get_source(synt.get_x())
def get_problem(self, event):
targets = self.get_targets(event)
......@@ -196,13 +196,13 @@ def get_best_x_and_gm(problem, xs, misfits):
def get_mean_source(problem, xs):
x_mean = get_mean_x(xs)
source = problem.unpack(x_mean)
source = problem.get_source(x_mean)
return source
def get_best_source(problem, xs, misfits):
x_best = get_best_x(problem, xs, misfits)
source = problem.unpack(x_best)
source = problem.get_source(x_best)
return source
......@@ -517,6 +517,11 @@ def solve(problem,
gxs = xhist[chains_i[0, :nlinks], :]
gms = chains_m[0, :nlinks]
col_width = 15
col_param_width = max([len(p) for p in problem.parameter_names])+2
console_output = '{:<{col_param_width}s}'
console_output += ''.join(['{:>{col_width}{type}}'] * 5)
if nlinks > (nlinks_cap-1)/2:
# mean and std of all bootstrap ensembles together
mbx = num.mean(bxs, axis=0)
......@@ -541,26 +546,41 @@ def solve(problem,
if 'state' in status:
lines.append(
'%-15s %15s %15s %15s %15s %15s' %
('parameter', 'B mean', 'B std', 'G mean', 'G std',
'G best'))
console_output.format(
'parameter', 'B mean', 'B std', 'G mean', 'G std',
'G best',
col_param_width=col_param_width,
col_width=col_width,
type='s'))
for (pname, mbv, sbv, mgv, sgv, bgv) in zip(
pnames, mbx, sbx, mgx, sgx, bgx):
lines.append(
'%-15s %15.4g %15.4g %15.4g %15.4g %15.4g' %
(pname, mbv, sbv, mgv, sgv, bgv))
console_output.format(
pname, mbv, sbv, mgv, sgv, bgv,
col_param_width=col_param_width,
col_width=col_width,
type='.4g'))
lines.append('%-15s %15s %15s %15.4g %15.4g %15.4g' % (
'misfit', '', '',
num.mean(gms), num.std(gms), num.min(gms)))
lines.append(
console_output.format(
'misfit', '', '',
'%.4g' % num.mean(gms),
'%.4g' % num.std(gms),
'%.4g' % num.min(gms),
col_param_width=col_param_width,
col_width=col_width,
type='s'))
if 'state' in status:
lines.append(
'%-15s %15i %-15s %15i %15i' % (
console_output.format(
'iteration', iiter+1, '(%s, %g)' % (phase, factor),
ntries_sample, ntries_preconstrain))
ntries_sample, ntries_preconstrain, '',
col_param_width=col_param_width,
col_width=col_width,
type=''))
if 'matrix' in status:
matrix = (chains_i[:, :30] % 94 + 32).T
......@@ -637,7 +657,7 @@ def forward(rundir_or_config_path, event_names):
ds.empty_cache()
ms, ns, results = problem.evaluate(x, result_mode='full')
event = problem.unpack(x).pyrocko_event()
event = problem.get_source(x).pyrocko_event()
events.append(event)
for result in results:
......@@ -766,7 +786,7 @@ def check(
else:
for i in xrange(n_random_synthetics):
x = problem.random_uniform(xbounds)
sources.append(problem.unpack(x))
sources.append(problem.get_source(x))
ms, ns, results = problem.evaluate(x, result_mode='full')
results_list.append(results)
......@@ -970,10 +990,6 @@ def go(config, event_names=None, force=False, nparallel=1, status=('state',)):
nevents = len(event_names)
from .baraddur import BaraddurProcess
baraddur = BaraddurProcess(project_dir=op.abspath(op.curdir))
baraddur.start()
for x in parimap.parimap(
process_event,
xrange(nevents),
......@@ -981,9 +997,6 @@ def go(config, event_names=None, force=False, nparallel=1, status=('state',)):
nprocs=nparallel):
pass
logger.info('Grond finished processing %d events. '
'Ctrl+C to kill Barad-dur')
baraddur.join()
def process_event(ievent, g_data_id):
......@@ -1007,7 +1020,7 @@ def process_event(ievent, g_data_id):
synt = ds.synthetic_test
if synt:
problem.base_source = problem.unpack(synt.get_x())
problem.base_source = problem.get_source(synt.get_x())
check_problem(problem)
......@@ -1161,11 +1174,11 @@ def export(what, rundirs, type=None, pnames=None, filename=None):
'%16.7g' % gm
elif type == 'source':
source = problem.unpack(x)
source = problem.get_source(x)
guts.dump(source, stream=out)
elif type == 'event':
ev = problem.unpack(x).pyrocko_event()
ev = problem.get_source(x).pyrocko_event()
model.dump_events([ev], stream=out)
else:
......
......@@ -9,6 +9,8 @@ from pyrocko import gf, util, guts, moment_tensor as mtm
from pyrocko.guts import (Object, String, Bool, List, Float, Dict, Int,
StringChoice)
from .targets import MisfitTarget, MisfitSatelliteTarget
guts_prefix = 'grond'
logger = logging.getLogger('grond')
......@@ -28,6 +30,7 @@ class Parameter(Object):
kwargs['name'] = args[0]
if len(args) >= 2:
kwargs['unit'] = args[1]
self.target = kwargs.pop('target', None)
Object.__init__(self, **kwargs)
......@@ -93,12 +96,14 @@ class Problem(Object):
targets = List.T()
def __init__(self, **kwargs):
Object.__init__(self, **kwargs)
self.init_satellite_target_leveling()
self._bootstrap_weights = None
self._target_weights = None
self._engine = None
self._group_mask = None
logger.name = self.__class__.__name__
def get_engine(self):
......@@ -110,8 +115,12 @@ class Problem(Object):
o._target_weights = None
return o
def parameter_dict(self, x):
return ADict((p.name, v) for (p, v) in zip(self.parameters, x))
def parameter_dict(self, x, target=None):
params = []
for ip, p in enumerate(self.parameters):
if p.target is target:
params.append((p.name, x[ip]))
return ADict(params)
def parameter_array(self, d):
return num.array([d[p.name] for p in self.parameters], dtype=num.float)
......@@ -161,9 +170,64 @@ class Problem(Object):
def combined(self):
return self.parameters + self.dependants
@property
def satellite_targets(self):
return [t for t in self.targets
if isinstance(t, MisfitSatelliteTarget)]
@property
def waveform_targets(self):
return [t for t in self.targets
if isinstance(t, MisfitTarget)]
@property
def has_statics(self):
if self.satellite_targets:
return True
return False
@property
def has_waveforms(self):
if self.waveform_targets:
return True
return False
def set_engine(self, engine):
self._engine = engine
def init_satellite_target_leveling(self):
for t in self.satellite_targets:
add = True
add_parameters = [
Parameter('%s:waterlevel' % t.scene_id.lower(), 'm',
label='Waterlevel (%s)' % t.scene_id,
target=t),
Parameter('%s:ramp_north' % t.scene_id.lower(), 'm/m',
label='Ramp North (%s)' % t.scene_id,
target=t),
Parameter('%s:ramp_east' % t.scene_id.lower(), 'm/m',
label='Ramp East (%s)' % t.scene_id,
target=t),
]
add_ranges = t.inner_misfit_config.leveling_ranges.copy()
for k in add_ranges.keys():
new_k = '%s:%s' % (t.scene_id.lower(), k)
if new_k in self.parameter_names:
add = False
add_ranges[new_k] = add_ranges.pop(k)
if add:
logger.info('Adding waterlevel parameters for %s' % t.scene_id)
self.parameters += add_parameters
self.ranges.update(add_ranges)
def set_satellite_scene_levels(self, x):
for target in self.satellite_targets:
levels = self.parameter_dict(x, target=target)
for k in levels.keys():
new_k = k.split(':')[-1]
levels[new_k] = levels.pop(k)
target.set_scene_levels(**levels)
def make_bootstrap_weights(self, nbootstrap):
ntargets = len(self.targets)
ws = num.zeros((nbootstrap, ntargets))
......@@ -283,6 +347,9 @@ class Problem(Object):
return self._group_mask
def evaluate(self, x, dump_data=False):
raise NotImplementedError()
class CMTProblem(Problem):
......@@ -308,7 +375,7 @@ class CMTProblem(Problem):
nbootstrap = Int.T(default=10)
mt_type = StringChoice.T(default='full', choices=['full', 'deviatoric'])
def unpack(self, x):
def get_source(self, x):
d = self.parameter_dict(x)
rm6 = num.array([d.rmnn, d.rmee, d.rmdd, d.rmne, d.rmnd, d.rmed],
dtype=num.float)
......@@ -334,7 +401,7 @@ class CMTProblem(Problem):
if pname in ('rel_moment_iso', 'rel_moment_clvd'):
y = num.zeros(xs.shape[0])
for i, x in enumerate(xs):
source = self.unpack(x)
source = self.get_source(x)
mt = source.pyrocko_moment_tensor()
res = mt.standard_decomposition()
......@@ -405,7 +472,7 @@ class CMTProblem(Problem):
d.rmnn, d.rmee, d.rmdd, d.rmne, d.rmnd, d.rmed = m6
x = self.parameter_array(d)
source = self.unpack(x)
source = self.get_source(x)
if any(self.distance_min > source.distance_to(t)
for t in self.targets):
raise Forbidden()
......@@ -420,7 +487,7 @@ class CMTProblem(Problem):