Commit d41ecb8b authored by Julien's avatar Julien
Browse files

Final PEP8 on estimators and xqml_utils

parent f5de1e0d
......@@ -75,7 +75,7 @@ def El(invCAA, invCBB, Pl, Bll=None, expend=False):
???
"""
# Tegmark B-matrix useless so far)
if Bll == None:
if Bll is None:
Bll = np.diagflat((np.ones(len(Pl))))
lmax = len(Pl) * 2**int(expend)
......@@ -83,9 +83,13 @@ def El(invCAA, invCBB, Pl, Bll=None, expend=False):
npix = len(invCAA)
# triangular shape ds_dcb
if expend:
El = np.array([np.dot(np.dot(invCAA, Expd(Pl,l)), invCBB) for l in lrange]).reshape((lmax,npix,npix))
El = np.array(
[np.dot(np.dot(invCAA, Expd(Pl, l)), invCBB) for l in lrange]
).reshape((lmax, npix, npix))
else:
El = np.array([np.dot(np.dot(invCAA, Pl[l]), invCBB) for l in lrange]).reshape((lmax,npix,npix))
El = np.array(
[np.dot(np.dot(invCAA, Pl[l]), invCBB) for l in lrange]
).reshape((lmax, npix, npix))
return El
def ElLong(invCAA, invCBB, Pl, Bll=None, expend=False):
......@@ -102,7 +106,7 @@ def ElLong(invCAA, invCBB, Pl, Bll=None, expend=False):
???
"""
# Tegmark B-matrix useless so far)
if Bll == None:
if Bll is None:
Bll = np.diagflat((np.ones(len(Pl))))
lmax = len(Pl) * 2**int(expend)
......@@ -111,14 +115,16 @@ def ElLong(invCAA, invCBB, Pl, Bll=None, expend=False):
# triangular shape ds_dcb
if expend:
for l in lrange:
Pl[l] = np.dot(np.dot(invCAA, Expd(Pl, l)), invCBB).reshape((npix, npix ))
Pl[l] = np.dot(
np.dot(invCAA, Expd(Pl, l)), invCBB).reshape((npix, npix))
else:
for l in lrange:
Pl[l] = np.dot(np.dot(invCAA, Pl[l]), invCBB).reshape((npix, npix ))
Pl[l] = np.dot(np.dot(invCAA, Pl[l]), invCBB).reshape((npix, npix))
def CrossFisherMatrix(El, Pl, expend=False):
"""
Compute cross (or auto) fisher matrix Fll = Trace[invCAA.Pl.invCBB.Pl] = Trace[El.Pl]
Compute cross (or auto) fisher matrix
Fll = Trace[invCAA.Pl.invCBB.Pl] = Trace[El.Pl]
Parameters
----------
......@@ -133,15 +139,20 @@ def CrossFisherMatrix(El, Pl, expend=False):
lrange = np.arange(lmax)
if expend:
# pas de transpose car symm
FAB = np.array([np.sum(El[il] * Expd(Pl,jl)) for il in lrange for jl in lrange]).reshape(lmax, lmax)
FAB = np.array(
[np.sum(El[il] * Expd(Pl, jl)) for il in lrange for jl in lrange]
).reshape(lmax, lmax)
else:
# pas de transpose car symm
FAB = np.array([np.sum(El[il] * Pl[jl]) for il in lrange for jl in lrange]).reshape(lmax, lmax)
FAB = np.array(
[np.sum(El[il] * Pl[jl]) for il in lrange for jl in lrange]
).reshape(lmax, lmax)
return FAB
def CrossWindowFunction(El, Pl):
"""
Compute Tegmark cross (or auto) window matrix Wll = Trace[invCAA.Pl.invCBB.Pl] = Trace[El.Pl]
Compute Tegmark cross (or auto) window matrix
Wll = Trace[invCAA.Pl.invCBB.Pl] = Trace[El.Pl]
Equivalent to Fisher matrix Fll when Tegmark B-matrix = 1
Parameters
......@@ -156,12 +167,15 @@ def CrossWindowFunction(El, Pl):
lmax = len(El)
lrange = np.arange((lmax))
# pas de transpose car symm
Wll = np.array([np.sum(El[il] * Pl[jl]) for il in lrange for jl in lrange]).reshape(lmax, lmax)
Wll = np.array(
[np.sum(El[il] * Pl[jl]) for il in lrange for jl in lrange]
).reshape(lmax, lmax)
return Wll
def CrossWindowFunctionLong(invCAA, invCBB, Pl):
"""
Compute Tegmark cross (or auto) window matrix Wll = Trace[invCAA.Pl.invCBB.Pl] = Trace[El.Pl]
Compute Tegmark cross (or auto) window matrix
Wll = Trace[invCAA.Pl.invCBB.Pl] = Trace[El.Pl]
Equivalent to Fisher matrix Fll when Tegmark B-matrix = 1
Parameters
......@@ -176,7 +190,9 @@ def CrossWindowFunctionLong(invCAA, invCBB, Pl):
lmax = len(Pl)
lrange = np.arange((lmax))
# Pas de transpose car symm
Wll = np.array([np.sum(np.dot(np.dot(invCAA, Pl[il]), invCBB) * Pl[jl]) for il in lrange for jl in lrange]).reshape(lmax, lmax)
Wll = np.array(
[np.sum(np.dot(np.dot(invCAA, Pl[il]), invCBB) * Pl[jl])
for il in lrange for jl in lrange]).reshape(lmax, lmax)
return Wll
def CrossMisherMatrix(El, CAA, CBB):
......@@ -194,9 +210,11 @@ def CrossMisherMatrix(El, CAA, CBB):
"""
lmax = len(El)
lrange = np.arange(lmax)
El_CAA = np.array([np.dot(CAA,El[il]) for il in lrange])
El_CBB = np.array([np.dot(CBB,El[il]) for il in lrange])
MAB = np.array([np.sum(El_CAA[il] * El_CBB[jl].T) for il in lrange for jl in lrange]).reshape(lmax, lmax)
El_CAA = np.array([np.dot(CAA, El[il]) for il in lrange])
El_CBB = np.array([np.dot(CBB, El[il]) for il in lrange])
MAB = np.array(
[np.sum(El_CAA[il] * El_CBB[jl].T) for il in lrange for jl in lrange]
).reshape(lmax, lmax)
return MAB
def CrossGisherMatrix(El, CAB):
......@@ -214,8 +232,10 @@ def CrossGisherMatrix(El, CAB):
"""
lmax = len(El)
lrange = np.arange(lmax)
El_CAB = np.array([np.dot(CAB,El[il]) for il in lrange])
GAB = np.array([np.sum(El_CAB[il] * El_CAB[jl].T) for il in lrange for jl in lrange]).reshape(lmax, lmax)
El_CAB = np.array([np.dot(CAB, El[il]) for il in lrange])
GAB = np.array(
[np.sum(El_CAB[il] * El_CAB[jl].T) for il in lrange for jl in lrange]
).reshape(lmax, lmax)
return GAB
def CrossGisherMatrixLong(El, CAB):
......@@ -233,7 +253,9 @@ def CrossGisherMatrixLong(El, CAB):
"""
lmax = len(El)
lrange = np.arange(lmax)
GAB = np.array([np.sum(np.dot(CAB, El[il]) * np.dot(CAB,El[jl]).T) for il in lrange for jl in lrange]).reshape(lmax, lmax)
GAB = np.array(
[np.sum(np.dot(CAB, El[il]) * np.dot(CAB, El[jl]).T)
for il in lrange for jl in lrange]).reshape(lmax, lmax)
return GAB
def yQuadEstimator(dA, dB, El):
......
......@@ -15,12 +15,15 @@ def GetSizeNumpies():
To be copy/past in ipython to get current memory size taken by numpy arrays
"""
import operator
np_arrays = {k:v for k,v in locals().items() if isinstance(v, np.ndarray)}
np_arrayssize = {key : round((np_arrays[key]).nbytes/1024./1024./1024., 4) for key in np_arrays }
np_arrays = {
k: v for k, v in locals().items() if isinstance(v, np.ndarray)}
np_arrayssize = {
key: round((np_arrays[key]).nbytes/1024./1024./1024., 4)
for key in np_arrays}
sorted_x = sorted(np_arrayssize.items(), key=operator.itemgetter(1))
for s in np.arange(len(sorted_x)):
print(sorted_x[s])
sumx = sum(np.array(sorted_x)[:,1].astype(np.float))
sumx = sum(np.array(sorted_x)[:, 1].astype(np.float))
print("total = " + str(sumx)+"Gb")
def ComputeSizeDs_dcb(nside, fsky, deltal=1):
......@@ -36,9 +39,10 @@ def ComputeSizeDs_dcb(nside, fsky, deltal=1):
----------
???
"""
sizeds_dcb = (2*12*nside**2*fsky)**2*8*2*(3.*nside/deltal)/1024./1024./1024.
print("size (Gb) = "+str(sizeds_dcb) )
print("possible reduced size (Gb) = "+str(sizeds_dcb/4) )
toGB = 1024. * 1024. * 1024.
sizeds_dcb = (2*12*nside**2*fsky)**2*8*2*(3.*nside/deltal) / toGB
print("size (Gb) = " + str(sizeds_dcb))
print("possible reduced size (Gb) = " + str(sizeds_dcb/4))
def get_colors(num_colors):
......@@ -55,16 +59,16 @@ def get_colors(num_colors):
???
"""
import colorsys
colors=[]
colors = []
ndeg = 250.
for i in np.arange(0., ndeg, ndeg / num_colors):
hue = i/360.
lightness = 0.5#(50 + np.random.rand() * 10)/100.
saturation = 0.7#(90 + np.random.rand() * 10)/100.
lightness = 0.5
saturation = 0.7
colors.append(colorsys.hls_to_rgb(hue, lightness, saturation))
return np.array(colors)
def progress_bar(i,n, dt):
def progress_bar(i, n, dt):
"""
???
......@@ -78,20 +82,22 @@ def progress_bar(i,n, dt):
???
"""
if n != 1:
ntot=50
ndone=ntot*i/(n-1)
a='\r|'
ntot = 50
ndone = ntot * i / (n - 1)
a = '\r|'
for k in np.arange(ndone):
a += '#'
for k in np.arange(ntot-ndone):
a += ' '
fra = i/(n-1.)
remain = round(dt/fra*(1-fra))
a += '| '+str(int(100.*fra))+'%'+" : "+str(remain)+" sec = " +str(round(remain/60.,1)) + " min"
min = str(round(remain/60., 1))
a += '| '+str(int(100.*fra))+'%'+" : "+str(remain)+" sec = "+min+" min"
sys.stdout.write(a)
sys.stdout.flush()
if i == n-1:
sys.stdout.write(' Done. Total time = '+str(np.ceil(dt))+" sec = " +str(round(remain/60.,1)) + " min \n")
sys.stdout.write(
' Done. Total time = '+str(np.ceil(dt))+" sec = "+min+" min\n")
sys.stdout.flush()
......@@ -126,11 +132,7 @@ def randomword(length):
"""
return ''.join(rd.choice(string.lowercase) for i in range(length))
################################# Not used #################################
def cov_from_maps(maps0,maps1):
def cov_from_maps(maps0, maps1):
"""
???
......@@ -143,26 +145,25 @@ def cov_from_maps(maps0,maps1):
----------
???
"""
sh=np.shape(maps0)
npix=sh[1]
nbmc=sh[0]
covmc=np.zeros((npix,npix))
mm0=np.mean(maps0,axis=0)
sh = np.shape(maps0)
npix = sh[1]
nbmc = sh[0]
covmc = np.zeros((npix, npix))
mm0 = np.mean(maps0, axis=0)
# print(mm0)
mm1=np.mean(maps1,axis=0)
mm1 = np.mean(maps1, axis=0)
# print(mm1)
themaps0=np.zeros((nbmc,npix))
themaps1=np.zeros((nbmc,npix))
themaps0 = np.zeros((nbmc, npix))
themaps1 = np.zeros((nbmc, npix))
start = timeit.default_timer()
for i in np.arange(npix):
progress_bar(i,npix,timeit.default_timer()-start)
themaps0[:,i]=maps0[:,i]-mm0[i]
themaps1[:,i]=maps1[:,i]-mm1[i]
progress_bar(i, npix, timeit.default_timer()-start)
themaps0[:, i] = maps0[:, i] - mm0[i]
themaps1[:, i] = maps1[:, i] - mm1[i]
print('hy')
start = timeit.default_timer()
for i in np.arange(npix):
progress_bar(i,npix,timeit.default_timer()-start)
progress_bar(i, npix, timeit.default_timer()-start)
for j in np.arange(npix):
covmc[i,j]=np.mean(themaps0[:,i]*themaps1[:,j])
#covmc[j,i]=covmc[i,j]
covmc[i, j] = np.mean(themaps0[:, i] * themaps1[:, j])
return(covmc)
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment