diff --git a/maoppy/psfmodel.py b/maoppy/psfmodel.py
index 85900c15cc2890aa493202829ec39c2b190e2308..ecfe25b57103c56af66dc6cf88bc8b6ae1ccf6ae 100644
--- a/maoppy/psfmodel.py
+++ b/maoppy/psfmodel.py
@@ -11,7 +11,6 @@ from scipy.optimize import least_squares
 from astropy.io import fits
 import time
 from numpy.fft import fft2, fftshift, ifft2
-from functools import lru_cache
 from maoppy.config import _EPSILON
 from maoppy.utils import binning
 
@@ -312,19 +311,15 @@ def moffat(XY, param, norm=None):
 class Moffat(ParametricPSF):
     
     def __init__(self,Npix,norm=np.inf):
-        #super(ParametricPSF,self).__init__(Npix)
         self.Npix = Npix
         self.norm = norm
         bounds_down = [_EPSILON,_EPSILON,-np.inf,1+_EPSILON]
         bounds_up = [np.inf for i in range(4)]
         self.bounds = (bounds_down,bounds_up)
-    
-    @lru_cache(maxsize=5)
-    def _XY(self,Npix):
-        YX = np.mgrid[0:Npix[0],0:Npix[1]]
-        YX[1] = YX[1] - Npix[0]/2
-        YX[0] = YX[0] - Npix[1]/2
-        return YX
+        
+        self._XY = np.mgrid[0:Npix[0],0:Npix[1]]
+        self._XY[1] -= Npix[0]/2
+        self._XY[0] -= Npix[1]/2
     
     def __call__(self,x,dx=0,dy=0):
         """
@@ -337,7 +332,7 @@ class Moffat(ParametricPSF):
             x[3] - Beta
         """
         y = np.concatenate(([1],x,[dx,dy]))
-        return moffat(self._XY(self.Npix),y,norm=self.norm)
+        return moffat(self._XY,y,norm=self.norm)
     
     def tofits(self,param,filename,*args,keys=["ALPHA_X","ALPHA_Y","THETA","BETA"],**kwargs):
         super(Moffat,self).tofits(param,filename,*args,keys=keys,**kwargs)
@@ -346,16 +341,19 @@ class Moffat(ParametricPSF):
 
 
 def oversample(samp):
-    """Oversample with an integer"""
-    if samp>2:
-        return (samp,1)
-    else:
-        k = int(np.ceil(2.0/samp))
-        return (k*samp,k)
+    """Find the minimal integer that allows oversampling"""
+    k = int(np.ceil(2.0/samp))
+    return (k*samp,k)
         
 
 
-class Psfao(ParametricPSF):    
+class Psfao(ParametricPSF):
+    """
+    Name:        PSFAO
+    Description: a long-exposure PSF model dedicated to AO correction
+    Reference:   Fétick et al., A&A, 2019 [Fétick2019b]
+    Note:        the PSF is parametrised through the PSD of the electromagnetic phase
+    """
     def __init__(self,Npix,system=None,Lext=10.,samp=None,symmetric=False):
         """
         Parameters
@@ -528,22 +526,15 @@ class Psfao(ParametricPSF):
         
     def tofits(self,param,filename,*args,keys=None,**kwargs):
         if keys is None:
-            if len(param)==5:
-                keys = ["R0","CST","SIGMA2","ALPHA","BETA"]
-                keys_comment = ["Fried parameter [m]",
-                                "PSD AO area constant C [rad2]",
-                                "PSD AO area Moffat variance A [rad2]",
-                                "PSD AO area Moffat alpha [1/m]",
-                                "PSD AO area Moffat beta"]
-            else: # if not 5, then equals 7
-                keys = ["R0","CST","SIGMA2","ALPHA_X","ALPHA_Y","THETA","BETA"]
-                keys_comment = ["Fried parameter [m]",
-                                "PSD AO area constant C [rad2]",
-                                "PSD AO area Moffat variance A [rad2]",
-                                "PSD AO area Moffat alpha X [1/m]",
-                                "PSD AO area Moffat alpha Y [1/m]",
-                                "PSD AO area Moffat theta [rad]",
-                                "PSD AO area Moffat beta"]
+            if len(param)==5: param = tuple(param[0:4])+(param[3],0,param[4])
+            keys = ["R0","CST","SIGMA2","ALPHA_X","ALPHA_Y","THETA","BETA"]
+            keys_comment = ["Fried parameter [m]",
+                            "PSD AO area constant C [rad2]",
+                            "PSD AO area Moffat variance A [rad2]",
+                            "PSD AO area Moffat alpha X [1/m]",
+                            "PSD AO area Moffat alpha Y [1/m]",
+                            "PSD AO area Moffat theta [rad]",
+                            "PSD AO area Moffat beta"]
         
         # redefine tofits() because extra hdr is required
         psf = self.__call__(param,*args,**kwargs)