# ICRS to galactocentric conversion using Astropy

**URL:** <https://community.openastronomy.org/t/icrs-to-galactocentric-conversion-using-astropy/724>\
**Category:** Astropy\
**Tags:** astropy, table, fits, question\
**Created:** [August 17, 2023, 3:08pm UTC](https://community.openastronomy.org/t/icrs-to-galactocentric-conversion-using-astropy/724 "2023-08-17T15:08:28Z")\
**Posts on this page:** 3\
**Page:** 1

<div class="post-metadata">

**Author:** ![vhpatel2000](https://dub1.discourse-cdn.com/flex005/user_avatar/community.openastronomy.org/vhpatel2000/32/402_2.png) [@vhpatel2000](https://community.openastronomy.org/u/vhpatel2000)\
**Post date:** [August 17, 2023, 3:08pm UTC](https://community.openastronomy.org/t/icrs-to-galactocentric-conversion-using-astropy/724/1 "2023-08-17T15:08:28Z")

</div>

Hi all ,  
I am using astropy package for converting RA, Dec and distance ( Line of Sight) to Galactocentric coordinates. The problem is the following:

I generated a sphere of randomly distributed 1 million points in RA , Dec and distance plot. It is like a synthetic cluster located at 4000 pc and having radius of 1 pc. Now when converted to Galactocentric frame, I am getting all points flattened to an ellipse. I intuitively expected it to still remain a sphere, as we are just looking from a different location , thus a sphere should, remain a sphere. Is this output expected , or am I doing something wrong ? Following is the code and corresponding output:

> from astropy.io.votable import parse, writeto  
> from astropy.table import Table, Column  
> import numpy as np  
> from [astropy.io](http://astropy.io) import fits  
> from astropy.coordinates import SkyCoord, Galactic, CartesianRepresentation  
> from astropy import units as u
> 
> **Generating Synthetic Cluster in ICRS**
> 
> #### Center coordinates
> 
> center\_ra = 343.81252 # degrees  
> center\_dec = 62.53608 # degrees  
> center\_distance = 4000 # parsecs
> 
> #### Number of stars
> 
> num\_stars = 10\*\*6
> 
> #### Radius of the sphere
> 
> radius = 1 # this is in pc
> 
> distances = (np.random.uniform(0, 1, num\_stars))\*\*(1/3) \* radius
> 
> #### Generating random spherical coordinates
> 
> theta = np.random.uniform(0, 2\*np.pi, num\_stars)  
> phi = np.arccos(2 \* np.random.uniform(0, 1, num\_stars) - 1)
> 
> #### Converting spherical coordinates to Cartesian coordinates
> 
> ra = center\_ra + np.arctan(distances \* np.sin(phi) \* np.cos(theta)/center\_distance)  
> dec = center\_dec + np.arctan(distances \* np.sin(phi) \* np.sin(theta)/center\_distance)  
> distance = center\_distance + distances \* np.cos(phi)
> 
> #### Creating a FITS table with the generated data (RA, Dec and Distance)
> 
> hdu = fits.BinTableHDU.from\_columns([  
> fits.Column(name=‘RA’, format=‘D’, array=ra),  
> fits.Column(name=‘Dec’, format=‘D’, array=dec), This text will be hidden  
> fits.Column(name=‘Distance’, format=‘D’, array=distance)  
> ])
> 
> #### Save the FITS table to a file
> 
> hdu.writeto(‘synthetic\_star\_cluster.fits’, overwrite=True)

Output in RA, Dec and Distance plot is as expected ( a perfect spherical synthetic cluster, refer the figure at the end.)

Now I am reading this fits file to convert coordinates to galactocentric frame :

> import astropy.units as u  
> import astropy.coordinates as coord
> 
> fits\_data = fits.getdata(‘synthetic\_star\_cluster.fits’)  
> table = Table(fits\_data)
> 
> #### Extracting data from the loaded data array
> 
> ra = table[‘RA’]  
> dec = table[‘Dec’]  
> distance = table[‘Distance’]
> 
> c = coord.SkyCoord(ra=ra \* u.degree,  
> dec=dec \* u.degree,  
> distance=(distance) \* u.pc,  
> frame=‘icrs’)  
> galactocentric\_cartesian = c.transform\_to(coord.Galactocentric(galcen\_distance = `8.3 * u.kpc` , z\_sun = `27 * u.pc`))
> 
> #### Define new columns for galactocentric coordinates
> 
> X\_column = Column(data=galactocentric\_cartesian.x.value, name=‘X’, dtype=‘float64’)  
> Y\_column = Column(data=galactocentric\_cartesian.y.value, name=‘Y’, dtype=‘float64’)  
> Z\_column = Column(data=galactocentric\_cartesian.z.value, name=‘Z’, dtype=‘float64’)
> 
> #### Append the new column to the table
> 
> table.add\_column(X\_column)  
> table.add\_column(Y\_column)  
> table.add\_column(Z\_column)
> 
> #### Write the modified Astropy Table to a FITS file
> 
> table.write(‘synthetic\_XYZ\_out’, format=‘fits’, overwrite=True)

The figure below clearly shows that a spherical synthetic cluster in ICRS frame is a flattened ellipse in galactocentric frame , shouldn’t it remain a sphere? . Why is this happening?

 ![image](https://europe1.discourse-cdn.com/flex005/uploads/openastronomy/original/1X/8782fddecffd77b850c91a3e14b92fd1a491092c.png)

---

<div class="post-metadata">

**Author:** ![ayshih](https://dub1.discourse-cdn.com/flex005/user_avatar/community.openastronomy.org/ayshih/32/84_2.png) [@ayshih](https://community.openastronomy.org/u/ayshih)\
**Post date:** [August 17, 2023, 4:52pm UTC](https://community.openastronomy.org/t/icrs-to-galactocentric-conversion-using-astropy/724/2 "2023-08-17T16:52:07Z")

</div>

> [@vhpatel2000](#):
>
> #### Converting spherical coordinates to Cartesian coordinates
> 
> ra = center\_ra + np.arctan(distances \* np.sin(phi) \* np.cos(theta)/center\_distance)  
> dec = center\_dec + np.arctan(distances \* np.sin(phi) \* np.sin(theta)/center\_distance)  
> distance = center\_distance + distances \* np.cos(phi)

Your equations for the modified RA and declination are not correct for adding two spherical vectors. RA and declination are coupled, and you’re missing the `cos(dec)` factor that converts between an angular offset in the RA direction and the corresponding offset in RA units.

I suggest you leverage Astropy’s spherical representations and being able to add them:

```python
from astropy.coordinates import SphericalRepresentation

center = SphericalRepresentation(center_ra, center_dec, center_distance)
offset = SphericalRepresentation(theta, phi - np.pi/2, distances)

c = coord.SkyCoord(center + offset, ...)

```

---

<div class="post-metadata">

**Author:** ![vhpatel2000](https://dub1.discourse-cdn.com/flex005/user_avatar/community.openastronomy.org/vhpatel2000/32/402_2.png) [@vhpatel2000](https://community.openastronomy.org/u/vhpatel2000)\
**Post date:** [August 18, 2023, 8:25pm UTC](https://community.openastronomy.org/t/icrs-to-galactocentric-conversion-using-astropy/724/3 "2023-08-18T20:25:12Z")

</div>

Thanks @ayshih , that was the right solution 🙂  
Sphere does transform to a sphere !