diff --git a/maoppy/psfmodel.py b/maoppy/psfmodel.py index e60ace53ce4a3868452c9653990428233fd8d27a..608d838a6e2ce387eb9f324c1413cb7e9286b0d1 100644 --- a/maoppy/psfmodel.py +++ b/maoppy/psfmodel.py @@ -42,7 +42,7 @@ class ParametricPSF: self.bounds = ((-np.inf,)*nparam, (np.inf,)*nparam) def __str__(self): - s = "ParametricPSF\n" + s = self.__class__.__name__ + "\n" s += "-------------\n" s += "pixels : (%u,%u)\n"%self.npix s += "nb param: %u"%self._nparam @@ -625,6 +625,25 @@ class Psfao(ParametricPSFfromPSD): return ["r0","bck","amp","alpha","ratio","theta","beta"] + def var_corr(self, parampsd): + """Return the numerical variance on the corrected area""" + psd,_ = self.psd(parampsd) + return np.sum(psd*self._maskin) * self.pix2freq**2 + + + def var_halo(self, parampsd, fovnorm=False): + """Return the numerical variance on the halo (fitting error)""" + psd,var_tot = self.psd(parampsd, fovnorm=fovnorm) + var_corr = np.sum(psd*self._maskin) * self.pix2freq**2 + return var_tot - var_corr + + + def var_total(self, parampsd, fovnorm=False): + """Return total variance halo + corrected area""" + _,var_tot = self.psd(parampsd, fovnorm=fovnorm) + return var_tot + + def psd(self, parampsd, grad=False, fovnorm=False): """Compute the PSD of the electromagnetic phase PSD is given in [rad²/f²] = [rad² m²] @@ -654,8 +673,9 @@ class Psfao(ParametricPSFfromPSD): r0, bck, amp, alpha, ratio, theta, beta = parampsd - if alpha < pix: - warnings.warn("`alpha < df`. Everything works fine but the value of `amp` might not correspond exactly to the PSD integral.") + if alpha < (2*pix): + # well, the factor 2 is really empirical... + warnings.warn("`alpha < 2*df`. Everything works fine but the value of `amp` might not correspond exactly to the PSD integral.") # Von-Karman psd = (r0**(-5./3.)) * self._vk