Files
Masterprojekt_V3/Funktionales_Modell.py
2026-01-19 14:48:24 +01:00

915 lines
45 KiB
Python

from typing import Any
from numpy import ndarray, dtype
from sympy import MutableDenseMatrix
from sympy.matrices.expressions.matexpr import MatrixElement
from Datenbank import Datenbankzugriff
import sympy as sp
from Export import Export
from Berechnungen import Berechnungen
import numpy as np
import importlib
from Koordinatentransformationen import Transformationen
from pathlib import Path
import pandas as pd
import numpy as np
import sympy as sp
class FunktionalesModell:
def __init__(self, pfad_datenbank: str, a: float, b: float, pfad_tif_quasigeoidundolation: str = None) -> None:
self.pfad_datenbank = pfad_datenbank
self.a = a
self.b = b
self.berechnungen = Berechnungen(self.a, self.b)
self.trafos = Transformationen(pfad_datenbank)
self.pfad_tif_quasigeoidundolation = pfad_tif_quasigeoidundolation
self.substitutionen_dict = self.dict_substitutionen_uebergeordnetes_system()
self.dict_punkt_symbole = {}
self.liste_symbole_lambdify = sorted(self.substitutionen_dict.keys(), key=lambda s: str(s))
self.func_beob0 = None
self.func_A0 = None
self.func_u0 = None
self.liste_beobachtungsvektor_symbolisch = None
def jacobi_matrix_symbolisch(self, datumsfestlegung: str = None, liste_unbekannte_datumsfestlegung: list = None) -> tuple[MutableDenseMatrix | MatrixElement | list[Any] | Any, list[Any], list[Any]] | None:
#liste_beobachtungsarten = ["tachymeter_distanz", "tachymeter_richtung", "tachymeter_zenitwinkel"]
liste_beobachtungsarten = ["tachymeter_distanz", "tachymeter_richtung", "tachymeter_zenitwinkel", "gnss_basislinien", "geometrisches_nivellement"]
#liste_beobachtungsarten = ["tachymeter_distanz", "tachymeter_richtung", "tachymeter_zenitwinkel", "gnss_basislinien"]
db_zugriff = Datenbankzugriff(self.pfad_datenbank)
liste_beobachtungen_rohdaten_gnssbasislinien = []
liste_beobachtungen_rohdaten_tachymeter = []
liste_beobachtungen_rohdaten_nivellement = []
liste_punktnummern =[]
liste_orientierungsunbekannte = []
for beobachtungsart in liste_beobachtungsarten:
#Tachymeter Block
if beobachtungsart.startswith("tachymeter"):
liste_id_standpunkt_zielpunkt = db_zugriff.get_beobachtungen_id_beobachtungsgruppe_standpunkt_zielpunkt(beobachtungsart)
for beobachtungenID, beobachtungsgruppeID, standpunkt, zielpunkt in liste_id_standpunkt_zielpunkt:
liste_beobachtungen_rohdaten_tachymeter.append(
(beobachtungsart, beobachtungenID, beobachtungsgruppeID, standpunkt, zielpunkt)
)
if standpunkt not in liste_punktnummern:
liste_punktnummern.append(standpunkt)
if zielpunkt not in liste_punktnummern:
liste_punktnummern.append(zielpunkt)
if beobachtungsart == "tachymeter_richtung":
if beobachtungsgruppeID not in liste_orientierungsunbekannte:
liste_orientierungsunbekannte.append(beobachtungsgruppeID)
#GNSS Block
#if beobachtungsart == "gnss_basislinien":
# liste_id_standpunkt_zielpunkt = db_zugriff.get_gnss_beobachtungen_punktnummern()
# for beobachtungenID, standpunkt, zielpunkt in liste_id_standpunkt_zielpunkt:
# standpunkt = str(standpunkt).strip()
# zielpunkt = str(zielpunkt).strip()
# liste_beobachtungen_rohdaten_gnssbasislinien.append((beobachtungsart, beobachtungenID, standpunkt, zielpunkt))
# if standpunkt not in liste_punktnummern:
# liste_punktnummern.append(standpunkt)
# if zielpunkt not in liste_punktnummern:
# liste_punktnummern.append(zielpunkt)
#GNSS Block
if beobachtungsart == "gnss_basislinien":
liste_id_standpunkt_zielpunkt = db_zugriff.get_gnss_beobachtungen_punktnummern("gnss_bx")
for beobachtungenID, standpunkt, zielpunkt in liste_id_standpunkt_zielpunkt:
standpunkt = str(standpunkt).strip()
zielpunkt = str(zielpunkt).strip()
liste_beobachtungen_rohdaten_gnssbasislinien.append(("gnssbx", beobachtungenID, standpunkt, zielpunkt))
if standpunkt not in liste_punktnummern:
liste_punktnummern.append(standpunkt)
if zielpunkt not in liste_punktnummern:
liste_punktnummern.append(zielpunkt)
liste_id_standpunkt_zielpunkt = db_zugriff.get_gnss_beobachtungen_punktnummern("gnss_by")
for beobachtungenID, standpunkt, zielpunkt in liste_id_standpunkt_zielpunkt:
standpunkt = str(standpunkt).strip()
zielpunkt = str(zielpunkt).strip()
liste_beobachtungen_rohdaten_gnssbasislinien.append(("gnssby", beobachtungenID, standpunkt, zielpunkt))
if standpunkt not in liste_punktnummern:
liste_punktnummern.append(standpunkt)
if zielpunkt not in liste_punktnummern:
liste_punktnummern.append(zielpunkt)
liste_id_standpunkt_zielpunkt = db_zugriff.get_gnss_beobachtungen_punktnummern("gnss_bz")
for beobachtungenID, standpunkt, zielpunkt in liste_id_standpunkt_zielpunkt:
standpunkt = str(standpunkt).strip()
zielpunkt = str(zielpunkt).strip()
liste_beobachtungen_rohdaten_gnssbasislinien.append(("gnssbz", beobachtungenID, standpunkt, zielpunkt))
if standpunkt not in liste_punktnummern:
liste_punktnummern.append(standpunkt)
if zielpunkt not in liste_punktnummern:
liste_punktnummern.append(zielpunkt)
if beobachtungsart == "geometrisches_nivellement":
liste_id_standpunkt_zielpunkt = db_zugriff.get_nivellement_beobachtungen_punktnummern()
for beobachtungenID, standpunkt, zielpunkt in liste_id_standpunkt_zielpunkt:
standpunkt = str(standpunkt).strip()
zielpunkt = str(zielpunkt).strip()
liste_beobachtungen_rohdaten_nivellement.append(
(beobachtungsart, beobachtungenID, standpunkt, zielpunkt))
if standpunkt not in liste_punktnummern:
liste_punktnummern.append(standpunkt)
if zielpunkt not in liste_punktnummern:
liste_punktnummern.append(zielpunkt)
#if liste_beobachtungen_rohdaten_tachymeter == []:
# return None
#dict_punkt_symbole = {}
liste_unbekannte = []
for punkt in liste_punktnummern:
X, Y, Z = sp.symbols(f"X{punkt} Y{punkt} Z{punkt}")
self.dict_punkt_symbole[punkt] = (X, Y, Z)
liste_unbekannte.append(X)
liste_unbekannte.append(Y)
liste_unbekannte.append(Z)
dict_orientierung_symbole = {}
for orientierungsunbekannte in liste_orientierungsunbekannte:
O = sp.symbols(f"O{orientierungsunbekannte}")
dict_orientierung_symbole[orientierungsunbekannte] = O
liste_unbekannte.append(O)
liste_beobachtungsgleichungen_distanz =[]
liste_zeilenbeschriftungen_distanz = []
liste_A_richtung_zeilen = []
liste_zeilenbeschriftungen_richtung = []
liste_A_zenitwinkel_zeilen = []
liste_zeilenbeschriftungen_zenitwinkel = []
liste_beobachtungsgleichungen_gnssbasislinien = []
liste_A_gnssbasislinien_zeilen = []
liste_zeilenbeschriftungen_gnssbasislinien = []
liste_beobachtungsgleichungen_nivellement = []
liste_A_nivellement_zeilen = []
liste_zeilenbeschriftungen_nivellement = []
if liste_beobachtungen_rohdaten_tachymeter != []:
for beobachtungsart, beobachtungenID, beobachtungsgruppeID, standpunkt, zielpunkt in liste_beobachtungen_rohdaten_tachymeter:
X_sp, Y_sp, Z_sp = self.dict_punkt_symbole[standpunkt]
X_zp, Y_zp, Z_zp = self.dict_punkt_symbole[zielpunkt]
B_sp, L_sp = sp.symbols(f"B{standpunkt} L{standpunkt}")
alpha = sp.symbols(f"{beobachtungenID}_R_{beobachtungsgruppeID}_{standpunkt}_{zielpunkt}")
zw = sp.symbols(f"{beobachtungenID}_ZW_{beobachtungsgruppeID}_{standpunkt}_{zielpunkt}")
s = sp.symbols(f"{beobachtungenID}_SD_{beobachtungsgruppeID}_{standpunkt}_{zielpunkt}")
azimut_berechnet = sp.symbols(f"azimut_berechnet_{beobachtungsgruppeID}_{standpunkt}_{zielpunkt}")
zw_berechnet = sp.symbols(f"zw_berechnet_{beobachtungsgruppeID}_{standpunkt}_{zielpunkt}")
s_berechnet = sp.symbols(f"strecke_berechnet_{beobachtungsgruppeID}_{standpunkt}_{zielpunkt}")
if beobachtungsart == "tachymeter_distanz":
beobachtungsgleichung = sp.sqrt((X_zp - X_sp) ** 2 + (Y_zp - Y_sp) ** 2 + (Z_zp - Z_sp) ** 2)
liste_beobachtungsgleichungen_distanz.append(beobachtungsgleichung)
liste_zeilenbeschriftungen_distanz.append(
f"{beobachtungenID}_SD_{beobachtungsgruppeID}_{standpunkt}_{zielpunkt}")
elif beobachtungsart == "tachymeter_richtung":
# for beobachtungenID, beobachtungsgruppeID, standpunkt, zielpunkt in liste_id_standpunkt_zielpunkt:
d_r_dX_zp = (
(sp.sin(B_sp) * sp.cos(L_sp) * sp.sin(azimut_berechnet) - sp.sin(L_sp) * sp.cos(azimut_berechnet)) / (
s_berechnet * sp.sin(zw_berechnet)))
d_r_dX_sp = - d_r_dX_zp
d_r_dY_zp = (
(sp.sin(B_sp) * sp.sin(L_sp) * sp.sin(azimut_berechnet) + sp.cos(L_sp) * sp.cos(azimut_berechnet)) / (
s_berechnet * sp.sin(zw_berechnet)))
d_r_dY_sp = - d_r_dY_zp
d_r_dZ_zp = ((-sp.cos(B_sp) * sp.sin(azimut_berechnet) / (s_berechnet * sp.sin(zw_berechnet))))
d_r_dZ_sp = - d_r_dZ_zp
d_r_dO_sp = -1
zeile_A_Matrix = []
for punkt in liste_punktnummern:
if punkt == standpunkt:
zeile_A_Matrix.extend([d_r_dX_sp, d_r_dY_sp, d_r_dZ_sp])
elif punkt == zielpunkt:
zeile_A_Matrix.extend([d_r_dX_zp, d_r_dY_zp, d_r_dZ_zp])
else:
zeile_A_Matrix.extend([0, 0, 0])
for orientierung in liste_orientierungsunbekannte:
if orientierung == beobachtungsgruppeID:
zeile_A_Matrix.append(d_r_dO_sp)
else:
zeile_A_Matrix.append(0)
liste_A_richtung_zeilen.append(zeile_A_Matrix)
liste_zeilenbeschriftungen_richtung.append(
f"{beobachtungenID}_R_{beobachtungsgruppeID}_{standpunkt}_{zielpunkt}"
)
elif beobachtungsart == "tachymeter_zenitwinkel":
d_r_dX_zp = ((X_zp - X_sp) * sp.cos(zw_berechnet) - s_berechnet * sp.cos(B_sp) * sp.cos(L_sp)) / (s_berechnet ** 2 * sp.sin(zw_berechnet))
d_r_dX_sp = - d_r_dX_zp
d_r_dY_zp = ((Y_zp - Y_sp) * sp.cos(zw_berechnet) - s_berechnet * sp.cos(B_sp) * sp.sin(L_sp)) / (s_berechnet ** 2 * sp.sin(zw_berechnet))
d_r_dY_sp = - d_r_dY_zp
d_r_dZ_zp = ((Z_zp - Z_sp) * sp.cos(zw_berechnet) - s_berechnet * sp.sin(B_sp)) / (s_berechnet ** 2 * sp.sin(zw_berechnet))
d_r_dZ_sp = - d_r_dZ_zp
zeile_A_Matrix = []
for punkt in liste_punktnummern:
if punkt == standpunkt:
zeile_A_Matrix.extend([d_r_dX_sp, d_r_dY_sp, d_r_dZ_sp])
elif punkt == zielpunkt:
zeile_A_Matrix.extend([d_r_dX_zp, d_r_dY_zp, d_r_dZ_zp])
else:
zeile_A_Matrix.extend([0, 0, 0])
for orientierung in liste_orientierungsunbekannte:
zeile_A_Matrix.append(0)
liste_A_zenitwinkel_zeilen.append(zeile_A_Matrix)
liste_zeilenbeschriftungen_zenitwinkel.append(
f"{beobachtungenID}_ZW_{beobachtungsgruppeID}_{standpunkt}_{zielpunkt}"
)
#if liste_beobachtungen_rohdaten_gnssbasislinien != []:
# for beobachtungsart, beobachtungenID, standpunkt, zielpunkt in liste_beobachtungen_rohdaten_gnssbasislinien:
# X_sp, Y_sp, Z_sp = self.dict_punkt_symbole[standpunkt]
# X_zp, Y_zp, Z_zp = self.dict_punkt_symbole[zielpunkt]
# if beobachtungsart == "gnss_basislinien":
# beobachtungsgleichung_bx = X_zp - X_sp
# beobachtungsgleichung_by = Y_zp - Y_sp
# beobachtungsgleichung_bz = Z_zp - Z_sp
# liste_beobachtungsgleichungen_gnssbasislinien.append(beobachtungsgleichung_bx)
# liste_beobachtungsgleichungen_gnssbasislinien.append(beobachtungsgleichung_by)
# liste_beobachtungsgleichungen_gnssbasislinien.append(beobachtungsgleichung_bz)
# liste_zeilenbeschriftungen_gnssbasislinien.append(
# f"{beobachtungenID}_gnssbx_{standpunkt}_{zielpunkt}")
# liste_zeilenbeschriftungen_gnssbasislinien.append(
# f"{beobachtungenID}_gnssby_{standpunkt}_{zielpunkt}")
# liste_zeilenbeschriftungen_gnssbasislinien.append(
# f"{beobachtungenID}_gnssbz_{standpunkt}_{zielpunkt}")
if liste_beobachtungen_rohdaten_gnssbasislinien != []:
for beobachtungsart, beobachtungenID, standpunkt, zielpunkt in liste_beobachtungen_rohdaten_gnssbasislinien:
X_sp, Y_sp, Z_sp = self.dict_punkt_symbole[standpunkt]
X_zp, Y_zp, Z_zp = self.dict_punkt_symbole[zielpunkt]
if beobachtungsart == "gnssbx":
beobachtungsgleichung_bx = X_zp - X_sp
liste_beobachtungsgleichungen_gnssbasislinien.append(beobachtungsgleichung_bx)
liste_zeilenbeschriftungen_gnssbasislinien.append(
f"{beobachtungenID}_gnssbx_{standpunkt}_{zielpunkt}")
if beobachtungsart == "gnssby":
beobachtungsgleichung_by = Y_zp - Y_sp
liste_beobachtungsgleichungen_gnssbasislinien.append(beobachtungsgleichung_by)
liste_zeilenbeschriftungen_gnssbasislinien.append(
f"{beobachtungenID}_gnssby_{standpunkt}_{zielpunkt}")
if beobachtungsart == "gnssbz":
beobachtungsgleichung_bz = Z_zp - Z_sp
liste_beobachtungsgleichungen_gnssbasislinien.append(beobachtungsgleichung_bz)
liste_zeilenbeschriftungen_gnssbasislinien.append(
f"{beobachtungenID}_gnssbz_{standpunkt}_{zielpunkt}")
if liste_beobachtungen_rohdaten_nivellement != []:
for beobachtungsart, beobachtungenID, standpunkt, zielpunkt in liste_beobachtungen_rohdaten_nivellement:
X_sp, Y_sp, Z_sp = self.dict_punkt_symbole[standpunkt]
X_zp, Y_zp, Z_zp = self.dict_punkt_symbole[zielpunkt]
B_sp, L_sp = sp.symbols(f"B{standpunkt} L{standpunkt}")
B_zp, L_zp = sp.symbols(f"B{zielpunkt} L{zielpunkt}")
if beobachtungsart == "geometrisches_nivellement":
d_r_dX_zp = sp.cos(B_zp) * sp.cos(L_zp)
d_r_dX_sp = -sp.cos(B_sp) * sp.cos(L_sp)
d_r_dY_zp = sp.cos(B_zp) * sp.sin(L_zp)
d_r_dY_sp = -sp.cos(B_sp) * sp.sin(L_sp)
d_r_dZ_zp = sp.sin(B_zp)
d_r_dZ_sp = -sp.sin(B_sp)
zeile_A_Matrix = []
for punkt in liste_punktnummern:
if punkt == standpunkt:
zeile_A_Matrix.extend([d_r_dX_sp, d_r_dY_sp, d_r_dZ_sp])
elif punkt == zielpunkt:
zeile_A_Matrix.extend([d_r_dX_zp, d_r_dY_zp, d_r_dZ_zp])
else:
zeile_A_Matrix.extend([0, 0, 0])
for orientierung in liste_orientierungsunbekannte:
zeile_A_Matrix.append(0)
liste_A_nivellement_zeilen.append(zeile_A_Matrix)
liste_zeilenbeschriftungen_nivellement.append(
f"{beobachtungenID}_niv_{standpunkt}_{zielpunkt}"
)
if liste_beobachtungsgleichungen_distanz:
f_matrix_dist = sp.Matrix(liste_beobachtungsgleichungen_distanz)
unbekanntenvektor = sp.Matrix(liste_unbekannte)
A_dist = f_matrix_dist.jacobian(unbekanntenvektor)
else:
A_dist = None
if liste_A_richtung_zeilen:
A_richtung = sp.Matrix(liste_A_richtung_zeilen)
else:
A_richtung = None
if liste_A_zenitwinkel_zeilen:
A_zenitwinkel = sp.Matrix(liste_A_zenitwinkel_zeilen)
else:
A_zenitwinkel = None
if liste_beobachtungsgleichungen_gnssbasislinien:
f_matrix_gnss = sp.Matrix(liste_beobachtungsgleichungen_gnssbasislinien)
unbekanntenvektor = sp.Matrix(liste_unbekannte)
A_gnssbasislinien = f_matrix_gnss.jacobian(unbekanntenvektor)
else:
A_gnssbasislinien = None
if liste_A_nivellement_zeilen:
#f_matrix_nivellement = sp.Matrix(liste_beobachtungsgleichungen_nivellement)
#unbekanntenvektor = sp.Matrix(liste_unbekannte)
#A_nivellement = f_matrix_nivellement.jacobian(unbekanntenvektor)
A_nivellement = sp.Matrix(liste_A_nivellement_zeilen)
else:
A_nivellement = None
A_gesamt = None
liste_zeilenbeschriftungen_gesamt = []
if A_dist is not None:
A_gesamt = A_dist
liste_zeilenbeschriftungen_gesamt.extend(liste_zeilenbeschriftungen_distanz)
if A_richtung is not None:
if A_gesamt is None:
A_gesamt = A_richtung
else:
A_gesamt = A_gesamt.col_join(A_richtung)
liste_zeilenbeschriftungen_gesamt.extend(liste_zeilenbeschriftungen_richtung)
if A_zenitwinkel is not None:
if A_gesamt is None:
A_gesamt = A_zenitwinkel
else:
A_gesamt = A_gesamt.col_join(A_zenitwinkel)
liste_zeilenbeschriftungen_gesamt.extend(liste_zeilenbeschriftungen_zenitwinkel)
if A_gnssbasislinien is not None:
if A_gesamt is None:
A_gesamt = A_gnssbasislinien
else:
A_gesamt = A_gesamt.col_join(A_gnssbasislinien)
liste_zeilenbeschriftungen_gesamt.extend(liste_zeilenbeschriftungen_gnssbasislinien)
if A_nivellement is not None:
if A_gesamt is None:
A_gesamt = A_nivellement
else:
A_gesamt = A_gesamt.col_join(A_nivellement)
liste_zeilenbeschriftungen_gesamt.extend(liste_zeilenbeschriftungen_nivellement)
if A_gesamt is None:
return None
if datumsfestlegung == "weiche Lagerung":
vertauschung = list(range(len(liste_unbekannte)))
if liste_unbekannte_datumsfestlegung is not None and liste_unbekannte_datumsfestlegung != []:
liste_unbekannte_alt = list(liste_unbekannte)
liste_unbekannte_datumsfestlegung = [str(u).strip() for u in liste_unbekannte_datumsfestlegung]
idx_rechts = []
for name in liste_unbekannte_datumsfestlegung:
for i, sym in enumerate(liste_unbekannte_alt):
if str(sym) == name and i not in idx_rechts:
idx_rechts.append(i)
break
idx_links = [i for i in range(len(liste_unbekannte_alt)) if i not in idx_rechts]
vertauschung = idx_links + idx_rechts
A_gesamt = A_gesamt[:, vertauschung]
liste_unbekannte = [liste_unbekannte_alt[i] for i in vertauschung]
# Zusatzgeleichungen der weichen Lagerung
anzhl_einheitsmatrix = len(liste_unbekannte_datumsfestlegung)
if anzhl_einheitsmatrix > 0:
nullenmatrix = sp.zeros(anzhl_einheitsmatrix, A_gesamt.shape[1] - anzhl_einheitsmatrix)
einheitsmatrix = sp.eye(anzhl_einheitsmatrix)
A_weiche_Lagerung = nullenmatrix.row_join(einheitsmatrix)
A_gesamt = A_gesamt.col_join(A_weiche_Lagerung)
for unbekannte_datumsfestlegung in liste_unbekannte_datumsfestlegung:
liste_zeilenbeschriftungen_gesamt.append(f"lA_{unbekannte_datumsfestlegung}")
self.liste_unbekanntenvektor_symbolisch = liste_unbekannte
Export.matrix_to_csv(r"Zwischenergebnisse\Jacobi_Matrix_Symbolisch.csv", liste_unbekannte,
liste_zeilenbeschriftungen_gesamt, A_gesamt, "Beobachtung")
return A_gesamt, liste_unbekannte, liste_zeilenbeschriftungen_gesamt
def jacobi_matrix_zahlen_iteration_0(self, A_symbolisch: sp.Matrix, koordinatenart: str,
liste_unbekannte: list = None,
liste_zeilenbeschriftungen_gesamt: list = None,
iterationsnummer: int = 0) -> ndarray[tuple[Any, ...], dtype[Any]] | None:
self.liste_beobachtungsvektor_symbolisch = [str(x) for x in liste_zeilenbeschriftungen_gesamt]
if koordinatenart == "naeherung_us":
#A_numerisch = A_symbolisch.xreplace(self.substitutionen_dict)
if self.func_A0 is None:
#self.liste_symbole_lambdify = sorted(self.substitutionen_dict.keys(), key=lambda s: str(s))
self.func_A0 = sp.lambdify(
self.liste_symbole_lambdify,
A_symbolisch,
modules="numpy",
cse=True
)
fehlend = [s for s in self.liste_symbole_lambdify if s not in self.substitutionen_dict]
if fehlend:
Export.matrix_to_csv(
r"Zwischenergebnisse\fehlende_substitutionen_A.csv",
[""],
[str(s) for s in fehlend],
sp.Matrix([[str(s)] for s in fehlend]),
"fehlend"
)
raise ValueError(f"Fehlende Substitutionen in A: {[str(s) for s in fehlend[:30]]}")
liste_werte = [self.substitutionen_dict[s] for s in self.liste_symbole_lambdify]
#A_numerisch = sp.Matrix(self.func_A0(*liste_werte))
A_numerisch = np.asarray(self.func_A0(*liste_werte), dtype=float)
Export.matrix_to_csv(fr"Zwischenergebnisse\{iterationsnummer}Jacobi_Matrix_Numerisch_Iteration0.csv", liste_unbekannte,
liste_zeilenbeschriftungen_gesamt, A_numerisch, "Beobachtung")
condA = float(np.linalg.cond(A_numerisch))
rankA = int(np.linalg.matrix_rank(A_numerisch))
Export.matrix_to_csv(
fr"Zwischenergebnisse\{iterationsnummer}_Jacobi_Matrix_Stats.csv",
[""],
["condA", "rankA"],
np.array([[condA], [rankA]], dtype=float),
"Wert"
)
return A_numerisch
else:
print("Koordinaten noch nicht implementiert!")
def beobachtungsvektor_numerisch(self, liste_beobachtungsvektor_symbolisch: list) -> MutableDenseMatrix:
liste_beobachtungsvektor_numerisch = []
for beobachtung_symbolisch in liste_beobachtungsvektor_symbolisch:
beobachtung_symbolisch = str(beobachtung_symbolisch).strip()
if beobachtung_symbolisch.startswith("lA_"):
beobachtung_symbolisch = str(beobachtung_symbolisch.split("_", 1)[1]).strip()
liste_beobachtungsvektor_numerisch.append(self.substitutionen_dict[sp.Symbol(beobachtung_symbolisch)])
beobachtungsvektor_numerisch = sp.Matrix(liste_beobachtungsvektor_numerisch)
Export.matrix_to_csv(r"Zwischenergebnisse\Beobachtungsvektor_Numerisch.csv", [""], liste_beobachtungsvektor_symbolisch, beobachtungsvektor_numerisch, "Beobachtungsvektor")
return beobachtungsvektor_numerisch
def beobachtungsvektor_naeherung_symbolisch(self, liste_beobachtungsvektor_symbolisch: list) -> sp.Matrix:
liste_beobachtungsgleichungen = []
self.dict_punkt_symbole = {}
liste_punktnummern = []
for beobachtung_symbolisch in liste_beobachtungsvektor_symbolisch:
aufgeteilt = str(beobachtung_symbolisch).strip().split("_")
if aufgeteilt[0] == "lA":
continue
if aufgeteilt [1] == "SD" or aufgeteilt [1] == "R" or aufgeteilt [1] == "ZW":
standpunkt = str(aufgeteilt[3])
zielpunkt = str(aufgeteilt[4])
if standpunkt not in liste_punktnummern:
liste_punktnummern.append(standpunkt)
if zielpunkt not in liste_punktnummern:
liste_punktnummern.append(zielpunkt)
if aufgeteilt[1] == "gnssbx" or aufgeteilt[1] == "gnssby" or aufgeteilt[1] == "gnssbz":
standpunkt = str(aufgeteilt[2])
zielpunkt = str(aufgeteilt[3])
if standpunkt not in liste_punktnummern:
liste_punktnummern.append(standpunkt)
if zielpunkt not in liste_punktnummern:
liste_punktnummern.append(zielpunkt)
for punkt in liste_punktnummern:
X, Y, Z = sp.symbols(f"X{punkt} Y{punkt} Z{punkt}")
self.dict_punkt_symbole[str(punkt)] = (X, Y, Z)
for beobachtung_symbolisch in liste_beobachtungsvektor_symbolisch:
beobachtung_symbolisch = str(beobachtung_symbolisch).strip()
aufgeteilt = beobachtung_symbolisch.split("_")
if aufgeteilt[0] == "lA":
anschlusspunkt = str(aufgeteilt[1])
liste_beobachtungsgleichungen.append(sp.Symbol(anschlusspunkt))
continue
if aufgeteilt[1] == "SD" or aufgeteilt[1] == "R" or aufgeteilt[1] == "ZW":
#beobachtungen_ID = aufgeteilt[0]
beobachtungsart = aufgeteilt[1] # "SD", "R", "ZW"
beobachtungsgruppeID = aufgeteilt[2]
standpunkt = str(aufgeteilt[3]).strip()
zielpunkt = str(aufgeteilt[4]).strip()
X_sp, Y_sp, Z_sp = self.dict_punkt_symbole[standpunkt]
X_zp, Y_zp, Z_zp = self.dict_punkt_symbole[zielpunkt]
dX = X_zp - X_sp
dY = Y_zp - Y_sp
dZ = Z_zp - Z_sp
s = sp.sqrt(dX ** 2 + dY ** 2 + dZ ** 2) # Schrägstrecke
B_sp = sp.Symbol(f"B{standpunkt}")
L_sp = sp.Symbol(f"L{standpunkt}")
if beobachtungsart == "SD":
s_geom = sp.sqrt(dX ** 2 + dY ** 2 + dZ ** 2)
liste_beobachtungsgleichungen.append(s_geom)
elif beobachtungsart == "R":
#O_sp = sp.Symbol(f"O_{beobachtungsgruppeID}")
r_sp_zp = sp.Symbol(f"richtung_berechnet_{beobachtungsgruppeID}_{standpunkt}_{zielpunkt}")
# Lokales System: x_loc = Nord, y_loc = Ost
#x_loc = (-sp.sin(B_sp) * sp.cos(L_sp)) * dX + (-sp.sin(B_sp) * sp.sin(L_sp)) * dY + (sp.cos(B_sp)) * dZ
#y_loc = (-sp.sin(L_sp)) * dX + (sp.cos(L_sp)) * dY
#a12 = sp.atan2(y_loc, x_loc)
# Richtung nach Otepka: r = a12 - O
liste_beobachtungsgleichungen.append(r_sp_zp)
elif beobachtungsart == "ZW":
#dX = X_zp - X_sp
#dY = Y_zp - Y_sp
#dZ = Z_zp - Z_sp
#s_geom = sp.sqrt(dX ** 2 + dY ** 2 + dZ ** 2)
#z_loc = (sp.cos(B_sp) * sp.cos(L_sp)) * dX + (sp.cos(B_sp) * sp.sin(L_sp)) * dY + (sp.sin(B_sp)) * dZ
#zw = sp.acos(z_loc / s_geom)
zw_sp_zp = sp.Symbol(f"zw_berechnet_{beobachtungsgruppeID}_{standpunkt}_{zielpunkt}")
liste_beobachtungsgleichungen.append(zw_sp_zp)
if aufgeteilt[1] == "gnssbx" or aufgeteilt[1] == "gnssby" or aufgeteilt[1] == "gnssbz":
beobachtungsart = aufgeteilt[1]
standpunkt = str(aufgeteilt[2]).strip()
zielpunkt = str(aufgeteilt[3]).strip()
X_sp, Y_sp, Z_sp = self.dict_punkt_symbole[standpunkt]
X_zp, Y_zp, Z_zp = self.dict_punkt_symbole[zielpunkt]
dX = X_zp - X_sp
dY = Y_zp - Y_sp
dZ = Z_zp - Z_sp
if beobachtungsart == "gnssbx":
liste_beobachtungsgleichungen.append(dX)
if beobachtungsart == "gnssby":
liste_beobachtungsgleichungen.append(dY)
if beobachtungsart == "gnssbz":
liste_beobachtungsgleichungen.append(dZ)
if aufgeteilt[1] == "niv":
beobachtungsart = aufgeteilt[1]
standpunkt = str(aufgeteilt[2]).strip()
zielpunkt = str(aufgeteilt[3]).strip()
nh_sp = sp.Symbol(f"NH{standpunkt}")
nh_zp = sp.Symbol(f"NH{zielpunkt}")
niv_sp_zp = nh_zp - nh_sp
liste_beobachtungsgleichungen.append(niv_sp_zp)
beobachtungsvektor_naeherung_symbolisch = sp.Matrix(liste_beobachtungsgleichungen)
Export.matrix_to_csv(r"Zwischenergebnisse\Beobachtungsvektor_Näherung_Symbolisch.csv", [""],
liste_beobachtungsvektor_symbolisch, beobachtungsvektor_naeherung_symbolisch, "Beobachtungsvektor")
return beobachtungsvektor_naeherung_symbolisch
def beobachtungsvektor_naeherung_numerisch_iteration0(self, liste_beobachtungsvektor_symbolisch: list,
beobachtungsvektor_naeherung_symbolisch: sp.Matrix,
iterationsnummer: int = 0) -> ndarray[tuple[int, int], Any]:
#beobachtungsvektor_naeherung_numerisch_iteration0 = beobachtungsvektor_naeherung_symbolisch.xreplace(self.substitutionen_dict)
if self.func_beob0 is None:
#self.liste_symbole_lambdify = sorted(self.substitutionen_dict.keys(), key=lambda s: str(s))
self.func_beob0 = sp.lambdify(
self.liste_symbole_lambdify,
beobachtungsvektor_naeherung_symbolisch,
modules="numpy",
cse=True
)
liste_werte = [self.substitutionen_dict[s] for s in self.liste_symbole_lambdify]
#beobachtungsvektor_naeherung_numerisch_iteration0 = sp.Matrix(self.func_beob0(*liste_werte))
beobachtungsvektor_naeherung_numerisch_iteration0 = np.asarray(self.func_beob0(*liste_werte),
dtype=float).reshape(-1, 1)
Export.matrix_to_csv(fr"Zwischenergebnisse\{iterationsnummer}_Beobachtungsvektor_Näherung_Numerisch_Iteration0.csv", [""],
liste_beobachtungsvektor_symbolisch, beobachtungsvektor_naeherung_numerisch_iteration0,
"Beobachtungsvektor")
return beobachtungsvektor_naeherung_numerisch_iteration0
def unbekanntenvektor_symbolisch(self, liste_unbekannte: list) -> sp.Matrix:
unbekanntenvektor_symbolisch = sp.Matrix(liste_unbekannte)
Export.matrix_to_csv(r"Zwischenergebnisse\Unbekanntenvektor_Symbolisch.csv", [""], liste_unbekannte, unbekanntenvektor_symbolisch,
"Unbekanntenvektor")
return(unbekanntenvektor_symbolisch)
def unbekanntenvektor_numerisch(self, liste_unbekanntenvektor_symbolisch: list,
unbekanntenvektor_symbolisch: sp.Matrix,
dX_Vektor: np.Matrix = None,
unbekanntenvektor_neumerisch_vorherige_Iteration: np.Matrix = None,
iterationsnummer: int = 0) -> ndarray[tuple[int, int], Any] | ndarray[tuple[Any, ...], dtype[Any]]:
self.liste_unbekanntenvektor_symbolisch = liste_unbekanntenvektor_symbolisch
#if not hasattr(self, "liste_unbekanntenvektor_symbolisch"):
# self.liste_unbekanntenvektor_symbolisch = liste_unbekanntenvektor_symbolisch
if dX_Vektor is None and unbekanntenvektor_neumerisch_vorherige_Iteration is None:
#unbekanntenvektor_numerisch = unbekanntenvektor_symbolisch.xreplace(self.substitutionen_dict)
#if self.func_u0 is None:
# self.func_u0 = sp.lambdify(
# self.liste_symbole_lambdify,
# unbekanntenvektor_symbolisch,
# modules="numpy",
# cse=True
# )
unbekanntenvektor_numerisch = np.asarray(
[[float(self.substitutionen_dict[sym])] for sym in self.liste_unbekanntenvektor_symbolisch],
dtype=float
).reshape(-1, 1)
#liste_werte = [self.substitutionen_dict[s] for s in self.liste_symbole_lambdify]
#unbekanntenvektor_numerisch = sp.Matrix(self.func_u0(*liste_werte))
#unbekanntenvektor_numerisch = np.asarray(self.func_u0(*liste_werte), dtype=float).reshape(-1, 1)
else:
#unbekanntenvektor_numerisch = unbekanntenvektor_neumerisch_vorherige_Iteration + dX_Vektor
unbekanntenvektor_neumerisch_vorherige_Iteration = np.asarray(
unbekanntenvektor_neumerisch_vorherige_Iteration, dtype=float).reshape(-1, 1)
dX_Vektor = np.asarray(dX_Vektor, dtype=float).reshape(-1, 1)
unbekanntenvektor_numerisch = unbekanntenvektor_neumerisch_vorherige_Iteration + dX_Vektor
self.substitutionen_dict = self.dict_substitutionen_uebergeordnetes_system(unbekanntenvektor_numerisch)
Export.matrix_to_csv(fr"Zwischenergebnisse\{iterationsnummer}_Unbekanntenvektor_Numerisch_Iteration0.csv", [""],
liste_unbekanntenvektor_symbolisch, unbekanntenvektor_numerisch,
"Unbekanntenvektor")
return unbekanntenvektor_numerisch
def unbekanntenvektor_numerisch_to_dict_unbekanntenvektor(self, liste_unbekanntenvektor_symbolisch: list, unbekanntenvektor_numerisch: np.Matrix) -> dict:
dict_unbekanntenvektor_numerisch = {}
#index = 0
unbekanntenvektor_numerisch = np.asarray(unbekanntenvektor_numerisch, dtype=float).reshape(-1, 1)
idx = {str(sym): i for i, sym in enumerate(liste_unbekanntenvektor_symbolisch)}
punktnummern = []
for sym in liste_unbekanntenvektor_symbolisch:
name = str(sym)
if name.startswith("X"):
pn = name[1:]
if pn not in punktnummern:
punktnummern.append(pn)
dict_koordinaten = {}
for pn in punktnummern:
iX = idx.get(f"X{pn}", None)
iY = idx.get(f"Y{pn}", None)
iZ = idx.get(f"Z{pn}", None)
if iX is None or iY is None or iZ is None:
continue
dict_koordinaten[pn] = sp.Matrix([
float(unbekanntenvektor_numerisch[iX, 0]),
float(unbekanntenvektor_numerisch[iY, 0]),
float(unbekanntenvektor_numerisch[iZ, 0]),
])
#dict_unbekanntenvektor_numerisch[punktnummer] = sp.Matrix([
# float(unbekanntenvektor_numerisch[index, 0]),
# float(unbekanntenvektor_numerisch[index + 1, 0]),
# float(unbekanntenvektor_numerisch[index + 2, 0])
#])
#index += 3
return dict_koordinaten
def berechnung_dl(self, beobachtungsvektor_numerisch: np.Matrix, beobachtungsvektor_naeherung_numerisch: sp.Matrix,
liste_beobachtungsvektor_symbolisch: list = None, iterationsnummer: int = 0) -> np.Matrix:
dl = beobachtungsvektor_numerisch - beobachtungsvektor_naeherung_numerisch
dl = np.asarray(dl, dtype=float)
if liste_beobachtungsvektor_symbolisch is None:
liste_beobachtungsvektor_symbolisch = self.liste_beobachtungsvektor_symbolisch
for i, name in enumerate(liste_beobachtungsvektor_symbolisch):
if "_R_" in str(name):
dl[i] = np.arctan2(np.sin(dl[i]), np.cos(dl[i]))
Export.matrix_to_csv(
fr"Zwischenergebnisse\{iterationsnummer}_dl.csv",
[""],
liste_beobachtungsvektor_symbolisch,
dl.reshape(-1, 1),
"dl"
)
return dl
def dict_substitutionen_uebergeordnetes_system(self,
unbekanntenvektor_aus_iteration: np.Matrix = None) -> dict[Any, Any]:
db_zugriff = Datenbankzugriff(self.pfad_datenbank)
berechnungen = Berechnungen(self.a, self.b)
if unbekanntenvektor_aus_iteration is None:
dict_koordinaten = db_zugriff.get_koordinaten("naeherung_us")
else:
dict_koordinaten = self.unbekanntenvektor_numerisch_to_dict_unbekanntenvektor(
self.liste_unbekanntenvektor_symbolisch,
unbekanntenvektor_aus_iteration
)
for punktnummer, matrix in dict_koordinaten.items():
dict_koordinaten[punktnummer] = [float(matrix[0]), float(matrix[1]), float(matrix[2])]
liste_beobachtungen_tachymeter = db_zugriff.get_beobachtungen_from_beobachtungenid()
liste_beobachtungen_gnssbasislinien = db_zugriff.get_beobachtungen_gnssbasislinien()
liste_beobachtungen_nivellemente = db_zugriff.get_beobachtungen_nivellement()
liste_azimut_richtungen, dict_orientierungen = berechnungen.berechnung_richtung_azimut_zenitwinkel(self.pfad_datenbank, dict_koordinaten)
dict_koordinaten_xyz_kopie = {pn: [v[0], v[1], v[2]] for pn, v in dict_koordinaten.items()}
dict_koordinaten_B_L = berechnungen.geometrische_breite_laenge(dict_koordinaten_xyz_kopie)
punktnummern_niv = set()
for beobachtungenID, pn_sp, pn_zp, niv_dh, niv_strecke, niv_anz_standpkte in liste_beobachtungen_nivellemente:
punktnummern_niv.add(str(pn_sp).strip())
punktnummern_niv.add(str(pn_zp).strip())
dict_koordinaten_niv = {}
for pn in punktnummern_niv:
if pn in dict_koordinaten:
dict_koordinaten_niv[pn] = dict_koordinaten[pn]
dict_koordinaten_utm = self.trafos.ecef_to_utm(
dict_koordinaten_niv,
self.pfad_tif_quasigeoidundolation)
substitutionen = {}
for punktnummer, vektor in dict_koordinaten_B_L.items():
X_sym, Y_sym, Z_sym, B_sym, L_Sym = sp.symbols(
f"X{punktnummer} Y{punktnummer} Z{punktnummer} B{punktnummer} L{punktnummer}")
substitutionen[X_sym] = float(vektor[0][0])
substitutionen[Y_sym] = float(vektor[0][1])
substitutionen[Z_sym] = float(vektor[0][2])
substitutionen[B_sym] = float(vektor[1])
substitutionen[L_Sym] = float(vektor[2])
#for beobachtungsgruppeID, orientierung in dict_orientierungen.items():
# O_sym = sp.symbols(f"O_{beobachtungsgruppeID}")
# substitutionen[O_sym] = float(orientierung)
for beobachtungsgruppeID, standpunkt, zielpunkt, azimut, richtung, zenitwinkel, schraegstrecke, orientierung in liste_azimut_richtungen:
richtung_sym = sp.symbols(f"richtung_berechnet_{beobachtungsgruppeID}_{standpunkt}_{zielpunkt}")
substitutionen[richtung_sym] = float(richtung)
azimut_sym = sp.symbols(f"azimut_berechnet_{beobachtungsgruppeID}_{standpunkt}_{zielpunkt}")
substitutionen[azimut_sym] = float(azimut)
zenitwinkel_sym = sp.symbols(f"zw_berechnet_{beobachtungsgruppeID}_{standpunkt}_{zielpunkt}")
substitutionen[zenitwinkel_sym] = float(zenitwinkel)
schraegstrecke_sym = sp.symbols(f"strecke_berechnet_{beobachtungsgruppeID}_{standpunkt}_{zielpunkt}")
substitutionen[schraegstrecke_sym] = float(schraegstrecke)
for punktnummer, koordinaten_utm in dict_koordinaten_utm.items():
normalhoehe_sym = sp.symbols(f"NH{punktnummer}")
substitutionen[normalhoehe_sym] = float(koordinaten_utm[2])
for standpunkt, zielpunkt, beobachtungenID, beobachtungsgruppeID, tachymeter_richtung, tachymeter_zenitwinkel, tachymeter_distanz in liste_beobachtungen_tachymeter:
alpha = sp.symbols(f"{beobachtungenID}_R_{beobachtungsgruppeID}_{standpunkt}_{zielpunkt}")
zw = sp.symbols(f"{beobachtungenID}_ZW_{beobachtungsgruppeID}_{standpunkt}_{zielpunkt}")
s = sp.symbols(f"{beobachtungenID}_SD_{beobachtungsgruppeID}_{standpunkt}_{zielpunkt}")
if tachymeter_richtung is None and tachymeter_zenitwinkel is None and tachymeter_distanz is None:
continue
substitutionen[alpha] = float(tachymeter_richtung)
substitutionen[zw] = float(tachymeter_zenitwinkel)
substitutionen[s] = float(tachymeter_distanz)
#substitutionen[sp.Symbol(f"O{beobachtungsgruppeID}")] = 0.0
for beobachtungenID, punktnummer_sp, punktnummer_zp, gnss_bx, gnss_by, gnss_bz, gnss_s0, gnss_cxx, gnss_cxy, gnss_cxz, gnss_cyy, gnss_cyz, gnss_czz in liste_beobachtungen_gnssbasislinien:
beobachtungenID = str(beobachtungenID).strip()
punktnummer_sp = str(punktnummer_sp).strip()
punktnummer_zp = str(punktnummer_zp).strip()
bx = sp.symbols(f"{beobachtungenID}_gnssbx_{punktnummer_sp}_{punktnummer_zp}")
by = sp.symbols(f"{beobachtungenID}_gnssby_{punktnummer_sp}_{punktnummer_zp}")
bz = sp.symbols(f"{beobachtungenID}_gnssbz_{punktnummer_sp}_{punktnummer_zp}")
if gnss_bx is None and gnss_by is None and gnss_bz is None:
continue
substitutionen[bx] = float(gnss_bx)
substitutionen[by] = float(gnss_by)
substitutionen[bz] = float(gnss_bz)
for beobachtungenID, punktnummer_sp, punktnummer_zp, niv_dh, niv_strecke, niv_anz_standpkte in liste_beobachtungen_nivellemente:
beobachtungenID = str(beobachtungenID).strip()
punktnummer_sp = str(punktnummer_sp).strip()
punktnummer_zp = str(punktnummer_zp).strip()
niv = sp.symbols(f"{beobachtungenID}_niv_{punktnummer_sp}_{punktnummer_zp}")
if niv_dh is None:
continue
substitutionen[niv] = float(niv_dh)
if unbekanntenvektor_aus_iteration is not None:
dict_O = self.unbekanntenvektor_numerisch_to_dict_orientierungen(
self.liste_unbekanntenvektor_symbolisch,
unbekanntenvektor_aus_iteration
)
for orientierungs_id, wert in dict_O.items():
substitutionen[sp.Symbol(f"O{orientierungs_id}")] = float(wert)
else:
#for standpunkt, zielpunkt, beobachtungenID, beobachtungsgruppeID, *_ in liste_beobachtungen_tachymeter:
# O_sym = sp.Symbol(f"O{beobachtungsgruppeID}")
# if O_sym not in substitutionen:
# substitutionen[O_sym] = 0
for beobachtungsgruppeID, standpunkt, zielpunkt, azimut, richtung, zenitwinkel, schraegstrecke, orientierung in liste_azimut_richtungen:
O_sym = sp.Symbol(f"O{beobachtungsgruppeID}")
if O_sym not in substitutionen:
substitutionen[O_sym] = orientierung
return substitutionen
def unbekanntenvektor_numerisch_to_dict_orientierungen(self, liste_unbekanntenvektor_symbolisch: list,
unbekanntenvektor_numerisch: np.Matrix) -> dict[Any, Any]:
dict_O = {}
unbekanntenvektor_numerisch = np.asarray(unbekanntenvektor_numerisch, dtype=float).reshape(-1, 1)
for i, symbol in enumerate(liste_unbekanntenvektor_symbolisch):
name = str(symbol)
if name.startswith("O"):
orientierungs_id = name[1:]
dict_O[orientierungs_id] = float(unbekanntenvektor_numerisch[i, 0])
return dict_O