import pyfits
import numpy

bl = pyfits.open('bllac.fits')
cts = bl[0].data
bl.close()

# Pozor na prehozeni sloupecku a radku
#print(cts[971,1126])
print("(Asi) Nejjasnejsi pixel hvezdy B = ",cts[1126,971])

# polovina velikosti ctverce, ve kterem se meri
# obloha i hvezda
box = 12
detail = cts[1126-box:1126+box,971-box:971+box]
s = numpy.sum(detail)
ns = detail.size
print("Suma fotonu (hvezda+obloha), pocet pixelu =",s,ns)

obloha = cts[1079-box:1079+box,1114-box:1114+box]
b = numpy.mean(obloha)
print("Uroven oblohy =",b)

ct = s - b*ns
print("Pocet pozorovanych fotonu od hvezdy B = ",ct)

m = 12.78
flam = 3.95e-11 * 10**(-m/2.5)
n = 0.3 * 90 * flam * 70 / (6.7e-34 * 3e8 / 555e-9)
eta = ct / n
r = n / ct
print("Magnituda, flam, ocekavane fotony, merene fotony, ucinnost,r",
      m,flam,n,ct,eta,r)