Add functions in instrument.py

paompy/instrument.py

@@ -6,3 +6,130 @@ Created on Mon May 27 17:30:51 2019
 @author: rfetick
 """
 
+from paompy.utils import circarr, airy, RAD2ARCSEC
+from paompy.config import _DEFAULT_RES
+
+#%% DETECTOR CLASS AND ITS SUBCLASSES
+
+class Detector(object):
+    """Represents a detector
+    
+    Attributes
+    ----------
+    resolution : float
+        Pixel scale [meter]
+    Npix : tuple, list, numpy.ndarray
+        Number of pixels on X and Y axis
+    gainADU : float
+        Detector gain [electron/ADU]
+    RON : float
+        Read-Out-Noise [electron]
+        
+    """
+    
+    def __init__(self,Npix,resolution=_DEFAULT_RES,gainADU=1.,RON=0.):
+        self.resolution_pixel = resolution
+        self.Npix = Npix
+        self.gainADU = gainADU
+        self.RON = RON
+        self.binning = 1
+    
+    def __repr__(self):
+        s  = "PAOMPY Detector\n"
+        s += "---------------\n"
+        s += "Pixels    : (%u,%u)\n" % (self.Npix[0],self.Npix[1])
+        s += "Resolution: %u um\n" % round(self.resolution*1e6)
+        s += "Binning   : %u\n" % self.binning
+        s += "Gain ADU  : %.2f e-/ADU\n" % self.gainADU
+        s += "RON       : %.2f e-" % self.RON
+        return s
+    
+    @property
+    def resolution(self):
+        return self.resolution_pixel * self.binning
+
+
+ZIMPOL_DETECTOR = Detector((1024,1024),30*1e-6,gainADU=10.5,RON=20.)
+
+# Equivalent detector resolution after image reduction pipeline
+MUSE_DETECTOR = Detector((200,200),237.14745*1e-6,gainADU=5.,RON=15.)
+
+
+
+
+#%% INSTRUMENT CLASS AND ITS SUBCLASSES
+
+class Instrument(object):
+    """Represents an optical system
+    
+    Attributes
+    ----------
+    D : float
+        Entrance pupil diameter [meter]
+    detector : Detector
+        Camera at the focal plane
+    filters : dict
+        Dictionary of available filters as tuples (central wvl, width) [meter]
+    AO_Nact : int
+        Linear number of actuators
+    """
+
+    def __init__(self,D,detector=None,occ=0.):
+        self.D = D
+        self.occ = occ # occultation ratio
+        self.detector = detector
+        self.filters = {}
+        self.AO_Nact = 0
+        self.focal_length = None
+        self.name = ""
+        
+    def __repr__(self):
+        s  = self.name+" OpticalSystem\n"
+        s += "-------------------------\n" 
+        s += "Diameter: %.2g m (occ=%u%%)\n" % (self.D,self.occ*100)
+        s += "AO_Nact : %u\n" % self.AO_Nact
+        s += "Focal   : %s m\n" % str(self.focal_length)
+        s += "Filters : %u" % len(self.filters)
+        return s
+    
+    @property
+    def resolution_mas(self):
+        if self.focal_length is None:
+            raise ValueError("Cannot compute `resolution_mas` if `focal_length` is not set")
+        return self.detector.resolution/self.focal_length * RAD2ARCSEC * 1e3
+    
+    def pupil(self,Npix,wvl=None,samp=None):
+        """Returns the 2D array of the pupil transmission function"""
+        Dpix = min(Npix)/2
+        pup = circarr(Npix)
+        return (pup < Dpix) * (pup >= Dpix*self.occ)
+    
+    def samp(self,wvl):
+        """Returns sampling value for the given wavelength"""
+        if self.detector is None:
+            raise ValueError("Cannot compute sampling if `detector` is not defined")
+        if self.focal_length is None:
+            raise ValueError("Cannot compute sampling if `focal_length` is not defined")
+        return wvl*self.focal_length/(self.detector.resolution*self.D)
+    
+    def wvl(self,samp):
+        """Returns wavelength for the given sampling"""
+        if self.detector is None:
+            raise ValueError("Cannot compute wavelength if `detector` is not defined")
+        if self.focal_length is None:
+            raise ValueError("Cannot compute wavelength if `focal_length` is not defined")
+        return samp*(self.detector.resolution*self.D)/self.focal_length
+    
+    
+    def PSFdl(self,Npix,wvl):
+        """Returns the diffraction limited PSF
+        
+        Parameters
+        ----------
+        Npix : tuple, list of 2 elements
+            Size of the output 2D array
+        wvl : float
+            Observation wavelength
+        """
+        return airy(Npix,self.samp(wvl),self.occ)
+