Commit 838e9678 authored by rfetick's avatar rfetick
Browse files

minor modif

parent e0fbfa4d
# -*- coding: utf-8 -*-
"""
Created on Tue Jun 1 10:56:24 2021
@author: rfetick
"""
import matplotlib.pyplot as plt
import numpy as np
import copy
from maoppy.psfmodel import Turbulent
from maoppy.instrument import muse_nfm
Npix = 128 # pixel size of PSF
wvl = 600*1e-9 # wavelength [m]
#%% Initialize PSF model
samp = muse_nfm.samp(wvl) # sampling (2.0 for Shannon-Nyquist)
Pmodel = Turbulent((Npix,Npix),system=muse_nfm,samp=samp)
#%% Choose parameters
r0 = 0.15 # Fried parameter [m]
L0 = 20 # external length [m]
param = [r0,L0]
#%% Compute PSD derivative
psd,_,g,_ = Pmodel.psd(param,grad=True)
gnum = 0*g
eps = 1e-5
for i in range(len(param)):
dp = copy.copy(param)
dp[i] += eps
psd2,_ = Pmodel.psd(dp)
gnum[i,...] = (psd2-psd)/eps
#%% Plot results
names = ['r0','L0']
plt.figure(1)
plt.clf()
for i in range(len(param)):
plt.subplot(2,len(param),i+1)
plt.pcolormesh(np.arcsinh(g[i,...]))
plt.axis('image')
plt.axis('off')
plt.title(names[i])
plt.colorbar(orientation='horizontal')
plt.subplot(2,len(param),i+1+len(param))
plt.pcolormesh(np.arcsinh(gnum[i,...]))
plt.axis('image')
plt.axis('off')
#plt.title(names[i])
plt.colorbar(orientation='horizontal')
#%% Compute PSF derivative
psf,g = Pmodel(param,grad=True)
gnum = 0*g
eps = 1e-5
for i in range(len(param)):
dp = copy.copy(param)
dp[i] += eps
psf2 = Pmodel(dp)
gnum[i,...] = (psf2-psf)/eps
#%% Plot results
names = ['r0','L0']
plt.figure(2)
plt.clf()
for i in range(len(param)):
plt.subplot(2,len(param),i+1)
plt.pcolormesh(g[i,...])
plt.axis('image')
plt.axis('off')
plt.title(names[i])
plt.colorbar(orientation='horizontal')
plt.subplot(2,len(param),i+1+len(param))
plt.pcolormesh(gnum[i,...])
plt.axis('image')
plt.axis('off')
#plt.title(names[i])
plt.colorbar(orientation='horizontal')
\ No newline at end of file
......@@ -699,9 +699,13 @@ class Psfao(ParametricPSFfromPSD):
super().__init__(7,npix,system=system,samp=samp)
self.Lext = Lext
# r0,C,A,alpha,ratio,theta,beta
bounds_down = [_EPSILON,0,0,_EPSILON,_EPSILON,-_np.inf,1+_EPSILON]
bounds_up = [_np.inf for i in range(7)]
# r0,C,A,alpha,ratio,theta,beta
### Mathematical bounds
#bounds_down = [_EPSILON,0,0,_EPSILON,_EPSILON,-_np.inf,1+_EPSILON]
#bounds_up = [_np.inf for i in range(7)]
### Physical bounds
bounds_down = [1e-3,0,0,1e-3,1e-2,-_np.inf,1.01]
bounds_up = [_np.inf]*4 + [1e2,_np.inf,5]
self.bounds = (bounds_down,bounds_up)
def psd(self,x0,grad=False):
......@@ -753,6 +757,7 @@ class Psfao(ParametricPSFfromPSD):
integral_out = 0.0229*6*_np.pi/5 * (r0*fmax)**(-5./3.) # analytical sum
integral = integral_in + integral_out
# TODO: finish gradient computation!
if grad:
g = _np.zeros((len(x0),)+F2.shape)
# derivative towards r0
......
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