|
On this page... (hide) Resolution degrading routine for polarization datasets, using 2D weighed average. 1. 2D Weighted mean estimator :map {$ m_i = [Q_i, U_i] $} mean estimator for big pixel {$ p $} : {$ \chi^2 = \sum_{i\in p} ( \mu - m_i) C_i^{-1} ( \mu - m_i)^T $} chi2 derivative : {$ \partial_\mu \chi^2 = 2 \sum_{i\in p} ( \mu - m_i)^T C_i^{-1} = 0 $} optimal estimator : {$ \mu = (\sum_{i\in p} C_i^{-1})^{-1} \cdot (\sum_{i\in p} C_i ^{-1} m_i) $} 2. implementation :
3. python code :Input map P = [Q, U] and C = [covQQ, covQU, covUU] unseenmap = (P[0] == hp.UNSEEN) + (P[1] == hp.UNSEEN) nside_in = hp.npix2nside(len(C[0])) dets = 1.0*(C[0]*C[2] - C[1]**2) dets[dets == 0] = np.inf dets[unseenmap] = np.inf # Compute Wud = ud(C^-1) W = np.array([C[2]/dets, -C[1]/dets, C[0]/dets]) W[:, unseenmap] = hp.UNSEEN Wud = np.array([hp.ud_grade(W[0], nside_out, order_in=order, order_out=order), hp.ud_grade(W[1], nside_out, order_in=order, order_out=order), hp.ud_grade(W[2], nside_out, order_in=order, order_out=order)]) Wud[Wud == hp.UNSEEN] = 0 # Compute ud_grade([Q,U].C^-1) pWP = np.array([(C[2]*P[0] - C[1]*P[1])/dets, (C[1]*P[0] + C[0]*P[1])/dets]) pWP[:, unseenmap] = hp.UNSEEN WP = np.array([hp.ud_grade(pWP[0], nside_out, order_in=order, order_out=order), hp.ud_grade(pWP[1], nside_out, order_in=order, order_out=order)]) #### WP[WP == hp.UNSEEN] = 0 unseenmapWP = (WP[0] == hp.UNSEEN) + (WP[1] == hp.UNSEEN) dets = 1.0*(Wud[0]*Wud[2] - Wud[1]**2) dets[dets == 0] = np.inf Q = (Wud[2]*WP[0] - Wud[1]*WP[1])/dets U = (-Wud[1]*WP[0] + Wud[0]*WP[1])/dets Pud = np.array([Q, U]) Pud[:, unseenmapWP] = hp.UNSEEN Wud = (nside_in/(1.0*nside_out))**2*Wud dets = 1.0*(Wud[0]*Wud[2] - Wud[1]**2) dets[dets == 0] = np.inf Cud = np.array([Wud[2], -Wud[1], Wud[0]])/dets 4. directory :/sps/planck/Users/vanneste/Work/Planck/maps/PR3/wudgraded/ |
|