VILNIAUSUNIVERSITETAS FIZIKOSFAKULTETAS TEORINĖSFIZIKOSKATEDRA JonasNarkeliūnas ANUMERICALAPPROACHTOESTIMATETHEBALLISTICCOEFFICIENTOFSPACE DEBRISFROMTLEORBITALDATA Mokslotiriamasisdarbas (studijųprograma–TEORINĖFIZIKAIRASTROFIZIKA) Studentas JonasNarkeliūnas Darbovadovas dr. JanStupl Katedrosvedėjas Vilnius2016 Contents Introduction ................................................................................................ 3 1 AtmosphericDrag ..................................................................................... 5 1.1 DragForceandDragCoefficient........................................................... 5 1.2 AtmosphericDensityModel................................................................ 7 2 NumericalApproach .................................................................................. 10 2.1 InputData .................................................................................... 10 2.2 BallisticCoefficient.......................................................................... 12 2.3 IntegrationUsingSGP4..................................................................... 13 2.4 Fast-FourierTransformation ................................................................ 14 3 Results .................................................................................................. 18 Conclusion ................................................................................................. 24 References .................................................................................................. 25 2 Introduction LowEarthOrbit(LEO)isfullofspacedebris,whichconsistofspentrocketstages,oldsatellitesand fragments from explosions and collisions. As of 2009, more than 21 000 orbital debris larger than 10cmareknowntoexist[1],andwhileitishardtotrackanythingsmallerthanthat,theestimated population of particles between 1 and 10 cm in diameter is approximately 500 000, whereas small as1cmexceeds100million. TheseobjectsorbitEarthwithhugekineticenergies–speedsusually exceed7km/s. Theshapeoftheirorbitvariesfromalmostcirculartohighlyellipticalandcoversall LEO,aregioninspacebetween160and2000kmabovesealevel. Unfortunately,LEOisalsothe placewheremostofouractivesatellitesaresituated1,aswellas,InternationalSpaceStation(ISS) andHubbleSpaceTelescope,whoseorbitsarearound400and550kmabovesealevel,respectively. This poses a real threat as debris can collide with satellites and deal substantial damage or even destroythem. Collisionsbetweentwoormoredebriscreatecloudsofsmallerdebris,whicharehardertotrack andincreaseoverallobjectdensityandcollisionprobability. Atsomepoint,thedebrisdensitycould thenreachacriticalvalue,whichwouldstartachainreactionandthenumberofspacedebriswould grow exponentially. This phenomenon was first described by Kessler in 1978 and he concluded that it would lead to creation of debris belt, which would vastly complicate satellite operations in LEO.[2]. The debris density is already relatively high, as seen from several necessary debris avoidance maneuvers donebyShuttle,beforeitwasdiscontinued,andISS[1,3,4]. Butnotallsatelliteshave apropulsionsystemtoavoidcollision,hencedifferentmethodsneedtobeapplied. Oneofthepro- posed collision avoidance concepts is called LightForce and it suggests using photon pressure to inducesmallorbitalcorrectionstodeflectdebrisfromcolliding[5,6]. Thismethodisveryefficient – as seen from theoretical simulations, even few continuous mode 10 kW ground-based lasers, fo- cused by 1.5 m telescopes with adaptive optics, were enough to prevent significant amount of the debriscollisions. SimulationsweredonebypropagatingallspaceobjectsinLEOby1yearintothe futureandcheckingwhethertheprobabilityofcollisionwashigh. Forthosespaceobjectsdifferent ground-based lasers were used to divert them, afterwards collision probabilities were reevaluated. However,theactualaccuracyoftheLightForcesoftware,whichhasbeendevelopedatNASAAmes ResearchCenter,dependsontheveracityoftheinputparameters,oneofwhichistheobject’sbal- listiccoefficient. Itisameasureofbody’sabilitytoovercomeairresistance,whichhasasignificant impact on the debris in LEO, and thus it is responsible for the shape of the trajectory of the de- bris. Having the exact values of the ballistic coefficient would make significantly better collision predictions,unfortunately,wedonotknowwhatarethevaluesformostoftheobjects. In this research, we were working with part of LightForce code, which estimates the ballis- 1Source: https://www.space-track.org. 3 tic coefficient from ephemerides2. Previously used method gave highly inaccurate values, when compared to known objects, and it needed to be changed. The goal of this work was to try out a differentmethodofestimatingtheballisticcoefficientandtocheckwhetherornotitgivesnoticeable improvementsinaccuracy. Herewepresentournewapproachandtheideasbehindit. 2Ephemeris(pluralephemerides)givesthepositionofastronomicalobjectsaswellasspacedebrisintheskyata giventime. 4 1 Atmospheric Drag Everything above Karman line, 100 km above sea level, is assumed to be outer space, however, atmosphere extends far beyond that point. It extends for at least 10 000 km more [7], although after 600 - 700 km it is not very significant (see Figure 1.1). Nevertheless, besides gravity, the atmosphericdrageffecthasthemostimpactontheobject’strajectoryinLEOandithastobetaken into account when modeling orbits of space objects in LEO. From here everything will eventually fall down to Earth, but the magnitude of the effect depends on many different parameters and we presenttheminthissection. 1.1 Drag Force and Drag Coefficient An object moving through medium always experiences drag – a resistance force acting in the op- posite direction of object’s movement. It is a result of three main effects: surface friction, flow interference and form drag. At relatively high velocities and low densities, like in our case, the lateral effect overshadows the first two, so they can be neglected [9]. In this case, the drag force is largelyproportionaltocharacteristicareaoftheobjectandtothesquareofitsvelocity. Figure1.1Variationsofairdensitywithheightfrom150to1000kmforhighandlowsolaractivity and diurnal variations (see chapter 1.2) based on COSPAR International Reference Atmosphere (CIRA1972). Courtesy: King-Hele[8] 5 Usingdimensionalanalysisfromfluiddynamics,wecandefinedragasfollows: F =−1ρ(cid:175)(cid:175)v 2(cid:175)(cid:175)C A, (1.1) d rel d 2 whereF isthedragforce;ρ isadensityofmedium,inourcase,thedensityofair;v isarelative d rel velocityofanobjectgoingthroughmedium;andAisacharacteristicareaofthatobject. C isadrag d coefficient; it has no units and depends on the object’s surface roughness and shape, as seen from Figure1.2–italsodependsonthemedium’sdensityandthemolecularvelocity[9]. Nevertheless, the dependence of different parameters for specific conditions are usually too small to take into account,thusitismostlyusedasaconstant. Whenmodelingaspaceobject,weusuallydonotknowitsshape,becausethisinformationisnot savedinephemerides,asseenfromchapter2.1,hencewecannotknowwhatistheexactvalueofthe dragcoefficient. Evenifwecouldcaptureanimageoftheobject’sprofileandmeasurethesurface area,C and A valueswouldn’tberight,asmostoftheobjectsarenon-sphericalandspinningwith d respect to their direction of motion, hence one would need multiple pictures of different sides, as wellas,toknowitsangularvelocity. Agoodwaytodealwiththisistoapproximatethemasspheres and then use drag coefficient of a sphere. The value we use isC =2.2, which was popularized by d theuseofJacchiamodelsforestimatingatmosphericdensityandotherparametersofupperlayers Figure1.2Dragcoefficientfordifferentshapeobjectscalculatedatlowvelocities[10] 6 of Earth’s atmosphere [11,12]. Note that this value is different from the one showed in Figure 1.2 asaresultofthefollowingassumptions[9]: • Themeanfreepathoftheairmoleculesismuchgreaterthanthesizeofthesatellite. • Thesatellitevelocityv ismuchgreaterthanthemeanrandomvelocityv¯oftheairmolecules. 0 • The mechanisms that deflect the impinging molecules range from elastic to inelastic impact anddiffuseevaporationataneffectivesurfacetemperatureT. 1.2 Atmospheric Density Model Asonecansee,thedragforcedependsontheatmosphericdensityindifferentways,butthedensity is not constant. Even though we assume that its changes do little effect to the drag coefficient, for thedragforcewehave touseaccuratevalues,asobjectsclosertothegroundwill interactwithita lotmoreanddeorbitquickerthanthosewhoarefurtheraway;itiseasytoseethisfromtheISSand Vanguard13 orbits. WhiletheISSisthebiggestman-madesatellite,ithasoneofthelowestorbits of all artificial satellites and thus it loses 2 or more kilometers of altitude per month, as seen from Figure 1.3. Without frequent orbital corrections it would deorbit within a year. Vanguard 1 is in a higher orbit, with a perigee4 of 656 km and apogee5 of 3842 km (as of December 2015), and it is expectedtostayinorbitfor240years6. Figure 1.3 Orbital height of the ISS throughout last year. Image source: www.heavens- above.com/ISSHeight.aspx 3Vanguard1istheoldestmanmadesatellitestillinorbit,launchedin1958;itwasusedtoobtaingeodeticmea- surements. 4ApointwherecelestialbodyisclosesttotheEarth. 5ApointwherecelestialbodyisfurthestfromtheEarth. 6Retrievedfromhttp://nssdc.gsfc.nasa.gov/nmc/spacecraftDisplay.do?id=1958-002B. 7 The Earth’s atmosphere is expanding and contracting in a complex manner all the time, not to mentionitisnotthesamethroughoutdifferentlatitudesandlongitudes. Manyeffectscontributeto that and a good atmospheric model should include at least most significant ones, if not all. Below, welisteffectsthatareconsideredtobemostimportant[13]: • Latitudinal variations: The effect comes from satellite passing over the Earth’s equatorial bulge,effectivelythischangestheactualaltitudeand,therefore,thedensity. • Diurnal variations: These variations come from the Earth’s rotation and the Sun’s irradia- tion. The day side of the Earth is warmer than the night side, hence the atmosphere expands and bulges over ecliptic plane in the direction of the Sun. Because of the diurnal variations, it’snecessarytoknowlatitude,localtime,andtimeoftheyear. • 27-day solar-rotation cycle: This effect arises from Sun’s rotation around its axis, which causesfluctuationsintheamountofradiationthatreachestheEarth. • 11-year cycle of sunspots: Solar cycle strongly varies the amount of incoming radiation reaching the Earth. As a result, the Earth’s atmosphere heats up different amounts and it expands different amounts throughout the cycle. At the solar maximum the temperature of upperlayersoftheatmosphereisabouttwotimeshigherthanatsolarminimum[14],meaning that the atmosphere has expanded a lot more and debris will be affected more by the drag. In fact, nine times more debris fall to the surface of the Earth during solar maximum than minimum [15]. Furthermore, depending on the altitude, this effect can even cause larger disturbancethroughsolar-radiationpressurethandrag[13]. • Seasonal variations: These variations come from solar flux fluctuations due to Earth’s ro- tation around the Sun. The Earth’s rotational axis is tilted and its orbit is elliptical, thus the Sun’sdeclinationandthedistancebetweentheEarthandtheSunvariesthroughouttheyear. • Rotating atmosphere: The atmosphere is rotating together with the Earth, however, the velocityoflowerlayerisusuallyhigherbecauseofthefrictionwiththeEarth’ssurface. • Winds: Atmosphericwindsandweathercreatelocalvariationsofthetemperatureand,there- fore, the density. Unfortunately, this is difficult to account for, as the weather itself is almost unpredictablebecauseofitssensitivitytochanges. • Magnetic-storm variations: Effect comes from fluctuations in the Earth’s magnetic field, whichaffectstheamountofchargedparticlesreachingtheEarth’satmospherefromtheSun. • Irregular short-periodic variations: Mostly solar flares: sudden increase in solar wind intensity,whichalsodisturbstheEarth’smagneticfield. 8 • Tides: Ocean tides and even atmospheric tides as a result of gravitational pull of the Moon cancausesmallvariationsintheatmosphericdensity. The first model of the Earth’s atmosphere above 120 km was created by Nicolet in 1961 [16]. It wasanempiricalmodelandusedasetofboundaryconditionsoftemperatureandpartialdensities of gases altogether with diffusion theory to get values for higher altitudes. Later it was improved uponseveraltimesbyJacchiawithempiricaldatafromthenewsatellites[11,17,18]. Itwasagood model and even today it is still widely used for spacecraft dynamics. However, for more accurate calculations like satellite and space debris orbital simulations, it is better to use a newer model, developed by Mike Picone, Alan Hedin, and Doug Drobof, which is called NRLMSISE-00 [19]. It isanabbreviationofUSNavalResearchLaboratoryandMassSpectrometerandIncoherentScatter Radar (00 stands for the year of release), the two primary data sources for development of earlier versionsofthemodel. NRLMSISE-00 is an empirical model of the Earth’s atmosphere from the ground through ex- osphere. Its outputs are atmospheric density, which we need, as well as temperature and other parameters, while all the necessary input parameters are: year and day time of day; local apparent solar time; geodetic altitude, latitude and longitude; daily and 81 day average of F solar flux7; 10.7 daily geomagnetic index. This model includes most of the effects listed above [19], it is imple- mentedinLightForcecodeandthusitistheperfectmodelforoursimulations. 710.7cmsolarradioflux,F ,isagreatindicatorofsolaractivity,asitcorrelateswithsunspotsandUVradiation 10.7 [14],aswellasitiseasytomeasureitonthegroundusingradioantenna. 9 2 Numerical Approach Fromthepreviouschapterwelearnedhowtheatmosphericdragaffectsthedebrisandwhyitvaries in time, but the problem we face that it is very impractical to experimentally measure it without knowing the shape and the size of the debris. However, we can find ratios between mass, drag coefficient and characteristic surface area of debris – the ballistic coefficient, and it can be done numericallyfromtheorbitaltrackingdata. Itisexplainedinthissection. 2.1 Input Data Satellitetrackingbeganin1957withthelaunchofthefirstman-madesatelliteintospace–Sputnik. Thefirstmethodsfortrackingusedradiocommunicationandweren’tveryreliableastheycouldn’t detectortrackobjectsinspacethatweren’ttransmittinginternationallyagreedfrequencies,assome did. At that time, the main reason for satellite tracking was to defend against missile attacks, as a resultbettermethodswereconceived. Telescopeswithcameraswereusedinstead,aswellasradars, although they had their own limitations. Optical tracking systems depended on the weather and only worked during night time, while radar was free from those problems, but it was costly [20]. Nowadaysthemethodsarethesame,onlytheaccuracyhasincreasedabit,albeitstillrelativelylow, andobjectsassmallas1cmcouldnowbedetectedandtrackedwithladars[21],whichisgreatfor estimatingvariousparametersoftheatmosphere,aswellaspredictingtrajectories. Themostpopularwayforacelestialbody’sdetectionandtrackingisusingaradarsystem. The detection is done by sending a powerful, yet short, impulse of radio signal into particular part of space and waiting for echoes. Computer analyzes the returned signal and from its strength, shape andtimeintooktocomebackitcalculatesaverageephemerisforthatobject[22]. Tracking,onthe other hand, uses a priori knowledge where the object is going to be at specific time, so the radar is prepared by pointing to that part of the sky. After detection, it then tracks that object to get new orbital data, which is then saved as one data set. But there is only a finite number of radar sensor sites around the world, and they have a narrow view of the sky, and limited capabilities, meaning thatsomespaceobjectaren’tinviewallthetime. Furthermore,theymightnotbeinagoodviewor trackedforcoupledays,thusmeasurementsarerare–frompublicdatabase8 weseethatitisabout onepertwodaysorless. Theformatforstoringephemeridesisthesameasitwas50yearsago. Backthen,asetoforbital elements were written in two lines on a punch card9; from here came the name Two-Line Element (TLE) set, meaning this data encoding format. Given a reference frame and an epoch (a reference time when ephemeriswere taken), actually onlysix orbital elements areneeded to unambiguously 8www.celestrak.com/NORAD/elements/ 9Punchcardisapieceofstiffpaperthatcontainedeithercommandsforcontrollingautomatedmachineryordata for data processing applications. Both commands and data were represented by the presence or absence of holes in predefinedpositions. 10