diff --git a/.gitlab-ci.yml b/.gitlab-ci.yml
index c33a8df2d35546de3f49d7c01bc6a7fe549a70a4..13e0e1a81c18f75aaf3592d1764ce81d28ac18f3 100644
--- a/.gitlab-ci.yml
+++ b/.gitlab-ci.yml
@@ -15,10 +15,10 @@ test:
- export PYTHONPATH=inspiral_range
- export MPLBACKEND=agg
- python3 -m gwinc aLIGO -s aLIGO.png
- - python3 -m gwinc A+ -s A+.png
+ - python3 -m gwinc Aplus -s Aplus.png
- python3 -m gwinc Voyager -s Voyager.png
- python3 -m gwinc.test aLIGO -t 10e-6 -k Seismic -k "Substrate Thermo-Elastic" -p -s aLIGO_test.png
- - python3 -m gwinc.test -t 100e-6 -k Seismic -k "Substrate Thermo-Elastic" A+ -p -s A+_test.png
+ - python3 -m gwinc.test -t 100e-6 -k Seismic -k "Substrate Thermo-Elastic" Aplus -p -s Aplus_test.png
after_script:
- rm gitID.txt
cache:
@@ -29,10 +29,10 @@ test:
expire_in: 4w
paths:
- aLIGO.png
- - A+.png
+ - Aplus.png
- Voyager.png
- aLIGO_test.png
- - A+_test.png
+ - Aplus_test.png
pages:
stage: deploy
@@ -41,10 +41,10 @@ pages:
script:
- mkdir public
- mv aLIGO.png public/
- - mv A+.png public/
+ - mv Aplus.png public/
- mv Voyager.png public/
- mv aLIGO_test.png public/ || true
- - mv A+_test.png public/ || true
+ - mv Aplus_test.png public/ || true
- apt-get install -y -qq python3-pip python3-dev make
- pip3 install sphinx sphinx-rtd-theme
- cd docs
diff --git a/HDF5_SCHEMATA b/HDF5_SCHEMATA
new file mode 100644
index 0000000000000000000000000000000000000000..b2ef34183ef687f19dc6bcb6c7f535c18cf3f531
--- /dev/null
+++ b/HDF5_SCHEMATA
@@ -0,0 +1,104 @@
+# HDF5 Schema for GWINC noise trace storage
+
+
+This file describes a schemata for HDF5 storage of noise trace data
+and plot styling GWINC noise budgets.
+
+HDF5 is a heirarchical, structured data storage format [0]. Content
+is organized into a heirarchical folder-like structure, with two
+types of named objects:
+
+ * groups: holder of other objects (like folders)
+ * datasets: holder of data arrays (like files)
+
+Objects can also have attributes as (almost) arbitrary key:value
+pairs.
+
+Bindings are available for most platforms including Python [1] and
+Matlab.
+
+[0] https://en.wikipedia.org/wiki/Hierarchical_Data_Format
+[1] http://www.h5py.org/
+
+
+## version history
+
+v1
+- first versioned schema release
+
+
+## schema
+
+The following describes the noise budget schema. Specific strings are
+enclosed in single quotes (''), and variables are described in
+brackets (<>). Group objects are indicated by a closing '/'
+separator, data set are indicated by a closing ':' followed by a
+specification of their length and type (e.g. "(N),float"), and
+attributes are specified in the .attrs[] dictionary format. Optional
+elements are enclosed in parentheses.
+
+A single trace is a length N array (with optional plot style specified
+in attributes:
+```
+/<trace>: (N),float
+(/<trace>.attrs['label'] = <label>)
+(/<trace>.attrs['color] = <color>)
+...
+```
+
+A budget item, i.e. a collection of noises is structured as follows:
+```
+/<budget>/
+ /'Total': (N),float
+ /<trace_0>: (N),float
+ (/<trace_1>: (N),float)
+```
+<!-- ``` -->
+<!-- /'noises'/ -->
+<!-- /*<noise_0>: (N),float -->
+<!-- ... -->
+<!-- /\*'references'/ -->
+<!-- /*<ref_0>: (N),float -->
+<!-- ... -->
+<!-- ``` -->
+
+
+## Top-level Budget
+
+The following two root attributes expected: a string describing the schema,
+and an int schema version number:
+```
+/.attrs['SCHEMA'] = 'GWINC Noise Budget'
+/.attrs['SCHEMA_VERSION'] = 1
+```
+
+Top-level attributes are generally used for specifying overall plot styling, but the
+following root attributes are typically defined:
+```
+/.attrs['title'] = <experiment description string (e.g. 'H1 Strain Budget')>
+/.attrs['date'] = <ISO-formatted string (e.g. '2015-10-24T20:30:00.000000Z')>
+```
+
+The budget frequency array is defined in a 'Freq' dataset:
+```
+/'Freq': (N), float
+```
+
+The budget traces are defined a traces group. The overall structure
+looks something like this:
+```
+/.attrs['SCHEMA'] = 'GWINC Noise Budget'
+/.attrs['SCHEMA_VERSION'] = 1
+/.attrs['title'] = <experiment description string (e.g. 'H1 Strain Budget')>
+/.attrs['date'] = <ISO-formatted string (e.g. '2015-10-24T20:30:00.000000Z')>
+/'Freq': (N), float
+/traces/
+ /'Total': (N),float
+ /<noise_0>: (N),float
+ /<noise_1>: (N),float
+ /<noise_2>/
+ /'Total': (N),float
+ /<noise_3>: (N),float
+ /<noise_4>: (N),float
+ ...
+```
diff --git a/README.md b/README.md
index f983a18fec3cbf381b10b54e98e220a158d53b87..adf5a057689297196c1a71085532c6aa65ee69bf 100644
--- a/README.md
+++ b/README.md
@@ -16,14 +16,11 @@ description is loaded, the noise budget can be calculated and plotted:
>>> import gwinc
>>> import numpy as np
>>> freq = np.logspace(1, 3, 1000)
->>> ifo = gwinc.load_ifo('aLIGO')
+>>> Budget, ifo, freq_, plot_style = gwinc.load_ifo('aLIGO')
>>> ifo = gwinc.precompIFO(freq, ifo)
->>> noises = gwinc.noise_calc(freq, ifo)
->>> gwinc.plot_noise(noises)
-```
-Or the `gwinc` convenience function can be used to handle it all:
-```
->>> score, data, ifo = gwinc.gwinc(freq, ifo, plot=True)
+>>> traces = Budget(freq, ifo=ifo).calc_trace()
+>>> fig = gwinc.plot_noise(freq, traces, **plot_style)
+>>> fig.show()
```
@@ -32,39 +29,7 @@ Or the `gwinc` convenience function can be used to handle it all:
You can make gwinc plots directly from the command line by executing
the package directly:
```shell
-$ python3 -m gwinc -h
-usage: gwinc [-h] [--flo FLO] [--fhi FHI] [--npoints NPOINTS] [--title TITLE]
- [--matlab] [--fom FOM] [--dump | --save SAVE | --interactive]
- [IFO]
-
-Plot GWINC noise budget for specified IFO.
-
-If the inspiral_range package is installed, various figures of merit
-can be calculated for the resultant spectrum with the --fom argument,
-e.g.:
-
- gwinc --fom horizon ...
- gwinc --fom range:m1=20,m2=20 ...
-
-See documentation for inspiral_range package for details.
-
-positional arguments:
- IFO IFO name or description file path (.yaml or .mat)
-
-optional arguments:
- -h, --help show this help message and exit
- --flo FLO, -fl FLO lower frequency bound in Hz [5]
- --fhi FHI, --fh FHI upper frequency bound in Hz [6000]
- --npoints NPOINTS, -n NPOINTS
- number of frequency points [3000]
- --title TITLE, -t TITLE
- plot title
- --matlab, -m use MATLAB gwinc engine to calculate noises
- --fom FOM calculate inspiral range for resultant spectrum
- ('func[:param=val,param=val]')
- --dump, -d print IFO parameters to stdout and exit
- --save SAVE, -s SAVE save figure to file
- --interactive, -i open interactive shell when plotting
+$ python3 -m gwinc aLIGO
```
@@ -80,6 +45,51 @@ various detectors:
* [Voyager.yaml](https://git.ligo.org/gwinc/pygwinc/blob/master/gwinc/ifo/Voyager.yaml)
+
+## noise budgets
+
+
+GWINC provides an `nb` package for defining arbitrary noise budgets:
+
+```python
+import numpy as np
+from gwinc import nb
+from gwinc import noise
+
+class ExcessGas(nb.Noise):
+ """Excess gas"""
+ style = dict(
+ label='Excess Gas',
+ color='#ad900d',
+ linestyle='--',
+ )
+
+ def calc(self):
+ return noise.residualgas.gas(self.freq, self.ifo)
+
+class MeasuredNoise(nb.Noise):
+ """My measured noise"""
+ style = dict(
+ label='Measured Noise,
+ color='#838209,
+ linestyle='-',
+ )
+
+ def load(self):
+ psd, freq = np.loadtxt('/path/to/measured/psd.txt')
+ self.data = self.interpolate(f, psd)
+
+ def calc(self):
+ return self.data
+
+class MyBudget(nb.Budget):
+ noises = [
+ ExcessGas,
+ MeasuredNoise,
+ ]
+```
+
+
## comparison with MATLAB gwinc
`pygwinc` includes the ability use MATLAB gwinc directly via the
diff --git a/gwinc/__init__.py b/gwinc/__init__.py
index abf5763fd5323f41731a5125b9b33246d2ee1b4e..640ed7365f593773a5f24691b2b2f909c42c943c 100644
--- a/gwinc/__init__.py
+++ b/gwinc/__init__.py
@@ -1,6 +1,5 @@
-from .struct import load_struct
from .ifo import available_ifos, load_ifo
+from .struct import load_struct
from .precomp import precompIFO
-from .gwinc import noise_calc
-from .gwinc import gwinc
from .plot import plot_noise
+from .io import load_hdf5, save_hdf5
diff --git a/gwinc/__main__.py b/gwinc/__main__.py
index 0005c5a72f9d4390fec4af9071408284b33ba24e..b20919da799a681221cac83e851e2ac67b4aedc0 100644
--- a/gwinc/__main__.py
+++ b/gwinc/__main__.py
@@ -11,9 +11,8 @@ logging.basicConfig(format='%(message)s',
from .ifo import available_ifos, load_ifo
from .precomp import precompIFO
-from . import gwinc
from . import plot_noise
-from . import util
+from . import io
##################################################
@@ -30,10 +29,11 @@ By default a GWINC noise budget of the specified IFO will calculated,
and plotted with an interactive plotter. If the --save option is
specified the plot will be saved directly to a file (without display)
(various formats are supported, indicated by file extension). If the
-requested extension is "hdf5" then the noise traces and IFO parameters
-will be saved to an HDF5 file (without plotting). The input file
-(IFO) can be an HDF5 file saved from a previous call, in which all
-noise traces and IFO parameters will be loaded from that file.
+requested extension is 'hdf5' or 'h5' then the noise traces and IFO
+parameters will be saved to an HDF5 file (without plotting). The
+input file (IFO) can be an HDF5 file saved from a previous call, in
+which case all noise traces and IFO parameters will be loaded from
+that file.
If the inspiral_range package is installed, various figures of merit
can be calculated for the resultant spectrum with the --fom argument,
@@ -51,9 +51,10 @@ FLO = 5
FHI = 6000
NPOINTS = 3000
-parser = argparse.ArgumentParser(prog='gwinc',
- description=description,
- formatter_class=argparse.RawDescriptionHelpFormatter)
+parser = argparse.ArgumentParser(
+ prog='gwinc',
+ description=description,
+ formatter_class=argparse.RawDescriptionHelpFormatter)
parser.add_argument('--flo', '-fl', default=FLO, type=float,
help="lower frequency bound in Hz [{}]".format(FLO))
parser.add_argument('--fhi', '--fh', default=FHI, type=float,
@@ -68,7 +69,7 @@ group = parser.add_mutually_exclusive_group()
group.add_argument('--interactive', '-i', action='store_true',
help="interactive plot with interactive shell")
group.add_argument('--save', '-s',
- help="save noise (hdf5) or figure (pdf/png/svg) to file")
+ help="save budget traces (.hdf5/.h5) or plot (.pdf/.png/.svg) to file")
group.add_argument('--yaml', '-y', action='store_true',
help="print IFO as yaml to stdout and exit")
group.add_argument('--text', '-x', action='store_true',
@@ -77,8 +78,8 @@ group.add_argument('--diff', '-d', metavar='IFO',
help="show differences table between another IFO description")
group.add_argument('--no-plot', '-np', action='store_false', dest='plot',
help="supress plotting")
-parser.add_argument('IFO', default=IFO,
- help="IFO name or description file path (.yaml, .mat, .m), or full HDF5 data file")
+parser.add_argument('IFO',
+ help="IFO name, description file path (.yaml, .mat, .m), budget module (.py), or HDF5 data file (.hdf5, .h5)")
def main():
@@ -89,26 +90,34 @@ def main():
##########
# initial arg processing
- if os.path.splitext(args.IFO)[1] in ['.hdf5', '.h5']:
- title, ifo, noises = util.load_hdf5(args.IFO)
+ if os.path.splitext(os.path.basename(args.IFO))[1] in ['.hdf5', '.h5']:
+ Budget = None
+ freq, traces, attrs = io.load_hdf5(args.IFO)
+ ifo = getattr(attrs, 'IFO', None)
+ plot_style = attrs
else:
- ifo = load_ifo(args.IFO)
- noises = None
- if args.title:
- title = args.title
- else:
- title = '{} GWINC Noise Budget'.format(args.IFO)
+ Budget, ifo, freq, plot_style = load_ifo(args.IFO)
+ # FIXME: this should be done only if specified, to allow for
+ # using any FREQ specified in budget module by default
+ freq = np.logspace(np.log10(args.flo), np.log10(args.fhi), args.npoints)
+ traces = None
if args.yaml:
+ if not ifo:
+ parser.exit(2, "no IFO structure available.")
print(ifo.to_yaml(), end='')
return
if args.text:
+ if not ifo:
+ parser.exit(2, "no IFO structure available.")
print(ifo.to_txt(), end='')
return
if args.diff:
+ if not ifo:
+ parser.exit(2, "no IFO structure available.")
fmt = '{:30} {:>20} {:>20}'
- ifoo = load_ifo(args.diff)
+ Budget, ifoo, freq, plot_style = load_ifo(args.diff)
diffs = ifo.diff(ifoo)
if diffs:
print(fmt.format('', args.IFO, args.diff))
@@ -118,6 +127,14 @@ def main():
ov = repr(p[2])
print(fmt.format(k, v, ov))
return
+
+ if args.title:
+ plot_style['title'] = args.title
+ elif Budget:
+ plot_style['title'] = "GWINC Noise Budget: {}".format(Budget.name)
+ else:
+ plot_style['title'] = "GWINC Noise Budget: {}".format(args.IFO)
+
if args.plot:
if args.save:
# FIXME: this silliness seems to be the only way to have
@@ -156,18 +173,20 @@ def main():
##########
# main calculations
- if noises:
- freq = noises['Freq']
-
- else:
- logging.info("calculating noises...")
- freq = np.logspace(np.log10(args.flo), np.log10(args.fhi), args.npoints)
- score, noises, ifo = gwinc(freq, ifo)
-
- logging.info('recycling factor: {: >0.3f}'.format(ifo.gwinc.prfactor))
- logging.info('BS power: {: >0.3f} W'.format(ifo.gwinc.pbs))
- logging.info('arm finesse: {: >0.3f}'.format(ifo.gwinc.finesse))
- logging.info('arm power: {: >0.3f} kW'.format(ifo.gwinc.parm/1000))
+ if not traces:
+ if ifo:
+ logging.info("precomputing ifo...")
+ ifo = precompIFO(freq, ifo)
+ logging.info("calculating budget...")
+ budget = Budget(freq=freq, ifo=ifo)
+ budget.load()
+ budget.update()
+ traces = budget.calc_trace()
+
+ # logging.info('recycling factor: {: >0.3f}'.format(ifo.gwinc.prfactor))
+ # logging.info('BS power: {: >0.3f} W'.format(ifo.gwinc.pbs))
+ # logging.info('arm finesse: {: >0.3f}'.format(ifo.gwinc.finesse))
+ # logging.info('arm power: {: >0.3f} kW'.format(ifo.gwinc.parm/1000))
if args.fom:
logging.info("calculating inspiral {}...".format(range_func))
@@ -181,30 +200,35 @@ def main():
m2=H.params['m2'],
fom=fom,
)
- title += '\n{}'.format(fom_title)
+ plot_style['title'] += '\n{}'.format(fom_title)
##########
# output
# save noise traces to HDF5 file
if args.save and os.path.splitext(args.save)[1] in ['.hdf5', '.h5']:
- util.save_hdf5(
- title,
- ifo,
- noises,
- args.save,
+ logging.info("saving budget traces {}...".format(args.save))
+ if ifo:
+ plot_style['IFO'] = ifo.to_yaml()
+ io.save_hdf5(
+ path=args.save,
+ freq=freq,
+ traces=traces,
+ **plot_style
)
# interactive shell plotting
elif args.interactive:
ipshell = InteractiveShellEmbed(
- user_ns={'freq': freq,
- 'noises': noises,
- 'ifo': ifo,
- 'plot_noise': plot_noise,
+ user_ns={
+ 'freq': freq,
+ 'traces': traces,
+ 'ifo': ifo,
+ 'plot_style': plot_style,
+ 'plot_noise': plot_noise,
},
banner1='''
-PYGWINC interactive plotter
+GWINC interactive plotter
You may interact with plot using "plt." methods, e.g.:
@@ -212,20 +236,21 @@ You may interact with plot using "plt." methods, e.g.:
>>> plt.savefig("foo.pdf")
''')
ipshell.enable_pylab()
- ipshell.run_code("plot_noise(noises)")
- ipshell.run_code("plt.title('{}')".format(title))
+ ipshell.run_code("plot_noise(freq, traces, **plot_style)")
+ ipshell.run_code("plt.title('{}')".format(plot_style['title']))
ipshell()
# standard plotting
elif args.plot:
- logging.info("plotting...")
+ logging.info("plotting noises...")
fig = plt.figure()
ax = fig.add_subplot(1, 1, 1)
plot_noise(
- noises,
+ freq,
+ traces,
ax=ax,
+ **plot_style
)
- ax.set_title(title)
fig.tight_layout()
if args.save:
fig.savefig(
diff --git a/gwinc/ifo/Aplus/__init__.py b/gwinc/ifo/Aplus/__init__.py
new file mode 100644
index 0000000000000000000000000000000000000000..81bb8b17cf0bbb4717ce34505d13ec5b248bc638
--- /dev/null
+++ b/gwinc/ifo/Aplus/__init__.py
@@ -0,0 +1,18 @@
+from gwinc.ifo.noises import *
+
+
+class Aplus(nb.Budget):
+
+ name = 'A+'
+
+ noises = [
+ QuantumVacuum,
+ Seismic,
+ Newtonian,
+ SuspensionThermal,
+ CoatingBrownian,
+ CoatingThermoOptic,
+ SubstrateBrownian,
+ SubstrateThermoElastic,
+ ExcessGas,
+ ]
diff --git a/gwinc/ifo/A+.yaml b/gwinc/ifo/Aplus/ifo.yaml
similarity index 100%
rename from gwinc/ifo/A+.yaml
rename to gwinc/ifo/Aplus/ifo.yaml
diff --git a/gwinc/ifo/CE/__init__.py b/gwinc/ifo/CE/__init__.py
new file mode 100644
index 0000000000000000000000000000000000000000..ca22f9e70a78ccac2cbafbe4edb3903fb0ba3487
--- /dev/null
+++ b/gwinc/ifo/CE/__init__.py
@@ -0,0 +1,20 @@
+from gwinc.ifo.noises import *
+
+
+class CE(nb.Budget):
+
+ name = 'Cosmic Explorer'
+
+ noises = [
+ QuantumVacuum,
+ Seismic,
+ Newtonian,
+ SuspensionThermal,
+ CoatingBrownian,
+ CoatingThermoOptic,
+ ITMThermoRefractive,
+ ITMCarrierDensity,
+ SubstrateBrownian,
+ SubstrateThermoElastic,
+ ExcessGas,
+ ]
diff --git a/gwinc/ifo/CE2.yaml b/gwinc/ifo/CE/ifo.yaml
similarity index 100%
rename from gwinc/ifo/CE2.yaml
rename to gwinc/ifo/CE/ifo.yaml
diff --git a/gwinc/ifo/Voyager/__init__.py b/gwinc/ifo/Voyager/__init__.py
new file mode 100644
index 0000000000000000000000000000000000000000..0839e3c2113765e734f5d1620dacf81832701d80
--- /dev/null
+++ b/gwinc/ifo/Voyager/__init__.py
@@ -0,0 +1,20 @@
+from gwinc.ifo.noises import *
+
+
+class Voyager(nb.Budget):
+
+ name = 'Voyager'
+
+ noises = [
+ QuantumVacuum,
+ Seismic,
+ Newtonian,
+ SuspensionThermal,
+ CoatingBrownian,
+ CoatingThermoOptic,
+ ITMThermoRefractive,
+ ITMCarrierDensity,
+ SubstrateBrownian,
+ SubstrateThermoElastic,
+ ExcessGas,
+ ]
diff --git a/gwinc/ifo/Voyager.yaml b/gwinc/ifo/Voyager/ifo.yaml
similarity index 100%
rename from gwinc/ifo/Voyager.yaml
rename to gwinc/ifo/Voyager/ifo.yaml
diff --git a/gwinc/ifo/__init__.py b/gwinc/ifo/__init__.py
index 7ea1c794bf7ae0d296e7ab38edc89d930fe6b6b0..a51ff4a2c365f6edc51021dfd932534666c1081a 100644
--- a/gwinc/ifo/__init__.py
+++ b/gwinc/ifo/__init__.py
@@ -1,48 +1,76 @@
import os
-import fnmatch
+import logging
-from ..struct import Struct
-from ..gwinc_matlab import Matlab
+from ..struct import load_struct, STRUCT_EXT
+from ..util import load_module
+
+
+PLOT_STYLE = dict(
+ ylabel=u"Strain [1/\u221AHz]",
+)
def available_ifos():
- """List available included IFO files"""
+ """List available included pre-defined IFOs
+
+ """
ifos = []
- for f in os.listdir(os.path.dirname(__file__)):
- if fnmatch.fnmatch(f, '*.yaml'):
- ifos.append(f.split('.')[0])
- return ifos
+ root = os.path.dirname(__file__)
+ for f in os.listdir(root):
+ if os.path.isdir(os.path.join(root, f)) and f[0] != '_':
+ ifos.append(f)
+ return sorted(ifos)
def load_ifo(name_or_path):
- """Load IFO by name or from file.
+ """Load GWINC IFO Budget by name or from file.
- Named IFOs should correspond to the basename of .yaml IFO
- definition files included with pygwinc (see available_ifos()
- above).
+ Named IFOs should correspond to one of the IFOs available in the
+ gwinc package (see gwinc.available_ifos()). If a path is provided
+ it should either be a budget package (directory) or module (ending
+ in .py), or an IFO struct (see gwinc.load_struct()). In the
+ latter case the base aLIGO budget definition will be used.
- When specifying by path files may be either .yaml, .mat or .m.
- For .m files, the file is expected to include either an object or
- function that corresponds to the basename of the file. The MATLAB
- engine will be invoked to execute the .m code and extract the
- resultant IFO data.
+ Returns primary Budget class, ifo structure, frequency array, and
+ plot style dictionary, with the last three being None if they are
+ not defined in the budget.
"""
+ ifo = None
+
if os.path.exists(name_or_path):
- path = name_or_path
- else:
- path = os.path.join(os.path.dirname(__file__),
- name_or_path+'.yaml')
+ path = name_or_path.rstrip('/')
+ bname, ext = os.path.splitext(os.path.basename(path))
- (root, ext) = os.path.splitext(path)
+ if ext in STRUCT_EXT:
+ logging.info("loading struct {}...".format(path))
+ ifo = load_struct(path)
+ bname = 'aLIGO'
+ modname = 'gwinc.ifo.aLIGO'
+ logging.info("loading budget {}...".format(modname))
- if ext == '.m':
- matlab = Matlab()
- matlab.addpath(os.path.dirname(path))
- func_name = os.path.basename(root)
- matlab.eval("ifo = {};".format(func_name), nargout=0)
- ifo = matlab.extract('ifo')
- return Struct.from_matstruct(ifo)
+ else:
+ modname = path
+ logging.info("loading module path {}...".format(modname))
else:
- return Struct.from_file(path)
+ if name_or_path not in available_ifos():
+ raise RuntimeError("Unknonw IFO '{}' (available IFOs: {}).".format(
+ name_or_path,
+ available_ifos(),
+ ))
+ bname = name_or_path
+ modname = 'gwinc.ifo.'+name_or_path
+ logging.info("loading module {}...".format(modname))
+
+ mod, modpath = load_module(modname)
+
+ Budget = getattr(mod, bname)
+ ifo = getattr(mod, 'IFO', ifo)
+ ifopath = os.path.join(modpath, 'ifo.yaml')
+ if not ifo and ifopath:
+ ifo = load_struct(ifopath)
+ freq = getattr(mod, 'FREQ', None)
+ plot_style = getattr(mod, 'PLOT_STYLE', PLOT_STYLE)
+
+ return Budget, ifo, freq, plot_style
diff --git a/gwinc/ifo/aLIGO/__init__.py b/gwinc/ifo/aLIGO/__init__.py
new file mode 100644
index 0000000000000000000000000000000000000000..1f1bf4efffb7c6c0c64b3b1ad14c86b5135733a5
--- /dev/null
+++ b/gwinc/ifo/aLIGO/__init__.py
@@ -0,0 +1,18 @@
+from gwinc.ifo.noises import *
+
+
+class aLIGO(nb.Budget):
+
+ name = 'Advanced LIGO'
+
+ noises = [
+ QuantumVacuum,
+ Seismic,
+ Newtonian,
+ SuspensionThermal,
+ CoatingBrownian,
+ CoatingThermoOptic,
+ SubstrateBrownian,
+ SubstrateThermoElastic,
+ ExcessGas,
+ ]
diff --git a/gwinc/ifo/aLIGO.yaml b/gwinc/ifo/aLIGO/ifo.yaml
similarity index 100%
rename from gwinc/ifo/aLIGO.yaml
rename to gwinc/ifo/aLIGO/ifo.yaml
diff --git a/gwinc/ifo/noises.py b/gwinc/ifo/noises.py
new file mode 100644
index 0000000000000000000000000000000000000000..2cea3c4c019af97523c36980a5e45e9a15301cbe
--- /dev/null
+++ b/gwinc/ifo/noises.py
@@ -0,0 +1,151 @@
+from .. import nb
+from .. import noise
+
+
+class QuantumVacuum(nb.Noise):
+ """Quantum Vacuum
+
+ """
+ style = dict(
+ label='Quantum',
+ color='#ad03de',
+ )
+
+ def calc(self):
+ return noise.quantum.shotrad(self.freq, self.ifo)
+
+
+class Seismic(nb.Noise):
+ """Seismic
+
+ """
+ style = dict(
+ label='Seismic',
+ color='#653700',
+ )
+
+ def calc(self):
+ return noise.seismic.seismic(self.freq, self.ifo)
+
+
+class Newtonian(nb.Noise):
+ """Newtonian Gravity
+
+ """
+ style = dict(
+ label='Newtonian Gravity',
+ color='#15b01a',
+ )
+
+ def calc(self):
+ return noise.newtonian.gravg(self.freq, self.ifo)
+
+
+class SuspensionThermal(nb.Noise):
+ """Suspension Thermal
+
+ """
+ style = dict(
+ label='Suspension Thermal',
+ color='#0d75f8',
+ )
+
+ def calc(self):
+ return noise.suspensionthermal.susptherm(self.freq, self.ifo)
+
+
+class CoatingBrownian(nb.Noise):
+ """Coating Brownian
+
+ """
+ style = dict(
+ label='Coating Brownian',
+ color='#fe0002',
+ )
+
+ def calc(self):
+ return noise.coatingthermal.coatbrownian(self.freq, self.ifo)
+
+
+class CoatingThermoOptic(nb.Noise):
+ """Coating Thermo-Optic
+
+ """
+ style = dict(
+ label='Coating Thermo-Optic',
+ color='#02ccfe',
+ linestyle='--',
+ )
+
+ def calc(self):
+ return noise.coatingthermal.thermooptic(self.freq, self.ifo)
+
+
+class ITMThermoRefractive(nb.Noise):
+ """ITM Thermo-Refractive
+
+ """
+ style = dict(
+ label='ITM Thermo-Refractive',
+ color='#448ee4',
+ linestyle='--',
+ )
+
+ def calc(self):
+ return noise.substratethermal.thermorefractiveITM(self.freq, self.ifo)
+
+
+class ITMCarrierDensity(nb.Noise):
+ """ITM Carrier Density
+
+ """
+ style = dict(
+ label='ITM Carrier Density',
+ color='#929591',
+ linestyle='--',
+ )
+
+ def calc(self):
+ return noise.substratethermal.carrierdensity(self.freq, self.ifo)
+
+
+class SubstrateBrownian(nb.Noise):
+ """Substrate Brownian
+
+ """
+ style = dict(
+ label='Substrate Brownian',
+ color='#fb7d07',
+ linestyle='--',
+ )
+
+ def calc(self):
+ return noise.substratethermal.subbrownian(self.freq, self.ifo)
+
+
+class SubstrateThermoElastic(nb.Noise):
+ """Substrate Thermo-Elastic
+
+ """
+ style = dict(
+ label='Substrate Thermo-Elastic',
+ color='#f5bf03',
+ linestyle='--',
+ )
+
+ def calc(self):
+ return noise.substratethermal.subtherm(self.freq, self.ifo)
+
+
+class ExcessGas(nb.Noise):
+ """Excess Gas
+
+ """
+ style = dict(
+ label='Excess Gas',
+ color='#ad900d',
+ linestyle='--',
+ )
+
+ def calc(self):
+ return noise.residualgas.gas(self.freq, self.ifo)
diff --git a/gwinc/io.py b/gwinc/io.py
new file mode 100644
index 0000000000000000000000000000000000000000..5184082732a6263d559e0483bf6056b6ea5a5bad
--- /dev/null
+++ b/gwinc/io.py
@@ -0,0 +1,60 @@
+import datetime
+import h5py
+
+
+SCHEMA = 'GWINC noise budget'
+SCHEMA_VERSION = 1
+
+
+def _write_trace_recursive(grp, traces):
+ for name, trace in traces.items():
+ if isinstance(trace, dict):
+ tgrp = grp.create_group(name)
+ _write_trace_recursive(tgrp, trace)
+ else:
+ data, style = trace
+ dset = grp.create_dataset(name, data=data)
+ for key, val in style.items():
+ dset.attrs[key] = val
+
+
+def save_hdf5(path, freq, traces, **kwargs):
+ """Save GWINC traces dict to HDF5 file.
+
+ See HDF5_SCHEMA.
+
+ """
+ with h5py.File(path, 'w') as f:
+ f.attrs['SCHEMA'] = SCHEMA
+ f.attrs['SCHEMA_VERSION'] = SCHEMA_VERSION
+ # FIXME: add budget code hash or something
+ f.attrs['date'] = datetime.datetime.now().isoformat()
+ for key, val in kwargs.items():
+ f.attrs[key] = val
+ f.create_dataset('Freq', data=freq)
+ tgrp = f.create_group('traces')
+ _write_trace_recursive(tgrp, traces)
+
+
+def _read_trace_recursive(element):
+ trace = {}
+ for name, item in element.items():
+ if isinstance(item, h5py.Group):
+ trace[name] = _read_trace_recursive(item)
+ else:
+ trace[name] = item.value, dict(item.attrs.items())
+ return trace
+
+
+def load_hdf5(path):
+ """Load traces from HDF5 file.
+
+ Returns a recursive traces dictionary. See HDF5_SCHEMA.
+
+ """
+ with h5py.File(path, 'r') as f:
+ # FIXME: check SCHEMA name/version
+ freq = f['Freq'].value
+ traces = _read_trace_recursive(f['/traces'])
+ attrs = dict(f.attrs.items())
+ return freq, traces, attrs
diff --git a/gwinc/nb.py b/gwinc/nb.py
new file mode 100644
index 0000000000000000000000000000000000000000..9434ef13388a5e49a6543c002cb1d7d46120ede4
--- /dev/null
+++ b/gwinc/nb.py
@@ -0,0 +1,420 @@
+import os
+import logging
+import itertools
+import importlib
+import importlib.util
+import collections
+import numpy as np
+import scipy.interpolate
+
+
+def quadsum(data):
+ """Calculate quadrature sum of list of data arrays.
+
+ Provided data are assumed to be power-referred, so this is a
+ simple point-by-point sum.
+
+ NaNs in sum elements do not contribute to sum.
+
+ """
+ return np.nansum(data, 0)
+
+
+
+class BudgetItem:
+ """GWINC BudgetItem class
+
+ """
+ def load(self):
+ """Overload method for initial loading of static data.
+
+ """
+ return None
+
+ def update(self, **kwargs):
+ """Overload method for updating data needed to calculate final PSD.
+
+ By default any keyword arguments provided are written directly
+ as attribute variables (as with __init__).
+
+ """
+ for key, val in kwargs.items():
+ setattr(self, key, val)
+
+ def calc(self):
+ """Overload method for calculation of final PSD.
+
+ Should return an array of power-referenced values evaluated at
+ all evaluation frequencies (self.freq).
+
+ """
+ return None
+
+ ##########
+
+ def __init__(self, freq, **kwargs):
+ """Initialize budget item.
+
+ Primary argument is the evaluation frequency array. Any
+ keyword arguments provided are simple written as attribute
+ variables in the initialized object.
+
+ """
+ self.__freq = freq
+ for key, val in kwargs.items():
+ setattr(self, key, val)
+
+ @property
+ def name(self):
+ """"Name of this BudgetItem class."""
+ return self.__class__.__name__
+
+ def __str__(self):
+ # FIXME: provide info on internal state (load/update/calc/etc.)
+ return '<{} {}>'.format(
+ ' '.join([c.__name__ for c in self.__class__.__bases__]),
+ self.name,
+ )
+
+ @property
+ def freq(self):
+ """Evaluation frequency array supplied at initialization."""
+ return self.__freq
+
+ def interpolate(self, freq, data):
+ """Interpolate data to the evaluation frequencies.
+
+ """
+ func = scipy.interpolate.interp1d(
+ freq, data,
+ kind='nearest',
+ copy=False,
+ assume_sorted=True,
+ bounds_error=False,
+ fill_value=np.nan,
+ )
+ return func(self.freq)
+
+
+class Calibration(BudgetItem):
+ """GWINC Calibration class
+
+ BudgetItem that represents a calibration transfer function for a
+ Noise. The calc() method should return a transfer function
+ amplitude array evaluated at the evaluation frequencies supplied
+ at initialization and available in the `freq` array attribute
+ (self.freq).
+
+ """
+ def __call__(self, data):
+ """Calibrate input data.
+
+ Returns calibrated version of input data array,
+ e.g. point-by-point product of data and calibration arrays.
+
+ """
+ cal = self.calc()
+ assert data.shape == cal.shape, \
+ "data shape does not match calibration ({} != {})".format(data.shape, cal.shape)
+ return data * cal
+
+
+class Noise(BudgetItem):
+ """GWINC Noise class
+
+ BudgetItem that represents a PSD noise calculation. The calc()
+ method should return the noise PSD spectrum array evaluated at the
+ evaluation frequencies supplied at initialization and available in
+ the `freq` array attribute (self.freq).
+
+ """
+
+ style = {}
+ """Trace plot style dictionary"""
+
+ def calc_trace(self, calibration=None, calc=True):
+ """Returns noise (PSD, style) tuple.
+
+ If `calibration` is not None it is assumed to be a
+ len(self.freq) array that will be multiplied to the output
+ PSD.
+
+ If calc=False, the noise will not be calculated and the PSD
+ will be None. This is useful for just getting style the
+ style.
+
+ """
+ if calc:
+ data = self.calc()
+ if calibration is not None:
+ data *= calibration
+ else:
+ data = None
+ return data, self.style
+
+
+class Budget(Noise):
+ """GWINC Budget class
+
+ This is a Noise that represents the budget of multiple sub noises.
+
+ The `noises` attribute of this class should be list constituent
+ Noise classes. Each element can be either a single Noise class,
+ or a tuple of (Noise, Calibration) classes, e.g.:
+
+ noises = [
+ Thermal,
+ (Shot, Sensing),
+ ]
+
+ When this object is initialized, all sub noises and calibrations
+ are initialized. Pre-defined load() and update() methods call the
+ load() and update() methods of all sub noises and calibrations.
+ When calc() is called, the PSD is calculated for all sub noises,
+ the relevant calibration is evaluated and applied, and the
+ quadrature sum of all calibrated consituent noises is returned.
+
+ Additionally a `references` attribute may be definied, similar to
+ the `noises` attribute described above except that the specified
+ noises do not contribute to the overall budget total.
+
+ NOTE: an `ifo` attribute is always passed as an initialization
+ argument to sub noises.
+
+ """
+
+ noises = []
+ """List of constituent noise classes"""
+
+ references = []
+ """List of reference nosie classes"""
+
+ def __init__(self, *args, noises=None, **kwargs):
+ """Initialize Budget object.
+
+ See BudgetItem for base initialization arguments.
+
+ If a `noises` keyword argument is provided it should be an
+ iterable of noise names (constituent or reference) which will
+ be used to filter the noises initialized in this budget.
+
+ """
+ super().__init__(*args, **kwargs)
+ # store args and kwargs for later use
+ self.args = args
+ self.kwargs = kwargs
+ # FIXME: special casing the IFO here, in case it's defined as
+ # a class attribute rather than passed at initialization. we
+ # do this because we're not defining a standard way to extract
+ # IFO variables that get passed around in a reasonable way.
+ # how can we clarify this?
+ if 'ifo' not in kwargs and getattr(self, 'ifo', None):
+ self.kwargs['ifo'] = getattr(self, 'ifo', None)
+ # all noise objects keyed by name
+ self._noise_objs = collections.OrderedDict()
+ # all cal objects keyed by name
+ self._cal_objs = {}
+ # noise to calibration mapping
+ self._noise_cal = {}
+ # set of all constituent budget noise names
+ self._budget_noises = set()
+ # initialize all noise objects
+ for nc in self.noises:
+ name, noise_obj, cal = self.__init_noise(nc)
+ if noises and name not in noises:
+ continue
+ self.__add_noise(noise_obj, cal)
+ self._budget_noises.add(name)
+ for nc in self.references:
+ name, noise_obj, cal = self.__init_noise(nc)
+ if noises and name not in noises:
+ continue
+ self.__add_noise(noise_obj, cal)
+ # error if requested noise is not present
+ if noises:
+ sset = set(noises)
+ nset = set([name for name in self._noise_objs.keys()])
+ if not sset <= nset:
+ raise AttributeError("unknown noise terms: {}".format(' '.join(sset-nset)))
+
+ def __init_noise(self, nc):
+ cal = None
+ if isinstance(nc, (list, tuple)):
+ noise, cal = nc[:2]
+ else:
+ noise = nc
+ noise_obj = noise(
+ *self.args,
+ **self.kwargs
+ )
+ return noise_obj.name, noise_obj, cal
+
+ def __add_noise(self, noise_obj, cal):
+ logging.debug("init {}".format(noise_obj))
+ # instantiate the noise object
+ name = noise_obj.name
+ self._noise_objs[name] = noise_obj
+ if cal is not None:
+ # if a cal object is specified, instantiate and store
+ cal_obj = cal(*self.args, **self.kwargs)
+ if cal_obj.name not in self._cal_objs:
+ self._cal_objs[cal_obj.name] = cal_obj
+ self._noise_cal[name] = cal_obj.name
+
+ def __getitem__(self, name):
+ try:
+ return self._noise_objs[name]
+ except KeyError:
+ try:
+ return self._cal_objs[name]
+ except KeyError:
+ raise KeyError("unknown noise or cal name '{}".format(name))
+
+ def keys(self):
+ """Iterate over budget noise names."""
+ return iter(self._noise_objs.keys())
+
+ def values(self):
+ """Iterate over budget noise objects."""
+ return iter(self._noise_objs.values())
+
+ def items(self):
+ """Iterate over budget noise (name, object) tuples."""
+ return iter(self._noise_objs.items())
+
+ def __iter__(self):
+ return iter(self.keys())
+
+ def walk(self):
+ """Walk recursively through every BudgetItem in the budget.
+
+ This includes Noise, Calibration and Budget objects, as well
+ as any decendents of Budget objects.
+
+ For each leaf item yields a tuple of all ancestor objects,
+ e.g.:
+
+ (self)
+ (self, BudgetItem)
+ (self, ChildBudget1)
+ (self, ChildBudget1, BudgetItem)
+ ...
+
+ """
+ yield (self,)
+ for item in itertools.chain(
+ self._cal_objs.values(),
+ self._noise_objs.values()):
+ if isinstance(item, Budget):
+ for i in item.walk():
+ yield (self,) + i
+ else:
+ yield (self, item)
+
+ def load(self):
+ """Load all noise and cal objects."""
+ for name, item in itertools.chain(
+ self._cal_objs.items(),
+ self._noise_objs.items()):
+ logging.debug("load {}".format(item))
+ item.load()
+
+ def update(self, **kwargs):
+ """Update all noise and cal objects with supplied kwargs."""
+ for name, item in itertools.chain(
+ self._cal_objs.items(),
+ self._noise_objs.items()):
+ logging.debug("update {}".format(item))
+ item.update(**kwargs)
+
+ def cal_for_noise(self, name):
+ """Return the calibration object for named noise."""
+ try:
+ return self._cal_objs[self._noise_cal[name]]
+ except KeyError:
+ return None
+
+ def calc_noise(self, name):
+ """Return calibrated individual noise term.
+
+ The noise PSD and calibration transfer functions are
+ calculated, and the calibrated noise array is returned.
+
+ """
+ noise = self._noise_objs[name]
+ nd = noise.calc()
+ cal = self.cal_for_noise(name)
+ if cal:
+ cd = cal.calc()
+ return nd * cd
+ else:
+ return nd
+
+ def calc(self):
+ """Calculate sum of all noises.
+
+ """
+ data = [self.calc_noise(name) for name in self._noise_objs.keys() if name in self._budget_noises]
+ return quadsum(data)
+
+ def calc_trace(self, calibration=None, calc=True):
+ """Returns a dictionary of noises traces, keyed by noise names.
+
+ Values are (data, style) trace tuples (see Noise.calc_trace).
+ The key of the budget sum total is 'Total'. The values of sub
+ budgets are themselves dictionaries returned from
+ calc_trace() of the sub budget.
+
+ If `calibration` is not None it is assumed to be a
+ len(self.freq) array that will be multiplied to the output
+ PSD of the budget and all sub noises.
+
+ If calc=False, the noise will not be calculated and the PSD
+ will be None. This is useful for just getting style the
+ style.
+
+ """
+ # start by creating an empty OrderedDict used for outputing trace data
+ # or style info with the following order:
+ # references
+ # total
+ # constituents
+ d = collections.OrderedDict()
+ # allocate references
+ for name, noise in self._noise_objs.items():
+ if name in self._budget_noises:
+ continue
+ d[name] = noise.calc_trace(calc=False)
+ # allocate total
+ if self._budget_noises:
+ d['Total'] = None, self.style
+ # allocate constituent
+ for name, noise in self._noise_objs.items():
+ if name not in self._budget_noises:
+ continue
+ d[name] = noise.calc_trace(calc=False)
+ # if we're not calc'ing, just return the dict now
+ if not calc:
+ return d
+ # calc all noises
+ for name, noise in self._noise_objs.items():
+ # extract/calc the budget-level calibration for this noise
+ cal_obj = self.cal_for_noise(name)
+ if cal_obj:
+ cal = cal_obj.calc()
+ else:
+ cal = np.ones_like(self.freq)
+ # then multiply by the supplied calibration
+ if calibration is not None:
+ cal *= calibration
+ d[name] = noise.calc_trace(calibration=cal, calc=True)
+ # calc budget total
+ constituent_data = []
+ for name in self._budget_noises:
+ if isinstance(d[name], dict):
+ data = d[name]['Total'][0]
+ else:
+ data = d[name][0]
+ constituent_data.append(data)
+ d['Total'] = quadsum(constituent_data), self.style
+ return d
diff --git a/gwinc/plot.py b/gwinc/plot.py
index 16f6da5ef67300e5685ce7aa02cb098bb5a59a17..4aec5cfa8b8dce10e60aff1c5c665ee7286488bc 100644
--- a/gwinc/plot.py
+++ b/gwinc/plot.py
@@ -1,91 +1,19 @@
from numpy import sqrt
-PRIORITY_MAP = {
- 'Quantum Vacuum': 0,
- 'Seismic': 1,
- 'Newtonian Gravity': 2,
- 'Suspension Thermal': 3,
- 'Coating Brownian': 4,
- 'Coating Thermo-Optic': 5,
- 'ITM Thermo-Refractive': 6,
- 'ITM Carrier Density': 7,
- 'Substrate Brownian': 8,
- 'Substrate Thermo-Elastic': 9,
- 'Excess Gas': 10,
-}
-
-STYLE_MAP = {
- 'Quantum Vacuum': dict(
- # color = 'xkcd:vibrant purple',
- color = '#ad03de',
- ),
- 'Seismic': dict(
- # color = 'xkcd:brown',
- color = '#653700',
- ),
- 'Newtonian Gravity': dict(
- # color = 'xkcd:green',
- color = '#15b01a',
- ),
- 'Atmospheric Infrasound': dict(
- # color = 'xkcd:neon green',
- color = '#0cff0c',
- ls = '--',
- ),
- 'Suspension Thermal': dict(
- # color = 'xkcd:deep sky blue',
- color = '#0d75f8',
- ),
- 'Coating Brownian': dict(
- # color = 'xkcd:fire engine red',
- color = '#fe0002',
- ),
- 'Coating Thermo-Optic': dict(
- # color = 'xkcd:bright sky blue',
- color = '#02ccfe',
- ls = '--',
- ),
- 'ITM Thermo-Refractive': dict(
- # color = 'xkcd:dark sky blue',
- color = '#448ee4',
- ls = '--',
- ),
- 'ITM Carrier Density': dict(
- # color = 'xkcd:grey',
- color = '#929591',
- ls = '--',
- ),
- 'Substrate Brownian': dict(
- # color = 'xkcd:pumpkin orange',
- color = '#fb7d07',
- ls = '--',
- ),
- 'Substrate Thermo-Elastic': dict(
- # color = 'xkcd:golden',
- color = '#f5bf03',
- ls = '--',
- ),
- 'Excess Gas': dict(
- # color = 'xkcd:baby shit brown',
- color = '#ad900d',
- ls = '--',
- ),
-}
-
def plot_noise(
- noises,
+ freq,
+ traces,
ax=None,
+ **kwargs
):
"""Plot a GWINC noise budget from calculated noises
If an axis handle is provided it will be used for the plot.
- Returns the figure handle used.
+ Returns the figure handle.
"""
- f = noises['Freq']
-
if ax is None:
import matplotlib.pyplot as plt
fig = plt.figure()
@@ -93,31 +21,27 @@ def plot_noise(
else:
fig = ax.figure
- def plot_dict(noises):
- #use sorted to force a consistent ordering
- #The key lambda in sorted gets the (name, noise) pair and so nn[0] returns the name
- #the return tuple causes sorting by priority, with a fallback to lexical sort on the name
- for name, noise in sorted(noises.items(), key = lambda nn : (PRIORITY_MAP.get(nn[0], 100), nn[0])):
- if name in ['Freq', 'Total']:
- continue
- if isinstance(noise, dict):
- plot_dict(noise)
- else:
- noise = noises[name]
- stylekw = dict(
- color = (0, 0, 0),
- label = name,
- lw = 3,
- )
- try:
- style = STYLE_MAP[name]
- stylekw.update(style)
- except KeyError:
- pass
- ax.loglog(f, sqrt(noise), **stylekw)
- plot_dict(noises)
-
- ax.loglog(f, sqrt(noises['Total']), color='#000000', alpha=0.6, label='Total', lw=4)
+ for name, trace in traces.items():
+ if isinstance(trace, dict):
+ data, style = trace['Total']
+ else:
+ data, style = trace
+ # assuming all data is PSD
+ data = sqrt(data)
+ if name == 'Total':
+ style = dict(
+ color='#000000',
+ alpha=0.6,
+ lw=4,
+ )
+ ylim = [min(data)/10, max(data)]
+ if 'label' not in style:
+ style['label'] = name
+ if 'linewidth' in style:
+ style['lw'] = style['linewidth']
+ elif 'lw' not in style:
+ style['lw'] = 3
+ ax.loglog(freq, data, **style)
ax.grid(
True,
@@ -132,9 +56,17 @@ def plot_noise(
fontsize='small',
)
+ ax.autoscale(enable=True, axis='y', tight=True)
+ if 'ylim' in kwargs:
+ ax.set_ylim(kwargs['ylim'])
+ else:
+ ax.set_ylim(ylim)
+ ax.set_xlim(freq[0], freq[-1])
+
ax.set_xlabel('Frequency [Hz]')
- ax.set_ylabel(u"Strain [1/\u221AHz]")
- ax.set_xlim([min(f), max(f)])
- ax.set_ylim([3e-25, 1e-21])
+ if 'ylabel' in kwargs:
+ ax.set_ylabel(kwargs['ylabel'])
+ if 'title' in kwargs:
+ ax.set_title(kwargs['title'])
return fig
diff --git a/gwinc/test/A+.pkl b/gwinc/test/A+.pkl
deleted file mode 100644
index 53308c4037b0eabca523951b605c60cac89049d4..0000000000000000000000000000000000000000
--- a/gwinc/test/A+.pkl
+++ /dev/null
@@ -1,3 +0,0 @@
-version https://git-lfs.github.com/spec/v1
-oid sha256:226246bdc5600e865629bfc56bae363a3a6db04d1acfdc665b8aa95edc7e3d38
-size 897625
diff --git a/gwinc/test/Aplus.pkl b/gwinc/test/Aplus.pkl
new file mode 100644
index 0000000000000000000000000000000000000000..1edf74d3b2833413ec2013fd1ea6a10dbe94b2a2
Binary files /dev/null and b/gwinc/test/Aplus.pkl differ
diff --git a/gwinc/test/__main__.py b/gwinc/test/__main__.py
index f648311bc3e817c511289e1a7debc966072cc014..dd46d6532d857b6dc916789389dd807cc99e8b5a 100644
--- a/gwinc/test/__main__.py
+++ b/gwinc/test/__main__.py
@@ -12,7 +12,8 @@ import logging
logging.basicConfig(format='%(message)s',
level=os.getenv('LOG_LEVEL', logging.INFO))
-from .. import load_ifo, gwinc
+from .. import load_ifo
+from ..gwinc import gwinc
from ..gwinc_matlab import gwinc_matlab
try:
import inspiral_range
@@ -57,7 +58,7 @@ def main():
args = parser.parse_args()
logging.info("loading IFO '{}'...".format(args.IFO))
- ifo = load_ifo(args.IFO)
+ Budget, ifo, freq, plot_style = load_ifo(args.IFO)
freq = np.logspace(np.log10(FLO), np.log10(FHI), NPOINTS)
diff --git a/gwinc/util.py b/gwinc/util.py
index b0a7618261a2dfddee7d1bf977f77db7ef4c7ac7..3f07dffef74ec53ad4a9ae86c2111dac0a4ac3a2 100644
--- a/gwinc/util.py
+++ b/gwinc/util.py
@@ -1,34 +1,39 @@
-import datetime
-import h5py
+import os
+import importlib
-from . import struct
+def lpath(file0, file1):
+ """Return path of file1 when expressed relative to file0.
-def load_hdf5(path):
- """Load IFO and noises from HDF5 file.
+ For instance, if file0 is "/path/to/file0" and file1 is
+ "../for/file1" then what is returned is "/path/for/file1".
- Returns (name, ifo, noises) tuple load from file.
+ This is useful for resolving paths within packages with e.g.:
+
+ rpath = lpath(__file__, '../resource.yaml')
"""
- with h5py.File(path, 'r') as f:
- name = f.attrs['name']
- ifo = struct.Struct.from_yaml(f.attrs['IFO'])
- noises = {}
- for name, data in f['/traces'].items():
- noises[name] = data.value
- return name, ifo, noises
+ return os.path.abspath(os.path.join(os.path.dirname(file0), file1))
+
+def load_module(name_or_path):
+ """Load module from name or path.
-def save_hdf5(name, ifo, noises, path):
- """Save IFO and noises to HDF5 file.
+ Return loaded module and module path.
"""
- with h5py.File(path, 'w') as f:
- f.attrs['schema'] = 'GWINC noise budget'
- # FIXME: add GWINC version
- f.attrs['date'] = datetime.datetime.now().isoformat()
- f.attrs['name'] = name
- f.attrs['IFO'] = ifo.to_yaml()
- tgrp = f.create_group('traces')
- for noise, data in noises.items():
- tgrp.create_dataset(noise, data=data)
+ if os.path.exists(name_or_path):
+ path = name_or_path.rstrip('/')
+ modname = os.path.splitext(os.path.basename(path))[0]
+ if os.path.isdir(path):
+ path = os.path.join(path, '__init__.py')
+ spec = importlib.util.spec_from_file_location(modname, path)
+ mod = importlib.util.module_from_spec(spec)
+ spec.loader.exec_module(mod)
+ else:
+ mod = importlib.import_module(name_or_path)
+ try:
+ path = mod.__path__[0]
+ except AttributeError:
+ path = mod.__file__
+ return mod, path