diff --git a/getting_started copy.ipynb b/getting_started copy.ipynb
new file mode 100644
index 000000000..8a67849bf
--- /dev/null
+++ b/getting_started copy.ipynb
@@ -0,0 +1,1149 @@
+{
+ "cells": [
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "# Getting Started"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Here we go through a simplified rocket trajectory simulation to get you started. Let's start by importing the rocketpy module."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 1,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "%load_ext autoreload\n",
+ "%autoreload 2"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 15,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "from rocketpy import Environment, SolidMotor, Rocket, Flight, Dispersion\n"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "If you are using Jupyter Notebooks, it is recommended to run the following line to make matplotlib plots which will be shown later interactive and higher quality."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 18,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "%matplotlib widget\n",
+ "%config InlineBackend.figure_formats = ['svg']\n",
+ "%matplotlib inline"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "## Setting Up a Simulation"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "### Creating an Environment for Spaceport America"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 4,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "Env = Environment(\n",
+ " railLength=5.2, latitude = -23.363611, longitude = -48.011389, elevation=1400)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "To get weather data from the GFS forecast, available online, we run the following lines.\n",
+ "\n",
+ "First, we set tomorrow's date."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 5,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "import datetime\n",
+ "\n",
+ "tomorrow = datetime.date.today() + datetime.timedelta(days=1)\n",
+ "\n",
+ "Env.setDate((2019, 8, 10, 12)) # Hour given in UTC time"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Then, we tell Env to use a GFS forecast to get the atmospheric conditions for flight.\n",
+ "\n",
+ "Don't mind the warning, it just means that not all variables, such as wind speed or atmospheric temperature, are available at all altitudes given by the forecast."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 6,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "Env.setAtmosphericModel(type='Ensemble', file='data\\weather\\LASC2019_TATUI_reanalysis_ensemble.nc', dictionary=\"ECMWF\")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "We can see what the weather will look like by calling the info method!"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "Env.info()"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "### Creating a Motor\n",
+ "\n",
+ "A solid rocket motor is used in this case. To create a motor, the SolidMotor class is used and the required arguments are given.\n",
+ "\n",
+ "The SolidMotor class requires the user to have a thrust curve ready. This can come either from a .eng file for a commercial motor, such as below, or a .csv file from a static test measurement.\n",
+ "\n",
+ "Besides the thrust curve, other parameters such as grain properties and nozzle dimensions must also be given."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 7,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "Pro75M1670 = SolidMotor(\n",
+ " thrustSource=\"data\\motors\\Cesaroni_M1670.eng\",\n",
+ " burnOut=3.9,\n",
+ " grainNumber=5,\n",
+ " grainSeparation=5 / 1000,\n",
+ " grainDensity=1815,\n",
+ " grainOuterRadius=33 / 1000,\n",
+ " grainInitialInnerRadius=15 / 1000,\n",
+ " grainInitialHeight=120 / 1000,\n",
+ " nozzleRadius=33 / 1000,\n",
+ " throatRadius=11 / 1000,\n",
+ " interpolationMethod=\"linear\",\n",
+ ")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "To see what our thrust curve looks like, along with other import properties, we invoke the info method yet again. You may try the allInfo method if you want more information all at once!"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "Pro75M1670.info()"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "### Creating a Rocket"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "A rocket is composed of several components. Namely, we must have a motor (good thing we have the Pro75M1670 ready), a couple of aerodynamic surfaces (nose cone, fins and tail) and parachutes (if we are not launching a missile).\n",
+ "\n",
+ "Let's start by initializing our rocket, named Calisto, supplying it with the Pro75M1670 engine, entering its inertia properties, some dimensions and also its drag curves."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 8,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "Calisto = Rocket(\n",
+ " motor=Pro75M1670,\n",
+ " radius=127 / 2000,\n",
+ " mass=19.197 - 2.956,\n",
+ " inertiaI=6.60,\n",
+ " inertiaZ=0.0351,\n",
+ " distanceRocketNozzle=-1.255,\n",
+ " distanceRocketPropellant=-0.85704,\n",
+ " powerOffDrag=\"data/calisto/powerOffDragCurve.csv\",\n",
+ " powerOnDrag=\"data/calisto/powerOnDragCurve.csv\",\n",
+ ")\n",
+ "\n",
+ "Calisto.setRailButtons([0.2, -0.5])"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "#### Adding Aerodynamic Surfaces"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Now we define the aerodynamic surfaces. They are really straight forward."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 9,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "NoseCone = Calisto.addNose(length=0.55829, kind=\"vonKarman\", distanceToCM=0.71971)\n",
+ "\n",
+ "FinSet = Calisto.addFins(\n",
+ " 4, span=0.100, rootChord=0.120, tipChord=0.040, distanceToCM=-1.04956\n",
+ ")\n",
+ "\n",
+ "Tail = Calisto.addTail(\n",
+ " topRadius=0.0635, bottomRadius=0.0435, length=0.060, distanceToCM=-1.194656\n",
+ ")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "#### Adding Parachutes"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Finally, we have parachutes! Calisto will have two parachutes, Drogue and Main.\n",
+ "\n",
+ "Both parachutes are activated by some special algorithm, which is usually really complex and a trade secret. Most algorithms are based on pressure sampling only, while some also use acceleration info.\n",
+ "\n",
+ "RocketPy allows you to define a trigger function which will decide when to activate the ejection event for each parachute. This trigger function is supplied with pressure measurement at a predefined sampling rate. This pressure signal is usually noisy, so artificial noise parameters can be given. Call help(Rocket.addParachute) for more details. Furthermore, the trigger function also receives the complete state vector of the rocket, allowing us to use velocity, acceleration or even attitude to decide when the parachute event should be triggered.\n",
+ "\n",
+ "Here, we define our trigger functions rather simply using Python. However, you can call the exact code which will fly inside your rocket as well."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 10,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "def drogueTrigger(p, y):\n",
+ " # p = pressure\n",
+ " # y = [x, y, z, vx, vy, vz, e0, e1, e2, e3, w1, w2, w3]\n",
+ " # activate drogue when vz < 0 m/s.\n",
+ " return True if y[5] < 0 else False\n",
+ "\n",
+ "\n",
+ "def mainTrigger(p, y):\n",
+ " # p = pressure\n",
+ " # y = [x, y, z, vx, vy, vz, e0, e1, e2, e3, w1, w2, w3]\n",
+ " # activate main when vz < 0 m/s and z < 800 + 1400 m (+1400 due to surface elevation).\n",
+ " return True if y[5] < 0 and y[2] < 800 + 1400 else False\n",
+ "\n",
+ "\n",
+ "Main = Calisto.addParachute(\n",
+ " \"Main\",\n",
+ " CdS=10.0,\n",
+ " trigger=mainTrigger,\n",
+ " samplingRate=105,\n",
+ " lag=1.5,\n",
+ " noise=(0, 8.3, 0.5),\n",
+ ")\n",
+ "\n",
+ "Drogue = Calisto.addParachute(\n",
+ " \"Drogue\",\n",
+ " CdS=1.0,\n",
+ " trigger=drogueTrigger,\n",
+ " samplingRate=105,\n",
+ " lag=1.5,\n",
+ " noise=(0, 8.3, 0.5),\n",
+ ")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Just be careful if you run this last cell multiple times! If you do so, your rocket will end up with lots of parachutes which activate together, which may cause problems during the flight simulation. We advise you to re-run all cells which define our rocket before running this, preventing unwanted old parachutes. Alternatively, you can run the following lines to remove parachutes.\n",
+ "\n",
+ "```python\n",
+ "Calisto.parachutes.remove(Drogue)\n",
+ "Calisto.parachutes.remove(Main)\n",
+ "```"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "## Simulating a Flight\n",
+ "\n",
+ "Simulating a flight trajectory is as simple as initializing a Flight class object givin the rocket and environnement set up above as inputs. The launch rail inclination and heading are also given here."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 11,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "TestFlight = Flight(rocket=Calisto, environment=Env, inclination=85, heading=0)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "## Analyzing the Results\n",
+ "\n",
+ "RocketPy gives you many plots, thats for sure! They are divided into sections to keep them organized. Alternatively, see the Flight class documentation to see how to get plots for specific variables only, instead of all of them at once."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {
+ "scrolled": false
+ },
+ "outputs": [],
+ "source": [
+ "TestFlight.allInfo()"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Export Flight Trajectory to a .kml file so it can be opened on Google Earth"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "TestFlight.exportKML(\n",
+ " fileName=\"trajectory.kml\",\n",
+ " extrude=True,\n",
+ " altitudeMode=\"relativetoground\",\n",
+ ")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "## Dispersion analysis\n",
+ "\n",
+ "Monte Carlo analysis to predict probability distributions of the rocket's \n",
+ "landing point, apogee and other relevant information."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 20,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "'Completed 10 iterations successfully. Total CPU time: 17.6875 s. Total wall time: 19.050102710723877 s'"
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ }
+ ],
+ "source": [
+ "D = Dispersion(number_of_simulations=10, flight=TestFlight, filename='teste', dispersionDict={'rocketMass':(20, 0.2), 'burnOut': (3.9,0.3), 'railLength': (5.2, 0.5)},environment=None, motor=None, rocket=None, distributionType='normal')"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 22,
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "image/svg+xml": "\n\n\n",
+ "text/plain": [
+ ""
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ },
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Apogee Altitude - Mean Value: 2555.689 m\n",
+ "Apogee Altitude - Standard Deviation: 18.910 m\n"
+ ]
+ },
+ {
+ "data": {
+ "image/svg+xml": "\n\n\n",
+ "text/plain": [
+ ""
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ },
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Out of Rail Velocity - Mean Value: 23.353 m/s\n",
+ "Out of Rail Velocity - Standard Deviation: 0.062 m/s\n"
+ ]
+ },
+ {
+ "data": {
+ "image/svg+xml": "\n\n\n",
+ "text/plain": [
+ ""
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ },
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Initial Static Margin - Mean Value: 2.220 c\n",
+ "Initial Static Margin - Standard Deviation: 0.004 c\n",
+ "Out of Rail Static Margin - Mean Value: 2.299 c\n",
+ "Out of Rail Static Margin - Standard Deviation: 0.004 c\n",
+ "Final Static Margin - Mean Value: 3.090 c\n",
+ "Final Static Margin - Standard Deviation: 0.000 c\n"
+ ]
+ },
+ {
+ "data": {
+ "image/svg+xml": "\n\n\n",
+ "text/plain": [
+ ""
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ },
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Lateral Surface Wind Speed - Mean Value: -3.754 m/s\n",
+ "Lateral Surface Wind Speed - Standard Deviation: 0.000 m/s\n"
+ ]
+ },
+ {
+ "data": {
+ "image/svg+xml": "\n\n\n",
+ "text/plain": [
+ ""
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ },
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Frontal Surface Wind Speed - Mean Value: -6.318 m/s\n",
+ "Frontal Surface Wind Speed - Standard Deviation: 0.000 m/s\n"
+ ]
+ },
+ {
+ "data": {
+ "image/svg+xml": "\n\n\n",
+ "text/plain": [
+ ""
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ },
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Out of Rail Time - Mean Value: 0.394 s\n",
+ "Out of Rail Time - Standard Deviation: 0.001 s\n"
+ ]
+ },
+ {
+ "data": {
+ "image/svg+xml": "\n\n\n",
+ "text/plain": [
+ ""
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ },
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Impact Time - Mean Value: 149.673 s\n",
+ "Impact Time - Standard Deviation: 1.276 s\n"
+ ]
+ },
+ {
+ "data": {
+ "image/svg+xml": "\n\n\n",
+ "text/plain": [
+ ""
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ },
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Apogee X Position - Mean Value: -241.404 m\n",
+ "Apogee X Position - Standard Deviation: 0.608 m\n"
+ ]
+ },
+ {
+ "data": {
+ "image/svg+xml": "\n\n\n",
+ "text/plain": [
+ ""
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ },
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Apogee Y Position - Mean Value: 875.081 m\n",
+ "Apogee Y Position - Standard Deviation: 4.025 m\n"
+ ]
+ },
+ {
+ "data": {
+ "image/svg+xml": "\n\n\n",
+ "text/plain": [
+ ""
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ },
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Impact Time - Mean Value: 149.673 s\n",
+ "Impact Time - Standard Deviation: 1.276 s\n"
+ ]
+ },
+ {
+ "data": {
+ "image/svg+xml": "\n\n\n",
+ "text/plain": [
+ ""
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ },
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Impact Velocity - Mean Value: -19.497 m/s\n",
+ "Impact Velocity - Standard Deviation: 0.056 m/s\n"
+ ]
+ },
+ {
+ "data": {
+ "image/svg+xml": "\n\n\n",
+ "text/plain": [
+ ""
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ },
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Impact X Position - Mean Value: 541.903 m\n",
+ "Impact X Position - Standard Deviation: 8.830 m\n"
+ ]
+ },
+ {
+ "data": {
+ "image/svg+xml": "\n\n\n",
+ "text/plain": [
+ ""
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ },
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Impact Y Position - Mean Value: 571.933 m\n",
+ "Impact Y Position - Standard Deviation: 0.985 m\n"
+ ]
+ },
+ {
+ "data": {
+ "image/svg+xml": "\n\n\n",
+ "text/plain": [
+ ""
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ },
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Maximum Velocity - Mean Value: 235.166 m/s\n",
+ "Maximum Velocity - Standard Deviation: 1.315 m/s\n"
+ ]
+ },
+ {
+ "data": {
+ "image/svg+xml": "\n\n\n",
+ "text/plain": [
+ ""
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ },
+ {
+ "data": {
+ "image/svg+xml": "\n\n\n",
+ "text/plain": [
+ ""
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ },
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Drogue Fully Inflated Time - Mean Value: 24.945 s\n",
+ "Drogue Fully Inflated Time - Standard Deviation: 0.067 s\n"
+ ]
+ },
+ {
+ "data": {
+ "image/svg+xml": "\n\n\n",
+ "text/plain": [
+ ""
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ },
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Drogue Parachute Fully Inflated Velocity - Mean Value: 39.770 m/s\n",
+ "Drogue Parachute Fully Inflated Velocity - Standard Deviation: 0.009 m/s\n"
+ ]
+ },
+ {
+ "data": {
+ "image/svg+xml": "\n\n\n",
+ "text/plain": [
+ ""
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ },
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Drogue Trigger Time - Mean Value: 23.445 s\n",
+ "Drogue Trigger Time - Standard Deviation: 0.067 s\n"
+ ]
+ },
+ {
+ "data": {
+ "image/svg+xml": "\n\n\n",
+ "text/plain": [
+ ""
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ }
+ ],
+ "source": [
+ "D = Dispersion(filename='teste.disp_outputs.txt', image=\"data/valetudo/valetudo_basemap.png\")\n",
+ "D.info()"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "## Using Simulation for Design\n",
+ "\n",
+ "Here, we go through a couple of examples which make use of RocketPy in cool ways to help us design our rocket."
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "### Dynamic Stability Analysis"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Ever wondered how static stability translates into dynamic stability? Different static margins result in different dynamic behavior, which also depends on the rocket's rotational inertial.\n",
+ "\n",
+ "Let's make use of RocketPy's helper class called Function to explore how the dynamic stability of Calisto varies if we change the fins span by a certain factor."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# Helper class\n",
+ "from rocketpy import Function\n",
+ "\n",
+ "# Prepare Rocket Class\n",
+ "Calisto = Rocket(\n",
+ " motor=Pro75M1670,\n",
+ " radius=127 / 2000,\n",
+ " mass=19.197 - 2.956,\n",
+ " inertiaI=6.60,\n",
+ " inertiaZ=0.0351,\n",
+ " distanceRocketNozzle=-1.255,\n",
+ " distanceRocketPropellant=-0.85704,\n",
+ " powerOffDrag=\"../../data/calisto/powerOffDragCurve.csv\",\n",
+ " powerOnDrag=\"../../data/calisto/powerOnDragCurve.csv\",\n",
+ ")\n",
+ "Calisto.setRailButtons([0.2, -0.5])\n",
+ "Nose = Calisto.addNose(length=0.55829, kind=\"vonKarman\", distanceToCM=0.71971)\n",
+ "FinSet = Calisto.addFins(\n",
+ " 4, span=0.1, rootChord=0.120, tipChord=0.040, distanceToCM=-1.04956\n",
+ ")\n",
+ "Tail = Calisto.addTail(\n",
+ " topRadius=0.0635, bottomRadius=0.0435, length=0.060, distanceToCM=-1.194656\n",
+ ")\n",
+ "\n",
+ "# Prepare Environment Class\n",
+ "Env = Environment(5.2, 9.8)\n",
+ "Env.setAtmosphericModel(type=\"CustomAtmosphere\", wind_v=-5)\n",
+ "\n",
+ "# Simulate Different Static Margins by Varying Fin Position\n",
+ "simulation_results = []\n",
+ "\n",
+ "for factor in [0.5, 0.7, 0.9, 1.1, 1.3]:\n",
+ " # Modify rocket fin set by removing previous one and adding new one\n",
+ " Calisto.aerodynamicSurfaces.remove(FinSet)\n",
+ " FinSet = Calisto.addFins(\n",
+ " 4, span=0.1, rootChord=0.120, tipChord=0.040, distanceToCM=-1.04956 * factor\n",
+ " )\n",
+ " # Simulate\n",
+ " print(\n",
+ " \"Simulating Rocket with Static Margin of {:1.3f}->{:1.3f} c\".format(\n",
+ " Calisto.staticMargin(0), Calisto.staticMargin(Calisto.motor.burnOutTime)\n",
+ " )\n",
+ " )\n",
+ " TestFlight = Flight(\n",
+ " rocket=Calisto,\n",
+ " environment=Env,\n",
+ " inclination=90,\n",
+ " heading=0,\n",
+ " maxTimeStep=0.01,\n",
+ " maxTime=5,\n",
+ " terminateOnApogee=True,\n",
+ " verbose=True,\n",
+ " )\n",
+ " # Post process flight data\n",
+ " TestFlight.postProcess()\n",
+ " # Store Results\n",
+ " staticMarginAtIgnition = Calisto.staticMargin(0)\n",
+ " staticMarginAtOutOfRail = Calisto.staticMargin(TestFlight.outOfRailTime)\n",
+ " staticMarginAtSteadyState = Calisto.staticMargin(TestFlight.tFinal)\n",
+ " simulation_results += [\n",
+ " (\n",
+ " TestFlight.attitudeAngle,\n",
+ " \"{:1.2f} c | {:1.2f} c | {:1.2f} c\".format(\n",
+ " staticMarginAtIgnition,\n",
+ " staticMarginAtOutOfRail,\n",
+ " staticMarginAtSteadyState,\n",
+ " ),\n",
+ " )\n",
+ " ]\n",
+ "\n",
+ "Function.comparePlots(\n",
+ " simulation_results,\n",
+ " lower=0,\n",
+ " upper=1.5,\n",
+ " xlabel=\"Time (s)\",\n",
+ " ylabel=\"Attitude Angle (deg)\",\n",
+ ")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "### Characteristic Frequency Calculation\n",
+ "\n",
+ "Here we analyze the characteristic frequency of oscillation of our rocket just as it leaves the launch rail. Note that when we ran TestFlight.allInfo(), one of the plots already showed us the frequency spectrum of our flight. Here, however, we have more control of what we are plotting."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "import numpy as np\n",
+ "import matplotlib.pyplot as plt\n",
+ "\n",
+ "Env = Environment(\n",
+ " railLength=5.2, latitude=32.990254, longitude=-106.974998, elevation=1400\n",
+ ")\n",
+ "\n",
+ "Env.setAtmosphericModel(type=\"CustomAtmosphere\", wind_v=-5)\n",
+ "\n",
+ "# Prepare Motor\n",
+ "Pro75M1670 = SolidMotor(\n",
+ " thrustSource=\"../../data/motors/Cesaroni_M1670.eng\",\n",
+ " burnOut=3.9,\n",
+ " grainNumber=5,\n",
+ " grainSeparation=5 / 1000,\n",
+ " grainDensity=1815,\n",
+ " grainOuterRadius=33 / 1000,\n",
+ " grainInitialInnerRadius=15 / 1000,\n",
+ " grainInitialHeight=120 / 1000,\n",
+ " nozzleRadius=33 / 1000,\n",
+ " throatRadius=11 / 1000,\n",
+ " interpolationMethod=\"linear\",\n",
+ ")\n",
+ "\n",
+ "# Prepare Rocket\n",
+ "Calisto = Rocket(\n",
+ " motor=Pro75M1670,\n",
+ " radius=127 / 2000,\n",
+ " mass=19.197 - 2.956,\n",
+ " inertiaI=6.60,\n",
+ " inertiaZ=0.0351,\n",
+ " distanceRocketNozzle=-1.255,\n",
+ " distanceRocketPropellant=-0.85704,\n",
+ " powerOffDrag=\"../../data/calisto/powerOffDragCurve.csv\",\n",
+ " powerOnDrag=\"../../data/calisto/powerOnDragCurve.csv\",\n",
+ ")\n",
+ "\n",
+ "Calisto.setRailButtons([0.2, -0.5])\n",
+ "\n",
+ "Nose = Calisto.addNose(length=0.55829, kind=\"vonKarman\", distanceToCM=0.71971)\n",
+ "FinSet = Calisto.addFins(\n",
+ " 4, span=0.1, rootChord=0.120, tipChord=0.040, distanceToCM=-1.04956\n",
+ ")\n",
+ "Tail = Calisto.addTail(\n",
+ " topRadius=0.0635, bottomRadius=0.0435, length=0.060, distanceToCM=-1.194656\n",
+ ")\n",
+ "\n",
+ "# Simulate first 5 seconds of Flight\n",
+ "TestFlight = Flight(\n",
+ " rocket=Calisto,\n",
+ " environment=Env,\n",
+ " inclination=90,\n",
+ " heading=0,\n",
+ " maxTimeStep=0.01,\n",
+ " maxTime=5,\n",
+ ")\n",
+ "TestFlight.postProcess()\n",
+ "\n",
+ "# Perform a Fourier Analysis\n",
+ "Fs = 100.0\n",
+ "# sampling rate\n",
+ "Ts = 1.0 / Fs\n",
+ "# sampling interval\n",
+ "t = np.arange(1, 400, Ts) # time vector\n",
+ "ff = 5\n",
+ "# frequency of the signal\n",
+ "y = TestFlight.attitudeAngle(t) - np.mean(TestFlight.attitudeAngle(t))\n",
+ "n = len(y) # length of the signal\n",
+ "k = np.arange(n)\n",
+ "T = n / Fs\n",
+ "frq = k / T # two sides frequency range\n",
+ "frq = frq[range(n // 2)] # one side frequency range\n",
+ "Y = np.fft.fft(y) / n # fft computing and normalization\n",
+ "Y = Y[range(n // 2)]\n",
+ "fig, ax = plt.subplots(2, 1)\n",
+ "ax[0].plot(t, y)\n",
+ "ax[0].set_xlabel(\"Time\")\n",
+ "ax[0].set_ylabel(\"Signal\")\n",
+ "ax[0].set_xlim((0, 5))\n",
+ "ax[1].plot(frq, abs(Y), \"r\") # plotting the spectrum\n",
+ "ax[1].set_xlabel(\"Freq (Hz)\")\n",
+ "ax[1].set_ylabel(\"|Y(freq)|\")\n",
+ "ax[1].set_xlim((0, 5))\n",
+ "plt.subplots_adjust(hspace=0.5)\n",
+ "plt.show()"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "### Apogee as a Function of Mass\n",
+ "\n",
+ "This one is a classic one! We always need to know how much our rocket's apogee will change when our payload gets heavier."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "def apogee(mass):\n",
+ " # Prepare Environment\n",
+ " Env = Environment(\n",
+ " railLength=5.2,\n",
+ " latitude=32.990254,\n",
+ " longitude=-106.974998,\n",
+ " elevation=1400,\n",
+ " date=(2018, 6, 20, 18),\n",
+ " )\n",
+ "\n",
+ " Env.setAtmosphericModel(type=\"CustomAtmosphere\", wind_v=-5)\n",
+ "\n",
+ " # Prepare Motor\n",
+ " Pro75M1670 = SolidMotor(\n",
+ " thrustSource=\"../../data/motors/Cesaroni_M1670.eng\",\n",
+ " burnOut=3.9,\n",
+ " grainNumber=5,\n",
+ " grainSeparation=5 / 1000,\n",
+ " grainDensity=1815,\n",
+ " grainOuterRadius=33 / 1000,\n",
+ " grainInitialInnerRadius=15 / 1000,\n",
+ " grainInitialHeight=120 / 1000,\n",
+ " nozzleRadius=33 / 1000,\n",
+ " throatRadius=11 / 1000,\n",
+ " interpolationMethod=\"linear\",\n",
+ " )\n",
+ "\n",
+ " # Prepare Rocket\n",
+ " Calisto = Rocket(\n",
+ " motor=Pro75M1670,\n",
+ " radius=127 / 2000,\n",
+ " mass=mass,\n",
+ " inertiaI=6.60,\n",
+ " inertiaZ=0.0351,\n",
+ " distanceRocketNozzle=-1.255,\n",
+ " distanceRocketPropellant=-0.85704,\n",
+ " powerOffDrag=\"../../data/calisto/powerOffDragCurve.csv\",\n",
+ " powerOnDrag=\"../../data/calisto/powerOnDragCurve.csv\",\n",
+ " )\n",
+ "\n",
+ " Calisto.setRailButtons([0.2, -0.5])\n",
+ " Nose = Calisto.addNose(length=0.55829, kind=\"vonKarman\", distanceToCM=0.71971)\n",
+ " FinSet = Calisto.addFins(\n",
+ " 4, span=0.1, rootChord=0.120, tipChord=0.040, distanceToCM=-1.04956\n",
+ " )\n",
+ " Tail = Calisto.addTail(\n",
+ " topRadius=0.0635, bottomRadius=0.0435, length=0.060, distanceToCM=-1.194656\n",
+ " )\n",
+ "\n",
+ " # Simulate Flight until Apogee\n",
+ " TestFlight = Flight(\n",
+ " rocket=Calisto,\n",
+ " environment=Env,\n",
+ " inclination=85,\n",
+ " heading=0,\n",
+ " terminateOnApogee=True,\n",
+ " )\n",
+ " return TestFlight.apogee\n",
+ "\n",
+ "\n",
+ "apogeebymass = Function(apogee, inputs=\"Mass (kg)\", outputs=\"Estimated Apogee (m)\")\n",
+ "apogeebymass.plot(8, 20, 20)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "### Out of Rail Speed as a Function of Mass\n",
+ "\n",
+ "To finish off, lets make a really important plot. Out of rail speed is the speed our rocket has when it is leaving the launch rail. This is crucial to make sure it can fly safely after leaving the rail. A common rule of thumb is that our rocket's out of rail speed should be 4 times the wind speed so that it does not stall and become unstable."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "def speed(mass):\n",
+ " # Prepare Environment\n",
+ " Env = Environment(\n",
+ " railLength=5.2,\n",
+ " latitude=32.990254,\n",
+ " longitude=-106.974998,\n",
+ " elevation=1400,\n",
+ " date=(2018, 6, 20, 18),\n",
+ " )\n",
+ "\n",
+ " Env.setAtmosphericModel(type=\"CustomAtmosphere\", wind_v=-5)\n",
+ "\n",
+ " # Prepare Motor\n",
+ " Pro75M1670 = SolidMotor(\n",
+ " thrustSource=\"../../data/motors/Cesaroni_M1670.eng\",\n",
+ " burnOut=3.9,\n",
+ " grainNumber=5,\n",
+ " grainSeparation=5 / 1000,\n",
+ " grainDensity=1815,\n",
+ " grainOuterRadius=33 / 1000,\n",
+ " grainInitialInnerRadius=15 / 1000,\n",
+ " grainInitialHeight=120 / 1000,\n",
+ " nozzleRadius=33 / 1000,\n",
+ " throatRadius=11 / 1000,\n",
+ " interpolationMethod=\"linear\",\n",
+ " )\n",
+ "\n",
+ " # Prepare Rocket\n",
+ " Calisto = Rocket(\n",
+ " motor=Pro75M1670,\n",
+ " radius=127 / 2000,\n",
+ " mass=mass,\n",
+ " inertiaI=6.60,\n",
+ " inertiaZ=0.0351,\n",
+ " distanceRocketNozzle=-1.255,\n",
+ " distanceRocketPropellant=-0.85704,\n",
+ " powerOffDrag=\"../../data/calisto/powerOffDragCurve.csv\",\n",
+ " powerOnDrag=\"../../data/calisto/powerOnDragCurve.csv\",\n",
+ " )\n",
+ "\n",
+ " Calisto.setRailButtons([0.2, -0.5])\n",
+ " Nose = Calisto.addNose(length=0.55829, kind=\"vonKarman\", distanceToCM=0.71971)\n",
+ " FinSet = Calisto.addFins(\n",
+ " 4, span=0.1, rootChord=0.120, tipChord=0.040, distanceToCM=-1.04956\n",
+ " )\n",
+ " Tail = Calisto.addTail(\n",
+ " topRadius=0.0635, bottomRadius=0.0435, length=0.060, distanceToCM=-1.194656\n",
+ " )\n",
+ "\n",
+ " # Simulate Flight until Apogee\n",
+ " TestFlight = Flight(\n",
+ " rocket=Calisto,\n",
+ " environment=Env,\n",
+ " inclination=85,\n",
+ " heading=0,\n",
+ " terminateOnApogee=True,\n",
+ " )\n",
+ " return TestFlight.outOfRailVelocity\n",
+ "\n",
+ "\n",
+ "speedbymass = Function(speed, inputs=\"Mass (kg)\", outputs=\"Out of Rail Speed (m/s)\")\n",
+ "speedbymass.plot(8, 20, 20)"
+ ]
+ }
+ ],
+ "metadata": {
+ "hide_input": false,
+ "kernelspec": {
+ "display_name": "Python 3.8.8 64-bit ('base': conda)",
+ "language": "python",
+ "name": "python3"
+ },
+ "language_info": {
+ "codemirror_mode": {
+ "name": "ipython",
+ "version": 3
+ },
+ "file_extension": ".py",
+ "mimetype": "text/x-python",
+ "name": "python",
+ "nbconvert_exporter": "python",
+ "pygments_lexer": "ipython3",
+ "version": "3.8.8"
+ },
+ "vscode": {
+ "interpreter": {
+ "hash": "ea183bb76f01b1c0f19d4faefaf72022d209702e86a95c61a348be375d9bcd4f"
+ }
+ }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 2
+}
diff --git a/rocketpy/Dispersion.py b/rocketpy/Dispersion.py
new file mode 100644
index 000000000..2ecdcab75
--- /dev/null
+++ b/rocketpy/Dispersion.py
@@ -0,0 +1,997 @@
+# -*- coding: utf-8 -*-
+
+__author__ = ""
+__copyright__ = "Copyright 20XX, Projeto Jupiter"
+__license__ = "MIT"
+
+
+from rocketpy import *
+import netCDF4
+
+
+
+from datetime import datetime
+from os import _Environ
+from time import process_time, perf_counter, time
+import glob
+import traceback
+
+
+import numpy as np
+from numpy.random import normal, uniform, choice, beta, binomial, chisquare, dirichlet, exponential, f, gamma, geometric, gumbel, hypergeometric, laplace, logistic, lognormal, logseries, multinomial, multivariate_normal, negative_binomial, noncentral_chisquare, noncentral_f, pareto, poisson, power, rayleigh, standard_cauchy, standard_exponential, standard_gamma, standard_normal, standard_t, triangular, vonmises, wald, weibull, zipf
+from IPython.display import display
+from rocketpy.Environment import Environment
+import matplotlib.pyplot as plt
+from imageio import imread
+from matplotlib.patches import Ellipse
+
+
+
+
+class Dispersion:
+
+ """Monte Carlo analysis to predict probability distributions of the rocket's
+ landing point, apogee and other relevant information.
+
+ Attributes
+ ----------
+
+ Parameters:
+ Dispersion.filename: string
+ When running a new simulation, this attribute represents the initial
+ part of the export filenames (e.g. 'filename.disp_outputs.txt').
+ When analyzing the results of a previous simulation, this attribute
+ shall be the filename containing the outputs of a dispersion calculation.
+ Dispersion.image: string
+ Launch site PNG file to be plotted along with the dispersion ellipses.
+ Attribute needed to run a new simulation.
+ Dispersion.realLandingPoint: tuple
+ Rocket's experimental landing point relative to launch point.
+ Dispersion.N: integer
+ Number of simulations in an output file.
+
+ Other classes:
+ Dispersion.environment: Environment
+ Launch environment.
+ Attribute needed to run a new simulation, when Dispersion.flight remains unchanged.
+ Dispersion.motor: Motor
+ Rocket's motor.
+ Attribute needed to run a new simulation, when Dispersion.flight remains unchanged.
+ Dispersion.rocket: Rocket
+ Rocket with nominal values.
+ Attribute needed to run a new simulation, when Dispersion.flight remains unchanged.
+
+ """
+
+ def __init__(
+ self,
+ filename,
+ number_of_simulations=0,
+ flight=None,
+ image=None,
+ dispersionDict={},
+ environment=None,
+ motor=None,
+ rocket=None,
+ distributionType='normal',
+ realLandingPoint=None
+ ):
+
+ '''
+ Parameters
+ ----------
+
+ filename: string
+ When running a new simulation, this parameter represents the initial
+ part of the export filenames (e.g. 'filename.disp_outputs.txt').
+ When analyzing the results of a previous simulation, this parameter
+ shall be the .txt filename containing the outputs of a dispersion calculation.
+
+ number_of_simulations: integer, needed when running a new simulation
+ Number of simulations desired, must be greater than zero.
+ Default is zero.
+
+ flight: Flight
+ Original rocket's flight with nominal values.
+ Parameter needed to run a new simulation, when environment,
+ motor and rocket remain unchanged.
+ Default is None.
+
+ image: string, needed when running a new simulation
+ Launch site PNG file to be plotted along with the dispersion ellipses.
+
+ dispersionDict: dictionary, optional
+ Contains the information of which environment, motor, rocket and flight variables
+ will vary according to its standard deviation.
+ Format {'parameter0': (nominal value, standard deviation), 'parameter1':
+ (nominal value, standard deviation), ...}
+ (e.g. {'rocketMass':(20, 0.2),
+ 'burnOut': (3.9, 0.3), 'railLength': (5.2, 0.05)})
+ Default is {}.
+
+ environment: Environment
+ Launch environment.
+ Parameter needed to run a new simulation, when Dispersion.flight remains unchanged.
+ Default is None.
+
+ motor: Motor, optional
+ Rocket's motor.
+ Parameter needed to run a new simulation, when Dispersion.flight remains unchanged.
+ Default is None.
+
+ rocket: Rocket, optional
+ Rocket with nominal values.
+ Parameter needed to run a new simulation, when Dispersion.flight remains unchanged.
+ Default is None.
+
+ distributionType: string, optional
+ Determines which type of distribution will be applied to variable parameters and
+ its respective standard deviation.
+ Default is 'normal'
+
+ realLandingPoint: tuple, optional
+ Rocket's experimental landing point relative to launch point.
+ Format (horizontal distance, vertical distance)
+
+ Returns
+ -------
+ None
+
+ '''
+
+ # Save parameters
+ self.filename = filename
+ self.image = image
+ self.realLandingPoint = realLandingPoint
+
+ # Run a new simulation
+ if number_of_simulations > 0:
+
+ self.flight = flight
+ self.environment = flight.env if not environment else environment
+ self.rocket = flight.rocket if not rocket else rocket
+ self.motor = flight.rocket.motor if not motor else motor
+
+ analysis_parameters = {i: j for i, j in dispersionDict.items()}
+
+ def flight_settings(analysis_parameters, number_of_simulations):
+ i = 0
+ while i < number_of_simulations:
+ # Generate a flight setting
+ flight_setting = {}
+ for parameter_key, parameter_value in analysis_parameters.items():
+ if type(parameter_value) is tuple:
+ if distributionType == "normal" or distributionType == None:
+ flight_setting[parameter_key] = normal(*parameter_value)
+ if distributionType == "beta":
+ flight_setting[parameter_key] = beta(*parameter_value)
+ if distributionType == "binomial":
+ flight_setting[parameter_key] = binomial(*parameter_value)
+ if distributionType == "chisquare":
+ flight_setting[parameter_key] = chisquare(*parameter_value)
+ if distributionType == "dirichlet":
+ flight_setting[parameter_key] = dirichlet(*parameter_value)
+ if distributionType == "exponential":
+ flight_setting[parameter_key] = exponential(*parameter_value)
+ if distributionType == "f":
+ flight_setting[parameter_key] = f(*parameter_value)
+ if distributionType == "gamma":
+ flight_setting[parameter_key] = gamma(*parameter_value)
+ if distributionType == "geometric":
+ flight_setting[parameter_key] = geometric(*parameter_value)
+ if distributionType == "gumbel":
+ flight_setting[parameter_key] = gumbel(*parameter_value)
+ if distributionType == "hypergeometric":
+ flight_setting[parameter_key] = hypergeometric(*parameter_value)
+ if distributionType == "laplace":
+ flight_setting[parameter_key] = laplace(*parameter_value)
+ if distributionType == "logistic":
+ flight_setting[parameter_key] = logistic(*parameter_value)
+ if distributionType == "lognormal":
+ flight_setting[parameter_key] = lognormal(*parameter_value)
+ if distributionType == "logseries":
+ flight_setting[parameter_key] = logseries(*parameter_value)
+ if distributionType == "multinomial":
+ flight_setting[parameter_key] = multinomial(*parameter_value)
+ if distributionType == "multivariate_normal":
+ flight_setting[parameter_key] = multivariate_normal(*parameter_value)
+ if distributionType == "negative_binomial":
+ flight_setting[parameter_key] = negative_binomial(*parameter_value)
+ if distributionType == "noncentral_chisquare":
+ flight_setting[parameter_key] = noncentral_chisquare(*parameter_value)
+ if distributionType == "noncentral_f":
+ flight_setting[parameter_key] = noncentral_f(*parameter_value)
+ if distributionType == "pareto":
+ flight_setting[parameter_key] = pareto(*parameter_value)
+ if distributionType == "poisson":
+ flight_setting[parameter_key] = poisson(*parameter_value)
+ if distributionType == "power":
+ flight_setting[parameter_key] = power(*parameter_value)
+ if distributionType == "rayleigh":
+ flight_setting[parameter_key] = rayleigh(*parameter_value)
+ if distributionType == "standard_cauchy":
+ flight_setting[parameter_key] = standard_cauchy(*parameter_value)
+ if distributionType == "standard_exponential":
+ flight_setting[parameter_key] = standard_exponential(*parameter_value)
+ if distributionType == "standard_gamma":
+ flight_setting[parameter_key] = standard_gamma(*parameter_value)
+ if distributionType == "standard_normal":
+ flight_setting[parameter_key] = standard_normal(*parameter_value)
+ if distributionType == "standard_t":
+ flight_setting[parameter_key] = standard_t(*parameter_value)
+ if distributionType == "triangular":
+ flight_setting[parameter_key] = triangular(*parameter_value)
+ if distributionType == "uniform":
+ flight_setting[parameter_key] = uniform(*parameter_value)
+ if distributionType == "vonmises":
+ flight_setting[parameter_key] = vonmises(*parameter_value)
+ if distributionType == "wald":
+ flight_setting[parameter_key] = wald(*parameter_value)
+ if distributionType == "weibull":
+ flight_setting[parameter_key] = weibull(*parameter_value)
+ if distributionType == "zipf":
+ flight_setting[parameter_key] = zipf(*parameter_value)
+ else:
+ flight_setting[parameter_key] = choice(parameter_value)
+
+ # Skip if certain values are negative, which happens due to the normal curve but isnt realistic
+ if "lag_rec" in analysis_parameters and flight_setting["lag_rec"] < 0:
+ continue
+ if "lag_se" in analysis_parameters and flight_setting["lag_se"] < 0:
+ continue
+ # Update counter
+ i += 1
+ # Yield a flight setting
+ yield flight_setting
+
+ def export_flight_data(flight_setting, flight_data, exec_time):
+ # Generate flight results
+ flight_result = {
+ "outOfRailTime": flight_data.outOfRailTime,
+ "outOfRailVelocity": flight_data.outOfRailVelocity,
+ "apogeeTime": flight_data.apogeeTime,
+ "apogeeAltitude": flight_data.apogee - flight_data.env.elevation,
+ "apogeeX": flight_data.apogeeX,
+ "apogeeY": flight_data.apogeeY,
+ "impactTime": flight_data.tFinal,
+ "impactX": flight_data.xImpact,
+ "impactY": flight_data.yImpact,
+ "impactVelocity": flight_data.impactVelocity,
+ "initialStaticMargin": flight_data.rocket.staticMargin(0),
+ "outOfRailStaticMargin": flight_data.rocket.staticMargin(flight_data.outOfRailTime),
+ "finalStaticMargin": flight_data.rocket.staticMargin(flight_data.rocket.motor.burnOutTime),
+ "numberOfEvents": len(flight_data.parachuteEvents),
+ "drogueTriggerTime": [],
+ "drogueInflatedTime": [],
+ "drogueInflatedVelocity": [],
+ "executionTime": exec_time,
+ 'lateralWind': flight_data.lateralSurfaceWind,
+ 'frontalWind': flight_data.frontalSurfaceWind}
+
+ # Calculate maximum reached velocity
+ sol = np.array(flight_data.solution)
+ flight_data.vx = Function(sol[:, [0, 4]], "Time (s)", "Vx (m/s)", "linear", extrapolation="natural")
+ flight_data.vy = Function(sol[:, [0, 5]], "Time (s)", "Vy (m/s)", "linear", extrapolation="natural")
+ flight_data.vz = Function(sol[:, [0, 6]], "Time (s)", "Vz (m/s)", "linear", extrapolation="natural")
+ flight_data.v = (flight_data.vx**2 + flight_data.vy**2 + flight_data.vz**2) ** 0.5
+ flight_data.maxVel = np.amax(flight_data.v.source[:, 1])
+ flight_result["maxVelocity"] = flight_data.maxVel
+
+ # Take care of parachute results
+ if len(flight_data.parachuteEvents) > 0:
+ flight_result["drogueTriggerTime"] = flight_data.parachuteEvents[0][0]
+ flight_result["drogueInflatedTime"] = (
+ flight_data.parachuteEvents[0][0] + flight_data.parachuteEvents[0][1].lag
+ )
+ flight_result["drogueInflatedVelocity"] = flight_data.v(
+ flight_data.parachuteEvents[0][0] + flight_data.parachuteEvents[0][1].lag
+ )
+ else:
+ flight_result["drogueTriggerTime"] = 0
+ flight_result["drogueInflatedTime"] = 0
+ flight_result["drogueInflatedVelocity"] = 0
+
+ # Write flight setting and results to file
+ dispersion_input_file.write(str(flight_setting) + "\n")
+ dispersion_output_file.write(str(flight_result) + "\n")
+
+ def export_flight_error(flight_setting):
+ dispersion_error_file.write(str(flight_setting) + "\n")
+
+ # Basic analysis info
+
+
+ # Create data files for inputs, outputs and error logging
+ dispersion_error_file = open(str(filename) + ".disp_errors.txt", "w")
+ dispersion_input_file = open(str(filename) + ".disp_inputs.txt", "w")
+ dispersion_output_file = open(str(filename) + ".disp_outputs.txt", "w")
+
+ # Initialize counter and timer
+ i = 0
+
+ initial_wall_time = time()
+ initial_cpu_time = process_time()
+
+ # Iterate over flight settings
+ out = display("Starting", display_id=True)
+ for setting in flight_settings(analysis_parameters, number_of_simulations):
+ start_time = process_time()
+ i += 1
+
+ # Creates copy of environment
+ envDispersion = self.environment
+
+ # Apply environment parameters variations on each iteration if possible
+ envDispersion.railLength = (
+ setting["railLength"]
+ if "railLength" in setting
+ else envDispersion.rL
+ )
+ envDispersion.gravity = (
+ setting["gravity"] if "gravity" in setting
+ else envDispersion.g
+ )
+ envDispersion.date = (
+ setting["date"] if "date" in setting
+ else envDispersion.date
+ )
+ envDispersion.latitude = (
+ setting["latitude"] if "latitude" in setting
+ else envDispersion.lat
+ )
+ envDispersion.longitude = (
+ setting["longitude"]
+ if "longitude" in setting
+ else envDispersion.lon
+ )
+ envDispersion.elevation = (
+ setting["elevation"]
+ if "elevation" in setting
+ else envDispersion.elevation
+ )
+ envDispersion.datum = (
+ setting["datum"] if "datum" in setting
+ else envDispersion.datum
+ )
+ if "ensembleMember" in setting:
+ envDispersion.selectEnsembleMember(setting["ensembleMember"])
+
+
+
+
+
+ # Creates copy of motor
+ motorDispersion = self.motor
+
+ # Apply motor parameters variations on each iteration if possible
+ motorDispersion = SolidMotor(
+ thrustSource = setting["thrustSource"] if "thrustSource" in setting else motorDispersion.thrustSource,
+ burnOut = setting["burnOut"] if "burnOut" in setting else motorDispersion.burnOut,
+ grainNumber = setting["grainNumber"] if "grainNumber" in setting else motorDispersion.grainNumber,
+ grainDensity = setting["grainDensity"] if "grainDensity" in setting else motorDispersion.grainDensity,
+ grainOuterRadius = setting["grainOuterRadius"] if "grainOuterRadius" in setting else motorDispersion.grainOuterRadius,
+ grainInitialInnerRadius = setting["grainInitialInnerRadius"] if "grainInitialInnerRadius" in setting else motorDispersion.grainInitialInnerRadius,
+ grainInitialHeight = setting["grainInitialHeight"] if "grainInitialHeight" in setting else motorDispersion.grainInitialHeight,
+ grainSeparation = setting["grainSeparation"] if "grainSeparation" in setting else motorDispersion.grainSeparation,
+ nozzleRadius = setting["nozzleRadius"] if "nozzleRadius" in setting else motorDispersion.nozzleRadius,
+ throatRadius = setting["throatRadius"] if "throatRadius" in setting else motorDispersion.throatRadius,
+ reshapeThrustCurve = setting["reshapeThrustCurve"] if "reshapeThrustCurve" in setting else motorDispersion.reshapeThrustCurve,
+ interpolationMethod = setting["interpolationMethod"] if "interpolationMethod" in setting else motorDispersion.interpolationMethod,
+ )
+
+
+
+ # Creates copy of rocket
+ rocketDispersion = self.rocket
+
+ # Apply rocket parameters variations on each iteration if possible
+ rocketDispersion = Rocket(
+ motor=motorDispersion,
+ mass = setting["rocketMass"] if "rocketMass" in setting else self.rocket.mass,
+ inertiaI = setting["inertiaI"] if "inertiaI" in setting else self.rocket.inertiaI,
+ inertiaZ = setting["inertiaZ"] if "inertiaZ" in setting else self.rocket.inertiaZ,
+ radius = setting["radius"] if "radius" in setting else self.rocket.radius,
+ distanceRocketNozzle = setting["distanceRocketNozzle"] if "distanceRocketNozzle" in setting else self.rocket.distanceRocketNozzle,
+ distanceRocketPropellant = setting["distanceRocketPropellant"] if "distanceRocketPropellant" in setting else self.rocket.distanceRocketPropellant,
+ powerOffDrag = setting["powerOffDrag"] if "powerOffDrag" in setting else self.rocket.powerOffDrag,
+ powerOnDrag = setting["powerOnDrag"] if "powerOnDrag" in setting else self.rocket.powerOnDrag,
+ )
+
+ # Add rocket nose, fins and tail
+ rocketDispersion.addNose(
+ length=setting["noseLength"] if "noseLength" in setting else self.rocket.noseLength,
+ kind=setting["noseKind"]if "noseKind" in setting else self.rocket.noseKind,
+ distanceToCM=setting["noseDistanceToCM"]if "noseDistanceToCM" in setting else self.rocket.noseDistanceToCM,
+ )
+ rocketDispersion.addFins(
+ n=setting["n"] if "n" in setting else self.rocket.numberOfFins,
+ rootChord=setting["rootChord"] if "rootChord" in setting else self.rocket.rootChord,
+ tipChord=setting["tipChord"] if "tipChord" in setting else self.rocket.tipChord,
+ span=setting["span"] if "span" in setting else self.rocket.span,
+ distanceToCM=setting["finDistanceToCM"] if "finDistanceToCM" in setting else self.rocket.distanceRocketFins,
+ radius=setting["radius"] if "radius" in setting else self.rocket.finRadius,
+ airfoil=setting["airfoil"] if "airfoil" in setting else self.rocket.finAirfoil
+ )
+ rocketDispersion.addTail(
+ topRadius =setting["topRadius"] if "topRadius" in setting else self.rocket.tailTopRadius,
+ bottomRadius=setting["bottomRadius"] if "bottomRadius" in setting else self.rocket.tailBottomRadius,
+ length=setting["length"] if "length" in setting else self.rocket.tailLength,
+ distanceToCM=setting["distanceToCM"] if "distanceToCM" in setting else self.rocket.tailDistanceToCM
+ )
+
+ # Add parachute
+ rocketDispersion.addParachute(
+ name=setting["name"] if "name" in setting else self.rocket.parachuteName,
+ CdS=setting["CdS"] if "CdS" in setting else self.rocket.parachuteCdS,
+ trigger=setting["trigger"] if "trigger" in setting else self.rocket.parachuteTrigger,
+ samplingRate=setting["samplingRate"] if "samplingRate" in setting else self.rocket.parachuteSamplingRate,
+ lag=setting["lag_rec"] if "lag_rec" in setting else self.rocket.lag_rec + setting["lag_se"] if "lag_se" in setting else self.rocket.parachuteLag,
+ noise=setting["noise"] if "noise" in setting else self.rocket.parachuteNoise,
+ )
+
+ rocketDispersion.setRailButtons(
+ distanceToCM = setting["RBdistanceToCM"] if "RBdistanceToCM" in setting else self.rocket.RBdistanceToCM,
+ angularPosition = setting["angularPosition"] if "angularPosition" in setting else self.rocket.angularPosition
+ )
+
+ # Run trajectory simulation
+ try:
+ TestFlight = Flight(
+ rocket=rocketDispersion,
+ environment=envDispersion,
+ inclination=setting["inclination"] if "inclination" in setting else self.flight.inclination,
+ heading= setting["heading"] if "heading" in setting else self.flight.heading,
+ #initialSolution=setting["initialSolution"] if "initialSolution" in setting else self.flight.initialSolution,
+ terminateOnApogee=setting["terminateOnApogee"] if "terminateOnApogee" in setting else self.flight.terminateOnApogee,
+ maxTime=setting["maxTime"] if "maxTime" in setting else self.flight.maxTime,
+ maxTimeStep=setting["maxTimeStep"] if "maxTimeStep" in setting else self.flight.maxTimeStep,
+ minTimeStep=setting["minTimeStep"] if "minTimeStep" in setting else self.flight.minTimeStep,
+ rtol=setting["rtol"] if "rtol" in setting else self.flight.rtol,
+ atol=setting["atol"] if "atol" in setting else self.flight.atol,
+ timeOvershoot=setting["timeOvershoot"] if "timeOvershoot" in setting else self.flight.timeOvershoot,
+ verbose=False,
+ )
+
+ export_flight_data(setting, TestFlight, process_time() - start_time)
+ except Exception as E:
+ print(E)
+ print(traceback.format_exc())
+ export_flight_error(setting)
+
+ # Register time
+ out.update(
+ f"Curent iteration: {i:06d} | Average Time per Iteration: {(process_time() - initial_cpu_time)/i:2.6f} s | Estimated time left: {int((number_of_simulations - i)*((process_time() - initial_cpu_time)/i))} s"
+ )
+
+ # Done
+
+ ## Print and save total time
+ final_string = f"Completed {i} iterations successfully. Total CPU time: {process_time() - initial_cpu_time} s. Total wall time: {time() - initial_wall_time} s"
+ out.update(final_string)
+ dispersion_input_file.write(final_string + "\n")
+ dispersion_output_file.write(final_string + "\n")
+ dispersion_error_file.write(final_string + "\n")
+
+ ## Close files
+ dispersion_input_file.close()
+ dispersion_output_file.close()
+ dispersion_error_file.close()
+
+ def importingDispersionResultsFromFile(self,dispersion_output_file):
+
+ # Initialize variable to store all results
+ dispersion_general_results = []
+
+ dispersion_results = {"outOfRailTime": [],
+ "outOfRailVelocity": [],
+ "apogeeTime": [],
+ "apogeeAltitude": [],
+ "apogeeX": [],
+ "apogeeY": [],
+ "impactTime": [],
+ "impactX": [],
+ "impactY": [],
+ "impactVelocity": [],
+ "initialStaticMargin": [],
+ "outOfRailStaticMargin": [],
+ "finalStaticMargin": [],
+ "numberOfEvents": [],
+ "maxVelocity": [],
+ "drogueTriggerTime": [],
+ "drogueInflatedTime": [],
+ "drogueInflatedVelocity": [],
+ "executionTime": [],
+ 'railDepartureAngleOfAttack': [],
+ 'lateralWind': [],
+ 'frontalWind': []}
+
+ # Get all dispersion results
+ # Get file
+ dispersion_output_file = open(dispersion_output_file, "r+")
+
+ # Read each line of the file and convert to dict
+ for line in dispersion_output_file:
+ # Skip comments lines
+ if line[0] != "{":
+ continue
+ # Eval results and store them
+ flight_result = eval(line)
+ dispersion_general_results.append(flight_result)
+ for parameter_key, parameter_value in flight_result.items():
+ dispersion_results[parameter_key].append(parameter_value)
+
+ # Close data file
+ dispersion_output_file.close()
+
+ #Number of flights simulated
+ self.N = len(dispersion_general_results)
+
+ return dispersion_results
+
+ def plotOutOfRailTime(self, dispersion_results):
+ print(
+ f'Out of Rail Time - Mean Value: {np.mean(dispersion_results["outOfRailTime"]):0.3f} s'
+ )
+ print(
+ f'Out of Rail Time - Standard Deviation: {np.std(dispersion_results["outOfRailTime"]):0.3f} s'
+ )
+
+ plt.figure()
+ plt.hist(dispersion_results["outOfRailTime"], bins=int(self.N**0.5))
+ plt.title("Out of Rail Time")
+ plt.xlabel("Time (s)")
+ plt.ylabel("Number of Occurences")
+ plt.show()
+
+ return None
+
+ def plotOutOfRailVelocity(self, dispersion_results):
+ print(
+ f'Out of Rail Velocity - Mean Value: {np.mean(dispersion_results["outOfRailVelocity"]):0.3f} m/s'
+ )
+ print(
+ f'Out of Rail Velocity - Standard Deviation: {np.std(dispersion_results["outOfRailVelocity"]):0.3f} m/s'
+ )
+
+ plt.figure()
+ plt.hist(dispersion_results["outOfRailVelocity"], bins=int(self.N**0.5))
+ plt.title("Out of Rail Velocity")
+ plt.xlabel("Velocity (m/s)")
+ plt.ylabel("Number of Occurences")
+ plt.show()
+
+ return None
+
+ def plotApogeeTime(self, dispersion_results):
+ print(
+ f'Impact Time - Mean Value: {np.mean(dispersion_results["impactTime"]):0.3f} s'
+ )
+ print(
+ f'Impact Time - Standard Deviation: {np.std(dispersion_results["impactTime"]):0.3f} s'
+ )
+
+ plt.figure()
+ plt.hist(dispersion_results["impactTime"], bins=int(self.N**0.5))
+ plt.title("Impact Time")
+ plt.xlabel("Time (s)")
+ plt.ylabel("Number of Occurences")
+ plt.show()
+
+ return None
+
+ def plotApogeeAltitude(self, dispersion_results):
+ print(
+ f'Apogee Altitude - Mean Value: {np.mean(dispersion_results["apogeeAltitude"]):0.3f} m'
+ )
+ print(
+ f'Apogee Altitude - Standard Deviation: {np.std(dispersion_results["apogeeAltitude"]):0.3f} m'
+ )
+
+ plt.figure()
+ plt.hist(dispersion_results["apogeeAltitude"], bins=int(self.N**0.5))
+ plt.title("Apogee Altitude")
+ plt.xlabel("Altitude (m)")
+ plt.ylabel("Number of Occurences")
+ plt.show()
+
+ return None
+
+ def plotApogeeXPosition(self, dispersion_results):
+ print(
+ f'Apogee X Position - Mean Value: {np.mean(dispersion_results["apogeeX"]):0.3f} m'
+ )
+ print(
+ f'Apogee X Position - Standard Deviation: {np.std(dispersion_results["apogeeX"]):0.3f} m'
+ )
+
+ plt.figure()
+ plt.hist(dispersion_results["apogeeX"], bins=int(self.N**0.5))
+ plt.title("Apogee X Position")
+ plt.xlabel("Apogee X Position (m)")
+ plt.ylabel("Number of Occurences")
+ plt.show()
+
+ return None
+
+ def plotApogeeYPosition(self, dispersion_results):
+ print(
+ f'Apogee Y Position - Mean Value: {np.mean(dispersion_results["apogeeY"]):0.3f} m'
+ )
+ print(
+ f'Apogee Y Position - Standard Deviation: {np.std(dispersion_results["apogeeY"]):0.3f} m'
+ )
+
+ plt.figure()
+ plt.hist(dispersion_results["apogeeY"], bins=int(self.N**0.5))
+ plt.title("Apogee Y Position")
+ plt.xlabel("Apogee Y Position (m)")
+ plt.ylabel("Number of Occurences")
+ plt.show()
+
+ return None
+
+ def plotImpactTime(self, dispersion_results):
+ print(
+ f'Impact Time - Mean Value: {np.mean(dispersion_results["impactTime"]):0.3f} s'
+ )
+ print(
+ f'Impact Time - Standard Deviation: {np.std(dispersion_results["impactTime"]):0.3f} s'
+ )
+
+ plt.figure()
+ plt.hist(dispersion_results["impactTime"], bins=int(self.N**0.5))
+ plt.title("Impact Time")
+ plt.xlabel("Time (s)")
+ plt.ylabel("Number of Occurences")
+ plt.show()
+
+ return None
+
+ def plotImpactXPosition(self, dispersion_results):
+ print(
+ f'Impact X Position - Mean Value: {np.mean(dispersion_results["impactX"]):0.3f} m'
+ )
+ print(
+ f'Impact X Position - Standard Deviation: {np.std(dispersion_results["impactX"]):0.3f} m'
+ )
+
+ plt.figure()
+ plt.hist(dispersion_results["impactX"], bins=int(self.N**0.5))
+ plt.title("Impact X Position")
+ plt.xlabel("Impact X Position (m)")
+ plt.ylabel("Number of Occurences")
+ plt.show()
+
+ return None
+
+ def plotImpactYPosition(self, dispersion_results):
+ print(
+ f'Impact Y Position - Mean Value: {np.mean(dispersion_results["impactY"]):0.3f} m'
+ )
+ print(
+ f'Impact Y Position - Standard Deviation: {np.std(dispersion_results["impactY"]):0.3f} m'
+ )
+
+ plt.figure()
+ plt.hist(dispersion_results["impactY"], bins=int(self.N**0.5))
+ plt.title("Impact Y Position")
+ plt.xlabel("Impact Y Position (m)")
+ plt.ylabel("Number of Occurences")
+ plt.show()
+
+ return None
+
+ def plotImpactVelocity(self, dispersion_results):
+ print(
+ f'Impact Velocity - Mean Value: {np.mean(dispersion_results["impactVelocity"]):0.3f} m/s'
+ )
+ print(
+ f'Impact Velocity - Standard Deviation: {np.std(dispersion_results["impactVelocity"]):0.3f} m/s'
+ )
+
+ plt.figure()
+ plt.hist(dispersion_results["impactVelocity"], bins=int(self.N**0.5))
+ plt.title("Impact Velocity")
+ plt.xlim(-35, 0)
+ plt.xlabel("Velocity (m/s)")
+ plt.ylabel("Number of Occurences")
+ plt.show()
+
+ return None
+
+ def plotStaticMargin(self, dispersion_results):
+ print(
+ f'Initial Static Margin - Mean Value: {np.mean(dispersion_results["initialStaticMargin"]):0.3f} c'
+ )
+ print(
+ f'Initial Static Margin - Standard Deviation: {np.std(dispersion_results["initialStaticMargin"]):0.3f} c'
+ )
+
+ print(
+ f'Out of Rail Static Margin - Mean Value: {np.mean(dispersion_results["outOfRailStaticMargin"]):0.3f} c'
+ )
+ print(
+ f'Out of Rail Static Margin - Standard Deviation: {np.std(dispersion_results["outOfRailStaticMargin"]):0.3f} c'
+ )
+
+ print(
+ f'Final Static Margin - Mean Value: {np.mean(dispersion_results["finalStaticMargin"]):0.3f} c'
+ )
+ print(
+ f'Final Static Margin - Standard Deviation: {np.std(dispersion_results["finalStaticMargin"]):0.3f} c'
+ )
+
+ plt.figure()
+ plt.hist(
+ dispersion_results["initialStaticMargin"],
+ label="Initial",
+ bins=int(self.N**0.5),
+ )
+ plt.hist(
+ dispersion_results["outOfRailStaticMargin"],
+ label="Out of Rail",
+ bins=int(self.N**0.5),
+ )
+ plt.hist(
+ dispersion_results["finalStaticMargin"], label="Final", bins=int(self.N**0.5)
+ )
+ plt.legend()
+ plt.title("Static Margin")
+ plt.xlabel("Static Margin (c)")
+ plt.ylabel("Number of Occurences")
+ plt.show()
+
+ return None
+
+ def plotMaximumVelocity(self, dispersion_results):
+ print(
+ f'Maximum Velocity - Mean Value: {np.mean(dispersion_results["maxVelocity"]):0.3f} m/s'
+ )
+ print(
+ f'Maximum Velocity - Standard Deviation: {np.std(dispersion_results["maxVelocity"]):0.3f} m/s'
+ )
+
+ plt.figure()
+ plt.hist(dispersion_results["maxVelocity"], bins=int(self.N**0.5))
+ plt.title("Maximum Velocity")
+ plt.xlabel("Velocity (m/s)")
+ plt.ylabel("Number of Occurences")
+ plt.show()
+
+ return None
+
+ def plotNumberOfParachuteEvents(self, dispersion_results):
+ plt.figure()
+ plt.hist(dispersion_results["numberOfEvents"])
+ plt.title("Parachute Events")
+ plt.xlabel("Number of Parachute Events")
+ plt.ylabel("Number of Occurences")
+ plt.show()
+
+ return None
+
+ def plotDrogueTriggerTime(self, dispersion_results):
+ print(
+ f'Drogue Trigger Time - Mean Value: {np.mean(dispersion_results["drogueTriggerTime"]):0.3f} s'
+ )
+ print(
+ f'Drogue Trigger Time - Standard Deviation: {np.std(dispersion_results["drogueTriggerTime"]):0.3f} s'
+ )
+
+ plt.figure()
+ plt.hist(dispersion_results["drogueTriggerTime"], bins=int(self.N**0.5))
+ plt.title("Drogue Trigger Time")
+ plt.xlabel("Time (s)")
+ plt.ylabel("Number of Occurences")
+ plt.show()
+
+ return None
+
+ def plotDrogueFullyInflatedTime(self, dispersion_results):
+ print(
+ f'Drogue Fully Inflated Time - Mean Value: {np.mean(dispersion_results["drogueInflatedTime"]):0.3f} s'
+ )
+ print(
+ f'Drogue Fully Inflated Time - Standard Deviation: {np.std(dispersion_results["drogueInflatedTime"]):0.3f} s'
+ )
+
+ plt.figure()
+ plt.hist(dispersion_results["drogueInflatedTime"], bins=int(self.N**0.5))
+ plt.title("Drogue Fully Inflated Time")
+ plt.xlabel("Time (s)")
+ plt.ylabel("Number of Occurences")
+ plt.show()
+
+ return None
+
+ def plotDrogueFullyVelocity(self, dispersion_results):
+ print(
+ f'Drogue Parachute Fully Inflated Velocity - Mean Value: {np.mean(dispersion_results["drogueInflatedVelocity"]):0.3f} m/s'
+ )
+ print(
+ f'Drogue Parachute Fully Inflated Velocity - Standard Deviation: {np.std(dispersion_results["drogueInflatedVelocity"]):0.3f} m/s'
+ )
+
+ plt.figure()
+ plt.hist(dispersion_results["drogueInflatedVelocity"], bins=int(self.N**0.5))
+ plt.title("Drogue Parachute Fully Inflated Velocity")
+ plt.xlabel("Velocity m/s)")
+ plt.ylabel("Number of Occurences")
+ plt.show()
+
+ return None
+
+ def plotEllipses(self, dispersion_results, image, realLandingPoint):
+
+ # Import background map
+ img = imread(image)
+
+ # Retrieve dispersion data por apogee and impact XY position
+ apogeeX = np.array(dispersion_results["apogeeX"])
+ apogeeY = np.array(dispersion_results["apogeeY"])
+ impactX = np.array(dispersion_results["impactX"])
+ impactY = np.array(dispersion_results["impactY"])
+
+ # Define function to calculate eigen values
+ def eigsorted(cov):
+ vals, vecs = np.linalg.eigh(cov)
+ order = vals.argsort()[::-1]
+ return vals[order], vecs[:, order]
+
+ # Create plot figure
+ plt.figure(num=None, figsize=(8, 6), dpi=150, facecolor="w", edgecolor="k")
+ ax = plt.subplot(111)
+
+ # Calculate error ellipses for impact
+ impactCov = np.cov(impactX, impactY)
+ impactVals, impactVecs = eigsorted(impactCov)
+ impactTheta = np.degrees(np.arctan2(*impactVecs[:, 0][::-1]))
+ impactW, impactH = 2 * np.sqrt(impactVals)
+
+ # Draw error ellipses for impact
+ impact_ellipses = []
+ for j in [1, 2, 3]:
+ impactEll = Ellipse(
+ xy=(np.mean(impactX), np.mean(impactY)),
+ width=impactW * j,
+ height=impactH * j,
+ angle=impactTheta,
+ color="black",
+ )
+ impactEll.set_facecolor((0, 0, 1, 0.2))
+ impact_ellipses.append(impactEll)
+ ax.add_artist(impactEll)
+
+ # Calculate error ellipses for apogee
+ apogeeCov = np.cov(apogeeX, apogeeY)
+ apogeeVals, apogeeVecs = eigsorted(apogeeCov)
+ apogeeTheta = np.degrees(np.arctan2(*apogeeVecs[:, 0][::-1]))
+ apogeeW, apogeeH = 2 * np.sqrt(apogeeVals)
+
+ # Draw error ellipses for apogee
+ for j in [1, 2, 3]:
+ apogeeEll = Ellipse(
+ xy=(np.mean(apogeeX), np.mean(apogeeY)),
+ width=apogeeW * j,
+ height=apogeeH * j,
+ angle=apogeeTheta,
+ color="black",
+ )
+ apogeeEll.set_facecolor((0, 1, 0, 0.2))
+ ax.add_artist(apogeeEll)
+
+ # Draw launch point
+ plt.scatter(0, 0, s=30, marker="*", color="black", label="Launch Point")
+ # Draw apogee points
+ plt.scatter(
+ apogeeX, apogeeY, s=5, marker="^", color="green", label="Simulated Apogee"
+ )
+ # Draw impact points
+ plt.scatter(
+ impactX,
+ impactY,
+ s=5,
+ marker="v",
+ color="blue",
+ label="Simulated Landing Point",
+ )
+ # Draw real landing point
+ if realLandingPoint != None:
+ plt.scatter(
+ realLandingPoint[0],
+ realLandingPoint[1],
+ s=20,
+ marker="X",
+ color="red",
+ label="Measured Landing Point",
+ )
+
+ plt.legend()
+
+ # Add title and labels to plot
+ ax.set_title(
+ "1$\sigma$, 2$\sigma$ and 3$\sigma$ Dispersion Ellipses: Apogee and Lading Points"
+ )
+ ax.set_ylabel("North (m)")
+ ax.set_xlabel("East (m)")
+
+ # Add background image to plot
+ # You can translate the basemap by changing dx and dy (in meters)
+ dx = 0
+ dy = 0
+ plt.imshow(img, zorder=0, extent=[-3000 - dx, 3000 - dx, -3000 - dy, 3000 - dy])
+ plt.axhline(0, color="black", linewidth=0.5)
+ plt.axvline(0, color="black", linewidth=0.5)
+ plt.xlim(-3000, 3000)
+ plt.ylim(-3000, 3000)
+
+ # Save plot and show result
+ plt.savefig(str(self.filename) + ".pdf", bbox_inches="tight", pad_inches=0)
+ plt.savefig(str(self.filename) + ".svg", bbox_inches="tight", pad_inches=0)
+ plt.show()
+
+ def plotLateralWindSpeed(self, dispersion_results):
+ print(f'Lateral Surface Wind Speed - Mean Value: {np.mean(dispersion_results["lateralWind"]):0.3f} m/s')
+ print(f'Lateral Surface Wind Speed - Standard Deviation: {np.std(dispersion_results["lateralWind"]):0.3f} m/s')
+
+ plt.figure()
+ plt.hist(dispersion_results["lateralWind"], bins=int(self.N**0.5))
+ plt.title('Lateral Surface Wind Speed')
+ plt.xlabel('Velocity (m/s)')
+ plt.ylabel('Number of Occurences')
+ plt.show()
+
+ def plotFrontalWindSpeed(self, dispersion_results):
+ print(f'Frontal Surface Wind Speed - Mean Value: {np.mean(dispersion_results["frontalWind"]):0.3f} m/s')
+ print(f'Frontal Surface Wind Speed - Standard Deviation: {np.std(dispersion_results["frontalWind"]):0.3f} m/s')
+
+ plt.figure()
+ plt.hist(dispersion_results["frontalWind"], bins=int(self.N**0.5))
+ plt.title('Frontal Surface Wind Speed')
+ plt.xlabel('Velocity (m/s)')
+ plt.ylabel('Number of Occurences')
+ plt.show()
+
+ def info(self):
+ dispersion_results = self.importingDispersionResultsFromFile(self.filename)
+
+ self.plotEllipses(dispersion_results, self.image, self.realLandingPoint)
+
+ self.plotApogeeAltitude(dispersion_results)
+
+ self.plotOutOfRailVelocity(dispersion_results)
+
+ self.plotStaticMargin(dispersion_results)
+
+ self.plotLateralWindSpeed(dispersion_results)
+
+ self.plotFrontalWindSpeed(dispersion_results)
+
+ self.plotOutOfRailTime(dispersion_results)
+
+ self.plotApogeeTime(dispersion_results)
+
+ self.plotApogeeXPosition(dispersion_results)
+
+ self.plotApogeeYPosition(dispersion_results)
+
+ self.plotImpactTime(dispersion_results)
+
+ self.plotImpactVelocity(dispersion_results)
+
+ self.plotImpactXPosition(dispersion_results)
+
+ self.plotImpactYPosition(dispersion_results)
+
+ self.plotMaximumVelocity(dispersion_results)
+
+ self.plotNumberOfParachuteEvents(dispersion_results)
+
+ self.plotDrogueFullyInflatedTime(dispersion_results)
+
+ self.plotDrogueFullyVelocity(dispersion_results)
+
+ self.plotDrogueTriggerTime(dispersion_results)
+
+
+
\ No newline at end of file
diff --git a/rocketpy/Environment.py b/rocketpy/Environment.py
index 0292b0e1c..9ae378b8b 100644
--- a/rocketpy/Environment.py
+++ b/rocketpy/Environment.py
@@ -358,6 +358,7 @@ def __init__(
self.datum = datum
# Save date
+ self.date = date
if date != None:
self.setDate(date, timeZone)
else:
@@ -373,6 +374,8 @@ def __init__(
self.setAtmosphericModel("StandardAtmosphere")
# Save latitude and longitude
+ self.lat = latitude
+ self.lon = longitude
if latitude != None and longitude != None:
self.setLocation(latitude, longitude)
else:
@@ -389,6 +392,7 @@ def __init__(
self.initialEW = convert[5]
# Save elevation
+ self.elevation = elevation
self.setElevation(elevation)
# Recalculate Earth Radius
diff --git a/rocketpy/Flight.py b/rocketpy/Flight.py
index 50ac5d97d..2e21dfb54 100644
--- a/rocketpy/Flight.py
+++ b/rocketpy/Flight.py
@@ -603,6 +603,7 @@ def __init__(
self.initialSolution = initialSolution
self.timeOvershoot = timeOvershoot
self.terminateOnApogee = terminateOnApogee
+ self.verbose = verbose
# Modifying Rail Length for a better out of rail condition
upperRButton = max(self.rocket.railButtons[0])
diff --git a/rocketpy/Motor.py b/rocketpy/Motor.py
index be9d62707..7fc605cdb 100644
--- a/rocketpy/Motor.py
+++ b/rocketpy/Motor.py
@@ -151,6 +151,16 @@ def __init__(
-------
None
"""
+
+ #Save parameters
+ self.thrustSource = thrustSource
+ self.burnOut = burnOut
+ self.nozzleRadius = nozzleRadius
+ self.throatRadius = throatRadius
+ self.reshapeThrustCurve = reshapeThrustCurve
+ self.interpolationMethod = interpolationMethod
+
+
# Thrust parameters
self.interpolate = interpolationMethod
self.burnOutTime = burnOut
diff --git a/rocketpy/Rocket.py b/rocketpy/Rocket.py
index 80d93fc89..b8181489b 100644
--- a/rocketpy/Rocket.py
+++ b/rocketpy/Rocket.py
@@ -431,6 +431,13 @@ def addNose(self, length, kind, distanceToCM):
self : Rocket
Object of the Rocket class.
"""
+
+ # Save parameters for Dispersion
+ self.noseLength = length
+ self.noseKind = kind
+ self.noseDistanceToCM = distanceToCM
+
+
# Analyze type
if kind == "conical":
k = 1 - 1 / 3
@@ -524,6 +531,11 @@ def addFins(
Object of the Rocket class.
"""
+ # Save parameters for Dispersion
+ self.numberOfFins = n
+ self.finRadius = radius
+ self.finAirfoil = airfoil
+
# Retrieve parameters for calculations
Cr = rootChord
Ct = tipChord
@@ -749,6 +761,15 @@ def addParachute(
noiseSignal and noisyPressureSignal which are filled in during
Flight simulation.
"""
+
+ # Save parameters for Dispersion
+ self.parachuteName = name
+ self.parachuteCdS = CdS
+ self.parachuteTrigger = trigger
+ self.parachuteSamplingRate = samplingRate
+ self.parachuteLag = lag
+ self.parachuteNoise = noise
+
# Create a parachute
parachute = Parachute(name, CdS, trigger, samplingRate, lag, noise)
@@ -787,11 +808,15 @@ def setRailButtons(self, distanceToCM, angularPosition=45):
-------
None
"""
+
+
# Order distance to CM
if distanceToCM[0] < distanceToCM[1]:
distanceToCM.reverse()
# Save
self.railButtons = self.railButtonPair(distanceToCM, angularPosition)
+ self.RBdistanceToCM = distanceToCM
+ self.angularPosition = angularPosition
return None
diff --git a/rocketpy/__init__.py b/rocketpy/__init__.py
index 88f9cba9f..6446fc604 100644
--- a/rocketpy/__init__.py
+++ b/rocketpy/__init__.py
@@ -27,4 +27,5 @@
from .Function import Function
from .Motor import HybridMotor, SolidMotor
from .Rocket import Rocket
+from .Dispersion import Dispersion
from .utilities import *
diff --git a/teste.disp_errors.txt b/teste.disp_errors.txt
new file mode 100644
index 000000000..63ae3e1a5
--- /dev/null
+++ b/teste.disp_errors.txt
@@ -0,0 +1 @@
+Completed 10 iterations successfully. Total CPU time: 17.6875 s. Total wall time: 19.050102710723877 s
diff --git a/teste.disp_inputs.txt b/teste.disp_inputs.txt
new file mode 100644
index 000000000..458975e50
--- /dev/null
+++ b/teste.disp_inputs.txt
@@ -0,0 +1,11 @@
+{'rocketMass': 19.90652137182852, 'burnOut': 4.122751149659427, 'railLength': 5.565123045049819}
+{'rocketMass': 20.040442434623223, 'burnOut': 4.040194754660794, 'railLength': 5.4460671407624535}
+{'rocketMass': 19.930837701395934, 'burnOut': 3.468467598491397, 'railLength': 5.150566738117259}
+{'rocketMass': 19.8506753734651, 'burnOut': 3.0803468079703005, 'railLength': 5.821128682799574}
+{'rocketMass': 20.040395390137192, 'burnOut': 3.4398519219680628, 'railLength': 5.0211721343161}
+{'rocketMass': 19.92409646924012, 'burnOut': 4.114787103323786, 'railLength': 5.2566996814270395}
+{'rocketMass': 19.83694037126801, 'burnOut': 4.363111003241498, 'railLength': 5.428369692561595}
+{'rocketMass': 20.170452293227846, 'burnOut': 4.235807754853636, 'railLength': 5.049820100552791}
+{'rocketMass': 20.168066054386124, 'burnOut': 4.369091258708742, 'railLength': 5.921815945979148}
+{'rocketMass': 19.91857803920068, 'burnOut': 4.342779134206488, 'railLength': 5.080259704139901}
+Completed 10 iterations successfully. Total CPU time: 17.6875 s. Total wall time: 19.050102710723877 s
diff --git a/teste.disp_outputs.txt b/teste.disp_outputs.txt
new file mode 100644
index 000000000..f4c3a0e35
--- /dev/null
+++ b/teste.disp_outputs.txt
@@ -0,0 +1,11 @@
+{'outOfRailTime': 0.39364752882342674, 'outOfRailVelocity': 23.393901564816144, 'apogeeTime': 23.48209630922026, 'apogeeAltitude': 2567.5769822795232, 'apogeeX': -241.77524247257443, 'apogeeY': 877.5895761489942, 'impactTime': 150.4764074480069, 'impactX': 547.4755585200927, 'impactY': 572.4876365926631, 'impactVelocity': -19.462199163839447, 'initialStaticMargin': 2.2172168363626756, 'outOfRailStaticMargin': 2.2959121643050797, 'finalStaticMargin': 3.0897188014956276, 'numberOfEvents': 1, 'drogueTriggerTime': 23.485714285714288, 'drogueInflatedTime': 24.985714285714288, 'drogueInflatedVelocity': 39.7631366423571, 'executionTime': 1.328125, 'lateralWind': -3.7535914518728144, 'frontalWind': -6.317898247256396, 'maxVelocity': 235.99104212868633}
+{'outOfRailTime': 0.3948759737597698, 'outOfRailVelocity': 23.319771149981346, 'apogeeTime': 23.4028291433462, 'apogeeAltitude': 2545.3078375989367, 'apogeeX': -241.06417203223654, 'apogeeY': 872.8585625767016, 'impactTime': 148.97167317343016, 'impactX': 537.0216691616343, 'impactY': 571.393383082964, 'impactVelocity': -19.52766516726175, 'initialStaticMargin': 2.2222979195834713, 'outOfRailStaticMargin': 2.300890920849731, 'finalStaticMargin': 3.0897188014956276, 'numberOfEvents': 1, 'drogueTriggerTime': 23.40952380952381, 'drogueInflatedTime': 24.90952380952381, 'drogueInflatedVelocity': 39.77449604815257, 'executionTime': 1.34375, 'lateralWind': -3.7535914518728144, 'frontalWind': -6.317898247256396, 'maxVelocity': 234.4425381986718}
+{'outOfRailTime': 0.39383938438037336, 'outOfRailVelocity': 23.379426763634843, 'apogeeTime': 23.46794207038415, 'apogeeAltitude': 2563.5792893018133, 'apogeeX': -241.66211213218298, 'apogeeY': 876.7672315674648, 'impactTime': 150.20226074219605, 'impactX': 545.5025868758161, 'impactY': 572.4551802719707, 'impactVelocity': -19.47410279394695, 'initialStaticMargin': 2.2181438376033484, 'outOfRailStaticMargin': 2.2968129802515764, 'finalStaticMargin': 3.0897188014956276, 'numberOfEvents': 1, 'drogueTriggerTime': 23.476190476190478, 'drogueInflatedTime': 24.976190476190478, 'drogueInflatedVelocity': 39.78032246823763, 'executionTime': 1.703125, 'lateralWind': -3.7535914518728144, 'frontalWind': -6.317898247256396, 'maxVelocity': 235.71407446451917}
+{'outOfRailTime': 0.3930704664495948, 'outOfRailVelocity': 23.418332470677253, 'apogeeTime': 23.515367492889425, 'apogeeAltitude': 2576.967280631924, 'apogeeX': -242.10277007283392, 'apogeeY': 879.6372495529988, 'impactTime': 151.10771800521675, 'impactX': 551.7575545784338, 'impactY': 573.2154731676559, 'impactVelocity': -19.43483438695117, 'initialStaticMargin': 2.215080359624552, 'outOfRailStaticMargin': 2.2938029229466124, 'finalStaticMargin': 3.0897188014956276, 'numberOfEvents': 1, 'drogueTriggerTime': 23.523809523809526, 'drogueInflatedTime': 25.023809523809526, 'drogueInflatedVelocity': 39.78158608378997, 'executionTime': 1.4375, 'lateralWind': -3.7535914518728144, 'frontalWind': -6.317898247256396, 'maxVelocity': 236.64866640993696}
+{'outOfRailTime': 0.39500696060662027, 'outOfRailVelocity': 23.325764409221105, 'apogeeTime': 23.402730892947652, 'apogeeAltitude': 2545.268387833093, 'apogeeX': -241.04337919100507, 'apogeeY': 872.8104354829267, 'impactTime': 148.96979322782306, 'impactX': 537.0293697505082, 'impactY': 571.3452280041232, 'impactVelocity': -19.527642212283883, 'initialStaticMargin': 2.2222961450651133, 'outOfRailStaticMargin': 2.3009208054590156, 'finalStaticMargin': 3.0897188014956276, 'numberOfEvents': 1, 'drogueTriggerTime': 23.40952380952381, 'drogueInflatedTime': 24.90952380952381, 'drogueInflatedVelocity': 39.77324414527358, 'executionTime': 2.15625, 'lateralWind': -3.7535914518728144, 'frontalWind': -6.317898247256396, 'maxVelocity': 234.43854100088942}
+{'outOfRailTime': 0.393786187885889, 'outOfRailVelocity': 23.383441067614577, 'apogeeTime': 23.47201858836762, 'apogeeAltitude': 2564.719672077175, 'apogeeX': -241.6949708415784, 'apogeeY': 877.0039280151173, 'impactTime': 150.28163849653603, 'impactX': 546.1156054888656, 'impactY': 572.3842677454285, 'impactVelocity': -19.47080301448961, 'initialStaticMargin': 2.2178870418425163, 'outOfRailStaticMargin': 2.2965634374473396, 'finalStaticMargin': 3.0897188014956276, 'numberOfEvents': 1, 'drogueTriggerTime': 23.476190476190478, 'drogueInflatedTime': 24.976190476190478, 'drogueInflatedVelocity': 39.76614816730068, 'executionTime': 2.546875, 'lateralWind': -3.7535914518728144, 'frontalWind': -6.317898247256396, 'maxVelocity': 235.79377306150013}
+{'outOfRailTime': 0.3929585800142951, 'outOfRailVelocity': 23.426281445945772, 'apogeeTime': 23.52347772692958, 'apogeeAltitude': 2579.2615570345256, 'apogeeX': -242.17290460742268, 'apogeeY': 880.1177661016592, 'impactTime': 151.26832676009747, 'impactX': 552.9942843408794, 'impactY': 573.080316911245, 'impactVelocity': -19.42810182800753, 'initialStaticMargin': 2.2145533013720162, 'outOfRailStaticMargin': 2.2932898618648787, 'finalStaticMargin': 3.0897188014956276, 'numberOfEvents': 1, 'drogueTriggerTime': 23.523809523809526, 'drogueInflatedTime': 25.023809523809526, 'drogueInflatedVelocity': 39.753529828392345, 'executionTime': 1.953125, 'lateralWind': -3.7535914518728144, 'frontalWind': -6.317898247256396, 'maxVelocity': 236.80796514901837}
+{'outOfRailTime': 0.3955845155910141, 'outOfRailVelocity': 23.239490423502964, 'apogeeTime': 23.32661064243672, 'apogeeAltitude': 2524.1052140855963, 'apogeeX': -240.41223767499687, 'apogeeY': 868.4084069471148, 'impactTime': 147.5418539513837, 'impactX': 527.1436049410614, 'impactY': 570.3680499438507, 'impactVelocity': -19.591010733074842, 'initialStaticMargin': 2.2271743136835758, 'outOfRailStaticMargin': 2.305552942215501, 'finalStaticMargin': 3.0897188014956276, 'numberOfEvents': 1, 'drogueTriggerTime': 23.333333333333336, 'drogueInflatedTime': 24.833333333333336, 'drogueInflatedVelocity': 39.778177727954976, 'executionTime': 1.1875, 'lateralWind': -3.7535914518728144, 'frontalWind': -6.317898247256396, 'maxVelocity': 232.9728405593581}
+{'outOfRailTime': 0.395702891244475, 'outOfRailVelocity': 23.25256596431927, 'apogeeTime': 23.328153215748333, 'apogeeAltitude': 2524.534808308825, 'apogeeX': -240.40043656304366, 'apogeeY': 868.4449324828357, 'impactTime': 147.5709274420821, 'impactX': 527.396121496818, 'impactY': 570.2834551160805, 'impactVelocity': -19.589849853526715, 'initialStaticMargin': 2.2270853048997874, 'outOfRailStaticMargin': 2.305499311444709, 'finalStaticMargin': 3.0897188014956276, 'numberOfEvents': 1, 'drogueTriggerTime': 23.333333333333336, 'drogueInflatedTime': 24.833333333333336, 'drogueInflatedVelocity': 39.770636265523834, 'executionTime': 1.71875, 'lateralWind': -3.7535914518728144, 'frontalWind': -6.317898247256396, 'maxVelocity': 233.001280524876}
+{'outOfRailTime': 0.3937427949045639, 'outOfRailVelocity': 23.386723682784606, 'apogeeTime': 23.47499691340987, 'apogeeAltitude': 2565.5736138519483, 'apogeeX': -241.71498339389854, 'apogeeY': 877.1714488168865, 'impactTime': 150.34230311624643, 'impactX': 546.5900688232151, 'impactY': 572.3173877624449, 'impactVelocity': -19.468105805968417, 'initialStaticMargin': 2.217676713980707, 'outOfRailStaticMargin': 2.2963590866966, 'finalStaticMargin': 3.0897188014956276, 'numberOfEvents': 1, 'drogueTriggerTime': 23.476190476190478, 'drogueInflatedTime': 24.976190476190478, 'drogueInflatedVelocity': 39.75527823501804, 'executionTime': 2.265625, 'lateralWind': -3.7535914518728144, 'frontalWind': -6.317898247256396, 'maxVelocity': 235.85313212739368}
+Completed 10 iterations successfully. Total CPU time: 17.6875 s. Total wall time: 19.050102710723877 s
diff --git a/teste.disp_outputs.txt.svg b/teste.disp_outputs.txt.svg
new file mode 100644
index 000000000..587db8574
--- /dev/null
+++ b/teste.disp_outputs.txt.svg
@@ -0,0 +1,1320 @@
+
+
+