# Lite Python-script for "Raketten",
# til seminaroppgave 5.2, seminar uke 05, 18.02.2025
# NB! Python 3, ikke python 2 som i boka...
# Ann-Cecilie Larsen, a.c.larsen@fys.uio.no

# Importerer nyttige Python-bibliotek
import numpy as np 
import matplotlib.pyplot as plt

# Definerer vektorer t (for tid) og a (for akselerasjon),
# leser input-filen, og putter verdiene i kolonne 0 inn i tidsvektoren
# og verdiene i kolonne 1 i akselerasjonsvektoren
t, a = np.loadtxt('therocketdata.txt', usecols = [0, 1], unpack = True)

# Finner antall datapunkter (for ? lage vektorer for hastighet og posisjon)
n = len(t)

# Lage vektorer for posisjon x(t) og v(t)
x = np.zeros(n, float)
v = np.zeros(n, float)

# Definerer initialbetingelser 
x[0] = 0.0
v[0] = 0.0

# Definerer tidsintervallet (her er det likt for alle m?lepunktene)
dt = t[1] - t[0]

# Setter opp l?kka for numerisk integrasjon med forward-Euler-metodene
for i in range(0, n-1):
	v[i+1] = v[i] + a[i]*dt
	# Hvis vi bytter ut v[i] med v[i+1] n?r vi beregner posisjonen, 
	# f?r vi Euler-Cromer-metoden i stedet for Euler
	x[i+1] = x[i] + v[i]*dt 

# N? skal vi plotte. Se https://matplotlib.org/3.1.1/api/_as_gen/matplotlib.pyplot.subplot.html
# F?rst posisjonen som funksjon av tid
plt.subplot(3,1,1)
plt.plot(t,x)
# Husk aksetitler og enheter!
plt.xlabel('t [s]')
plt.ylabel('x [m]')

# S? plotter vi hastighet som funksjon av tid
plt.subplot(3,1,2)
plt.plot(t,v)
plt.xlabel('t [s]')
plt.ylabel('v [m/s]')

# Tilslutt plotter vi akselerasjonen som funksjon av tid
plt.subplot(3,1,3)
plt.plot(t,a)
plt.xlabel('t [s]')
plt.ylabel('a [m/s$^2$]') 

plt.show()