773 lines
30 KiB
Python
Executable File
773 lines
30 KiB
Python
Executable File
#!/usr/bin/env python
|
|
"""
|
|
@package pmsco.reports.rfactor
|
|
graphics rendering module to show R-factor plots.
|
|
|
|
the module can be used in several different ways:
|
|
|
|
1. via the command line on a pmsco database or .dat results file.
|
|
this is the most simple but least flexible way.
|
|
2. via python functions on given population arrays or database queries.
|
|
this is the most flexible way but requires understanding of the required data formats.
|
|
3. as a _report_ function during an optimization.
|
|
reporting is configurable in the run file.
|
|
|
|
@author Matthias Muntwiler, matthias.muntwiler@psi.ch
|
|
|
|
@copyright (c) 2021 by Paul Scherrer Institut @n
|
|
Licensed under the Apache License, Version 2.0 (the "License"); @n
|
|
you may not use this file except in compliance with the License.
|
|
You may obtain a copy of the License at
|
|
http://www.apache.org/licenses/LICENSE-2.0
|
|
"""
|
|
|
|
import argparse
|
|
import itertools
|
|
import logging
|
|
import numpy as np
|
|
from pathlib import Path
|
|
from scipy.interpolate import griddata
|
|
import sys
|
|
|
|
if __name__ == "__main__":
|
|
pmsco_root = Path(__file__).resolve().parent.parent.parent
|
|
if str(pmsco_root) not in sys.path:
|
|
sys.path.insert(0, str(pmsco_root))
|
|
|
|
import pmsco.reports.results as rp_results
|
|
import pmsco.database.util as db_util
|
|
import pmsco.database.query as db_query
|
|
from pmsco.reports.base import ProjectReport
|
|
|
|
logger = logging.getLogger(__name__)
|
|
|
|
try:
|
|
from matplotlib.figure import Figure
|
|
from matplotlib.ticker import MaxNLocator
|
|
from matplotlib.backends.backend_agg import FigureCanvasAgg as FigureCanvas
|
|
# from matplotlib.backends.backend_pdf import FigureCanvasPdf
|
|
# from matplotlib.backends.backend_svg import FigureCanvasSVG
|
|
except ImportError:
|
|
Figure = None
|
|
FigureCanvas = None
|
|
MaxNLocator = None
|
|
logger.warning("error importing matplotlib. graphics rendering disabled.")
|
|
|
|
|
|
def plot_rfactor_1d(filename, data, params, title=None, canvas=None):
|
|
"""
|
|
one-dimensional R-factor versus parameter plot
|
|
|
|
this is a low-level function containing just the plotting commands from numpy arrays.
|
|
|
|
the graphics file format can be changed by providing a specific canvas. default is PNG.
|
|
|
|
this function requires the matplotlib module.
|
|
if it is not available, the function raises an error.
|
|
|
|
@param filename: path and base name of the output file without extension.
|
|
a generation index and the file extension according to the file format are appended.
|
|
@param data: structured ndarray containing parameter and R-factor values.
|
|
the '_rfac' column is required and must contain R-factor values.
|
|
the column specified by `param_name` is required and must contain the parameter values.
|
|
other columns may be present and are ignored.
|
|
@param params: dictionary of parameter and range.
|
|
the key corresponds to the parameter column of the data array.
|
|
the value is a list [minimum, maximum] that defines the axis range.
|
|
if the value is None, the axes are auto-scaled.
|
|
the dictionary should contain one parameter exactly.
|
|
only the first one is plotted.
|
|
@param title: (str) title of the chart.
|
|
default: derived from parameter names.
|
|
@param canvas: a FigureCanvas class reference from a matplotlib backend.
|
|
if None, the default FigureCanvasAgg is used which produces a bitmap file in PNG format.
|
|
some other options are:
|
|
matplotlib.backends.backend_pdf.FigureCanvasPdf or
|
|
matplotlib.backends.backend_svg.FigureCanvasSVG.
|
|
|
|
@return (str) path and name of the generated graphics file.
|
|
None if no file was generated due to an error.
|
|
"""
|
|
pname = list(params.keys())[0]
|
|
if canvas is None:
|
|
canvas = FigureCanvas
|
|
if canvas is None or Figure is None:
|
|
return None
|
|
if title is None:
|
|
title = f'{pname} R-factor map'
|
|
|
|
fig = Figure()
|
|
canvas(fig)
|
|
ax = fig.add_subplot(111)
|
|
|
|
ax.scatter(data[pname], data['_rfac'], c='b', marker='o', s=4.0)
|
|
# ax.grid(True)
|
|
rng = params[pname]
|
|
if rng is not None:
|
|
ax.set_xlim(rng)
|
|
ax.set_ylim([0., 1.])
|
|
ax.set_xlabel(pname)
|
|
ax.set_ylabel('R-factor')
|
|
ax.set_title(title)
|
|
|
|
out_filename = "{base}.{ext}".format(base=filename, ext=canvas.get_default_filetype())
|
|
try:
|
|
fig.savefig(out_filename)
|
|
except OSError:
|
|
logger.exception(f"exception while saving figure {out_filename}")
|
|
out_filename = None
|
|
|
|
return out_filename
|
|
|
|
|
|
def plot_rfactor_2d(filename, data, params, title=None, cmap=None, canvas=None):
|
|
"""
|
|
plot a two-dimensional R-factor scatter plot.
|
|
|
|
the plot consists of a pseudo-color scatter plot of R-factors projected on two parameter dimensions.
|
|
|
|
this is a low-level function containing just the plotting commands from numpy arrays.
|
|
|
|
the graphics file format can be changed by providing a specific canvas. default is PNG.
|
|
|
|
this function requires the matplotlib module.
|
|
if it is not available, the function raises an error.
|
|
|
|
@param filename: path and base name of the output file without extension.
|
|
a generation index and the file extension according to the file format are appended.
|
|
@param data: structured ndarray containing parameter and R-factor values.
|
|
@param params: dictionary of two parameters to be plotted.
|
|
the keys correspond to columns of the data array.
|
|
the values are lists [minimum, maximum] that define the axis range.
|
|
if the values are None, the axes are auto-scaled.
|
|
the dictionary must contain at least two parameters and should not contain more than two.
|
|
only the first two are plotted.
|
|
@param title: (str) title of the chart.
|
|
default: derived from parameter names.
|
|
@param cmap: (str) name of colour map supported by matplotlib.
|
|
default is 'plasma'.
|
|
other good-looking options are 'viridis', 'plasma', 'inferno', 'magma', 'cividis'.
|
|
@param canvas: a FigureCanvas class reference from a matplotlib backend.
|
|
if None, the default FigureCanvasAgg is used which produces a bitmap file in PNG format.
|
|
some other options are:
|
|
matplotlib.backends.backend_pdf.FigureCanvasPdf or
|
|
matplotlib.backends.backend_svg.FigureCanvasSVG.
|
|
|
|
@return (str) path and name of the generated graphics file.
|
|
None if no file was generated due to an error.
|
|
"""
|
|
pnames = list(params.keys())
|
|
if canvas is None:
|
|
canvas = FigureCanvas
|
|
if canvas is None or Figure is None:
|
|
return None
|
|
if cmap is None:
|
|
cmap = 'plasma'
|
|
if title is None:
|
|
title = f'{pnames[0]} - {pnames[1]} R-factor map'
|
|
|
|
fig = Figure()
|
|
canvas(fig)
|
|
ax = fig.add_subplot(111)
|
|
|
|
s = ax.scatter(data[pnames[0]], data[pnames[1]], s=5, c=data['_rfac'], cmap=cmap, vmin=0, vmax=1)
|
|
cb = ax.figure.colorbar(s, ax=ax)
|
|
cb.ax.set_ylabel("R-factor", rotation=-90, va="bottom")
|
|
rng = params[pnames[0]]
|
|
if rng is not None:
|
|
ax.set_xlim(rng)
|
|
rng = params[pnames[1]]
|
|
if rng is not None:
|
|
ax.set_ylim(rng)
|
|
ax.set_xlabel(pnames[0])
|
|
ax.set_ylabel(pnames[1])
|
|
ax.set_title(title)
|
|
|
|
out_filename = "{base}.{ext}".format(base=filename, ext=canvas.get_default_filetype())
|
|
try:
|
|
fig.savefig(out_filename)
|
|
except OSError:
|
|
logger.exception(f"exception while saving figure {out_filename}")
|
|
out_filename = None
|
|
|
|
return out_filename
|
|
|
|
|
|
def triplet_to_grid(x, y, z):
|
|
"""
|
|
reshape and interpolate three flat arrays to a two-dimensional grid array.
|
|
|
|
input is two coordinate arrays and one value array.
|
|
the coordinates arrays must contain coordinates of a rectangular grid with (per dimension) equidistant steps.
|
|
some defects (missing values, imprecise coordinate values) are accepted.
|
|
|
|
the function currently has a very basic implementation which relies on the following strong assumptions:
|
|
- the x and y arrays contain coordinates of a two-dimensional, rectangular grid
|
|
- jitter up to one third of step size is accepted.
|
|
- missing data points are interpolated (nearest neighbour interpolation).
|
|
a warning message is sent to the logger in such a case.
|
|
|
|
@param x: 1D array of x-coordinates. x, y and z must have the same shape.
|
|
@param y: 1D array of y-coordinates. x, y and z must have the same shape.
|
|
@param z: 1D array of values. x, y and z must have the same shape.
|
|
@return: three 2D arrays of same shape containing: x-coordinates, y-coordinates, values.
|
|
"""
|
|
x_sort = np.sort(x)
|
|
x_diff = np.diff(x_sort)
|
|
x_big_step = np.amax(x_diff)
|
|
x_steps = x_diff[np.where(x_diff > x_big_step / 3)]
|
|
x_step = np.amin(x_steps)
|
|
x_min = x_sort[0]
|
|
x_max = x_sort[-1]
|
|
nx = round((x_max - x_min) / x_step + 1)
|
|
|
|
y_sort = np.sort(y)
|
|
y_diff = np.diff(y_sort)
|
|
y_big_step = np.amax(y_diff)
|
|
y_steps = y_diff[np.where(y_diff > y_big_step / 3)]
|
|
y_step = np.amin(y_steps)
|
|
y_min = y_sort[0]
|
|
y_max = y_sort[-1]
|
|
ny = round((y_max - y_min) / y_step + 1)
|
|
|
|
if nx * ny != len(z):
|
|
logger.warning(f"grid plot: input array size ({len(z)}) does not match grid size ({nx}*{ny}). "
|
|
f"some array elements may be interpolated!")
|
|
|
|
grid_x, grid_y = np.mgrid[x_min:x_max:nx * 1j, y_min:y_max:ny * 1j]
|
|
grid_z = griddata((x, y), z, (grid_x, grid_y), method='nearest')
|
|
|
|
return grid_x, grid_y, grid_z
|
|
|
|
|
|
def plot_rfactor_grid(filename, data, params, title=None, cmap=None, canvas=None):
|
|
"""
|
|
plot a two-dimensional R-factor map.
|
|
|
|
the plot consists of a pseudo-color R-factor map on a Cartesian grid.
|
|
|
|
this is a low-level function containing just the plotting commands from numpy arrays.
|
|
the input data can be provided in a matrix or structured array.
|
|
in the former case, the coordinate scales are defined by the params argument.
|
|
in the latter case, each row of the array contains coordinates and value,
|
|
and the coordinate columns to be used are defined by the params argument.
|
|
the array is reshaped to a matrix internally.
|
|
|
|
the graphics file format can be changed by providing a specific canvas. default is PNG.
|
|
|
|
this function requires the matplotlib module.
|
|
if it is not available, the function raises an error.
|
|
|
|
@param filename: path and base name of the output file without extension.
|
|
a generation index and the file extension according to the file format are appended.
|
|
@param data: R-factor data.
|
|
this can be a two-dimensional numpy array of R-factor values
|
|
or a one-dimensional numpy structured array containing parameter R-factor values.
|
|
in the first case, for an array of shape (M, N),
|
|
the M (N) index corresponds to params[0] (params[1]), respectively.
|
|
M/params[0] is plotted on the horizontal axis.
|
|
in the latter case, the `_rfac` column must be a flat version of a rectangular array
|
|
with axes `params[0]` (horizontal) and `params[1]` (vertical).
|
|
@param params: dictionary of two parameters to be plotted.
|
|
the keys contain the parameter names,
|
|
the values are lists [minimum, maximum] that define the axis range.
|
|
if the values are None, the axes are auto-scaled.
|
|
the dictionary must contain at least two parameters and should not contain more than two.
|
|
only the first two are plotted.
|
|
params.keys()[0] runs along the horizontal axis.
|
|
@param title: (str) title of the chart.
|
|
default: derived from parameter names.
|
|
@param cmap: (str) name of colour map supported by matplotlib.
|
|
default is 'plasma'.
|
|
other good-looking options are 'viridis', 'plasma', 'inferno', 'magma', 'cividis'.
|
|
@param canvas: a FigureCanvas class reference from a matplotlib backend.
|
|
if None, the default FigureCanvasAgg is used which produces a bitmap file in PNG format.
|
|
some other options are:
|
|
matplotlib.backends.backend_pdf.FigureCanvasPdf or
|
|
matplotlib.backends.backend_svg.FigureCanvasSVG.
|
|
|
|
@return (str) path and name of the generated graphics file.
|
|
None if no file was generated due to an error.
|
|
"""
|
|
pnames = list(params.keys())
|
|
if canvas is None:
|
|
canvas = FigureCanvas
|
|
if canvas is None or Figure is None:
|
|
return None
|
|
if cmap is None:
|
|
cmap = 'plasma'
|
|
if title is None:
|
|
title = f'{pnames[0]} - {pnames[1]} R-factor map'
|
|
|
|
try:
|
|
scales = [params[pnames[0]][0], params[pnames[0]][1], params[pnames[1]][0], params[pnames[1]][1]]
|
|
extent = tuple(scales)
|
|
except IndexError:
|
|
scales = None
|
|
extent = None
|
|
|
|
if len(data.shape) == 1:
|
|
try:
|
|
grid_x, grid_y, data = triplet_to_grid(data[pnames[0]], data[pnames[1]], data['_rfac'])
|
|
except ValueError:
|
|
logger.exception("error reshaping data array")
|
|
return None
|
|
extent = (grid_x.min(), grid_x.max(), grid_y.min(), grid_y.max())
|
|
scales = None
|
|
else:
|
|
scales = None
|
|
|
|
nx, ny = data.shape
|
|
if extent is None:
|
|
extent = (0, nx - 1, 0, ny - 1)
|
|
xstep = (extent[1] - extent[0]) / (nx - 1)
|
|
ystep = (extent[3] - extent[2]) / (ny - 1)
|
|
extent = (extent[0] - xstep / 2, extent[1] + xstep / 2, extent[2] - ystep / 2, extent[3] + ystep / 2)
|
|
|
|
fig = Figure()
|
|
canvas(fig)
|
|
ax = fig.add_subplot(111)
|
|
im = ax.imshow(data.T, aspect='auto', cmap=cmap, vmin=0.0, vmax=1.0, extent=extent, origin='lower')
|
|
|
|
if scales is not None:
|
|
ax.set_xlim(scales[0:2])
|
|
if scales is not None:
|
|
ax.set_ylim(scales[2:4])
|
|
ax.set_xlabel(pnames[0])
|
|
ax.set_ylabel(pnames[1])
|
|
ax.set_title(title)
|
|
|
|
cb = ax.figure.colorbar(im, ax=ax)
|
|
cb.ax.set_ylabel("R-factor", rotation=-90, va="bottom")
|
|
|
|
out_filename = "{base}.{ext}".format(base=filename, ext=canvas.get_default_filetype())
|
|
try:
|
|
fig.savefig(out_filename)
|
|
except OSError:
|
|
logger.exception(f"exception while saving figure {out_filename}")
|
|
out_filename = None
|
|
|
|
return out_filename
|
|
|
|
|
|
class Rfactor2DPlot(ProjectReport):
|
|
"""
|
|
produce two-dimensional pseudo-colour R-factor scatter plots
|
|
|
|
this class collects and validates all parameters and data for generating a 2D R-factor plot.
|
|
|
|
the plots consist of a pseudo-colour scatter or image plot of R-factors versus two parameters.
|
|
scatter plots can contain an arbitrary number of results from any optimization algorithm.
|
|
image plots require results from a rectangular grid and should be used in 2D grid optimization mode only.
|
|
|
|
see the render_rfactor_2d function code for an example of how to use this class.
|
|
"""
|
|
|
|
## @var params
|
|
# parameter that should be plotted
|
|
#
|
|
# this should be a list of pairs (sequence or tuple of two strings) of parameter names.
|
|
# for each pair, a plot is generated.
|
|
# by default, all non-degenerate parameters are plotted in all combinations.
|
|
|
|
## @var result_data
|
|
# R-factor data for scatter plot
|
|
#
|
|
# pmsco.reports.results.ResultData object holding the data or filter criteria
|
|
# for the scatter plot in the background.
|
|
# the scatter plot shows a map of the model space in terms of R-factor.
|
|
|
|
## @var plot_type
|
|
# type of the generated plot (scatter or image)
|
|
#
|
|
# - 'auto': (default) choose automatically based on project mode and number of parameters.
|
|
# 'auto' will choose 'image' if the data fit on a rectangular grid and 'scatter' otherwise.
|
|
# - 'scatter': most general plot, suitable in any situation.
|
|
# - 'image': works in two-dimensional grid mode only.
|
|
|
|
def __init__(self):
|
|
super().__init__()
|
|
self._modes = ['genetic', 'grid', 'swarm']
|
|
self.result_data = rp_results.ResultData()
|
|
self.filename_format = "${base}-rfactor2d-${param0}-${param1}"
|
|
self.title_format = "${param0} ${param1} R-factor"
|
|
self.cmap = None
|
|
self.params = []
|
|
self.plot_type = 'auto'
|
|
|
|
def select_data(self, jobs=-1, calcs=None):
|
|
"""
|
|
query data from the database
|
|
|
|
@param jobs: filter by job.
|
|
the argument can be a singleton or sequence of orm.Job objects or numeric id.
|
|
if None, results from all jobs are loaded.
|
|
if -1 (default), results from the most recent job (by datetime field) are loaded.
|
|
|
|
@param calcs: the calcs argument is ignored.
|
|
|
|
@return: None
|
|
"""
|
|
|
|
with self.get_session() as session:
|
|
if jobs == -1:
|
|
jobs = db_query.query_newest_job(session)
|
|
|
|
self.result_data.reset_filters()
|
|
self.result_data.levels = {'scan': -1}
|
|
self.result_data.load_from_db(session, jobs=jobs)
|
|
|
|
if self._project:
|
|
self.result_data.set_model_space(self._project.model_space)
|
|
|
|
def create_report(self):
|
|
"""
|
|
generate the plots based on the stored attributes.
|
|
|
|
this method essentially loops over parameter combinations,
|
|
and compiles the input for plot_rfactor_1d.
|
|
|
|
combinations of parameter names are either set by the user in self.params
|
|
or constructed from all non-degenerate (not constant) parameters.
|
|
|
|
@return: list of created files
|
|
"""
|
|
|
|
vmin = self.result_data.model_space.min
|
|
vmax = self.result_data.model_space.max
|
|
nd_params = self.result_data.non_degenerate_params()
|
|
|
|
if self.params:
|
|
ppairs = [p for p in self.params if len(p) == 2 and p[0] in nd_params and p[1] in nd_params]
|
|
else:
|
|
pnames = sorted(list(nd_params), key=str.lower)
|
|
ppairs = [p for p in itertools.combinations(pnames, 2)]
|
|
|
|
plot_type = self.plot_type
|
|
if plot_type == 'auto':
|
|
if self._project and self._project.mode == 'grid' and len(nd_params) == 2:
|
|
plot_type = 'image'
|
|
else:
|
|
plot_type = 'scatter'
|
|
|
|
kwargs = {}
|
|
if self.cmap is not None:
|
|
kwargs['cmap'] = self.cmap
|
|
if self.canvas is not None:
|
|
kwargs['canvas'] = self.canvas
|
|
|
|
files = []
|
|
fdict = {'base': self.base_filename}
|
|
|
|
for ppair in ppairs:
|
|
fdict['param0'] = ppair[0]
|
|
fdict['param1'] = ppair[1]
|
|
filename = Path(self.report_dir, self.filename_format)
|
|
filename = Path(self.resolve_template(filename, fdict))
|
|
kwargs['title'] = self.resolve_template(self.title_format, fdict)
|
|
params = {ppair[0]: [vmin[ppair[0]], vmax[ppair[0]]],
|
|
ppair[1]: [vmin[ppair[1]], vmax[ppair[1]]]}
|
|
if plot_type == 'image':
|
|
of = plot_rfactor_grid(filename, self.result_data.values, params, **kwargs)
|
|
else:
|
|
of = plot_rfactor_2d(filename, self.result_data.values, params, **kwargs)
|
|
if of:
|
|
files.append(of)
|
|
|
|
return files
|
|
|
|
|
|
class Rfactor1DPlot(ProjectReport):
|
|
"""
|
|
produce one-dimensional R-factor plots
|
|
|
|
this class collects and validates all parameters and data for generating 1D R-factor plots.
|
|
|
|
the plots consist of scatter plots of R-factor versus one model parameter.
|
|
this can contain an arbitrary number of results.
|
|
|
|
see the render_rfactor_1d function code for an example of how to use this class.
|
|
"""
|
|
|
|
## @var params
|
|
# parameter that should be plotted
|
|
#
|
|
# this should be a list of parameter names.
|
|
# for each parameter, a plot is generated.
|
|
# by default, all non-degenerate parameters are plotted.
|
|
|
|
## @var result_data
|
|
# R-factor data for scatter plot
|
|
#
|
|
# pmsco.reports.results.ResultData object holding the data or filter criteria.
|
|
|
|
def __init__(self):
|
|
super().__init__()
|
|
self._modes = ['bayesian', 'genetic', 'grid', 'swarm']
|
|
self.result_data = rp_results.ResultData()
|
|
self.filename_format = "${base}-rfactor1d-${param}"
|
|
self.title_format = "${param} R-factor"
|
|
self.cmap = None
|
|
self.params = []
|
|
|
|
def select_data(self, jobs=-1, calcs=None):
|
|
"""
|
|
query data from the database
|
|
|
|
@param jobs: filter by job.
|
|
the argument can be a singleton or sequence of orm.Job objects or numeric id.
|
|
if None, results from all jobs are loaded.
|
|
if -1 (default), results from the most recent job (by datetime field) are loaded.
|
|
|
|
@param calcs: the calcs argument is ignored.
|
|
|
|
@return: None
|
|
"""
|
|
|
|
with self.get_session() as session:
|
|
if jobs == -1:
|
|
jobs = db_query.query_newest_job(session)
|
|
|
|
self.result_data.reset_filters()
|
|
self.result_data.levels = {'scan': -1}
|
|
self.result_data.load_from_db(session, jobs=jobs)
|
|
|
|
if self._project:
|
|
self.result_data.set_model_space(self._project.model_space)
|
|
|
|
def create_report(self):
|
|
"""
|
|
generate the plots based on the stored attributes.
|
|
|
|
this method essentially loops over parameter combinations,
|
|
and compiles the input for plot_rfactor_1d.
|
|
|
|
the list of parameter names is either set by the user in self.params
|
|
or defaults to all non-degenerate (not constant) parameters.
|
|
|
|
@return: list of created files
|
|
"""
|
|
|
|
vmin = self.result_data.model_space.min
|
|
vmax = self.result_data.model_space.max
|
|
nd_params = self.result_data.non_degenerate_params()
|
|
|
|
if self.params:
|
|
pnames = self.params
|
|
else:
|
|
pnames = sorted(list(nd_params), key=str.lower)
|
|
|
|
kwargs = {}
|
|
if self.cmap is not None:
|
|
kwargs['cmap'] = self.cmap
|
|
if self.canvas is not None:
|
|
kwargs['canvas'] = self.canvas
|
|
|
|
files = []
|
|
fdict = {'base': self.base_filename}
|
|
|
|
for param in pnames:
|
|
fdict['param'] = param
|
|
filename = Path(self.report_dir, self.filename_format)
|
|
filename = Path(self.resolve_template(filename, fdict))
|
|
kwargs['title'] = self.resolve_template(self.title_format, fdict)
|
|
params = {param: [vmin[param], vmax[param]]}
|
|
of = plot_rfactor_1d(filename, self.result_data.values, params, **kwargs)
|
|
if of:
|
|
files.append(of)
|
|
|
|
return files
|
|
|
|
|
|
def render_rfactor_2d(output_file, values, model_space=None, title=None, plot_type=None, cmap=None, canvas=None):
|
|
"""
|
|
render two-dimensional R-factor scatter plots from database or result file.
|
|
|
|
the function generates one plot for each combination of non-degenerate parameters.
|
|
|
|
the function accepts input in one of the following forms:
|
|
- an open pmsco database session. the most recent job results are loaded.
|
|
- a result (.dat) file or numpy structured array.
|
|
the array must contain columns of regular model parameters and the _rfac column.
|
|
- a pmsco.optimizers.population.Population object with valid data.
|
|
|
|
the graphics file format can be changed by providing a specific canvas. default is PNG.
|
|
|
|
this is a high-level function with minimal input requirements.
|
|
some assumptions are made about which data to plot.
|
|
more fine grained control is provided by the Rfactor2DPlot class or the plot_rfactor_2d function.
|
|
the function code shows how to use the Rfactor2DPlot class.
|
|
|
|
@param output_file: path and base name of the output file without extension.
|
|
parameter names and the file extension according to the file format are appended.
|
|
@param values: result data to plot.
|
|
the input can be specified in several ways - see the main description for details.
|
|
@param model_space: model space can be a pmsco.project.ModelSpace object,
|
|
any object that contains the same min and max attributes as pmsco.project.ModelSpace,
|
|
or a dictionary with to keys 'min' and 'max' that provides the corresponding ModelSpace dictionaries.
|
|
by default, the model space boundaries are derived from the input data.
|
|
if a model_space is specified, only the parameters listed in it are plotted.
|
|
@param title: (str) title of the chart.
|
|
the title is a {}-style format string, where {base} is the output file name and {gen} is the generation.
|
|
default: derived from file name.
|
|
@param plot_type: (str) override plot type: 'scatter' or 'image'.
|
|
by default, @ref Rfactor2DPlot chooses a scatter plot, which works in all situations.
|
|
for data on a rectangular grid, an image type may be suitable and can be forced using this option.
|
|
@param cmap: (str) name of colour map supported by matplotlib.
|
|
default is 'plasma'.
|
|
@param canvas: a FigureCanvas class reference from a matplotlib backend.
|
|
if None, the default FigureCanvasAgg is used which produces a bitmap file in PNG format.
|
|
some other options are:
|
|
matplotlib.backends.backend_pdf.FigureCanvasPdf or
|
|
matplotlib.backends.backend_svg.FigureCanvasSVG.
|
|
|
|
@return (list of str) paths of the generated graphics files.
|
|
empty if an error occurred.
|
|
|
|
@raise TypeError if matplotlib is not available.
|
|
"""
|
|
data = rp_results.ResultData()
|
|
data.levels = {'scan': -1}
|
|
data.load_any(values)
|
|
if model_space is not None:
|
|
data.set_model_space(model_space)
|
|
|
|
plot = Rfactor2DPlot()
|
|
if plot_type:
|
|
plot.plot_type = plot_type
|
|
plot.canvas = canvas
|
|
plot.cmap = cmap
|
|
if title:
|
|
plot.title_format = title
|
|
else:
|
|
plot.title_format = "R-factor"
|
|
plot.report_dir = Path(output_file).parent
|
|
plot.filename_format = Path(output_file).name + "-${param0}-${param1}"
|
|
plot.validate(None)
|
|
plot.result_data = data
|
|
files = plot.create_report()
|
|
|
|
return files
|
|
|
|
|
|
def render_rfactor_1d(output_file, values, model_space=None, title=None, cmap=None, canvas=None):
|
|
"""
|
|
render one-dimensional R-factor scatter plots from database or result file.
|
|
|
|
the function generates one plot for each non-degenerate parameter.
|
|
|
|
the function accepts input in one of the following forms:
|
|
- an open pmsco database session. the most recent job results are loaded.
|
|
- a result (.dat) file or numpy structured array.
|
|
the array must contain columns of regular model parameters and the _rfac column.
|
|
- a pmsco.optimizers.population.Population object with valid data.
|
|
|
|
the graphics file format can be changed by providing a specific canvas. default is PNG.
|
|
|
|
this is a high-level function with minimal input requirements.
|
|
some assumptions are made about which data to plot.
|
|
more fine grained control is provided by the Rfactor1DPlot class or the plot_rfactor_1d function.
|
|
the function code shows how to use the Rfactor1DPlot class.
|
|
|
|
@param output_file: path and base name of the output file without extension.
|
|
parameter name and the file extension according to the file format are appended.
|
|
@param values: result data to plot.
|
|
the input can be specified in several ways - see the main description for details.
|
|
@param model_space: model space can be a pmsco.project.ModelSpace object,
|
|
any object that contains the same min and max attributes as pmsco.project.ModelSpace,
|
|
or a dictionary with to keys 'min' and 'max' that provides the corresponding ModelSpace dictionaries.
|
|
by default, the model space boundaries are derived from the input data.
|
|
if a model_space is specified, only the parameters listed in it are plotted.
|
|
@param title: (str) title of the chart.
|
|
the title is a {}-style format string, where {base} is the output file name and {gen} is the generation.
|
|
default: derived from file name.
|
|
@param cmap: (str) name of colour map supported by matplotlib.
|
|
default is 'plasma'.
|
|
@param canvas: a FigureCanvas class reference from a matplotlib backend.
|
|
if None, the default FigureCanvasAgg is used which produces a bitmap file in PNG format.
|
|
some other options are:
|
|
matplotlib.backends.backend_pdf.FigureCanvasPdf or
|
|
matplotlib.backends.backend_svg.FigureCanvasSVG.
|
|
|
|
@return (list of str) paths of the generated graphics files.
|
|
empty if an error occurred.
|
|
|
|
@raise TypeError if matplotlib is not available.
|
|
"""
|
|
data = rp_results.ResultData()
|
|
data.levels = {'scan': -1}
|
|
data.load_any(values)
|
|
if model_space is not None:
|
|
data.set_model_space(model_space)
|
|
|
|
plot = Rfactor1DPlot()
|
|
plot.canvas = canvas
|
|
plot.cmap = cmap
|
|
if title:
|
|
plot.title_format = title
|
|
else:
|
|
plot.title_format = "R-factor"
|
|
plot.report_dir = Path(output_file).parent
|
|
plot.filename_format = Path(output_file).name + "-${param}"
|
|
plot.validate(None)
|
|
plot.result_data = data
|
|
files = plot.create_report()
|
|
|
|
return files
|
|
|
|
|
|
def main(args=None):
|
|
parser = argparse.ArgumentParser(
|
|
formatter_class=argparse.ArgumentDefaultsHelpFormatter,
|
|
description="""
|
|
R-factor result plots for multiple-scattering optimization results
|
|
|
|
this module operates on result or database files and produces graphics files.
|
|
|
|
database files contain the complete information for all plot types.
|
|
data from the most recent job stored in the database is used.
|
|
|
|
the plot type is selected by the mode argument.
|
|
""")
|
|
parser.add_argument('results_file',
|
|
help="path to results file (.dat) or sqlite3 database file.")
|
|
parser.add_argument('output_file',
|
|
help="base name of output file. generation and extension will be appended.")
|
|
parser.add_argument('-m', '--mode', choices=['1d', '2d', 'grid'], default='1d',
|
|
help="plot type. "
|
|
"1d: R-factor versus one parameter; "
|
|
"2d: pseudo-color scatter plot versus two parameters; "
|
|
"grid: pseudo-color image plot (2d grid data only)")
|
|
parser.add_argument('-t', '--title', default=None,
|
|
help='graph title. may contain {gen} as a placeholder for the generation number.')
|
|
|
|
args, unknown_args = parser.parse_known_args(args=args)
|
|
|
|
kwargs = {}
|
|
if args.title is not None:
|
|
kwargs['title'] = args.title
|
|
|
|
if args.mode == '1d':
|
|
render_func = render_rfactor_1d
|
|
elif args.mode == '2d':
|
|
render_func = render_rfactor_2d
|
|
elif args.mode == 'grid':
|
|
render_func = render_rfactor_2d
|
|
kwargs['plot_type'] = 'image'
|
|
else:
|
|
raise ValueError(f"unknown graph mode {args.mode}")
|
|
|
|
if db_util.is_sqlite3_file(args.results_file):
|
|
import pmsco.database.access as db_access
|
|
db = db_access.DatabaseAccess()
|
|
db.connect(args.results_file)
|
|
with db.session() as session:
|
|
render_func(args.output_file, session, **kwargs)
|
|
else:
|
|
render_func(args.output_file, args.results_file, **kwargs)
|
|
|
|
|
|
if __name__ == '__main__':
|
|
main()
|
|
sys.exit(0)
|