Calculează tu: stereograma, înclinarea aparentă, cele trei puncte

6.6. Calculează tu: stereograma, înclinarea aparentă, cele trei puncte#

Trei calcule scurte din acest capitol, pe care le poți rula și modifica. Schimbă numerele din prima celulă a fiecărui exemplu și vezi ce se întâmplă.

6.6.1. Un plan pe rețeaua stereografică#

Planul este dat prin sens de înclinare/înclinare. Codul desenează urma lui ciclografică și polul, în proiecție stereografică pe emisfera inferioară, cu relația (6.1).

sens, inclinare = 120, 35     # poziția planului: sens de înclinare / înclinare [grade]

Hide code cell source

 
import numpy as np
import matplotlib.pyplot as plt

def proiectie(azimut, plonjare):
    """linie (azimut, plonjare) -> (x, y) pe rețeaua de rază 1; nordul în sus"""
    r = np.tan(np.radians(90 - plonjare) / 2)
    a = np.radians(azimut)
    return r * np.sin(a), r * np.cos(a)

def urma_ciclografica(sens, inclinare, n=181):
    """liniile conținute în plan: unghiul t se măsoară în plan, de la direcție"""
    t = np.radians(np.linspace(0, 180, n))
    d, s = np.radians(inclinare), np.radians(sens)
    # vector în plan: combinație între direcție (orizontală) și linia de cea mai mare pantă
    e = np.cos(t) * np.sin(s - np.pi / 2) + np.sin(t) * np.cos(d) * np.sin(s)   # est
    n_ = np.cos(t) * np.cos(s - np.pi / 2) + np.sin(t) * np.cos(d) * np.cos(s)  # nord
    jos = np.sin(t) * np.sin(d)                                                 # în jos
    az = np.degrees(np.arctan2(e, n_)) % 360
    pl = np.degrees(np.arcsin(np.clip(jos, 0, 1)))
    return proiectie(az, pl)

fig, ax = plt.subplots(figsize=(5, 5))
u = np.linspace(0, 2 * np.pi, 361)
ax.plot(np.sin(u), np.cos(u), "k", lw=1)
for az in range(0, 360, 30):
    ax.plot(*proiectie(np.array([az, az]), np.array([0, 3])), "k", lw=0.8)
for txt, az in (("N", 0), ("E", 90), ("S", 180), ("V", 270)):
    ax.text(*(1.09 * np.array([np.sin(np.radians(az)), np.cos(np.radians(az))])), txt, ha="center", va="center")

ax.plot(*urma_ciclografica(sens, inclinare), color="tab:red", lw=2, label=f"planul {sens:03d}/{inclinare}")
pol_az, pol_pl = (sens + 180) % 360, 90 - inclinare
ax.plot(*proiectie(pol_az, pol_pl), "o", color="tab:red", label=f"polul: azimut {pol_az:03d}°, plonjare {pol_pl}°")
ax.plot(*proiectie(sens, inclinare), "s", color="tab:blue", label="linia de cea mai mare pantă")
ax.plot(0, 0, "k+")
ax.set_aspect("equal"); ax.axis("off")
ax.legend(loc="upper center", bbox_to_anchor=(0.5, -0.02), frameon=False)
plt.show()
print(f"polul se află la r = {np.tan(np.radians(inclinare) / 2):.3f} R de centru")
../_images/10a2f10bf58dbd81f496ff46dab2f0af5ae9bdb5a5a79e821d40a52d605046fa.png 
polul se află la r = 0.315 R de centru

Încearcă inclinare = 0 (plan orizontal: urma este cercul primitiv, polul este în centru) și inclinare = 90 (plan vertical: urma este un diametru, polul este pe margine).

6.6.2. Înclinarea aparentă#

Cât de înclinat pare un strat într-o secțiune care nu este trasată pe sensul lui de înclinare? Relația (6.2).

alpha = 40        # înclinarea reală [grade]

Hide code cell source

 
beta = np.arange(0, 91, 15)
ap = np.degrees(np.arctan(np.tan(np.radians(alpha)) * np.cos(np.radians(beta))))
print("unghiul dintre secțiune și sensul de înclinare  ->  înclinarea aparentă")
for b, a in zip(beta, ap):
    print(f"   {b:3d}°   ->   {a:5.1f}°")

b = np.linspace(0, 90, 91)
plt.figure(figsize=(6, 3.2))
plt.plot(b, np.degrees(np.arctan(np.tan(np.radians(alpha)) * np.cos(np.radians(b)))))
plt.xlabel("β – unghiul dintre traseul secțiunii și sensul de înclinare [°]")
plt.ylabel("înclinarea aparentă [°]"); plt.grid(alpha=0.3)
plt.show()

unghiul dintre secțiune și sensul de înclinare  ->  înclinarea aparentă
     0°   ->    40.0°
    15°   ->    39.0°
    30°   ->    36.0°
    45°   ->    30.7°
    60°   ->    22.8°
    75°   ->    12.3°
    90°   ->     0.0°
../_images/14571d68145e1b391b2fbce71a05a5da782e36cc257b7813185aefc132f6f774.png

Până la β ≈ 30° eroarea este mică; dincolo de 60° stratul pare aproape de două ori mai puțin înclinat decât este.

6.6.3. Problema celor trei puncte#

Trei puncte în care a fost întâlnită aceeași suprafață (de exemplu baza unui strat, în trei foraje) definesc un plan. Codul îl găsește pentru orice trei puncte necoliniare; aplicația A2 din Aplicații (laborator) este același calcul, făcut grafic.

# (x spre est [m], y spre nord [m], cota suprafeței [m])
A = (100.0, 400.0, 400.0)
B = (500.0, 450.0, 250.0)
C = (350.0, 100.0, 250.0)

Hide code cell source

 
p1, p2, p3 = map(np.array, (A, B, C))
n = np.cross(p2 - p1, p3 - p1)          # normala la plan
if abs(n[2]) < 1e-9:
    raise ValueError("puncte coliniare (în plan) sau plan vertical: schimbă datele")
if n[2] < 0:
    n = -n                              # normala orientată în sus
# planul coboară în sensul componentei orizontale a normalei
sens = np.degrees(np.arctan2(n[0], n[1])) % 360
inclinare = np.degrees(np.arctan(np.hypot(n[0], n[1]) / n[2]))
print(f"poziția planului: {sens:03.0f}/{inclinare:.0f}   (sens de înclinare / înclinare)")
print(f"direcția: {(sens - 90) % 180:03.0f}° – {(sens - 90) % 180 + 180:03.0f}°")
grad = np.tan(np.radians(inclinare))
print(f"izobate din 50 în 50 m: distanța orizontală dintre ele = {50 / grad:.0f} m")

poziția planului: 113/23   (sens de înclinare / înclinare)
direcția: 023° – 203°
izobate din 50 în 50 m: distanța orizontală dintre ele = 116 m

Mută punctul C pe dreapta AB și rulează din nou: calculul refuză, pentru că trei puncte coliniare nu definesc un plan (vezi și exercițiul 8 al capitolului).

Cuprins