Update scenarios to use new logging syntax

examples/scenarioCSSFilters.py

@@ -483,6 +483,7 @@ def setupCSS(CSS):
     # Setup filter
     #
     numStates = 6
+    bVecLogger = None
     if FilterType == 'EKF':
         moduleConfig = sunlineEKF.sunlineEKFConfig()
         moduleWrap = scSim.setModelDataWrap(moduleConfig)
@@ -522,7 +523,8 @@ def setupCSS(CSS):
 
         # Add test module to runtime call list
         scSim.AddModelToTask(simTaskName, moduleWrap, moduleConfig)
-        scSim.AddVariableForLogging('SunlineSEKF.bVec_B', simulationTimeStep, 0, 2)
+        bVecLogger = moduleWrap.logger('bVec_B', simulationTimeStep)
+        scSim.AddModelToTask(simTaskName, bVecLogger)
 
     if FilterType == 'SuKF':
         numStates = 6
@@ -533,7 +535,8 @@ def setupCSS(CSS):
 
         # Add test module to runtime call list
         scSim.AddModelToTask(simTaskName, moduleWrap, moduleConfig)
-        scSim.AddVariableForLogging('SunlineSuKF.bVec_B', simulationTimeStep, 0, 2)
+        bVecLogger = moduleWrap.logger('bVec_B', simulationTimeStep)
+        scSim.AddModelToTask(simTaskName, bVecLogger)
 
     moduleConfig.cssDataInMsg.subscribeTo(cssConstelation.constellationOutMsg)
     moduleConfig.cssConfigInMsg.subscribeTo(cssConstMsg)
@@ -575,13 +578,15 @@ def addTimeColumn(time, data):
     covarLog = addTimeColumn(timeAxis, filtLog.covar[:, range(numStates*numStates)])
     obsLog = addTimeColumn(timeAxis, filtLog.numObs)
 
+    # Get bVec_B through the variable logger
+    bVecLog = None if bVecLogger is None else addTimeColumn(timeAxis, bVecLogger.bVec_B)
+
     dcmLog = np.zeros([len(stateLog[:,0]),3,3])
     omegaExp = np.zeros([len(stateLog[:,0]),3])
     if FilterType == 'SEKF':
         dcm = sunlineSEKF.new_doubleArray(3 * 3)
         for j in range(9):
             sunlineSEKF.doubleArray_setitem(dcm, j, 0)
-        bVecLog = scSim.GetLogVariableData('SunlineSEKF.bVec_B')
         for i in range(len(stateLog[:,0])):
             sunlineSEKF.sunlineSEKFComputeDCM_BS(stateLog[i,1:4].tolist(), bVecLog[i, 1:4].tolist(), dcm)
             dcmOut = []
@@ -593,7 +598,6 @@ def addTimeColumn(time, data):
         dcm = sunlineSuKF.new_doubleArray(3 * 3)
         for j in range(9):
             sunlineSuKF.doubleArray_setitem(dcm, j, 0)
-        bVecLog = scSim.GetLogVariableData('SunlineSuKF.bVec_B')
         for i in range(len(stateLog[:,0])):
             sunlineSuKF.sunlineSuKFComputeDCM_BS(stateLog[i,1:4].tolist(), bVecLog[i, 1:4].tolist(), dcm)
             dcmOut = []
@@ -650,6 +654,6 @@ def addTimeColumn(time, data):
 if __name__ == "__main__":
     run(False,       # save figures to file
         True,      # show_plots
-        'SuKF',
+        'SEKF',
          400
        )

examples/scenarioDragDeorbit.py

@@ -250,9 +250,10 @@ def run(show_plots, initialAlt=250, deorbitAlt=100, model="exponential"):
     # Setup data logging before the simulation is initialized
     dataRec = scObject.scStateOutMsg.recorder(samplingTime)
     dataAtmoLog = atmo.envOutMsgs[0].recorder(samplingTime)
+    forceLog = dragEffector.logger("forceExternal_B", samplingTime)
     scSim.AddModelToTask(simTaskName, dataRec)
     scSim.AddModelToTask(simTaskName, dataAtmoLog)
-    scSim.AddVariableForLogging('DragEff.forceExternal_B', samplingTime, StartIndex=0, StopIndex=2)
+    scSim.AddModelToTask(simTaskName, forceLog)
 
     # Event to terminate the simulation
     scSim.createNewEvent("Deorbited", simulationTimeStep, True,
@@ -282,7 +283,7 @@ def run(show_plots, initialAlt=250, deorbitAlt=100, model="exponential"):
     # retrieve the logged data
     posData = dataRec.r_BN_N
     velData = dataRec.v_BN_N
-    dragForce = scSim.GetLogVariableData('DragEff.forceExternal_B')
+    dragForce = forceLog.forceExternal_B
     denseData = dataAtmoLog.neutralDensity
 
     figureList = plotOrbits(dataRec.times(), posData, velData, dragForce, denseData, oe, mu, planet, model)
@@ -348,7 +349,7 @@ def register_fig(i):
 
     # draw drag as a function of time
     fig, ax = register_fig(4)
-    plt.semilogy(timeAxis * macros.NANO2HOUR, np.linalg.norm(dragForce[:, 1:], 2, 1))
+    plt.semilogy(timeAxis * macros.NANO2HOUR, np.linalg.norm(dragForce, 2, 1))
     plt.xlabel('$t$ [hr]')
     plt.ylabel('$|F_drag|$ [N]')
 

examples/scenarioDragRendezvous.py

@@ -252,10 +252,6 @@ def drag_simulator(altOffset, trueAnomOffset, densMultiplier, ctrlType='lqr', us
     scSim.AddModelToTask(dynTaskName, chiefNav,ModelPriority=800)
     scSim.AddModelToTask(dynTaskName, depNav,ModelPriority=801)
     # scSim.AddModelToTask(dynTaskName, gravFactory.spiceObject, 999)
-    scSim.AddVariableForLogging(chiefDrag.ModelTag + ".forceExternal_B",
-                                      dynTimeStep, 0, 2, 'double')
-    scSim.AddVariableForLogging(depDrag.ModelTag + ".forceExternal_B",
-                                      dynTimeStep, 0, 2, 'double')
 
     ##  FSW setup
     #   Chief S/C
@@ -307,12 +303,14 @@ def drag_simulator(altOffset, trueAnomOffset, densMultiplier, ctrlType='lqr', us
     numDataPoints = 21384
     samplingTime = simulationTime // (numDataPoints - 1)
 
-    #   Set up recorders
+    #   Set up recorders and loggers
     chiefStateRec = chiefSc.scStateOutMsg.recorder()
     depStateRec = depSc.scStateOutMsg.recorder()
     hillStateRec = hillStateNavData.hillStateOutMsg.recorder()
     depAttRec = depAttRefData.attRefOutMsg.recorder()
     chiefAttRec = chiefAttRefData.attRefOutMsg.recorder()
+    chiefDragForceLog = chiefDrag.logger("forceExternal_B")
+    depDragForceLog = depDrag.logger("forceExternal_B")
     atmoRecs = []
     for msg in atmosphere.envOutMsgs:
         atmoRec = msg.recorder()
@@ -323,8 +321,10 @@ def drag_simulator(altOffset, trueAnomOffset, densMultiplier, ctrlType='lqr', us
     scSim.AddModelToTask(dynTaskName, hillStateRec, ModelPriority=702)
     scSim.AddModelToTask(dynTaskName, depAttRec, ModelPriority=703)
     scSim.AddModelToTask(dynTaskName, chiefAttRec, ModelPriority=704)
+    scSim.AddModelToTask(dynTaskName, chiefDragForceLog, ModelPriority=705)
+    scSim.AddModelToTask(dynTaskName, depDragForceLog, ModelPriority=706)
     for ind,rec in enumerate(atmoRecs):
-        scSim.AddModelToTask(dynTaskName, rec, 705+ind)
+        scSim.AddModelToTask(dynTaskName, rec, 707+ind)
 
     # if this scenario is to interface with the BSK Viz, uncomment the following lines
     # to save the BSK data to a file, uncomment the saveFile line below
@@ -346,8 +346,8 @@ def drag_simulator(altOffset, trueAnomOffset, densMultiplier, ctrlType='lqr', us
     print(f"Sim complete in {execTime} seconds, pulling data...")
     # ----- pull ----- #
     results_dict = {}
-    results_dict['chiefDrag_B'] = scSim.GetLogVariableData(chiefDrag.ModelTag + ".forceExternal_B")
-    results_dict['depDrag_B'] = scSim.GetLogVariableData(depDrag.ModelTag + ".forceExternal_B")
+    results_dict['chiefDrag_B'] = chiefDragForceLog.forceExternal_B
+    results_dict['depDrag_B'] = depDragForceLog.forceExternal_B
     results_dict['dynTimeData'] = chiefStateRec.times()
     results_dict['fswTimeData'] = depAttRec.times()
     results_dict['wiggum.r_BN_N'] = chiefStateRec.r_BN_N 
@@ -454,27 +454,27 @@ def run(show_plots, altOffset, trueAnomOffset, densMultiplier, ctrlType='lqr', u
 
     #   Debug plots
     plt.figure()
-    plt.plot(timeData[1:], depDrag[1:,1]-chiefDrag[1:,1],label="delta a_1")
-    plt.plot(timeData[1:], depDrag[1:,2]-chiefDrag[1:,2],label="delta a_2")
-    plt.plot(timeData[1:], depDrag[1:,3]-chiefDrag[1:,3],label="delta a_3")
+    plt.plot(timeData[1:], depDrag[1:,0]-chiefDrag[1:,0],label="delta a_1")
+    plt.plot(timeData[1:], depDrag[1:,1]-chiefDrag[1:,1],label="delta a_2")
+    plt.plot(timeData[1:], depDrag[1:,2]-chiefDrag[1:,2],label="delta a_3")
     plt.grid()
     plt.legend()
     plt.xlabel('Time')   
     plt.ylabel('Relative acceleration due to drag, body frame (m/s)')
 
     plt.figure()
-    plt.plot(timeData[1:], chiefDrag[1:,1],label="chief a_1")
-    plt.plot(timeData[1:], chiefDrag[1:,2],label="chief a_2")
-    plt.plot(timeData[1:], chiefDrag[1:,3],label="chief a_3")
+    plt.plot(timeData[1:], chiefDrag[1:,0],label="chief a_1")
+    plt.plot(timeData[1:], chiefDrag[1:,1],label="chief a_2")
+    plt.plot(timeData[1:], chiefDrag[1:,2],label="chief a_3")
     plt.grid()
     plt.legend()
     plt.xlabel('Time')   
     plt.ylabel('Relative acceleration due to drag, body frame (m/s)')
 
     plt.figure()
-    plt.plot(timeData[1:], depDrag[1:,1],label="dep a_1")
-    plt.plot(timeData[1:], depDrag[1:,2],label="dep a_2")
-    plt.plot(timeData[1:], depDrag[1:,3],label="dep a_3")
+    plt.plot(timeData[1:], depDrag[1:,0],label="dep a_1")
+    plt.plot(timeData[1:], depDrag[1:,1],label="dep a_2")
+    plt.plot(timeData[1:], depDrag[1:,2],label="dep a_3")
     plt.grid()
     plt.legend()
     plt.xlabel('Time')   

examples/scenarioFuelSlosh.py

@@ -184,6 +184,7 @@
 from Basilisk.utilities import orbitalMotion
 from Basilisk.utilities import simIncludeGravBody
 from Basilisk.utilities import unitTestSupport
+from Basilisk.utilities import variable_logger
 
 filename = inspect.getframeinfo(inspect.currentframe()).filename
 path = os.path.dirname(os.path.abspath(filename))
@@ -319,46 +320,52 @@ def run(show_plots, damping_parameter, timeStep):
     dataLog = scObject.scStateOutMsg.recorder(samplingTime)
     scSim.AddModelToTask(simTaskName, dataLog)
 
-    scSim.InitializeSimulation()
+    scLog = scObject.logger(
+        ["totOrbEnergy", "totOrbAngMomPntN_N", "totRotAngMomPntC_N", "totRotEnergy"], 
+        simulationTimeStep)
+    scSim.AddModelToTask(simTaskName, scLog)
 
-    scSim.AddVariableForLogging(scObject.ModelTag + ".totOrbEnergy", simulationTimeStep, 0, 0, 'double')
-    scSim.AddVariableForLogging(scObject.ModelTag + ".totOrbAngMomPntN_N", simulationTimeStep, 0, 2, 'double')
-    scSim.AddVariableForLogging(scObject.ModelTag + ".totRotAngMomPntC_N", simulationTimeStep, 0, 2, 'double')
-    scSim.AddVariableForLogging(scObject.ModelTag + ".totRotEnergy", simulationTimeStep, 0, 0, 'double')
+    damperRhoLog = None
     if damping_parameter != 0.0:
-        scSim.AddVariableForLogging(
-            scObject.ModelTag +
-            ".dynManager.getStateObject('linearSpringMassDamperRho1').getState()", simulationTimeStep, 0, 0, 'double')
-        scSim.AddVariableForLogging(
-            scObject.ModelTag +
-            ".dynManager.getStateObject('linearSpringMassDamperRho2').getState()", simulationTimeStep, 0, 0, 'double')
-        scSim.AddVariableForLogging(
-            scObject.ModelTag +
-            ".dynManager.getStateObject('linearSpringMassDamperRho3').getState()", simulationTimeStep, 0, 0, 'double')
+        def get_rho(CurrentSimNanos, i):
+            stateName = f'linearSpringMassDamperRho{i}'
+            return scObject.dynManager.getStateObject(stateName).getState()[0][0]
+
+        damperRhoLog = variable_logger.VariableLogger({
+                "damperRho1": lambda CurrentSimNanos: get_rho(CurrentSimNanos, 1),
+                "damperRho2": lambda CurrentSimNanos: get_rho(CurrentSimNanos, 2),
+                "damperRho3": lambda CurrentSimNanos: get_rho(CurrentSimNanos, 3),
+            }, simulationTimeStep)
+
+        scSim.AddModelToTask(simTaskName, damperRhoLog)
 
     #
     #   configure a simulation stop time and execute the simulation run
     #
+    scSim.InitializeSimulation()
+
     scSim.ConfigureStopTime(simulationTime)
     scSim.ExecuteSimulation()
 
     # request states to the simulation
     posData = dataLog.r_CN_N
     velData = dataLog.v_CN_N
 
-    orbEnergy = scSim.GetLogVariableData(scObject.ModelTag + ".totOrbEnergy")
-    orbAngMom_N = scSim.GetLogVariableData(scObject.ModelTag + ".totOrbAngMomPntN_N")
-    rotAngMom_N = scSim.GetLogVariableData(scObject.ModelTag + ".totRotAngMomPntC_N")
-    rotEnergy = scSim.GetLogVariableData(scObject.ModelTag + ".totRotEnergy")
+    time = scLog.times() * 1e-9
+    orbEnergy = scLog.totOrbEnergy
+    orbAngMom_N = scLog.totOrbAngMomPntN_N
+    rotAngMom_N = scLog.totRotAngMomPntC_N
+    rotEnergy = scLog.totRotEnergy
+
+    print("orbEnergy", orbEnergy)
+    print("orbAngMom_N", orbAngMom_N)
+    print("rotAngMom_N", rotAngMom_N)
+    print("rotEnergy", rotEnergy)
 
-    rhoj1Out = rhoj2Out = rhoj3Out = []
+    rhojOuts = None
     if damping_parameter != 0.0:
-        rhoj1Out = scSim.GetLogVariableData(
-            scObject.ModelTag + ".dynManager.getStateObject('linearSpringMassDamperRho1').getState()")
-        rhoj2Out = scSim.GetLogVariableData(
-            scObject.ModelTag + ".dynManager.getStateObject('linearSpringMassDamperRho2').getState()")
-        rhoj3Out = scSim.GetLogVariableData(
-            scObject.ModelTag + ".dynManager.getStateObject('linearSpringMassDamperRho3').getState()")
+        rhojOuts = [
+            damperRhoLog.damperRho1, damperRhoLog.damperRho2, damperRhoLog.damperRho3]
 
     fileName = os.path.basename(os.path.splitext(__file__)[0])
     if damping_parameter == 0.0 and timeStep == 0.75:
@@ -401,43 +408,40 @@ def run(show_plots, damping_parameter, timeStep):
     figureList[pltName] = plt.figure(1)
 
     plt.figure(2, figsize=(5, 4))
-    plt.plot(orbAngMom_N[:, 0] * 1e-9, (orbAngMom_N[:, 1] - orbAngMom_N[0, 1]) / orbAngMom_N[0, 1],
-             orbAngMom_N[:, 0] * 1e-9,
-             (orbAngMom_N[:, 2] - orbAngMom_N[0, 2]) / orbAngMom_N[0, 2],
-             orbAngMom_N[:, 0] * 1e-9, (orbAngMom_N[:, 3] - orbAngMom_N[0, 3]) / orbAngMom_N[0, 3])
+    for i in range(3):
+        plt.plot(time, (orbAngMom_N[:, i] - orbAngMom_N[0, i]) / orbAngMom_N[0, i])
     plt.xlabel('Time (s)')
     plt.ylabel('Relative Orbital Angular Momentum Variation')
     pltName = fileName + "OAM" + str(setupNo)
     figureList[pltName] = plt.figure(2)
 
     plt.figure(3, figsize=(5, 4))
-    plt.plot(orbEnergy[:, 0] * 1e-9, (orbEnergy[:, 1] - orbEnergy[0, 1]) / orbEnergy[0, 1])
+    plt.plot(time, (orbEnergy - orbEnergy[0]) / orbEnergy[0])
     plt.xlabel('Time (s)')
     plt.ylabel('Relative Orbital Energy Variation')
     pltName = fileName + "OE" + str(setupNo)
     figureList[pltName] = plt.figure(3)
 
     plt.figure(4, figsize=(5, 4))
-    plt.plot(rotAngMom_N[:, 0] * 1e-9,
-             (rotAngMom_N[:, 1] - rotAngMom_N[0, 1]) / rotAngMom_N[0, 1],
-             rotAngMom_N[:, 0] * 1e-9, (rotAngMom_N[:, 2] - rotAngMom_N[0, 2]) / rotAngMom_N[0, 2],
-             rotAngMom_N[:, 0] * 1e-9, (rotAngMom_N[:, 3] - rotAngMom_N[0, 3]) / rotAngMom_N[0, 3])
+    for i in range(3):
+        plt.plot(time, (rotAngMom_N[:, i] - rotAngMom_N[0, i]) / rotAngMom_N[0, i])
     plt.xlabel('Time (s)')
     plt.ylabel('Relative Rotational Angular Momentum Variation')
     pltName = fileName + "RAM" + str(setupNo)
     figureList[pltName] = plt.figure(4)
 
     plt.figure(5, figsize=(5, 4))
-    plt.plot(rotEnergy[:, 0] * 1e-9, (rotEnergy[:, 1] - rotEnergy[0, 1]) / rotEnergy[0, 1])
+    plt.plot(time, (rotEnergy - rotEnergy[0]) / rotEnergy[0])
     plt.xlabel('Time (s)')
     plt.ylabel('Relative Rotational Energy Variation')
     pltName = fileName + "RE" + str(setupNo)
     figureList[pltName] = plt.figure(5)
 
     if damping_parameter != 0.0:
         plt.figure(6, figsize=(5, 4))
-        plt.plot(rhoj1Out[:, 0] * 1e-9, rhoj1Out[:, 1], rhoj2Out[:, 0] *
-                 1e-9, rhoj2Out[:, 1], rhoj3Out[:, 0] * 1e-9, rhoj3Out[:, 1])
+        for rhojOut in rhojOuts:
+            plt.plot(time, rhojOut)
+
         plt.legend(['Particle 1', 'Particle 2', 'Particle 3'], loc='lower right')
         plt.xlabel('Time (s)')
         plt.ylabel('Displacement (m)')
@@ -451,12 +455,12 @@ def run(show_plots, damping_parameter, timeStep):
     plt.close("all")
 
 
-    return rhoj1Out, rhoj2Out, rhoj3Out, figureList
+    return time, rhojOuts, figureList
 
 
 if __name__ == "__main__":
     run(
         True,               # show_plots
-        0.0,				 # damping_parameter
+        15.0,				 # damping_parameter
         0.75,				 # timeStep
     )

src/tests/test_scenarioFuelSlosh.py

@@ -46,36 +46,23 @@ def test_scenarioFuelSlosh(show_plots, damping_parameter, timeStep):
     testFailCount = 0  # zero unit test result counter
     testMessages = []  # create empty array to store test log messages
 
-    rhoj1Out, rhoj2Out, rhoj3Out, figureList = \
-        scenarioFuelSlosh.run(show_plots, damping_parameter, timeStep)
+    time, rhojOuts, figureList = scenarioFuelSlosh.run(
+        show_plots, damping_parameter, timeStep)
 
     if damping_parameter != 0:
-
-        zita1 = damping_parameter / (2 * np.sqrt(1500.0 * 1.0))
-        omegan1 = np.sqrt(1.0 / 1500.0)
-        settling_time1 = -1.0 / (zita1 * omegan1) * np.log(0.05 * np.sqrt(1 - zita1**2))
-        index_settling_time1 = np.argmax(rhoj1Out[:, 0] * 1e-9 > settling_time1)
-
-        zita2 = damping_parameter / (2 * np.sqrt(1400.0 * 1.0))
-        omegan2 = np.sqrt(1.0 / 1400.0)
-        settling_time2 = -1.0 / (zita2 * omegan2) * np.log(0.05 * np.sqrt(1 - zita2**2))
-        index_settling_time2 = np.argmax(rhoj2Out[:, 0] * 1e-9 > settling_time2)
-
-        zita3 = damping_parameter / (2 * np.sqrt(1300.0 * 1.0))
-        omegan3 = np.sqrt(1.0 / 1300.0)
-        settling_time3 = -1.0 / (zita3 * omegan3) * np.log(0.05 * np.sqrt(1 - zita3**2))
-        index_settling_time3 = np.argmax(rhoj3Out[:, 0] * 1e-9 > settling_time3)
-
         accuracy = 0.05
-        if abs(rhoj1Out[index_settling_time1, 1] - rhoj1Out[-1, 1]) > accuracy:
-            testFailCount = testFailCount + 1
-            testMessages = [testMessages, "Particle 1 settling time does not match second order systems theories"]
-        if abs(rhoj2Out[index_settling_time2, 1] - rhoj2Out[-1, 1]) > accuracy:
-            testFailCount = testFailCount + 1
-            testMessages = [testMessages, "Particle 1 settling time does not match second order systems theories"]
-        if abs(rhoj3Out[index_settling_time3, 1] - rhoj3Out[-1, 1]) > accuracy:
-            testFailCount = testFailCount + 1
-            testMessages = [testMessages, "Particle 1 settling time does not match second order systems theories"]
+        mass = (1500, 1400, 1300)
+
+        for i in range(3):
+            zita = damping_parameter / (2 * np.sqrt(mass[i] * 1.0))
+            omegan = np.sqrt(1.0 / mass[i])
+            settling_time = -1.0 / (zita * omegan) * np.log(0.05 * np.sqrt(1 - zita**2))
+            index_settling_time = np.argmax(time > settling_time)
+
+            if abs(rhojOuts[i][index_settling_time] - rhojOuts[i][-1]) > accuracy:
+                testFailCount = testFailCount + 1
+                testMessages = [testMessages, f"Particle {i+1} settling time does "
+                                "not match second order systems theories"]
 
     # save the figures to the Doxygen scenario images folder
     for pltName, plt in list(figureList.items()):