Commit babfa76c authored by Sean Leavey's avatar Sean Leavey

Merge branch 'develop'

parents aec40581 6b2557ae
Pipeline #112688 passed with stage
in 26 minutes and 22 seconds
......@@ -95,7 +95,7 @@ test validation:3.7:
# Generate the documentation only on creation of new tags.
deploy pages:
image: python:3.7
image: python:3.8
stage: deploy
only:
refs:
......@@ -125,7 +125,7 @@ deploy pages:
# Generate PyPI release only on creation of new tags.
deploy pypi:
image: python:3.7
image: python:3.8
stage: deploy
only:
refs:
......
......@@ -19,9 +19,6 @@ STATICDIR = _static
### Dynamic plot dependencies
# Static files
ZEROSTATICDEPS = $(STATICDIR)/liso-input-node-graph.svg
# Zero Python files.
ZEROPYTHONDEPS = $(STATICDIR)/liso-input-response.svg $(STATICDIR)/solution-combination.svg \
$(STATICDIR)/resistor-current-noise.svg
# Native plots with LISO files.
ZEROLISODEPS = $(STATICDIR)/liso-two-noises.svg
# Native/LISO script comparisons.
......@@ -83,11 +80,6 @@ $(CLIOPAMPGAINDEP):
@echo "Generating $@"
zero library search "gbw > 800M & ((vnoise < 10n & inoise < 10p) | (vnoise < 100n & inoise < 1p)) & model != OP00" --no-plot-gain --save-gain-figure $@ --fstop 1M
# Generate SVG plots of Zero output using Python script.
$(STATICDIR)/%.svg: $(STATICDIR)/zero-python/%.py
@echo "Generating $@ from $<"
python $< $@
# Generate SVG plot of two noise example.
$(STATICDIR)/liso-two-noises.svg:
@echo "Generating $@"
......
"""Documentation example. Requires target figure filename as argument."""
import sys
from zero.liso import LisoInputParser
parser = LisoInputParser()
parser.parse("""
c c1 10u gnd n1
r r1 430 n1 nm
r r2 43k nm nout
c c2 47p nm nout
op o1 lt1124 nin nm nout
freq log 1 100k 100
uinput nin 0
uoutput nout:db:deg
""")
solution = parser.solution()
plotter = solution.plot_responses()
plotter.save(sys.argv[1])
"""Documentation example. Requires target figure filename as argument."""
import sys
import numpy as np
from zero import Circuit
from zero.analysis import AcNoiseAnalysis
from zero.noise import VoltageNoise
# Create a new noise type.
class ResistorCurrentNoise(VoltageNoise):
"""Resistor current noise source.
This models resistor current noise. See e.g. https://dcc.ligo.org/LIGO-T0900200/public
for more details. This noise depends on resistor composition and on its current. Be
careful when using this noise - it generally does not transfer to different circuits
with identical resistors as it depends on the voltage drop across the resistor.
Parameters
----------
vnoise : :class:`float`
The voltage noise at the specified frequency (V/sqrt(Hz)).
frequency : :class:`float`
The frequency at which the specified voltage noise is defined (Hz).
exponent : :class:`float`
The frequency exponent to use for calculating the frequency response.
"""
def __init__(self, vnoise, frequency=1.0, exponent=0.5, **kwargs):
super().__init__(**kwargs)
self.vnoise = vnoise
self.frequency = frequency
self.exponent = exponent
def noise_voltage(self, frequencies, **kwargs):
return self.vnoise * self.frequency / frequencies ** self.exponent
@property
def label(self):
return f"RE({self.component.name})"
# 1000 frequencies between 0.1 Hz to 10 kHz
frequencies = np.logspace(-1, 4, 1000)
# Create circuit object.
circuit = Circuit()
# Add components.
circuit.add_capacitor(value="10u", node1="gnd", node2="n1")
circuit.add_resistor(value="430", node1="n1", node2="nm", name="r1")
circuit.add_resistor(value="43k", node1="nm", node2="nout")
circuit.add_capacitor(value="47p", node1="nm", node2="nout")
circuit.add_library_opamp(model="LT1124", node1="gnd", node2="nm", node3="nout")
# Add resistor current noise to r1 with 10 nV/sqrt(Hz) at 1 Hz, with 1/f^2 drop-off.
r1 = circuit["r1"]
r1.add_noise(ResistorCurrentNoise(vnoise=1e-8, frequency=1.0, exponent=0.5))
# Solve circuit.
analysis = AcNoiseAnalysis(circuit=circuit)
solution = analysis.calculate(frequencies=frequencies, input_type="voltage", node="n1",
sink="nout", incoherent_sum=True)
# Plot.
plotter = solution.plot_noise(sink="nout")
plotter.save(sys.argv[1])
"""Documentation example. Requires target figure filename as argument."""
import sys
from zero.liso import LisoInputParser
# Create parser.
parser = LisoInputParser()
base_circuit = """
l l2 420n nlf nout
c c4 47p nlf nout
c c1 1n nrf gnd
r r1 1k nrf gnd
l l1 600n nrf n_l1_c2
c c2 330p n_l1_c2 n_c2_c3
c c3 33p n_c2_c3 nout
c load 20p nout gnd
freq log 100k 100M 1000
uoutput nout
"""
# Parse the base circuit.
parser.parse(base_circuit)
# Set the circuit input to the low frequency port.
parser.parse("uinput nlf 50")
# Ground the unused input.
parser.parse("r nrfsrc 5 nrf gnd")
# Calculate the solution.
solutionlf = parser.solution()
solutionlf.name = "LF Circuit"
# Reset the parser's state.
parser.reset()
# Parse the base circuit.
parser.parse(base_circuit)
# Set the input to the radio frequency port.
parser.parse("uinput nrf 50")
# Ground the unused input.
parser.parse("r nlfsrc 5 nlf gnd")
# Calculate the solution.
solutionrf = parser.solution()
solutionrf.name = "RF Circuit"
# Combine the solutions. By default, this keeps the functions from each source solution in different
# groups in the resulting solution. This makes the plot show the functions with different styles and
# shows the source solution's name as a suffix on each legend label.
solution = solutionlf.combine(solutionrf)
# Plot.
plotter = solution.plot_responses()
plotter.save(sys.argv[1])
......@@ -65,7 +65,8 @@ return the corresponding noise.
Here is an example of defining a resistor current noise source and using it in a circuit:
.. code-block:: python
.. plot::
:include-source:
import numpy as np
from zero import Circuit
......@@ -129,6 +130,3 @@ Here is an example of defining a resistor current noise source and using it in a
# Plot.
solution.plot_noise(sink="nout")
solution.show()
.. image:: /_static/resistor-current-noise.svg
......@@ -43,7 +43,8 @@ extensions = [
'sphinx.ext.coverage',
'sphinx.ext.viewcode',
'sphinx.ext.napoleon',
'sphinx_click.ext'
'sphinx_click.ext',
'matplotlib.sphinxext.plot_directive',
]
# Add any paths that contain templates here, relative to this directory.
......@@ -378,3 +379,11 @@ napoleon_use_param = False
# show return type on same line as description
napoleon_use_rtype = False
# -- Options for Matplotlib plot directive --------------------------------
# Code to run before plot scripts.
plot_pre_code = None
# Plot formats to show under each generated plot.
plot_formats = [('png', 200), 'pdf']
This diff is collapsed.
.. include:: /defs.txt
Input file parsing
==================
.. code:: python
>>> from zero.liso import LisoInputParser
|Zero| is capable of parsing :ref:`most <liso/input:Known incompatibilities>` LISO input files.
To start, create a new parser:
.. code:: python
>>> parser = LisoInputParser()
To parse a LISO circuit, either call the :meth:`~.LisoParser.parse` method with text:
.. code:: python
>>> parser.parse("""
c c1 10u gnd n1
r r1 430 n1 nm
r r2 43k nm nout
c c2 47p nm nout
op o1 lt1124 nin nm nout
freq log 1 100k 100
uinput nin 0
uoutput nout:db:deg
""")
Or point it to a file using the :code:`path` parameter:
.. code:: python
>>> parser.parse(path="/path/to/liso/script.fil")
Get the solution with :meth:`~.LisoParser.solution` and plot and show it with
:meth:`.Solution.plot` and :meth:`.Solution.show`:
.. code:: python
>>> solution = parser.solution()
>>> solution.plot()
>>> solution.show()
.. image:: /_static/liso-input-response.svg
You can at any time list the circuit's constituent components:
.. code-block:: python
>>> parser.circuit
Circuit with 6 components and 5 nodes
1. c1 [in=gnd, out=n1, C=1e-05]
2. c2 [in=nm, out=nout, C=4.7e-11]
3. input [in=gnd, out=nin, Z=default]
4. o1 [in+=nin, in-=nm, out=nout, model=LT1124]
5. r1 [in=n1, out=nm, R=430.0]
6. r2 [in=nm, out=nout, R=43000.0]
......@@ -212,7 +212,8 @@ The resulting solution's group names can be changed using :meth:`.rename_group`.
Here is an example of solution combination using a :ref:`LISO model <liso/index:LISO compatibility>`
of an RF summing box with two inputs and one output:
.. code-block:: python
.. plot::
:include-source:
from zero.liso import LisoInputParser
......@@ -265,8 +266,6 @@ of an RF summing box with two inputs and one output:
solution.plot()
solution.show()
.. image:: /_static/solution-combination.svg
.. hint::
The above example makes a call to :meth:`~.Solution.plot`. This relies on :ref:`default
......
......@@ -32,7 +32,17 @@ plot:
sum_greyscale_cycle_start: 1 # black
sum_greyscale_cycle_stop: 0.2 # light grey
sum_greyscale_cycle_count: 4 # number of steps
grid_zorder: 1
# Plot grid options.
grid:
zorder: 1
alpha_major: 1
alpha_minor: 0.3
# BodePlotter options.
bode:
# If not limits are specified, show ±180° (plus a little extra margin).
show_full_phase_limits: true
# Matplotlib configuration overrides.
matplotlib:
......
......@@ -566,6 +566,76 @@ class MatplotlibPlotter(BasePlotter, metaclass=abc.ABCMeta):
self.legend = legend
self.legend_loc = legend_loc
def _expand_linear_axis_limits(self, dmin, dmax, step, margin):
"""Intelligently expand specified limits for a linearly scaled axis.
This is similar in behaviour to Matplotlib when axes.autolimit_mode = round_numbers, except
it applies the autoscale margin rcparam to the current y-axis limits before deciding whether
to add another major tick step to the limits.
The reason why this method might be used instead of the default Matplotlib behaviour is to
ensure that when axis data is close to the upper or lower limit of its axis, the axis limits
are expanded by a full `step` on top of the closest higher/lower step (for upper/lower
limits, respectively). By default, Matplotlib will simply expand by the margin settings,
which is by default 5%. Expanding by a full step is in some contexts neater.
"""
# closeto, le and ge adapted from :class:`matplotlib.ticker._Edge_integer` (3.2.0).
def closeto(ms, edge):
tol = 1e-10
return abs(ms - edge) < tol
def le(x):
d, m = divmod(x, step)
if closeto(m / step, 1):
return (d + 1)
return d
def ge(x):
d, m = divmod(x, step)
if closeto(m / step, 0):
return d
return (d + 1)
factor = 1 + margin
vmin = le(dmin * factor) * step
vmax = ge(dmax * factor) * step
if vmin == vmax:
vmin -= 1
vmax += 1
return vmin, vmax
def _expand_log_axis_limits(self, dmin, dmax, base, margin):
"""Intelligently expand specified limits for a log scaled axis.
This is similar in behaviour to Matplotlib when axes.autolimit_mode = round_numbers, except
it applies the autoscale margin rcparam to the current y-axis limits before deciding whether
to add another major tick step to the limits.
The reason why this method might be used instead of the default Matplotlib behaviour is to
ensure that when axis data is close to the upper or lower limit of its axis, the axis limits
are expanded by a full `base` on top of the closest higher/lower multiple of `base` (for
upper/lower limits, respectively). By default, Matplotlib will simply expand by the margin
settings, which is by default 5%. Expanding by a full step is in some contexts neater.
"""
def _decade_less_equal(x, base):
return (x if x == 0 else
-_decade_greater_equal(-x, base) if x < 0 else
base ** np.floor(np.log(x) / np.log(base)))
def _decade_greater_equal(x, base):
return (x if x == 0 else
-_decade_less_equal(-x, base) if x < 0 else
base ** np.ceil(np.log(x) / np.log(base)))
factor = 1 + margin
vmin = _decade_less_equal(dmin * factor, base)
vmax = _decade_greater_equal(dmax * factor, base)
return vmin, vmax
@property
def figure(self):
if self._figure is None:
......@@ -670,7 +740,7 @@ class MplGroupPlotter(MatplotlibPlotter, BaseGroupPlotter, metaclass=abc.ABCMeta
class BodePlotter(MplGroupPlotter):
def __init__(self, scale_db=True, xlim=None, mag_ylim=None, phase_ylim=None, xlabel=None,
ylabel_mag=None, ylabel_phase=None, db_tick_major_step=20, db_tick_minor_step=10,
phase_tick_major_step=30, phase_tick_minor_step=15, **kwargs):
phase_tick_major_step=45, phase_tick_minor_step=15, **kwargs):
super().__init__(**kwargs)
self.scale_db = scale_db
self.xlim = xlim
......@@ -719,8 +789,13 @@ class BodePlotter(MplGroupPlotter):
ax1.set_ylabel(self.ylabel_mag)
if self.ylabel_phase is not None:
ax2.set_ylabel(self.ylabel_phase)
ax1.grid(zorder=CONF["plot"]["grid_zorder"])
ax2.grid(zorder=CONF["plot"]["grid_zorder"])
gridconf = CONF["plot"]["grid"]
ax1.grid(which="major", alpha=gridconf["alpha_major"], zorder=gridconf["zorder"])
ax1.grid(which="minor", alpha=gridconf["alpha_minor"], zorder=gridconf["zorder"])
ax2.grid(which="major", alpha=gridconf["alpha_major"], zorder=gridconf["zorder"])
ax2.grid(which="minor", alpha=gridconf["alpha_minor"], zorder=gridconf["zorder"])
# Magnitude and phase tick locators.
if self.scale_db:
ax1.yaxis.set_major_locator(MultipleLocator(base=self.db_tick_major_step))
......@@ -746,14 +821,48 @@ class BodePlotter(MplGroupPlotter):
def _finalise_plot(self):
# Update the legend.
self.ax1.legend(loc=self.legend_loc)
# Set limits.
self._set_limits()
def _set_limits(self):
"""Set appropriate plot limits, if not explicitly specified."""
if self.xlim is not None:
self.ax1.set_xlim(self.xlim)
self.ax2.set_xlim(self.xlim)
if self.mag_ylim is not None:
self.ax1.set_ylim(self.mag_ylim)
if self.phase_ylim is not None:
self.ax2.set_ylim(self.phase_ylim)
mag_ylim = self.mag_ylim
if mag_ylim is None:
LOGGER.info("No magnitude y-axis limits specified; attempting to use reasonable values.")
yinterval = self.ax1.yaxis.get_data_interval()
_, mag_ymargin = self.ax1.margins()
if self.scale_db:
# Round up/down to nearest multiple.
# Note: you might be tempted to turn on axes.autolimit_mode = round_numbers here,
# which also does this rounding up/down, but also sets the log x-axis limits a full
# decade lower/above.
mag_ylim = self._expand_linear_axis_limits(*yinterval, self.db_tick_minor_step,
mag_ymargin)
else:
mag_ylim = self._expand_log_axis_limits(*yinterval, 10, mag_ymargin)
self.ax1.set_ylim(mag_ylim)
phase_ylim = self.phase_ylim
if phase_ylim is None:
_, phase_ymargin = self.ax2.margins()
if CONF["plot"]["bode"]["show_full_phase_limits"]:
LOGGER.info("No phase y-axis limits specified; defaulting to full span due to "
"show_full_phase_limits setting.")
yinterval = (-180, 180)
else:
LOGGER.info("No phase y-axis limits specified; attempting to use reasonable values.")
yinterval = self.ax2.yaxis.get_data_interval()
phase_ylim = self._expand_linear_axis_limits(*yinterval, self.phase_tick_minor_step,
phase_ymargin)
self.ax2.set_ylim(phase_ylim)
class SpectralDensityPlotter(MplGroupPlotter):
......@@ -774,19 +883,23 @@ class SpectralDensityPlotter(MplGroupPlotter):
def _create_figure(self):
figure = super()._create_figure()
axis = figure.add_subplot(111)
ax = figure.add_subplot(111)
# Draw labels etc.
if self.title is not None:
# Use ax1 since it's at the top. We could use figure.suptitle but this doesn't
# behave with tight_layout.
axis.set_title(self.title)
ax.set_title(self.title)
if self.legend:
axis.legend(loc=self.legend_loc)
ax.legend(loc=self.legend_loc)
if self.xlabel is not None:
axis.set_xlabel(self.xlabel)
ax.set_xlabel(self.xlabel)
if self.ylabel is not None:
axis.set_ylabel(self.ylabel)
axis.grid(zorder=CONF["plot"]["grid_zorder"])
ax.set_ylabel(self.ylabel)
gridconf = CONF["plot"]["grid"]
ax.grid(which="major", alpha=gridconf["alpha_major"], zorder=gridconf["zorder"])
ax.grid(which="minor", alpha=gridconf["alpha_minor"], zorder=gridconf["zorder"])
# Magnitude and phase tick locators.
return figure
......@@ -820,11 +933,23 @@ class SpectralDensityPlotter(MplGroupPlotter):
def _finalise_plot(self):
# Update the legend.
self.axis.legend()
# Set limits.
self._set_limits()
def _set_limits(self):
"""Set appropriate plot limits, if not explicitly specified."""
if self.xlim is not None:
self.axis.set_xlim(self.xlim)
if self.ylim is not None:
self.axis.set_ylim(self.ylim)
ylim = self.ylim
if ylim is None:
LOGGER.info("No y-axis limits specified; attempting to use reasonable values.")
yinterval = self.axis.yaxis.get_data_interval()
_, ymargin = self.axis.margins()
ylim = self._expand_log_axis_limits(*yinterval, 10, ymargin)
self.axis.set_ylim(ylim)
class OpAmpGainPlotter(BodePlotter):
......
......@@ -958,6 +958,7 @@ class Solution:
default_none_param_names = ("group", "groups", "source", "sources", "sink", "sinks")
default_none_params = (group, groups, source, sources, sink, sinks)
if all([param is None for param in default_none_params]):
# No filters have been applied.
groups = self.default_responses
else:
groups = self.filter_responses(source=source, sources=sources, sink=sink, sinks=sinks,
......@@ -976,8 +977,8 @@ class Solution:
return plotter
def plot_noise(self, figure=None, group=None, groups=None, source=None, sources=None, sink=None,
sinks=None, type=None, types=None, show_sums=True, xlabel=None, ylabel=None,
**kwargs):
sinks=None, type=None, types=None, show_individual=True, show_sums=True,
xlabel=None, ylabel=None, **kwargs):
"""Plot noise.
Note: if only some of "groups", "sources", "sinks", "types" are specified, the others
......@@ -996,6 +997,8 @@ class Solution:
The sink(s) to plot noise at. If None, all matched sinks are plotted.
type, types : :class:`str` or list of :class:`str`, optional
The noise type(s) to plot. If None, all matched noise types are plotted.
show_individual : :class:`bool`, optional
Plot any individual noise spectra contained in this solution.
show_sums : :class:`bool`, optional
Plot any sums contained in this solution.
xlabel, ylabel : :class:`str`, optional
......@@ -1023,24 +1026,32 @@ class Solution:
default_none_param_names = ("group", "groups", "source", "sources", "sink", "sinks", "type",
"types")
default_none_params = (group, groups, source, sources, sink, sinks, type, types)
grouped_noise = {}
if all([param is None for param in default_none_params]):
# Filter against sum flag.
groups = self._filter_default_noise()
# No filters have been applied.
if show_individual:
grouped_noise = self._merge_groupsets(grouped_noise, self._filter_default_noise())
if show_sums:
groups = self._merge_groupsets(groups, self.default_noise_sums)
grouped_noise = self._merge_groupsets(grouped_noise, self.default_noise_sums)
else:
groups = self.filter_noise(source=source, sources=sources, sink=sink, sinks=sinks,
group=group, groups=groups, type=type, types=types)
if show_individual:
individual = self.filter_noise(source=source, sources=sources, sink=sink,
sinks=sinks, group=group, groups=groups, type=type,
types=types)
grouped_noise = self._merge_groupsets(grouped_noise, individual)
if show_sums:
groups = self._merge_groupsets(groups, self.noise_sums)
if not groups:
sums = self.filter_noise_sums(sink=sink, sinks=sinks, group=group, groups=groups,
type=type, types=types)
grouped_noise = self._merge_groupsets(grouped_noise, sums)
if not grouped_noise:
filters = ", ".join([f"'{param}'" for param in default_none_param_names])
raise NoDataException(f"No noise spectra found. Consider setting one of {filters}.")
if ylabel is None:
# Show all plotted noise units.
unit_tex = []
noise_units = set()
for spectra in groups.values():
for spectra in grouped_noise.values():
noise_units.update([spectral_density.sink_unit for spectral_density in spectra])
if noise_units:
unit_tex = [r"$\frac{\mathrm{" + unit + r"}}{\sqrt{\mathrm{Hz}}}$"
......@@ -1048,11 +1059,11 @@ class Solution:
unit = ", ".join(unit_tex)
ylabel = r"$\bf{Noise}$" + f" ({unit})"
# Add reference functions.
groups[self.DEFAULT_REF_GROUP_NAME] = self.noise_references
grouped_noise[self.DEFAULT_REF_GROUP_NAME] = self.noise_references
# Draw plot.
plotter = self.noise_plotter(figure=figure, xlabel=xlabel, ylabel=ylabel,
hidden_group_names=[self.DEFAULT_GROUP_NAME], **kwargs)
plotter.plot_groups(groups)
plotter.plot_groups(grouped_noise)
self._last_plotter = plotter
return plotter
......@@ -1069,27 +1080,39 @@ class Solution:
if self.has_responses:
data += "\n\tResponses:"
for response in self.responses:
data += f"\n\t\t{response}"
if self.is_default_response(response):
data += " (default)"
for group, responses in self.responses.items():
if group == self.DEFAULT_GROUP_NAME:
data += f"\n\t\tDefault group:"
else:
data += f"\n\t\tGroup '{group}':"
for response in responses:
data += f"\n\t\t\t{response}"
if self.is_default_response(response):
<