unbeam cov ok

python/spherelib/__init__.py

@@ -29,3 +29,4 @@ from fitsfunc2 import *
 from bin import *
 from alm import *
 from utils import *
+from cov import *

python/spherelib/cov.py

@@ -17,7 +17,42 @@
 # You should have received a copy of the GNU General Public License
 # along with this program; if not, see http://www.gnu.org/licenses/gpl.html
 
+from pylab import *
+import bin
+
 def unbeamcov (cov, beams, max_inv=10000):
+    """ Correct covariance matrices from beam effect
+
+    Parameters
+    ----------
+    cov: 3D array, symetric covariance matrices
+    beams: vector, beam fwhm in arcmin
+    max_inv: upper limit for inversion
     """
+    
+    nin  = shape(cov)[0]
+    lmax = shape(cov)[2]-1;
+
+    ib = zeros((nin, lmax+1))
+    for m in range(nin):
+        ib[m,:] = squeeze(1/bin.gaussbeam(beams[m], lmax))
+        ib[m,[ib[m,:]>max_inv]] = max_inv
+        
+    for ell in range(lmax+1):
+        cov[:,:,ell] = dot(dot(diag(ib[:,ell]),cov[:,:,ell]),diag(ib[:,ell]))
+
+
+def plotautocov (cov, figfile=None, leg=None):
+    """ 
     """
-    pass
+    nmap = shape(cov)[0]
+    figure()
+    for m in range(nmap):
+        semilogy(squeeze(cov[m,m,:]))
+
+    if leg is not None:
+        legend(leg)
+    if figfile is None:
+        show()
+    else:
+        savefig(figfile)