""" TODO : document this module """
import math
import numpy
import mpi4py.MPI as MPI
import h5py
import sys
myrank = MPI.COMM_WORLD.Get_rank()
def Check_must_have_keys(dizio, keys):
if( set(dizio.keys())!=set(keys) ):
B = keys
A = dizio.keys()
for a in A:
print a , a in B
print " ################ "
raise Exception, " set(dizio.keys())!=set(keys) "
def writeh5(map_complex, nome):
h5=h5py.File( nome+".h5",'w')
h5["/map/map_squareamplitude"]= map_complex.real*map_complex.real + map_complex.imag*map_complex.imag
h5["/map/map_phase"] = numpy.arctan2(map_complex.imag,map_complex.real)
h5["/map/map_real"] =map_complex.real
h5["/map/map_imag"] =map_complex.imag
h5.flush()
h5.close()
h5=None
def writeh5R(map_complex, nome):
h5=h5py.File( nome+".h5",'w')
h5["/map/map_real"] =map_complex
h5.flush()
h5.close()
h5=None
[docs]def main() :
""" TODO : document this function """
import getCpuSet
ncpus, mygpus=getCpuSet.getCpuSet()
print " process " , myrank, " mygpus = " , mygpus
print " process " , myrank, " ncpus= " ,ncpus
import ParametersR_module
if mygpus!=[]:
ParametersR_module.Parameters.MYGPU= mygpus[0]
else:
ParametersR_module.Parameters.MYGPU= -1
arrays_space_in_python = ParametersR_module.create_arrays_space_as_dictionary()
Check_must_have_keys( arrays_space_in_python, ["datas","solution"] )
datas = arrays_space_in_python["datas"]
#datas[:]=numpy.sqrt(abs(datas))
solution = arrays_space_in_python["solution"]
derived_pars=ParametersR_module.Parameters.get_derived_pars()
print " IL FOCUS EST " , derived_pars.focus
datas[:]=numpy.fft.ifftshift(datas)
writeh5R( datas, "datoiniziale2")
import Cspace
cspace=Cspace.Cspace( arrays_space_in_python , ParametersR_module.Parameters, derived_pars )
setup_number_of_real_arrays=3
setup_number_of_cmplx_arrays=6
cspace.initialise_gpu( setup_number_of_real_arrays, setup_number_of_cmplx_arrays )
realA_phase_a = 0
realA_field_s0 = 1
realA_field_support = 2
cmplxA_field_s =0
cmplxA_field_dumA = 1
cmplxA_obj_fl_new = 2
cmplxA_obj_fl_last = 3
cmplxA_obj_fl = 4
cmplxA_obj_mem = 5
params = ParametersR_module.Parameters
writeh5R( datas, "initial_field_s0")
cspace.initialise_realA ( realA_field_s0, datas )
data_dum=datas*0
cspace.initialise_realA( realA_phase_a, (0.0*math.pi *numpy.random.random([params.N_pxl,params.N_pxl]) ).astype(numpy.float32) )
if(0):
h5=h5py.File("phase0.h5",'r')
data1 = h5["/map/map_real"][:]
phase = data1
cspace.initialise_realA( realA_phase_a, phase )
h5=None
cspace.get_realA(data_dum , realA_phase_a)
writeh5R(data_dum , "realA_phase_a")
######### Support in direct space for full illumination ############
supp1=numpy.ones((params.N_pxl,params.N_pxl),dtype=numpy.float32)
X=numpy.linspace(-(derived_pars.pxl_FZP* params.N_pxl)/2.0, (derived_pars.pxl_FZP*params.N_pxl)/2.0, params.N_pxl) # da rimettere dopo confronto , endpoint=False)
test = X*X<(params.supp_radius)**2
print "derived_pars.pxl_FZP " , derived_pars.pxl_FZP
supporto_iniziale=((X*X +(X*X)[:,None]) < (params.supp_radius)**2).astype(numpy.float32)
cspace.get_realA(data_dum ,realA_field_support )
writeh5R(numpy.fft.fftshift(data_dum) , "support")
# supporto_iniziale= (supporto_iniziale*( (X*X +(X*X)[:,None]) > ( params.CS/2 )**2)).astype(numpy.float32)
####################################################################
cspace.initialise_realA( realA_field_support ,numpy.fft.ifftshift( supporto_iniziale) )
errori=[]
def from_detector_through_focus_to_support( params,derived_pars , target = None, source= None, dumA = None, focus=None):
# prima del backward lo zero dell' array e' il centro del detector
cspace.go_backward(dumA, source )
cspace.apply_fresnel_phase( dumA , -derived_pars.z_fd , derived_pars.size_focal_pixel , dumA )
# nel nel mezzo lo zero e' un angolo , ed e' invece il centro dell' array che e' il centro della zona
if focus is not None: cspace.copy( focus , dumA )
cspace.apply_fresnel_phase( dumA , -derived_pars.focus *params.OSA_prox_factor, derived_pars.size_focal_pixel , dumA )
cspace.go_backward(target, dumA )
# si ritorna qui nella disposizione zero=centro
def from_support_through_focus_to_detector( params,derived_pars , target = None, source= None, dumA = None, focus=None):
cspace.go_forward( dumA, source )
cspace.apply_fresnel_phase( dumA , derived_pars.focus*params.OSA_prox_factor ,derived_pars.size_focal_pixel , dumA )
if focus is not None: cspace.copy( focus , dumA )
cspace.apply_fresnel_phase( dumA , derived_pars.z_fd , derived_pars.size_focal_pixel , dumA )
cspace.go_forward( target , dumA )
renorm1 = 1.0/params.N_pxl/params.N_pxl
renorm=renorm1*renorm1
for p in range(params.cycles_TOT):
for n in range(params.cycles_HIO):
if p==0:
break
else:
# il target (se c' e' un target ) di tutte le operazioni e' sempre a sinistra
cspace.set_phase(cmplxA_field_s, realA_field_s0, realA_phase_a )
from_detector_through_focus_to_support( params,derived_pars , target = cmplxA_obj_fl_new, source= cmplxA_field_s, dumA = cmplxA_field_dumA)
cspace.Fienup( cmplxA_obj_fl_new , cmplxA_obj_fl_new, cmplxA_obj_fl_last , realA_field_support, params.Beta ) ## ----------------
cspace.copy(cmplxA_obj_fl_last, cmplxA_obj_fl_new )
from_support_through_focus_to_detector( params,derived_pars , target =cmplxA_field_s , source=cmplxA_obj_fl_new ,
dumA = cmplxA_field_dumA,
focus = cmplxA_obj_mem)
cspace.get_phase( realA_phase_a ,cmplxA_field_s )
cspace.apply_factor( cmplxA_field_s , cmplxA_field_s ,renorm )
newerr = cspace.error(cmplxA_field_s, realA_field_s0 )
print " fienup " , newerr
errori.append(newerr)
for l in range(params.cycles_ER):
# cspace.initialise_realA( realA_phase_a, phase )
# cspace.get_realA(data_dum , realA_phase_a )
# writeh5R(numpy.fft.fftshift(data_dum ), "phase")
cspace.set_phase(cmplxA_field_s, realA_field_s0, realA_phase_a )
from_detector_through_focus_to_support( params,derived_pars , target = cmplxA_obj_fl_new, source= cmplxA_field_s, dumA = cmplxA_field_dumA,
focus = cmplxA_obj_mem)
cspace.apply_support( cmplxA_obj_fl, cmplxA_obj_fl_new ,realA_field_support ) ## -----------------------------
# cspace.copy( cmplxA_obj_fl, cmplxA_obj_fl_new ) ;
from_support_through_focus_to_detector( params,derived_pars , target =cmplxA_field_s , source=cmplxA_obj_fl,
dumA = cmplxA_field_dumA,
focus = cmplxA_obj_mem)
# cspace.get_cmplxA( solution,cmplxA_obj_mem )
# writeh5(numpy.fft.fftshift(solution*renorm ), "focus")
cspace.apply_factor( cmplxA_field_s , cmplxA_field_s ,renorm )
cspace.get_phase( realA_phase_a ,cmplxA_field_s )
newerr = cspace.error(cmplxA_field_s, realA_field_s0 )
print " error reduction ", l, newerr
# cspace.get_realA(data_dum ,realA_field_s0 )
# writeh5R(numpy.fft.fftshift(data_dum) , "data")
# cspace.get_cmplxA( solution,cmplxA_field_s )
# writeh5(numpy.fft.fftshift(solution), "redetector")
errori.append(newerr)
cspace.get_realA(data_dum ,realA_field_support )
writeh5R( numpy.fft.fftshift (data_dum) , "support")
cspace.get_cmplxA( solution, cmplxA_obj_mem )
writeh5( numpy.fft.fftshift (solution), "field_focus_dopoerr")
cspace.shrinkwrap_support( realA_field_support , cmplxA_obj_fl, params.thl_sw )
for l in range(params.cycles_CF ):
cspace.set_phase(cmplxA_field_s, realA_field_s0, realA_phase_a )
from_detector_through_focus_to_support( params,derived_pars , target = cmplxA_obj_fl_new, source= cmplxA_field_s, dumA = cmplxA_field_dumA)
cspace.apply_chargeflip( cmplxA_obj_fl_new , cmplxA_obj_fl_new ,realA_field_support ) ## -----------------------------
from_support_through_focus_to_detector( params,derived_pars , target =cmplxA_field_s , source=cmplxA_obj_fl_new , dumA = cmplxA_field_dumA)
cspace.get_phase( realA_phase_a ,cmplxA_field_s )
cspace.apply_factor( cmplxA_field_s , cmplxA_field_s ,renorm )
newerr = cspace.error(cmplxA_field_s, realA_field_s0 )
errori.append(newerr)
print " ----- cahrge flip " , newerr
cspace.copy(cmplxA_obj_fl_last, cmplxA_obj_fl_new )
cspace.get_cmplxA( solution,cmplxA_field_s )
writeh5(numpy.fft.fftshift(solution), "field_detector")
# numpy.save("field_detector",numpy.fft.fftshift(solution))
cspace.set_phase(cmplxA_field_s, realA_field_s0, realA_phase_a )
cspace.go_backward( cmplxA_field_dumA , cmplxA_field_s )
cspace.apply_fresnel_phase( cmplxA_field_dumA , -derived_pars.z_fd , derived_pars.size_focal_pixel , cmplxA_field_dumA )
cspace.get_cmplxA( solution, cmplxA_field_dumA)
writeh5(numpy.fft.fftshift(solution), "field_focus")
#numpy.save("field_focus",numpy.fft.fftshift(solution))
cspace.apply_fresnel_phase( cmplxA_field_dumA , -derived_pars.focus*params.OSA_prox_factor, derived_pars.size_focal_pixel , cmplxA_field_dumA )
cspace.go_backward( cmplxA_obj_fl_new , cmplxA_field_dumA )
cspace.get_cmplxA( solution,cmplxA_obj_fl_new )
writeh5(numpy.fft.fftshift(solution), "field_support")
cspace.apply_fresnel_phase( cmplxA_field_dumA , -derived_pars.focus*(1-params.OSA_prox_factor),
derived_pars.size_focal_pixel , cmplxA_field_dumA )
cspace.go_backward( cmplxA_obj_fl_new , cmplxA_field_dumA )
cspace.get_cmplxA( solution,cmplxA_obj_fl_new )
writeh5(numpy.fft.fftshift(solution), "field_support")
print " FINE "
if(sys.argv[0][-12:]!="sphinx-build"):
main()