Skip to content
Snippets Groups Projects

Compare revisions

Changes are shown as if the source revision was being merged into the target revision. Learn more about comparing revisions.

Source

Select target project
No results found

Target

Select target project
  • wenxuan.jia/pygwinc
  • sean-leavey/pygwinc
  • sebastian.steinlechner/pygwinc
  • nicholas.demos/pygwinc
  • chris.whittle/pygwinc
  • raymond.robie/pygwinc
  • mateusz.bawaj/pygwinc
  • anchal.gupta/pygwinc
  • 40m/pygwinc
  • evan.hall/pygwinc
  • kevin.kuns/pygwinc
  • geoffrey-lovelace/pygwinc
  • brittany.kamai/pygwinc
  • daniel-brown/pygwinc
  • lee-mcculler/pygwinc
  • jameson.rollins/pygwinc
  • gwinc/pygwinc
17 results
Show changes
Commits on Source (25)
Hide whitespace changes
Inline Side-by-side
Showing
with 539 additions and 257 deletions
.gitignore
View file @ b69d4267
#for docs and setup.py outputs
build/
# test cache
gwinc/test/cache
# Byte-compiled / optimized / DLL files
__pycache__/
......
.gitlab-ci.yml
View file @ b69d4267
image: igwn/base:buster
stages:
- test
- deploy
- dist
- test
- review
- deploy
# have to specify this so that all jobs execute for all commits
# including merge requests
workflow:
rules:
- if: $CI_MERGE_REQUEST_ID
- if: $CI_COMMIT_BRANCH
variables:
GIT_STRATEGY: clone
test:
# build the docker image we will use in all the jobs, with all
# dependencies pre-installed/configured.
gwinc/base:
stage: dist
variables:
IMAGE_TAG: $CI_REGISTRY_IMAGE/$CI_JOB_NAME:$CI_COMMIT_REF_NAME
GIT_STRATEGY: none
script:
- docker login -u gitlab-ci-token -p $CI_JOB_TOKEN $CI_REGISTRY
- |
cat <<EOF > Dockerfile
FROM igwn/base:buster
RUN apt-get update -qq
RUN apt-get -y install --no-install-recommends git python3-gitlab python3 python3-yaml python3-scipy python3-matplotlib python3-ipython python3-lalsimulation python3-pypdf2 python3-h5py
RUN git clone https://gitlab-ci-token:ci_token@git.ligo.org/gwinc/inspiral_range.git
EOF
- docker build -t $IMAGE_TAG .
- docker push $IMAGE_TAG
# create plots for the canonical IFOs
generate_budgets:
stage: test
before_script:
- echo $CI_COMMIT_SHA | cut -b1-8 > gitID.txt
image: $CI_REGISTRY_IMAGE/gwinc/base:$CI_COMMIT_REF_NAME
script:
- rm -rf ifo gwinc_test_report.pdf
- mkdir ifo
- apt-get update -qq
- apt-get install -y -qq git python3-yaml python3-scipy python3-matplotlib python3-ipython lalsimulation-python3 python3-pypdf2 python3-h5py
- git clone https://gitlab-ci-token:ci_token@git.ligo.org/gwinc/inspiral_range.git
- export PYTHONPATH=inspiral_range
- export MPLBACKEND=agg
- python3 -m gwinc.test -r gwinc_test_report.pdf
- for ifo in aLIGO Aplus Voyager CE1 CE2; do
- python3 -m gwinc $ifo -s ifo/$ifo.png
- python3 -m gwinc $ifo -s ifo/$ifo.h5
- done
- python3 -m gwinc.ifo -s ifo/all_compare.png
after_script:
- rm gitID.txt
cache:
key: "$CI_PROJECT_NAMESPACE:$CI_PROJECT_NAME:$CI_JOB_NAME"
untracked: true
- mkdir -p ifo
- export PYTHONPATH=/inspiral_range
- for ifo in $(python3 -c "import gwinc; print(' '.join(gwinc.IFOS))"); do
- python3 -m gwinc $ifo -s ifo/$ifo.png
- python3 -m gwinc $ifo -s ifo/$ifo.h5
- done
- python3 -m gwinc.ifo -s ifo/all_compare.png
artifacts:
when: always
expire_in: 4w
paths:
- ifo
- gwinc_test_report.pdf
- ifo
# this is a special job intended to run only for merge requests.
# budgets are compared against those from the target branch. if the
# merge request has not yet been approved and noise changes are found,
# the job will fail. once the merge request is approved the job can
# be re-run, at which point the pipeline should succeed allowing the
# merge to be merged.
noise_change_approval:
stage: review
rules:
# - if: '$CI_MERGE_REQUEST_TARGET_BRANCH_NAME == "master"'
- if: $CI_MERGE_REQUEST_ID
image: $CI_REGISTRY_IMAGE/gwinc/base:$CI_COMMIT_REF_NAME
script:
- |
cat <<EOF > check_approved.py
import sys
import gitlab
project_id = sys.argv[1]
mr_iid = sys.argv[2]
# this only works for public repos, otherwise need to specify
# private_token=
gl = gitlab.Gitlab('https://git.ligo.org')
project = gl.projects.get(project_id)
mr = project.mergerequests.get(mr_iid)
approvals = mr.approvals.get()
print(approvals.approved)
EOF
- echo $CI_MERGE_REQUEST_PROJECT_ID, $CI_MERGE_REQUEST_IID, $CI_MERGE_REQUEST_TARGET_BRANCH_NAME,
- approved=$(python3 check_approved.py $CI_MERGE_REQUEST_PROJECT_ID $CI_MERGE_REQUEST_IID )
- if [[ $approved != True ]] ; then
- target=origin/$CI_MERGE_REQUEST_TARGET_BRANCH_NAME
- if ! python3 -m gwinc.test --git-rev $target -r gwinc_test_report.pdf ; then
- echo "NOISE CHANGES RELATIVE TO $CI_MERGE_REQUEST_TARGET_BRANCH_NAME."
- echo "Approval required to merge this branch."
- /bin/false
- else
- echo "No noise changes detected."
- fi
- else
- echo "Merge request approved, noise change accepted."
- fi
artifacts:
when: on_failure
paths:
- gwinc_test_report.pdf
expose_as: 'noise changes relative to target branch'
# generate the html doc web pages. the "pages" job has special
# meaning, as it's "public" artifact becomes the directory served
# through gitlab static pages
pages:
stage: deploy
dependencies:
- test
only:
- master
needs:
- job: generate_budgets
artifacts: true
image: $CI_REGISTRY_IMAGE/gwinc/base:$CI_COMMIT_REF_NAME
script:
- rm -rf public
- mv ifo public
- apt-get install -y -qq python3-sphinx-rtd-theme
- cd docs
- make html
- cd ..
- mv ./build/sphinx/html/* public/
- rm -rf public
- apt-get install -y -qq python3-sphinx-rtd-theme
- cd docs
- make html
- cd ..
- mv ./build/sphinx/html public
- mv ifo public/
artifacts:
when: always
paths:
- public
expire_in: 4w
only:
- master
- public
CONTRIBUTING.md
View file @ b69d4267
......@@ -2,25 +2,48 @@
The pygwinc project welcomes your contributions. Our policy is that
all contributions should be peer-reviewed. To facilitate the review,
please do not commit your work to this repository yourself. Instead,
fork the repository and send a merge request.
please do not commit your work to this repository directly. Instead,
please fork the repository and create a [merge
request](https://git.ligo.org/gwinc/pygwinc/-/merge_requests/new)
against the main pygwinc master branch.
`pygwinc` comes with a test suite that will calculate all canonical
IFO budgets with the current state of the code, and compare them
against cached hdf5 traces (stored in the repo with
[git-lfs](https://git-lfs.github.com/) ( [see
also](https://wiki.ligo.org/Computing/GitLFS)). Contributors should
run the tests before submitting any merge requests:
When submitting code for merge, please follow good coding practice.
Respect the existing coding style, which for `pygwinc` is standard
[PEP8](https://www.python.org/dev/peps/pep-0008/) (with some
exceptions). Make individual commits as logically distinct and atomic
as possible, and provide a complete, descriptive log of the changes
(including a top summary line). Review efficiency follows code
legibility.
`pygwinc` comes with a validation command that can compare budgets
from the current code against those produced from different versions
in git (by default it compares against the current HEAD). The command
can be run with:
```shell
$ python3 -m gwinc.test
```
The test suite will be run automatically upon push to git.ligo.org as
part of the continuous integration, and code that fails the CI will
not be merged, unless failures are well justified.
Use the '--plot' or '--report' options to produce visual comparisons
of the noise differences. The comparison can be done against an
arbitrary commit using the '--git-ref' option. Traces for referenced
commits are cached, which speeds up subsequent comparison runs
significantly.
Once you submit your merge request a special CI job will determine if
there are budgets differences between your code and the master branch.
If there are, explicit approval from reviewers will be required before
your changes can be merged (see "approving noise" below).
## For reviewers: approving noise curve changes
If your change affects the shape of a noise curve, your commit message
should make note of that, and provide a justification. It will also
be necessary to update the reference curves stored in cache, but that
should be done in a separate commit not a part of the original merge
request, so that all reviewers can see the changes introduced in the
CI test failure report.
As discussed above, merge requests that generate noise changes will
cause a pipeline failure in the `review:noise_change_approval` CI job.
The job will generate a report comparing the new noise traces against
those from master, which can be found under the 'View exposed
artifacts' menu item in the pipeline report. Once you have reviewed
the report and the code, and understand and accept the noise changes,
click the 'Approve' button in the MR. Once sufficient approval has
been given, `review:noise_change_approval` job can be re-run, which
should now pick up that approval has been given and allow the pipeline
to succeed. Once the pipeline succeeds the merge request can be
merged. Click the 'Merge' button to finally merge the code.
IFO.md
View file @ b69d4267
......@@ -2,43 +2,43 @@
CI-generated plots and data for all IFOs included in pygwinc.
![IFO compare](https://gwinc.docs.ligo.org/pygwinc/all_compare.png)
![IFO compare](https://gwinc.docs.ligo.org/pygwinc/ifo/all_compare.png)
## aLIGO
* [ifo.yaml](gwinc/ifo/aLIGO/ifo.yaml)
* [aLIGO.h5](https://gwinc.docs.ligo.org/pygwinc/aLIGO.h5)
* [aLIGO.h5](https://gwinc.docs.ligo.org/pygwinc/ifo/aLIGO.h5)
![aLIGO](https://gwinc.docs.ligo.org/pygwinc/aLIGO.png)
![aLIGO](https://gwinc.docs.ligo.org/pygwinc/ifo/aLIGO.png)
## A+
* [ifo.yaml](gwinc/ifo/Aplus/ifo.yaml)
* [Aplus.h5](https://gwinc.docs.ligo.org/pygwinc/Aplus.h5)
* [Aplus.h5](https://gwinc.docs.ligo.org/pygwinc/ifo/Aplus.h5)
![Aplus](https://gwinc.docs.ligo.org/pygwinc/Aplus.png)
![Aplus](https://gwinc.docs.ligo.org/pygwinc/ifo/Aplus.png)
## Voyager
* [ifo.yaml](gwinc/ifo/Voyager/ifo.yaml)
* [Voyager.h5](https://gwinc.docs.ligo.org/pygwinc/Voyager.h5)
* [Voyager.h5](https://gwinc.docs.ligo.org/pygwinc/ifo/Voyager.h5)
![Voyager](https://gwinc.docs.ligo.org/pygwinc/Voyager.png)
![Voyager](https://gwinc.docs.ligo.org/pygwinc/ifo/Voyager.png)
## Cosmic Explorer 1
* [ifo.yaml](gwinc/ifo/CE1/ifo.yaml)
* [CE1.h5](https://gwinc.docs.ligo.org/pygwinc/CE1.h5)
* [CE1.h5](https://gwinc.docs.ligo.org/pygwinc/ifo/CE1.h5)
![CE1](https://gwinc.docs.ligo.org/pygwinc/CE1.png)
![CE1](https://gwinc.docs.ligo.org/pygwinc/ifo/CE1.png)
## Cosmic Explorer 2
* [ifo.yaml](gwinc/ifo/CE2/ifo.yaml)
* [CE2.h5](https://gwinc.docs.ligo.org/pygwinc/CE2.h5)
* [CE2.h5](https://gwinc.docs.ligo.org/pygwinc/ifo/CE2.h5)
![CE2](https://gwinc.docs.ligo.org/pygwinc/CE2.png)
![CE2](https://gwinc.docs.ligo.org/pygwinc/ifo/CE2.png)
README.md
View file @ b69d4267
......@@ -2,12 +2,14 @@
# Python Gravitational Wave Interferometer Noise Calculator
![gwinc](https://gwinc.docs.ligo.org/pygwinc/aLIGO.png)
[![aLIGO](https://gwinc.docs.ligo.org/pygwinc/ifo/aLIGO.png "Canonical
IFOs")](IFO.md)
`pygwinc` is a multi-faceted tool for processing and plotting noise
budgets for ground-based gravitational wave detectors. It's primary
feature is a collection of mostly analytic noise calculation functions
for various sources of noise affecting detectors (`gwinc.noise`):
feature is a collection of mostly analytic [noise calculation
functions](#noise-functions) for various sources of noise affecting
detectors (`gwinc.noise`):
* quantum noise
* mirror coating thermal noise
......@@ -17,22 +19,20 @@ for various sources of noise affecting detectors (`gwinc.noise`):
* Newtonian/gravity-gradient noise
* residual gas noise
See [noise functions](#noise-functions) below.
`pygwinc` also includes a generalized noise budgeting tool
(`gwinc.nb`) that allows users to create arbitrary noise budgets (for
any experiment, not just ground-based GW detectors) using measured or
analytically calculated data. See the [budget
interface](#Budget-interface) section below.
`pygwinc` is also a generalized noise budgeting tool (`gwinc.nb`) that
allows users to create arbitrary noise budgets (for any experiment,
not just ground-based GW detectors) using measured or analytically
calculated data. See the [budget interface](#Budget-interface)
section below.
`pygwinc` includes canonical budgets for various well-known current
and future detectors (`gwinc.ifo`):
and future GW detectors (`gwinc.ifo`):
* [aLIGO](https://gwinc.docs.ligo.org/pygwinc/aLIGO.png)
* [A+](https://gwinc.docs.ligo.org/pygwinc/Aplus.png)
* [Voyager](https://gwinc.docs.ligo.org/pygwinc/Voyager.png)
* [Cosmic Explorer 1](https://gwinc.docs.ligo.org/pygwinc/CE1.png)
* [Cosmic Explorer 2](https://gwinc.docs.ligo.org/pygwinc/CE2.png)
* [aLIGO](https://gwinc.docs.ligo.org/pygwinc/ifo/aLIGO.png)
* [A+](https://gwinc.docs.ligo.org/pygwinc/ifo/Aplus.png)
* [Voyager](https://gwinc.docs.ligo.org/pygwinc/ifo/Voyager.png)
* [Cosmic Explorer 1](https://gwinc.docs.ligo.org/pygwinc/ifo/CE1.png)
* [Cosmic Explorer 2](https://gwinc.docs.ligo.org/pygwinc/ifo/CE2.png)
See [IFO.md](IFO.md) for the latest CI-generated plots and hdf5 cached
data.
......@@ -68,6 +68,12 @@ Materials.Coating.Philown 5e-05 3e-05
$ python3 -m gwinc my_aLIGO.yaml
```
You can also use the `--ifo` option to change parameters from the
command line:
```shell
$ python3 -m gwinc aLIGO --ifo Optics.SRM.Tunephase=3.14
```
Stand-alone YAML files will always assume the nominal ['aLIGO' budget
description](gwinc/ifo/aLIGO).
......@@ -325,13 +331,7 @@ traces = budget.run()
```
The IFOs included in `gwinc.ifo` provide examples of the use of the
budget interface:
* [aLIGO](gwinc/ifo/aLIGO)
* [Aplus](gwinc/ifo/Aplus)
* [Voyager](gwinc/ifo/Voyager)
* [CE1](master/gwinc/ifo/CE1)
* [CE2](master/gwinc/ifo/CE2)
budget interface (e.g. [gwinc.ifo.aLIGO](gwinc/ifo/aLIGO)).
### extracting single noise terms
......
gwinc/__main__.py
View file @ b69d4267
......@@ -34,8 +34,12 @@ 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.
Individual IFO parameters can be overriden with the --ifo option:
gwinc --ifo Optics.SRM.Tunephase=3.14 ...
If the inspiral_range package is installed, various figures of merit
can be calculated for the resultant spectrum with the --fom argument,
can be calculated for the resultant spectrum with the --fom option,
e.g.:
gwinc --fom horizon ...
......@@ -46,39 +50,59 @@ See documentation for inspiral_range package for details.
"""
IFO = 'aLIGO'
FLO = 5
FHI = 6000
NPOINTS = 3000
FREQ = '5:3000:6000'
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,
help="upper frequency bound in Hz [{}]".format(FHI))
parser.add_argument('--npoints', '-n', default=NPOINTS,
help="number of frequency points [{}]".format(NPOINTS))
parser.add_argument('--title', '-t',
help="plot title")
parser.add_argument('--fom',
help="calculate inspiral range for resultant spectrum ('func[:param=val,param=val]')")
parser.add_argument(
'--freq', '-f', metavar='FSPEC',
help="frequency array specification in Hz, either as 'flo:fhi' or 'flo:npoints:fhi' [{}]".format(FREQ))
parser.add_argument(
'--ifo', '-o',
#nargs='+', action='extend',
action='append',
help="override budget IFO parameter (may be specified multiple times)")
parser.add_argument(
'--title', '-t',
help="plot title")
parser.add_argument(
'--fom',
help="calculate inspiral range for resultant spectrum ('func[:param=val,param=val]')")
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 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',
help="print IFO as text table to stdout and exit")
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',
help="IFO name, description file path (.yaml, .mat, .m), budget module (.py), or HDF5 data file (.hdf5, .h5)")
group.add_argument(
'--interactive', '-i', action='store_true',
help="interactive plot with interactive shell")
group.add_argument(
'--save', '-s',
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',
help="print IFO as text table to stdout and exit")
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',
help="IFO name, description file path (.yaml, .mat, .m), budget module (.py), or HDF5 data file (.hdf5, .h5)")
def freq_from_spec(spec):
fspec = spec.split(':')
if len(fspec) == 2:
fspec = fspec[0], FREQ.split(':')[1], fspec[1]
return np.logspace(
np.log10(float(fspec[0])),
np.log10(float(fspec[2])),
int(fspec[1]),
)
def main():
......@@ -95,16 +119,27 @@ def main():
freq, traces, attrs = load_hdf5(args.IFO)
ifo = getattr(attrs, 'IFO', None)
plot_style = attrs
if args.freq:
logging.warning("ignoring frequency specification for frequencies defined in HDF5...")
else:
Budget = load_budget(args.IFO)
ifo = Budget.ifo
# FIXME: this should be done only if specified, to allow for
# using any FREQ specified in the Budget
freq = np.logspace(np.log10(args.flo), np.log10(args.fhi), args.npoints)
if args.freq:
try:
freq = freq_from_spec(args.freq)
except IndexError:
parser.error("improper frequency specification '{}'".format(args.freq))
else:
freq = getattr(Budget, 'freq', freq_from_spec(FREQ))
plot_style = getattr(Budget, 'plot_style', {})
traces = None
if args.ifo:
for paramval in args.ifo:
param, val = paramval.split('=', 1)
ifo[param] = float(val)
if args.yaml:
if not ifo:
parser.exit(2, "no IFO structure available.")
......
gwinc/noise/substratethermal.py
View file @ b69d4267
......@@ -61,8 +61,6 @@ def substrate_carrierdensity(f, materials, wBeam, exact=False):
def substrate_thermorefractive(f, materials, wBeam, exact=False):
"""Substrate thermal displacement noise spectrum from thermorefractive fluctuations
For semiconductor substrates.
:f: frequency array in Hz
:materials: gwinc optic materials structure
:wBeam: beam radius (at 1 / e^2 power)
......
gwinc/suspension.py
View file @ b69d4267
......@@ -100,7 +100,6 @@ def getJointParams(sus, n):
TempUpper = stages[n-1].Temp
else:
TempUpper = sus.Temp
TempUpper = Temp # FIXME: reproduces the old calculation
##############################
# material parameters
......@@ -114,10 +113,6 @@ def getJointParams(sus, n):
elif 'WireMaterial' in stage:
wireMatUpper = stage.WireMaterial
wireMatLower = stage.WireMaterial
# FIXME: reproduces the old calculation
elif n == last_stage and sus.Type == 'BQuad':
wireMatUpper = 'Silicon'
wireMatLower = 'Silicon'
elif n == last_stage:
wireMatUpper = 'Silica'
wireMatLower = 'Silica'
......@@ -132,9 +127,6 @@ def getJointParams(sus, n):
# support blade (upper joint only)
if 'BladeMaterial' in stage:
bladeMat = stage.BladeMaterial
# FIXME: reproduces the old calculation
elif n == last_stage and sus.Type == 'BQuad':
bladeMat = 'Silicon'
else:
bladeMat = 'MaragingSteel'
......@@ -151,11 +143,13 @@ def wireGeometry(r, N, RibbonThickness=None, TaperedEndRadius=None, **kwargs):
"""Compute wire geometry-dependent factors
r is the wire radius, or ribbon width.
N is the number of wires or ribbons
RibbonThickness must be set when ribbons are used, or
TaperedEndRadius when tapered fibers are used.
Other kwargs are ignored.
Returns the cross-sectional area and moment of inertia,
and the modified surface to volume ratios (vertical and horizontal).
Returns cross-sectional areas (central and end) and moment of inertia,
and modified surface to volume ratios (vertical and horizontal).
"""
# Usual case: round wire/fiber
......@@ -206,25 +200,32 @@ def wireGeometry(r, N, RibbonThickness=None, TaperedEndRadius=None, **kwargs):
return (xsect, xsectEnd, xII, mu_v, mu_h)
def wireTELoss(w, tension, xsectEnd, xII, Temp, alpha, beta, rho, C, K, Y, xsect,
def wireTELoss(w, tension, xsectEnd, xII, Temp, alpha, beta, rho, C, K, Y,
RibbonThickness=None, **kwargs):
"""Thermoelastic calculations for wires
Repeated for upper and lower joint of each stage.
w = angular frequency
tension = weight supported per wire
xsectEnd = cross sectional area of wire end
xII = cross sectional moment of inertia
Temp = temperature
alpha = coeff of thermal expansion
beta = temp dependence of Young's modulus
rho = mass density
C = heat capacity
K = thermal conductivity W/(m K)
Y = Young's modulus
RibbonThickness must be set when ribbons are used
Other kwargs are ignored
Returns the loss angle associated with thermoelastic damping
(wire horizontal)
"""
# horizontal TE time constant, wires
# The constant 7.37e-2 is 1/(4*q0^2) from eq 12, C. Zener 10.1103/PhysRev.53.90 (1938)
tau = 7.37e-2 * 4 * (rho * C * xsectEnd) / (pi * K)
# FIXME: reproduces the old calculation
# then xsect can be dropped from the argument list
tau = 7.37e-2 * 4 * (rho * C * xsect) / (pi * K)
# deal with ribbon geometry
if RibbonThickness is not None:
......@@ -248,6 +249,9 @@ def bladeTELoss(w, t, Temp, alpha, beta, rho, C, K, Y):
"""Thermoelastic calculations for blades
Invoked for upper joint only (there is no lower blade)
w = angular frequency
t = blade thickness
Temp = temperature
alpha = coeff of thermal expansion
beta = temp dependence of Young's modulus
rho = mass density
......@@ -255,6 +259,9 @@ def bladeTELoss(w, t, Temp, alpha, beta, rho, C, K, Y):
K = thermal conductivity W/(m K)
Y = Young's modulus
Returns the loss angle associated with thermoelastic damping
(blade vertical)
"""
# vertical TE time constant, blades
tau = (rho * C * t**2) / (K * pi**2)
......@@ -285,6 +292,17 @@ def bladeTELoss(w, t, Temp, alpha, beta, rho, C, K, Y):
def continuumWireKh(w, N, length, tension, xsect, xII, rho, Y, phi):
"""Horizontal spring constant, including violin modes
w = angular frequency
N = number of wires
length = wire length
tension = weight supported per wire
xsect = wire cross sectional area
xII = cross sectional moment of inertia
rho = mass density
Y = Young's modulus
phi = loss angle
Returns the spring constant (wire horizontal)
"""
Y = Y * (1 + 1j * phi) # complex Young's modulus
......@@ -306,16 +324,12 @@ def continuumWireKh(w, N, length, tension, xsect, xII, rho, Y, phi):
# = T k (cos(k L) + k delta sin(k L))
# for w -> 0, this reduces to N_w * T * k
khnum = N * tension * k * (coskl + dk * sinkl)
# FIXME: reproduces the old calculation
khnum = N * tension * k * (1 + dk**2) * (coskl + dk * sinkl)
# denominator, horiz spring constant
# D after equation 8 in GG
# D = sin(k L) - 2 k delta cos(k L)
# for w -> 0, this reduces to k (L - 2 delta)
khden = sinkl - 2 * dk * coskl
# FIXME: reproduces the old calculation
khden = (1 - dk**2) * sinkl - 2 * dk * coskl
return khnum/khden
......@@ -324,6 +338,18 @@ def continuumWireKv(w, N, length, xsect, xsectEnd, rho, Y, phi,
TaperedEndLength=None, **kwargs):
"""Vertical spring constant, including bounce mode.
w = angular frequency
N = number of wires
length = wire length
xsect = wire cross sectional area
xsectEnd = cross sectional area of wire end
rho = mass density
Y = Young's modulus
phi = loss angle
TaperedEndLength must be set when tapered fibers are used
Other kwargs are ignored
Returns the spring constant (wire vertical)
"""
Y = Y * (1 + 1j * phi) # complex Young's modulus
k = sqrt(rho / Y) * w
......@@ -336,9 +362,6 @@ def continuumWireKv(w, N, length, xsect, xsectEnd, rho, Y, phi,
kv_mid = N * xsect * Y * k / tan(k * l_mid)
kv_end = N * xsectEnd * Y * k / tan(k * l_end)
kv = 1/(2/kv_end + 1/kv_mid)
# FIXME: reproduces the old calculation
kv_end = N * xsectEnd * Y * k / tan(k * 2*l_end)
kv = 1/(1/kv_end + 1/kv_mid)
return kv
......@@ -491,7 +514,7 @@ def suspQuad(f, sus):
# nominally this term becomes zero for steel wires as ds_w is zero
# The last term is for TE loss
phih_TE = wireTELoss(w, tension, xsectEnd, xII, Temp,
alpha_w, beta_w, rho_w, C_w, K_w, Y_w, xsect,
alpha_w, beta_w, rho_w, C_w, K_w, Y_w,
**wireShape)
phih = phi_w * (1 + mu_h * ds_w) + phih_TE
......
gwinc/test/__main__.py
View file @ b69d4267
import os
import sys
import glob
import shutil
import signal
import logging
import tempfile
import argparse
import subprocess
import numpy as np
import matplotlib.pyplot as plt
from collections import OrderedDict
from collections.abc import Mapping
from PyPDF2 import PdfFileReader, PdfFileWriter
logging.basicConfig(format='%(message)s',
level=os.getenv('LOG_LEVEL', logging.INFO))
from .. import IFOS, load_budget
from ..io import load_hdf5, save_hdf5
from ..io import load_hdf5
try:
import inspiral_range
......@@ -23,8 +22,103 @@ except ImportError:
inspiral_range = None
logging.basicConfig(
format='%(message)s',
level=os.getenv('LOG_LEVEL', logging.INFO))
TOLERANCE = 1e-6
CACHE_PATH = os.path.join(os.path.dirname(__file__), 'cache')
CACHE_LIMIT = 5
def test_path(*args):
"""Return path to package file."""
return os.path.join(os.path.dirname(__file__), *args)
def git_find_upstream_name():
try:
remotes = subprocess.run(
['git', 'remote', '-v'],
capture_output=True, universal_newlines=True,
check=True,
).stdout
except subprocess.CalledProcessError as e:
logging.error(e.stderr.split('\n')[0])
for remote in remotes.split('\n'):
name, url, fp = remote.split()
if 'gwinc/pygwinc.git' in url:
return name
def git_rev_resolve_hash(git_rev):
"""Resolve a git revision into its hash string."""
try:
return subprocess.run(
['git', 'show', '-s', '--format=format:%H', git_rev],
capture_output=True, universal_newlines=True,
check=True,
).stdout
except subprocess.CalledProcessError as e:
logging.error(e.stderr.split('\n')[0])
def prune_cache_dir():
"""Prune all but the N most recently accessed caches.
"""
cache_dir = test_path('cache')
if not os.path.exists(cache_dir):
return
expired_paths = sorted(
[os.path.join(cache_dir, path) for path in os.listdir(cache_dir)],
key=lambda path: os.stat(path).st_atime, reverse=True,
)[CACHE_LIMIT:]
if not expired_paths:
return
for path in expired_paths:
logging.info("pruning old cache: {}".format(path))
shutil.rmtree(path)
def gen_cache(git_hash, path):
"""generate cache for specified git hash at the specified path
The included shell script is used to extract the gwinc code from
the appropriate git commit, and invoke a new python instance to
generate the noise curves.
"""
logging.info("creating new cache for hash {}...".format(git_hash))
subprocess.run(
[test_path('gen_cache.sh'), git_hash, path],
check=True,
)
def load_cache(path):
"""load a cache from the specified path
returns a "cache" dictionary with 'git_hash' and 'ifos' keys.
"""
logging.info("loading cache {}...".format(path))
cache = {}
git_hash_path = os.path.join(path, 'git_hash')
if os.path.exists(git_hash_path):
with open(git_hash_path) as f:
git_hash = f.read().strip()
else:
git_hash = None
logging.debug("cache hash: {}".format(git_hash))
cache['git_hash'] = git_hash
cache['ifos'] = {}
for f in sorted(os.listdir(path)):
name, ext = os.path.splitext(f)
if ext != '.h5':
continue
cache['ifos'][name] = os.path.join(path, f)
return cache
def walk_traces(traces, root=()):
......@@ -105,17 +199,15 @@ def compare_traces(tracesA, tracesB, tolerance=TOLERANCE, skip=None):
return diffs
def plot_diffs(freq, diffs, tolerance,
name, labelA, labelB, fom_title='',
save=None):
def plot_diffs(freq, diffs, styleA, styleB):
spec = (len(diffs)+1, 2)
sharex = None
for i, nname in enumerate(diffs):
noiseA, noiseB, frac = diffs[nname]
axl = plt.subplot2grid(spec, (i, 0), sharex=None)
axl.loglog(freq, np.sqrt(noiseA), label=labelA)
axl.loglog(freq, np.sqrt(noiseB), label=labelB)
axl.loglog(freq, np.sqrt(noiseA), **styleA)
axl.loglog(freq, np.sqrt(noiseB), **styleB)
axl.grid()
axl.legend(loc='upper right')
axl.set_ylabel(nname)
......@@ -134,133 +226,164 @@ def plot_diffs(freq, diffs, tolerance,
if i == 0:
axr.set_title("fractional difference")
plt.suptitle('''{} {}/{} noise comparison
(noises that differ by more than {} ppm)
{}'''.format(name, labelA, labelB, tolerance*1e6, fom_title))
axl.set_xlabel("frequency [Hz]")
axr.set_xlabel("frequency [Hz]")
plt.subplots_adjust(top=0.8, right=0.85, wspace=0.3)
if save:
pwidth = 10
pheight = (len(diffs) * 5) + 2
plt.gcf().set_size_inches(pwidth, pheight)
plt.savefig(save)
else:
plt.show()
##################################################
def main():
parser = argparse.ArgumentParser()
parser.add_argument('--tolerance', '-t', type=float, default=TOLERANCE,
help='fractional tolerance [{}]'.format(TOLERANCE))
parser.add_argument('--skip', '-k', metavar='NOISE', action='append',
help='traces to skip in comparison (multiple may be specified)')
parser.add_argument('--cache', '-c', metavar='PATH', default=CACHE_PATH,
help='specify alternate IFO traces cache path')
parser = argparse.ArgumentParser(
description="""GWINC noise validation
This command calculates the canonical noise budgets with the current
code and compares them against those calculated with code from a
specified git revision. You must be running from a git checkout of
the source for this to work. The command will fail if it detects any
noise differences. Plots or a PDF report of differences can be
generated with the '--plot' or '--report' commands respectively.
By default it will attempt to determine the git reference for upstream
master for your current configuration (usually 'origin/master' or
'upstream/master'). You may specify an arbitrary git revision with
the --git-rev command. For example, to compare against another
remote/branch use:
$ python3 -m gwinc.test --git-rev remote/dev-branch
or if you have uncommitted changes compare against the current head
with:
$ python3 -m gwinc.test -g HEAD
See gitrevisions(7) for various ways to refer to git revisions.
A cache of traces from reference git revisions will be stored in
gwinc/test/cache/<SHA1>. Old caches are automatically pruned.""",
formatter_class=argparse.RawDescriptionHelpFormatter,
)
parser.add_argument(
'--tolerance', '-t', type=float, default=TOLERANCE,
help="fractional tolerance of comparison [{}]".format(TOLERANCE))
parser.add_argument(
'--skip', '-k', metavar='NOISE', action='append',
help="traces to skip in comparison (multiple may be specified)")
rgroup = parser.add_mutually_exclusive_group()
rgroup.add_argument('--plot', '-p', action='store_true',
help='plot differences')
rgroup.add_argument('--report', '-r', metavar='REPORT.pdf',
help='create PDF report of test results (only created if differences found)')
rgroup.add_argument('--gen-cache', action='store_true',
help='update/create IFO traces cache directory')
parser.add_argument('ifo', metavar='IFO', nargs='*',
help='specific ifos to test (default all)')
rgroup.add_argument(
'--git-rev', '-g', metavar='REV',
help="specify specific git revision to compare against")
ogroup = parser.add_mutually_exclusive_group()
ogroup.add_argument(
'--plot', '-p', action='store_true',
help="show interactive plot differences")
ogroup.add_argument(
'--report', '-r', metavar='REPORT.pdf',
help="create PDF report of test results (only created if differences found)")
parser.add_argument(
'ifo', metavar='IFO', nargs='*',
help="specific ifos to test (default all)")
args = parser.parse_args()
if args.gen_cache:
try:
os.makedirs(args.cache)
except FileExistsError:
pass
freq = np.logspace(np.log10(5), np.log10(6000), 3000)
for name in IFOS:
Budget = load_budget(name)
traces = Budget(freq).run()
path = os.path.join(args.cache, name+'.h5')
save_hdf5(path, freq, traces)
return
# get the reference hash
if args.git_rev:
git_rev = args.git_rev
else:
remote = git_find_upstream_name()
if not remote:
sys.exit("Could not resolve upstream remote name")
git_rev = '{}/master'.format(remote)
logging.info("presumed upstream git reference: {}".format(git_rev))
git_hash = git_rev_resolve_hash(git_rev)
if not git_hash:
sys.exit("Could not resolve reference, could not run test.")
logging.info("git hash: {}".format(git_hash))
# load the cache
cache_path = test_path('cache', git_hash)
if not os.path.exists(cache_path):
prune_cache_dir()
gen_cache(git_hash, cache_path)
cache = load_cache(cache_path)
if args.report:
base, ext = os.path.splitext(args.report)
if ext != '.pdf':
parser.error("Test reports only support PDF format.")
outdir = tempfile.TemporaryDirectory()
# find all cached IFOs
logging.info("loading cache {}...".format(args.cache))
cached_ifos = {}
for f in sorted(os.listdir(args.cache)):
name, ext = os.path.splitext(f)
if ext != '.h5':
continue
cached_ifos[name] = os.path.join(args.cache, f)
# select
if args.ifo:
ifos = {name:cached_ifos[name] for name in args.ifo}
ifos = args.ifo
else:
ifos = cached_ifos
ifos = IFOS
labelA = 'cache'
labelB = 'head'
style_cache = dict(label='reference', linestyle='-')
style_head = dict(label='head', linestyle='--')
fail = False
# compare
for name, path in ifos.items():
for name in ifos:
logging.info("{} tests...".format(name))
freq, tracesA, attrs = load_hdf5(path)
try:
path = cache['ifos'][name]
except KeyError:
logging.warning("IFO {} not found in cache")
fail |= True
continue
freq, traces_cache, attrs = load_hdf5(path)
Budget = load_budget(name)
tracesB = Budget(freq).run()
traces_head = Budget(freq).run()
if inspiral_range:
totalA = tracesA['Total'][0]
totalB = tracesB['Total'][0]
total_cache = traces_cache['Total'][0]
total_head = traces_head['Total'][0]
range_func = inspiral_range.range
H = inspiral_range.waveform.CBCWaveform(freq)
fomA = range_func(freq, totalA, H=H)
tracesA['int73'] = inspiral_range.int73(freq, totalA)[1], None
fomB = range_func(freq, totalB, H=H)
tracesB['int73'] = inspiral_range.int73(freq, totalB)[1], None
fom_cache = range_func(freq, total_cache, H=H)
traces_cache['int73'] = inspiral_range.int73(freq, total_cache)[1], None
fom_head = range_func(freq, total_head, H=H)
traces_head['int73'] = inspiral_range.int73(freq, total_head)[1], None
fom_summary = """
inspiral {func} {m1}/{m2} Msol:
{labelA}: {fomA:.2f} Mpc
{labelB}: {fomB:.2f} Mpc
{label_cache}: {fom_cache:.2f} Mpc
{label_head}: {fom_head:.2f} Mpc
""".format(
func=range_func.__name__,
m1=H.params['m1'],
m2=H.params['m2'],
labelA=labelA,
fomA=fomA,
labelB=labelB,
fomB=fomB,
label_cache=style_cache['label'],
fom_cache=fom_cache,
label_head=style_head['label'],
fom_head=fom_head,
)
else:
fom_summary = ''
diffs = compare_traces(tracesA, tracesB, args.tolerance, args.skip)
diffs = compare_traces(traces_cache, traces_head, args.tolerance, args.skip)
if diffs:
logging.warning("{} tests FAIL".format(name))
fail |= True
if args.plot or args.report:
if args.report:
save = os.path.join(outdir.name, name+'.pdf')
else:
save = None
plot_diffs(
freq, diffs, args.tolerance,
name, labelA, labelB, fom_summary,
save=save,
)
else:
if not diffs:
logging.info("{} tests pass.".format(name))
continue
logging.warning("{} tests FAIL".format(name))
fail |= True
if args.plot or args.report:
plot_diffs(freq, diffs, style_cache, style_head)
plt.suptitle('''{} {}/{} noise comparison
(noises that differ by more than {} ppm)
reference git hash: {}
{}'''.format(name, style_cache['label'], style_head['label'],
args.tolerance*1e6, cache['git_hash'], fom_summary))
if args.report:
pwidth = 10
pheight = (len(diffs) * 5) + 2
plt.gcf().set_size_inches(pwidth, pheight)
plt.savefig(os.path.join(outdir.name, name+'.pdf'))
else:
plt.show()
if not fail:
logging.info("all tests pass.")
......
gwinc/test/cache/Aplus.h5 deleted 100644 → 0 LFS
View file @ 9ff4ba54
File deleted
gwinc/test/cache/CE1.h5 deleted 100644 → 0 LFS
View file @ 9ff4ba54
File deleted
gwinc/test/cache/CE2.h5 deleted 100644 → 0 LFS
View file @ 9ff4ba54
File deleted
gwinc/test/cache/Voyager.h5 deleted 100644 → 0 LFS
View file @ 9ff4ba54
File deleted
gwinc/test/cache/aLIGO.h5 deleted 100644 → 0 LFS
View file @ 9ff4ba54
File deleted
gwinc/test/gen_cache.sh 0 → 100755
View file @ b69d4267
#!/bin/bash -e
if [ -z "$1" ] || [ -z "$2" ] ; then
echo "usage: $(basename $0) git_hash cache_dir_path"
echo "generate a cache of IFO budget traces from a particular git commit"
exit 1
fi
git_hash="$1"
cache_dir="$2"
mkdir -p $cache_dir
cache_dir=$(cd $cache_dir && pwd)
gwinc_dir=$cache_dir/gwinc
mkdir -p $gwinc_dir
git archive $git_hash | tar -x -C $gwinc_dir
cd $gwinc_dir
export LOG_LEVEL=INFO
for ifo in $(python3 -c "import gwinc; print(' '.join(gwinc.IFOS))") ; do
python3 -m gwinc --save $cache_dir/${ifo}.h5 $ifo
done
echo $git_hash > $cache_dir/git_hash