Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Inconsistency in position angles ? #613

Closed
mariorpd13 opened this issue Dec 14, 2023 · 5 comments · Fixed by #624
Closed

Inconsistency in position angles ? #613

mariorpd13 opened this issue Dec 14, 2023 · 5 comments · Fixed by #624
Assignees

Comments

@mariorpd13
Copy link

mariorpd13 commented Dec 14, 2023

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]

@hsinfan1996
Copy link
Collaborator

hsinfan1996 commented Dec 15, 2023

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").

@mariorpd13
Copy link
Author

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).

@hsinfan1996
Copy link
Collaborator

hsinfan1996 commented Dec 16, 2023

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.

@mariorpd13
Copy link
Author

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.

@caioolivv
Copy link
Collaborator

caioolivv commented Jul 1, 2024

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.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
4 participants