Apodisation window functions

[lukashergt]

Hi David,

two questions regarding the apodisation options. The documentation for mask_apodization lists the following two functions:

• 'C1' apodisation: $f = x - \sin(2 \pi x) / (2 \pi)$
• 'C2' apodisation: $f = (1/2) \cdot [1 - \cos(\pi x)]$

for $x < 1$ and where $x = \sqrt{(1 - \cos(\theta)) / (1 - \cos(\theta_\ast))}$.

This does not match the cited source (Grain et al. 2009) in two ways:

1. The naming of C1 and C2 in pymaster appears to be the opposite of the naming in the paper, i.e.:

   • pymaster 'C1' apodisation resembles more equation (31) in the paper which is referred to as the $C^2$-window therein.
   • pymaster 'C2' apodisation resembles more equation (30) in the paper which is referred to as the $C^1$-window therein.

   Is this just a typo? Or am I maybe missing something subtle?

2. The equations above match equations (30) and (31) in the paper almost exactly, except for the definition of x. The paper effectively uses $x=\theta/\theta_\ast$.
   The difference can be seen in the plot below. The paper definition results in symmetric transitions with a mask value of 0.5 at the half way point. The pymaster doc definition results in an asymmetric transition.
   Why is that? Was this shown to perform better?

import numpy as np
from scipy import stats
import matplotlib.pyplot as plt
fw, fh = plt.rcParams['figure.figsize']

sigma = 1
theta_star = 2.5 * sigma
theta = np.linspace(0, theta_star, 101)
norm = stats.norm(loc=theta_star, scale=sigma)

fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(2*fw, fh))

x = np.sqrt((1-np.cos(theta)) / (1-np.cos(theta_star)))
ax1.set_title(r"pymaster definition of $x = \sqrt{(1-\cos\theta)~/~(1-\cos\theta_\ast))}$")
ax1.plot(theta, norm.pdf(theta)/norm.pdf(theta).max(), c='k', label="pymaster: apotype='Smooth'")
ax1.plot(theta, x - np.sin(2 * np.pi * x) / (2 * np.pi),      label="pymaster: apotype='C1'")
ax1.plot(theta, 1 / 2 * (1 - np.cos(np.pi * x)),              label="pymaster: apotype='C2'")
ax1.axvline(theta_star/2, c='k', ls=':', lw=0.6)

x = theta / theta_star
ax2.set_title(r"Grain et al. (2009) definition of $x = \theta~/~\theta_\ast$")
ax2.plot(theta, norm.pdf(theta)/norm.pdf(theta).max(), c='k', label="pymaster: apotype='Smooth'")
ax2.plot(theta, 1 / 2 * (1 - np.cos(np.pi * x)),              label="Grain et al.: $C^1$-window, eq. (30)")
ax2.plot(theta, x - np.sin(2 * np.pi * x) / (2 * np.pi),      label="Grain et al.: $C^2$-window, eq. (31)")
ax2.axvline(theta_star/2, c='k', ls=':', lw=0.6)

ax1.set_xlabel(r"Distance $\theta$ from masked pixel")
ax2.set_xlabel(r"Distance $\theta$ from masked pixel")
ax1.set_ylabel(r"Mask value")
ax2.set_ylabel(r"Mask value")
ax1.legend(loc='lower right');
ax2.legend(loc='lower right');

Cheers,
Lukas

[damonge]

Hi @lukashergt

Ah, good catch on the C1/C2 thing. That is indeed a typo, and I'm very surprised that it hasn't been caught before.

Indeed, the definition of the apodisation kernels is not exactly the same as used in Grain et al. 2009 (although they should coincide for small angles -- not sure what value of theta_* your plot corresponds to). The main reason for this is efficiency: calculating angles between two pixels requires taking the arccosine of the dot product between their unit vectors, which is a slow operation. The x variable requires only the cosine of that angle, which is more efficient while preserving the differentiability of the kernels. If there is a strong reason for preferring something that adheres exactly to the Grain et al. 2009 definition, it could be implemented.

[lukashergt]

Thanks for the quick reply!

not sure what value of theta_* your plot corresponds to

I was using theta_star = 2.5 * sigma with sigma = 1 (in degrees). The value of 2.5 comes from the apotype='Smooth' recipe, which sets all pixels within 2.5*sigma to zero before smooting. I figured that value was chosen to kinda match the C1 and C2 methods (see right hand plot above).

That said, the problem that you are hinting at is that I did not convert my input units from degree to radians in this little mock example...^^'
Minimal change to sigma makes all the difference:

sigma_deg = 1  # in degree
sigma = np.deg2rad(sigma_deg)

Now we have a more realistically small apodisation scale (yes, an apodisation scale of over 120 degrees would be rather absurd...) and the plots look the same (except for the C1/C2 confusion...):

The main reason for this is efficiency

I see, I was just a little surprised and got curious (after fooling myself into thinking the difference was bigger than it is...).

Thanks! Apologies for the needless confusion... Feel free to close (or keep open if intending to fix the C1/C2 thing).

[damonge]

OK, great. I'll leave it open, since we should fix the C1/C2 confusion. Thanks for catching that!