Commit 7de3670d authored by Daniel Brown's avatar Daniel Brown
Browse files

Added in the beginnings of the aplus model and code tools. This has the BHD LO...

Added in the beginnings of the aplus model and code tools. This has the BHD LO pick off from PR2 added but no filter cavity yet. Mode matching is way off. The ifo classes also had to be changed to accomodate Outputs being able to use two node homodyne detectors, before hand it was just single node detectors. Multiple nodes are specified with a list for node names now. ALternate node names are set using the new keyword arg
parent eb03fce5
......@@ -22,6 +22,7 @@ import six
########################
# Global helper functions
isContainer = lambda c: (not isinstance(c, six.string_types)) and hasattr(c, "__iter__")
is_container = lambda c: (not isinstance(c, six.string_types)) and hasattr(c, "__iter__")
########################
import imp
......
......@@ -1452,21 +1452,65 @@ class Output(object):
You can add many detectors at a given port, which readout different quadratures or types of transfer functions
or signals.
"""
def __init__(self, IFO, name, nodeNames, f_property_name=None, phase=0, block=None):
def __init__(self, IFO, name, nodeNames, f_property_name=None, phase=0, block=None,
alternate_node_name=None):
"""
An Output object is a generic description of an output port of an interferometer.
At an output port, such as the AS port, we can do DC or RF demodulations. Each
output can then compute the transfer functions or or photodiode error signals,
or optical amplitudes.
Parameters
----------
IFO : pykat.ifo.IFO
Interferometer class to associate with
name : str
Name of the output
nodeNames : list(str)
List of nodes that this output uses, one for normal detectors, multiple for
balanced homodyne like outputs.
f_property_name : str
Name of the frequency attribute in IFO. This should be a name of a frequency
attribute in the IFO class, this is used when computing what frequency to
demodulate at. None means a DC output
phase : foat
Demodulation phase
block : str
Which block in the kat object to add the commands to
"""
if f_property_name is not None and not isinstance(f_property_name, six.string_types):
raise Exception("f_property_name should be a property name of a frequency in class {}".format(IFO.__class__))
self.__IFO = IFO
self.name = name
self.nodeNames = make_list_copy(nodeNames)
self.__nodeNames = make_list_copy(nodeNames)
if alternate_node_name is None:
self.__alternate_node_name = None
else:
self.__alternate_node_name = make_list_copy(alternate_node_name)
self.__f_property_name = f_property_name
self.phase = phase # demodulation phase for I quadrature, float
self._block = block
if pykat.is_container(nodeNames):
self.num_nodes = len(nodeNames)
self.__nodeNames = [self.__nodeNames]
else:
self.num_nodes = 1
self.__nodeNames = [self.__nodeNames]
if self.num_nodes > 1 and self.__f_property_name is not None:
raise Exception("Can't use multiple node outputs and RF calculations")
@property
def kat(self):
"""For referencing the kat object this DOF is associated with"""
"""For referencing the actual kat object this output is associated with"""
return self.__IFO.kat
@property
def f(self):
......@@ -1475,20 +1519,23 @@ class Output(object):
return getattr(self.__IFO, self.__f_property_name)
def check_nodeName(self):
self.nodeName = None
for node in self.nodeNames:
_node = node
if _node[-1] == "*":
_node = _node[:-1]
if _node in self.__IFO.kat.nodes:
self.nodeName=node
if self.__alternate_node_name is not None:
_nodes = zip(self.__nodeNames, self.__alternate_node_name)
else:
_nodes = self.__nodeNames
for nodes in _nodes:
if all([_.strip('* ') in self.__IFO.kat.nodes for _ in nodes]):
self.nodeName = nodes
break
if self.nodeName==None:
raise pkex.BasePyKatException("port {}: cannot find any of these nodes: '{}'".format(self.name,self.nodeNames))
raise pkex.BasePyKatException("port {}: cannot find any of these nodes: '{}'".format(self.name, self.__nodeNames))
def get_amplitude_name(self, f, n=None, m=None):
''''
......@@ -1513,7 +1560,11 @@ class Output(object):
rtn = "{}_{}_{}{}_ad".format(self.name,fstr,n,m)
return rtn
def get_amplitude_cmds(self, f, n=None, m=None):
if self.num_nodes != 1:
raise Exception("Amplitude output doesn't work for multiple node outputs")
rtn = []
self.check_nodeName()
......@@ -1521,11 +1572,12 @@ class Output(object):
name = self.get_amplitude_name(f, n=n, m=m)
if n==None and m==None:
rtn.append("ad {} {} {}".format(name, f, self.nodeName))
rtn.append("ad {} {} {}".format(name, f, self.nodeName[0]))
else:
rtn.append("ad {} {} {} {} {}".format(name, f, n, m, self.nodeName))
rtn.append("ad {} {} {} {} {}".format(name, f, n, m, self.nodeName[0]))
return rtn
def _pdtype(self, name, sigtype):
rtn = []
......@@ -1537,6 +1589,7 @@ class Output(object):
return rtn
def add_signal(self, quad=None, sigtype=None):
"""
Adds a photodiode detector to the kat object at this port. Must
......@@ -1561,6 +1614,7 @@ class Output(object):
self.__IFO.kat.parse(cmds, addToBlock=self._block)
return self.get_signal_name(quad, sigtype)
def get_signal_name(self, quad="I", sigtype='z'):
name = self.name
......@@ -1575,6 +1629,7 @@ class Output(object):
return name
def get_signal_cmds(self, dof=None, **kwargs):
"""
Returns the Finesse commands for a detector added to this port's location.
......@@ -1620,19 +1675,20 @@ class Output(object):
name = self.get_signal_name(quad=quad, sigtype=sigtype)
if self.f==None:
rtn.append("pd {} {}".format(name, self.nodeName))
rtn.append("pd {} {}".format(name, self.nodeName[0]))
else:
if quad !="I" and quad != "Q":
raise pkex.BasePyKatException("quadrature must be 'I' or 'Q'")
phase = self.IQ_phase(quad, self.phase)
rtn.append("pd1 {} {} {} {}".format(name, self.f, phase, self.nodeName))
rtn.append("pd1 {} {} {} {}".format(name, self.f, phase, self.nodeName[0]))
rtn.extend(self._pdtype(name, sigtype))
return rtn
def IQ_phase(self, quad, phase):
if phase is None:
raise pkex.BasePyKatException("Phase cannot be None")
......@@ -1645,11 +1701,13 @@ class Output(object):
return phase
def add_transfer(self, quad=None, sigtype="z"):
cmds = self.get_transfer_cmds(quad=quad, sigtype=sigtype)
self.__IFO.kat.parse(cmds, addToBlock=self._block)
return self.get_transfer_name(quad, sigtype)
def get_transfer_name(self, quad="I", sigtype="z"):
name = self.name
......@@ -1663,32 +1721,49 @@ class Output(object):
return name + "_TF"
def get_transfer_cmds(self, quad="I", sigtype="z", phase2=None):
self.check_nodeName()
name = self.get_transfer_name(quad=quad, sigtype=sigtype)
rtn = []
if self.num_nodes == 2:
cmd = 'hd'
elif self.num_nodes == 1:
cmd = 'pd1'
else:
raise Exception("Unexpected number of nodes")
if self.f is None:
if phase2 is None:
rtn.append("pd1 {} {} {}".format(name, "$fs", self.nodeName))
if cmd == 'pd1':
rtn.append("{} {} {} {}".format(cmd, name, "$fs", *self.nodeName))
elif cmd == 'hd':
rtn.append("{} {} 180 {} {}".format(cmd, name, *self.nodeName))
else:
raise Expection("Unexpected command")
else:
rtn.append("pd1 {} {} {} {}".format(name, "$fs", phase2, self.nodeName))
if self.num_nodes != 1: raise Exception("Can't do RF homodyne output")
rtn.append("pd1 {} {} {} {}".format(name, "$fs", phase2, self.nodeName[0]))
else:
if self.num_nodes != 1: raise Exception("Can't do RF homodyne output")
if quad not in ("I", "Q"):
raise pkex.BasePyKatException("quadrature must be 'I' or 'Q'")
phase = self.IQ_phase(quad, self.phase)
if phase2 is None:
rtn.append("pd2 {} {} {} {} {}".format(name, self.f , phase, "$fs", self.nodeName))
rtn.append("pd2 {} {} {} {} {}".format(name, self.f , phase, "$fs", self.nodeName[0]))
else:
rtn.append("pd2 {} {} {} {} {} {}".format(name, self.f, phase, "$fs", phase2, self.nodeName))
rtn.append("pd2 {} {} {} {} {} {}".format(name, self.f, phase, "$fs", phase2, self.nodeName[0]))
rtn.extend(self._pdtype(name, sigtype))
return rtn
def get_diff_cmds(self, quad='I', sigtype='z'):
'''
......@@ -1700,14 +1775,14 @@ class Output(object):
name = self.get_transfer_name(quad=quad, sigtype=sigtype)
rtn = []
if self.f is None:
rtn.append("pd {} {}".format(name, self.nodeName))
rtn.append("pd {} {}".format(name, self.nodeName[0]))
else:
if quad not in ("I", "Q"):
raise pkex.BasePyKatException("quadrature must be 'I' or 'Q'")
phase = self.IQ_phase(quad, self.phase)
rtn.append("pd1 {} {} {} {}".format(name, self.f, phase, self.nodeName))
rtn.append("pd1 {} {} {} {}".format(name, self.f, phase, self.nodeName[0]))
rtn.extend(self._pdtype(name, sigtype))
......
from .aplus import *
\ No newline at end of file
This diff is collapsed.
import matplotlib.pyplot as plt
from . import assert_aplus_ifo_kat
from .. import scan_optics_string
from pykat.ifo.plot import *
import pykat.ifo
import numpy as np
import six
def f1_PRC_resonance(_kat, ax=None, show=True):
"""
Plot the sideband amplitudes for modulation frequecy
f1 (~ 9MHz) in the PRC, to check the resonance
condition.
"""
assert_aplus_ifo_kat(_kat)
kat = _kat.deepcopy()
# Don't need locks for this plot so remove if present
kat.removeBlock('locks', False)
if ax is None:
fig = plt.figure()
ax = fig.add_subplot(1,1,1)
startf = kat.IFO.f1 - 4000.0
stopf = kat.IFO.f1 + 4000.0
if _kat.maxtem == "off":
nmStr = None
else:
nmStr = "0 0"
code = """
ad f1p {0} {1} nPRM2
ad f1m {0} -{1} nPRM2
xaxis mod1 f lin {2} {3} 200
put f1p f $x1
put f1m f $mx1
""".format(nmStr, kat.IFO.f1, startf, stopf)
kat.parse(code)
out = kat.run()
ax.plot(out.x-kat.IFO.f1,np.abs(out["f1p"]), label=" f1")
ax.plot(out.x-kat.IFO.f1,np.abs(out["f1m"]), label="-f1", ls="--")
ax.set_xlim([np.min(out.x-kat.IFO.f1), np.max(out.x-kat.IFO.f1)])
ax.set_xlabel("delta_f1 [Hz]")
ax.set_ylabel('sqrt(W) ')
ax.grid(True)
ax.legend()
ax.figure.set_tight_layout(True)
if show: plt.show()
def pretuning_powers(self, _kat, xlimits=[-10,10]):
assert_aplus_ifo_kat(_kat)
kat = _kat.deepcopy()
kat.verbose = False
kat.noxaxis = True
dofs = [self.preARMX, self.preARMY, self.preMICH, self.prePRCL]
idx=1
fig = plt.figure()
for d in dofs:
ax = fig.add_subplot(2,2,idx, label=d.name)
idx+=1
out = scan_DOF(kat, d, xlimits = np.multiply(d.scale, xlimits), relative = True)
ax.semilogy(out.x,out[d.signal_name(kat)])
ax.set_xlim([np.min(out.x), np.max(out.x)])
ax.set_xlabel("phi [deg] {}".format(d.optics[0]))
ax.set_ylabel('{} [W] '.format(d.signal_name(kat)))
ax.grid(True)
plt.tight_layout()
plt.show(block=0)
def error_signals(_kat, xlimits=[-1,1], DOFs=None, plotDOFs=None,
replaceDOFSignals=False, block=True, fig=None, label=None, steps=100):
"""
Displays error signals for a given kat file. Can also be used to plot multiple
DOF's error signals against each other for visualising any cross coupling.
_kat: LIGO-like kat object.
xlimits: Range of DOF to plot in degrees
DOFs: list, DOF names to compute. Default: DARM, CARM, PRCL, SRCL, MICH
plotDOFs: list, DOF names to plot against each DOF. If None the same DOF as in DOFs is plotted.
block: Boolean, for plot blocking terminal or not if being shown
replaceDOFSignals: Bool, replaces already present signals for any DOF if already defined in kat. Regardless of this value, it will add default signals if none found.
fig: figure, uses predefined figure, when defined it won't be shown automatically
label: string, if no plotDOFs is defined this legend is shown
Example:
import pykat
from pykat.gw_detectors import ifo
ligo = ifo.aLIGO()
# Plot default
ligo.plot_error_signals(ligo.kat, block=True)
# plot MICH and CARM against themselves
ligo.plot_error_signals(ligo.kat, DOFs=["MICH", "CARM"], block=True)
# plot DARM and CARM against MICH
ligo.plot_error_signals(ligo.kat, DOFs=["MICH"], plotDOFs=["DARM", "CARM"], block=True)
"""
kat = _kat.deepcopy()
kat.verbose = False
kat.noxaxis = True
kat.removeBlock("locks", False)
if DOFs is None:
dofs = [kat.IFO.DARM, kat.IFO.CARM, kat.IFO.PRCL, kat.IFO.SRCL, kat.IFO.MICH]
else:
dofs = kat.IFO.strsToDOFs(DOFs)
# add in signals for those DOF to plot
for _ in dofs:
if replaceDOFSignals and not hasattr(kat, _.signal_name()):
kat.parse(_.signal())
toShow = None
if plotDOFs is not None:
toShow = self._strToDOFs(plotDOFs)
# Check if other DOF signals we need to include for plotting
for _ in toShow:
if not (not replaceDOFSignals and hasattr(kat, _.signal_name())):
kat.parse(_.signal())
if fig is not None:
_fig = fig
else:
_fig = plt.figure()
nrows = 2
ncols = 3
if DOFs is not None:
n = len(DOFs)
if n < 3:
nrows = 1
ncols = n
ax_list = _fig.axes
for d, idx in zip(dofs, range(1, len(dofs)+1)):
### TODO this seems fragile, needs more testing (adf 25032018)
if idx < len(ax_list)+1:
ax = ax_list[idx-1]
#print("reuse axis")
else:
ax = _fig.add_subplot(nrows, ncols, idx)
#print("new axis")
kat.removeBlock("SCAN", False)
scan_cmd = scan_optics_string(d.optics, d.factors, "scan", linlog="lin",
xlimits=np.multiply(d.scale, xlimits), steps=steps,
axis=1, relative=True)
kat.parse(scan_cmd, addToBlock="SCAN")
out = kat.run(cmd_args=['-cr=on'])
DC_Offset = None
# Get a lock offset if used
if (d.name + '_lock') in _kat.commands:
DC_Offset = _kat.commands[d.name + '_lock'].offset.value
if DC_Offset is None:
DC_Offset = 0
else:
DC_Offset = float(DC_Offset)
if toShow is None:
ax.plot(out.x, out[d.signal_name()] + DC_Offset, label=label)
else:
for _ in toShow:
ax.plot(out.x, out[_.signal_name()] + DC_Offset, label=label)
ax.set_xlim([np.min(out.x), np.max(out.x)])
ax.set_xlabel("{} [deg]".format(d.name))
if plotDOFs is None:
ax.set_ylabel('{} [W] '.format(d.port.name))
else:
ax.set_ylabel('Error signal [W]')
ax.grid(True)
if toShow is not None or label is not None:
plt.legend(loc=0)
plt.tight_layout()
if fig is None:
plt.show(block=block)
return(fig)
def amps_vs_dof(kat, DoF, f, n=None, m=None, xaxis = [-10,10,100], noplot=False):
'''
Plotting amplitude vs tuning for one LSC DoF.
DoF - Degree of freedom to sweep.
f - Frequency component relative to default
n, m - Mode numbers
xaxis - range to plot over [min, max, steps]
'''
_kat = kat.deepcopy()
if isinstance(DoF, six.string_types):
DoF = _kat.IFO.DOFs[DoF]
# Adding detectors
code = ""
names = []
for o in _kat.IFO.CAV_POWs:
code += "{}\n".format(o.get_amplitude_cmds(f,n,m)[0])
names.append(o.get_amplitude_name(f,n,m))
# Adding simulation instructions
code += pykat.ifo.scan_DOF_cmds(DoF, xlimits=[xaxis[0], xaxis[1]], steps=xaxis[2], relative=True)
_kat.parse(code)
out = _kat.run()
if noplot:
rtn = {'x': out.x}
for n in names:
rtn[n] = out[n]
return rtn
else:
FS = 13
LS = 12
fig = plt.figure()
ax = fig.add_subplot(111)
for n in names:
ax.semilogy(out.x, out[n], label = n)
ax.set_ylabel('$\mathrm{Amplitude}\ [\sqrt{\mathrm{W}}]$', fontsize=FS)
ax.set_xlabel('$\mathrm{{ {}\ tuning\ [deg] }} $'.format(DoF.name), fontsize=FS)
ax.set_xlim(out.x.min(), out.x.max())
ax.grid()
ax.legend(loc=3, fontsize=LS)
plt.show(fig)
return fig, ax
def amps_vs_dofs(kat, f, n=None, m=None, xaxis = [-1,1,100]):
'''
Plotting amplitude vs tuning for all LSC DoFs.
f - Frequency component relative to default
n, m - Mode numbers
xaxis - range to plot over [min, max, steps]
'''
dic = {}
for d in kat.IFO.LSC_DOFs:
xax = [d.scale*xaxis[0], d.scale*xaxis[1], xaxis[2]]
dic[d.name] = amps_vs_dof(kat, d, f, n=n, m=m, xaxis = xax, noplot=True)
N = len(dic)
FS = 13
LS = 11
fig = plt.figure(figsize=(17,8))
axs = []
for k in range(N):
axs.append(fig.add_subplot(2,3,k+1), label = 'axis{}'.format(k))
for ax, v in zip(axs,dic.keys()):
for n in dic[v].keys():
if n != 'x':
ax.semilogy(dic[v]['x'], dic[v][n], label = n)
ax.set_ylabel('$\mathrm{Amplitude}\ [\sqrt{\mathrm{W}}]$', fontsize=FS)
ax.set_xlabel('$\mathrm{{ {}\ tuning\ [deg] }} $'.format(v), fontsize=FS)
ax.set_xlim(dic[v]['x'].min(), dic[v]['x'].max())
ax.grid()
ax.legend(loc=3, fontsize=LS)
plt.show(fig)
def pows_vs_dof(kat, DoF, xaxis = [-10,10,100], noplot=False):
'''
Plotting amplitude vs tuning for one LSC DoF.
DoF - Degree of freedom to sweep.
xaxis - range to plot over [min, max, steps]
'''
_kat = kat.deepcopy()
if isinstance(DoF, six.string_types):
DoF = _kat.IFO.DOFs[DoF]
# Adding detectors
code = ""
names = []
for o in _kat.IFO.CAV_POWs:
code += "{}\n".format(o.get_signal_cmds()[0])
names.append(o.get_signal_name())
# Adding simulation instructions
code += pykat.ifo.scan_DOF_cmds(DoF, xlimits=[xaxis[0], xaxis[1]], steps=xaxis[2], relative=True)
_kat.parse(code)
out = _kat.run()
if noplot:
rtn = {'x': out.x}
for n in names:
rtn[n] = out[n]
return rtn
else:
FS = 13
LS = 12
fig = plt.figure()
ax = fig.add_subplot(111)
for n in names:
ax.semilogy(out.x, out[n], label = n)
ax.set_ylabel('$\mathrm{Power\ [W]}$', fontsize=FS)
ax.set_xlabel('$\mathrm{{ {}\ tuning\ [deg] }} $'.format(DoF.name), fontsize=FS)
ax.set_xlim(out.x.min(), out.x.max())
ax.grid()
ax.legend(loc=3, fontsize=LS)
plt.show(fig)
return fig, ax