Skip to content
GitLab
Projects
Groups
Snippets
Help
Loading...
Help
Help
Support
Community forum
Keyboard shortcuts
?
Submit feedback
Contribute to GitLab
Sign in / Register
Toggle navigation
Open sidebar
xQML
xQML
Commits
e7960521
Commit
e7960521
authored
Sep 26, 2018
by
Syl
Browse files
use dlss for P02 instead of F10
parent
3325d8d4
Changes
2
Hide whitespace changes
Inline
Side-by-side
Showing
2 changed files
with
13 additions
and
16 deletions
+13
-16
TestEsti.py
TestEsti.py
+8
-8
libcov.py
libcov.py
+5
-8
No files found.
TestEsti.py
View file @
e7960521
...
...
@@ -30,16 +30,16 @@ clth = np.array(hp.read_cl('planck_base_planck_2015_TTlowP.fits'))
Clthshape
=
zeros
(((
6
,)
+
shape
(
clth
)[
1
:]))
Clthshape
[:
4
]
=
clth
clth
=
Clthshape
EB
=
0.
1
EB
=
0.
5
clth
[
4
]
=
EB
*
sqrt
(
clth
[
1
]
*
clth
[
2
])
TB
=
0.
0
TB
=
0.
5
clth
[
5
]
=
TB
*
sqrt
(
clth
[
0
]
*
clth
[
2
])
lth
=
arange
(
2
,
lmax
+
1
)
spec
=
None
#['E
B']
#
spec =
['TE', 'EB', 'T
B']
temp
=
True
polar
=
Fals
e
corr
=
Fals
e
polar
=
Tru
e
corr
=
Tru
e
pixwin
=
True
ellbins
=
arange
(
2
,
lmax
+
2
,
deltal
)
...
...
@@ -51,7 +51,7 @@ nbins = len(ellbins) - 1
# Create mask
t
,
p
=
hp
.
pix2ang
(
nside
,
range
(
hp
.
nside2npix
(
nside
)))
mask
=
np
.
ones
(
hp
.
nside2npix
(
nside
),
bool
)
#
mask[abs(90 - rad2deg(t)) < 30] = False
mask
[
abs
(
90
-
rad2deg
(
t
))
<
30
]
=
False
npix
=
sum
(
mask
)
fwhm
=
0.5
...
...
@@ -93,7 +93,7 @@ esti.construct_esti(NoiseVar, NoiseVar)
fpixw
=
extrapolpixwin
(
nside
,
lmax
+
2
,
pixwin
=
pixwin
)
start
=
timeit
.
default_timer
()
for
n
in
np
.
arange
(
nsimu
):
progress_bar
(
n
,
nsimu
,
timeit
.
default_timer
()
-
start
)
#
progress_bar(n, nsimu, timeit.default_timer() - start)
cmb
=
np
.
array
(
hp
.
synfast
(
clth
[:,
:
len
(
fpixw
)]
*
(
fpixw
*
bl
)
**
2
,
nside
,
pixwin
=
False
,
lmax
=
Slmax
,
fwhm
=
deg2rad
(
fwhm
),
new
=
True
,
verbose
=
False
))
...
...
@@ -164,7 +164,7 @@ grid()
show
()
savefig
(
"../Plots/Git/"
+
"test0.png"
)
#
savefig("../Plots/Git/"+"test0.png")
if
__name__
==
"__main__"
:
"""
...
...
libcov.py
View file @
e7960521
...
...
@@ -12,9 +12,6 @@ import healpy as hp
import
math
from
scipy
import
special
# from spin_functions import dlss, pl0
# from spin_functions import F1l2, F2l2
# from libangles import polrotangle
from
simulation
import
getstokes
from
simulation
import
extrapolpixwin
from
xqml_utils
import
progress_bar
...
...
@@ -215,8 +212,9 @@ def compute_PlS(
if
timing
:
progress_bar
(
i
,
npi
,
-
0.5
*
(
start
-
timeit
.
default_timer
()))
for
j
in
np
.
arange
(
i
,
npi
):
cos_chi
=
allcosang
[
i
,
j
]
if
temp
:
pl
=
norm
*
pl0
(
allcosang
[
i
,
j
]
,
Slmax
)[
2
:]
pl
=
norm
*
pl0
(
cos_chi
,
Slmax
)[
2
:]
elem
=
np
.
sum
((
pl
*
clthn
[
0
]))
S
[
i
,
j
]
=
elem
S
[
j
,
i
]
=
elem
...
...
@@ -230,7 +228,6 @@ def compute_PlS(
jj
=
j
+
ponpi
cij
,
sij
=
polrotangle
(
rpix
[:,
i
],
rpix
[:,
j
])
cji
,
sji
=
polrotangle
(
rpix
[:,
j
],
rpix
[:,
i
])
cos_chi
=
allcosang
[
i
,
j
]
# Tegmark version
Q22
=
norm
*
F1l2
(
cos_chi
,
Slmax
)
...
...
@@ -245,9 +242,9 @@ def compute_PlS(
if
TE
or
TB
:
# # Matt version
#
d20 = -dlss(cos_chi, 2, 0, Slmax)[2:]
#
P02 = norm * d20
P02
=
-
norm
*
F1l0
(
cos_chi
,
Slmax
)
d20
=
-
dlss
(
cos_chi
,
2
,
0
,
Slmax
)[
2
:]
P02
=
norm
*
d20
#
P02 = -norm * F1l0(cos_chi, Slmax)
elemA
=
0
elemB
=
0
elemC
=
0
...
...
Write
Preview
Markdown
is supported
0%
Try again
or
attach a new file
.
Attach a file
Cancel
You are about to add
0
people
to the discussion. Proceed with caution.
Finish editing this message first!
Cancel
Please
register
or
sign in
to comment