import numpy as np
import matplotlib.pyplot as plt

# a)
# Vid odvodenie ukol5_1.jpg.

# b)
# Z teoretickeho odvodenia v casti a) pouzijeme logaritmicky vztah:
# m = b1 + b2*log(N), kde b2 = 5/3 = 1.67.
matrix = np.loadtxt('list_of_stars.dat')
x = matrix[:, 0]   # hodnoty poradia hviezd
y = matrix[:, 1]   # hodnoty zdanlivej hviezdnej velkosti
n = len(x)
X = np.ones((n, 2))   # matica parcialnych derivacii
X[:, 1] = np.log10(x)
V = np.dot(X.T, X)
U = np.dot(X.T, y)
# Hodnoty parametrov b = [b1, b2] pre najpresnejsi fit.
b = np.dot(np.linalg.inv(V), U)
f = np.dot(X, b)   # funkcne hodnoty predpovedane fitom
print('Modelova funkcia ma tvar: m = {} + {}*log(N)'
      .format(np.round(b[0], 2), np.round(b[1], 2)))
plt.plot(x, y, '*')
plt.plot(x, f)
plt.xlabel('poradie')
plt.ylabel('zdanliva hviezdna velkost [mag]')

# V skutocnosti vidime volnym okom menej hviezd, nez bolo teoreticky
# predpovedane. Dovodmi mozu byt medzihviezdna extinkcia, nehomogenne
# rozlozenie hviezd v priestore, ktore sa sustreduju hlavne v oblasti disku a
# tiez rozna absolutna hviezdna velkost hviezd s preferenciou tych jasnejsich.