From 882b7820fc4b59da36ace51189901696558b4a7a Mon Sep 17 00:00:00 2001
From: Romain Fetick <rfetick@OSU.PYTHEAS>
Date: Tue, 9 Jul 2024 15:47:12 +0200
Subject: [PATCH] simplify file headers, update FITS functions

---
 maoppy/__init__.py   | 26 +++----------
 maoppy/instrument.py |  4 +-
 maoppy/psffit.py     |  3 +-
 maoppy/psfmodel.py   | 90 +++++++++++++++++++++++++++-----------------
 maoppy/psfutils.py   | 23 +----------
 maoppy/utils.py      |  4 +-
 setup.py             |  9 +----
 7 files changed, 67 insertions(+), 92 deletions(-)

diff --git a/maoppy/__init__.py b/maoppy/__init__.py
index dfc7ad7..eed7660 100644
--- a/maoppy/__init__.py
+++ b/maoppy/__init__.py
@@ -1,21 +1,7 @@
-#!/usr/bin/env python3
-# -*- coding: utf-8 -*-
-"""
-Created on Mon May 27 17:15:07 2019
+"""MAOPPY library"""
 
-@author: rfetick
-"""
-
-import sys as _sys
-import warnings as _warnings
-
-if _sys.version_info[0]<3:
-    _warnings.warn("MAOPPY was developped on Python 3, but your Python version is anterior. We hope everything will be fine")
-
-# MAOPPY release version
-__version__ = "1.6.0"
-__author__ = "Romain JL. Fétick (LAM, France)"
-__date__ = "January 18th, 2022"
-
-from . import utils, instrument, psfutils, psfmodel, psffit
-# from maoppy import utils, instrument, psfmodel
+from . import utils
+from . import instrument
+from . import psfutils
+from . import psfmodel
+from . import psffit
diff --git a/maoppy/instrument.py b/maoppy/instrument.py
index 6888d7e..0ff45a8 100644
--- a/maoppy/instrument.py
+++ b/maoppy/instrument.py
@@ -1,7 +1,5 @@
-#!/usr/bin/env python3
-# -*- coding: utf-8 -*-
 """
-Created on Mon May 27 17:30:51 2019
+List of all available instruments loaded from the <data> folder.
 
 @author: rfetick
 """
diff --git a/maoppy/psffit.py b/maoppy/psffit.py
index a723afd..d3c2b1b 100644
--- a/maoppy/psffit.py
+++ b/maoppy/psffit.py
@@ -1,6 +1,5 @@
-# -*- coding: utf-8 -*-
 """
-Created on Wed Nov 17 19:39:53 2021
+Methods to fit a PSF model on data.
 
 @author: rfetick
 """
diff --git a/maoppy/psfmodel.py b/maoppy/psfmodel.py
index 96ece54..ca4345d 100644
--- a/maoppy/psfmodel.py
+++ b/maoppy/psfmodel.py
@@ -1,7 +1,5 @@
-#!/usr/bin/env python3
-# -*- coding: utf-8 -*-
 """
-Created on Mon May 27 17:31:18 2019
+All PSF models available in this library.
 
 @author: rfetick
 """
@@ -11,7 +9,7 @@ import numpy as np
 from scipy.fft import fft2, ifft2, fftshift, ifftshift
 from astropy.io import fits
 from maoppy.utils import binning, random_sample
-from maoppy.psfutils import moffat, gauss, oversample, moffat_center, make_fits_hdr
+from maoppy.psfutils import moffat, gauss, oversample, moffat_center
 import warnings
 
 __all__ = ["ParametricPSF", "ConstantPSF", "Moffat", "Gaussian",
@@ -53,6 +51,11 @@ class ParametricPSF:
         """Ordered list of the parameters names"""
         return ["param_%u"%i for i in range(self._nparam)]
     
+    @property
+    def param_comment(self):
+        """Comment on the parameters (self.param_name)"""
+        return ["default name" for i in range(self._nparam)]
+    
     def dict2list(self, dctnry):
         """Transform a dictionary of parameters into a list"""
         return [dctnry[k] for k in self.param_name]
@@ -73,12 +76,32 @@ class ParametricPSF:
         """Return the Modulation Transfer Function (MTF)"""
         return np.abs(self.otf(args,kwargs))
     
-    def tofits(self, param, filename, *args, **kwargs):
-        """Save PSF as FITS data"""
+    def fits_header_read(self, hdr):
+        """Read a FITS header and return PSF parameters.
+        Also accepts a path to the FITS file."""
+        if type(hdr)==str:
+            fts = fits.open(hdr)
+            hdr = fts[0].header
+            fts.close()
+        param = [0]*len(self.param_name)
+        for i in range(len(self.param_name)):
+            param[i] = float(hdr[self.param_name[i]])
+        return param
+    
+    def fits_header_make(self, param):
+        """Generate a FITS header and fill it with the given PSF parameters"""
+        hdr = fits.Header()
+        hdr['ORIGIN'] = 'MAOPPY automatic header'
+        for i in range(len(param)):
+            hdr[self.param_name[i]] = (param[i], self.param_comment[i])
+        return hdr
+    
+    def fits_psf_write(self, param, filename, *args, **kwargs):
+        """Write PSF into a FITS file"""
         psf = self.__call__(param, *args, **kwargs)
-        hdr = make_fits_hdr(param, keys=self.param_name)
+        hdr = self.fits_header_make(param)
         hdu = fits.PrimaryHDU(psf, hdr)
-        hdu.writeto(filename) # overwrite ?
+        hdu.writeto(filename)
 
 
 class ConstantPSF(ParametricPSF):
@@ -650,8 +673,8 @@ class Psfao(ParametricPSFfromPSD):
         #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]
+        bounds_down = [1e-3, 0, 0, 1e-5, 1e-2, -np.inf, 1.01]
+        bounds_up = [np.inf]*4 + [1e2, np.inf, 5]
         self.bounds = (bounds_down,bounds_up)
     
     
@@ -678,6 +701,17 @@ class Psfao(ParametricPSFfromPSD):
         """Ordered list of the parameters names"""
         return ["r0","bck","amp","alpha","ratio","theta","beta"]
     
+    @property
+    def param_comment(self):
+        """Comment on the parameters (self.param_name)"""
+        comments = ["Fried parameter [m]",
+                    "AO background [rad2 m2]",
+                    "AO Moffat variance [rad2]",
+                    "AO Moffat alpha [1/m]",
+                    "AO Moffat sqrt(ax/ay) ratio",
+                    "AO Moffat theta [rad]",
+                    "AO Moffat beta"]
+        return comments
     
     def var_corr(self, parampsd):
         """Return the numerical variance on the corrected area"""
@@ -802,31 +836,17 @@ class Psfao(ParametricPSFfromPSD):
         return psd, integral
         
     
-    def tofits(self, param, filename, *args, **kwargs):
-        keys_comment = ["Fried parameter [m]",
-                        "AO background [rad2 m2]",
-                        "AO Moffat variance A [rad2]",
-                        "AO Moffat alpha [1/m]",
-                        "AO Moffat sqrt(ax/ay) ratio",
-                        "AO Moffat theta [rad]",
-                        "AO Moffat beta"]
-        
-        # redefine tofits() because extra hdr is required
-        psf = self.__call__(param, *args, **kwargs)
-        hdr = make_fits_hdr(param, keys=self.param_name, keys_comment=keys_comment)
-        
-        hdr['CDELT1'] = (self.system.resolution_mas,"pixel size")
-        hdr['CUNIT1'] = ("mas","pixel size unit unit")
-        hdr['CDELT2'] = (self.system.resolution_mas,"pixel size")
-        hdr['CUNIT2'] = ("mas","pixel size unit unit")
-        hdr["HIERARCH SYSTEM"] = (self.system.name,"System name")
-        hdr["HIERARCH SYSTEM D"] = (self.system.D,"Primary mirror diameter")
-        hdr["HIERARCH SYSTEM NACT"] = (self.system.Nact,"Linear number of AO actuators")
-        hdr["HIERARCH SAMP"] = (self.samp,"Sampling (eg. 2 for Shannon)")
-        hdr["HIERARCH lext"] = (self.lext,"Von-Karman outer scale")
-        
-        hdu = fits.PrimaryHDU(psf, hdr)
-        hdu.writeto(filename)
+    def fits_header_make(self, param):
+        """Generate a FITS header and fill it with the given PSF parameters"""
+        hdr = super().fits_header_make(param)
+        # hdr['CDELT1'] = (self.system.resolution_mas,"pixel size")
+        # hdr['CUNIT1'] = ("mas","pixel size unit unit")
+        # hdr['CDELT2'] = (self.system.resolution_mas,"pixel size")
+        # hdr['CUNIT2'] = ("mas","pixel size unit unit")
+        hdr["SYSTEM"] = (self.system.name,"System name")
+        hdr["SAMP"] = (self.samp,"Sampling (eg. 2 for Shannon)")
+        hdr["LEXT"] = (self.lext,"Von-Karman outer scale")
+        return hdr
 
 
 
diff --git a/maoppy/psfutils.py b/maoppy/psfutils.py
index 20015d0..c5216f4 100644
--- a/maoppy/psfutils.py
+++ b/maoppy/psfutils.py
@@ -1,15 +1,13 @@
-# -*- coding: utf-8 -*-
 """
-Created on Tue Jan 18 12:02:10 2022
+List of useful functions to generate PSF.
 
 @author: rfetick
 """
 
 import numpy as np
-from astropy.io import fits
 
 __all__ = ["oversample", "reduced_coord", "reduced_center_coord", "moffat",
-           "moffat_center", "gauss", "make_fits_hdr"]
+           "moffat_center", "gauss"]
 
 
 def oversample(samp, fixed_k = None):
@@ -216,20 +214,3 @@ def gauss(xxyy,param):
     uu = reduced_coord(xxyy,ax,ay,param[2],param[3],param[4])
     return np.exp(-uu) / (2*np.pi*param[0]*param[1])
 
-
-def make_fits_hdr(param, keys=None, keys_comment=None):
-    """Define a header to save a .fits file"""
-    if len(keys)!=len(param):
-        raise ValueError("`keys` must be same size as `param`")
-    if keys_comment is not None:
-        if len(keys_comment)!=len(param):
-            raise ValueError("When defined, `keys_comment` must be same size as `param`")
-    hdr = fits.Header()
-    
-    hdr["HIERARCH ORIGIN"] = "MAOPPY automatic header"
-    for i in range(len(param)):
-        if keys_comment is None:
-            hdr["HIERARCH PARAM "+keys[i]] = param[i]
-        else:
-            hdr["HIERARCH PARAM "+keys[i]] = (param[i],keys_comment[i])
-    return hdr
diff --git a/maoppy/utils.py b/maoppy/utils.py
index e9c22aa..d096b11 100644
--- a/maoppy/utils.py
+++ b/maoppy/utils.py
@@ -1,7 +1,5 @@
-#!/usr/bin/env python3
-# -*- coding: utf-8 -*-
 """
-Created on Mon May 27 17:27:44 2019
+List of useful functions.
 
 @author: rfetick
 """
diff --git a/setup.py b/setup.py
index a35bfbd..d1b052f 100644
--- a/setup.py
+++ b/setup.py
@@ -1,13 +1,6 @@
-#!/usr/bin/env python3
-# -*- coding: utf-8 -*-
-"""
-Created on Mon May 27 17:11:01 2019
-
-@author: rfetick
-"""
 
 from setuptools import setup, find_packages
- 
+
 setup(name='maoppy',
       version='1.6.0',
       url='https://gitlab.lam.fr/lam-grd-public/maoppy.git',
-- 
GitLab