LGS tip-tillt correction

soapy/confParse.py

@@ -871,6 +871,10 @@ class LgsConfig(ConfigObj):
         ``naProfile``        list: The relative sodium layer
                              strength for each elongation
                              layer. If None, all equal.          ``None``
+        ``correctLgsTT``     bool: correct the LGS tip-tilt on   ``False``
+                                the sensor. Required 'removeTT'
+                                == True and 'uplinkgain>0'
+        ``uplinkgain``       float: LGS tip-tilt gain            0
         ==================== =================================   ===========
 
     """
@@ -888,6 +892,8 @@ class LgsConfig(ConfigObj):
                         ("elongationLayers", 10),
                         ("launchPosition",  numpy.array([0,0])),
                         ("naProfile", None),
+                        ("correctLgsTT", False),
+                        ("uplinkgain", 0)
                         ]
     calculatedParams = ["position"]
 
@@ -907,6 +913,8 @@ def calcParams(self):
 
         self.wavelength = float(self.wavelength)
         self.height = float(self.height)
+        self.launchPosition = numpy.array(self.launchPosition)
+        self.uplinkgain = float(self.uplinkgain)
 
 class DmConfig(ConfigObj):
     """

soapy/wfs/base.py

@@ -83,7 +83,7 @@
 
 import aotools
 
-from .. import AOFFT, LGS, logger, lineofsight_legacy
+from .. import AOFFT, LGS, logger, lineofsight_legacy, interp
 
 # xrange now just "range" in python3.
 # Following code means fastest implementation used in 2 and 3
@@ -161,6 +161,29 @@ def __init__(
         self.calcTiltCorrect()
         self.getStatic()
 
+        # Set up array for tip-tilt uplink compensation
+        if self.config.lgs:
+            if self.lgsConfig.correctLgsTT:
+                self._upTTCommand = numpy.array([0, 0], dtype='float')
+                self._tiptilt = numpy.array([0, 0], dtype='float')
+                self._nx_elements = int(round(self.sim_size))
+                fact = self.sim_size / self.pupil_size
+                coords = numpy.linspace(
+                            -1 * fact, 1 * fact, self._nx_elements) / 2
+                self.tip1arcsec, self.tilt1arcsec = numpy.meshgrid(coords,
+                                                                   coords)
+
+                tt_amp = -(ASEC2RAD * self.soapy_config.tel.telDiam/2.) * 1e9
+                self.tip1arcsec *= tt_amp
+                self.tilt1arcsec *= tt_amp
+
+                if self.config.removeTT == 0:
+                    logger.warning("LGS tiptilt correction will not work: "
+                                   "correctLgs==1 but removeTT==0")
+                if self.lgsConfig.uplinkgain == 0:
+                    logger.warning("LGS tiptilt correction will not work: "
+                                   "correctLgs==1 but uplinkgain==0")
+
 ############################################################
 # Initialisation routines
 
@@ -287,7 +310,7 @@ def calcElongPhaseAddition(self, elongLayer):
         h = self.elongHeights[elongLayer]
         dh = h - self.config.GSHeight
         H = float(self.lgsConfig.height)
-        d = numpy.array(self.lgsLaunchPos).astype('float32') * self.los.telDiam/2.
+        d = numpy.array(self.lgsLaunchPos).astype('float32') * self.los.telescope_diameter/2.
         D = self.telescope_diameter
         theta = (d.astype("float")/H) - self.config.GSPosition
 
@@ -410,6 +433,8 @@ def frame(self, scrns, phase_correction=None, read=True, iMatFrame=False):
             self.iMat = True
             removeTT = self.config.removeTT
             self.config.removeTT = False
+            lgsDic = self.config.lgs
+            self.config.lgs = None
             photonNoise = self.config.photonNoise
             self.config.photonNoise = False
             eReadNoise = self.config.eReadNoise
@@ -418,6 +443,9 @@ def frame(self, scrns, phase_correction=None, read=True, iMatFrame=False):
         self.zeroData(detector=read, FP=False)
 
         self.los.frame(scrns)
+        if self.config.lgs:
+            if self.lgsConfig.correctLgsTT:
+                self.correctUplinkTilt()
 
         # If LGS elongation simulated
         if self.config.lgs and self.elong!=0:
@@ -429,7 +457,7 @@ def frame(self, scrns, phase_correction=None, read=True, iMatFrame=False):
             self.uncorrectedPhase = self.los.phase.copy()/self.los.phs2Rad
             if phase_correction is not None:
                 self.los.performCorrection(phase_correction)
-                
+
             self.calcFocalPlane()
 
         if read:
@@ -441,6 +469,7 @@ def frame(self, scrns, phase_correction=None, read=True, iMatFrame=False):
         if iMatFrame:
             self.iMat=False
             self.config.removeTT = removeTT
+            self.config.lgs = lgsDic
             self.config.photonNoise = photonNoise
             self.config.eReadNoise = eReadNoise
 
@@ -481,6 +510,12 @@ def makeDetectorPlane(self):
     def LGSUplink(self):
         pass
 
+    def correctUplinkTilt(self):
+        self._upTTCommand += self.lgsConfig.uplinkgain * self._tiptilt
+
+        self.los.phase -= self._upTTCommand[0] * self.tip1arcsec
+        self.los.phase -= self._upTTCommand[1] * self.tilt1arcsec
+
     def calculateSlopes(self):
         self.slopes = self.los.EField
 

soapy/wfs/shackhartmann.py

@@ -74,7 +74,7 @@ def calcInitParams(self):
             self.SUBAP_OVERSIZE = 1
         else:
             self.SUBAP_OVERSIZE = 2
-        
+
         self.nx_detector_pixels = self.nx_subaps * (self.nx_subap_pixels + self.nx_guard_pixels) + self.nx_guard_pixels
 
         self.nx_subap_interp *= self.SUBAP_OVERSIZE
@@ -200,8 +200,7 @@ def initFFTs(self):
             self.iFFT = pyfftw.FFTW(
                 self.ifft_input_data, self.ifft_output_data, axes=(-2, -1),
                 threads=self.threads, flags=(self.config.fftwFlag, "FFTW_DESTROY_INPUT"),
-                direction="FFTW_BACKWARD"
-                )
+                direction="FFTW_BACKWARD")
 
             self.lgs_ifft_input_data = pyfftw.empty_aligned(
                 (self.subapFFTPadding, self.subapFFTPadding), dtype=CDTYPE)
@@ -210,8 +209,7 @@ def initFFTs(self):
             self.lgs_iFFT = pyfftw.FFTW(
                 self.lgs_ifft_input_data, self.lgs_ifft_output_data, axes=(0, 1),
                 threads=self.threads, flags=(self.config.fftwFlag, "FFTW_DESTROY_INPUT"),
-                direction="FFTW_BACKWARD"
-                )
+                direction="FFTW_BACKWARD")
 
     def initLGS(self):
         super(ShackHartmann, self).initLGS()
@@ -399,7 +397,7 @@ def applyLgsUplink(self):
 
         self.lgs.getLgsPsf(self.los.scrns)
 
-        self.lgs_ifft_input_data[:] = self.lgs.psf[::-1, ::-1]
+        self.lgs_ifft_input_data[:] = self.lgs.psf  #self.lgs.psf[::-1, ::-1]
         self.lgs_iFFT()
 
         self.ifft_input_data[:] = self.subap_focus_intensity
@@ -442,8 +440,13 @@ def calculateSlopes(self):
         self.slopes[:] = slopes.reshape(self.n_subaps * 2)
 
         if self.config.removeTT == True:
-            self.slopes[:self.n_subaps] -= self.slopes[:self.n_subaps].mean()
-            self.slopes[self.n_subaps:] -= self.slopes[self.n_subaps:].mean()
+            _tip = self.slopes[:self.n_subaps].mean()
+            _tilt = self.slopes[self.n_subaps:].mean()
+            self.slopes[:self.n_subaps] -= _tip
+            # self.slopes[:self.n_subaps].mean()
+            self.slopes[self.n_subaps:] -= _tilt
+            # self.slopes[self.n_subaps:].mean()
+            self._tiptilt = numpy.array([_tip, _tilt]) * self.pixel_scale
 
         if self.config.angleEquivNoise and not self.iMat:
             pxlEquivNoise = (
@@ -579,4 +582,4 @@ def photons_per_mag(mag, mask, phase_scale, exposureTime, zeropoint):
     # N photons for mag and exposure time
     n_photons *= (10**(-float(mag)/2.5)) * exposureTime
 
-    return n_photons
+    return n_photons