From 03d046cda475da44873bddbee9125cc30e29bbed Mon Sep 17 00:00:00 2001
From: Jesse Mapel <jmapel@usgs.gov>
Date: Mon, 24 Feb 2020 14:26:25 -0700
Subject: [PATCH] Fixed covariance clalculations
---
plio/utils/covariance.py | 14 +++++++++-----
plio/utils/tests/test_covariance.py | 8 ++++----
2 files changed, 13 insertions(+), 9 deletions(-)
diff --git a/plio/utils/covariance.py b/plio/utils/covariance.py
index d4cbb56..ca12526 100644
--- a/plio/utils/covariance.py
+++ b/plio/utils/covariance.py
@@ -1,7 +1,7 @@
import numpy as np
import math
-def compute_covariance(lat, lon, rad, latsigma=10., lonsigma=10., radsigma=15., semimajor_axis=1737400.0):
+def compute_covariance(lat, lon, rad, latsigma=10., lonsigma=10., radsigma=15., semimajor_axis=None):
"""
Given geospatial coordinates, desired accuracy sigmas, and an equitorial radius, compute a 2x3
sigma covariange matrix.
@@ -27,7 +27,8 @@ def compute_covariance(lat, lon, rad, latsigma=10., lonsigma=10., radsigma=15.,
The desired radius accuracy in meters (Defualt: 15.0)
semimajor_axis : float
- The semi-major or equitorial radius in meters (Default: 1737400.0 - Moon)
+ The semi-major or equitorial radius in meters. If not entered,
+ the local radius will be used.
Returns
-------
@@ -37,12 +38,15 @@ def compute_covariance(lat, lon, rad, latsigma=10., lonsigma=10., radsigma=15.,
"""
lat = math.radians(lat)
lon = math.radians(lon)
-
+
+ if semimajor_axis is None:
+ semimajor_axis = rad
+
# SetSphericalSigmasDistance
scaled_lat_sigma = latsigma / semimajor_axis
# This is specific to each lon.
- scaled_lon_sigma = lonsigma * math.cos(lat) / semimajor_axis
+ scaled_lon_sigma = lonsigma / (math.cos(lat) * semimajor_axis)
# SetSphericalSigmas
cov = np.eye(3,3)
@@ -77,4 +81,4 @@ def compute_covariance(lat, lon, rad, latsigma=10., lonsigma=10., radsigma=15.,
rectcov[1,0] = mat[1,1]
rectcov[1,1] = mat[1,2]
rectcov[1,2] = mat[2,2]
- return rectcov
\ No newline at end of file
+ return rectcov
diff --git a/plio/utils/tests/test_covariance.py b/plio/utils/tests/test_covariance.py
index 35bc270..a2c921d 100644
--- a/plio/utils/tests/test_covariance.py
+++ b/plio/utils/tests/test_covariance.py
@@ -8,12 +8,12 @@ def test_compute_covariance():
# This test is using values from an ISIS control network
lat = 86.235596
lon = 140.006195
- rad = 1736777.625
+ rad = 1736777.625
sigmalat = 15.0
sigmalon = 15.0
sigmarad = 25.0
semimajor_rad = 1737400.0
cov = covariance.compute_covariance(lat, lon, rad, sigmalat, sigmalon, sigmarad, semimajor_rad)
- expected = np.array([[132.97888775695, -111.55453178747, -20.08405416233],
- [93.588991760138, 16.848822261184, 623.27512692323]])
- np.testing.assert_almost_equal(cov, expected)
\ No newline at end of file
+ expected = np.array([[225.85120968, -0.84930178, -20.08405416233],
+ [225.55132129, 16.848822261184, 623.27512692323]])
+ np.testing.assert_almost_equal(cov, expected)
--
GitLab