In [1]:
%matplotlib inline
In [5]:
import numpy as np
import matplotlib.pyplot as plt
from astropy.io import ascii
In [8]:
def normalization(data):
'''
Normalize the magnitude of the star
'''
norm = np.ma.copy(data)
norm[:,1] = (norm[:,1] - np.min(norm[:,1])) \
/ (np.max(norm[:,1]) - np.min(norm[:,1]))
return norm
In [45]:
def cond_entropy(period, data, p_bins=10, m_bins=5):
'''
Compute the conditional entropy for the normalized observations
'''
if period <= 0:
return np.PINF
r = rephase(data, period)
bins, *_ = np.histogram2d(r[:,0], r[:,1], [p_bins, m_bins],
[[0,1], [0,1]])
size = r.shape[0]
if size > 0:
divided_bins = bins / size
arg_positive = divided_bins > 0
column_sums = np.sum(divided_bins, axis=1) #change 0 by 1
column_sums = np.repeat(np.reshape(column_sums,(p_bins,1)),
m_bins, axis=1)
#column_sums = np.repeat(np.reshape(column_sums, (1,-1)),
# p_bins, axis=0)
select_divided_bins = divided_bins[arg_positive]
select_column_sums = column_sums[arg_positive]
A = np.empty((p_bins, m_bins), dtype=float)
# bins[i,j]/size * log(bins[i,:] / size / (bins[i,j]/size))
A[ arg_positive] = select_divided_bins \
* np.log(select_column_sums / select_divided_bins)
A[~arg_positive] = 0
return np.sum(A)
else:
return np.PINF
In [15]:
def rephase(data, period, col=0, copy=True):
'''
transform the time of observations to the phase space
'''
rephased = np.ma.array(data, copy=copy)
rephased[:, col] = get_phase(rephased[:, col], period)
return rephased
def get_phase(time, period):
'''
divide the time of observations by the period
'''
return (time / period) % 1
In [52]:
#data = ascii.read('/home/glauffer/Dropbox/FURG/final_project/data/OGLE-LMC-CEP-0018.dat')
data = np.array(np.loadtxt('/home/glauffer/Dropbox/FURG/final_project/data/OGLE-LMC-CEP-0018.dat'))
norm = normalization(data)
phased = rephase(norm, 4.0478)
ent = cond_entropy(4.0478, phased)
plt.plot(phased.T[0], phased.T[1], 'bo')
plt.xticks(np.linspace(0,1,11))
plt.grid()
plt.xlabel('Fase' , fontsize = 14)
plt.ylabel('Magnitude Normalizada', fontsize = 14 )
plt.title('OGLE-LMC-CEP-0018', fontsize = 18)
ent
Out[52]:
In [53]:
phased1 = rephase(norm, 1.0)
ent1 = cond_entropy(1.0, phased1)
plt.plot(phased1.T[0], phased1.T[1], 'bo')
plt.xticks(np.linspace(0,1,11))
plt.grid()
plt.xlabel('Fase' , fontsize = 14)
plt.ylabel('Magnitude Normalizada', fontsize = 14 )
plt.title('OGLE-LMC-CEP-0018', fontsize = 18)
ent1
Out[53]:
In [55]:
phased3 = rephase(norm, 3.0123)
ent3 = cond_entropy(3.0123, phased3)
plt.plot(phased3.T[0], phased3.T[1], 'bo')
plt.xticks(np.linspace(0,1,11))
plt.grid()
plt.xlabel('Fase' , fontsize = 14)
plt.ylabel('Magnitude Normalizada', fontsize = 14 )
plt.title('OGLE-LMC-CEP-0018', fontsize = 18)
ent3
Out[55]:
In [ ]: