Released v2.0.6

README.md

@@ -37,7 +37,7 @@ Installation
 Development
 -----------
 
-1. Download and extract <https://github.com/ut-astria/orbdetpy/releases/download/2.0.5/orekit-data.tar.gz>
+1. Download and extract <https://github.com/ut-astria/orbdetpy/releases/download/2.0.6/orekit-data.tar.gz>
    under the `orbdetpy/` sub-folder.
 
 2. Developers will need Apache Maven 3+ to build the Java library. Build

examples/test_estimation.py

@@ -81,11 +81,11 @@
 # Propagate orbits and generate measurements
 meas = list(m for m in propagate_orbits([cfg])[0].array)
 
-# Zero valued measurements are ignored but they force the filters to output
-# propagated states/covariances at arbitrary intermediate times
+# Empty measurements are ignored in orbit estimation but they force the filters
+# to output propagated states/covariances at arbitrary intermediate times
 t0 = get_J2000_epoch_offset("2019-05-01T00:26:02.853Z")
 for i in range(5):
-    meas.append(build_measurement(t0 + i*60.0, "Maui", [0.0, 0.0]))
+    meas.append(build_measurement(t0 + i*60.0, "", []))
 
 # Measurements must be sorted in ascending chronological order
 meas.sort(key=lambda m: m.time)
@@ -99,8 +99,7 @@
 
 # Perturb truth initial state before running OD
 cfg.prop_initial_state[:] = multivariate_normal(cfg.prop_initial_state, diag([1E6, 1E6, 1E6, 1E2, 1E2, 1E2]))
-# Specify non-zero step to get propagated states/covariances between measurement updates
-cfg.prop_step = 0.0
+# Use Filter.EXTENDED_KALMAN for the EKF
 cfg.estm_filter = Filter.UNSCENTED_KALMAN
 
 # Run OD on simulated measurements

orbdetpy/plotting/estimation.py

@@ -71,8 +71,8 @@ def plot(cfg: Settings, measurements, orbit_fit, interactive: bool=False,
                 parnames.append(f"{sk}-{meas_names[m]}")
 
     diag_pos = [0, 2, 5, 9, 14, 20]
-    has_truth = len(inp) > 0 and hasattr(inp[0], "true_state") and len(inp[0].true_state) >= 6
     tstamp,prefit,posfit,inocov,params,estmacc,estmcov,colors,state_err,state_cov = [],[],[],[],[],[],[],[],[],[]
+
     for i, o in zip(inp, out):
         tstamp.append(i.time)
         colors.append(cmap[o.station if len(o.station) else None])
@@ -83,7 +83,7 @@ def plot(cfg: Settings, measurements, orbit_fit, interactive: bool=False,
         else:
             inocov.append([0.0]*len(prefit[0]))
 
-        if (has_truth):
+        if (hasattr(i, "true_state") and len(i.true_state) >= 6):
             state_err.append([ix - ox for ix, ox in zip(i.true_state, o.estimated_state)])
             state_cov.append([3.0*numpy.sqrt(o.estimated_covariance[m]) for m in diag_pos if (m < len(o.estimated_covariance))])
 
@@ -103,9 +103,9 @@ def plot(cfg: Settings, measurements, orbit_fit, interactive: bool=False,
 
     if (len(tstamp) > 0):
         start = tstamp[0] if (tstamp[0] < tstamp[-1]) else tstamp[-1]
-        tim = [(t - start)/3600 for t in tstamp]
+        times = [(t - start)/3600 for t in tstamp]
     else:
-        tim = []
+        times = []
 
     angles = [MeasurementType.AZIMUTH, MeasurementType.ELEVATION, MeasurementType.RIGHT_ASCENSION, MeasurementType.DECLINATION]
     if (MeasurementType.POSITION in key):
@@ -136,7 +136,7 @@ def plot(cfg: Settings, measurements, orbit_fit, interactive: bool=False,
             plt.subplot(3, 2, order[i])
         else:
             plt.subplot(pre.shape[-1], 1, i + 1)
-        plt.scatter(tim, pre[:,i], color=colors, marker="o", s=7)
+        plt.scatter(times, pre[:,i], color=colors, marker="o", s=7)
         plt.legend(handles=patches, loc="best")
         plt.xlabel("Time [hr]")
         plt.ylabel(f"{ylabs[i]} [{units[i]}]")
@@ -156,10 +156,10 @@ def plot(cfg: Settings, measurements, orbit_fit, interactive: bool=False,
             plt.subplot(3, 2, order[i])
         else:
             plt.subplot(pos.shape[-1], 1, i + 1)
-        plt.scatter(tim, pos[:,i], color=colors, marker="o", s=7)
+        plt.scatter(times, pos[:,i], color=colors, marker="o", s=7)
         plt.legend(handles=patches, loc="best")
-        plt.plot(tim, -cov[:,i], "-r")
-        plt.plot(tim,  cov[:,i], "-r")
+        plt.plot(times, -cov[:,i], "-r")
+        plt.plot(times,  cov[:,i], "-r")
         plt.xlabel("Time [hr]")
         plt.ylabel(f"{ylabs[i]} [{units[i]}]")
         plt.grid(True)
@@ -175,7 +175,7 @@ def plot(cfg: Settings, measurements, orbit_fit, interactive: bool=False,
                 plt.figure(2)
                 plt.suptitle("Estimated Parameters")
             plt.subplot(par.shape[1], 1, i + 1)
-            plt.scatter(tim, par[:,i], marker="o", s=7)
+            plt.scatter(times, par[:,i], marker="o", s=7)
             plt.xlabel("Time [hr]")
             plt.ylabel(parnames[i])
             plt.grid(True)
@@ -192,7 +192,7 @@ def plot(cfg: Settings, measurements, orbit_fit, interactive: bool=False,
             for i in range(3):
                 plt.subplot(3, 1, i+1)
                 if (len(estmacc) > 0):
-                    plt.plot(tim, estmacc[:,i])
+                    plt.plot(times, estmacc[:,i])
                 plt.xlabel("Time [hr]")
                 plt.ylabel(lab[i])
                 plt.grid(True)
@@ -202,7 +202,7 @@ def plot(cfg: Settings, measurements, orbit_fit, interactive: bool=False,
                 outfiles.append(output_file_path + "_estacc.png")
                 plt.savefig(outfiles[-1], format="png")
 
-    if (has_truth):
+    if (len(times) == len(state_err)):
         state_err = numpy.array(state_err)
         state_cov = numpy.array(state_cov)
         order = [1, 3, 5, 2, 4, 6]
@@ -212,9 +212,9 @@ def plot(cfg: Settings, measurements, orbit_fit, interactive: bool=False,
         plt.suptitle("Position and Velocity Errors")
         for i in range(6):
             plt.subplot(3, 2, order[i])
-            plt.scatter(tim, state_err[:,i], color=colors, marker="o", s=7)
-            plt.plot(tim, -state_cov[:,i], "-r")
-            plt.plot(tim,  state_cov[:,i], "-r")
+            plt.scatter(times, state_err[:,i], color=colors, marker="o", s=7)
+            plt.plot(times, -state_cov[:,i], "-r")
+            plt.plot(times,  state_cov[:,i], "-r")
             plt.xlabel("Time [hr]")
             plt.ylabel(f"{ylabs[i]} [{units[i]}]")
             plt.grid(True)

orbdetpy/pom.xml

@@ -2,19 +2,17 @@
   <modelVersion>4.0.0</modelVersion>
   <groupId>org.astria</groupId>
   <artifactId>orbdetpy-server</artifactId>
-  <version>2.0.5</version>
+  <version>2.0.6</version>
   <packaging>jar</packaging>
 
   <properties>
     <application.main.class>org.astria.rpc.RPCServer</application.main.class>
     <os.detected.classifier>linux-x86_64</os.detected.classifier>
     <project.build.sourceEncoding>UTF-8</project.build.sourceEncoding>
     <project.reporting.outputEncoding>UTF-8</project.reporting.outputEncoding>
-    <javax.version>1.3.2</javax.version>
-    <grpc.version>1.35.0</grpc.version>
+    <grpc.version>1.36.0</grpc.version>
+    <javax.annotation.version>1.3.2</javax.annotation.version>
     <maven.compiler.version>3.8.1</maven.compiler.version>
-    <maven.compiler.source>1.8</maven.compiler.source>
-    <maven.compiler.target>1.8</maven.compiler.target>
     <maven.javadoc.version>3.2.0</maven.javadoc.version>
     <maven.shade.version>3.2.4</maven.shade.version>
     <maven.surefire.version>2.22.2</maven.surefire.version>
@@ -36,7 +34,7 @@
     <dependency>
       <groupId>javax.annotation</groupId>
       <artifactId>javax.annotation-api</artifactId>
-      <version>${javax.version}</version>
+      <version>${javax.annotation.version}</version>
     </dependency>
     <dependency>
       <groupId>io.grpc</groupId>
@@ -89,6 +87,7 @@
         <version>${maven.compiler.version}</version>
         <configuration>
 	  <compilerArgument>-Xlint:unchecked</compilerArgument>
+	  <release>11</release>
 	  <showDeprecation>true</showDeprecation>
 	  <showWarnings>true</showWarnings>
         </configuration>

orbdetpy/rpc/server.py

@@ -14,9 +14,9 @@
 # You should have received a copy of the GNU General Public License
 # along with this program.  If not, see <http://www.gnu.org/licenses/>.
 
-from os import path
 from socket import socket
 from atexit import register
+from os import environ, path
 from psutil import Popen, process_iter
 from grpc import channel_ready_future, insecure_channel
 
@@ -41,7 +41,15 @@ def connect(cls, data_dir, jar_file):
             sock.bind((rpc_host, 0))
             rpc_port = f"{sock.getsockname()[1]}"
             sock.close()
-            cls.rpc_server = Popen(["java", "-Xmx2G", "-XX:+UseG1GC", "-jar", jar_file, rpc_port, data_dir])
+
+            java_home = environ.get("JAVA_HOME")
+            if (java_home and path.isdir(path.join(java_home, "bin"))):
+                java_exec = path.join(java_home, "bin", "java")
+            else:
+                java_exec = "java"
+
+            cls.rpc_server = Popen([java_exec, "-Xmx2G", "-XX:+UseG1GC", "-jar", jar_file, rpc_port, data_dir])
+
         rpc_uri = f"{rpc_host}:{rpc_port}"
         cls.rpc_channel = insecure_channel(rpc_uri, options=[("grpc.max_send_message_length", 2147483647),
                                                              ("grpc.max_receive_message_length", 2147483647)])

orbdetpy/src/main/java/org/astria/DataManager.java

@@ -1,6 +1,6 @@
 /*
  * DataManager.java - Functions for handling data files.
- * Copyright (C) 2018-2020 University of Texas
+ * Copyright (C) 2018-2021 University of Texas
  *
  * This program is free software: you can redistribute it and/or modify
  * it under the terms of the GNU General Public License as published by
@@ -19,38 +19,21 @@
 package org.astria;
 
 import java.io.File;
-import java.nio.file.Files;
-import java.nio.file.Path;
-import java.nio.file.Paths;
 import java.util.HashMap;
-import java.util.Scanner;
-import java.util.TreeMap;
 import java.util.concurrent.Executors;
 import java.util.concurrent.ExecutorService;
-import java.util.regex.Pattern;
 import org.orekit.bodies.OneAxisEllipsoid;
 import org.orekit.data.DataContext;
 import org.orekit.data.DirectoryCrawler;
 import org.orekit.frames.FramesFactory;
 import org.orekit.frames.Predefined;
 import org.orekit.time.AbsoluteDate;
-import org.orekit.time.DateTimeComponents;
-import org.orekit.time.TimeScalesFactory;
 import org.orekit.utils.Constants;
 
 public final class DataManager
 {
-    public static class MSISEData
-    {
-	public AbsoluteDate minDate;
-	public AbsoluteDate maxDate;
-	public HashMap<String, double[]> data;
-    }
-
-    private static String dataPath;
+    protected static String dataPath;
     protected static ExecutorService threadPool;
-    protected static MSISEData msiseData;
-    protected static TreeMap<Double, String> wamFileMap;
     protected static OneAxisEllipsoid earthShape;
     protected static HashMap<Integer, AbsoluteDate> epochs;
 
@@ -66,8 +49,6 @@ public static synchronized void initialize(String dataPath) throws Exception
 	DataManager.earthShape = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS, Constants.WGS84_EARTH_FLATTENING,
 						      FramesFactory.getFrame(Predefined.ITRF_CIO_CONV_2010_ACCURATE_EOP));
 	loadEpochs();
-	loadMSISEData();
-	loadWAMData();
     }
 
     public static synchronized void shutdown()
@@ -96,58 +77,6 @@ private static void loadEpochs()
 	epochs.put(11, AbsoluteDate.QZSS_EPOCH);
     }
 
-    private static void loadMSISEData() throws Exception
-    {
-	String[] toks = null;
-	DataManager.msiseData = new MSISEData();
-	DataManager.msiseData.data = new HashMap<String, double[]>();
-	Scanner scan = new Scanner(new File(dataPath, "SpaceWeather.dat"));
-	while (scan.hasNextLine())
-	{
-	    toks = scan.nextLine().split(",");
-	    double[] vals = new double[toks.length];
-	    for (int i = 0; i < toks.length; i++)
-		toks[i] = toks[i].trim();
-
-	    if (DataManager.msiseData.minDate == null)
-		DataManager.msiseData.minDate = new AbsoluteDate(Integer.parseInt(toks[0]), Integer.parseInt(toks[1]),
-								 Integer.parseInt(toks[2]), TimeScalesFactory.getUTC());
-	    for (int i = 0; i < toks.length; i++)
-	    {
-		if (toks[i].length() > 0)
-		    vals[i] = Double.parseDouble(toks[i]);
-	    }
-	    DataManager.msiseData.data.put(String.format("%s%s%s", toks[0], toks[1], toks[2]), vals);
-	}
-
-	scan.close();
-	if (toks != null)
-	    DataManager.msiseData.maxDate = new AbsoluteDate(Integer.parseInt(toks[0]), Integer.parseInt(toks[1]),
-							     Integer.parseInt(toks[2]), TimeScalesFactory.getUTC());
-    }
-
-    private static void loadWAMData() throws Exception
-    {
-	Path wamPath = Paths.get(dataPath, "WAM");
-	DataManager.wamFileMap = new TreeMap<Double, String>();
-	if (!wamPath.toFile().isDirectory())
-	    return;
-	Files.find(wamPath, Integer.MAX_VALUE, (p, a) -> a.isRegularFile() && p.toString().endsWith(".nc"))
-	    .forEach(cdf -> {
-		    for (String s: cdf.getFileName().toString().split(Pattern.quote(".")))
-		    {
-			if (!s.startsWith("20") || s.length() != 15)
-			    continue;
-			s = s.replace("_", "");
-			String utc = String.format("%s-%s-%sT%s:%s:%sZ", s.substring(0, 4), s.substring(4, 6), s.substring(6, 8),
-						   s.substring(8, 10), s.substring(10, 12), s.substring(12, 14));
-			DataManager.wamFileMap.put(new AbsoluteDate(DateTimeComponents.parseDateTime(utc), TimeScalesFactory.getUTC())
-						   .durationFrom(AbsoluteDate.J2000_EPOCH), cdf.toString());
-			break;
-		    }
-		});
-    }
-
     public static AbsoluteDate getEpoch(int epoch)
     {
 	return(epochs.get(epoch));