COMPASS  5.0.0
End-to-end AO simulation tool using GPU acceleration
psfMap.py
1 
37 
38 import numpy as np
39 import astropy.io.fits as fits
40 import matplotlib.pyplot as plt
41 import os
42 import shutil
43 
44 
45 class PSF_map:
46 
47  def __init__(self, SR=np.zeros((0, 0)), radius=0, wl=0, NGS=None, LGS=None, TS=None,
48  AOtype="AO", sup=None, filename=None):
49  """ SR_map class
50 
51  SR : np.ndarray : array of strehl (organised by position)
52  radius : float : radius of the map (in arcsec)
53  wl : float : wavelength of observation (in nm)
54  NGS : np.ndarray : position (arcsec) of the NGS (1st line:x axis, 2nd:y axis)
55  LGS : np.ndarray : position (arcsec) of the LGS (1st line:x axis, 2nd:y axis)
56  TS : np.ndarray : position (arcsec) of the TS (1st line:x axis, 2nd:y axis)
57  sup : shesha.supervisor : shesha supervisor to retrieve all data from
58  filename: str : filename to read the psf map from
59 
60  warning:
61  using 'sup' will take only the supervisor into consideration
62  using 'filename' will overide previous values of the psf map (unless called with sup)
63 
64  """
65  self.map = SR
66 
67  self._Rtar = radius
68  self._Rngs = 0
69  self._Rlgs = 0
70  self._Rlgs = 0
71  if (NGS is not None):
72  self.NGSx = NGS[0]
73  self.NGSy = NGS[1]
74  self._Rngs = max((self.NGSx.max(), self.NGSy.max()))
75  else:
76  self.NGSx = np.zeros((0))
77  self.NGSy = np.zeros((0))
78 
79  if (LGS is not None):
80  self.LGSx = LGS[0]
81  self.LGSy = LGS[1]
82  self._Rlgs = max(self.LGSx.max(), self.LGSy.max())
83  else:
84  self.LGSx = np.zeros((0))
85  self.LGSy = np.zeros((0))
86 
87  if (TS is not None):
88  self.TSx = TS[0]
89  self.TSy = TS[1]
90  self._Rts = max(self.TSx.max(), self.TSy.max())
91  else:
92  self.TSx = np.zeros((0))
93  self.TSy = np.zeros((0))
94 
95  self.wavelength = wl
96  self.type = AOtype
97 
98  if (sup is not None):
99  self.NGSx = np.array([t.xpos for t in sup.config.p_wfs_ngs[:-1]])
100  self.NGSy = np.array([t.ypos for t in sup.config.p_wfs_ngs[:-1]])
101  self._Rngs = max((self.NGSx.max(), self.NGSy.max()))
102  self.LGSx = np.array([t.xpos for t in sup.config.p_wfs_lgs])
103  self.LGSy = np.array([t.ypos for t in sup.config.p_wfs_lgs])
104  self._Rlgs = max(self.LGSx.max(), self.LGSy.max())
105  self.TSx = np.array([t.xpos for t in sup.config.p_wfs_ts])
106  self.TSy = np.array([t.ypos for t in sup.config.p_wfs_ts])
107  self._Rts = max(self.TSx.max(), self.TSy.max())
108  self.wavelength = sup.config.p_targets[0].Lambda
109  NTAR = len(sup.config.p_targets)
110  NTAR_side = int(np.sqrt(NTAR))
111  if (NTAR != NTAR_side**2):
112  raise ValueError("not a square nb of targets")
113  self.map = np.zeros((NTAR_side, NTAR_side))
114  for i in range(NTAR):
115  #self.map.itemset(i,sup._sim.get_strehl(i)[1])
116  self.map.itemset(i, sup._sim.get_strehl(i)[0])
117  tar = sup._sim.tar.d_targets[i]
118  self._Rtar = max(self._Rtar, tar.posx, tar.posy)
119 
120  elif (filename is not None):
121  self.read(filename)
122 
123  def setNGS(self, NGS, NGSy=None):
124  if (NGSy is None):
125  self.NGSx = NGS[0]
126  self.NGSy = NGS[1]
127  else:
128  self.NGS = NGS
129  self.NGSy = NGS
130  self._Rngs = max((self.NGSx.max(), self.NGSy.max()))
131 
132  def setLGS(self, LGS, LGSy=None):
133  if (LGSy is None):
134  self.LGSx = LGS[0]
135  self.LGSy = LGS[1]
136  else:
137  self.LGS = LGS
138  self.LGSy = LGS
139  self._Rlgs = max(self.LGSx.max(), self.LGSy.max())
140 
141  def setTS(self, TS, TSy=None):
142  if (TSy is None):
143  self.TSx = LGS[0]
144  self.TSy = LGS[1]
145  else:
146  self.TS = LGS
147  self.TSy = LGS
148  self._Rts = max(self.TSx.max(), self.TSy.max())
149 
150  def setWaveLength(self, wl):
151  self.wavelength = wl
152 
153  def plot(self, title=False, GS=False, WFS=False, LGS=False, NGS=False, TS=False):
154  if (self.map.shape[0] > 0 and self.map.shape[1] > 0):
155  plt.matshow(self.map,
156  extent=[-self._Rtar, self._Rtar, -self._Rtar, self._Rtar])
157  plt.colorbar()
158  if (GS or WFS or LGS):
159  plt.scatter(self.LGSy, self.LGSx, color="red")
160  if (GS or WFS or NGS):
161  plt.scatter(self.NGSy, self.NGSx, color="blue")
162  if (GS or TS):
163  plt.scatter(self.TSy, self.TSx, color="yellow", s=1)
164 
165  t = self.type + " Strehl"
166  if (self.wavelength > 0):
167  t += " @ {:.3f} nm".format(self.wavelength)
168 
169  if (self._Rts > 0):
170  t += ", {:.2f}".format(self._Rts) + " arcsec ring optim"
171  if (title):
172  plt.title(t)
173 
174  def save(self, name=""):
175  hdu_map = fits.PrimaryHDU(self.map)
176  hdu_map.header["TYPE"] = self.type
177  hdu_map.header["LAMBDA"] = self.wavelength
178  hdu_map.header["RTAR"] = self._Rtar
179  hdu_map.header["RLGS"] = self._Rlgs
180  hdu_map.header["RNGS"] = self._Rngs
181  hdu_map.header["RTS"] = self._Rts
182  hdu_LGSX = fits.ImageHDU(self.LGSx, name="LGSX")
183  hdu_LGSY = fits.ImageHDU(self.LGSy, name="LGSY")
184  hdu_NGSX = fits.ImageHDU(self.NGSx, name="NGSX")
185  hdu_NGSY = fits.ImageHDU(self.NGSy, name="NGSY")
186  hdu_TSX = fits.ImageHDU(self.TSx, name="TSX")
187  hdu_TSY = fits.ImageHDU(self.TSy, name="TSY")
188 
189  hdul = fits.HDUList([
190  hdu_map, hdu_LGSX, hdu_LGSY, hdu_NGSX, hdu_NGSY, hdu_TSX, hdu_TSY
191  ])
192  t = self.type + "_StrehlMap"
193  if (self.wavelength > 0):
194  t += "_{:.3f}nm".format(self.wavelength)
195  if (self._Rts > 0):
196  t += "_{:.2f}arcsec".format(self._Rts)
197  t += ".fits"
198  if (name == ""):
199  name = t
200  hdul.writeto(name, overwrite=1)
201 
202  def read(self, name):
203  hdu_map = fits.open(name)
204  self.type = hdu_map[0].header["TYPE"]
205  self.wavelength = hdu_map[0].header["LAMBDA"]
206  self._Rtar = hdu_map[0].header["RTAR"]
207  self._Rlgs = hdu_map[0].header["RLGS"]
208  self._Rngs = hdu_map[0].header["RNGS"]
209  self._Rts = hdu_map[0].header["RTS"]
210  self.map = hdu_map[0].data
211  self.LGSx = hdu_map["LGSx"].data
212  self.LGSy = hdu_map["LGSy"].data
213  self.NGSx = hdu_map["NGSx"].data
214  self.NGSy = hdu_map["NGSy"].data
215  self.TSx = hdu_map["TSx"].data
216  self.TSy = hdu_map["TSy"].data
shesha.util.psfMap.PSF_map.read
def read(self, name)
Definition: psfMap.py:202
shesha.util.psfMap.PSF_map.TSx
TSx
Definition: psfMap.py:87
shesha.util.psfMap.PSF_map.save
def save(self, name="")
Definition: psfMap.py:174
shesha.util.psfMap.PSF_map.setLGS
def setLGS(self, LGS, LGSy=None)
Definition: psfMap.py:132
shesha.util.psfMap.PSF_map.LGSy
LGSy
Definition: psfMap.py:80
shesha.util.psfMap.PSF_map.map
map
Definition: psfMap.py:64
shesha.util.psfMap.PSF_map._Rtar
_Rtar
Definition: psfMap.py:66
shesha.util.psfMap.PSF_map._Rngs
_Rngs
Definition: psfMap.py:67
shesha.util.psfMap.PSF_map
Definition: psfMap.py:45
shesha.util.psfMap.PSF_map.NGSx
NGSx
Definition: psfMap.py:71
shesha.util.psfMap.PSF_map._Rts
_Rts
Definition: psfMap.py:89
shesha.util.psfMap.PSF_map.NGSy
NGSy
Definition: psfMap.py:72
shesha.util.psfMap.PSF_map.__init__
def __init__(self, SR=np.zeros((0, 0)), radius=0, wl=0, NGS=None, LGS=None, TS=None, AOtype="AO", sup=None, filename=None)
SR_map class.
Definition: psfMap.py:63
shesha.util.psfMap.PSF_map.LGS
LGS
Definition: psfMap.py:137
shesha.util.psfMap.PSF_map.wavelength
wavelength
Definition: psfMap.py:94
shesha.util.psfMap.PSF_map.setWaveLength
def setWaveLength(self, wl)
Definition: psfMap.py:150
shesha.util.psfMap.PSF_map._Rlgs
_Rlgs
Definition: psfMap.py:68
shesha.util.psfMap.PSF_map.LGSx
LGSx
Definition: psfMap.py:79
shesha.util.psfMap.PSF_map.TS
TS
Definition: psfMap.py:146
shesha.util.psfMap.PSF_map.NGS
NGS
Definition: psfMap.py:128
shesha.util.psfMap.PSF_map.TSy
TSy
Definition: psfMap.py:88
shesha.util.psfMap.PSF_map.setNGS
def setNGS(self, NGS, NGSy=None)
Definition: psfMap.py:123
shesha.util.psfMap.PSF_map.type
type
Definition: psfMap.py:95
shesha.util.psfMap.PSF_map.plot
def plot(self, title=False, GS=False, WFS=False, LGS=False, NGS=False, TS=False)
Definition: psfMap.py:153
shesha.util.psfMap.PSF_map.setTS
def setTS(self, TS, TSy=None)
Definition: psfMap.py:141