From c7fa8197e3160ca591ebe23b5268b35a701fe293 Mon Sep 17 00:00:00 2001
From: Jameson Graef Rollins <jrollins@finestructure.net>
Date: Tue, 22 May 2018 17:30:46 -0700
Subject: [PATCH] move MATLAB constants into separate const module

This is mostly used for MATLAB translation, but is useful for a couple of
other things as well (R_earth for instance).
---
 gwinc/const.py        | 27 +++++++++++++++
 gwinc/gwinc_matlab.py | 80 +++++++++----------------------------------
 gwinc/ifo/A+.yaml     |  2 --
 gwinc/precomp.py      |  3 +-
 4 files changed, 45 insertions(+), 67 deletions(-)
 create mode 100644 gwinc/const.py

diff --git a/gwinc/const.py b/gwinc/const.py
new file mode 100644
index 0000000..8c01e84
--- /dev/null
+++ b/gwinc/const.py
@@ -0,0 +1,27 @@
+import scipy.constants
+
+
+CONSTANTS = {
+    'c': scipy.constants.c,
+    'G': scipy.constants.G,
+    'g': scipy.constants.g,
+    'E0': scipy.constants.epsilon_0,
+    'hbar': scipy.constants.hbar,
+    'kB': scipy.constants.k,
+
+    'yr': scipy.constants.year,
+    'AU': scipy.constants.astronomical_unit,
+    'parsec': scipy.constants.parsec,
+    'Mpc': scipy.constants.parsec * 1e6,
+
+    # FIXME: use astropy for the following astrophysical/cosmological
+    # constants
+    'R_earth': 6.3781e6,
+    'SolarMassParameter': 1.3271244e20,
+    'MSol': 1.3271244e20 / scipy.constants.G,
+    # http://physics.nist.gov/cuu/Constants/
+    'H0': 67110,
+    # http://arxiv.org/pdf/1303.5076v3.pdf
+    'omegaM': 0.3175,
+    'omegaLambda': 1 - 0.3175,
+}
diff --git a/gwinc/gwinc_matlab.py b/gwinc/gwinc_matlab.py
index f419337..a315bbf 100644
--- a/gwinc/gwinc_matlab.py
+++ b/gwinc/gwinc_matlab.py
@@ -1,10 +1,10 @@
 import os
 import tempfile
 import scipy.io
-import scipy.constants
 import numpy as np
 import logging
 
+from . import const
 from .struct import Struct
 
 ##################################################
@@ -96,68 +96,6 @@ class Matlab:
 
 ##################################################
 
-def ifo_add_constants(ifo):
-    """Add "constants" sub-Struct to ifo Struct
-
-    This is required by MATLAB gwinc.
-
-    """
-    # permittivity of free space [F / m]
-    #ifo.Constants.E0      = 8.8541878176e-12
-    ifo.Constants.E0      = scipy.constants.epsilon_0
-    # Plancks constant [J s]
-    #ifo.Constants.hbar    = 1.054572e-34
-    ifo.Constants.hbar    = scipy.constants.hbar
-    # Planks constant [J s]
-    #ifo.Constants.h       = ifo.Constants.hbar * 2 * pi
-    # Boltzman constant [J / K]
-    #ifo.Constants.kB      = 1.380658e-23
-    ifo.Constants.kB      = scipy.constants.k
-    # gas constant [J / (K * mol)]
-    #ifo.Constants.R       = 8.31447215
-    # electron mass [kg]
-    #ifo.Constants.m_e     = 9.10938291e-31
-    # speed of light in vacuum [m / s]
-    #ifo.Constants.c       = 2.99792458e8
-    ifo.Constants.c       = scipy.constants.c
-    # seconds in a year [s]
-    ifo.Constants.yr      = 365.2422 * 86400
-    # mass of Earth [kg]
-    #ifo.Constants.M_earth = 5.972e24
-    # radius of Earth [m]
-    ifo.Constants.R_earth = 6.3781e6
-    # sampling frequency [Hz]
-    #ifo.Constants.fs      = 16384
-    # Astronomical unit, IAU 2012 Resolution B2 [m]
-    ifo.Constants.AU      = 149597870700
-    # IAU 2015 Resolution B2 [m]
-    ifo.Constants.parsec  = ifo.Constants.AU * (648000 / np.pi)
-    # IAU 2015 Resolution B2 [m]
-    ifo.Constants.Mpc     = ifo.Constants.parsec * 1e6
-    # IAU 2015 Resolution B3 [m^3 / s^2; G * MSol]
-    ifo.Constants.SolarMassParameter = 1.3271244e20
-    # gravitational const [m^3 / (kg  s^2)]
-    #ifo.Constants.G       = 6.67408e-11
-    ifo.Constants.G       = scipy.constants.G
-    # solar mass [kg]
-    # http://arxiv.org/abs/1507.07956
-    ifo.Constants.MSol    = ifo.Constants.SolarMassParameter / ifo.Constants.G
-    # gravitational acceleration [m / s^2]
-    #ifo.Constants.g       = 9.806
-    ifo.Constants.g       = scipy.constants.g
-    # Hubble constant [ms^( - 1)]
-    # http://physics.nist.gov/cuu/Constants/
-    ifo.Constants.H0      = 67110
-    # http://arxiv.org/pdf/1303.5076v3.pdf
-    # Mass density parameter
-    ifo.Constants.omegaM  = 0.3175
-    # http://arxiv.org/pdf/1303.5076v3.pdf
-    # Cosmological constant density parameter
-    # omegaK = 0 (flat universe) is assumed
-    ifo.Constants.omegaLambda = 1 - ifo.Constants.omegaM
-    return ifo
-
-
 NOISE_NAME_MAP = {
     'Quantum': 'Quantum Vacuum',
     'Newtonian': 'Newtonian Gravity',
@@ -168,6 +106,20 @@ NOISE_NAME_MAP = {
     'SuspThermal': 'Suspension Thermal',
     'ResGas': 'Excess Gas',
 }
+
+
+def ifo_matlab_transform(ifo):
+    """Prep the ifo structure for use with MATLAB gwinc
+
+    * add "constants" sub-Struct
+
+    """
+    # add constants
+    ifo.Constants = Struct.from_dict(const.CONSTANTS)
+
+    return ifo
+
+
 def _rename_noises(d):
     nd = {}
     for k,v in d.items():
@@ -199,7 +151,7 @@ def gwinc_matlab(f, ifo, plot=False):
     matlab = Matlab()
 
     # add Constants attribute to ifo structure
-    ifo_add_constants(ifo)
+    ifo_matlab_transform(ifo)
 
     matlab.load_array('f', f)
     matlab.load_struct('ifo', ifo)
diff --git a/gwinc/ifo/A+.yaml b/gwinc/ifo/A+.yaml
index 536e47d..8dc1cef 100644
--- a/gwinc/ifo/A+.yaml
+++ b/gwinc/ifo/A+.yaml
@@ -30,8 +30,6 @@
 # Updated numbers March 2018: LIGO-T1800044
 
 Constants:
-  # earth radius
-  R_earth: 6.3781e6
   # Temperature of the Vacuum
   Temp: 290                       # K
 
diff --git a/gwinc/precomp.py b/gwinc/precomp.py
index 1c501e8..fae365d 100644
--- a/gwinc/precomp.py
+++ b/gwinc/precomp.py
@@ -4,6 +4,7 @@ from scipy.io import loadmat
 import scipy.special
 import logging
 
+from .const import CONSTANTS
 from .struct import Struct
 from .noise.coatingthermal import getCoatDopt
 
@@ -33,7 +34,7 @@ def precompIFO(ifo, PRfixed=True):
 
     if 'VHCoupling' not in ifo.Suspension:
         ifo.Suspension.VHCoupling = Struct()
-        ifo.Suspension.VHCoupling.theta = ifo.Infrastructure.Length / ifo.Constants.R_earth
+        ifo.Suspension.VHCoupling.theta = ifo.Infrastructure.Length / CONSTANTS['R_earth']
 
     ##############################
     # optics values
-- 
GitLab