Untitled

 avatar
unknown
plain_text
2 years ago
23 kB
6
Indexable
import argparse
import datetime
import gzip
import matplotlib.pyplot as plt
import matplotlib.image as mpimg
import numpy as np
import os
from ftplib import FTP
from matplotlib.animation import FuncAnimation

# Pfad zum Verzeichnis des Skripts
skript_verzeichnis = os.path.dirname(os.path.abspath(__file__))

# Name des neuen Ordners
neuer_ordner = "proj02SatelliteVisibility_DATA_2008"

# Pfad zum neuen Ordner
neuer_ordner_pfad = os.path.join(skript_verzeichnis, neuer_ordner)

# Überprüfen, ob der Ordner bereits existiert
if not os.path.exists(neuer_ordner_pfad):
    # Ordner erstellen
    os.mkdir(neuer_ordner_pfad)
    print("Der Ordner wurde erfolgreich erstellt.")
else:
    print("Der Ordner existiert bereits.")

#%% Parse all

parser = argparse.ArgumentParser(description='Generates an animated visualization of GRACE satellite visibility to the GPS constellation for a single day.')

# Argumente definieren
parser.add_argument('date', type=str, help='date to visualize (format: YYYY-MM-DD)')
parser.add_argument("-n", "--novisibility", help="don't show visibility lines", action="store_true")
parser.add_argument("-o", "--outfile", type=str, help="save animation to this file", nargs='*')
group = parser.add_mutually_exclusive_group()
group.add_argument("-t", "--time", type=str, help="start and end time in hours to limit animation", nargs='*')

# Argumente analysieren
args = parser.parse_args()

# Analysierte Argumente abrufen
date = args.date
novisibility = args.novisibility
outfile = args.outfile
time = args.time

# Überprüfen der Eingaben Date
if date[:4] == "2008":
    pass
else:
    print("Fehler: Daten sind nur verfügbar für das Jahr 2008.")
    exit(1)

# Überprüfung und Anpassung der Start- und Endzeiten
if time is not None:
    start_time_parts = time[0].split('.')
    end_time_parts = time[1].split('.')

    # Konvertierung der Stunden und Minuten in MJD
    start_hour = int(start_time_parts[0])
    start_minute = int(start_time_parts[1]) if len(start_time_parts) > 1 else 0
    end_hour = int(end_time_parts[0])
    end_minute = int(end_time_parts[1]) if len(end_time_parts) > 1 else 0

    # Überprüfung, ob Startzeit kleiner als Endzeit ist
    if start_hour > end_hour or (start_hour == end_hour and start_minute >= end_minute):
        print("Fehler: Die Startzeit muss kleiner als die Endzeit sein.")
        exit(1)

    # Überprüfung, ob Startzeit und Endzeit im Bereich von 0 bis 24 Stunden liegen
    if start_hour < 0 or start_hour >= 24 or end_hour < 0 or end_hour >= 24:
        print("Fehler: Die Startzeit und die Endzeit müssen im Bereich von 0 bis 24 Stunden liegen.")
        exit(1)

    # Umrechnung der Start- und Endzeiten in MJD
    current_date = datetime.datetime.strptime(date, "%Y-%m-%d")
    start_datetime = current_date.replace(hour=start_hour, minute=start_minute)
    end_datetime = current_date.replace(hour=end_hour, minute=end_minute)

    # Umrechnung in MJD (nur die Nachkommastellen)
    mjd_start = start_datetime.toordinal() - 678576 + start_datetime.hour/24 + start_datetime.minute/(24*60)
    mjd_end = end_datetime.toordinal() - 678576 + end_datetime.hour/24 + end_datetime.minute/(24*60)

    time_mjd = [str(mjd_start).split('.')[1], str(mjd_end).split('.')[1]]
    
    # Umrechnung in Minuten
    start_minute = str(start_minute * 6).zfill(2)
    end_minute = str(end_minute * 6).zfill(2)
    
    # Erstellen der Liste time mit den formatierten Start- und Endzeiten
    time = [f"{start_hour}:{start_minute}", f"{end_hour}:{end_minute}"]

# Die Werte der Argumente ausgeben
print("Date:", date)
print(f"Time: {time[0]} - {time[1]}")
print("Novisibility:", novisibility)
print("Outfile:", outfile)


#%% Daten herunterladen

# Extract the year, month, and day from the date
year = date[:4]
month = date[5:7]
day = date[8:10]

def download_blue_marble_image(date, target_folder):
    # FTP-Server Informationen
    server = "ftp.tugraz.at"

    # FTP-Ordnerpfad für Blue Marble-Bilder
    folder_path = "/outgoing/ITSG/teaching/2023SS_Informatik2/bluemarble/"

    # Dateiname des Blue Marble-Bildes
    image_filename = f"bluemarble{month}.jpg"

    # Überprüfen, ob die Bilddatei bereits im lokalen Ordner vorhanden ist
    local_file_path = os.path.join(target_folder, month, image_filename)
    if os.path.exists(local_file_path):
        print(f"Blue Marble-Bild für Monat {month} ist bereits vorhanden. Überspringe den Download.")
        return

    # Lokalen Ordnerpfad erstellen, falls er nicht vorhanden ist
    new_folder = os.path.join(target_folder, month)
    if not os.path.exists(new_folder):
        os.makedirs(new_folder)

    # Lokaler Dateipfad zum Speichern des Blue Marble-Bildes
    local_file_path = os.path.join(new_folder, image_filename)

    # Mit dem FTP-Server verbinden
    ftp = FTP(server)
    ftp.login()

    # Zum FTP-Ordner wechseln
    ftp.cwd(folder_path)

    # Das Blue Marble-Bild herunterladen
    with open(local_file_path, "wb") as local_file:
        ftp.retrbinary("RETR " + image_filename, local_file.write)

    print(f"Blue Marble-Bild für Monat {month} heruntergeladen.")

    # Die Verbindung zum FTP-Server trennen
    ftp.quit()

# Festlegen des Zielordners, in dem das Bild gespeichert werden soll
target_folder1 = os.path.join(neuer_ordner_pfad)

# Aufruf der Funktion, um das Blue Marble-Bild herunterzuladen
download_blue_marble_image(date, target_folder1)

def download_files_from_ftp(date, target_folder):
    # FTP-Server Informationen
    server = "ftp.tugraz.at"

    # FTP-Ordnerpfad
    folder_path = f"/outgoing/ITSG/teaching/2023SS_Informatik2/orbit/{date}/"

    # Erstellen der Ordnerstruktur
    folder_structure = os.path.join(target_folder1, month, day)

    # Überprüfen, ob die Dateien bereits im lokalen Ordner vorhanden sind
    if os.path.exists(folder_structure):
        print(f"Dateien für das Datum {date} sind bereits vorhanden. Überspringe den Download.")
        return

    # Lokalen Ordner erstellen, falls er nicht vorhanden ist
    os.makedirs(folder_structure)

    # Mit dem FTP-Server verbinden
    ftp = FTP(server)
    ftp.login()

    # Zum FTP-Ordner wechseln
    ftp.cwd(folder_path)

    # Liste der Dateien auf dem FTP-Server abrufen
    files = ftp.nlst()

    # Jede Datei herunterladen
    for file in files:
        local_file_path = os.path.join(folder_structure, file)

        # Die Datei herunterladen
        with open(local_file_path, "wb") as local_file:
            ftp.retrbinary("RETR " + file, local_file.write)

        print(f"Heruntergeladen: {file}")

    # Die Verbindung zum FTP-Server trennen
    ftp.quit()

# Zielordner für die Dateien festlegen
target_folder2 = os.path.join(target_folder1, month, day)

# Funktion aufrufen, um die Dateien herunterzuladen und in der angegebenen Ordnerstruktur zu speichern
download_files_from_ftp(date, target_folder2)


#%% Daten einlesen

# Pfad zum Ordner mit den heruntergeladenen Dateien
folder_path = target_folder2

# Dateinamen im Ordner auflisten
files = os.listdir(folder_path)
# print(files)

# Liste für Satellitennamen erstellen
satellite_names = []

# Satellitennamen aus Dateinamen extrahieren und zur Liste hinzufügen
for file in files:
    parts = file.split(".")
    satellite_name = parts[1]  # Den Satellitennamen extrahieren
    satellite_names.append(satellite_name)

# print(satellite_names)
# print(np.shape(satellite_names))

# Liste für den Inhalt der Textdateien erstellen
daten = []

for file in files:
    try:
        with gzip.open(f"{folder_path}/{file}", "rt") as f:
            content = np.loadtxt(f, skiprows=2)
            daten.append(content)
    except ValueError as e:
        print(f"Fehler beim Einlesen der Datei {file}: {e}")

#%% Daten vorbereiten

# Liste für die umgerechneten Koordinaten erstellen
umgerechnete_daten = []

for content in daten:
    # Spalten auswählen
    mjd = content[:, 0]
    x = content[:, 1]
    y = content[:, 2]
    z = content[:, 3]
    
    # MJD in gregorianische Datum umrechnen
    startdatum = datetime.datetime(1858, 11, 17)
    datum = [startdatum + datetime.timedelta(days=int(m), seconds=int((m % 1) * 24 * 60 * 60)) for m in mjd]

    # Formatieren des Datums im gewünschten Format
    formatiertes_datum = [d.strftime('%Y-%m-%d %H:%M') for d in datum]

    # Koordinaten umrechnen
    r = np.sqrt(x**2 + y**2 + z**2)
    phi = np.arcsin(z / r)
    lam = np.arctan2(y , x)

    # Umgerechnete Koordinaten zur Liste hinzufügen
    umgerechnete_daten.append(np.column_stack((formatiertes_datum, phi, lam, z)))

# Funktion zum Lesen der heruntergeladenen Dateien und Berechnen der Norm und des Winkels
def norm_angle_calculation(list_of_satellite_array):
    normed_satellite_arrays = []
    list_with_XYZandAngle = []

    for array in list_of_satellite_array:
        x = array[:, 1].astype(float)
        y = array[:, 2].astype(float)
        z = array[:, 3].astype(float)

        # Norm berechnen
        norm = np.sqrt(x**2 + y**2 + z**2)
        normed_satellite_arrays.append(norm)

        # Winkel berechnen und in Grad umrechnen
        angle_rad = np.arccos(z / norm)
        angle_deg = np.degrees(angle_rad)
        list_with_XYZandAngle.append(np.column_stack((x, y, z, angle_deg)))

    return normed_satellite_arrays, list_with_XYZandAngle

# Norm und Winkel berechnen
normed_data, data_with_angle = norm_angle_calculation(umgerechnete_daten)

# print(np.shape(normed_data))
# print(np.shape(data_with_angle))

# Liste für die Ausgabe erstellen
all_data = []

for i in range(len(umgerechnete_daten)):
    date = umgerechnete_daten[i][:, 0]
    phi = umgerechnete_daten[i][:, 1]
    lam = umgerechnete_daten[i][:, 2]
    alpha = data_with_angle[i][:, 3]
    
    all_data.append(np.column_stack((date, phi, lam, alpha)))

# Extrahieren der Start- und Endzeiten
if time:
    start_time = time[0]
    end_time = time[1]
else:
    start_time = "00:00"
    end_time = "23:59"

# Speichern des Bereichs zwischen Start- und Endzeiten
time_data = []

for content in all_data:
    dates = content[:, 0]
    times = [date.split()[1] for date in dates]
    start_index = times.index(start_time)
    end_index = times.index(end_time)
    filtered_content = content[start_index:end_index+1]
    time_data.append(filtered_content)

# TODO: print_data so struckturieren, dass aus der Liste(time_data) aus den ganzen zeilen jeweils die ersten zeile werte in einer unterliste sind usw

# print_data = [[] for _ in range(len(all_data[0][0]))]

# for data in time_data:
#     for i, content in enumerate(data.T):
#         print_data[i].extend(content)






#%% Visualisieren

# Laden des Blue Marble-Bildes
blue_marble_path = os.path.join(target_folder1, month, f"bluemarble{month}.jpg")
blue_marble_img = mpimg.imread(blue_marble_path)

# Anzeigen des Bildes
plt.imshow(blue_marble_img)

# # Iteriere über jede Zeile in print_data
# for row in print_data:
#     # Extrahiere Phi- und Lambda-Werte aus der aktuellen Zeile
#     phi_values = row[:, 1]
#     lambda_values = row[:, 2]

#     # Plotten der Phi- und Lambda-Werte als Scatter-Plot
#     plt.scatter(lambda_values, phi_values, marker='o')

# Anzeigen des Plots
plt.legend()
plt.show()

import datetime
import gzip
import matplotlib.pyplot as plt
import matplotlib.image as mpimg
import numpy as np
import os
from ftplib import FTP
from matplotlib.animation import FuncAnimation

# Pfad zum Verzeichnis des Skripts
skript_verzeichnis = os.path.dirname(os.path.abspath(__file__))

# Name des neuen Ordners
neuer_ordner = "proj02SatelliteVisibility_DATA_2008"

# Pfad zum neuen Ordner
neuer_ordner_pfad = os.path.join(skript_verzeichnis, neuer_ordner)

# Überprüfen, ob der Ordner bereits existiert
if not os.path.exists(neuer_ordner_pfad):
    # Ordner erstellen
    os.mkdir(neuer_ordner_pfad)
    print("Der Ordner wurde erfolgreich erstellt.")
else:
    print("Der Ordner existiert bereits.")

#%% Parse all

parser = argparse.ArgumentParser(description='Generates an animated visualization of GRACE satellite visibility to the GPS constellation for a single day.')

# Argumente definieren
parser.add_argument('date', type=str, help='date to visualize (format: YYYY-MM-DD)')
parser.add_argument("-n", "--novisibility", help="don't show visibility lines", action="store_true")
parser.add_argument("-o", "--outfile", type=str, help="save animation to this file", nargs='*')
group = parser.add_mutually_exclusive_group()
group.add_argument("-t", "--time", type=str, help="start and end time in hours to limit animation", nargs='*')

# Argumente analysieren
args = parser.parse_args()

# Analysierte Argumente abrufen
date = args.date
novisibility = args.novisibility
outfile = args.outfile
time = args.time

# Überprüfen der Eingaben Date
if date[:4] == "2008":
    pass
else:
    print("Fehler: Daten sind nur verfügbar für das Jahr 2008.")
    exit(1)

# Überprüfung und Anpassung der Start- und Endzeiten
if time is not None:
    start_time_parts = time[0].split('.')
    end_time_parts = time[1].split('.')

    # Konvertierung der Stunden und Minuten in MJD
    start_hour = int(start_time_parts[0])
    start_minute = int(start_time_parts[1]) if len(start_time_parts) > 1 else 0
    end_hour = int(end_time_parts[0])
    end_minute = int(end_time_parts[1]) if len(end_time_parts) > 1 else 0

    # Überprüfung, ob Startzeit kleiner als Endzeit ist
    if start_hour > end_hour or (start_hour == end_hour and start_minute >= end_minute):
        print("Fehler: Die Startzeit muss kleiner als die Endzeit sein.")
        exit(1)

    # Überprüfung, ob Startzeit und Endzeit im Bereich von 0 bis 24 Stunden liegen
    if start_hour < 0 or start_hour >= 24 or end_hour < 0 or end_hour >= 24:
        print("Fehler: Die Startzeit und die Endzeit müssen im Bereich von 0 bis 24 Stunden liegen.")
        exit(1)

    # Umrechnung der Start- und Endzeiten in MJD
    current_date = datetime.datetime.strptime(date, "%Y-%m-%d")
    start_datetime = current_date.replace(hour=start_hour, minute=start_minute)
    end_datetime = current_date.replace(hour=end_hour, minute=end_minute)

    # Umrechnung in MJD (nur die Nachkommastellen)
    mjd_start = start_datetime.toordinal() - 678576 + start_datetime.hour/24 + start_datetime.minute/(24*60)
    mjd_end = end_datetime.toordinal() - 678576 + end_datetime.hour/24 + end_datetime.minute/(24*60)

    time_mjd = [str(mjd_start).split('.')[1], str(mjd_end).split('.')[1]]
    
    # Umrechnung in Minuten
    start_minute = str(start_minute * 6).zfill(2)
    end_minute = str(end_minute * 6).zfill(2)
    
    # Erstellen der Liste time mit den formatierten Start- und Endzeiten
    time = [f"{start_hour}:{start_minute}", f"{end_hour}:{end_minute}"]

# Die Werte der Argumente ausgeben
print("Date:", date)
print(f"Time: {time[0]} - {time[1]}")
print("Novisibility:", novisibility)
print("Outfile:", outfile)


#%% Daten herunterladen

# Extract the year, month, and day from the date
year = date[:4]
month = date[5:7]
day = date[8:10]

def download_blue_marble_image(date, target_folder):
    # FTP-Server Informationen
    server = "ftp.tugraz.at"

    # FTP-Ordnerpfad für Blue Marble-Bilder
    folder_path = "/outgoing/ITSG/teaching/2023SS_Informatik2/bluemarble/"

    # Dateiname des Blue Marble-Bildes
    image_filename = f"bluemarble{month}.jpg"

    # Überprüfen, ob die Bilddatei bereits im lokalen Ordner vorhanden ist
    local_file_path = os.path.join(target_folder, month, image_filename)
    if os.path.exists(local_file_path):
        print(f"Blue Marble-Bild für Monat {month} ist bereits vorhanden. Überspringe den Download.")
        return

    # Lokalen Ordnerpfad erstellen, falls er nicht vorhanden ist
    new_folder = os.path.join(target_folder, month)
    if not os.path.exists(new_folder):
        os.makedirs(new_folder)

    # Lokaler Dateipfad zum Speichern des Blue Marble-Bildes
    local_file_path = os.path.join(new_folder, image_filename)

    # Mit dem FTP-Server verbinden
    ftp = FTP(server)
    ftp.login()

    # Zum FTP-Ordner wechseln
    ftp.cwd(folder_path)

    # Das Blue Marble-Bild herunterladen
    with open(local_file_path, "wb") as local_file:
        ftp.retrbinary("RETR " + image_filename, local_file.write)

    print(f"Blue Marble-Bild für Monat {month} heruntergeladen.")

    # Die Verbindung zum FTP-Server trennen
    ftp.quit()

# Festlegen des Zielordners, in dem das Bild gespeichert werden soll
target_folder1 = os.path.join(neuer_ordner_pfad)

# Aufruf der Funktion, um das Blue Marble-Bild herunterzuladen
download_blue_marble_image(date, target_folder1)

def download_files_from_ftp(date, target_folder):
    # FTP-Server Informationen
    server = "ftp.tugraz.at"

    # FTP-Ordnerpfad
    folder_path = f"/outgoing/ITSG/teaching/2023SS_Informatik2/orbit/{date}/"

    # Erstellen der Ordnerstruktur
    folder_structure = os.path.join(target_folder1, month, day)

    # Überprüfen, ob die Dateien bereits im lokalen Ordner vorhanden sind
    if os.path.exists(folder_structure):
        print(f"Dateien für das Datum {date} sind bereits vorhanden. Überspringe den Download.")
        return

    # Lokalen Ordner erstellen, falls er nicht vorhanden ist
    os.makedirs(folder_structure)

    # Mit dem FTP-Server verbinden
    ftp = FTP(server)
    ftp.login()

    # Zum FTP-Ordner wechseln
    ftp.cwd(folder_path)

    # Liste der Dateien auf dem FTP-Server abrufen
    files = ftp.nlst()

    # Jede Datei herunterladen
    for file in files:
        local_file_path = os.path.join(folder_structure, file)

        # Die Datei herunterladen
        with open(local_file_path, "wb") as local_file:
            ftp.retrbinary("RETR " + file, local_file.write)

        print(f"Heruntergeladen: {file}")

    # Die Verbindung zum FTP-Server trennen
    ftp.quit()

# Zielordner für die Dateien festlegen
target_folder2 = os.path.join(target_folder1, month, day)

# Funktion aufrufen, um die Dateien herunterzuladen und in der angegebenen Ordnerstruktur zu speichern
download_files_from_ftp(date, target_folder2)


#%% Daten einlesen

# Pfad zum Ordner mit den heruntergeladenen Dateien
folder_path = target_folder2

# Dateinamen im Ordner auflisten
files = os.listdir(folder_path)
# print(files)

# Liste für Satellitennamen erstellen
satellite_names = []

# Satellitennamen aus Dateinamen extrahieren und zur Liste hinzufügen
for file in files:
    parts = file.split(".")
    satellite_name = parts[1]  # Den Satellitennamen extrahieren
    satellite_names.append(satellite_name)

# print(satellite_names)
# print(np.shape(satellite_names))

# Liste für den Inhalt der Textdateien erstellen
daten = []

for file in files:
    try:
        with gzip.open(f"{folder_path}/{file}", "rt") as f:
            content = np.loadtxt(f, skiprows=2)
            daten.append(content)
    except ValueError as e:
        print(f"Fehler beim Einlesen der Datei {file}: {e}")

#%% Daten vorbereiten

# Liste für die umgerechneten Koordinaten erstellen
umgerechnete_daten = []

for content in daten:
    # Spalten auswählen
    mjd = content[:, 0]
    x = content[:, 1]
    y = content[:, 2]
    z = content[:, 3]
    
    # MJD in gregorianische Datum umrechnen
    startdatum = datetime.datetime(1858, 11, 17)
    datum = [startdatum + datetime.timedelta(days=int(m), seconds=int((m % 1) * 24 * 60 * 60)) for m in mjd]

    # Formatieren des Datums im gewünschten Format
    formatiertes_datum = [d.strftime('%Y-%m-%d %H:%M') for d in datum]

    # Koordinaten umrechnen
    r = np.sqrt(x**2 + y**2 + z**2)
    phi = np.arcsin(z / r)
    lam = np.arctan2(y , x)

    # Umgerechnete Koordinaten zur Liste hinzufügen
    umgerechnete_daten.append(np.column_stack((formatiertes_datum, phi, lam, z)))

# Funktion zum Lesen der heruntergeladenen Dateien und Berechnen der Norm und des Winkels
def norm_angle_calculation(list_of_satellite_array):
    normed_satellite_arrays = []
    list_with_XYZandAngle = []

    for array in list_of_satellite_array:
        x = array[:, 1].astype(float)
        y = array[:, 2].astype(float)
        z = array[:, 3].astype(float)

        # Norm berechnen
        norm = np.sqrt(x**2 + y**2 + z**2)
        normed_satellite_arrays.append(norm)

        # Winkel berechnen und in Grad umrechnen
        angle_rad = np.arccos(z / norm)
        angle_deg = np.degrees(angle_rad)
        list_with_XYZandAngle.append(np.column_stack((x, y, z, angle_deg)))

    return normed_satellite_arrays, list_with_XYZandAngle

# Norm und Winkel berechnen
normed_data, data_with_angle = norm_angle_calculation(umgerechnete_daten)

# print(np.shape(normed_data))
# print(np.shape(data_with_angle))

# Liste für die Ausgabe erstellen
all_data = []

for i in range(len(umgerechnete_daten)):
    date = umgerechnete_daten[i][:, 0]
    phi = umgerechnete_daten[i][:, 1]
    lam = umgerechnete_daten[i][:, 2]
    alpha = data_with_angle[i][:, 3]
    
    all_data.append(np.column_stack((date, phi, lam, alpha)))

# Extrahieren der Start- und Endzeiten
if time:
    start_time = time[0]
    end_time = time[1]
else:
    start_time = "00:00"
    end_time = "23:59"

# Speichern des Bereichs zwischen Start- und Endzeiten
time_data = []

for content in all_data:
    dates = content[:, 0]
    times = [date.split()[1] for date in dates]
    start_index = times.index(start_time)
    end_index = times.index(end_time)
    filtered_content = content[start_index:end_index+1]
    time_data.append(filtered_content)

# TODO: print_data so struckturieren, dass aus der Liste(time_data) aus den ganzen zeilen jeweils die ersten zeile werte in einer unterliste sind usw

# print_data = [[] for _ in range(len(all_data[0][0]))]

# for data in time_data:
#     for i, content in enumerate(data.T):
#         print_data[i].extend(content)






#%% Visualisieren

# Laden des Blue Marble-Bildes
blue_marble_path = os.path.join(target_folder1, month, f"bluemarble{month}.jpg")
blue_marble_img = mpimg.imread(blue_marble_path)

# Anzeigen des Bildes
plt.imshow(blue_marble_img)

# # Iteriere über jede Zeile in print_data
# for row in print_data:
#     # Extrahiere Phi- und Lambda-Werte aus der aktuellen Zeile
#     phi_values = row[:, 1]
#     lambda_values = row[:, 2]

#     # Plotten der Phi- und Lambda-Werte als Scatter-Plot
#     plt.scatter(lambda_values, phi_values, marker='o')

# Anzeigen des Plots
plt.legend()
plt.show()
Editor is loading...