import sympy as sp from dataclasses import dataclass, field from typing import Dict, Tuple, Iterable from Export import Export from Datenbank import Datenbankzugriff @dataclass class StochastischesModell: n_beob: int sigma_beob: Iterable[float] =None #σ a priori der einzelnen Beobachtung gruppe_beob: Iterable[int] =None #Gruppenzugehörigkeit jeder Beobachtung (Distanz, Richtung, GNSS, Nivellement,...,) sigma0_gruppe: Dict[int, float] = field(default_factory=dict) #σ0² für jede Gruppe def __post_init__(self): # Defaults setzen if self.sigma_beob is None: self.sigma_beob = [1.0] * int(self.n_beob) if self.gruppe_beob is None: self.gruppe_beob = [1] * int(self.n_beob) # In SymPy-Spaltenvektoren umwandeln self.sigma_beob = sp.Matrix(list(self.sigma_beob)) self.gruppe_beob = sp.Matrix(list(self.gruppe_beob)) # Dimension prüfen if self.sigma_beob.rows != self.gruppe_beob.rows: raise ValueError("sigma_beob und gruppe_beob müssen gleich viele Einträge haben.") if self.sigma_beob.rows != int(self.n_beob): raise ValueError("n_beob passt nicht zur Länge von sigma_beob / gruppe_beob.") # Fehlende Gruppen mit sigma0_sq = 1.0 ergänzen unique_groups = sorted({int(g) for g in self.gruppe_beob}) for g in unique_groups: if g not in self.sigma0_gruppe: self.sigma0_gruppe[g] = 1.0 def berechne_Qll(self) -> Tuple[sp.Matrix, sp.Matrix]: n = self.n_beob Q_ll = sp.zeros(n, n) P = sp.zeros(n, n) for i in range(self.n_beob): sigma_i = self.sigma_beob[i, 0] #σ-Wert der i-ten Beobachtung holen g = int(self.gruppe_beob[i, 0]) #Gruppenzugehörigkeit der Beobachtung bestimmen sigma0_sq = self.sigma0_gruppe[g] #Den Varianzfaktor der Gruppe holen q_ii = sigma_i**2 #σ² berechnen Q_ll[i, i] = q_ii #Diagonale return Q_ll def Qll_symbolisch(self, pfad_datenbank, liste_beobachtungen_symbolisch): liste_standardabweichungen_symbole = [] liste_beobachtungen_symbolisch = [str(b) for b in liste_beobachtungen_symbolisch] Qll = sp.zeros(len(liste_beobachtungen_symbolisch), len(liste_beobachtungen_symbolisch)) db_zugriff = Datenbankzugriff(pfad_datenbank) dict_beobachtungenID_instrumenteID = db_zugriff.get_instrumenteID_beobachtungenID_dict() for i, beobachtung_symbolisch_i in enumerate(liste_beobachtungen_symbolisch): aufgeteilt_i = beobachtung_symbolisch_i.split("_") beobachtungenID_i = int(aufgeteilt_i[0]) instrumenteID_i = dict_beobachtungenID_instrumenteID[beobachtungenID_i] beobachtungsart_i = str(aufgeteilt_i[1]) if beobachtungsart_i == "SD": stabw_apriori_konstant = sp.Symbol(f"stabw_apriori_konstant_{beobachtungsart_i}_{instrumenteID_i}") stabw_apriori_streckenprop = sp.Symbol(f"stabw_apriori_streckenprop_{beobachtungsart_i}_{instrumenteID_i}") tachymeter_distanz = sp.Symbol(f"SD_{beobachtungenID_i}") sigma = sp.sqrt(stabw_apriori_konstant ** 2 + (stabw_apriori_streckenprop * tachymeter_distanz / 1000000) ** 2) liste_standardabweichungen_symbole.append(sigma) Qll[i, i] = sigma ** 2 elif beobachtungsart_i == "R" or beobachtungsart_i == "ZW": stabw_apriori_konstant = sp.Symbol(f"stabw_apriori_konstant_{beobachtungsart_i}_{instrumenteID_i}") stabw_apriori_konstant_distanz = sp.Symbol(f"stabw_apriori_konstant_SD_{instrumenteID_i}") tachymeter_distanz = sp.Symbol(f"SD_{beobachtungenID_i}") sigma = sp.sqrt( stabw_apriori_konstant ** 2 + (stabw_apriori_konstant_distanz / tachymeter_distanz) ** 2) liste_standardabweichungen_symbole.append(sigma) Qll[i, i] = sigma ** 2 for j in range(i + 1, len(liste_beobachtungen_symbolisch)): beobachtung_symbolisch_j = liste_beobachtungen_symbolisch[j] aufgeteilt_j = beobachtung_symbolisch_j.split("_") beobachtungsart_j = str(aufgeteilt_j[1]) if beobachtungsart_i == "SD" and beobachtungsart_j == "SD": Qll[i, j] = 0 Qll[j, i] = 0 Export.matrix_to_csv(r"Zwischenergebnisse\Qll_Symbolisch.csv", liste_beobachtungen_symbolisch, liste_beobachtungen_symbolisch, Qll, "Qll") return Qll def Qll_numerisch(self, pfad_datenbank, Qll_Matrix_Symbolisch, liste_beobachtungen_symbolisch): db_zugriff = Datenbankzugriff(pfad_datenbank) dict_genauigkeiten = db_zugriff.get_genauigkeiten_dict() dict_beobachtungenID_instrumenteID = db_zugriff.get_instrumenteID_beobachtungenID_dict() liste_beobachtungen = db_zugriff.get_beobachtungen_from_beobachtungenid() dict_beobachtungenID_distanz = {} for standpunkt, zielpunkt, beobachtungenID, beobachtungsgruppeID, tachymeter_richtung, tachymeter_zenitwinkel, tachymeter_distanz in liste_beobachtungen: dict_beobachtungenID_distanz[int(beobachtungenID)] = tachymeter_distanz dict_genauigkeiten_neu = {} for genauigkeitenID, eintrag in dict_genauigkeiten.items(): instrumenteID = int(eintrag[0]) beobachtungsart = str(eintrag[1]) stabw_apriori_konstant = eintrag[2] stabw_apriori_streckenprop = eintrag[3] dict_genauigkeiten_neu[(instrumenteID, beobachtungsart)] = (stabw_apriori_konstant, stabw_apriori_streckenprop) substitutionen = {} dict_konstante_sd = {} for (instrumenteID, beobachtungsart), (stabw_apriori_konstant, stabw_apriori_streckenprop) in dict_genauigkeiten_neu.items(): if beobachtungsart == "Tachymeter_Strecke": if stabw_apriori_konstant is not None: dict_konstante_sd[instrumenteID] = float(stabw_apriori_konstant) for (instrumenteID, beobachtungsart), (stabw_apriori_konstant, stabw_apriori_streckenprop) in dict_genauigkeiten_neu.items(): if beobachtungsart == "Tachymeter_Strecke": beobachtungsart_kurz = "SD" elif beobachtungsart == "Tachymeter_Richtung": beobachtungsart_kurz = "R" elif beobachtungsart == "Tachymeter_Zenitwinkel": beobachtungsart_kurz = "ZW" if stabw_apriori_konstant is not None: substitutionen[sp.Symbol(f"stabw_apriori_konstant_{beobachtungsart_kurz}_{instrumenteID}")] = float(stabw_apriori_konstant) if stabw_apriori_streckenprop is not None: substitutionen[sp.Symbol(f"stabw_apriori_streckenprop_{beobachtungsart_kurz}_{instrumenteID}")] = float(stabw_apriori_streckenprop) for instrumenteID, wert in dict_konstante_sd.items(): substitutionen[sp.Symbol(f"stabw_apriori_konstant_SD_{instrumenteID}")] = float(wert) liste_beobachtungen_symbolisch = [str(b) for b in liste_beobachtungen_symbolisch] for beobachtung_symbolisch in liste_beobachtungen_symbolisch: aufgeteilt = beobachtung_symbolisch.split("_") beobachtungenID = int(aufgeteilt[0]) distanz = dict_beobachtungenID_distanz.get(beobachtungenID, None) if distanz is not None: substitutionen[sp.Symbol(f"SD_{beobachtungenID}")] = float(distanz) Qll_numerisch = Qll_Matrix_Symbolisch.xreplace(substitutionen) Export.matrix_to_csv( r"Zwischenergebnisse\Qll_Numerisch.csv", liste_beobachtungen_symbolisch, liste_beobachtungen_symbolisch, Qll_numerisch, "Qll" ) return Qll_numerisch def berechne_P(Q_ll): P = Q_ll.inv() return P def berechne_Qvv(self, A: sp.Matrix, P: sp.Matrix, Q_xx: sp.Matrix) -> sp.Matrix: Q_vv = P.inv() - A * Q_xx * A.T return Q_vv #Kofaktormatrix der Beobachtungsresiduen def berechne_R(self, Q_vv: sp.Matrix, P: sp.Matrix) -> sp.Matrix: R = Q_vv * P return R #Redundanzmatrix def berechne_r(self, R: sp.Matrix) -> sp.Matrix: n = R.rows r = sp.zeros(n, 1) for i in range(n): r[i, 0] = R[i, i] return r #Redundanzanteile