From 20d424bd2ba576212f9162f3af0e7b316e3ea94b Mon Sep 17 00:00:00 2001 From: Kartik Kumar Date: Tue, 11 Oct 2016 16:38:51 +0200 Subject: [PATCH 1/8] Add missing include statements for SQLite library to access macros (due to change in SQLiteCpp). --- src/atomScanner.cpp | 2 ++ src/j2Analysis.cpp | 2 ++ src/lambertFetch.cpp | 2 ++ src/lambertScanner.cpp | 8 +++++--- src/sgp4Scanner.cpp | 2 ++ 5 files changed, 13 insertions(+), 3 deletions(-) diff --git a/src/atomScanner.cpp b/src/atomScanner.cpp index 5984f15..b10156b 100644 --- a/src/atomScanner.cpp +++ b/src/atomScanner.cpp @@ -23,6 +23,8 @@ #include #include +#include + #include #include diff --git a/src/j2Analysis.cpp b/src/j2Analysis.cpp index 8ed768c..6288763 100644 --- a/src/j2Analysis.cpp +++ b/src/j2Analysis.cpp @@ -17,6 +17,8 @@ #include #include +#include + #include #include diff --git a/src/lambertFetch.cpp b/src/lambertFetch.cpp index d7e153c..187027a 100644 --- a/src/lambertFetch.cpp +++ b/src/lambertFetch.cpp @@ -10,6 +10,8 @@ #include +#include + #include #include diff --git a/src/lambertScanner.cpp b/src/lambertScanner.cpp index 073e8e0..13bcb76 100644 --- a/src/lambertScanner.cpp +++ b/src/lambertScanner.cpp @@ -21,6 +21,8 @@ #include #include +#include + #include #include @@ -191,7 +193,7 @@ void executeLambertScanner( const rapidjson::Document& config ) { DateTime departureEpoch = input.departureEpochInitial; departureEpoch = departureEpoch.AddSeconds( input.departureEpochStepSize * m ); - + const Eci tleDepartureState = sgp4Departure.FindPosition( departureEpoch ); const Vector6 departureState = getStateVector( tleDepartureState ); @@ -209,7 +211,7 @@ void executeLambertScanner( const rapidjson::Document& config ) = astro::convertCartesianToKeplerianElements( departureState, earthGravitationalParameter ); const int departureObjectId = static_cast< int >( departureObject.NoradNumber( ) ); - + // Loop over time-of-flight grid. for ( int k = 0; k < input.timeOfFlightSteps; k++ ) { @@ -354,7 +356,7 @@ void executeLambertScanner( const rapidjson::Document& config ) // Reset SQL insert query. query.reset( ); } - } + } } ++showProgress; diff --git a/src/sgp4Scanner.cpp b/src/sgp4Scanner.cpp index dc53987..1738731 100644 --- a/src/sgp4Scanner.cpp +++ b/src/sgp4Scanner.cpp @@ -27,6 +27,8 @@ #include #include +#include + #include #include From b28aec8621110b56e3c8434fc400087ebc4aa610 Mon Sep 17 00:00:00 2001 From: Enne Hekma Date: Tue, 29 Mar 2016 11:36:06 +0200 Subject: [PATCH 2/8] Add tle plotting scripts --- python/plotcatalog.py | 149 +++++++++++++++++++++++++++++++++++++++ python/twoBodyMethods.py | 33 +++++++++ 2 files changed, 182 insertions(+) create mode 100644 python/plotcatalog.py create mode 100644 python/twoBodyMethods.py diff --git a/python/plotcatalog.py b/python/plotcatalog.py new file mode 100644 index 0000000..81fc954 --- /dev/null +++ b/python/plotcatalog.py @@ -0,0 +1,149 @@ +''' +Copyright (c) 2014, K. Kumar (me@kartikkumar.com) +All rights reserved. +''' + +################################################################################################### +# Set up input deck +################################################################################################### + +# Set path to TLE catalog file. +tleCatalogFilePath = "../data/catalogs/catalog_rocketbodiesLEOfull.txt" + +# Set number of lines per entry in TLE catalog (2 or 3). +tleEntryNumberOfLines = 3 + +# Set path to output directory. +outputPath = "../data/napacatalog/" + +# Set figure DPI. +figureDPI = 300 + +# Set font size for axes labels. +fontSize = 24 + +################################################################################################### + +''' + DO NOT EDIT PARAMETERS BEYOND THIS POINT!!! +''' + +################################################################################################### +# Set up modules and packages +################################################################################################### + +from sgp4.earth_gravity import wgs72 +from sgp4.io import twoline2rv +from sgp4.propagation import getgravconst + +from matplotlib import rcParams +import matplotlib.pyplot as plt +import numpy as np + +from twoBodyMethods import convertMeanMotionToSemiMajorAxis + +################################################################################################### + +################################################################################################### +# Read and store TLE catalog +################################################################################################### + +# Read in catalog and store lines in list. +fileHandle = open(tleCatalogFilePath) +catalogLines = fileHandle.readlines() +fileHandle.close() + +# Strip newline and return carriage characters. +for i in xrange(len(catalogLines)): + catalogLines[i] = catalogLines[i].strip('\r\n') + +# Parse TLE entries and store debris objects. +debrisObjects = [] +for tleEntry in xrange(0,len(catalogLines),tleEntryNumberOfLines): + debrisObjects.append(twoline2rv(catalogLines[tleEntry+1], catalogLines[tleEntry+2], wgs72)) + +# Sort list of debris objects based on inclination. +inclinationSortedObjects = sorted(debrisObjects, key=lambda x: x.inclo, reverse=False) + +inclinations = [] +raan = [] +ecc = [] +aop = [] +for i in xrange(len(inclinationSortedObjects)-1): + inclinations.append(inclinationSortedObjects[i].inclo) + raan.append(inclinationSortedObjects[i].nodeo) + ecc.append(inclinationSortedObjects[i].ecco) + aop.append(inclinationSortedObjects[i].argpo) + +################################################################################################### + +################################################################################################### +# Generate plots +################################################################################################### + +# Set font size for plot labels. +rcParams.update({'font.size': fontSize}) + +# Plot distribution of eccentricity [-] against semi-major axis [km]. +figure = plt.figure() +axis = figure.add_subplot(111) +plt.xlabel("Semi-major axis [km]") +plt.ylabel("Eccentricity [-]") +plt.ticklabel_format(style='sci', axis='x', scilimits=(0,0)) +plt.plot([convertMeanMotionToSemiMajorAxis(debrisObject.no/60.0, getgravconst('wgs72')[1]) \ + for debrisObject in debrisObjects], \ + [debrisObject.ecco for debrisObject in debrisObjects], \ + marker='o', markersize=5, linestyle='none') +axis.set_xlim(xmax=0.8e4) +figure.set_tight_layout(True) +plt.savefig(outputPath + "/figure1_debrisPopulation_eccentricityVsSemiMajorAxis.pdf", \ + dpi = figureDPI) +plt.close() + +# Plot components of eccentricity vector [-]. +figure = plt.figure() +axis = figure.add_subplot(111) +plt.xlabel("$e \cos{\omega}$ [-]") +plt.ylabel("$e \sin{\omega}$ [-]") +plt.plot(ecc*np.cos(aop),ecc*np.sin(aop), marker='o', markersize=5, linestyle='none') +plt.axis('equal') +axis.set_xlim(xmin=-.21, xmax=.21) +axis.set_ylim(ymin=-.21, ymax=.21) +axis.set(xticks=[-.2,-.1,0,.1,.2])#, xticklabels=datelabels) #Same as plt.xticks +figure.set_tight_layout(True) +plt.savefig(outputPath + "/figure2_debrisPopulation_eccentricityVector.pdf", dpi = figureDPI) +plt.savefig(outputPath + "/figure2_debrisPopulation_eccentricityVector.pdf", dpi = figureDPI) +plt.close() + +# Plot distribution of inclination [deg] against semi-major axis [km]. +figure = plt.figure() +axis = figure.add_subplot(111) +plt.xlabel("Semi-major axis [km]") +plt.ylabel("Inclination [deg]") +plt.ticklabel_format(style='sci', axis='x', scilimits=(0,0)) +plt.plot([convertMeanMotionToSemiMajorAxis(debrisObject.no/60.0, getgravconst('wgs72')[1]) \ + for debrisObject in debrisObjects], \ + [np.rad2deg(debrisObject.inclo) for debrisObject in debrisObjects], \ + marker='o', markersize=5, linestyle='none') +axis.set_xlim(xmax=0.8e4) +figure.set_tight_layout(True) +plt.savefig(outputPath + "/figure3_debrisPopulation_inclinationVsSemiMajorAxis.pdf", \ + dpi = figureDPI) +plt.close() + +# Plot components of inclination vector [deg]. +figure = plt.figure() +axis = figure.add_subplot(111) +plt.xlabel("$i \cos{\Omega}$ [rad]") +plt.ylabel("$i \sin{\Omega}$ [rad]") +plt.plot(np.rad2deg(inclinations)*np.cos(raan),np.rad2deg(inclinations)*np.sin(raan), \ + marker='o', markersize=5, linestyle='none') +plt.axis('equal') +axis.set_xlim(xmin=-180.0, xmax=180.0) +axis.set_ylim(ymin=-180.0, ymax=180.0) +figure.set_tight_layout(True) +plt.savefig(outputPath + "/figure4_debrisPopulation_inclinationVector.pdf", dpi = figureDPI) +plt.close() + +################################################################################################### +################################################################################################### \ No newline at end of file diff --git a/python/twoBodyMethods.py b/python/twoBodyMethods.py new file mode 100644 index 0000000..748532e --- /dev/null +++ b/python/twoBodyMethods.py @@ -0,0 +1,33 @@ +''' +Copyright (c) 2014, K. Kumar (me@kartikkumar.com) +All rights reserved. +''' + +################################################################################################### +# Set up modules and packages +################################################################################################### + +import numpy as np + +################################################################################################### + +################################################################################################### + +################################################################################################### +# Convert mean-motion to semi-major axis. +################################################################################################### + +def convertMeanMotionToSemiMajorAxis(meanMotion, gravitionalParameter): + return ( gravitionalParameter / ( meanMotion ** 2.0 ) ) ** ( 1.0 / 3.0 ) + +################################################################################################### + +################################################################################################### +# Convert mean-motion to semi-major axis. +################################################################################################### + +# Convert semi-major axis to mean-motion. +def convertSemiMajorAxisToMeanMotion(semiMajorAxis, gravitionalParameter): + return np.sqrt( gravitionalParameter / ( semiMajorAxis ** 3.0 ) ) + +################################################################################################### \ No newline at end of file From cb2a7485bc9b33258638e32b295c25037c6b53bf Mon Sep 17 00:00:00 2001 From: Enne Hekma Date: Tue, 29 Mar 2016 17:45:49 +0200 Subject: [PATCH 3/8] Add possibility to plot several TLE sets on top of eachother --- python/plotcatalog_asr.py | 311 ++++++++++++++++++++++++++++++++++++++ 1 file changed, 311 insertions(+) create mode 100644 python/plotcatalog_asr.py diff --git a/python/plotcatalog_asr.py b/python/plotcatalog_asr.py new file mode 100644 index 0000000..8304c76 --- /dev/null +++ b/python/plotcatalog_asr.py @@ -0,0 +1,311 @@ +''' +Copyright (c) 2014, K. Kumar (me@kartikkumar.com) +All rights reserved. +''' + +################################################################################################### +# Set up input deck +################################################################################################### + +# Set path to TLE catalog file. +tleCatalogFilePathall = "../data/catalog/3le.txt" +tleCatalogFilePathSSO = "../data/SSO/SSO_tle.txt" +tleCatalogFilePathGEO = "../data/GEO/GEO_tle.txt" +tleCatalogFilePathHEO = "../data/HEO/HEO_tle.txt" + +# Set number of lines per entry in TLE catalog (2 or 3). +tleEntryNumberOfLines = 3 + +# Set path to output directory. +outputPath = "../data/SSO/plots/" + +# Set figure DPI. +figureDPI = 300 + +# Set font size for axes labels. +fontSize = 22 + +################################################################################################### + +''' + DO NOT EDIT PARAMETERS BEYOND THIS POINT!!! +''' + +################################################################################################### +# Set up modules and packages +################################################################################################### + +from sgp4.earth_gravity import wgs72 +from sgp4.io import twoline2rv +from sgp4.propagation import getgravconst + +from matplotlib import rcParams +import matplotlib.pyplot as plt +import numpy as np +import pandas as pd +from twoBodyMethods import convertMeanMotionToSemiMajorAxis + +################################################################################################### + +################################################################################################### +# Read and store TLE catalog +################################################################################################### + + +# years = [2013, 2014, 2015] +# dn = [] +# for year in years: +# df1 = pd.DataFrame({'Incidents': [ 'C', 'B','A'], +# year: [1, 1, 1 ], +# }).set_index('Incidents') +# dn.append(df1) +# dn = pd.concat(dn, axis=1) +# print(dn) + +incl5 = [] +raan3 = [] +ecc3 = [] +ecc5 = [] +eccentricity = [] +sma = [] +inclinations3 = [] +aop3 = [] +order = ['all','SSO','GEO','HEO'] +markers = ['.','s','+','D'] +colors = ['k','b','g','r'] +for x in order: + # print eval(str('tleCatalogFilePath' + order[x])) + # print x + # time.sleep(10) + tleCatalogFilePathNew = eval(str('tleCatalogFilePath' + str(x))) + # print tleCatalogFilePathNew + # Read in catalog and store lines in list. + fileHandle = open(tleCatalogFilePathNew) + catalogLines = fileHandle.readlines() + fileHandle.close() + + # Strip newline and return carriage characters. + for i in xrange(len(catalogLines)): + catalogLines[i] = catalogLines[i].strip('\r\n') + + # Parse TLE entries and store debris objects. + debrisObjects = [] + for tleEntry in xrange(0,len(catalogLines),tleEntryNumberOfLines): + debrisObjects.append(twoline2rv(catalogLines[tleEntry+1], catalogLines[tleEntry+2], wgs72)) + + # Sort list of debris objects based on inclination. + inclinationSortedObjects = sorted(debrisObjects, key=lambda x: x.inclo, reverse=False) + + inclinations = [] + raan = [] + ecc = [] + aop = [] + for i in xrange(len(inclinationSortedObjects)): + inclinations.append(inclinationSortedObjects[i].inclo) + raan.append(inclinationSortedObjects[i].nodeo) + ecc.append(inclinationSortedObjects[i].ecco) + aop.append(inclinationSortedObjects[i].argpo) + + smatemp = [] + smatemp = [convertMeanMotionToSemiMajorAxis(debrisObject.no/60.0, getgravconst('wgs72')[1])-6373 \ + for debrisObject in debrisObjects] + smatemp2 = [] + smatemp2 = pd.DataFrame(smatemp, columns=[str(x)]) + sma.append(smatemp2) + + ecctemp = [] + ecctemp = [debrisObject.ecco for debrisObject in debrisObjects] + ecc4 = [] + ecc4 = pd.DataFrame(ecctemp, columns=[str(x)]) + ecc5.append(ecc4) + + incltemp = [] + incltemp = [np.rad2deg(debrisObject.inclo) for debrisObject in debrisObjects] + incl4 = [] + incl4 = pd.DataFrame(incltemp, columns=[str(x)]) + incl5.append(incl4) + + raan2 = [] + raan2 = pd.DataFrame(raan, columns=[str(x)]) + raan3.append(raan2) + + ecc2 = [] + ecc2 = pd.DataFrame(ecc, columns=[str(x)]) + ecc3.append(ecc2) + + aop2 = [] + aop2 = pd.DataFrame(aop, columns=[str(x)]) + aop3.append(aop2) + + inclinations2 = [] + inclinations2 = pd.DataFrame(inclinations, columns=[str(x)]) + inclinations3.append(inclinations2) +# print 'exit' +sma = pd.concat(sma, axis=1) +ecc5 = pd.concat(ecc5, axis=1) +incl5 = pd.concat(incl5, axis=1) + +raan3 = pd.concat(raan3, axis=1) +ecc3 = pd.concat(ecc3, axis=1) +aop3 = pd.concat(aop3, axis=1) +inclinations3 = pd.concat(inclinations3, axis=1) +# print raan3['GEO'] + + # print str('inclinations' + order[x] + 'new' ) + # str('inclinations' + order[x] + 'new' ) = [] + # str('inclinations' + order[x] + 'new' ) = inclinations + # str('raan' + order[x] + 'new' ) = [] + # str('raan' + order[x] + 'new' ) = raan + # str('ecc' + order[x] + 'new' ) = [] + # str('ecc' + order[x] + 'new' ) = ecc + # str('aop' + order[x] + 'new' ) = [] + # str('aop' + order[x] + 'new' ) = aop + # print raan[1] +# print raanSSOnew[1] +################################################################################################### + +################################################################################################### +# Generate plots +################################################################################################### + +# Set font size for plot labels. +rcParams.update({'font.size': fontSize}) +markers = ['.','s','+','D'] +colors = ['k','b','g','r'] + +# Plot distribution of eccentricity [-] against semi-major axis [km]. +figure = plt.figure() +axis = figure.add_subplot(111) +plt.xlabel("Semi-major axis altitude [km]") +plt.ylabel("Eccentricity [-]") +plt.ticklabel_format(style='sci', axis='x', scilimits=(0,0)) +# print +# print inclinationSortedObjects[1].no +plt.plot([convertMeanMotionToSemiMajorAxis(debrisObject.no/60.0, getgravconst('wgs72')[1]) \ + for debrisObject in debrisObjects], \ + [debrisObject.ecco for debrisObject in debrisObjects], \ + marker='o', markersize=1, color='k', linestyle='none') +# axis.set_xlim(xmax=0.8e4) +figure.set_tight_layout(True) +plt.savefig(outputPath + "/figure1_debrisPopulation_eccentricityVsSemiMajorAxis.pdf", \ + dpi = figureDPI) +plt.close() + +figure = plt.figure() +axis = figure.add_subplot(111) +plt.xlabel("Semi-major axis altitude [km]") +plt.ylabel("Eccentricity [-]") +plt.ticklabel_format(style='sci', axis='x', scilimits=(0,0)) +plt.plot(sma['all'],ecc3['all'], \ + marker='.', markersize=1, color='k', linestyle='none') +plt.plot(sma['SSO'],ecc3['SSO'], \ + marker='s', markersize=10, color='c', linestyle='none') +plt.plot(sma['GEO'],ecc3['GEO'], \ + marker='^', markersize=10, color='g', linestyle='none') +plt.plot(sma['HEO'],ecc3['HEO'], \ + marker='D', markersize=6, color='r', linestyle='none') +axis.set_xlim(xmax=0.5e5) +figure.set_tight_layout(True) +plt.savefig(outputPath + "/figure1_debrisPopulation_eccentricityVsSemiMajorAxisNew.pdf", \ + dpi = figureDPI) +plt.close() + +# Plot components of eccentricity vector [-]. +figure = plt.figure() +axis = figure.add_subplot(111) +plt.xlabel("$e \cos{\omega}$ [-]") +plt.ylabel("$e \sin{\omega}$ [-]") +plt.plot(ecc*np.cos(aop),ecc*np.sin(aop), marker='o', markersize=1, color='k', linestyle='none') +plt.axis('equal') +axis.set_xlim(xmin=-.21, xmax=.21) +axis.set_ylim(ymin=-.21, ymax=.21) +axis.set(xticks=[-.2,-.1,0,.1,.2])#, xticklabels=datelabels) #Same as plt.xticks +figure.set_tight_layout(True) +plt.savefig(outputPath + "/figure2_debrisPopulation_eccentricityVector.pdf", dpi = figureDPI) +# plt.savefig(outputPath + "/figure2_debrisPopulation_eccentricityVector.pdf", dpi = figureDPI) +plt.close() + +figure = plt.figure() +axis = figure.add_subplot(111) +plt.xlabel("$e \cos{\omega}$ [-]") +plt.ylabel("$e \sin{\omega}$ [-]") +# plt.plot(ecc*np.cos(aop),ecc*np.sin(aop), marker='o', markersize=1, color='k', linestyle='none') +plt.plot(ecc3['all']*np.cos(aop3['all']),ecc3['all']*np.sin(aop3['all']), marker='.', markersize=1, color='k', linestyle='none') +plt.plot(ecc3['SSO']*np.cos(aop3['SSO']),ecc3['SSO']*np.sin(aop3['SSO']), marker='s', markersize=10, color='c', linestyle='none') +plt.plot(ecc3['GEO']*np.cos(aop3['GEO']),ecc3['GEO']*np.sin(aop3['GEO']), marker='^', markersize=10, color='g', linestyle='none') +plt.plot(ecc3['HEO']*np.cos(aop3['HEO']),ecc3['HEO']*np.sin(aop3['HEO']), marker='D', markersize=6, color='r', linestyle='none') +plt.axis('equal') +# axis.set_xlim(xmin=-.21, xmax=.21) +# axis.set_ylim(ymin=-.21, ymax=.21) +# axis.set(xticks=[-.2,-.1,0,.1,.2])#, xticklabels=datelabels) #Same as plt.xticks +figure.set_tight_layout(True) +plt.savefig(outputPath + "/figure2_debrisPopulation_eccentricityVectorNew.pdf", dpi = figureDPI) +# plt.savefig(outputPath + "/figure2_debrisPopulation_eccentricityVector.pdf", dpi = figureDPI) +plt.close() + + +# Plot distribution of inclination [deg] against semi-major axis [km]. +figure = plt.figure() +axis = figure.add_subplot(111) +plt.xlabel("Semi-major axis altitude [km]") +plt.ylabel("Inclination [deg]") +plt.ticklabel_format(style='sci', axis='x', scilimits=(0,0)) +plt.plot(sma['all'],incl5['all'], \ + marker='.', markersize=1, color='k', linestyle='none') +plt.plot(sma['SSO'],incl5['SSO'], \ + marker='s', markersize=10, color='c', linestyle='none') +plt.plot(sma['GEO'],incl5['GEO'], \ + marker='^', markersize=10, color='g', linestyle='none') +plt.plot(sma['HEO'],incl5['HEO'], \ + marker='D', markersize=6, color='r', linestyle='none') +# plt.plot([convertMeanMotionToSemiMajorAxis(debrisObject.no/60.0, getgravconst('wgs72')[1]) \ +# for debrisObject in debrisObjects], \ +# [np.rad2deg(debrisObject.inclo) for debrisObject in debrisObjects], \ +# marker='o', markersize=1, color='k', linestyle='none') +axis.set_xlim(xmax=0.5e5) +figure.set_tight_layout(True) +plt.savefig(outputPath + "/figure3_debrisPopulation_inclinationVsSemiMajorAxisNew.pdf", \ + dpi = figureDPI) +plt.close() + +# Plot components of inclination vector [deg]. +figure = plt.figure() +axis = figure.add_subplot(111) +plt.xlabel("$i \cos{\Omega}$ [deg]") +plt.ylabel("$i \sin{\Omega}$ [deg]") +plt.plot(np.rad2deg(inclinations)*np.cos(raan),np.rad2deg(inclinations)*np.sin(raan), \ + marker='o', markersize=1, color='k', linestyle='none') +plt.axis('equal') +axis.set_xlim(xmin=-180.0, xmax=180.0) +axis.set_ylim(ymin=-180.0, ymax=180.0) +figure.set_tight_layout(True) +plt.savefig(outputPath + "/figure4_debrisPopulation_inclinationVector.pdf", dpi = figureDPI) +plt.close() + +# Plot components of inclination vector [deg]. +figure = plt.figure() +axis = figure.add_subplot(111) +plt.xlabel("$i \cos{\Omega}$ [deg]") +plt.ylabel("$i \sin{\Omega}$ [deg]") +# plt.plot(np.rad2deg(inclinations)*np.cos(raan),np.rad2deg(inclinations)*np.sin(raan), \ +# marker='o', markersize=1, color='k', linestyle='none') +plt.plot(np.rad2deg(inclinations3['all'])*np.cos(raan3['all']),np.rad2deg(inclinations3['all'])*np.sin(raan3['all']), \ + marker='.', markersize=1, color='k', linestyle='none') +plt.plot(np.rad2deg(inclinations3['SSO'])*np.cos(raan3['SSO']),np.rad2deg(inclinations3['SSO'])*np.sin(raan3['SSO']), \ + marker='s', markersize=10, color='c', linestyle='none') +plt.plot(np.rad2deg(inclinations3['GEO'])*np.cos(raan3['GEO']),np.rad2deg(inclinations3['GEO'])*np.sin(raan3['GEO']), \ + marker='^', markersize=10, color='g', linestyle='none') +plt.plot(np.rad2deg(inclinations3['HEO'])*np.cos(raan3['HEO']),np.rad2deg(inclinations3['HEO'])*np.sin(raan3['HEO']), \ + marker='D', markersize=6, color='r', linestyle='none') +plt.axis('equal') +axis.set_xlim(xmin=-180.0, xmax=180.0) +axis.set_ylim(ymin=-180.0, ymax=180.0) +figure.set_tight_layout(True) +plt.savefig(outputPath + "/figure4_debrisPopulation_inclinationVectorNew.pdf", dpi = figureDPI) +plt.close() + + + + +################################################################################################### +################################################################################################### \ No newline at end of file From 49483c4f05e92852647f8127669fd0495325b2b8 Mon Sep 17 00:00:00 2001 From: Enne Hekma Date: Tue, 29 Mar 2016 17:49:24 +0200 Subject: [PATCH 4/8] Removed unnecessary plot script --- python/plotcatalog.py | 149 ------------------------------------------ 1 file changed, 149 deletions(-) delete mode 100644 python/plotcatalog.py diff --git a/python/plotcatalog.py b/python/plotcatalog.py deleted file mode 100644 index 81fc954..0000000 --- a/python/plotcatalog.py +++ /dev/null @@ -1,149 +0,0 @@ -''' -Copyright (c) 2014, K. Kumar (me@kartikkumar.com) -All rights reserved. -''' - -################################################################################################### -# Set up input deck -################################################################################################### - -# Set path to TLE catalog file. -tleCatalogFilePath = "../data/catalogs/catalog_rocketbodiesLEOfull.txt" - -# Set number of lines per entry in TLE catalog (2 or 3). -tleEntryNumberOfLines = 3 - -# Set path to output directory. -outputPath = "../data/napacatalog/" - -# Set figure DPI. -figureDPI = 300 - -# Set font size for axes labels. -fontSize = 24 - -################################################################################################### - -''' - DO NOT EDIT PARAMETERS BEYOND THIS POINT!!! -''' - -################################################################################################### -# Set up modules and packages -################################################################################################### - -from sgp4.earth_gravity import wgs72 -from sgp4.io import twoline2rv -from sgp4.propagation import getgravconst - -from matplotlib import rcParams -import matplotlib.pyplot as plt -import numpy as np - -from twoBodyMethods import convertMeanMotionToSemiMajorAxis - -################################################################################################### - -################################################################################################### -# Read and store TLE catalog -################################################################################################### - -# Read in catalog and store lines in list. -fileHandle = open(tleCatalogFilePath) -catalogLines = fileHandle.readlines() -fileHandle.close() - -# Strip newline and return carriage characters. -for i in xrange(len(catalogLines)): - catalogLines[i] = catalogLines[i].strip('\r\n') - -# Parse TLE entries and store debris objects. -debrisObjects = [] -for tleEntry in xrange(0,len(catalogLines),tleEntryNumberOfLines): - debrisObjects.append(twoline2rv(catalogLines[tleEntry+1], catalogLines[tleEntry+2], wgs72)) - -# Sort list of debris objects based on inclination. -inclinationSortedObjects = sorted(debrisObjects, key=lambda x: x.inclo, reverse=False) - -inclinations = [] -raan = [] -ecc = [] -aop = [] -for i in xrange(len(inclinationSortedObjects)-1): - inclinations.append(inclinationSortedObjects[i].inclo) - raan.append(inclinationSortedObjects[i].nodeo) - ecc.append(inclinationSortedObjects[i].ecco) - aop.append(inclinationSortedObjects[i].argpo) - -################################################################################################### - -################################################################################################### -# Generate plots -################################################################################################### - -# Set font size for plot labels. -rcParams.update({'font.size': fontSize}) - -# Plot distribution of eccentricity [-] against semi-major axis [km]. -figure = plt.figure() -axis = figure.add_subplot(111) -plt.xlabel("Semi-major axis [km]") -plt.ylabel("Eccentricity [-]") -plt.ticklabel_format(style='sci', axis='x', scilimits=(0,0)) -plt.plot([convertMeanMotionToSemiMajorAxis(debrisObject.no/60.0, getgravconst('wgs72')[1]) \ - for debrisObject in debrisObjects], \ - [debrisObject.ecco for debrisObject in debrisObjects], \ - marker='o', markersize=5, linestyle='none') -axis.set_xlim(xmax=0.8e4) -figure.set_tight_layout(True) -plt.savefig(outputPath + "/figure1_debrisPopulation_eccentricityVsSemiMajorAxis.pdf", \ - dpi = figureDPI) -plt.close() - -# Plot components of eccentricity vector [-]. -figure = plt.figure() -axis = figure.add_subplot(111) -plt.xlabel("$e \cos{\omega}$ [-]") -plt.ylabel("$e \sin{\omega}$ [-]") -plt.plot(ecc*np.cos(aop),ecc*np.sin(aop), marker='o', markersize=5, linestyle='none') -plt.axis('equal') -axis.set_xlim(xmin=-.21, xmax=.21) -axis.set_ylim(ymin=-.21, ymax=.21) -axis.set(xticks=[-.2,-.1,0,.1,.2])#, xticklabels=datelabels) #Same as plt.xticks -figure.set_tight_layout(True) -plt.savefig(outputPath + "/figure2_debrisPopulation_eccentricityVector.pdf", dpi = figureDPI) -plt.savefig(outputPath + "/figure2_debrisPopulation_eccentricityVector.pdf", dpi = figureDPI) -plt.close() - -# Plot distribution of inclination [deg] against semi-major axis [km]. -figure = plt.figure() -axis = figure.add_subplot(111) -plt.xlabel("Semi-major axis [km]") -plt.ylabel("Inclination [deg]") -plt.ticklabel_format(style='sci', axis='x', scilimits=(0,0)) -plt.plot([convertMeanMotionToSemiMajorAxis(debrisObject.no/60.0, getgravconst('wgs72')[1]) \ - for debrisObject in debrisObjects], \ - [np.rad2deg(debrisObject.inclo) for debrisObject in debrisObjects], \ - marker='o', markersize=5, linestyle='none') -axis.set_xlim(xmax=0.8e4) -figure.set_tight_layout(True) -plt.savefig(outputPath + "/figure3_debrisPopulation_inclinationVsSemiMajorAxis.pdf", \ - dpi = figureDPI) -plt.close() - -# Plot components of inclination vector [deg]. -figure = plt.figure() -axis = figure.add_subplot(111) -plt.xlabel("$i \cos{\Omega}$ [rad]") -plt.ylabel("$i \sin{\Omega}$ [rad]") -plt.plot(np.rad2deg(inclinations)*np.cos(raan),np.rad2deg(inclinations)*np.sin(raan), \ - marker='o', markersize=5, linestyle='none') -plt.axis('equal') -axis.set_xlim(xmin=-180.0, xmax=180.0) -axis.set_ylim(ymin=-180.0, ymax=180.0) -figure.set_tight_layout(True) -plt.savefig(outputPath + "/figure4_debrisPopulation_inclinationVector.pdf", dpi = figureDPI) -plt.close() - -################################################################################################### -################################################################################################### \ No newline at end of file From 03b98b9bc87ee8058790c5d00cb2f4b521e8b9f2 Mon Sep 17 00:00:00 2001 From: Enne Hekma Date: Tue, 29 Mar 2016 18:05:39 +0200 Subject: [PATCH 5/8] Clean up code and break at 100 columns --- python/plotcatalog_asr.py | 192 ++++++++++++-------------------------- 1 file changed, 62 insertions(+), 130 deletions(-) diff --git a/python/plotcatalog_asr.py b/python/plotcatalog_asr.py index 8304c76..2050840 100644 --- a/python/plotcatalog_asr.py +++ b/python/plotcatalog_asr.py @@ -1,5 +1,6 @@ ''' -Copyright (c) 2014, K. Kumar (me@kartikkumar.com) +Copyright (c) 2014-2016, K. Kumar (me@kartikkumar.com) +Copyright (c) 2016, E.J. Hekma (ennehekma@gmail.com) All rights reserved. ''' @@ -51,34 +52,22 @@ # Read and store TLE catalog ################################################################################################### - -# years = [2013, 2014, 2015] -# dn = [] -# for year in years: -# df1 = pd.DataFrame({'Incidents': [ 'C', 'B','A'], -# year: [1, 1, 1 ], -# }).set_index('Incidents') -# dn.append(df1) -# dn = pd.concat(dn, axis=1) -# print(dn) - +sma = [] +ecc5 = [] incl5 = [] + raan3 = [] ecc3 = [] -ecc5 = [] -eccentricity = [] -sma = [] -inclinations3 = [] aop3 = [] +inclinations3 = [] + +# colors = ['k','b','g','r'] +# markers = ['.','s','+','D'] + order = ['all','SSO','GEO','HEO'] -markers = ['.','s','+','D'] -colors = ['k','b','g','r'] for x in order: - # print eval(str('tleCatalogFilePath' + order[x])) - # print x - # time.sleep(10) tleCatalogFilePathNew = eval(str('tleCatalogFilePath' + str(x))) - # print tleCatalogFilePathNew + # Read in catalog and store lines in list. fileHandle = open(tleCatalogFilePathNew) catalogLines = fileHandle.readlines() @@ -107,8 +96,9 @@ aop.append(inclinationSortedObjects[i].argpo) smatemp = [] - smatemp = [convertMeanMotionToSemiMajorAxis(debrisObject.no/60.0, getgravconst('wgs72')[1])-6373 \ - for debrisObject in debrisObjects] + smatemp = [convertMeanMotionToSemiMajorAxis(debrisObject.no/60.0, \ + getgravconst('wgs72')[1])-6373 \ + for debrisObject in debrisObjects] smatemp2 = [] smatemp2 = pd.DataFrame(smatemp, columns=[str(x)]) sma.append(smatemp2) @@ -140,7 +130,7 @@ inclinations2 = [] inclinations2 = pd.DataFrame(inclinations, columns=[str(x)]) inclinations3.append(inclinations2) -# print 'exit' + sma = pd.concat(sma, axis=1) ecc5 = pd.concat(ecc5, axis=1) incl5 = pd.concat(incl5, axis=1) @@ -149,19 +139,7 @@ ecc3 = pd.concat(ecc3, axis=1) aop3 = pd.concat(aop3, axis=1) inclinations3 = pd.concat(inclinations3, axis=1) -# print raan3['GEO'] - - # print str('inclinations' + order[x] + 'new' ) - # str('inclinations' + order[x] + 'new' ) = [] - # str('inclinations' + order[x] + 'new' ) = inclinations - # str('raan' + order[x] + 'new' ) = [] - # str('raan' + order[x] + 'new' ) = raan - # str('ecc' + order[x] + 'new' ) = [] - # str('ecc' + order[x] + 'new' ) = ecc - # str('aop' + order[x] + 'new' ) = [] - # str('aop' + order[x] + 'new' ) = aop - # print raan[1] -# print raanSSOnew[1] + ################################################################################################### ################################################################################################### @@ -170,8 +148,6 @@ # Set font size for plot labels. rcParams.update({'font.size': fontSize}) -markers = ['.','s','+','D'] -colors = ['k','b','g','r'] # Plot distribution of eccentricity [-] against semi-major axis [km]. figure = plt.figure() @@ -179,34 +155,17 @@ plt.xlabel("Semi-major axis altitude [km]") plt.ylabel("Eccentricity [-]") plt.ticklabel_format(style='sci', axis='x', scilimits=(0,0)) -# print -# print inclinationSortedObjects[1].no -plt.plot([convertMeanMotionToSemiMajorAxis(debrisObject.no/60.0, getgravconst('wgs72')[1]) \ - for debrisObject in debrisObjects], \ - [debrisObject.ecco for debrisObject in debrisObjects], \ - marker='o', markersize=1, color='k', linestyle='none') -# axis.set_xlim(xmax=0.8e4) -figure.set_tight_layout(True) -plt.savefig(outputPath + "/figure1_debrisPopulation_eccentricityVsSemiMajorAxis.pdf", \ - dpi = figureDPI) -plt.close() - -figure = plt.figure() -axis = figure.add_subplot(111) -plt.xlabel("Semi-major axis altitude [km]") -plt.ylabel("Eccentricity [-]") -plt.ticklabel_format(style='sci', axis='x', scilimits=(0,0)) -plt.plot(sma['all'],ecc3['all'], \ - marker='.', markersize=1, color='k', linestyle='none') -plt.plot(sma['SSO'],ecc3['SSO'], \ - marker='s', markersize=10, color='c', linestyle='none') -plt.plot(sma['GEO'],ecc3['GEO'], \ - marker='^', markersize=10, color='g', linestyle='none') -plt.plot(sma['HEO'],ecc3['HEO'], \ - marker='D', markersize=6, color='r', linestyle='none') -axis.set_xlim(xmax=0.5e5) +plt.plot(sma['all'],ecc3['all'], \ + marker='.', markersize=1, color='k', linestyle='none') +plt.plot(sma['SSO'],ecc3['SSO'], \ + marker='s', markersize=10, color='c', linestyle='none') +plt.plot(sma['GEO'],ecc3['GEO'], \ + marker='^', markersize=10, color='g', linestyle='none') +plt.plot(sma['HEO'],ecc3['HEO'], \ + marker='D', markersize=6, color='r', linestyle='none') +axis.set_xlim(xmax=0.4e5) figure.set_tight_layout(True) -plt.savefig(outputPath + "/figure1_debrisPopulation_eccentricityVsSemiMajorAxisNew.pdf", \ +plt.savefig(outputPath + "/figure1_debrisPopulation_eccentricityVsSemiMajorAxis.pdf", \ dpi = figureDPI) plt.close() @@ -215,32 +174,21 @@ axis = figure.add_subplot(111) plt.xlabel("$e \cos{\omega}$ [-]") plt.ylabel("$e \sin{\omega}$ [-]") -plt.plot(ecc*np.cos(aop),ecc*np.sin(aop), marker='o', markersize=1, color='k', linestyle='none') +plt.plot(ecc3['all']*np.cos(aop3['all']),ecc3['all']*np.sin(aop3['all']), \ + marker='.', markersize=1, color='k', linestyle='none') +plt.plot(ecc3['SSO']*np.cos(aop3['SSO']),ecc3['SSO']*np.sin(aop3['SSO']), \ + marker='s', markersize=10, color='c', linestyle='none') +plt.plot(ecc3['GEO']*np.cos(aop3['GEO']),ecc3['GEO']*np.sin(aop3['GEO']), \ + marker='^', markersize=10, color='g', linestyle='none') +plt.plot(ecc3['HEO']*np.cos(aop3['HEO']),ecc3['HEO']*np.sin(aop3['HEO']), \ + marker='D', markersize=6, color='r', linestyle='none') plt.axis('equal') -axis.set_xlim(xmin=-.21, xmax=.21) -axis.set_ylim(ymin=-.21, ymax=.21) -axis.set(xticks=[-.2,-.1,0,.1,.2])#, xticklabels=datelabels) #Same as plt.xticks +axis.set_xlim(xmin=-.82, xmax=.82) +axis.set_ylim(ymin=-.82, ymax=.82) +axis.set(xticks=[-.8,-.4,0,.4,.8]) +axis.set(yticks=[-.8,-.4,0,.4,.8]) figure.set_tight_layout(True) plt.savefig(outputPath + "/figure2_debrisPopulation_eccentricityVector.pdf", dpi = figureDPI) -# plt.savefig(outputPath + "/figure2_debrisPopulation_eccentricityVector.pdf", dpi = figureDPI) -plt.close() - -figure = plt.figure() -axis = figure.add_subplot(111) -plt.xlabel("$e \cos{\omega}$ [-]") -plt.ylabel("$e \sin{\omega}$ [-]") -# plt.plot(ecc*np.cos(aop),ecc*np.sin(aop), marker='o', markersize=1, color='k', linestyle='none') -plt.plot(ecc3['all']*np.cos(aop3['all']),ecc3['all']*np.sin(aop3['all']), marker='.', markersize=1, color='k', linestyle='none') -plt.plot(ecc3['SSO']*np.cos(aop3['SSO']),ecc3['SSO']*np.sin(aop3['SSO']), marker='s', markersize=10, color='c', linestyle='none') -plt.plot(ecc3['GEO']*np.cos(aop3['GEO']),ecc3['GEO']*np.sin(aop3['GEO']), marker='^', markersize=10, color='g', linestyle='none') -plt.plot(ecc3['HEO']*np.cos(aop3['HEO']),ecc3['HEO']*np.sin(aop3['HEO']), marker='D', markersize=6, color='r', linestyle='none') -plt.axis('equal') -# axis.set_xlim(xmin=-.21, xmax=.21) -# axis.set_ylim(ymin=-.21, ymax=.21) -# axis.set(xticks=[-.2,-.1,0,.1,.2])#, xticklabels=datelabels) #Same as plt.xticks -figure.set_tight_layout(True) -plt.savefig(outputPath + "/figure2_debrisPopulation_eccentricityVectorNew.pdf", dpi = figureDPI) -# plt.savefig(outputPath + "/figure2_debrisPopulation_eccentricityVector.pdf", dpi = figureDPI) plt.close() @@ -250,21 +198,17 @@ plt.xlabel("Semi-major axis altitude [km]") plt.ylabel("Inclination [deg]") plt.ticklabel_format(style='sci', axis='x', scilimits=(0,0)) -plt.plot(sma['all'],incl5['all'], \ - marker='.', markersize=1, color='k', linestyle='none') -plt.plot(sma['SSO'],incl5['SSO'], \ - marker='s', markersize=10, color='c', linestyle='none') -plt.plot(sma['GEO'],incl5['GEO'], \ - marker='^', markersize=10, color='g', linestyle='none') -plt.plot(sma['HEO'],incl5['HEO'], \ - marker='D', markersize=6, color='r', linestyle='none') -# plt.plot([convertMeanMotionToSemiMajorAxis(debrisObject.no/60.0, getgravconst('wgs72')[1]) \ -# for debrisObject in debrisObjects], \ -# [np.rad2deg(debrisObject.inclo) for debrisObject in debrisObjects], \ -# marker='o', markersize=1, color='k', linestyle='none') -axis.set_xlim(xmax=0.5e5) +plt.plot(sma['all'],incl5['all'], \ + marker='.', markersize=1, color='k', linestyle='none') +plt.plot(sma['SSO'],incl5['SSO'], \ + marker='s', markersize=10, color='c', linestyle='none') +plt.plot(sma['GEO'],incl5['GEO'], \ + marker='^', markersize=10, color='g', linestyle='none') +plt.plot(sma['HEO'],incl5['HEO'], \ + marker='D', markersize=6, color='r', linestyle='none') +axis.set_xlim(xmax=0.4e5) figure.set_tight_layout(True) -plt.savefig(outputPath + "/figure3_debrisPopulation_inclinationVsSemiMajorAxisNew.pdf", \ +plt.savefig(outputPath + "/figure3_debrisPopulation_inclinationVsSemiMajorAxis.pdf", \ dpi = figureDPI) plt.close() @@ -273,37 +217,25 @@ axis = figure.add_subplot(111) plt.xlabel("$i \cos{\Omega}$ [deg]") plt.ylabel("$i \sin{\Omega}$ [deg]") -plt.plot(np.rad2deg(inclinations)*np.cos(raan),np.rad2deg(inclinations)*np.sin(raan), \ - marker='o', markersize=1, color='k', linestyle='none') +plt.plot(np.rad2deg(inclinations3['all'])*np.cos(raan3['all']), \ + np.rad2deg(inclinations3['all'])*np.sin(raan3['all']), \ + marker='.', markersize=1, color='k', linestyle='none') +plt.plot(np.rad2deg(inclinations3['SSO'])*np.cos(raan3['SSO']), \ + np.rad2deg(inclinations3['SSO'])*np.sin(raan3['SSO']), \ + marker='s', markersize=10, color='c', linestyle='none') +plt.plot(np.rad2deg(inclinations3['GEO'])*np.cos(raan3['GEO']), \ + np.rad2deg(inclinations3['GEO'])*np.sin(raan3['GEO']), \ + marker='^', markersize=10, color='g', linestyle='none') +plt.plot(np.rad2deg(inclinations3['HEO'])*np.cos(raan3['HEO']), \ + np.rad2deg(inclinations3['HEO'])*np.sin(raan3['HEO']), \ + marker='D', markersize=6, color='r', linestyle='none') plt.axis('equal') -axis.set_xlim(xmin=-180.0, xmax=180.0) -axis.set_ylim(ymin=-180.0, ymax=180.0) +axis.set_xlim(xmin=-110.0, xmax=110.0) +axis.set_ylim(ymin=-110.0, ymax=110.0) figure.set_tight_layout(True) plt.savefig(outputPath + "/figure4_debrisPopulation_inclinationVector.pdf", dpi = figureDPI) plt.close() -# Plot components of inclination vector [deg]. -figure = plt.figure() -axis = figure.add_subplot(111) -plt.xlabel("$i \cos{\Omega}$ [deg]") -plt.ylabel("$i \sin{\Omega}$ [deg]") -# plt.plot(np.rad2deg(inclinations)*np.cos(raan),np.rad2deg(inclinations)*np.sin(raan), \ -# marker='o', markersize=1, color='k', linestyle='none') -plt.plot(np.rad2deg(inclinations3['all'])*np.cos(raan3['all']),np.rad2deg(inclinations3['all'])*np.sin(raan3['all']), \ - marker='.', markersize=1, color='k', linestyle='none') -plt.plot(np.rad2deg(inclinations3['SSO'])*np.cos(raan3['SSO']),np.rad2deg(inclinations3['SSO'])*np.sin(raan3['SSO']), \ - marker='s', markersize=10, color='c', linestyle='none') -plt.plot(np.rad2deg(inclinations3['GEO'])*np.cos(raan3['GEO']),np.rad2deg(inclinations3['GEO'])*np.sin(raan3['GEO']), \ - marker='^', markersize=10, color='g', linestyle='none') -plt.plot(np.rad2deg(inclinations3['HEO'])*np.cos(raan3['HEO']),np.rad2deg(inclinations3['HEO'])*np.sin(raan3['HEO']), \ - marker='D', markersize=6, color='r', linestyle='none') -plt.axis('equal') -axis.set_xlim(xmin=-180.0, xmax=180.0) -axis.set_ylim(ymin=-180.0, ymax=180.0) -figure.set_tight_layout(True) -plt.savefig(outputPath + "/figure4_debrisPopulation_inclinationVectorNew.pdf", dpi = figureDPI) -plt.close() - From 44dcdd753a4163ead9dd504c4524d0ded60cc5cf Mon Sep 17 00:00:00 2001 From: Enne Hekma Date: Mon, 27 Jun 2016 11:26:44 +0200 Subject: [PATCH 6/8] Add plotting script for literature --- python/plotcatalog_literature.py | 244 +++++++++++++++++++++++++++++++ 1 file changed, 244 insertions(+) create mode 100644 python/plotcatalog_literature.py diff --git a/python/plotcatalog_literature.py b/python/plotcatalog_literature.py new file mode 100644 index 0000000..63b54df --- /dev/null +++ b/python/plotcatalog_literature.py @@ -0,0 +1,244 @@ +''' +Copyright (c) 2014-2016, K. Kumar (me@kartikkumar.com) +Copyright (c) 2016, E.J. Hekma (ennehekma@gmail.com) +All rights reserved. +''' + +################################################################################################### +# Set up input deck +################################################################################################### + +# Set path to TLE catalog file. +tleCatalogFilePathall = "../data/SSO/leo_all.txt" +tleCatalogFilePathSSO = "../data/SSO/ARIANE.txt" +tleCatalogFilePathGEO = "../data/SSO/SSO_tle.txt" +tleCatalogFilePathHEO = "../data/SSO/ARIANE_RB.txt" + +# Set number of lines per entry in TLE catalog (2 or 3)self. +tleEntryNumberOfLines = 3 + +# Set path to output directory. +outputPath = "../data/plots_literature/" + +# Set figure DPI. +figureDPI = 300 + +# Set font size for axes labels. +fontSize = 22 + +################################################################################################### + +''' + DO NOT EDIT PARAMETERS BEYOND THIS POINT!!! +''' + +################################################################################################### +# Set up modules and packages +################################################################################################### + +from sgp4.earth_gravity import wgs72 +from sgp4.io import twoline2rv +from sgp4.propagation import getgravconst + +from matplotlib import rcParams +import matplotlib.pyplot as plt +import numpy as np +import pandas as pd +from twoBodyMethods import convertMeanMotionToSemiMajorAxis + +################################################################################################### + +################################################################################################### +# Read and store TLE catalog +################################################################################################### + +sma = [] +ecc5 = [] +incl5 = [] + +raan3 = [] +ecc3 = [] +aop3 = [] +inclinations3 = [] + +# colors = ['k','b','g','r'] +# markers = ['.','s','+','D'] + +order = ['all','SSO','GEO','HEO'] +for x in order: + tleCatalogFilePathNew = eval(str('tleCatalogFilePath' + str(x))) + + # Read in catalog and store lines in list. + fileHandle = open(tleCatalogFilePathNew) + catalogLines = fileHandle.readlines() + fileHandle.close() + + # Strip newline and return carriage characters. + for i in xrange(len(catalogLines)): + catalogLines[i] = catalogLines[i].strip('\r\n') + + # Parse TLE entries and store debris objects. + debrisObjects = [] + for tleEntry in xrange(0,len(catalogLines),tleEntryNumberOfLines): + debrisObjects.append(twoline2rv(catalogLines[tleEntry+1], catalogLines[tleEntry+2], wgs72)) + + # Sort list of debris objects based on inclination. + inclinationSortedObjects = sorted(debrisObjects, key=lambda x: x.inclo, reverse=False) + + inclinations = [] + raan = [] + ecc = [] + aop = [] + for i in xrange(len(inclinationSortedObjects)): + inclinations.append(inclinationSortedObjects[i].inclo) + raan.append(inclinationSortedObjects[i].nodeo) + ecc.append(inclinationSortedObjects[i].ecco) + aop.append(inclinationSortedObjects[i].argpo) + + smatemp = [] + smatemp = [convertMeanMotionToSemiMajorAxis(debrisObject.no/60.0, \ + getgravconst('wgs72')[1])-6373 \ + for debrisObject in debrisObjects] + smatemp2 = [] + smatemp2 = pd.DataFrame(smatemp, columns=[str(x)]) + sma.append(smatemp2) + + ecctemp = [] + ecctemp = [debrisObject.ecco for debrisObject in debrisObjects] + ecc4 = [] + ecc4 = pd.DataFrame(ecctemp, columns=[str(x)]) + ecc5.append(ecc4) + + incltemp = [] + incltemp = [np.rad2deg(debrisObject.inclo) for debrisObject in debrisObjects] + incl4 = [] + incl4 = pd.DataFrame(incltemp, columns=[str(x)]) + incl5.append(incl4) + + raan2 = [] + raan2 = pd.DataFrame(raan, columns=[str(x)]) + raan3.append(raan2) + + ecc2 = [] + ecc2 = pd.DataFrame(ecc, columns=[str(x)]) + ecc3.append(ecc2) + + aop2 = [] + aop2 = pd.DataFrame(aop, columns=[str(x)]) + aop3.append(aop2) + + inclinations2 = [] + inclinations2 = pd.DataFrame(inclinations, columns=[str(x)]) + inclinations3.append(inclinations2) + +sma = pd.concat(sma, axis=1) +ecc5 = pd.concat(ecc5, axis=1) +incl5 = pd.concat(incl5, axis=1) + +raan3 = pd.concat(raan3, axis=1) +ecc3 = pd.concat(ecc3, axis=1) +aop3 = pd.concat(aop3, axis=1) +inclinations3 = pd.concat(inclinations3, axis=1) + +################################################################################################### + +################################################################################################### +# Generate plots +################################################################################################### + +# Set font size for plot labels. +rcParams.update({'font.size': fontSize}) + +# Plot distribution of eccentricity [-] against semi-major axis [km]. +figure = plt.figure() +axis = figure.add_subplot(111) +plt.xlabel("Semi-major axis altitude [km]") +plt.ylabel("Eccentricity [-]") +plt.ticklabel_format(style='sci', axis='x', scilimits=(0,0)) +plt.plot(sma['all'],ecc3['all'], \ + marker='.', markersize=1, color='k', linestyle='none') +# plt.plot(sma['SSO'],ecc3['SSO'], \ +# marker='s', markersize=10, color='c', linestyle='none') +# plt.plot(sma['GEO'],ecc3['GEO'], \ +# marker='^', markersize=10, color='g', linestyle='none') +plt.plot(sma['HEO'],ecc3['HEO'], \ + marker='D', markersize=6, color='r', linestyle='none') +axis.set_xlim(xmax=1500) +axis.set_ylim(ymax=0.04) +figure.set_tight_layout(True) +plt.savefig(outputPath + "/figure1_debrisPopulation_eccentricityVsSemiMajorAxis.pdf", \ + dpi = figureDPI) +plt.close() + +# Plot components of eccentricity vector [-]. +figure = plt.figure() +axis = figure.add_subplot(111) +plt.xlabel("$e \cos{\omega}$ [-]") +plt.ylabel("$e \sin{\omega}$ [-]") +plt.plot(ecc3['all']*np.cos(aop3['all']),ecc3['all']*np.sin(aop3['all']), \ + marker='.', markersize=1, color='k', linestyle='none') +# plt.plot(ecc3['SSO']*np.cos(aop3['SSO']),ecc3['SSO']*np.sin(aop3['SSO']), \ +# marker='s', markersize=10, color='c', linestyle='none') +# plt.plot(ecc3['GEO']*np.cos(aop3['GEO']),ecc3['GEO']*np.sin(aop3['GEO']), \ +# marker='^', markersize=10, color='g', linestyle='none') +plt.plot(ecc3['HEO']*np.cos(aop3['HEO']),ecc3['HEO']*np.sin(aop3['HEO']), \ + marker='D', markersize=6, color='r', linestyle='none') +plt.axis('equal') +axis.set_xlim(xmin=-.021, xmax=.021) +axis.set_ylim(ymin=-.021, ymax=.021) +axis.set(xticks=[-.02,-.01,0,.01,.02]) +axis.set(yticks=[-.02,-.01,0,.01,.02]) +figure.set_tight_layout(True) +plt.savefig(outputPath + "/figure2_debrisPopulation_eccentricityVector.pdf", dpi = figureDPI) +plt.close() + + +# Plot distribution of inclination [deg] against semi-major axis [km]. +figure = plt.figure() +axis = figure.add_subplot(111) +plt.xlabel("Semi-major axis altitude [km]") +plt.ylabel("Inclination [deg]") +plt.ticklabel_format(style='sci', axis='x', scilimits=(0,0)) +plt.plot(sma['all'],incl5['all'], \ + marker='.', markersize=1, color='k', linestyle='none') +# plt.plot(sma['SSO'],incl5['SSO'], \ +# marker='s', markersize=10, color='c', linestyle='none') +# plt.plot(sma['GEO'],incl5['GEO'], \ +# marker='^', markersize=10, color='g', linestyle='none') +plt.plot(sma['HEO'],incl5['HEO'], \ + marker='D', markersize=6, color='r', linestyle='none') +axis.set_xlim(xmax=1500) +figure.set_tight_layout(True) +plt.savefig(outputPath + "/figure3_debrisPopulation_inclinationVsSemiMajorAxis.pdf", \ + dpi = figureDPI) +plt.close() + +# Plot components of inclination vector [deg]. +figure = plt.figure() +axis = figure.add_subplot(111) +plt.xlabel("$i \cos{\Omega}$ [deg]") +plt.ylabel("$i \sin{\Omega}$ [deg]") +plt.plot(np.rad2deg(inclinations3['all'])*np.cos(raan3['all']), \ + np.rad2deg(inclinations3['all'])*np.sin(raan3['all']), \ + marker='.', markersize=1, color='k', linestyle='none') +# plt.plot(np.rad2deg(inclinations3['SSO'])*np.cos(raan3['SSO']), \ +# np.rad2deg(inclinations3['SSO'])*np.sin(raan3['SSO']), \ +# marker='s', markersize=10, color='c', linestyle='none') +# plt.plot(np.rad2deg(inclinations3['GEO'])*np.cos(raan3['GEO']), \ +# np.rad2deg(inclinations3['GEO'])*np.sin(raan3['GEO']), \ +# marker='^', markersize=10, color='g', linestyle='none') +plt.plot(np.rad2deg(inclinations3['HEO'])*np.cos(raan3['HEO']), \ + np.rad2deg(inclinations3['HEO'])*np.sin(raan3['HEO']), \ + marker='D', markersize=6, color='r', linestyle='none') +plt.axis('equal') +axis.set_xlim(xmin=-110.0, xmax=110.0) +axis.set_ylim(ymin=-110.0, ymax=110.0) +figure.set_tight_layout(True) +plt.savefig(outputPath + "/figure4_debrisPopulation_inclinationVector.pdf", dpi = figureDPI) +plt.close() + + + + +################################################################################################### +################################################################################################### From 99cada7e788e18e57c2b356a8bb09c2fa667aadd Mon Sep 17 00:00:00 2001 From: Enne Hekma Date: Thu, 14 Jul 2016 10:11:54 +0200 Subject: [PATCH 7/8] minor alterations --- python/plotcatalog_literature.py | 16 +++++++++------- 1 file changed, 9 insertions(+), 7 deletions(-) diff --git a/python/plotcatalog_literature.py b/python/plotcatalog_literature.py index 63b54df..cdb5fcb 100644 --- a/python/plotcatalog_literature.py +++ b/python/plotcatalog_literature.py @@ -12,7 +12,8 @@ tleCatalogFilePathall = "../data/SSO/leo_all.txt" tleCatalogFilePathSSO = "../data/SSO/ARIANE.txt" tleCatalogFilePathGEO = "../data/SSO/SSO_tle.txt" -tleCatalogFilePathHEO = "../data/SSO/ARIANE_RB.txt" +# tleCatalogFilePathHEO = "../data/SSO/ARIANE_RB.txt" +tleCatalogFilePathHEO = "../data/SSO/catalog_800_99_01.txt" # Set number of lines per entry in TLE catalog (2 or 3)self. tleEntryNumberOfLines = 3 @@ -198,16 +199,17 @@ axis = figure.add_subplot(111) plt.xlabel("Semi-major axis altitude [km]") plt.ylabel("Inclination [deg]") -plt.ticklabel_format(style='sci', axis='x', scilimits=(0,0)) -plt.plot(sma['all'],incl5['all'], \ - marker='.', markersize=1, color='k', linestyle='none') -# plt.plot(sma['SSO'],incl5['SSO'], \ +# plt.ticklabel_format(style='sci', axis='x', scilimits=(0,0)) +# plt.plot(sma['all'],incl5['all'], \ +# marker='.', markersize=1, color='k', linestyle='none') +# # plt.plot(sma['SSO'],incl5['SSO'], \ # marker='s', markersize=10, color='c', linestyle='none') # plt.plot(sma['GEO'],incl5['GEO'], \ # marker='^', markersize=10, color='g', linestyle='none') plt.plot(sma['HEO'],incl5['HEO'], \ - marker='D', markersize=6, color='r', linestyle='none') -axis.set_xlim(xmax=1500) + marker='D', markersize=2, color='r', linestyle='none') +axis.set_xlim(xmax=900,xmin=700) +axis.set_ylim(ymin = 98, ymax=100) figure.set_tight_layout(True) plt.savefig(outputPath + "/figure3_debrisPopulation_inclinationVsSemiMajorAxis.pdf", \ dpi = figureDPI) From 4e07975ed88c0a0006b910066b5cf2e3e240c530 Mon Sep 17 00:00:00 2001 From: Enne Hekma Date: Tue, 31 Jan 2017 13:28:19 +0100 Subject: [PATCH 8/8] Update script with structured I/O and empty JSON file. --- config/plot_tle_catalog.json.empty | 22 ++ ...plotcatalog_asr.py => plot_tle_catalog.py} | 164 +++++++----- python/plotcatalog_literature.py | 246 ------------------ 3 files changed, 115 insertions(+), 317 deletions(-) create mode 100644 config/plot_tle_catalog.json.empty rename python/{plotcatalog_asr.py => plot_tle_catalog.py} (56%) delete mode 100644 python/plotcatalog_literature.py diff --git a/config/plot_tle_catalog.json.empty b/config/plot_tle_catalog.json.empty new file mode 100644 index 0000000..3d56e16 --- /dev/null +++ b/config/plot_tle_catalog.json.empty @@ -0,0 +1,22 @@ +// Copyright (c) 2014-2016 Kartik Kumar, Dinamica Srl (me@kartikkumar.com) +// Copyright (c) 2016 Enne Hekma, Delft University of Technology (ennehekma@gmail.com) +// Distributed under the MIT License. +// See accompanying file LICENSE.md or copy at http://opensource.org/licenses/MIT +{ + // Path to SQLite database containing scan data. Maximum of + "databases" : [ "", + "", + "", + ""], + // Directory where scan figure is stored. + "output_directory" : ".", + + // Set parameters for figures. + "figure_dpi" : 300, + + // # Set number of lines per entry in TLE catalog (2 or 3). + "tleEntryNumberOfLines" : 3, + + // # Set font size for axes labels. + "fontSize" : 22 +} diff --git a/python/plotcatalog_asr.py b/python/plot_tle_catalog.py similarity index 56% rename from python/plotcatalog_asr.py rename to python/plot_tle_catalog.py index 2050840..0decc43 100644 --- a/python/plotcatalog_asr.py +++ b/python/plot_tle_catalog.py @@ -9,22 +9,7 @@ ################################################################################################### # Set path to TLE catalog file. -tleCatalogFilePathall = "../data/catalog/3le.txt" -tleCatalogFilePathSSO = "../data/SSO/SSO_tle.txt" -tleCatalogFilePathGEO = "../data/GEO/GEO_tle.txt" -tleCatalogFilePathHEO = "../data/HEO/HEO_tle.txt" -# Set number of lines per entry in TLE catalog (2 or 3). -tleEntryNumberOfLines = 3 - -# Set path to output directory. -outputPath = "../data/SSO/plots/" - -# Set figure DPI. -figureDPI = 300 - -# Set font size for axes labels. -fontSize = 22 ################################################################################################### @@ -45,6 +30,28 @@ import numpy as np import pandas as pd from twoBodyMethods import convertMeanMotionToSemiMajorAxis +import sys + +# I/O +import commentjson +import json +from pprint import pprint +import sqlite3 + + +# Parse JSON configuration file +# Raise exception if wrong number of inputs are provided to script +if len(sys.argv) != 2: + raise Exception("Only provide a JSON config file as input!") + +json_data = open(sys.argv[1]) +config = commentjson.load(json_data) +json_data.close() +pprint(config) + +fontSize = config['fontSize'] +figureDPI = config['figure_dpi'] +# exit() ################################################################################################### @@ -61,13 +68,10 @@ aop3 = [] inclinations3 = [] -# colors = ['k','b','g','r'] -# markers = ['.','s','+','D'] - -order = ['all','SSO','GEO','HEO'] -for x in order: - tleCatalogFilePathNew = eval(str('tleCatalogFilePath' + str(x))) - +# order = ['all','SSO','GEO','HEO'] +for x in xrange(0,len(config['databases'])): + # tleCatalogFilePathNew = eval(str('tleCatalogFilePath' + str(x))) + tleCatalogFilePathNew = config['databases'][x] # Read in catalog and store lines in list. fileHandle = open(tleCatalogFilePathNew) catalogLines = fileHandle.readlines() @@ -79,7 +83,7 @@ # Parse TLE entries and store debris objects. debrisObjects = [] - for tleEntry in xrange(0,len(catalogLines),tleEntryNumberOfLines): + for tleEntry in xrange(0,len(catalogLines),config['tleEntryNumberOfLines']): debrisObjects.append(twoline2rv(catalogLines[tleEntry+1], catalogLines[tleEntry+2], wgs72)) # Sort list of debris objects based on inclination. @@ -155,17 +159,21 @@ plt.xlabel("Semi-major axis altitude [km]") plt.ylabel("Eccentricity [-]") plt.ticklabel_format(style='sci', axis='x', scilimits=(0,0)) -plt.plot(sma['all'],ecc3['all'], \ - marker='.', markersize=1, color='k', linestyle='none') -plt.plot(sma['SSO'],ecc3['SSO'], \ - marker='s', markersize=10, color='c', linestyle='none') -plt.plot(sma['GEO'],ecc3['GEO'], \ - marker='^', markersize=10, color='g', linestyle='none') -plt.plot(sma['HEO'],ecc3['HEO'], \ - marker='D', markersize=6, color='r', linestyle='none') -axis.set_xlim(xmax=0.4e5) +plt.plot(sma['0'],ecc3['0'], \ + marker='.', markersize=1, color='k', linestyle='none') +if len(config['databases']) > 1: + plt.plot(sma['1'],ecc3['1'], \ + marker='s', markersize=10, color='c', linestyle='none') + if len(config['databases']) > 2: + plt.plot(sma['2'],ecc3['2'], \ + marker='^', markersize=10, color='g', linestyle='none') + if len(config['databases']) > 3: + plt.plot(sma['3'],ecc3['3'], \ + marker='D', markersize=6, color='r', linestyle='none') +# axis.set_xlim(xmax=0.4e5) figure.set_tight_layout(True) -plt.savefig(outputPath + "/figure1_debrisPopulation_eccentricityVsSemiMajorAxis.pdf", \ +plt.savefig(config['output_directory'] + + "/figure1_debrisPopulation_eccentricityVsSemiMajorAxis.pdf", dpi = figureDPI) plt.close() @@ -174,21 +182,26 @@ axis = figure.add_subplot(111) plt.xlabel("$e \cos{\omega}$ [-]") plt.ylabel("$e \sin{\omega}$ [-]") -plt.plot(ecc3['all']*np.cos(aop3['all']),ecc3['all']*np.sin(aop3['all']), \ - marker='.', markersize=1, color='k', linestyle='none') -plt.plot(ecc3['SSO']*np.cos(aop3['SSO']),ecc3['SSO']*np.sin(aop3['SSO']), \ - marker='s', markersize=10, color='c', linestyle='none') -plt.plot(ecc3['GEO']*np.cos(aop3['GEO']),ecc3['GEO']*np.sin(aop3['GEO']), \ - marker='^', markersize=10, color='g', linestyle='none') -plt.plot(ecc3['HEO']*np.cos(aop3['HEO']),ecc3['HEO']*np.sin(aop3['HEO']), \ - marker='D', markersize=6, color='r', linestyle='none') +plt.plot(ecc3['0']*np.cos(aop3['0']),ecc3['0']*np.sin(aop3['0']), \ + marker='.', markersize=1, color='k', linestyle='none') +if len(config['databases']) > 1: + plt.plot(ecc3['1']*np.cos(aop3['1']),ecc3['1']*np.sin(aop3['1']), \ + marker='s', markersize=10, color='c', linestyle='none') + if len(config['databases']) > 2: + plt.plot(ecc3['2']*np.cos(aop3['2']),ecc3['2']*np.sin(aop3['2']), \ + marker='^', markersize=10, color='g', linestyle='none') + if len(config['databases']) > 3: + plt.plot(ecc3['3']*np.cos(aop3['3']),ecc3['3']*np.sin(aop3['3']), \ + marker='D', markersize=6, color='r', linestyle='none') plt.axis('equal') -axis.set_xlim(xmin=-.82, xmax=.82) -axis.set_ylim(ymin=-.82, ymax=.82) -axis.set(xticks=[-.8,-.4,0,.4,.8]) -axis.set(yticks=[-.8,-.4,0,.4,.8]) +# axis.set_xlim(xmin=-.82, xmax=.82) +# axis.set_ylim(ymin=-.82, ymax=.82) +# axis.set(xticks=[-.8,-.4,0,.4,.8]) +# axis.set(yticks=[-.8,-.4,0,.4,.8]) figure.set_tight_layout(True) -plt.savefig(outputPath + "/figure2_debrisPopulation_eccentricityVector.pdf", dpi = figureDPI) +plt.savefig(config['output_directory'] + + "/figure2_debrisPopulation_eccentricityVector.pdf", + dpi = figureDPI) plt.close() @@ -198,17 +211,21 @@ plt.xlabel("Semi-major axis altitude [km]") plt.ylabel("Inclination [deg]") plt.ticklabel_format(style='sci', axis='x', scilimits=(0,0)) -plt.plot(sma['all'],incl5['all'], \ - marker='.', markersize=1, color='k', linestyle='none') -plt.plot(sma['SSO'],incl5['SSO'], \ - marker='s', markersize=10, color='c', linestyle='none') -plt.plot(sma['GEO'],incl5['GEO'], \ - marker='^', markersize=10, color='g', linestyle='none') -plt.plot(sma['HEO'],incl5['HEO'], \ - marker='D', markersize=6, color='r', linestyle='none') -axis.set_xlim(xmax=0.4e5) +plt.plot(sma['0'],incl5['0'], \ + marker='.', markersize=1, color='k', linestyle='none') +if len(config['databases']) > 1: + plt.plot(sma['1'],incl5['1'], \ + marker='s', markersize=10, color='c', linestyle='none') + if len(config['databases']) > 2: + plt.plot(sma['2'],incl5['2'], \ + marker='^', markersize=10, color='g', linestyle='none') + if len(config['databases']) > 3: + plt.plot(sma['3'],incl5['3'], \ + marker='D', markersize=6, color='r', linestyle='none') +# axis.set_xlim(xmax=0.4e5) figure.set_tight_layout(True) -plt.savefig(outputPath + "/figure3_debrisPopulation_inclinationVsSemiMajorAxis.pdf", \ +plt.savefig(config['output_directory'] + + "/figure3_debrisPopulation_inclinationVsSemiMajorAxis.pdf", dpi = figureDPI) plt.close() @@ -217,23 +234,28 @@ axis = figure.add_subplot(111) plt.xlabel("$i \cos{\Omega}$ [deg]") plt.ylabel("$i \sin{\Omega}$ [deg]") -plt.plot(np.rad2deg(inclinations3['all'])*np.cos(raan3['all']), \ - np.rad2deg(inclinations3['all'])*np.sin(raan3['all']), \ - marker='.', markersize=1, color='k', linestyle='none') -plt.plot(np.rad2deg(inclinations3['SSO'])*np.cos(raan3['SSO']), \ - np.rad2deg(inclinations3['SSO'])*np.sin(raan3['SSO']), \ - marker='s', markersize=10, color='c', linestyle='none') -plt.plot(np.rad2deg(inclinations3['GEO'])*np.cos(raan3['GEO']), \ - np.rad2deg(inclinations3['GEO'])*np.sin(raan3['GEO']), \ - marker='^', markersize=10, color='g', linestyle='none') -plt.plot(np.rad2deg(inclinations3['HEO'])*np.cos(raan3['HEO']), \ - np.rad2deg(inclinations3['HEO'])*np.sin(raan3['HEO']), \ - marker='D', markersize=6, color='r', linestyle='none') +plt.plot(np.rad2deg(inclinations3['0'])*np.cos(raan3['0']), + np.rad2deg(inclinations3['0'])*np.sin(raan3['0']), + marker='.', markersize=1, color='k', linestyle='none') +if len(config['databases']) > 1: + plt.plot(np.rad2deg(inclinations3['1'])*np.cos(raan3['1']), \ + np.rad2deg(inclinations3['1'])*np.sin(raan3['1']), \ + marker='s', markersize=10, color='c', linestyle='none') + if len(config['databases']) > 2: + plt.plot(np.rad2deg(inclinations3['2'])*np.cos(raan3['2']), \ + np.rad2deg(inclinations3['2'])*np.sin(raan3['2']), \ + marker='^', markersize=10, color='g', linestyle='none') + if len(config['databases']) > 3: + plt.plot(np.rad2deg(inclinations3['3'])*np.cos(raan3['3']), \ + np.rad2deg(inclinations3['3'])*np.sin(raan3['3']), \ + marker='D', markersize=6, color='r', linestyle='none') plt.axis('equal') -axis.set_xlim(xmin=-110.0, xmax=110.0) -axis.set_ylim(ymin=-110.0, ymax=110.0) +# axis.set_xlim(xmin=-110.0, xmax=110.0) +# axis.set_ylim(ymin=-110.0, ymax=110.0) figure.set_tight_layout(True) -plt.savefig(outputPath + "/figure4_debrisPopulation_inclinationVector.pdf", dpi = figureDPI) +plt.savefig(config['output_directory'] + + "/figure4_debrisPopulation_inclinationVector.pdf", + dpi = figureDPI) plt.close() diff --git a/python/plotcatalog_literature.py b/python/plotcatalog_literature.py deleted file mode 100644 index cdb5fcb..0000000 --- a/python/plotcatalog_literature.py +++ /dev/null @@ -1,246 +0,0 @@ -''' -Copyright (c) 2014-2016, K. Kumar (me@kartikkumar.com) -Copyright (c) 2016, E.J. Hekma (ennehekma@gmail.com) -All rights reserved. -''' - -################################################################################################### -# Set up input deck -################################################################################################### - -# Set path to TLE catalog file. -tleCatalogFilePathall = "../data/SSO/leo_all.txt" -tleCatalogFilePathSSO = "../data/SSO/ARIANE.txt" -tleCatalogFilePathGEO = "../data/SSO/SSO_tle.txt" -# tleCatalogFilePathHEO = "../data/SSO/ARIANE_RB.txt" -tleCatalogFilePathHEO = "../data/SSO/catalog_800_99_01.txt" - -# Set number of lines per entry in TLE catalog (2 or 3)self. -tleEntryNumberOfLines = 3 - -# Set path to output directory. -outputPath = "../data/plots_literature/" - -# Set figure DPI. -figureDPI = 300 - -# Set font size for axes labels. -fontSize = 22 - -################################################################################################### - -''' - DO NOT EDIT PARAMETERS BEYOND THIS POINT!!! -''' - -################################################################################################### -# Set up modules and packages -################################################################################################### - -from sgp4.earth_gravity import wgs72 -from sgp4.io import twoline2rv -from sgp4.propagation import getgravconst - -from matplotlib import rcParams -import matplotlib.pyplot as plt -import numpy as np -import pandas as pd -from twoBodyMethods import convertMeanMotionToSemiMajorAxis - -################################################################################################### - -################################################################################################### -# Read and store TLE catalog -################################################################################################### - -sma = [] -ecc5 = [] -incl5 = [] - -raan3 = [] -ecc3 = [] -aop3 = [] -inclinations3 = [] - -# colors = ['k','b','g','r'] -# markers = ['.','s','+','D'] - -order = ['all','SSO','GEO','HEO'] -for x in order: - tleCatalogFilePathNew = eval(str('tleCatalogFilePath' + str(x))) - - # Read in catalog and store lines in list. - fileHandle = open(tleCatalogFilePathNew) - catalogLines = fileHandle.readlines() - fileHandle.close() - - # Strip newline and return carriage characters. - for i in xrange(len(catalogLines)): - catalogLines[i] = catalogLines[i].strip('\r\n') - - # Parse TLE entries and store debris objects. - debrisObjects = [] - for tleEntry in xrange(0,len(catalogLines),tleEntryNumberOfLines): - debrisObjects.append(twoline2rv(catalogLines[tleEntry+1], catalogLines[tleEntry+2], wgs72)) - - # Sort list of debris objects based on inclination. - inclinationSortedObjects = sorted(debrisObjects, key=lambda x: x.inclo, reverse=False) - - inclinations = [] - raan = [] - ecc = [] - aop = [] - for i in xrange(len(inclinationSortedObjects)): - inclinations.append(inclinationSortedObjects[i].inclo) - raan.append(inclinationSortedObjects[i].nodeo) - ecc.append(inclinationSortedObjects[i].ecco) - aop.append(inclinationSortedObjects[i].argpo) - - smatemp = [] - smatemp = [convertMeanMotionToSemiMajorAxis(debrisObject.no/60.0, \ - getgravconst('wgs72')[1])-6373 \ - for debrisObject in debrisObjects] - smatemp2 = [] - smatemp2 = pd.DataFrame(smatemp, columns=[str(x)]) - sma.append(smatemp2) - - ecctemp = [] - ecctemp = [debrisObject.ecco for debrisObject in debrisObjects] - ecc4 = [] - ecc4 = pd.DataFrame(ecctemp, columns=[str(x)]) - ecc5.append(ecc4) - - incltemp = [] - incltemp = [np.rad2deg(debrisObject.inclo) for debrisObject in debrisObjects] - incl4 = [] - incl4 = pd.DataFrame(incltemp, columns=[str(x)]) - incl5.append(incl4) - - raan2 = [] - raan2 = pd.DataFrame(raan, columns=[str(x)]) - raan3.append(raan2) - - ecc2 = [] - ecc2 = pd.DataFrame(ecc, columns=[str(x)]) - ecc3.append(ecc2) - - aop2 = [] - aop2 = pd.DataFrame(aop, columns=[str(x)]) - aop3.append(aop2) - - inclinations2 = [] - inclinations2 = pd.DataFrame(inclinations, columns=[str(x)]) - inclinations3.append(inclinations2) - -sma = pd.concat(sma, axis=1) -ecc5 = pd.concat(ecc5, axis=1) -incl5 = pd.concat(incl5, axis=1) - -raan3 = pd.concat(raan3, axis=1) -ecc3 = pd.concat(ecc3, axis=1) -aop3 = pd.concat(aop3, axis=1) -inclinations3 = pd.concat(inclinations3, axis=1) - -################################################################################################### - -################################################################################################### -# Generate plots -################################################################################################### - -# Set font size for plot labels. -rcParams.update({'font.size': fontSize}) - -# Plot distribution of eccentricity [-] against semi-major axis [km]. -figure = plt.figure() -axis = figure.add_subplot(111) -plt.xlabel("Semi-major axis altitude [km]") -plt.ylabel("Eccentricity [-]") -plt.ticklabel_format(style='sci', axis='x', scilimits=(0,0)) -plt.plot(sma['all'],ecc3['all'], \ - marker='.', markersize=1, color='k', linestyle='none') -# plt.plot(sma['SSO'],ecc3['SSO'], \ -# marker='s', markersize=10, color='c', linestyle='none') -# plt.plot(sma['GEO'],ecc3['GEO'], \ -# marker='^', markersize=10, color='g', linestyle='none') -plt.plot(sma['HEO'],ecc3['HEO'], \ - marker='D', markersize=6, color='r', linestyle='none') -axis.set_xlim(xmax=1500) -axis.set_ylim(ymax=0.04) -figure.set_tight_layout(True) -plt.savefig(outputPath + "/figure1_debrisPopulation_eccentricityVsSemiMajorAxis.pdf", \ - dpi = figureDPI) -plt.close() - -# Plot components of eccentricity vector [-]. -figure = plt.figure() -axis = figure.add_subplot(111) -plt.xlabel("$e \cos{\omega}$ [-]") -plt.ylabel("$e \sin{\omega}$ [-]") -plt.plot(ecc3['all']*np.cos(aop3['all']),ecc3['all']*np.sin(aop3['all']), \ - marker='.', markersize=1, color='k', linestyle='none') -# plt.plot(ecc3['SSO']*np.cos(aop3['SSO']),ecc3['SSO']*np.sin(aop3['SSO']), \ -# marker='s', markersize=10, color='c', linestyle='none') -# plt.plot(ecc3['GEO']*np.cos(aop3['GEO']),ecc3['GEO']*np.sin(aop3['GEO']), \ -# marker='^', markersize=10, color='g', linestyle='none') -plt.plot(ecc3['HEO']*np.cos(aop3['HEO']),ecc3['HEO']*np.sin(aop3['HEO']), \ - marker='D', markersize=6, color='r', linestyle='none') -plt.axis('equal') -axis.set_xlim(xmin=-.021, xmax=.021) -axis.set_ylim(ymin=-.021, ymax=.021) -axis.set(xticks=[-.02,-.01,0,.01,.02]) -axis.set(yticks=[-.02,-.01,0,.01,.02]) -figure.set_tight_layout(True) -plt.savefig(outputPath + "/figure2_debrisPopulation_eccentricityVector.pdf", dpi = figureDPI) -plt.close() - - -# Plot distribution of inclination [deg] against semi-major axis [km]. -figure = plt.figure() -axis = figure.add_subplot(111) -plt.xlabel("Semi-major axis altitude [km]") -plt.ylabel("Inclination [deg]") -# plt.ticklabel_format(style='sci', axis='x', scilimits=(0,0)) -# plt.plot(sma['all'],incl5['all'], \ -# marker='.', markersize=1, color='k', linestyle='none') -# # plt.plot(sma['SSO'],incl5['SSO'], \ -# marker='s', markersize=10, color='c', linestyle='none') -# plt.plot(sma['GEO'],incl5['GEO'], \ -# marker='^', markersize=10, color='g', linestyle='none') -plt.plot(sma['HEO'],incl5['HEO'], \ - marker='D', markersize=2, color='r', linestyle='none') -axis.set_xlim(xmax=900,xmin=700) -axis.set_ylim(ymin = 98, ymax=100) -figure.set_tight_layout(True) -plt.savefig(outputPath + "/figure3_debrisPopulation_inclinationVsSemiMajorAxis.pdf", \ - dpi = figureDPI) -plt.close() - -# Plot components of inclination vector [deg]. -figure = plt.figure() -axis = figure.add_subplot(111) -plt.xlabel("$i \cos{\Omega}$ [deg]") -plt.ylabel("$i \sin{\Omega}$ [deg]") -plt.plot(np.rad2deg(inclinations3['all'])*np.cos(raan3['all']), \ - np.rad2deg(inclinations3['all'])*np.sin(raan3['all']), \ - marker='.', markersize=1, color='k', linestyle='none') -# plt.plot(np.rad2deg(inclinations3['SSO'])*np.cos(raan3['SSO']), \ -# np.rad2deg(inclinations3['SSO'])*np.sin(raan3['SSO']), \ -# marker='s', markersize=10, color='c', linestyle='none') -# plt.plot(np.rad2deg(inclinations3['GEO'])*np.cos(raan3['GEO']), \ -# np.rad2deg(inclinations3['GEO'])*np.sin(raan3['GEO']), \ -# marker='^', markersize=10, color='g', linestyle='none') -plt.plot(np.rad2deg(inclinations3['HEO'])*np.cos(raan3['HEO']), \ - np.rad2deg(inclinations3['HEO'])*np.sin(raan3['HEO']), \ - marker='D', markersize=6, color='r', linestyle='none') -plt.axis('equal') -axis.set_xlim(xmin=-110.0, xmax=110.0) -axis.set_ylim(ymin=-110.0, ymax=110.0) -figure.set_tight_layout(True) -plt.savefig(outputPath + "/figure4_debrisPopulation_inclinationVector.pdf", dpi = figureDPI) -plt.close() - - - - -################################################################################################### -###################################################################################################