Skip to content

Commit

Permalink
Update scenarios for the new logging syntax
Browse files Browse the repository at this point in the history
  • Loading branch information
juan-g-bonilla committed Sep 12, 2023
1 parent d212729 commit e8eed16
Show file tree
Hide file tree
Showing 5 changed files with 71 additions and 65 deletions.
Original file line number Diff line number Diff line change
Expand Up @@ -143,10 +143,10 @@ def log_outputs(self):
self.rwLogs.append(DynModel.rwStateEffector.rwOutMsgs[item].recorder(samplingTime))
self.masterSim.AddModelToTask(DynModel.taskName, self.rwLogs[item])

self.masterSim.AddVariableForLogging('headingUKF.bVec_B', samplingTime, 0, 2)

self.headingBVecLog = FswModel.headingUKF.logger("bVec_B")
self.filtRec = FswModel.headingUKF.filtDataOutMsg.recorder(samplingTime)
self.opNavFiltRec = FswModel.headingUKF.opnavDataOutMsg.recorder(samplingTime)
self.masterSim.AddModelToTask(DynModel.taskName, self.headingBVecLog)
self.masterSim.AddModelToTask(DynModel.taskName, self.filtRec)
self.masterSim.AddModelToTask(DynModel.taskName, self.opNavFiltRec)

Expand All @@ -167,7 +167,7 @@ def pull_outputs(self, showPlots):
circleRadii = unitTestSupport.addTimeColumn(self.circlesRec.times(), self.circlesRec.circlesRadii)
validCircle = unitTestSupport.addTimeColumn(self.circlesRec.times(), self.circlesRec.valid)

frame = self.masterSim.GetLogVariableData('headingUKF.bVec_B')
frame = unitTestSupport.addTimeColumn(self.headingBVecLog.times(), self.headingBVecLog.bVec_B)

numRW = 4
dataRW = []
Expand Down
12 changes: 8 additions & 4 deletions examples/scenarioCSSFilters.py
Original file line number Diff line number Diff line change
Expand Up @@ -483,6 +483,7 @@ def setupCSS(CSS):
# Setup filter
#
numStates = 6
bVecLogger = None
if FilterType == 'EKF':
module = sunlineEKF.sunlineEKF()
module.ModelTag = "SunlineEKF"
Expand Down Expand Up @@ -518,7 +519,8 @@ def setupCSS(CSS):

# Add test module to runtime call list
scSim.AddModelToTask(simTaskName, module)
scSim.AddVariableForLogging('SunlineSEKF.bVec_B', simulationTimeStep, 0, 2)
bVecLogger = module.logger('bVec_B', simulationTimeStep)
scSim.AddModelToTask(simTaskName, bVecLogger)

if FilterType == 'SuKF':
numStates = 6
Expand All @@ -528,7 +530,8 @@ def setupCSS(CSS):

# Add test module to runtime call list
scSim.AddModelToTask(simTaskName, module)
scSim.AddVariableForLogging('SunlineSuKF.bVec_B', simulationTimeStep, 0, 2)
bVecLogger = module.logger('bVec_B', simulationTimeStep)
scSim.AddModelToTask(simTaskName, bVecLogger)

module.cssDataInMsg.subscribeTo(cssConstelation.constellationOutMsg)
module.cssConfigInMsg.subscribeTo(cssConstMsg)
Expand Down Expand Up @@ -570,13 +573,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 = []
Expand All @@ -588,7 +593,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 = []
Expand Down
8 changes: 5 additions & 3 deletions examples/scenarioDragDeorbit.py
Original file line number Diff line number Diff line change
Expand Up @@ -250,9 +250,11 @@ 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,
Expand Down Expand Up @@ -282,7 +284,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)
Expand Down Expand Up @@ -348,7 +350,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]')

Expand Down
35 changes: 18 additions & 17 deletions examples/scenarioDragRendezvous.py
Original file line number Diff line number Diff line change
Expand Up @@ -255,10 +255,6 @@ def drag_simulator(altOffset, trueAnomOffset, densMultiplier, ctrlType='lqr', us
scSim.AddModelToTask(dynTaskName, chiefNav,800)
scSim.AddModelToTask(dynTaskName, depNav,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
Expand Down Expand Up @@ -306,12 +302,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 = hillStateNavObj.hillStateOutMsg.recorder()
depAttRec = depAttRef.attRefOutMsg.recorder()
chiefAttRec = chiefAttRef.attRefOutMsg.recorder()
chiefDragForceLog = chiefDrag.logger("forceExternal_B")
depDragForceLog = depDrag.logger("forceExternal_B")
atmoRecs = []
for msg in atmosphere.envOutMsgs:
atmoRec = msg.recorder()
Expand All @@ -322,8 +320,11 @@ def drag_simulator(altOffset, trueAnomOffset, densMultiplier, ctrlType='lqr', us
scSim.AddModelToTask(dynTaskName, hillStateRec, 702)
scSim.AddModelToTask(dynTaskName, depAttRec, 703)
scSim.AddModelToTask(dynTaskName, chiefAttRec, 704)
scSim.AddModelToTask(dynTaskName, chiefDragForceLog, 705)
scSim.AddModelToTask(dynTaskName, depDragForceLog, 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
Expand All @@ -345,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
Expand Down Expand Up @@ -453,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')
Expand Down
75 changes: 37 additions & 38 deletions examples/scenarioFuelSlosh.py
Original file line number Diff line number Diff line change
Expand Up @@ -184,6 +184,7 @@
from Basilisk.utilities import orbitalMotion
from Basilisk.utilities import simIncludeGravBody
from Basilisk.utilities import unitTestSupport
from Basilisk.utilities import pythonVariableLogger

filename = inspect.getframeinfo(inspect.currentframe()).filename
path = os.path.dirname(os.path.abspath(filename))
Expand Down Expand Up @@ -319,46 +320,47 @@ 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 = pythonVariableLogger.PythonVariableLogger({
"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

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:
Expand Down Expand Up @@ -401,43 +403,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)')
Expand All @@ -451,7 +450,7 @@ def run(show_plots, damping_parameter, timeStep):
plt.close("all")


return rhoj1Out, rhoj2Out, rhoj3Out, figureList
return time, rhojOuts, figureList


if __name__ == "__main__":
Expand Down

0 comments on commit e8eed16

Please sign in to comment.