Commit 4ca48457 authored by jansen31's avatar jansen31
Browse files

P_Sparsity_real created as child of P_Sparsity

P_Sparsity updated to work correctly with numpy array indexing methods
OrbitalTomog/phase updated to accept sparse real and sparse complex
Created jupyter notebook for quick and easy code running
parent 1eff1f4f
This diff is collapsed.
......@@ -146,7 +146,7 @@ new_config = {
# options are 0=no, 1=yes, 2=make a movie
# default is 1.
'graphics_display': 'Phase_graphics', # unless specified, a default
'dataprocessor_plotting':True
'dataprocessor_plotting': True
# plotting subroutine will generate
# the graphics. Otherwise, the user
# can write their own plotting subroutine
......
......@@ -62,6 +62,10 @@ class Phase(Problem):
used_proxoperators[0] = 'P_SP'
elif self.config['constraint'] == 'amplitude only':
used_proxoperators[0] = 'P_amp'
elif self.config['constraint'] == 'sparse real':
used_proxoperators[0] = 'P_Sparsity_real'
elif self.config['constraint'] == 'sparse complex':
used_proxoperators[0] = 'P_Sparsity'
# Projector 2 (k / Fourier space)
# used_proxoperators[1] = 'P_M' # 'Approx_P_FreFra_Poisson'
......@@ -156,5 +160,5 @@ class Phase(Problem):
print("Calculation time:")
print(self.elapsed_time)
_graphics = getattr(Graphics, self.config['graphics_display'])
print("Retrieved graphics function")
# print("Retrieved graphics function")
_graphics(self.config, self.output)
......@@ -242,7 +242,7 @@
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.6.8"
"version": "3.7.3"
}
},
"nbformat": 4,
......
from numpy import zeros_like
import numpy as np
from .proxoperators import ProxOperator
......@@ -30,6 +31,23 @@ class P_Sparsity(ProxOperator):
'sparsity_parameter' : int
"""
self.sparsity_parameter = config['sparsity_parameter']
self.ny = config['Ny']
if self.sparsity_parameter > 30:
def value_selection(original, indices, sparsity_parameter):
idx_for_threshold = divmod(indices[-sparsity_parameter], self.ny)
threshold_val = abs(original[idx_for_threshold].get())
return (abs(original) >= threshold_val) * original
else:
def value_selection(original, indices, sparsity_parameter):
out = zeros_like(original)
hits = indices[-sparsity_parameter:].get()
hit_idx = [divmod(hit, self.ny) for hit in hits]
for _idx in hit_idx:
out[_idx[0], _idx[1]] = original[_idx[0], _idx[1]]
return out
self.value_selection = value_selection
def work(self, u):
"""
......@@ -43,8 +61,16 @@ class P_Sparsity(ProxOperator):
"""
sparsity_parameter = self.sparsity_parameter
p_Sparsity = 0 * u
I = np.argsort(abs(u)) # gives indices of sorted array in ascending order
I = np.flip(I, axis=0) # reverses the array so that I gives the indices of sorted array in descending order
p_Sparsity[I[0:sparsity_parameter]] = u[I[0:sparsity_parameter]]
sorting = np.argsort(abs(u), axis=None) # gives indices of sorted array in ascending order
indices = np.asarray([divmod(sorting[i], self.ny) for i in range(-1 * sparsity_parameter, 0)])
p_Sparsity[indices[:, 0], indices[:, 1]] = u[indices[:, 0], indices[:, 1]]
return p_Sparsity
class P_Sparsity_real(P_Sparsity):
def __init__(self, config):
super(P_Sparsity_real, self).__init__(config)
def work(self, u):
return super(P_Sparsity_real, self).work(u.real)
from proxoperators import Proxoperator
import numpy as np
#original matlab comment
# original matlab comment
# P_nonneg_sparsity.m
# written on August 20, 2012 by
# Russell Luke
......@@ -19,27 +21,27 @@ import numpy as np
# USAGE: p_Sparsity = P_Sparsity(input,u)
class P_nonneg_sparsity(ProxOperator):
def __init__(self, config):
"""
def __init__(self, config):
"""
Initialization
Parameters
----------
config : dict - Dictionary containing the problem configuration. It must contain the following mappings:
'm' : int -row-dim
'n' : int -Column-dim
'p' : int -Depth-dim
'q' : int - with [m,n,p,q]=size(u)
sparsity_parameter : int
"""
self.m = config['m']
self.n = config['n']
self.p = config['p']
self.q = config['q']
self.sparsity_parameter = config['sparsity_parameter']
'm' : int -row-dim
'n' : int -Column-dim
'p' : int -Depth-dim
'q' : int - with [m,n,p,q]=size(u)
sparsity_parameter : int
"""
self.m = config['m']
self.n = config['n']
self.p = config['p']
self.q = config['q']
self.sparsity_parameter = config['sparsity_parameter']
def work(self,u):
"""
def work(self, u):
"""
Applies the prox operator to some input data
Parameters
......@@ -50,13 +52,12 @@ class P_nonneg_sparsity(ProxOperator):
-------
p_Sparsity : array_like - Result of the operation
"""
sparsity_parameter = self.sparsity_parameter
tmp = max(0,np.real(u))
tmp = np.reshape(tmp,m*n*p*q)
I = np.argsort(tmp) #gives indices of sorted array in ascending order
I = np.flip(I,axis=0)#reverses the array so that I gives the indices of sorted array in descending order
p_Sparsity = np.zeros(m*n*p*q, dtype= tmp.dtype)
p_Sparsity[I[0:sparsity_parameter+1]] = tmp[I[0:sparsity_parameter+1]]
p_Sparsity = np.reshape(p_Sparsity,(m,n,p,q))
return p_Sparsity
sparsity_parameter = self.sparsity_parameter
tmp = max(0, np.real(u))
tmp = np.reshape(tmp, m * n * p * q)
I = np.argsort(tmp) # gives indices of sorted array in ascending order
I = np.flip(I, axis=0) # reverses the array so that I gives the indices of sorted array in descending order
p_Sparsity = np.zeros(m * n * p * q, dtype=tmp.dtype)
p_Sparsity[I[0:sparsity_parameter + 1]] = tmp[I[0:sparsity_parameter + 1]]
p_Sparsity = np.reshape(p_Sparsity, (m, n, p, q))
return p_Sparsity
......@@ -8,8 +8,18 @@ The "ProxOperators"-module contains all proximal operators that are already impl
"""
from .proxoperators import *
from .P_Sparsity import *
__all__ = ["P_diag","P_parallel","magproj", "Approx_P_JWST_Poisson", "P_amp", "P_SP", "Approx_PM_Gaussian",
"Approx_PM_Poisson", "P_S", "P_S_real", "P_sequential_hyperplane_odd", "P_sequential_hyperplane_even",
"P_parallel_hyperplane", "P_block_parallel_hyperplane", "P_block_sequential_hyperplane", "Q_Heau",
"P_Amod", "Approx_P_FreFra_Poisson"]
__all__ = ["P_diag",
"P_parallel",
"magproj",
"Approx_P_JWST_Poisson",
"P_amp",
"Approx_PM_Gaussian", "Approx_PM_Poisson", "Approx_P_FreFra_Poisson",
"P_S", "P_S_real", "P_SP",
'P_Sparsity', 'P_Sparsity_real',
"P_sequential_hyperplane_odd", "P_sequential_hyperplane_even",
"P_parallel_hyperplane", "P_block_parallel_hyperplane",
"P_block_sequential_hyperplane",
"Q_Heau",
"P_Amod"]
Supports Markdown
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