CLMM
CLMM copied to clipboard
LSSTDESC
Inconsistency in position angles ?
I found a possible inconsistency in how azimuthal angles (phi) are computed in clmm.dataops (e.g. by _compute_lensing_angles_astropy).
I would compute phi as (sk0 and sk1 being the coordinates of lens and sources): phi=(np.pi/2-sk0.position_angle(sk1).rad) phi=Angle(phi, u.rad).wrap_at(360 * u.deg).rad but this does not give the same result as CLMM
To give an example, the script below compares the values computed as described above and by CLMM, and those given in https://github.com/CosmoStat/lenspack/blob/master/lenspack/geometry/measures.py (printed in square brackets in the script)
from astropy.coordinates import Angle
from astropy.coordinates import SkyCoord
import numpy as np
from astropy import units as u
from clmm.dataops import _compute_lensing_angles_astropy
sk0=SkyCoord(0,0, unit='deg')
for sk1 ,refval in zip([SkyCoord(-30,-30, unit='deg'),SkyCoord(-1,-1, unit='deg') ],[229.10660535086907, 225.00436354465515]):
phi=(np.pi/2-sk0.position_angle(sk1).rad)
phi1=Angle(phi, u.rad).wrap_at(360 * u.deg).rad
_,phi2=_compute_lensing_angles_astropy(sk0.ra.deg, sk0.dec.deg, sk1.ra.deg, sk1.dec.deg)
print ("CLMM: {:.4f}, Proposed: {:.4f} [{:.4f}]".format(np.degrees(phi2), np.degrees(phi1), refval))
Giving as output: CLMM: -49.1066, Proposed: 229.1066 [229.1066] CLMM: -45.0044, Proposed: 225.0044 [225.0044]
CLMM measures the angle the same way as the polar coordinate in math, which is similar to starting at negative RA axis and increasing in the counter-clockwise direction. And then it puts the angle in between -180 to 180 degrees. Astropy position_angle measures the angle East-of-North on the celestial sphere, which starts from positive Dec axis and similarly increases in CCW direction (see positive RA axis as the negative x axis and positive Dec axis as the positive y axis in the Cartesian coordinate). The proposed calculation flips the direction from CCW to CW and then adds 90 degrees. I am wondering how you calculated the refval. I got 220.89 degree, instead of 229.1066, using Astropy position_angle for the case of (RA, Dec)=("-30d", "-30d").
Hi the refval are the values given as example in the link I gave (function spherical_polar_angle). I put that just as comparison. 220.89 is East-of-North, whereas 229.1066 is the angle staring from RA. It is just a matter of convention, but it may produce inconsistencies depending on the ellipticity definition (e.g. in sky or pixel coordinates).
Indeed the proposed calculation is the same as the example you posted. Part of the reason why CLMM performs such transformation is that we wanted the lensing angle calculation from astropy to be compatible with the older flat-sky approximation. I think @m-aguena can provide more detail on this part as well as opinions on the inconsistency you talk about.
Just another thing: the two angles are such that phi1+phi2=180. Since the main effect of the different definition is on how the tangential/cross ellipticity components are computed, maybe the solution could be to add an option in the function where this is done (something like: if "option" is True, then use phi=180-phi). Then it would be up to the user to decide whether to use this option, according to how ellipticities are defined in the his/her data.
Just to add to this conversation: HSC shear catalogs provide e1 and e2 in sky coordinates, which means all tangential and cross componentes of ellipticity calculated by CLMM from HSC data are wrong. In specific, the notebook CLMM/examples/mass_fitting/Example4_Fit_Halo_mass_to_HSC_data.ipynb provides a wrong analysis of this data. I think the user should be able to tell the GalaxyCluster object which coordinate system it should use to avoid problems such as this.
