Add Walker Delta constellation setup scenario

[squash] fixed copyright year

[squash] fix typos

[squash] remove un-needed import

[squash] removed un-needed comments

[squash] make simulation run faster

needed to avoid having pytest take too long

[squash] make sure the vizInterface code only runs if vizInterface is built

also, comment out the saveFile argument as we don't want to safe off sim data by default

[squash] fix a matplotlib warning about the color sequence

[squash] The mercator map is not part of the repo.  I'm removing including it in the figure.

This code was not running as committed.

added unit test for scenario script

[squash] make sure the plots are closed first before plotting

This is important if these scripts are called in a squential manner as prior plotting can leak into this script

[squash] make sure the plotted figures are return

They must be returned as a dictionary list so the pytest can save them off.

[squash] add the missing Vizard image (small jpg) to make the documentation fully build

[squash] Expand documentation to start by defining Walker numbers

add Walker Delta constellation setup scenario#

docs/source/Support/bskReleaseNotes.rst

@@ -32,6 +32,7 @@ Basilisk Release Notes
 
 Version |release|
 -----------------
+- Created a new example scenario :ref:`scenarioSatelliteConstellation` demonstrating setup of a Walker-Delta constellation
 - Created a new :ref:`pinholeCamera` module to support generation of landmarks-based measurements around a
   small body.
 - Corrected a memory leak in the ``swig`` access to standard vectors inside messages.

examples/_default.rst

@@ -294,6 +294,7 @@ Stand Alone Architecture
    Electrostatic Tractor Debris Reorbiting <scenarioDebrisReorbitET>
    Attitude-Driven differential drag control <scenarioDragRendezvous>
    Servicer approaching a debris object with 3 flight modes <scenarioRendezVous>
+   Walker-Delta Satellite Constellation <scenarioSatelliteConstellation>
 
 
 ``FormationBskSim`` Architecture

examples/scenarioSatelliteConstellation.py

@@ -0,0 +1,296 @@
+#
+#  ISC License
+#
+#  Copyright (c) 2023, Autonomous Vehicle Systems Lab, University of Colorado at Boulder
+#
+#  Permission to use, copy, modify, and/or distribute this software for any
+#  purpose with or without fee is hereby granted, provided that the above
+#  copyright notice and this permission notice appear in all copies.
+#
+#  THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES
+#  WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF
+#  MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR
+#  ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
+#  WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN AN
+#  ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF
+#  OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.
+#
+
+r"""
+
+Overview
+--------
+
+This scenario demonstrates how to set up a Walker-Delta constellation of satellites. Note that this scenario 
+uses the stand-alone Basilisk architecture rather than using the ''examples/FormationBskSim`` or ``examples/MultiSatBskSim``
+architectures for simultaneously simulating multiple spacecraft.
+
+The script is found in the folder ``basilisk/examples`` and executed by using::
+
+    python3 scenarioSatelliteConstellation.py
+
+When designing a satellite constellation, symmetry provides repeatable performance and makes design and operations simpler. One 
+such symmetric constellation design methodology is the Walker constellation which consists of circular orbits in evenly spaced 
+planes at a chosen inclination. Walker constellation layouts are fully specified by four parameters written as "i:T/P/F" where:
+    * i = inclination angle
+    * T = total number of satellites in the constellation
+    * P = numer of orbit planes evenly spaced in node angle
+    * F = relative spacing between satellites in adjacent orbit planes
+
+In order to simulate the constellation a semi-major axis 'a' is also specified as an input. The total number of satellites T must 
+be specified as an integer multiple of the number of planes P. The relative spacing F must be an integer between 0 and (P-1).
+
+Illustration of Simulation Results
+----------------------------------
+
+::
+
+    show_plots = True, a = 29994000, i = 56, T = 24, P = 3, F = 1
+
+The default scenario is modelled after the Galileo constellation with Walker numbers 56:24/3/1 and a semi-major axis of 29,994 km.
+
+.. image:: /_images/Scenarios/scenarioSatelliteConstellation.svg
+   :align: center
+
+Each line color corresponds to a different orbital plane, in this case 3 planes. Green circles mark the start of each of the 24 
+individual satellite's ground tracks and a red circle marks their end.
+
+Simulation Visualization In Vizard
+----------------------------------
+
+To add a visualization of each spacecraft's Earth nadir pointing, a cone aligned with the pointing direction
+is created through the ``vizInterface``::
+
+    vizSupport.createConeInOut(viz, fromBodyName=scList[i].ModelTag, toBodyName='earth', coneColor='teal',
+                               normalVector_B=[1, 0, 0], incidenceAngle=30/macros.D2R, isKeepIn=False,
+                               coneHeight=10.0, coneName='pointingCone')
+
+.. image:: /_images/static/scenarioSatelliteConstellationVizard.jpg
+   :scale: 50%
+   :align: center
+
+"""
+
+#
+# Basilisk Scenario Script and Integrated Test
+#
+# Purpose:  Demonstrates how to use the stand alone Basilisk architecture to automate 
+#           the setup of a Walker constellation.
+# Author:   Andrew Morell
+# Creation Date:  Dec. 13, 2023
+#
+
+import os
+
+import matplotlib.pyplot as plt
+import numpy as np
+from Basilisk import __path__
+
+bskPath = __path__[0]
+fileName = os.path.basename(os.path.splitext(__file__)[0])
+
+from Basilisk.simulation import (spacecraft, ephemerisConverter)
+from Basilisk.utilities import (SimulationBaseClass, macros, orbitalMotion,
+                                simIncludeGravBody, unitTestSupport, vizSupport)
+from Basilisk.fswAlgorithms import locationPointing
+from Basilisk.simulation import (simpleNav, planetEphemeris)
+
+# always import the Basilisk messaging support
+
+def run(show_plots, a, i, T, P, F):
+    """
+    The scenario can be run with the followings setups parameters:
+
+    Args:
+        show_plots (bool): Determines if the script should display plots
+        a (float): semi-major axis of all satellites [m]
+        i (float): orbital inclination [rad]
+        T (int): total number of satellites
+        P (int): number of equally spaced orbital planes
+        F (int): phasing between satellites in adjacent places
+
+    """
+
+    # Create simulation variable names
+    simTaskName = "simTask"
+    simProcessName = "simProcess"
+
+    #  Create a sim module as an empty container
+    scSim = SimulationBaseClass.SimBaseClass()
+
+    # (Optional) If you want to see a simulation progress bar in the terminal window, the
+    # use the following SetProgressBar(True) statement
+    scSim.SetProgressBar(True)
+
+    #  create the simulation process
+    dynProcess = scSim.CreateNewProcess(simProcessName)
+
+    # create the dynamics task and specify the integration update time
+    simulationTimeStep = macros.sec2nano(30.)
+    dynProcess.addTask(scSim.CreateNewTask(simTaskName, simulationTimeStep))
+
+    # clear prior gravitational body and SPICE setup definitions
+    gravFactory = simIncludeGravBody.gravBodyFactory()
+
+    # setup Earth Gravity Body
+    earth = gravFactory.createEarth()
+    earth.isCentralBody = True  # ensure this is the central gravitational body
+    mu = earth.mu
+
+    # Create the ephemeris converter module
+    spiceObject = gravFactory.createSpiceInterface(time="2029 June 12 5:30:30.0", epochInMsg=True)
+    spiceObject.zeroBase = 'Earth'
+    scSim.AddModelToTask(simTaskName, spiceObject)
+    ephemObject = ephemerisConverter.EphemerisConverter()
+    ephemObject.ModelTag = "ephemData"
+    ephemObject.addSpiceInputMsg(spiceObject.planetStateOutMsgs[0]) # Earth
+    scSim.AddModelToTask(simTaskName, ephemObject)
+
+    # spacecraft inertia
+    I = [900., 0., 0.,
+        0., 800., 0.,
+        0., 0., 600.]
+
+    # set the simulation time to be two orbit periods
+    n = np.sqrt(mu / a / a / a)
+    Pr = 2. * np.pi / n
+    simulationTime = macros.sec2nano(1.5 * Pr)
+
+    # create a logging task object of the spacecraft output message at the desired down sampling ratio
+    numDataPoints = 500
+    samplingTime = unitTestSupport.samplingTime(simulationTime, simulationTimeStep, numDataPoints)
+
+    if (np.mod(T,1) !=0 or np.mod(P,1) !=0 or np.mod(F,1) != 0 or T < 1 or P < 1 or F < 1):
+        raise Exception('number of satellites T, number of planes P, and relative spaceing F must be positive integer numbers')
+    if (np.mod(T,P) !=0):
+        raise Exception('number of satellites T must be an integer multiple of P')
+
+    # pre-compute Walker constellation parameters
+    PU = 360/T # patter unit in degrees
+    S = T/P # number of satellites in each plane
+    delta_RAAN = S*PU # delta RAAN in degrees
+    delta_nu = P*PU # in plane spacing of satellites
+    oe = orbitalMotion.ClassicElements()
+    oe.a = a
+    oe.e = 0.01
+    oe.i = i
+
+    # initialize T spacecraft objects and set properties
+    scList = []
+    navList = []
+    guideList = []
+    dataRec = []
+    for i in range(T):
+        # initialize spacecraft object
+        scList.append(spacecraft.Spacecraft())
+        scList[i].ModelTag = "bsk-Sat-"+str(i)
+        scSim.AddModelToTask(simTaskName, scList[i])
+
+        # attach gravity model to spacecraft
+        scList[i].gravField.gravBodies = spacecraft.GravBodyVector(list(gravFactory.gravBodies.values()))
+
+        # setup the orbit using classical orbit elements
+        RAAN = delta_RAAN * np.floor(i/S)
+        oe.Omega = RAAN * macros.D2R
+        oe.omega = F * np.floor(i/S) * macros.D2R
+        oe.f = delta_nu * np.mod(i,S) * macros.D2R
+        rN, vN = orbitalMotion.elem2rv(mu, oe)
+        scList[i].hub.r_CN_NInit = rN  # [m]
+        scList[i].hub.v_CN_NInit = vN  # [m/s]
+
+        # set the spacecraft mass and inertia
+        scList[i].hub.mHub = 750.0  # kg - spacecraft mass
+        scList[i].hub.IHubPntBc_B = unitTestSupport.np2EigenMatrix3d(I)
+
+        # set up prescribed attitude guidance
+        navList.append(simpleNav.SimpleNav())
+        navList[i].ModelTag = "SimpleNav"+str(i)
+        scSim.AddModelToTask(simTaskName, navList[i])
+        navList[i].scStateInMsg.subscribeTo(scList[i].scStateOutMsg)
+
+        # setup nadir pointing guidance module
+        guideList.append(locationPointing.locationPointing())
+        guideList[i].ModelTag = "nadirPoint"+str(i)
+        guideList[i].pHat_B = [1, 0, 0]
+        guideList[i].scTransInMsg.subscribeTo(navList[i].transOutMsg)
+        guideList[i].scAttInMsg.subscribeTo(navList[i].attOutMsg)
+        guideList[i].celBodyInMsg.subscribeTo(ephemObject.ephemOutMsgs[0])
+        scSim.AddModelToTask(simTaskName, guideList[i])
+        scList[i].attRefInMsg.subscribeTo(guideList[i].attRefOutMsg)
+
+        # record spacecraft states
+        dataRec.append(scList[i].scStateOutMsg.recorder(samplingTime))
+        scSim.AddModelToTask(simTaskName, dataRec[i])
+
+    # If you wish to transmit the simulation data to the United based Vizard Visualization application,
+    # then uncomment the following saveFile line
+    if vizSupport.vizFound:
+        viz = vizSupport.enableUnityVisualization(scSim, simTaskName, scList,
+                                                  # saveFile=__file__
+                                                  )
+        for i in range(T):
+            vizSupport.createConeInOut(viz, fromBodyName=scList[i].ModelTag, toBodyName='earth', coneColor='teal',
+                           normalVector_B=[1, 0, 0], incidenceAngle=30/macros.D2R, isKeepIn=False,
+                           coneHeight=10.0, coneName='pointingCone')
+
+    # initialize simulation
+    scSim.InitializeSimulation()
+
+    # configure a simulation stop time and execute the simulation run
+    scSim.ConfigureStopTime(simulationTime)
+    scSim.ExecuteSimulation()
+
+    # plot the ground track
+    plt.close("all")  # clears out plots from earlier test runs
+    figureList = {}
+    plt.figure(1)
+    for n in range(T):
+        # calculate latitude and longitude
+        [lon, lat] = rv2latlon(dataRec[n])
+        color = unitTestSupport.getLineColor(int(np.floor(n/S)), 3)
+        plt.scatter(lat, lon, s=5, c=[color])
+        plt.scatter(lat[0], lon[0], s=20, c='#00A300')
+        plt.scatter(lat[-1], lon[-1], s=20, c='#A30000')
+    plt.xlim([-180,180])
+    plt.ylim([-90,90])
+    plt.xticks(range(-180,181,30))
+    plt.yticks(range(-90,91,15))
+    plt.xlabel('longitude (degrees)')
+    plt.ylabel('lattitude (degrees)')
+
+    pltName = fileName
+    figureList[pltName] = plt.figure(1)
+
+    if show_plots:
+        plt.show()
+
+    # close the plots being saved off to avoid over-writing old and new figures
+    plt.close("all")
+
+    return figureList
+
+def rv2latlon(dataRec):
+    times = dataRec.times() * macros.NANO2SEC
+    rECI = dataRec.r_BN_N
+    lat = np.zeros(times.shape)
+    lon = np.zeros(times.shape)
+    for i in range(len(times)):
+        theta = times[i]*planetEphemeris.OMEGA_EARTH
+        dcm_ECI2ECEF = [[np.cos(theta), np.sin(theta), 0],
+                        [-np.sin(theta), np.cos(theta), 0],
+                        [0, 0, 1]]
+        rECEF = dcm_ECI2ECEF@rECI[i,:]
+        lat[i] = np.arcsin(rECEF[2]/np.linalg.norm(rECEF)) * macros.R2D
+        lon[i] = np.arctan2(rECEF[1],rECEF[0]) * macros.R2D
+
+    return lat, lon
+
+if __name__ == "__main__":
+    run(
+        True,       # show_plots
+        29994000,   # semi-major axis [m]
+        56,         # orbit inclination [deg]
+        24,         # total number of satellites (int)
+        3,          # number of orbit planes (int)
+        1           # phasing (int)
+    )

src/tests/test_scenarioSatelliteConstellation.py

@@ -0,0 +1,76 @@
+#
+#  ISC License
+#
+#  Copyright (c) 2023, Autonomous Vehicle Systems Lab, University of Colorado at Boulder
+#
+#  Permission to use, copy, modify, and/or distribute this software for any
+#  purpose with or without fee is hereby granted, provided that the above
+#  copyright notice and this permission notice appear in all copies.
+#
+#  THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES
+#  WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF
+#  MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR
+#  ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
+#  WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN AN
+#  ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF
+#  OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.
+#
+
+
+#
+# Basilisk Scenario Script and Integrated Test
+#
+# Author:   Hanspeter Schaub
+# Creation Date:  Dec. 13, 2023
+#
+
+import inspect
+import os
+import sys
+
+import pytest
+from Basilisk.utilities import unitTestSupport
+
+# Get current file path
+filename = inspect.getframeinfo(inspect.currentframe()).filename
+path = os.path.dirname(os.path.abspath(filename))
+
+sys.path.append(path + '/../../examples')
+import scenarioSatelliteConstellation
+
+@pytest.mark.scenarioTest
+
+def test_bskAttitudePointing(show_plots):
+    """This function is called by the py.test environment."""
+    # each test method requires a single assert method to be called
+
+    testFailCount = 0  # zero unit test result counter
+    testMessages = []  # create empty array to store test log messages
+
+    try:
+        figureList = scenarioSatelliteConstellation.run(show_plots,
+                                                        29994000,   # semi-major axis [m]
+                                                        56,         # orbit inclination [deg]
+                                                        24,         # total number of satellites (int)
+                                                        3,          # number of orbit planes (int)
+                                                        1)          # phasing (int)
+
+        # save the figures to the Doxygen scenario images folder
+        for pltName, plt in list(figureList.items()):
+            unitTestSupport.saveScenarioFigure(pltName, plt, path)
+
+    except OSError as err:
+        testFailCount += 1
+        testMessages.append("scenarioSatelliteConstellation  test are failed.")
+
+    #   print out success message if no error were found
+    if testFailCount == 0:
+        print("PASSED ")
+    else:
+        print(testFailCount)
+        print(testMessages)
+
+    # each test method requires a single assert method to be called
+    # this check below just makes sure no sub-test failures were found
+    assert testFailCount < 1, testMessages
+