Commit c39b1b5f authored by rfetick's avatar rfetick
Browse files

minor modif

parent c0143f51
......@@ -726,30 +726,27 @@ class Turbulent(ParametricPSFfromPSD):
"""
self.check_parameters(x0)
nx0,ny0 = self._nullFreqIndex
f2D = self._fxfy
F2 = f2D[0]**2. + f2D[1]**2.
r0,Lext = x0
PSD = 0.0229* r0**(-5./3.) * ((1./Lext**2.) + F2)**(-11./6.)
PSD = 0.0229* r0**(-5./3.) * ((1./Lext**2.) + self._f2)**(-11./6.)
# Set PSD = 0 at null frequency
PSD[nx0,ny0] = 0.0
# Compute PSD integral up to infinity
fmax = _np.min([nx0,ny0])*self._pix2freq
integral_in = _np.sum(PSD*(F2<(fmax**2))) * self._pix2freq**2 # numerical sum (in the array's tangent circle)
integral_in = _np.sum(PSD*(self._f2<(fmax**2))) * self._pix2freq**2 # numerical sum (in the array's tangent circle)
integral_out = 0.0229*6*_np.pi/5 * (r0*fmax)**(-5./3.) # analytical sum (outside the array's tangent circle)
integral = integral_in + integral_out
if grad:
g = _np.zeros((len(x0),)+F2.shape)
g = _np.zeros((len(x0),)+self._f2.shape)
g[0,...] = PSD*(-5./3)/r0
g[1,...] = PSD*(-11./6)/((1./Lext**2.) + F2)*(-2/Lext**3)
g[1,...] = PSD*(-11./6)/((1./Lext**2.) + self._f2)*(-2/Lext**3)
# compute integral gradient
integral_g = [0,0]
for i in range(2): integral_g[i] = _np.sum(g[i,...]*(F2<(fmax**2))) * self._pix2freq**2 # numerical sum (in the array's tangent circle)
for i in range(2): integral_g[i] = _np.sum(g[i,...]*(self._f2<(fmax**2))) * self._pix2freq**2 # numerical sum (in the array's tangent circle)
integral_g[0] += integral_out*(-5./3)/r0
if grad: return PSD,integral,g,integral_g
......@@ -798,9 +795,9 @@ class Psfao(ParametricPSFfromPSD):
fixed_k : int
Define a fixed oversampling factor
"""
self.Lext = Lext
super().__init__(7,npix,system=system,samp=samp,fixed_k=fixed_k)
self.Lext = Lext
# r0,C,A,alpha,ratio,theta,beta
### Mathematical bounds
#bounds_down = [_EPSILON,0,0,_EPSILON,_EPSILON,-_np.inf,1+_EPSILON]
......@@ -817,6 +814,7 @@ class Psfao(ParametricPSFfromPSD):
def _computeXYarray(self):
super()._computeXYarray()
self._maskin = (self._f2 < self.cutoffAOfreq**2.)
self._vk = (1-self._maskin) * 0.0229 * ((1. / self.Lext**2.) + self._f2)**(-11./6.)
def psd(self,x0,grad=False):
"""Compute the PSD model from parameters
......@@ -838,16 +836,14 @@ class Psfao(ParametricPSFfromPSD):
"""
self.check_parameters(x0)
nx0,ny0 = self._nullFreqIndex
#f2D = self._fxfy
pix = self._pix2freq
#F2 = self._fxfy[0]**2. + self._fxfy[1]**2.
#Fao = self.system.Nact/(2.0*self.system.D)
#maskin = (self._f2 < self.cutoffAOfreq**2.)
r0,C,A,alpha,ratio,theta,beta = x0
PSD = (1-self._maskin) * 0.0229 * r0**(-5./3.) * ((1. / self.Lext**2.) + self._f2)**(-11./6.)
# Von-Karman
PSD = (r0**(-5./3.)) * self._vk
# Moffat
ax = alpha*ratio
ay = alpha/ratio
param = (ax,ay,theta,beta,0,0)
......@@ -866,13 +862,13 @@ class Psfao(ParametricPSFfromPSD):
if grad: raise ValueError("PSFAO analytical gradient computation is not compatible with numerical Moffat normalization")
# Warning: Moffat numerical normalization generates strehlOTF jump when "_k" is changed
PSD += self._maskin * (C + A*moff) # before I wrote "moffat(self._fxfy,param,norm=Fao)" instead of "moff", so no numerical normalization
PSD += self._maskin * (C + A*moff)
PSD[nx0,ny0] = 0.0 # Set PSD = 0 at null frequency
# Compute PSD integral up to infinity
fmax = _np.min([nx0,ny0])*pix
integral_in = _np.sum(PSD*(self._f2<(fmax**2))) * pix**2 # numerical sum
integral_out = 0.0229*6*_np.pi/5 * (r0*fmax)**(-5./3.) # analytical sum
integral_in = _np.sum(PSD*(self._f2<(fmax**2))) * pix**2 # numerical sum on the array
integral_out = 0.0229*6*_np.pi/5 * (r0*fmax)**(-5./3.) # analytical sum outside
integral = integral_in + integral_out
if grad:
......
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