In [1]:
%load_ext autoreload
%autoreload 2
%matplotlib inline

execfile ("_ImportScript.py")


1

In [2]:
MOCK = 0

In [3]:
execfile ("_ReconstructionScript.py")


/Users/LaurencePeanuts/Documents/Travail/Stanford/Music/Music/beatbox/universe.py:573: RuntimeWarning: divide by zero encountered in power
  self.Power_Spectrum = self.PSnorm*10000*np.power((self.k/self.kstar) ,(-3+(self.n_s-1)))
/Users/LaurencePeanuts/Documents/Travail/Stanford/Music/Music/beatbox/universe.py:407: ComplexWarning: Casting complex values to real discards the imaginary part
  ay_real[zero_ind] = value[zero_ind].astype(np.float)
Generated  122  potential Fourier coefficients
 with phases uniformly distributed between 0 and  3.14159265359
Built potential grid, with dimensions  (41, 41, 41)  and mean value  0.0 +/- 0.9633692
NSIDE = 256
ORDERING = RING in fits file
INDXSCHM = IMPLICIT
Sigma is 0.000000 arcmin (0.000000 rad) 
-> fwhm is 0.000000 arcmin
Sigma is 0.000000 arcmin (0.000000 rad) 
-> fwhm is 0.000000 arcmin
Built potential grid, with dimensions  (41, 41, 41)  and mean value  -0.0 +/- 39.5908795
Sigma is 0.000000 arcmin (0.000000 rad) 
-> fwhm is 0.000000 arcmin
/Users/LaurencePeanuts/Documents/Travail/Stanford/Music/Music/beatbox/multiverse.py:359: ComplexWarning: Casting complex values to real discards the imaginary part
  R_real[zero_ind,:] = beatbox.Universe.R[zero_ind,:].astype(np.float)
<matplotlib.figure.Figure at 0x108fcb9d0>
<matplotlib.figure.Figure at 0x11140d1d0>

In [6]:
from matplotlib import cm
import matplotlib.pyplot as plt
cmap = cm.RdBu_r
cmap.set_under('w')
dpi = 300
figsize_inch = 60, 40
fig = plt.figure(figsize=figsize_inch, dpi=dpi)
hp.mollview(WeRes.Tmap,  rot=(-90,0,0),min=-50, max=50,title="Residuals, $\ell_{max}$=%d" % We.truncated_lmax, cmap=cmap, unit="$\mu$K") 
plt.savefig("Residuals.png", dpi=dpi, bbox_inches="tight")


<matplotlib.figure.Figure at 0x117010e50>

In [ ]:
We.Pdist=1
We.Pmax=2*np.pi
We.show_potential_with_yt(angle=np.pi/3,  N_layer=4, cmap='jet', show3D=1, Slice=1, alpha_norm=50)

In [ ]:
beatbox.You.all_simulated_universes[-1].show_potential_with_yt(angle=np.pi/3,  N_layer=4, cmap='jet', show3D=1, Slice=0, alpha_norm=50)

In [ ]:
import matplotlib.pyplot as plt
fig = plt.figure()
plt.imshow(We.phi[:,:, 20], clim=(-0.25, 2), cmap=cmap)
plt.colorbar()

circ = plt.Circle((20, 20), radius=10, color='k', fill=False)
fig.add_subplot(111).add_artist(circ)

plt.savefig('../doc/proposals/2016NASA/figures/Mock_rec'+str(beatbox.Universe.truncated_lmax)+'_lmin'+str(beatbox.Universe.truncated_lmin)+'_nmax'+str(beatbox.Universe.truncated_nmax)+'_nmin'+str(beatbox.Universe.truncated_nmin)+'.pdf')


plt.show()

In [ ]:
fig = plt.figure()
imgplot = plt.imshow(beatbox.You.all_simulated_universes[-1].phi[:,:,20], clim=(-0.25, 2), cmap=cmap)
plt.colorbar()

circ = plt.Circle((20, 20), radius=10, color='k', fill=False)
fig.add_subplot(111).add_artist(circ)

plt.savefig('../doc/proposals/2016NASA/figures/Mock_true'+str(beatbox.Universe.truncated_lmax)+'_lmin'+str(beatbox.Universe.truncated_lmin)+'_nmax'+str(beatbox.Universe.truncated_nmax)+'_nmin'+str(beatbox.Universe.truncated_nmin)+'.pdf')


plt.show()

In [ ]:
They=beatbox.Universe()
They.read_in_CMB_T_map(from_this = '../data/commander_32band_Clsamples100/cmb_Cl_c0001_k00031.fits')
They.decompose_T_map_into_spherical_harmonics()
They.show_CMB_T_map(Tmap=They.Tmap, max=250, cmap=cmap, title='$\ell_{max}=$'+str(They.lmax))

In [ ]:
from matplotlib import cm
cmap = cm.RdBu_r
cmap.set_under('w')

They.show_lowest_spherical_harmonics_of_CMB_T_map(lmax=10,max=110/3, cmap=cmap, title='$\ell_{max}=10$')

In [ ]:
They.show_lowest_spherical_harmonics_of_CMB_T_map(lmax=2,max=25/3, cmap=cmap, title='$\ell_{max}=2$')

In [ ]:
They.show_lowest_spherical_harmonics_of_CMB_T_map(lmax=5,max=75/3, cmap=cmap, title='$\ell_{max}=5$')

In [ ]:
10*1e-5

In [ ]:
0.8*1e-5

In [ ]:
plt.plot(beatbox.You.all_simulated_universes[-1].fn, 'ro')
plt.plot(We.fn, 'k-')

error=np.zeros((len(We.fn),1))
error[:,0]= np.diag(beatbox.You.inv_A)
x=np.zeros((len(We.fn),1))
x[:,0]=range(len(We.fn))

plt.fill_between(x, We.fn-error, We.fn+error, where=((We.fn+error)>(We.fn-error)))
plt.show()

In [ ]:
print x.shape, (We.fn-error).shape, (We.fn+error).shape, ((We.fn+error)>(We.fn-error)).shape

In [ ]:
x=np.zeros((len(We.fn),1))
x[:,0]=range(len(We.fn))

In [ ]:
np.savetxt( "scratch/foryashar/errorplot-true_fnimagsorted.txt", beatbox.You.all_simulated_universes[-1].fn[index+len(index)])
np.savetxt( "scratch/foryashar/errorplot-bestfit_fnimagsorted.txt", We.fn[index+len(index)])
np.savetxt( "scratch/foryashar/errorplot-Covimagsorted.txt", np.diag(beatbox.You.inv_A)[index+len(index)])
np.savetxt( "scratch/foryashar/errorplot-x.txt",x)

In [ ]:
index=np.argsort(beatbox.Universe.k, axis=None)

In [ ]:
((beatbox.You.all_simulated_universes[-1].fngrid).reshape(-1))[index]

In [ ]:
ind = np.where(We.kfilter>0)
a=(We.k[ind])
print a[:len(ind[1])/2]

In [ ]:
index=np.argsort(a[:len(ind[1])/2], axis=None)

In [ ]:
print beatbox.You.all_simulated_universes[-1].fn[index+len(index)]

In [ ]:
index[0:2]

In [ ]:
x=np.zeros(len(We.fn)/2)
x=np.sort(a[:len(ind[1])/2], axis=None)
print x

In [ ]: