Commit 30c1c4d2 authored by Samuel Rowlinson's avatar Samuel Rowlinson

Trying new tracing structure formats for mode-matching tools (WIP)

Want to allow for the use of arbitrary node beam parameter properties as target
variables in minimisation function routines. This requires a change to how tracing
results are initialised and stored as the beam parameter associated with the target
variable must remain as a valid reference during the beam tracing - i.e. no overwriting
of BeamParam instances in Model.__last_trace.data, only resetting of the q property.
parent babda885
......@@ -20,8 +20,8 @@ from .element import ModelElement
from .gaussian import BeamParam
from .model import Model
from .simulation import Simulation
from datetime import datetime
import pkg_resources
from datetime import datetime
import pkg_resources
setup_version=pkg_resources.get_distribution("finesse").version
[version, githash] = setup_version.split('+')
......@@ -30,14 +30,14 @@ print("""-----------------------------------------------------------------------
FINESSE {} {} (build {})
o_.-=. Frequency domain INterferomEter Simulation SoftwarE
(\\'\".\\| http://www.gwoptics.org/finesse/
.>' (_--.
_=/d ,^\\
~~ \\)-' '
/ |
.>' (_--.
_=/d ,^\\
~~ \\)-' '
/ |
' ' {}
------------------------------------------------------------------------""".format(
version, ' '*num_spaces, githash, datetime.now().strftime("%d %B %Y, %H:%M:%S")))
if __doc__:
def collect_submodules():
import os
......@@ -70,14 +70,11 @@ if __doc__:
# by "import finesse" - remove the ones which aren't from
# the toc (anything not exposed by this module import should
# not appear in the toc below)
from finesse.beamtrace import BeamTrace
#from finesse.beamtrace import BeamTrace
from finesse.utilities import check_name
from finesse.components import Node, Connector
toc = (
('Beam Tracing', (
BeamTrace,
)),
('Components', (
Node, Connector,
)),
......
......@@ -428,7 +428,8 @@ class Cavity(ModelElement):
self.__FSR = self._model._CLIGHT/self.__length
def _update_power(self):
from finesse.beamtrace import TraceType, Tracer
#from finesse.beamtrace import TraceType, Tracer
from finesse.tracing import Tracer, TraceType
tracer = Tracer(self._model)
power = tracer.run(TraceType.CAVPOWER, cavity_path=self.__path.data)
......@@ -466,7 +467,8 @@ class Cavity(ModelElement):
self.__resolution_y = self.__finesse*0.5*self.__rt_gouy_y/np.pi
def _update_ABCD(self):
from finesse.beamtrace import TraceType, Tracer
#from finesse.beamtrace import TraceType, Tracer
from finesse.tracing import Tracer, TraceType
tracer = Tracer(self._model)
self.__rt_ABCD_x, self.__rt_ABCD_y = tracer.run(
......
......@@ -442,7 +442,7 @@ class OpticalNode(Node):
"""
qx = None
if self._model.last_trace is not None:
qx, _, _ = self._model.last_trace.get(self, (None, None, None))
qx, _, _ = self._model.last_trace.data.get(self, (None, None, None))
# if beam tracing not run on the model yet then grab qx
# value from gauss_commands dict of model
if qx is None:
......@@ -466,7 +466,7 @@ class OpticalNode(Node):
"""
qy = None
if self._model.last_trace is not None:
_, qy, _ = self._model.last_trace.get(self, (None, None, None))
_, qy, _ = self._model.last_trace.data.get(self, (None, None, None))
# if beam tracing not run on the model yet then grab qy
# value from gauss_commands dict of model
if qy is None:
......
......@@ -43,8 +43,12 @@ class BeamParam(object):
self.__q = None
self.__lambda = wavelength
self.__nr = nr
self.__is_null = False
if len(args) == 1:
if not args and not kwargs:
self.__is_null = True
elif len(args) == 1:
self.__q = complex(args[0])
elif len(kwargs) == 1:
......@@ -84,6 +88,10 @@ class BeamParam(object):
def __repr__(self):
return "<%s (w0=%s, w=%s, z=%s) at %s>" % (self.__class__.__name__, self.w0, self.w, self.z, hex(id(self)))
@property
def is_null(self):
return self.__is_null
@property
def nr(self):
"""The refractive index associated with
......@@ -103,6 +111,11 @@ class BeamParam(object):
"""
return self.__q
@q.setter
def q(self, value):
self.__q = value
self.__is_null = False
@property
def z(self):
"""The relative distance to the waist of
......
......@@ -13,11 +13,14 @@ from copy import copy
from collections import Iterable
from collections import OrderedDict
from .beamtrace import BeamTrace, TraceType, Tracer, TracerOut
#from .beamtrace import BeamTrace, TraceType, Tracer, TracerOut
from .components import Port, Connector, Beamsplitter, Space, FrequencyGenerator, Cavity, Wire
from .components.node import QSetMode, NodeType, OpticalNode
from .components.general import CouplingType
from .detectors import Detector
from .tracing import TraceType, Tracer, TracerOut
from .tracing.beamtrace import BeamTrace
from .element import ModelElement
from .exceptions import NodeException
from .freeze import canFreeze, Freezable
......@@ -133,7 +136,7 @@ class Model(object):
self.__fsig = None
self.__tracer = Tracer(self)
self.__gauss_commands = OrderedDict()
self.__last_trace = None
self.__last_trace = TracerOut()
self.__lambda0 = 1064e-9
# Arbitrary factor epsilon, for conversion between e.g. mechanical
......@@ -543,6 +546,7 @@ class Model(object):
if node.type == NodeType.OPTICAL:
self.network.node[node]['optical'] = True
self.__last_trace.reset(node)
elif node.type == NodeType.MECHANICAL:
self.network.node[node]['mechanical'] = True
elif node.type == NodeType.ELECTRICAL:
......@@ -1135,25 +1139,25 @@ class Model(object):
raise NodeException("Specified node does not exist within the model", node)
if hasattr(q, "__getitem__") and not isinstance(q, complex):
self.__gauss_commands[node] = (q[0], q[1], QSetMode.USER)
self.__gauss_commands[node] = [q[0], q[1], QSetMode.USER]
else:
if direction is None:
raise ValueError("direction arg must not be None when setting the beam "
"parameter for a node in a given plane.")
if direction == "both":
self.__gauss_commands[node] = (q, q, QSetMode.USER)
self.__gauss_commands[node] = [q, q, QSetMode.USER]
elif direction == 'x':
if node in self.__gauss_commands:
self.__gauss_commands[node][0] = q
self.__gauss_commands[node][2] = QSetMode.USER
else:
self.__gauss_commands[node] = (q, None, QSetMode.USER)
self.__gauss_commands[node] = [q, None, QSetMode.USER]
elif direction == 'y':
if node in self.__gauss_commands:
self.__gauss_commands[node][1] = q
self.__gauss_commands[node][2] = QSetMode.USER
else:
self.__gauss_commands[node] = (None, q, QSetMode.USER)
self.__gauss_commands[node] = [None, q, QSetMode.USER]
def ABCD(self, from_node, to_node, via_node=None, direction=None):
"""Computes the ABCD matrix between `from_node` and `to_node`,
......@@ -1237,15 +1241,17 @@ class Model(object):
"must add a Cavity to the model or manually "
"set the q at any node in the model.")
trace_out = TracerOut()
for node in self.optical_nodes:
self.__last_trace.reset(node)
for node, qtuple in self.__gauss_commands.items():
trace_out.data[node] = (qtuple[0], qtuple[1], qtuple[2])
for i in range(2):
self.__last_trace.data[node][i].q = qtuple[i].q
self.__last_trace.data[node][2] = QSetMode.USER
if set_symmetric:
trace_out.data[node.opposite] = (
-qtuple[0].conjugate(),
-qtuple[1].conjugate(),
qtuple[2]
)
for i in range(2):
self.__last_trace.data[node.opposite][i].q = -qtuple[i].conjugate().q
self.__last_trace.data[node.opposite][2] = QSetMode.USER
for cav in self.__cavities.keys(): # compute eigenmodes first
cav.update()
......@@ -1255,20 +1261,20 @@ class Model(object):
# don't set the node to cavities' eigenmode if this node
# is in gauss_commands => use user set param preferentially
if cav.source in self.__gauss_commands.keys(): continue
trace_out.data[cav.source] = (cav.eigenmode_x, cav.eigenmode_y, QSetMode.AUTO)
self.__last_trace.data[cav.source][0].q = cav.eigenmode_x.q
self.__last_trace.data[cav.source][1].q = cav.eigenmode_y.q
self.__last_trace.data[cav.source][2] = QSetMode.AUTO
if set_symmetric:
trace_out.data[cav.source.opposite] = (
-cav.eigenmode_x.conjugate(),
-cav.eigenmode_y.conjugate(),
QSetMode.AUTO
)
self.__last_trace.data[cav.source.opposite][0].q = -cav.eigenmode_x.conjugate().q
self.__last_trace.data[cav.source.opposite][1].q = -cav.eigenmode_y.conjugate().q
self.__last_trace.data[cav.source.opposite][2] = QSetMode.AUTO
# start trace from explicitly user-defined node...
if startnode is not None:
self.__tracer.run(
TraceType.SIMULATION,
current=startnode,
trace_out=trace_out,
trace_out=self.__last_trace,
from_user=True,
set_symmetric=set_symmetric
)
......@@ -1284,7 +1290,7 @@ class Model(object):
self.__tracer.run(
TraceType.SIMULATION,
current=cav.source,
trace_out=trace_out,
trace_out=self.__last_trace,
from_cav=cav,
traced_cavs=traced_cavities,
set_symmetric=set_symmetric
......@@ -1296,7 +1302,7 @@ class Model(object):
self.__tracer.run(
TraceType.SIMULATION,
current=cav.source,
trace_out=trace_out,
trace_out=self.__last_trace,
from_cav=cav,
traced_cavs=traced_cavities,
set_symmetric=set_symmetric
......@@ -1307,25 +1313,25 @@ class Model(object):
self.__tracer.run(
TraceType.SIMULATION,
current=node,
trace_out=trace_out,
trace_out=self.__last_trace,
from_user=True,
set_symmetric=set_symmetric
)
if len(trace_out.data) != len(self.optical_nodes):
diff = list(set(self.optical_nodes).difference(trace_out.data.keys()))
if len(self.__last_trace.data) != len(self.optical_nodes):
diff = list(set(self.optical_nodes).difference(self.__last_trace.data.keys()))
# if we're not setting opposite node qs to -q* ...
if not set_symmetric:
for node in diff:
# ... but we've encountered a node which has no successors,
# then set node.opposite.q to -q* as it's the only option
if not list(self.optical_network.successors(node.opposite)):
qx, qy, _ = trace_out.data[node.opposite]
trace_out.data[node] = (
qx, qy, _ = self.__last_trace.data[node.opposite]
self.__last_trace.data[node] = (
-qx.conjugate(), -qy.conjugate(), QSetMode.AUTO
)
diff.remove(node)
for node, qs in trace_out.data.items():
for node, qs in self.__last_trace.data.items():
if qs[0] is None or qs[1] is None:
diff.append(node)
if diff:
......@@ -1333,8 +1339,8 @@ class Model(object):
"node q values have been set! The following "
"nodes have not been set: {}".format(diff))
self.__last_trace = trace_out.data
return trace_out
#self.__last_trace = trace_out.data
return self.__last_trace
def beam_trace_path(self, q_in, from_node, to_node, via_node=None, direction=None):
"""Finds the nodes and components along the path `from_node`
......
from .tracer import Tracer
from .tracer import TracerOut
from .tracer import TraceType
"""
Tools for tracing Gaussian beams through a Finesse :class:`.Model`
instance.
"""
from collections import OrderedDict
import numpy as np
from ..utilities import SI, SI_VALUE, SI_LABEL
class BeamTrace(object):
"""A class for storing and displaying traced beam data from
a Model. Calling :meth:`.Model.beam_trace_path` produces an
object of this type, containing the data for the trace performed.
"""
def __init__(self):
"""Constructs a new, empty instance of `BeamTrace`.
.. note::
Construction of a `BeamTrace` instance happens in :meth:`.Model.beam_trace_path`
where this object is also modified with the trace performed. Initialisation of
`BeamTrace` objects should not typically happen directly.
"""
self.data = OrderedDict()
self.Mabcd = np.matrix([])
def print(self):
"""Prints a tabular form of the beam trace, showing values for the
distance traversed, waist-size of the beam, beam radius, RoC, accumulated
Gouy phase and the beam parameter at each component.
"""
from tabulate import tabulate
nodes = self.data["nodes"]
comps = self.data["components"]
data = [
[
comp.name,
self.data[comp]['z'],
self.data[nodes[comps.index(comp)]]['q'].w0/1e-3,
self.data[nodes[comps.index(comp)]]['q'].w/1e-3,
self.data[nodes[comps.index(comp)]]['q'].Rc,
self.data[comp]['gouy'],
self.data[comp]['qout']
] for comp in comps if not self.data[comp]['is_space']
]
first_node = nodes[1] if nodes[0].is_input else nodes[0]
data.insert(
0,
[
first_node,
0,
self.data[first_node]['q'].w0/1e-3,
self.data[first_node]['q'].w/1e-3,
self.data[first_node]['q'].Rc,
0,
self.data[first_node]['q'],
]
)
last_space_node = nodes[-2] if nodes[-1].is_input else nodes[-1]
data.append([
last_space_node,
self.data[last_space_node]['z'],
self.data[last_space_node]['q'].w0/1e-3,
self.data[last_space_node]['q'].w/1e-3,
self.data[last_space_node]['q'].Rc,
self.data[last_space_node]['gouy'],
self.data[last_space_node]['q'],
])
print(
tabulate(data,
["Name", "z [m]", "w0 [mm]", "Beam radius [mm]", "RoC [m]", "Acc. Gouy [deg]", "q"],
tablefmt='psql'
)
)
def plot(self, filename=None, show=True, w_scale=SI.MILLI, markers=None):
"""Plots the beam size and Gouy phase as a function of distance
along the optical axis which the beam was traced through. The
positions of components along the traced path are marked on the
figure.
Parameters
----------
filename : str, optional
Optional name of file for saving the figure to.
show : bool, optional
Flag determining whether to show the figure, `True` by default.
w_scale : :class:`.SI`, optional
Distance scale for beam sizes, `SI.MILLI` by default.
markers : list, optional
Optional list of markers to place on the figure.
Returns
-------
fig, axes : :class:`matplotlib.figure.Figure`, tuple of :class:`matplotlib.axis.Axis`
A tuple of the figure and another tuple containing the two axis handles.
"""
import matplotlib.pyplot as plt
data = self.data
fig = plt.figure()
ax1 = plt.subplot(211)
ax2 = plt.subplot(212)
w_max = 0
g_max = 0
gouy = 0
for comp, from_node in zip(data['components'], data['nodes']):
gouy = data[comp]['gouy_i']
gouy_ref = data[comp]['gouy_ref']
if data[comp]['is_space']:
L = data[comp]['L']
q = data[from_node]['q']
z = data[from_node]['z']
_z = np.linspace(0, L, 1000)
g = np.rad2deg(q.gouy(_z+q.z))
# want to plot accumulated gouy phase so need to use
# a reference from where it started
_g = gouy + g - gouy_ref
w = q.beamsize(_z + q.z)/SI_VALUE[w_scale]
w_max = max(w_max, w.max())
g_max = max(g_max, _g.max())
ax1.plot(z+_z, w)
ax2.plot(z+_z, _g)
else:
z = data[comp]['z']
ax1.scatter(z, 0, marker='x', color='k')
ax1.text(z, 0, comp.name+"\n", ha="center", va='bottom', zorder=100)
ax2.scatter(z, 0, marker='x', color='k')
ax2.text(z, 0, comp.name+"\n", ha="center", va='bottom', zorder=100)
if markers is not None:
for _, z in markers:
ax1.scatter(z, 0, marker='x', color='r')
ax1.text(z, 0, _+"\n", ha="center", va='bottom', zorder=100)
ax2.scatter(z, 0, marker='x', color='r')
ax2.text(z, 0, _+"\n", ha="center", va='bottom', zorder=100)
ax1.grid(True, zorder=-10)
ax1.set_xlim(0, None)
if w_scale is None:
ax1.set_ylabel("Beam size [m]")
else:
ax1.set_ylabel("Beam size [{}m]".format(SI_LABEL[w_scale]))
ax1.set_xlabel("Distance [m]")
ax1.set_ylim(0, w_max)
ax2.set_xlim(0, None)
ax2.grid(True, zorder=-10)
ax2.set_ylabel("Gouy phase [deg]")
ax2.set_xlabel("Distance [m]")
ax2.set_ylim(0, g_max)
plt.tight_layout()
if filename is not None: plt.savefig(filename)
if show: plt.show()
return z + _z, w, _g
"""
Tools for tracing Gaussian beams through a Finesse :class:`.Model`
instance.
"""
import enum
from collections import OrderedDict
from functools import reduce
import networkx as nx
import numpy as np
from .components import Beamsplitter, Space
from .components.general import InteractionType, Surface
from .components.node import QSetMode
from .exceptions import NodeException, TotalReflectionError
from .paths import OpticalPath
from .utilities import SI, SI_VALUE, SI_LABEL, apply_ABCD, attached_refractive_indices
from ..components import Beamsplitter, Space
from ..components.general import InteractionType, Surface
from ..components.node import QSetMode
from ..exceptions import NodeException, TotalReflectionError
from ..gaussian import BeamParam
from ..paths import OpticalPath
from ..utilities import apply_ABCD, attached_refractive_indices
from .beamtrace import BeamTrace
class TraceType(enum.Enum):
"""Enum determining the type of tracing algorithm
......@@ -84,6 +84,9 @@ class TracerOut(object):
qys[node] = qy
return qys
def reset(self, node):
self.__data[node] = [BeamParam(), BeamParam(), QSetMode.NOTSET]
class Tracer(object):
"""
An object handling any type of tracing to be performed on a model
......@@ -457,19 +460,18 @@ class Tracer(object):
curr_node, next_node, to_comp
)
if (not next_node in data.keys() or
data[next_node][0] is None or data[next_node][1] is None):
if (data[next_node][0].is_null or data[next_node][1].is_null):
qx = apply_ABCD(Mabcd_x, data[curr_node][0], curr_node_nr, next_node_nr)
qy = apply_ABCD(Mabcd_y, data[curr_node][1], curr_node_nr, next_node_nr)
data[next_node] = (qx, qy, QSetMode.USER_TRACED if from_user else QSetMode.AUTO)
data[next_node][0].q = qx.q
data[next_node][1].q = qy.q
data[next_node][2] = QSetMode.USER_TRACED if from_user else QSetMode.AUTO
# set the opposite direction node to -q.conjugate
if set_symmetric:
data[next_node.opposite] = (
-qx.conjugate(),
-qy.conjugate(),
QSetMode.USER_TRACED if from_user else QSetMode.AUTO
)
data[next_node.opposite][0].q = -qx.conjugate().q
data[next_node.opposite][1].q = -qy.conjugate().q
data[next_node.opposite][2] = QSetMode.USER_TRACED if from_user else QSetMode.AUTO
def __run_utiltrace(self, path, q_in, datax, datay):
"""Runs a utility beam trace using a user-specified input
......@@ -588,158 +590,4 @@ class Tracer(object):
if datax is not None and datay is not None:
return trace_x, trace_y
if datax is not None: return trace_x
return trace_y
class BeamTrace(object):
"""A class for storing and displaying traced beam data from
a Model. Calling :meth:`.Model.beam_trace_path` produces an
object of this type, containing the data for the trace performed.
"""
def __init__(self):
"""Constructs a new, empty instance of `BeamTrace`.
.. note::
Construction of a `BeamTrace` instance happens in :meth:`.Model.beam_trace_path`
where this object is also modified with the trace performed. Initialisation of
`BeamTrace` objects should not typically happen directly.
"""
self.data = OrderedDict()
self.Mabcd = np.matrix([])
def print(self):
"""Prints a tabular form of the beam trace, showing values for the
distance traversed, waist-size of the beam, beam radius, RoC, accumulated
Gouy phase and the beam parameter at each component.
"""
from tabulate import tabulate
nodes = self.data["nodes"]
comps = self.data["components"]
data = [
[
comp.name,
self.data[comp]['z'],
self.data[nodes[comps.index(comp)]]['q'].w0/1e-3,
self.data[nodes[comps.index(comp)]]['q'].w/1e-3,
self.data[nodes[comps.index(comp)]]['q'].Rc,
self.data[comp]['gouy'],
self.data[comp]['qout']
] for comp in comps if not self.data[comp]['is_space']
]
first_node = nodes[1] if nodes[0].is_input else nodes[0]
data.insert(
0,
[
first_node,
0,
self.data[first_node]['q'].w0/1e-3,
self.data[first_node]['q'].w/1e-3,
self.data[first_node]['q'].Rc,
0,
self.data[first_node]['q'],
]
)
last_space_node = nodes[-2] if nodes[-1].is_input else nodes[-1]
data.append([
last_space_node,
self.data[last_space_node]['z'],
self.data[last_space_node]['q'].w0/1e-3,
self.data[last_space_node]['q'].w/1e-3,
self.data[last_space_node]['q'].Rc,
self.data[last_space_node]['gouy'],
self.data[last_space_node]['q'],
])