import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
from matplotlib import cm

import cvicenie5 as c5


#a)
def my_modus(x):
    x = np.sort(x)
    N = 1   #najvyssi pocet vyskytov 
    i = 0
    result = x[0]   #aktualny modus
    while i <= len(x):
        n = 1
        while i+n < len(x) and x[i] == x[i+n]:
            n += 1
        if n > N:   #ak najdeme novy modus s castejsim vyskytom, menime N a result
            N = n
            result = x[i]
        i += n
    return result, N   #vracia modus a jeho pocetnost

m44 = np.loadtxt("m44.dat")
print(my_modus(m44.flatten()))

#b)
m44_smooth = m44
background = 1721
sigma_r = c5.sigma_r
N = len(m44)
for i in range(len(m44)):
    for j in range(len(m44)):
        if np.abs(m44_smooth[i][j] - background) < 3.5*sigma_r:
            m44_smooth[i][j] = background
x, y = c5.x, c5.y
c5.surface(x, y, m44_smooth, 1)   #posledne 4 riadky su skopirovane z predchadzajuceho cvicenia

#c)
def inst_mag(matrix, background, x, y, r, sigma):
    n = (2*r + 1)**2.0   #pocet pixelov vyrezu z okolia hviezdy
    I = np.sum(matrix[x - r : x + r + 1, y - r : y + r + 1]) - n*background
    dI = np.sqrt(n)*sigma
    i_mag = -2.5*np.log10(I)
    di_mag1 = -2.5*np.log10(I + dI) + 2.5*np.log10(I)   #vypocet hornej chyby
    di_mag2 = -2.5*np.log10(I) + 2.5*np.log10(I - dI)   #vypocet dolnej chyby
    return i_mag, di_mag1, -di_mag2

#d)
def find_stars(matrix, n):   #najde n najjasnejsich pixelov spolu s ich poziciou
    sequence = np.sort(matrix.flatten())   #zoradi prvky matice vzostupne
    sequence = sequence[::-1]   #zoradi vektor zostupne
    result = []
    for i in range(n):
        x, y = np.where(matrix == sequence[i])   #vrati pozicie najjasnejsich pixelov
        result.append([sequence[i], x, y])
    return np.array(result)

print(find_stars(m44,20))
imag_gCnc = inst_mag(m44, 1721, 193, 19, 3, sigma_r)
imag_M44 = inst_mag(m44, 1721, 125, 65, 30, sigma_r)
imag_lim = -2.5*np.log10(5*sigma_r)
mag_corr = 4.65 - imag_gCnc[0]   #korekcia medzi instrumentalnou a zdanlivou velkostou
print('instrumentalna hviezdna velkost gamma Cnc:', imag_gCnc)
print('instrumentalna hviezdna velkost M44:', imag_M44)
print('zdanliva hviezdna velkost M44:', imag_M44[0] + mag_corr)
print('limitna hviezdna velkost:', imag_lim + mag_corr)