Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
51 changes: 19 additions & 32 deletions apps/beat
Original file line number Diff line number Diff line change
Expand Up @@ -12,13 +12,14 @@ import shutil

from optparse import OptionParser

from beat import heart, config, utility, inputf, plotting
from beat import heart, config, utility, inputf
from beat.models import load_model, Stage, estimate_hypers, sample
from beat.backend import TextChain
from beat.backend import TextChain
from beat.sampler import SamplingHistory
from beat.sources import MTSourceWithMagnitude
from beat.utility import list2string
from numpy import savez, atleast_2d
from beat import backend

from pyrocko import model, util
from pyrocko.trace import snuffle
Expand Down Expand Up @@ -445,7 +446,7 @@ def command_import(args):
stage.load_results(
model=problem.model, stage_number=-1, load='full')

point = plotting.get_result_point(stage, problem.config, 'max')
point = backend.get_result_point(stage, problem.config, 'max')
n_sources = problem.config.problem_config.n_sources

source_params = problem.config.problem_config.priors.keys()
Expand Down Expand Up @@ -482,10 +483,10 @@ def command_update(args):

parser.add_option(
'--parameters',
default=['structure'], type='string',
default=[], type='string',
action='callback', callback=list_callback,
help='Parameters to update; "structure, hypers, hierarchicals". '
'Default: ["structure"] (config file-structure only)')
help='Parameters to update; "hypers, hierarchicals". '
'Default: [] (config file-structure only)')

parser.add_option(
'--mode', dest='mode',
Expand Down Expand Up @@ -1059,6 +1060,7 @@ def command_build_gfs(args):

def command_plot(args):

from beat import plotting
command_str = 'plot'

def setup(parser):
Expand Down Expand Up @@ -1151,7 +1153,6 @@ def command_plot(args):
dest='build',
action='store_true',
help='Build models during problem loading.')

plots_avail = plotting.available_plots()

details = '''Available <plot types> are: %s or "all". Multiple plots can be
Expand Down Expand Up @@ -1387,17 +1388,13 @@ def command_export(args):
model=problem.model, stage_number=options.stage_number,
load='trace', chains=[-1])

trace_name = 'chain--1.csv'
results_path = pjoin(problem.outfolder, config.results_dir_name)
logger.info('Saving results to %s' % results_path)
util.ensuredir(results_path)

results_trace = pjoin(stage.handler.stage_path(-1), trace_name)
shutil.copy(results_trace, pjoin(results_path, trace_name))

point = plotting.get_result_point(
point = backend.get_result_point(
stage, problem.config, point_llk=options.post_llk)

mpoint = dict()
for var in problem.varnames:
mpoint[var] = atleast_2d(stage.mtrace.get_values(var).mean())

for datatype, composite in problem.composites.items():
logger.info(
'Exporting "%s" synthetics for "%s" likelihood parameters:' % (
Expand All @@ -1407,6 +1404,9 @@ def command_export(args):
varname, list2string(value.ravel().tolist())))

results = composite.assemble_results(point)
results_path = pjoin(problem.outfolder, config.results_dir_name)
logger.info('Saving results to %s' % results_path)
util.ensuredir(results_path)

if datatype == 'seismic':
from pyrocko import io
Expand Down Expand Up @@ -1435,7 +1435,7 @@ def command_export(args):
if hasattr(sc.parameters, 'update_covariances'):
if sc.parameters.update_covariances:
logger.info('Saving velocity model covariance matrixes...')
composite.update_weights(point)
composite.update_weights(mpoint)
for wmap in composite.wavemaps:
pcovs = {
list2string(dataset.nslc_id):
Expand All @@ -1444,23 +1444,10 @@ def command_export(args):

outname = pjoin(
results_path, '%s_C_vm_%s' % (
datatype, wmap._mapid))
logger.info('"%s" to: %s' % (wmap._mapid, outname))
datatype, wmap.name))
logger.info('"%s" to: %s' % (wmap.name, outname))
savez(outname, **pcovs)

logger.info('Saving data covariance matrixes...')
for wmap in composite.wavemaps:
dcovs = {
list2string(dataset.nslc_id):
dataset.covariance.data
for dataset in wmap.datasets}

outname = pjoin(
results_path, '%s_C_d_%s' % (
datatype, wmap._mapid))
logger.info('"%s" to: %s' % (wmap._mapid, outname))
savez(outname, **dcovs)

elif datatype == 'geodetic':
for ifgs, attribute in heart.results_for_export(
results, datatype=datatype):
Expand Down
53 changes: 53 additions & 0 deletions src/backend.py
Original file line number Diff line number Diff line change
Expand Up @@ -38,6 +38,7 @@
from beat.utility import load_objects, dump_objects, \
ListArrayOrdering, ListToArrayBijection
from beat.covariance import calc_sample_covariance
from beat import utility

from pyrocko import util
from time import time
Expand Down Expand Up @@ -830,6 +831,58 @@ def load_sampler_params(project_dir, stage_number, mode):
sample_p_outname)
return load_objects(stage_path)

def get_result_point(stage, config, point_llk='max'):
"""
Return point of a given stage result.

Parameters
----------
stage : :class:`models.Stage`
config : :class:`config.BEATConfig`
point_llk : str
with specified llk(max, mean, min).

Returns
-------
dict
"""
if config.sampler_config.name == 'Metropolis':
if stage.step is None:
raise AttributeError(
'Loading Metropolis results requires'
' sampler parameters to be loaded!')

sc = config.sampler_config.parameters
from beat.sampler.metropolis import get_trace_stats
pdict, _ = get_trace_stats(
stage.mtrace, stage.step, sc.burn, sc.thin)
point = pdict[point_llk]

elif config.sampler_config.name == 'SMC':
llk = stage.mtrace.get_values(
varname='like',
combine=True)

posterior_idxs = utility.get_fit_indexes(llk)

point = stage.mtrace.point(idx=posterior_idxs[point_llk])

elif config.sampler_config.name == 'PT':
params = config.sampler_config.parameters
llk = stage.mtrace.get_values(
varname='like',
burn=int(params.n_samples * params.burn),
thin=params.thin)

posterior_idxs = utility.get_fit_indexes(llk)

point = stage.mtrace.point(idx=posterior_idxs[point_llk])

else:
raise NotImplementedError(
'Sampler "%s" is not supported!' % config.sampler_config.name)

return point

def concatenate_traces(mtraces):
"""
Expand Down
55 changes: 0 additions & 55 deletions src/plotting.py
Original file line number Diff line number Diff line change
Expand Up @@ -10,7 +10,6 @@

from beat import utility
from beat.models import Stage
from beat.sampler.metropolis import get_trace_stats
from beat.heart import init_seismic_targets, init_geodetic_targets
from beat.colormap import slip_colormap

Expand Down Expand Up @@ -473,60 +472,6 @@ def plot_log_cov(cov_mat):
plt.colorbar(im)
plt.show()


def get_result_point(stage, config, point_llk='max'):
"""
Return point of a given stage result.

Parameters
----------
stage : :class:`models.Stage`
config : :class:`config.BEATConfig`
point_llk : str
with specified llk(max, mean, min).

Returns
-------
dict
"""
if config.sampler_config.name == 'Metropolis':
if stage.step is None:
raise AttributeError(
'Loading Metropolis results requires'
' sampler parameters to be loaded!')

sc = config.sampler_config.parameters
pdict, _ = get_trace_stats(
stage.mtrace, stage.step, sc.burn, sc.thin)
point = pdict[point_llk]

elif config.sampler_config.name == 'SMC':
llk = stage.mtrace.get_values(
varname='like',
combine=True)

posterior_idxs = utility.get_fit_indexes(llk)

point = stage.mtrace.point(idx=posterior_idxs[point_llk])

elif config.sampler_config.name == 'PT':
params = config.sampler_config.parameters
llk = stage.mtrace.get_values(
varname='like',
burn=int(params.n_samples * params.burn),
thin=params.thin)

posterior_idxs = utility.get_fit_indexes(llk)

point = stage.mtrace.point(idx=posterior_idxs[point_llk])

else:
raise NotImplementedError(
'Sampler "%s" is not supported!' % config.sampler_config.name)

return point


def plot_quadtree(ax, data, target, cmap, colim, alpha=0.8):
"""
Plot UnwrappedIFG displacements on the respective quadtree rectangle.
Expand Down