Wie schreibe ich Daten in eine .dat Datei Python neue Zeile mehrere Listen Satellite Orbit

Question

Ich habe eine Funktion zum Erzeugen und Verbreiten einer Satellitenumlaufbahn erstellt. Jetzt möchte ich alle meine Daten in einer .dat-Datei speichern. Ich bin mir nicht sicher, wie viele For-Loops ich brauche oder wie ich das machen soll.

Ich möchte Zeit, Breite, Länge und Höhe alle auf einer Zeile für jeden Ausbreitungsschritt.

-Code für Daten:

myOrbitJ2000Time = [1085.0, 2170.0, 3255.0, 4340.1, 5425.1] 

lat = [48.5, 26.5, -28.8, -48.1, 0.1] 

lon = [-144.1, -50.4, -1.6, 91.5, 152.8] 

alt = [264779.5, 262446.1, 319661.8, 355717.3, 306129.0] 

gewünschte Ausgabe in .dat-Datei:

J2000 Zeit, Lat, Lon, Alt

1085.0, 48.6, -144.2, 264779.5 

2170.0, 26.5, -50.4, 262446.1 

3255.0, -28.7, -1.6, 319661.8 

4340.1, -48.0, 91.5, 355717.2 

5425.1, 0.1, 152.8, 06129.0 

Voll Code:

import orbital 
from orbital import earth, KeplerianElements, plot 
import matplotlib.pyplot as plt 
import numpy as np 
from astropy import time 
from astropy.time import TimeDelta, Time 
from astropy import units as u 
from astropy import coordinates as coord 


#def orbitPropandcoordTrans(orbitPineapple_J2000time, _ecc, _inc, _raan, _arg_pe, _meananom, meanMotion): 
def orbitPropandcoordTrans(propNum, orbitPineapple_J2000time, _ecc, _inc, _raan, _arg_pe, _meananom, meanMotion): 
     ''' 
     Create original orbit and run for 100 propagations (in total one whole orbit) 
     in order to get xyz and time for each propagation step. 
     The end goal is to plot the lat, lon, & alt data to see if it matches ISS groundtrace. 
     ''' 
    import orbital 
    from orbital import earth, KeplerianElements, plot 
    import matplotlib.pyplot as plt 
    import numpy as np 
    from astropy import time 
    from astropy.time import TimeDelta, Time 
    from astropy import units as u 
    from astropy import coordinates as coord 

    'Calculate Avg. Period from Mean Motion' 
    _avgPeriod = 86400/meanMotion 
    #print('_avgPeriod', _avgPeriod) 

    ###### 
    #propNum = int(_avgPeriod) 

    'Generate Orbit' 
    #orbitPineapple = KeplerianElements.with_period(5560, body=earth, e=0.0004525, i=(np.deg2rad(51.6414)), raan=(np.deg2rad(247.1662))) 
    orbitPineapple = KeplerianElements.with_period(_avgPeriod, body=earth, e=_ecc, i=(np.deg2rad(_inc)), raan=(np.deg2rad(_raan)), arg_pe=(np.deg2rad(_arg_pe)), M0=(np.deg2rad(_meananom))) #ref_epoch= 
    #plot(orbitPineapple) 
    #plt.show() 
    #print('') 
    #print('') 

    'Propagate Orbit and retrieve xyz' 
    myOrbitX = []   #X Coordinate for propagated orbit step 
    myOrbitY = []   #Y Coordinate for propagated orbit step 
    myOrbitZ = []   #Z Coordinate for propagated orbit step 
    myOrbitTime = []  #Time for each propagated orbit step 
    myOrbitJ2000Time = [] #J2000 times 
    #propNum = 100  #Number of propagations and Mean Anomaly size (one orbit 2pi/propNum) 

    for i in range(propNum): 
     orbitPineapple.propagate_anomaly_by(M=(2.0*np.pi/propNum)) #Propagate the orbit by the Mean Anomaly 
     myOrbitX.append(orbitPineapple.r.x)      #x vals 
     myOrbitY.append(orbitPineapple.r.y)      #y vals 
     myOrbitZ.append(orbitPineapple.r.z)      #z vals 
     myOrbitTime.append(orbitPineapple_J2000time)    #J2000 time vals 
     myOrbitJ2000Time.append(orbitPineapple.t) 

     #plot(orbitPineapple) 
    #print('X:', 'Length:', len(myOrbitX)) 
    #print(myOrbitX) 
    #print('Y:', 'Length:',len(myOrbitY)) 
    #print(myOrbitY) 
    #print('Z:', 'Length:', len(myOrbitZ)) 
    #print(myOrbitZ) 
    #print('J2000 Time:', 'Length:',len(myOrbitTime)) 
    #print(myOrbitTime) 


    '''Because the myOrbitTime is only the time between each step to be the sum of itself plus 
    all the previous times. And then I need to convert that time from seconds after J2000 to UTC.''' 
    myT = [] #UTC time list 

    for i in range(propNum): 
     myT.append((Time(2000, format='jyear') + TimeDelta(myOrbitTime[i]*u.s)).iso) #Convert time from J2000 to UTC 
    #print('UTC Time List Length:', len(myT)) 
    #print('UTC Times:', myT) 

    '''Now I have xyz and time for each propagation step and need to convert the coordinates from 
    ECI to Lat, Lon, & Alt''' 
    now = []  #UTC time at each propagation step 
    xyz =[]  #Xyz coordinates from OrbitalPy initial orbit propagation 
    cartrep = [] #Cartesian Representation 
    gcrs = [] #Geocentric Celestial Reference System/Geocentric Equatorial Inertial, the default coord system of OrbitalPy 
    itrs =[]  #International Terrestrial Reference System coordinates 
    lat = []  #Longitude of the location, for the default ellipsoid 
    lon = []  #Longitude of the location, for the default ellipsoid 
    alt = []  #Height of the location, for the default ellipsoid 


    for i in range(propNum): 
     xyz = (myOrbitX[i], myOrbitY[i], myOrbitZ[i])     #Xyz coord for each prop. step 
     now = time.Time(myT[i])           #UTC time at each propagation step 
     cartrep = coord.CartesianRepresentation(*xyz, unit=u.m)   #Add units of [m] to xyz 
     gcrs = coord.GCRS(cartrep, obstime=time.Time(myT[i]))   #Let AstroPy know xyz is in GCRS 
     itrs = gcrs.transform_to(coord.ITRS(obstime=time.Time(myT[i]))) #Convert GCRS to ITRS 
     loc = coord.EarthLocation(*itrs.cartesian.xyz)     #Get lat/lon/height from ITRS 
     lat.append(loc.lat.value)          #Create latitude list 
     lon.append(loc.lon.value)          #Create longitude list 
     alt.append(loc.height.value)   

    #print('Current Time:', now) 
    #print('') 
    #print('UTC Time:') 
    #print(myT) 
    print('myOrbitJ2000Time', myOrbitJ2000Time) 
    print('') 
    print('Lat:') 
    print('') 
    print(lat) 
    print('') 
    print('Lon:') 
    print(lon) 
    print('') 
    print('Alt:') 
    print(alt) 

orbitPropandcoordTrans (5, -34963095, 0,0073662, 51,5946, 154,8079, 103,6176, 257,3038, 15,92610159)

Answer 1

Um Ihre allgemeine Frage zu beantworten, können Sie, anstatt eine Reihe von Python-Listen zu definieren (an die Sie nur langsam anhängen und mit denen Sie arbeiten, insbesondere beim Umgang mit einer großen Anzahl von Werten beim Skalieren Ihrer Auflösung) Erstellen von Numpy-Arrays anstelle von Listen Bei der Initialisierung von Numpy-Arrays ist es normalerweise notwendig, die Größe des Arrays anzugeben, aber in diesem Fall ist das einfach, da Sie wissen, wie viele Propagierungen Sie wollen. Zum Beispiel:

>>> orbitX = np.empty(propNum, dtype=float) 

schafft eine leere Numpy Array von propNum Gleitkommawerte (hier „leer“ bedeutet, wird das Array einfach nicht mit irgendwelchen Werten initialisiert - dies ist der schnellste Weg, ein neues Array, wie wir erstellen ist geht später sowieso darin alle Werte füllen

dann in der Schleife statt orbitX.append(x) verwenden Sie auf den Wert in den Array zuweisen würden den aktuellen Tick entsprechen:.. orbitX[i] = x das Gleiche gilt für die anderen Fälle

.

Dann gibt es mehrere Möglichkeiten, wie man herausführt Setzen Sie Ihre Daten, aber mit der Astropy Table bietet viel Flexibilität. Sie können eine Tabelle erstellen, die Spalten enthält, die Sie leicht wie wollen:

>>> from astropy.table import Table 
>>> table = Table([J2000, lat, lon, alt], names=('J2000', 'lat', 'lon', 'alt')) 

Die schöne an mit dem Table-Objekt ist, dass es eine Tonne von Optionen für das Ausgabeformat ist. Z.B. an dem Python-Prompt:

>>> table 
    J2000   lat   lon   alt 
    float64  float64  float64  float64 
------------- -------------- -------------- ------------- 
1085.01128806 48.5487129749 -144.175054697 264779.500823 
2170.02257613 26.5068122883 -50.3805485685 262446.085716 
3255.03386419 -28.7915478311 -1.6090285674 319661.817451 
4340.04515225 -48.0536526356 91.5416828221 355717.274021 
5425.05644032 0.084422392655 152.811717713 306129.02576 

Zur Ausgabe in eine Datei zuerst müssen Sie überlegen, wie Sie die Daten formatiert werden sollen. Es gibt bereits viele gängige Datenformate, die Sie in Betracht ziehen können, aber es hängt davon ab, wofür die Daten stehen und wer sie verbraucht (".dat file "bedeutet nichts in Bezug auf Dateiformate, oder vielmehr könnte es bedeuten irgendetwas. Aber in Ihrer Frage haben Sie angegeben, dass was Sie wollen ist, was" kommagetrennte Werte "(CSV) heißt Daten werden durch ein Komma getrennt in einer Reihe mit jedem Wert mit Spalt nach unten geht, formatieren Das Table Klasse ausgeben kann CSV (und jede Variante) ganz einfach:..

>>> import sys 
>>> table.write(sys.stdout, format='ascii.csv') # Note: I'm just using sys.stdout for demonstration purposes; normally you would give a filename 
J2000,lat,lon,alt 
1085.011288063746,48.54871297493748,-144.17505469691633,264779.5008225624 
2170.022576127492,26.506812288280788,-50.38054856853237,262446.0857159357 
3255.0338641912376,-28.79154783108773,-1.6090285674024463,319661.81745081506 
4340.045152254984,-48.05365263557444,91.54168282208444,355717.2740210131 
5425.05644031873,0.08442239265500713,152.81171771323176,306129.0257600865 

Es gibt viele andere Möglichkeiten, obwohl zum Beispiel Wenn Sie möchten, dass die Daten in berechtigten Spalten formatiert sind, können Sie das auch tun.Sie können mehr darüber lesen here. (Ich würde auch vorschlagen, dass wenn Sie ein einfaches Textformat wünschen, versuchen Sie ascii.ecsv, die den Vorteil hat, dass es Ausgänge addit ional Metadaten, die es leicht wieder in eine Astropy Tabelle) gelesen werden kann:

>>> table.write(sys.stdout, format='ascii.ecsv') 
# %ECSV 0.9 
# --- 
# datatype: 
# - {name: J2000, datatype: float64} 
# - {name: lat, datatype: float64} 
# - {name: lon, datatype: float64} 
# - {name: alt, datatype: float64} 
# schema: astropy-2.0 
J2000 lat lon alt 
1085.01128806 48.5487129749 -144.175054697 264779.500823 
2170.02257613 26.5068122883 -50.3805485685 262446.085716 
3255.03386419 -28.7915478311 -1.6090285674 319661.817451 
4340.04515225 -48.0536526356 91.5416828221 355717.274021 
5425.05644032 0.084422392655 152.811717713 306129.02576 

Ein weiteres unabhängiges, was ich ist zu beachten Sie werde feststellen, dass viele Objekte in Astropy über ganzen Arrays arbeiten können, zusätzlich zu den Einzelwerten und dieses kann oft viel effizienter sein, besonders für viele Werte. Insbesondere haben Sie diese Python-Schleife:

for i in range(propNum): 
    xyz = (myOrbitX[i], myOrbitY[i], myOrbitZ[i])     #Xyz coord for each prop. step 
    now = time.Time(myT[i])           #UTC time at each propagation step 
    cartrep = coord.CartesianRepresentation(*xyz, unit=u.m)   #Add units of [m] to xyz 
    gcrs = coord.GCRS(cartrep, obstime=time.Time(myT[i]))   #Let AstroPy know xyz is in GCRS 
    itrs = gcrs.transform_to(coord.ITRS(obstime=time.Time(myT[i]))) #Convert GCRS to ITRS 
    loc = coord.EarthLocation(*itrs.cartesian.xyz)     #Get lat/lon/height from ITRS 
    lat.append(loc.lat.value)          #Create latitude list 
    lon.append(loc.lon.value)          #Create longitude list 
    alt.append(loc.height.value) 

Diese neu geschrieben, ohne eine explizite Schleife vollständig werden kann, indem sie als Arrays von Koordinaten statt zu behandeln. Zum Beispiel:

>>> times = time.Time(myT) # All the times, not just a single one 
>>> cartrep = coord.CartesianRepresentation(orbitX, orbitY, orbitZ, unit=u.m) # Again, an array of coordinates 
>>> gcrs = coord.GCRS(cartrep, obstime=times) 
>>> itrs = gcrs.transform_to(coord.ITRS(obstime=ts)) 
>>> loc = coord.EarthLocation(*itrs.cartesian.xyz) # I'm not sure this is the most efficient way to do this but I'm not an expert on the coordinates package 

Mit diesem können wir alle Koordinaten auch als Arrays erhalten. Beispiel:

>>> loc.lat 
<Latitude [ 48.54871297, 26.50681229,-28.79154783,-48.05365264, 
      0.08442239] deg> 

Dies ist im Allgemeinen eine effizientere Möglichkeit, mit einer großen Anzahl von Koordinatenwerten zu arbeiten. Ähnlich, wenn Sie die myT in Ihrem ursprünglichen Code konvertieren, anstatt alle Zeitversätze zu durchlaufen, können Sie ein einzelnes TimeDelta Array erstellen und das zu Ihrer Anfangszeit hinzufügen.

Leider bin ich kein Experte für das orbital Paket, aber es scheint nicht so einfach zu sein, wie man die Koordinaten an verschiedenen Punkten auf einer Umlaufbahn bekommen möchte. ISTM, wie es sein sollte.

Answer 2

Wahrscheinlich ist der einfachste Weg zip() zu verwenden:

data = zip(myOrbitJ2000Time, lat, lon, alt)

So haben "Daten" das Format, das Sie serialisieren möchten. Vergessen Sie nicht, die gewünschte Kopfzeile zu serialisieren.