"""
#------------------------------------------------------------------------------
# Name:        Seiloptimierungstool
# Purpose:
#
# Author:      Patricia Moll
#
# Created:     14.05.2013
# Copyright:   (c) mollpa 2012
# Licence:     <your licence>
#------------------------------------------------------------------------------
"""

from math import cos, inf, pi, sin

import numpy as np


def xfrange(start, stop, step):
    while start < stop:
        yield start
        start += step


def preciseCable(IS, poles, STA):
    """Berechnet die exakte Seilmechanik für ein vorgegebenes Seilsystem:
    Seilstatik (gemäss dem Paper von Zweifel 1960)
    Fall: Stützen unbeweglich, Seile auf Zwischenstützen gleitend
         Inputdaten:
    Felddimensionen:
    b: Breite der Seilfelder [m] Bsp.: [98.2700, 184.8300]
    h: Höhe der Seilfelder   [m] Bsp.: [56.211, 55.22]
    HM: Höhe der Stützen     [m] Bsp.: [12, 12, 12]
    IS: Eigenschaften des Seilsystems
    AnkerIncluded: Ist der Anker bereits in den Seilfeldern inbegriffen?
    (True: ja, False: nein)
         Output:
    Der Output erfolgt als array in EigSeillinie
    Bem.: Die Ausgabe der x Koordinaten hat als Ursprung den Anfangspunkt
    (Auf! der ersten Stütze)
    .laengenprofil_y: y-Koordinaten des Laengenprofils
    .x_leer:          x-Koordinaten des Leerseils
    .x_last_zweifel   x-Koordinaten des Lastseils nach Zweifel (Zw)
    .x_last_pestal    x-Koordinaten des Lastseils nach Pestal (Pe)
    .Kraft_Mast       Einwirkende Kraft auf den Masten
    .STfm_Last        Seilzugkraft in Feldmitte unter Last (Zw)
    .y_LeerSeil       Durchhang in Feldmitte des Leerseils
    .d_ST_out         Erhöhung der Seilzugkraft bei Belastung
    .Hm   Horizontalkomponente der Seilzugkraft in Feldmitte unter Last (Zw)
    .Hs   Horizontalkomp. d. Seilzugkraft an den Stützen (Leerseilverhältn.)
    .ym   Durchhang in Feldmitte unter Last (Zw)"""

    Q = IS["Q"]  # [kN]
    qT = IS["qT"]  # [kN/m]
    F = IS["A"]  # [mm^2]
    E = IS["E"]  # [kN/mm^2]
    # beta = IS["beta"]     # [1/°C] Ausdehnungskoeffizient von Stahl
    Federkonstante = inf  # Federkonstante der Verankerung
    qz1 = IS["qz1"]  # [kN/m] Zugseilgewicht links
    qz2 = IS["qz2"]  # [kN/m] Zugseilgewicht rechts

    # Field dimension
    b, h = poles.getCableFieldDimension()
    # Round the field dimensions do dm because this is going to be
    # the resolution of the cable calculation
    b = np.round(b, 1)
    h = np.round(h, 1)

    if any(b == 0):
        raise ValueError("At least one field has a length of zero")

    anzFelder = b.size
    anzStue = anzFelder + 1
    bquad = b**2
    seil_possible = True

    if np.where(b == 0)[0].size == 0:
        korrektur = 0
    else:
        korrektur = h[0]

    z = np.zeros(anzStue)
    laengenprofil_y = np.zeros(anzStue)
    z[0] = korrektur
    z[1:] = np.cumsum(h)
    zfm = z[:-1] + h * 0.5
    b_cum = np.cumsum(b)
    laengenprofil_y[1:] = b_cum

    # Berechnung der Sehne
    c = (bquad + h**2) ** 0.5
    c_sum = np.sum(c)
    # Mit Berücksichtigung des Zugseilgewichtes
    q_strich = qT + (qz1 + qz2) / 2

    # Leerseil
    # --------
    zfmqT = zfm * qT
    STfm = STA + zfmqT  # Seilzugkraft in Feldmitte
    HT = (STfm * b) / c  # Horizontalkomponente der Seilzugkraft
    # Überlänge der Lastlänge gegenüber dem Sehnenzug (Leerfeld)
    delta_s = ((bquad**2 * qT**2) / (24 * c * HT**2)) * (
        1 + (3 * bquad + 8 * h**2) / (240 * c**2) * (b * qT / HT) ** 2
    )

    delta_s_sum = np.sum(delta_s)
    UeberLaenge_Leerseil = delta_s_sum
    Laenge_Leerseil = c_sum + delta_s_sum
    Laenge_Leerseil_bei_0_kN = Laenge_Leerseil + Laenge_Leerseil * STA / (F * E)

    # Vollseil
    # --------
    # Iteration zur Berechnung des d_ST, welches unter dem Gewicht von Q
    # entsteht. Berechnung könnte noch beschleunigt werden, falls nicht alle
    # Felder sondern nur das betreffende berechnet werden!
    CableImpossible = False
    d_ST_out = np.zeros(anzFelder)  # Änderung der Seilzugkraft unter Last
    gen = 6  # Genauigkeit der Iteration: 6 sollte ausreichen
    basis = 2  # Hilfsgrösse für die Iteration
    for n in range(anzFelder):
        expon = -5  # -5 ist Standard
        d_ST = 0.0  # Änderung der Zugkraft
        d_Laenge = 1  # Startgrösse
        change_dir = False  # Hilfsgrösse für die Iteration
        za = 0  # Hilfsgrösse für die Iteration
        d_ST_alt = 0  # Startgrösse
        genauig = basis**-gen
        # Diese Prüfung führt in wenigen Einzelfällen zu anderen Resultaten
        # als das Matlab-Programm weil die Berechnung von d_Laenge nach
        # einigen Iterationen zu Abweichungen im Submilimeter-Bereich führen
        while (abs(d_Laenge) > genauig) and CableImpossible is False:
            STfm_neu = STfm + d_ST
            ym = c / (4 * STfm_neu) * (Q + c * qT / 2)  # 8a
            # Funktion mit Berücksichtigung des Zugseilgewichtes:
            # ym = c / (4 * STfm_neu) * (Q + c * q_strich/2) #(8)
            d_c = 2 * bquad / c**3 * ym**2  # 10a
            # 11a Leerseil
            delta_s_leer = (bquad * c * qT**2) / (24 * STfm_neu**2)
            delta_s_last = delta_s_leer * 0.25
            # delta s für die Belastung im betrachteten Abschnitt
            delta_s = np.copy(delta_s_leer)
            delta_s[n] = delta_s_last[n]
            UeberLaenge_Vollseil = d_c[n] + np.sum(delta_s)

            d_Laenge = (
                UeberLaenge_Vollseil
                - UeberLaenge_Leerseil
                - d_ST / (F * E) * Laenge_Leerseil
                - d_ST / (F * E) * poles.anchor["len"]
                - d_ST / Federkonstante
            )
            # print "Laenge ={}".format(d_Laenge)
            d_ST_out[n] = d_ST

            # Abfangen von falschen Berechnungen
            if (d_ST == 0) and (d_Laenge < 0):
                CableImpossible = True
            # Steuerung delta ST (Seilzugkraft)
            if d_Laenge > genauig:
                d_ST_alt = d_ST
                d_ST += basis**-expon
                za += 1
            else:
                expon += 1
                d_ST = d_ST_alt + basis**-expon
                change_dir = True
                za = 0
            if za == (basis - 1) and change_dir:
                expon += 1
                d_ST = d_ST_alt + basis**-expon
                za = 0

    # Eigenschaften des Seils
    # -----------------------
    # Seilzugkraft in Feldmitte unter Last [kN]
    STfm_Last = STfm + d_ST_out
    # Horizontalanteil der Seilzugkraft in Feldmitte unter Last [kN]
    Hm = b / c * STfm_Last
    # Durchbiegung in Feldmitte unter Last [m]
    ym = c / (4 * STfm_Last) * (Q + c * q_strich / 2)  # 8
    # Horizontalkräfte an den Stützen (Leerseilverhältnisse) [kN]
    Hs = b / c * STfm

    # Genaue Seil-Koordinaten bestimmen
    # ---------------------------------
    # Achtung, komplett anders implementier als in matlab!

    # Um genauere Seildaten zu erhalten, muss Schrittgrösse von defaul=1 auf
    # Submeter reduziert werden. Siehe matlab: DP = 0.01 m, also multiplikation mit 10
    Q_Null = 0
    multipl = 10
    # Round the decimeter values and then convert to int to use as index
    lenSeil = int(round(b_cum[-1] * multipl, 0)) + 1

    b1_base = np.zeros(lenSeil)  # 1 cm Schritte zwischen den Stützen
    # Leere Koordinaten Arrays initialisieren
    y_leer = np.copy(b1_base)
    H = np.copy(b1_base)
    y_last_zweifel = np.copy(b1_base)
    z_coord_leer = np.copy(b1_base)
    z_coord_zweifel = np.copy(b1_base)
    l_coord = np.copy(b1_base)

    start = 0
    for n in range(anzFelder):
        bn = b[n]
        Hmn = Hm[n]
        s_small = int(round(bn * multipl, 0))
        end = start + s_small + 1
        b1 = np.arange(0, end - start, 1) / 10
        b2 = bn - b1
        try:
            l_coord[start:end] = laengenprofil_y[n] + b1
        except ValueError as e:
            raise ValueError(f"Seilfeld nicht berechenbar - {e}") from e
        y_leer[start:end] = b1 * b2 / (bn * HT[n]) * (Q_Null + c[n] * qT / 2)
        # Interpolationsbeziehung für die Berechnung der Lastwegkurve
        H[start:end] = (
            Hmn * (1 - (1 - (Hs[n] / Hmn) ** 2) * (1 - 2 * b1 / bn) ** 2) ** 0.5
        )
        y_last_zweifel[start:end] = (
            ym[n] * Hmn / H[start:end] * (1 - (1 - 2 * b1 / bn) ** 2)
        )
        # Y-Koordinate (Länge) der Lastwegkurve unter Zweifel
        h_sehne = b1 / bn * h[n]
        # Leerseil
        z_coord_leer[start:end] = h_sehne - y_leer[start:end] + z[n]
        # X-Koordinate der Lastkurve unter Zweifel
        z_coord_zweifel[start:end] = h_sehne - y_last_zweifel[start:end] + z[n]

        start = end - 1

        # for b1 in xfrange(0, int(bn)+0.1, 0.1):
        # #for b1 in range(int(b[n])):
        #     di_ex.append(b1)
        #     b2 = bn-b1
        #     y_leer = b1*b2/(bn * HT[n]) * (Q_Null + c[n]*qT/2)
        #     # Interpolationsbeziehung für die Berechnung der Lastwegkurve
        #     H.append(Hmn * (1-(1-(Hs[n] / Hmn)**2) * (1-2*b1 / bn)**2)**0.5)
        #     y_last_zweifel = ym[n] * Hmn / H[-1] * (1-(1-2*b1/bn)**2)
        #     # Y-Koordinate (Länge) der Lastwegkurve unter Zweifel
        #     h_sehne = b1 / bn * h[n]
        #     # Leerseil
        #     z_coord_leer.append(h_sehne - y_leer + z[n])
        #     # X-Koordinate der Lastkurve unter Zweifel
        #     z_coord_zweifel.append(h_sehne - y_last_zweifel + z[n])
        #     l_coord.append(l[n]+b1)

    k = anzFelder - 1
    lastElement = h[k] + z[k]
    z_coord_leer[-1] = lastElement
    z_coord_zweifel[-1] = lastElement
    l_coord[-1] = laengenprofil_y[k] + b[k]
    H[-1] = Hm[k] * (1 - (1 - (Hs[k] / Hm[k]) ** 2) * (1 - 2 * b[k] / b[k]) ** 2) ** 0.5

    # Seilparameter berechnen
    # -----------------------

    kraft = {}
    # Suche der Indizes zur Berechnung der Seilwinkel
    idxStuetze = np.zeros(anzStue).astype("int")

    idxStuetze[1:] = b_cum * multipl
    idxNachher = idxStuetze[:-1] + 1
    idxVorher = idxStuetze[1:] - 1

    # Last nicht auf Stütze, sondern in unmittelbarer Nähe
    # Berechnung der Sattelkräfte

    # Berechnung der Lastseilwinkel wie bei Leerseil
    tg_phi_ob = (h / 2 + 2 * ym) / (b / 2)
    tg_phi_un = (h / 2 - 2 * ym) / (b / 2)

    # Winkel bei Last unmittelbar in der Nähe der Stütze
    phi = np.array([[np.nan] * anzStue] * 2)
    phi_ob = np.array([np.nan] * anzStue)
    phi_un = np.array([np.nan] * anzStue)
    phi_ob[1:] = np.arctan(tg_phi_ob)
    phi_un[:-1] = np.arctan(tg_phi_un)
    phi[0] = phi_ob / pi * 180
    phi[1] = phi_un / pi * 180
    kraft["Anlegewinkel_Lastseil"] = phi

    # Seilzugkräfte
    Hn = H[idxNachher]
    Hv = H[idxVorher]
    dln = l_coord[idxNachher + 1] - l_coord[idxNachher - 1]
    dxn = z_coord_zweifel[idxNachher + 1] - z_coord_zweifel[idxNachher - 1]
    sn = np.array([np.nan] * anzStue)
    sn[:-1] = Hn / dln * (dln**2 + dxn**2) ** 0.5

    dlv = l_coord[idxVorher + 1] - l_coord[idxVorher - 1]
    dxv = z_coord_zweifel[idxVorher + 1] - z_coord_zweifel[idxVorher - 1]
    sv = np.array([np.nan] * anzStue)
    sv[1:] = Hv / dlv * (dlv**2 + dxv**2) ** 0.5
    # Werte richtig anordnen (links/rechts)
    seilzugkr = np.array([np.nan] * anzStue * 2)
    seilzugkr[::2] = sv
    seilzugkr[1::2] = sn
    # Seilzugkraft unmittelbar bei Stütze
    kraft["Seilzugkraft_beiStuetze"] = seilzugkr

    # Untenstehend wird der Nachweis für die Sattelkraft und die Änderung
    # des Winkels berechnet: (gemäss Skript von Oplatka)
    # Muss noch verifiziert werden
    y_LeerSeil = b / (4 * HT) * (c * qT / 2)  # Leerseildurchhang

    tg_phi_o = (h / 2 + 2 * y_LeerSeil) / (b / 2)
    tg_phi_u = (h / 2 - 2 * y_LeerSeil) / (b / 2)
    phi_o = np.array([np.nan] * anzStue)
    phi_u = np.array([np.nan] * anzStue)
    phi_o[1:] = np.arctan(tg_phi_o)
    phi_u[:-1] = np.arctan(tg_phi_u)

    dAnkerA = poles.anchor["field"][0]
    zAnkerA = poles.anchor["field"][1]
    dAnkerE = poles.anchor["field"][2]
    zAnkerE = poles.anchor["field"][3]

    if dAnkerA != 0:
        phi_oA = np.arctan(zAnkerA / dAnkerA)
    else:
        phi_oA = np.nan
    if dAnkerE != 0:
        phi_oE = np.arctan((-1 * zAnkerE) / dAnkerE)
    else:
        phi_oE = np.nan
    oldsettings = np.geterr()
    j = np.seterr(all="ignore")

    phi_o[0] = phi_oA
    phi_u[-1] = phi_oE
    phi_o = phi_o
    phi_u = phi_u
    phi_leer = np.array([[np.nan] * anzStue] * 2)
    phi_leer[0] = phi_o / pi * 180
    phi_leer[1] = phi_u / pi * 180
    kraft["Anlegewinkel_Leerseil"] = phi_leer
    phi_leer_knick = (phi_o - phi_u) / pi * 180
    kraft["Leerseilknickwinkel"] = phi_leer_knick

    # Berechnung Lastseilknickwinkel mit der Last rechts bzw. links direkt neben der Stütze
    # Hinzufügen des Leerseilwinkels vor dem Anfangspunkt resp. nach dem Endpunkt
    phi_ob[0] = phi_oA
    phi_un[-1] = phi_oE

    # Der grössere Knickwinkel wird als Lastseilknickwinkel ausgegeben
    phi_last_knick_ob = (phi_ob - phi_u) / pi * 180
    phi_last_knick_un = (phi_o - phi_un) / pi * 180
    phi_last_knick = np.fmax(phi_last_knick_un, phi_last_knick_ob)
    kraft["Lastseilknickwinkel"] = phi_last_knick

    ST = STA + z * qT
    Vi = ST * (np.sin(phi_o) - np.sin(phi_u))  # Vi darf nicht negativ sein
    # Nachweis basierend auf Leerseilknickwinkel
    kraft["Nachweis"] = np.where(phi_leer_knick >= IS["LeerKnickMit"], ["Ja"], ["Nein"])
    # Do not test for np.nan values (crane, anchor)
    kraft["Nachweis"][np.isnan(phi_leer_knick)] = "-"
    if "Nein" in kraft["Nachweis"][1:-1]:  # Test für die Zwischenstützen
        seil_possible = False

    Hi = ST * (np.cos(phi_o) - np.cos(phi_u))
    kraft["Sattelkraft_ausSeil"] = np.array([(Vi**2 + Hi**2) ** 0.5, Vi, Hi])

    # Tragseilkraft
    # Ueberlange Lastlänge gegenüber dem Sehnenzug (Leerfeld):
    Gew_Seil = (c + delta_s) * qT
    KraftTragseil = np.append(Gew_Seil, 0) / 2 + np.append(0, Gew_Seil) / 2

    # Zugseilkraft
    c_cum = np.append(0, np.cumsum(c))
    Gew_Seil_li = c_cum * IS["qz1"]
    Gew_Seil_re = (c_cum[-1] - c_cum) * IS["qz2"]
    KraftZugseil = (Gew_Seil_li + Gew_Seil_re) / 2

    sattelkrTotV = Vi + IS["Q"] + KraftZugseil + KraftTragseil
    kraft["Sattelkraft_Total"] = np.array(
        [(sattelkrTotV**2 + Hi**2) ** 0.5, sattelkrTotV, Hi]
    )

    kraft["UebrigeKraft_befStuetze"] = IS["Q"] + KraftZugseil + KraftTragseil

    # Sattelkraft für nicht befahrbare Stuetzen
    ###########################################
    # Kraft in unmittelbarer Nähe der Stützen
    # 1. Fall: Lastseil auf der linken Seite
    Vi_nbS_L = sv * (np.sin(phi_ob) - np.sin(phi_u))
    Hi_nbS_L = sv * (np.cos(phi_ob) - np.cos(phi_u))

    # Wird das benötigt?
    # kraft['Sattelkraft_beiStuetze_L'] = np.concatenate((
    #     (Vi_nbS_L**2 + Hi_nbS_L**2)**0.5, Vi_nbS_L, Hi_nbS_L))

    # 2. Fall: Lastseil auf der rechten Seite
    Vi_nbS_R = sn * (np.sin(phi_o) - np.sin(phi_un))
    Hi_nbS_R = sn * (np.cos(phi_o) - np.cos(phi_un))

    # Wird das benötigt?
    # kraft['Sattelkraft_beiStuetze_R'] = np.concatenate((
    #     (Vi_nbS_R**2 + Hi_nbS_R**2)**0.5, Vi_nbS_R, Hi_nbS_R))

    # Sattelkraft bei Stützen berechnen

    # sSt = np.array([[0]*anzSt*2]*3)
    # sSt[0][:anzSt] = (Vi_nbS_L**2 + Hi_nbS_L**2)**0.5
    # sSt[0][anzSt:] = (Vi_nbS_R**2 + Hi_nbS_R**2)**0.5
    # sSt[1][:anzSt] = Vi_nbS_L
    # sSt[1][anzSt:] = Vi_nbS_R
    # sSt[2][:anzSt] = Hi_nbS_L
    # sSt[2][anzSt:] = Hi_nbS_R

    sSt = np.array(
        [
            np.array(
                [
                    (Vi_nbS_L**2 + Hi_nbS_L**2) ** 0.5,
                    (Vi_nbS_R**2 + Hi_nbS_R**2) ** 0.5,
                ]
            ),
            np.array([Vi_nbS_L, Vi_nbS_R]),
            np.array([Hi_nbS_L, Hi_nbS_R]),
        ]
    )
    # richtige Reihenfolge herstellen (immer links, rechts, links, ...)
    sattelkraft = []
    kompCount = 0
    while kompCount < 3:  # für resultierende, horizontale und vertikale Komponente
        werte = [0] * anzStue * 2
        j = 0
        for i in range(anzStue):  # es wird die i-te Stütze angesprochen
            komp = sSt[kompCount]
            werte[j : j + 2] = [
                komp[0][i],
                komp[1][i],
            ]  # 0 = links, 1 = rechts
            j += 2
        sattelkraft.append(werte)
        kompCount += 1
    kraft["Sattelkraft_beiStuetze"] = np.array(sattelkraft)

    # Max auftretende Seilzugkraft
    STfm_Last_max = np.max(STfm_Last)
    feld_max = np.argmax(STfm_Last)

    # Auftretende Maximalseilzugkraft...
    kraft["MaxSeilzugkraft_L"] = np.array([np.nan] * 5)
    #   am höchsten Punkt
    # TODO Leo: Die hier berechnete max. Seilzugkraft (siehe Kurzbericht) berücksichtigt die toten Anker nicht --> anpassen?
    # Comment Leo (24.8.2020), Nein so lassen ist OK, genau genommen, wird die Kraft auf Sattelhöhe ausgegeben, aber
    # welche Kraft dann am Anker wirkt ist unrealistisch zu berechnen (mehrere Ankerseile)
    kraft["MaxSeilzugkraft_L"][0] = STfm_Last_max + (np.max(z) - zfm[feld_max]) * qT
    #   am Startpunkt der Optimierung (1. Stütze oder befahrbare Verankerung --> pole oder anchor_pole)
    kraft["MaxSeilzugkraft_L"][1] = STfm_Last_max + (z[0] - zfm[feld_max]) * qT
    #   am Endpunkt der Optimierung (letzte Stütze oder befahrbare Verankerung --> pole oder anchor_pole)
    kraft["MaxSeilzugkraft_L"][2] = STfm_Last_max + (z[-1] - zfm[feld_max]) * qT

    # Falls tote Anker vorhanden sind
    # TODO Leo: folgende Zeilen berechnen TmaxA/E für tote Anker, bitte überprüfen
    # genau genommen, wird die Kraft auf Sattelhöhe ausgegeben, aber
    # für die Praxis ist das genügend genau
    if poles.hasAnchorA:
        kraft["MaxSeilzugkraft_L"][3] = STfm_Last_max + (z[0] - zfm[feld_max]) * qT
    if poles.hasAnchorE:
        kraft["MaxSeilzugkraft_L"][4] = STfm_Last_max + (z[-1] - zfm[feld_max]) * qT

    # Auftretende Maximalseilzugkraft...
    kraft["MaxSeilzugkraft"] = np.array([[np.nan] * anzFelder] * 3)
    #   in Feldmitte
    kraft["MaxSeilzugkraft"][0] = STfm_Last
    #   horiz. Komponente
    kraft["MaxSeilzugkraft"][1] = Hm
    #   am höchsten Punkt (mit mit Last in Feldmitte)
    kraft["MaxSeilzugkraft"][2] = STfm_Last + (np.max(z) - zfm) * qT

    j = np.seterr(**oldsettings)

    kraft["Durchhang"] = [y_LeerSeil, ym]
    kraft["LaengeSeil"] = [Laenge_Leerseil, Laenge_Leerseil_bei_0_kN, c]

    kraft["Spannkraft"] = [ST[0], ST[-1]]
    kraft["Seilzugkraft"] = [ST, Hs]

    poles.calculateAdvancedProperties(kraft, IS["Bundstelle"])

    anchorField = poles.getAnchorCable()
    cableline = {
        "xaxis": l_coord + poles.firstPole["dtop"],  # X-data starts at first pole
        "empty": z_coord_leer
        + poles.firstPole["ztop"],  # Y-data is calculated relative
        "load": z_coord_zweifel + poles.firstPole["ztop"],
        "anchorA": anchorField["A"],
        "anchorE": anchorField["E"],
        "groundclear_di": [],
        "groundclear": [],
        "groundclear_under": [],
        "groundclear_rel": [],
    }

    return cableline, kraft, seil_possible


def preciseCableLight(zi, di, IS, STA, HM, LP):
    qT = IS["qT"]  # [kN/m]

    # Seilfeld berechnen, b = Breite, h = Höhe
    b = di[LP[1:]] - di[LP[:-1]]
    h = zi[LP[1:]] * 0.1 + HM[1:] - zi[LP[:-1]] * 0.1 - HM[:-1]

    anzFelder = b.size
    anzStue = anzFelder + 1
    bquad = b**2

    if np.where(b == 0)[0].size == 0:
        korrektur = 0
    else:
        korrektur = h[0]

    z = np.zeros(anzStue)
    laengenprofil_y = np.zeros(anzStue)
    z[0] = korrektur
    z[1:] = np.cumsum(h)
    zfm = z[:-1] + h * 0.5
    b_cum = np.cumsum(b)
    laengenprofil_y[1:] = b_cum

    # Berechnung der Sehne
    c = (bquad + h**2) ** 0.5

    # Leerseil
    # --------
    zfmqT = zfm * qT
    STfm = STA + zfmqT  # Seilzugkraft in Feldmitte
    HT = (STfm * b) / c  # Horizontalkomponente der Seilzugkraft

    y_LeerSeil = b / (4 * HT) * (c * qT / 2)  # Leerseildurchhang

    tg_phi_o = (h / 2 + 2 * y_LeerSeil) / (b / 2)
    tg_phi_u = (h / 2 - 2 * y_LeerSeil) / (b / 2)
    phi_o = np.array([np.nan] * anzStue)
    phi_u = np.array([np.nan] * anzStue)
    phi_o[1:] = np.arctan(tg_phi_o)
    phi_u[:-1] = np.arctan(tg_phi_u)

    dAnkerA = IS["Ank"][0][0]
    zAnkerA = IS["Ank"][0][1]
    dAnkerE = IS["Ank"][0][2]
    zAnkerE = IS["Ank"][0][3]
    if dAnkerA != 0:
        phi_oA = np.arctan(zAnkerA / dAnkerA)
    else:
        phi_oA = np.nan
    if dAnkerE != 0:
        phi_oE = np.arctan((-1 * zAnkerE) / dAnkerE)
    else:
        phi_oE = np.nan
    # oldsettings = np.geterr()
    # j = np.seterr(all='ignore')

    phi_o[0] = phi_oA
    phi_u[-1] = phi_oE

    # TODO: Werte sind +- 2 Grad anders als in matlab
    leerseilknickwinkel = (phi_o - phi_u) / pi * 180

    zqT = z * qT
    ST = STA + zqT
    Vi = ST * (np.sin(phi_o) - np.sin(phi_u))  # Vi darf nicht negativ sein
    seilHebtAb = np.prod(Vi[1:-1] >= 0.01) == 0

    return leerseilknickwinkel, seilHebtAb


def updateWithCableCoordinates(cableline, pointA, azimut):
    cableline["coordx"] = pointA[0] + cableline["xaxis"] * sin(azimut)
    cableline["coordy"] = pointA[1] + cableline["xaxis"] * cos(azimut)
