Skip to content
GitLab
Menu
Projects
Groups
Snippets
Help
Help
Support
Community forum
Keyboard shortcuts
?
Submit feedback
Contribute to GitLab
Sign in
Toggle navigation
Menu
Open sidebar
ift
NIFTy
Commits
dd87b630
Commit
dd87b630
authored
Jul 14, 2017
by
Theo Steininger
Browse files
Merge branch 'Branch_master' into 'master'
Branch master See merge request
!164
parents
997bbf68
e7db670d
Pipeline
#14878
passed with stages
in 12 minutes and 30 seconds
Changes
6
Pipelines
1
Hide whitespace changes
Inline
Side-by-side
demos/wiener_filter_easy.py
deleted
100644 → 0
View file @
997bbf68
import
numpy
as
np
from
nifty
import
RGSpace
,
PowerSpace
,
Field
,
FFTOperator
,
ComposedOperator
,
\
SmoothingOperator
,
DiagonalOperator
,
create_power_operator
from
nifty.library
import
WienerFilterCurvature
#import plotly.offline as pl
#import plotly.graph_objs as go
from
mpi4py
import
MPI
comm
=
MPI
.
COMM_WORLD
rank
=
comm
.
rank
if
__name__
==
"__main__"
:
distribution_strategy
=
'fftw'
# Setting up physical constants
# total length of Interval or Volume the field lives on, e.g. in meters
L
=
2.
# typical distance over which the field is correlated (in same unit as L)
correlation_length
=
0.1
# variance of field in position space sqrt(<|s_x|^2>) (in unit of s)
field_variance
=
2.
# smoothing length of response (in same unit as L)
response_sigma
=
0.1
# defining resolution (pixels per dimension)
N_pixels
=
512
# Setting up derived constants
k_0
=
1.
/
correlation_length
# note that field_variance**2 = a*k_0/4. for this analytic form of power
# spectrum
a
=
field_variance
**
2
/
k_0
*
4.
pow_spec
=
(
lambda
k
:
a
/
(
1
+
k
/
k_0
)
**
4
)
pixel_length
=
L
/
N_pixels
# Setting up the geometry
s_space
=
RGSpace
([
N_pixels
,
N_pixels
],
distances
=
pixel_length
)
fft
=
FFTOperator
(
s_space
,
domain_dtype
=
np
.
float
,
target_dtype
=
np
.
complex
)
h_space
=
fft
.
target
[
0
]
inverse_fft
=
FFTOperator
(
h_space
,
target
=
s_space
,
domain_dtype
=
np
.
complex
,
target_dtype
=
np
.
float
)
p_space
=
PowerSpace
(
h_space
,
distribution_strategy
=
distribution_strategy
)
# Creating the mock data
S
=
create_power_operator
(
h_space
,
power_spectrum
=
pow_spec
,
distribution_strategy
=
distribution_strategy
)
sp
=
Field
(
p_space
,
val
=
pow_spec
,
distribution_strategy
=
distribution_strategy
)
sh
=
sp
.
power_synthesize
(
real_signal
=
True
)
ss
=
fft
.
inverse_times
(
sh
)
R
=
SmoothingOperator
(
s_space
,
sigma
=
response_sigma
)
R_harmonic
=
ComposedOperator
([
inverse_fft
,
R
],
default_spaces
=
[
0
,
0
])
signal_to_noise
=
1
N
=
DiagonalOperator
(
s_space
,
diagonal
=
ss
.
var
()
/
signal_to_noise
,
bare
=
True
)
n
=
Field
.
from_random
(
domain
=
s_space
,
random_type
=
'normal'
,
std
=
ss
.
std
()
/
np
.
sqrt
(
signal_to_noise
),
mean
=
0
)
d
=
R
(
ss
)
+
n
# Wiener filter
j
=
R_harmonic
.
adjoint_times
(
N
.
inverse_times
(
d
))
wiener_curvature
=
WienerFilterCurvature
(
S
=
S
,
N
=
N
,
R
=
R_harmonic
)
m
=
wiener_curvature
.
inverse_times
(
j
)
m_s
=
inverse_fft
(
m
)
demos/wiener_filter_via_curvature.py
0 → 100644
View file @
dd87b630
import
numpy
as
np
from
nifty
import
RGSpace
,
PowerSpace
,
Field
,
FFTOperator
,
ComposedOperator
,
\
DiagonalOperator
,
ResponseOperator
,
plotting
,
\
create_power_operator
from
nifty.library
import
WienerFilterCurvature
if
__name__
==
"__main__"
:
distribution_strategy
=
'not'
# Setting up variable parameters
# Typical distance over which the field is correlated
correlation_length
=
0.01
# Variance of field in position space sqrt(<|s_x|^2>)
field_variance
=
2.
# smoothing length of response (in same unit as L)
response_sigma
=
0.1
# The signal to noise ratio
signal_to_noise
=
0.7
# note that field_variance**2 = a*k_0/4. for this analytic form of power
# spectrum
def
power_spectrum
(
k
):
a
=
4
*
correlation_length
*
field_variance
**
2
return
a
/
(
1
+
k
*
correlation_length
)
**
4
# Setting up the geometry
# Total side-length of the domain
L
=
2.
# Grid resolution (pixels per axis)
N_pixels
=
512
signal_space
=
RGSpace
([
N_pixels
,
N_pixels
],
distances
=
L
/
N_pixels
)
harmonic_space
=
FFTOperator
.
get_default_codomain
(
signal_space
)
fft
=
FFTOperator
(
harmonic_space
,
target
=
signal_space
,
domain_dtype
=
np
.
complex
,
target_dtype
=
np
.
float
)
power_space
=
PowerSpace
(
harmonic_space
,
distribution_strategy
=
distribution_strategy
)
# Creating the mock data
S
=
create_power_operator
(
harmonic_space
,
power_spectrum
=
power_spectrum
,
distribution_strategy
=
distribution_strategy
)
mock_power
=
Field
(
power_space
,
val
=
power_spectrum
,
distribution_strategy
=
distribution_strategy
)
np
.
random
.
seed
(
43
)
mock_harmonic
=
mock_power
.
power_synthesize
(
real_signal
=
True
)
mock_signal
=
fft
(
mock_harmonic
)
R
=
ResponseOperator
(
signal_space
,
sigma
=
(
response_sigma
,))
data_domain
=
R
.
target
[
0
]
R_harmonic
=
ComposedOperator
([
fft
,
R
],
default_spaces
=
[
0
,
0
])
N
=
DiagonalOperator
(
data_domain
,
diagonal
=
mock_signal
.
var
()
/
signal_to_noise
,
bare
=
True
)
noise
=
Field
.
from_random
(
domain
=
data_domain
,
random_type
=
'normal'
,
std
=
mock_signal
.
std
()
/
np
.
sqrt
(
signal_to_noise
),
mean
=
0
)
data
=
R
(
mock_signal
)
+
noise
# Wiener filter
j
=
R_harmonic
.
adjoint_times
(
N
.
inverse_times
(
data
))
wiener_curvature
=
WienerFilterCurvature
(
S
=
S
,
N
=
N
,
R
=
R_harmonic
)
m
=
wiener_curvature
.
inverse_times
(
j
)
m_s
=
fft
(
m
)
plotter
=
plotting
.
RG2DPlotter
()
plotter
.
title
=
'mock_signal.html'
;
plotter
(
mock_signal
)
plotter
.
title
=
'data.html'
plotter
(
Field
(
signal_space
,
val
=
data
.
val
.
get_full_data
().
reshape
(
signal_space
.
shape
)))
plotter
.
title
=
'map.html'
;
plotter
(
m_s
)
\ No newline at end of file
demos/wiener_filter_
advanced
.py
→
demos/wiener_filter_
via_hamiltonian
.py
View file @
dd87b630
File moved
nifty/operators/smoothing_operator/fft_smoothing_operator.py
View file @
dd87b630
...
...
@@ -21,7 +21,6 @@ class FFTSmoothingOperator(SmoothingOperator):
# transform to the (global-)default codomain and perform all remaining
# steps therein
transformator
=
self
.
_get_transformator
(
x
.
dtype
)
transformed_x
=
transformator
(
x
,
spaces
=
spaces
)
codomain
=
transformed_x
.
domain
[
spaces
[
0
]]
coaxes
=
transformed_x
.
domain_axes
[
spaces
[
0
]]
...
...
nifty/plotting/plotter/healpix_plotter.py
View file @
dd87b630
...
...
@@ -21,3 +21,6 @@ class HealpixPlotter(PlotterBase):
def
_initialize_figure
(
self
):
return
Figure2D
(
plots
=
None
)
def
_parse_data
(
self
,
data
,
field
,
spaces
):
return
data
nifty/spaces/lm_space/lm_space.py
View file @
dd87b630
...
...
@@ -22,7 +22,7 @@ import numpy as np
from
nifty.spaces.space
import
Space
from
d2o
import
arange
from
d2o
import
arange
,
distributed_data_object
class
LMSpace
(
Space
):
...
...
@@ -138,8 +138,10 @@ class LMSpace(Space):
for
l
in
range
(
1
,
lmax
+
1
):
ldist
[
idx
:
idx
+
2
*
(
lmax
+
1
-
l
)]
=
tmp
[
2
*
l
:]
idx
+=
2
*
(
lmax
+
1
-
l
)
dists
=
arange
(
start
=
0
,
stop
=
self
.
shape
[
0
],
distribution_strategy
=
distribution_strategy
)
dists
=
distributed_data_object
(
global_shape
=
self
.
shape
,
dtype
=
np
.
float
,
distribution_strategy
=
distribution_strategy
)
dists
.
set_local_data
(
ldist
)
return
dists
...
...
Write
Preview
Supports
Markdown
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