Skip to content
GitLab
Projects
Groups
Snippets
/
Help
Help
Support
Community forum
Keyboard shortcuts
?
Submit feedback
Contribute to GitLab
Sign in
Toggle navigation
Menu
Open sidebar
ift
starblade
Commits
c8ff7067
Commit
c8ff7067
authored
Apr 10, 2018
by
Jakob Knollmueller
Browse files
documentation, cleanup
parent
8563ad22
Changes
7
Hide whitespace changes
Inline
Side-by-side
1d_separation.py
View file @
c8ff7067
from
point_separation
import
build_problem
,
problem_iteration
from
sugar
import
build_problem
,
problem_iteration
import
nifty4
as
ift
import
numpy
as
np
from
matplotlib
import
rc
...
...
demo.py
0 → 100644
View file @
c8ff7067
from
sugar
import
build_starblade
,
starblade_iteration
from
matplotlib
import
pyplot
as
plt
from
astropy.io
import
fits
import
numpy
as
np
if
__name__
==
'__main__'
:
#specifying location of the input file:
path
=
'hst_05195_01_wfpc2_f702w_pc_sci.fits'
data
=
fits
.
open
(
path
)[
1
].
data
data
=
data
.
clip
(
min
=
0.001
)
data
=
np
.
ndarray
.
astype
(
data
,
float
)
vmin
=
np
.
log
(
data
.
min
()
+
0.01
)
vmax
=
np
.
log
(
data
.
max
())
alpha
=
1.25
Starblade
=
build_starblade
(
data
,
alpha
=
alpha
)
for
i
in
range
(
10
):
Starblade
=
starblade_iteration
(
Starblade
)
#plotting on logarithmic scale
plt
.
imsave
(
'diffuse_component.png'
,
Starblade
.
s
.
val
,
vmin
=
vmin
,
vmax
=
vmax
)
plt
.
imsave
(
'pointlike_component.png'
,
Starblade
.
u
.
val
,
vmin
=
vmin
,
vmax
=
vmax
)
\ No newline at end of file
gui_app.py
View file @
c8ff7067
...
...
@@ -2,7 +2,7 @@ import matplotlib
matplotlib
.
use
(
'agg'
)
# matplotlib.use('module://kivy.garden.matplotlib.backend_kivy')
from
point_separation
import
build_multi_problem
,
multi_problem_iteration
,
load_data
from
sugar
import
build_multi_problem
,
multi_problem_iteration
,
load_data
from
kivy.app
import
App
from
kivy.uix.widget
import
Widget
...
...
hubble_separation.py
View file @
c8ff7067
from
point_separation
import
build_problem
,
problem_iteration
,
load_data
from
sugar
import
build_problem
,
problem_iteration
,
load_data
from
nifty4
import
*
import
numpy
as
np
from
matplotlib
import
rc
...
...
@@ -8,14 +8,20 @@ from matplotlib import pyplot as plt
from
matplotlib.colors
import
LogNorm
from
mpl_toolkits.axes_grid1
import
make_axes_locatable
from
mpl_toolkits.axes_grid1
import
AxesGrid
from
astropy.io
import
fits
np
.
random
.
seed
(
42
)
if
__name__
==
'__main__'
:
path
=
'hst_05195_01_wfpc2_f702w_pc_sci.fits'
data
=
load_data
(
path
)
alpha
=
1.3
data
=
fits
.
open
(
path
)[
1
].
data
data
=
data
.
clip
(
min
=
0.001
)
data
=
np
.
ndarray
.
astype
(
data
,
float
)
alpha
=
1.25
...
...
rgb_separation
.py
→
multichannel_demo
.py
View file @
c8ff7067
from
point_separation
import
build_multi_
problem
,
multi_problem
_iteration
from
sugar
import
build_multi_
starblade
,
multi_starblade
_iteration
from
matplotlib
import
pyplot
as
plt
import
numpy
as
np
if
__name__
==
'__main__'
:
# data = plt.imread('eso1242a.jpg')
data
=
plt
.
imread
(
'10Keso1242a.tif'
)
# data = plt.imread('10Keso1242a.tif')
data
=
plt
.
imread
(
'eso1242a.jpg'
)
data
=
data
.
astype
(
float
)
data
=
data
.
clip
(
0.0001
)
energy_list
=
build_multi_problem
(
data
,
1.2
)
for
i
in
range
(
10
):
energy_list
=
multi_problem_iteration
(
energy_list
)
alpha
=
1.25
MultiStarblade
=
build_multi_starblade
(
data
,
alpha
)
for
i
in
range
(
2
):
MultiStarblade
=
multi_starblade_iteration
(
MultiStarblade
,
multiprocessing
=
True
)
#plotting a three channel RGB image in each iteration
diffuse
=
np
.
empty_like
(
data
)
point
=
np
.
empty_like
(
data
)
for
i
in
range
(
len
(
energy_list
)):
diffuse
[...,
i
]
=
np
.
exp
(
energy_list
[
i
].
s
.
val
)
point
[...,
i
]
=
np
.
exp
(
energy_list
[
i
].
u
.
val
)
for
i
in
range
(
len
(
MultiStarblade
)):
diffuse
[...,
i
]
=
np
.
exp
(
MultiStarblade
[
i
].
s
.
val
)
point
[...,
i
]
=
np
.
exp
(
MultiStarblade
[
i
].
u
.
val
)
plt
.
imsave
(
'rgb_diffuse.jpg'
,
diffuse
/
255.
)
plt
.
imsave
(
'rgb_point.jpg'
,
point
/
255.
)
...
...
separation_energy.py
View file @
c8ff7067
from
nifty4
import
Energy
,
Field
,
log
,
exp
,
DiagonalOperator
from
nifty4.library
import
WienerFilterCurvature
from
nifty4.library.nonlinearities
import
PositiveTanh
class
SeparationEnergy
(
Energy
):
class
StarbladeEnergy
(
Energy
):
"""The Energy for the starblade problem.
It implements the Information Hamiltonian of the separation of d
Parameters
----------
position : Field
The current position of the separation.
parameters : Dictionary
Dictionary containing all relevant quantities for the inference,
data : Field
The image data.
alpha : Field
Slope parameter of the point-source prior
q : Field
Cutoff parameter of the point-source prior
correlation : Field
A field in the Fourier space which encodes the diagonal of the prior
correlation structure of the diffuse component
FFT : FFTOperator
An operator performing the Fourier transform
inverter : ConjugateGradient
the minimization strategy to use for operator inversion
"""
def
__init__
(
self
,
position
,
parameters
):
x
=
position
.
val
.
clip
(
-
9
,
9
)
position
=
Field
(
position
.
domain
,
val
=
x
)
super
(
S
eparation
Energy
,
self
).
__init__
(
position
=
position
)
super
(
S
tarblade
Energy
,
self
).
__init__
(
position
=
position
)
self
.
parameters
=
parameters
self
.
inverter
=
parameters
[
'inverter'
]
...
...
@@ -17,8 +42,7 @@ class SeparationEnergy(Energy):
self
.
correlation
=
parameters
[
'correlation'
]
self
.
alpha
=
parameters
[
'alpha'
]
self
.
q
=
parameters
[
'q'
]
pos_tanh
=
parameters
[
'pos_tanh'
]
pos_tanh
=
PositiveTanh
()
self
.
S
=
self
.
FFT
*
self
.
correlation
*
self
.
FFT
.
adjoint
self
.
a
=
pos_tanh
(
self
.
position
)
self
.
a_p
=
pos_tanh
.
derivative
(
self
.
position
)
...
...
point_separation
.py
→
sugar
.py
View file @
c8ff7067
import
nifty4
as
ift
import
numpy
as
np
from
matplotlib
import
pyplot
as
plt
from
astropy.io
import
fits
from
separation_energy
import
SeparationEnergy
from
nifty4.library.nonlinearities
import
PositiveTanh
from
multiprocessing
import
Pool
from
separation_energy
import
StarbladeEnergy
def
load_data
(
path
):
if
path
[
-
5
:]
==
'.fits'
:
data
=
fits
.
open
(
path
)[
1
].
data
else
:
data
=
plt
.
imread
(
path
)[:,:,
0
]
data
=
data
.
clip
(
min
=
0.001
)
data
=
np
.
ndarray
.
astype
(
data
,
float
)
return
data
def
build_starblade
(
data
,
alpha
=
1.5
,
q
=
1e-40
,
cg_iterations
=
500
):
""" Setting up the StarbladeEnergy for the given data and parameters
Parameters
----------
data : array
The data in a numpy array
alpha : float
The slope parameter of the point source prior (default: 1.5).
q : float
The cutoff parameter of the point source prior (default: 1e-40).
cg_iterations : int
Maximum number of conjugate gradient iterations for numerical operator inversion (default: 500).
"""
def
build_problem
(
data
,
alpha
):
s_space
=
ift
.
RGSpace
(
data
.
shape
,
distances
=
len
(
data
.
shape
)
*
[
1
])
h_space
=
s_space
.
get_default_codomain
()
data
=
ift
.
Field
(
s_space
,
val
=
data
)
FFT
=
ift
.
FFTOperator
(
h_space
)
FFT
=
ift
.
FFTOperator
(
h_space
,
target
=
s_space
)
binbounds
=
ift
.
PowerSpace
.
useful_binbounds
(
h_space
,
logarithmic
=
False
)
p_space
=
ift
.
PowerSpace
(
h_space
,
binbounds
=
binbounds
)
initial_spectrum
=
ift
.
power_analyze
(
FFT
.
inverse_times
(
ift
.
log
(
data
)),
binbounds
=
p_space
.
binbounds
)
initial_spectrum
/=
(
p_space
.
k_lengths
+
1.
)
**
2
initial_correlation
=
ift
.
create_power_operator
(
h_space
,
initial_spectrum
)
initial_x
=
ift
.
Field
(
s_space
,
val
=-
1.
)
alpha
=
ift
.
Field
(
s_space
,
val
=
alpha
)
q
=
ift
.
Field
(
s_space
,
val
=
10e-40
)
pos_tanh
=
PositiveTanh
()
ICI
=
ift
.
GradientNormController
(
iteration_limit
=
500
,
q
=
ift
.
Field
(
s_space
,
val
=
q
)
ICI
=
ift
.
GradientNormController
(
iteration_limit
=
cg_iterations
,
tol_abs_gradnorm
=
1e-5
)
inverter
=
ift
.
ConjugateGradient
(
controller
=
ICI
)
parameters
=
dict
(
data
=
data
,
correlation
=
initial_correlation
,
alpha
=
alpha
,
q
=
q
,
inverter
=
inverter
,
FFT
=
FFT
,
pos_tanh
=
pos_tanh
)
separationEnergy
=
Separation
Energy
(
position
=
initial_x
,
parameters
=
parameters
)
return
separationEnergy
inverter
=
inverter
,
FFT
=
FFT
)
Starblade
=
Starblade
Energy
(
position
=
initial_x
,
parameters
=
parameters
)
return
Starblade
def
problem_iteration
(
energy
,
iterations
=
3
):
controller
=
ift
.
GradientNormController
(
name
=
"test1"
,
tol_abs_gradnorm
=
0.00000001
,
iteration_limit
=
iterations
)
def
starblade_iteration
(
starblade
,
iterations
=
3
):
""" Performing one Newton minimization step
Parameters
----------
starblade : StarbladeEnergy
An instance of an Starblade Energy
iterations : int
The number of steps with the Newton scheme (default: 3).
"""
controller
=
ift
.
GradientNormController
(
name
=
"Newton"
,
tol_abs_gradnorm
=
1e-8
,
iteration_limit
=
iterations
)
minimizer
=
ift
.
RelaxedNewton
(
controller
=
controller
)
energy
,
convergence
=
minimizer
(
energy
)
energy
,
convergence
=
minimizer
(
starblade
)
new_position
=
energy
.
position
h_space
=
energy
.
correlation
.
domain
[
0
]
FFT
=
energy
.
FFT
binbounds
=
ift
.
PowerSpace
.
useful_binbounds
(
h_space
,
logarithmic
=
False
)
new_power
=
ift
.
power_analyze
(
FFT
.
inverse_times
(
energy
.
s
),
binbounds
=
binbounds
)
# new_power /= (new_power.domain[0].k_lengths+1.)**2
new_correlation
=
ift
.
create_power_operator
(
h_space
,
new_power
)
new_parameters
=
energy
.
parameters
new_parameters
[
'correlation'
]
=
new_correlation
n
ew
_energy
=
Separation
Energy
(
new_position
,
new_parameters
)
return
n
ew
_energy
#
new_parameters['correlation'] = new_correlation
N
ew
Starblade
=
Starblade
Energy
(
new_position
,
new_parameters
)
return
N
ew
Starblade
def
build_multi_problem
(
data
,
alpha
):
energy_list
=
[]
def
build_multi_starblade
(
data
,
alpha
=
1.5
,
q
=
1e-40
,
cg_iterations
=
500
):
""" Builds a list of StarbladeEnergies for the given multi-channel dataset
Parameters
----------
data : array
The data in a numpy array of the multi-channel dataset with channel axis data[-1].
alpha : float
The slope parameter of the point source prior (default: 1.5).
q : float
The cutoff parameter of the point source prior (default: 1e-40).
cg_iterations : int
Maximum number of conjugate gradient iterations for numerical operator inversion (default: 500).
"""
MultiStarblade
=
[]
for
i
in
range
(
data
.
shape
[
-
1
]):
energy
=
build_problem
(
data
[...,
i
],
alpha
)
energy_list
.
append
(
energy
)
return
energy_list
starblade
=
build_starblade
(
data
[...,
i
],
alpha
=
alpha
,
q
=
q
,
cg_iterations
=
cg_iterations
)
MultiStarblade
.
append
(
starblade
)
return
MultiStarblade
def
multi_problem_iteration
(
energy_list
):
new_energy
=
[]
for
energy
in
energy_list
:
new_energy
.
append
(
problem_iteration
(
energy
))
return
new_energy
def
multi_starblade_iteration
(
MultiStarblade
,
multiprocessing
=
False
):
""" Performing one Newton minimization step for all entries of the MultiStarblade list.
Parameters
----------
MultiStarblade : list of StarbladeEnergy
A list of instances of an Starblade Energy
iterations : int
The number of steps with the Newton scheme (default: 3).
"""
if
multiprocessing
:
NewStarblades
=
list
(
Pool
(
processes
=
3
).
map
(
starblade_iteration
,
MultiStarblade
))
else
:
NewStarblades
=
[]
for
starblade
in
MultiStarblade
:
NewStarblades
.
append
(
starblade_iteration
(
starblade
))
return
NewStarblades
if
__name__
==
'__main__'
:
pass
...
...
Write
Preview
Supports
Markdown
0%
Try again
or
attach a new 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