Commit 613577ae authored by Russell Luke's avatar Russell Luke
Browse files

Merge branch 'master' into 'release'

Master

See merge request !3
parents 17239173 cb0297bc
......@@ -98,11 +98,13 @@ don't forget to type
exit()
1''. For this example, we already done the hard work for you. This time we'll
run the ProxToolbox from using a "DRprotein_demo.py" script. You can peek inside it
run the ProxToolbox from using a "DRprotein_demo.py" script inside the TestSuite
subdirectory. You can peek inside it
by openning it by your favourite editor, or even simpler type
cat DRprotein_demo.py
from the prompt in the TestSuite directory.
See how its got all the algorithmic parameters already in the file. Pretty neat huh?
2''. Ok. Time to run the script! It going to take a bit longer... actually a lot longer
......@@ -110,6 +112,7 @@ than the first two examples (hours maybe). We can work on speeding it up later.
python3 DRprotein_demo.py
(again, from in the TestSuite subdirectory).
You should see some (boring) information about the problem instance, etc. While the algorithm
is running you a report of its progress even 10 iterations, and will stop when the relative
error drops below 10e-5.
......
# a hack to get python to look in directories one branch up
import sys
sys.path.append("..")
#Input, Output and Logfiles
problemName = '1PTQ'
inputFile = 'inputdata/1PTQ.pdb' #in InputFiles directory
inputFile = '../InputData/Protein/1PTQ.pdb' #in InputFiles directory
outputFile = '1PTQrecon.pdb'
#Parameters used by all variants of the algorithm
......@@ -38,4 +42,4 @@ for rep in range(replications):
#problem.saveJPG()
problem.show()
del problem
del problem
\ No newline at end of file
List of files in ProxMatlab that have changed since Ziehe&Nahme version of ProxPython:
README
Algorithms:
ADMMPlus.m
DRlPHeBIE.m
PALM.m
RCAAR.m
ADMMPlus.m.old
HPR.m
RAAR.m
AP.m
Nesterov_PALM.m
RAAR_PALM.m
Wirttinger_grad.m
ProxOperators(top directory):
P_RN.m
P_amp.m
P_Diag.m
P_parallel_hyperplane.m
P_S.m
P_SP.m
P_parallel_scaled_hyperplane.m
P_S_real.m
Utilities:
FFT.m
IFFT.m
RANDsrc.m
drivers:
main_ProxToolbox.m
Custom.m
Ptychography/Ptychography_data:
Ptychography_data_processor.m
Sparsity
Phase
Misc
CT
Sudoku/Custom.m
Sensor_location
# -*- coding: utf-8 -*-
new_config = {
#We start very general.
#What type of problem is being solved? Classification
#is according to the geometry: Affine, Cone, Convex,
#Phase, Affine-sparsity, Nonlinear-sparsity, Sudoku'''
'problem_family' : 'Phase',
#==========================================
#Problem parameters
#==========================================
#What is the name of the data file?'''
'data_filename' : 'JWST_data_processor',
#What type of object are we working with?
#Options are: 'phase', 'real', 'nonnegative', 'complex' '''
'object' : 'complex',
#What type of constraints do we have?
#Options are: 'support only', 'real and support', 'nonnegative and support',
# 'amplitude only', 'sparse real', 'sparse complex', and 'hybrid' '''
'constraint' : 'amplitude only',
#What type of measurements are we working with?
#Options are: 'single diffraction', 'diversity diffraction',
# 'ptychography', and 'complex' '''
'experiment' : 'JWST',
#Next we move to things that most of our users will know
#better than we will. Some of these may be overwritten in the
#data processor file which the user will most likely write.
#Are the measurements in the far field or near field?
#Options are: 'far field' or 'near field' '''
'distance' : 'far field',
# What are the dimensions of the measurements?
'Nx' : 128,
'Ny' : 128,
'Nz' : 1, # do not formulate in the product space
'dim' : 4, #size of the product space
#if 'distance' =='near field':
# 'fresnel_nr' : 1*2*pi*config['Nx'],
# 'use_farfield_formula' : 0,
#else:
'fresnel_nr' : 0, #1*2*pi*prbl.Nx;
'use_farfield_formula' : 1,
# What are the noise characteristics (Poisson or Gaussian or none)?
'noise' : 'none', #'Poisson',
#==========================================
# Algorithm parameters
#==========================================
# Now set some algorithm parameters that the user should be
# able to control (without too much damage)'''
# Algorithm:
'algorithm' : 'RAAR', #'Accelerated_AP_product_space';
'numruns' : 1, # the only time this parameter will
# be different than 1 is when we are
# benchmarking...not something a normal user
# would be doing.
# The following are parameters specific to RAAR, HPR, and HAAR that the
# user should be able to set/modify. Surely
# there will be other algorithm specific parameters that a user might
# want to play with. Don't know how best
# to do this. Thinking of a GUI interface, we could hard code all the
# parameters the user might encounter and have the menu options change
# depending on the value of the prbl.method field.
# do different things depending on the chosen algorithm:
#maximum number of iterations and tolerances
'MAXIT' : 500,
'TOL' : 1e-8,
# relaxaton parameters in RAAR, HPR and HAAR
'beta0' : 1.0, # starting relaxation prameter (only used with
# HAAR, HPR and RAAR)
'beta_max' : 1.0, # maximum relaxation prameter (only used with
# HAAR, RAAR, and HPR)
'beta_switch' : 20, # iteration at which beta moves from beta_0 -> beta_max
# parameter for the data regularization
# need to discuss how/whether the user should
# put in information about the noise
'data_ball' : 1e-15,
# the above is the percentage of the gap
# between the measured data and the
# initial guess satisfying the
# qualitative constraints. For a number
# very close to one, the gap is not expected
# to improve much. For a number closer to 0
# the gap is expected to improve a lot.
# Ultimately the size of the gap depends
# on the inconsistency of the measurement model
# with the qualitative constraints.
#==========================================
# parameters for plotting and diagnostics
#==========================================
'verbose' : 1, # options are 0 or 1
'graphics' : 1, # whether or not to display figures, options are 0 or 1.
# default is 1.
'anim' : 2, # whether or not to disaply ``real time" reconstructions
# options are 0=no, 1=yes, 2=make a movie
# default is 1.
'graphics_display' : 'JWST_graphics' # unless specified, a default
# plotting subroutine will generate
# the graphics. Otherwise, the user
# can write their own plotting subroutine
#======================================================================
# Technical/software specific parameters
#======================================================================
# Given the parameter values above, the following technical/algorithmic
# parameters are automatically set. The user does not need to know
# about these details, and so probably these parameters should be set in
# a module one level below this one.
}
......@@ -8,5 +8,5 @@ The "ProxOperators"-module contains all proximal operators that are already impl
"""
from .proxoperators import *
__all__ = ["P_diag","P_parallel","magproj"]
from .breman import *
__all__ = ["P_diag","P_parallel","magproj","bregProj"]
# -*- coding: utf-8 -*-
"""
Created on Sun May 7 12:18:52 2017
@author: Theodor
"""
import numpy as np
__all__=["bregProj"]
def J(x,p,norm=np.linalg.norm):
"""
Duality function
Parameters
----------
x : array_like - The point which dual is searched
p : value - The .....
"""
if norm(x,ord=2)>0:
return x*norm(x,ord=2)**(p-2)
else:
return x
def quasi_newton(args):
"""
Quasi Newtonian method for finding roots.
starts at x=1 and assumes as gradient the secant going througn 1,f(1)
and 0,f(0). Afterwards it assumes as gradient the secant going though the
points defined by the last and current iteration.
Parameters
----------
args : dict - The dictionary containing all of the below
f : function - The function for which we want to find the root
epsilon : float - The wished precision
maxiter : int - The maximum allowed iterationcount
Returns
-------
float - The root
"""
std_args={"epsilon":10**(-9),"maxiter":10**3}
std_args.update(args)
f=args["f"]
epsilon=args["epsilon"]
maxiter=args["maxiter"]
last_x=0
f_x=f(last_x)
f_prime_guess=f(1)-f_x
last_x=1
next_x=last_x-f_x/f_prime_guess
f_x=f(last_x)
iteration=0
while iteration<maxiter and abs(next_x-last_x)>epsilon:
f_nx=f(next_x)
f_prime_guess=(f_x-f_nx)/(last_x-next_x)
tmp=next_x
next_x=last_x-f_x/f_prime_guess
last_x=tmp
f_x=f_nx
iteration+=1
print(abs(last_x-next_x))
return next_x
def bregProj(args):
"""
Bregman projection operator onto a hyperplane given defined by
all x such that u*x=a.
The bregman projection is to the standard bregman function ||.||**p
Parameters
----------
args : dict - The dictionary containing all of the below
u : array_like - The orthogonal vector defining the hyperplane
a : value - The offset needed to define the hyperplane
p : value - The parameter needed to define the bregman function
x : array_like - The point which we want to project onto the hyperplane
method : function- The rootfinding algorithm
norm : function - The underlying norm
Returns
-------
array_like - The projection
"""
#parse the arguments
std_args={"x":np.zeros(1),
"u":np.zeros(1),
"a":0,
"p":2,
"method":quasi_newton,
"norm":np.linalg.norm}
std_args.update(args)
x=std_args["x"]
u=std_args["u"]
a=std_args["a"]
p=std_args["p"]
method=std_args["method"]
norm=std_args["norm"]
hp= lambda t: -np.dot(u,J(J(x,p)-t*u,p/(p-1)))+a
if not np.any(u):
t0=0
elif a==0:
t0=np.dot(J(x,p),u)/norm(u,ord=2)**2
else:
arguments=dict(f=hp)
t0=method(arguments)
return J(J(x,p)-t0*u,p/(p-1))
Markdown is supported
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