<p>You don’t say what language your code is in, but there are Python functions to do this in the <a href="https://docs.astropy.org/en/stable/wcs/index.html" rel="nofollow noreferrer">WCS module in astropy</a>.</p>
<p>In particular, look at the <code>pixel_to_world</code> function.</p>
<p>If you want to look at AstroImageJ’s implementation, the code is open source, <a href="https://github.com/AstroImageJ/astroimagej" rel="nofollow noreferrer">here</a>.  The file <a href="https://github.com/AstroImageJ/astroimagej/blob/master/Astronomy_/src/main/java/astroj/WCS.java" rel="nofollow noreferrer">WCS.java</a> has the main world coordinate routines.</p>

CC BY-SA 4.0
2022-07-22T04:17:00.820Z
Thanks for the links, looks like good info! I'm working in Matlab right now, so I was askng to see if there were any algorithms I was missing beyond those I'd already implemented. I didn't think programming language would matter, but I can make a change if it's well solved in another language. I'll check out the links tomorrow.
Conor W
CC BY-SA 4.0
2022-07-22T15:20:24.850Z
In principle you can do it in any language and write your own implementation, but there are enough subtleties and potential gotchas in things like this that it's always better to use a well-tested library if one already exists.
Eric Jensen
CC BY-SA 4.0
2022-07-23T12:28:41.257Z
Thinking more about this, I’d say that the easiest mistake to make in writing your own code for this is to forget the cos(Dec) factor that converts distance on the sky in the RA direction into actual change in RA to apply to your coordinates.  The linear FITS transformations will give you the angular distance components from the reference pixel to your desired pixel.  The change in RA is then that true angular distance (in the RA direction) divided by cos(Dec), accounting for convergence of lines of constant RA as you go from equator to poles.
Eric Jensen
CC BY-SA 4.0
2022-07-24T03:33:16.103Z
I'll admit I haven't had time to check the links, but that sounds useful. While I was comparing the solutions I got from the linear FITS transformation I noticed that my RA would be off from AIJ's value much more than my Dec for a given point, and that the offset differed between images of different star fields.
I looked into correcting for the gnomonic projection because I figured I had an issue with correcting between linear/nonlinear spaces, so this advice looks very promising. Thanks!
Conor W
CC BY-SA 4.0
2022-07-24T20:23:57.193Z
Yup, I think that's what I was looking for, thanks! My error in RA is now at least an order of magnitude smaller than it was without that correction. Let me do a bit more playing around with my data to verify it's repeatable, if it is, I'll edit my question to take this answer into account.
Conor W
CC BY-SA 4.0
2022-07-21T13:11:46.920Z
2022-07-21T13:11:46.920Z
Eric Jensen
<p><strong>EDIT: solved, thanks to Eric Jensen's suggestion (in a comment on his answer) that I include the correction to my right ascension value as a function of my declination. All other comments were helpful but included information I'd already referenced, it was the the correction to RA that I was missing.</strong></p>
<p>As the title states, I'm trying to use the WCS data added to a FITS header by a plate solve to calculate RA and Dec at any given pixel for an image <em>in decimal form</em>. I am aware that programs like AstroImageJ can display this sort of info, but only in HHMSS/DDMMSS format and only by using the cursor. I'd like to automate this with my own code.</p>
<p>Problem summary:</p>
<p>I have a number of plate solved images of star fields with satellite streaks, and I'd like to automate the extraction of satellite RA/Dec from the images based on their positions within the FITS pixel space with a function that allows me to input FITS header data and pixel location to generate RA/Dec. I know that this is a solved problem because AstroImageJ can give me truth data from a plate solved image when I mouse over certain pixels, but I have so many satellite streaks to measure that using the cursor for RA/Dec readout for each point along them would be prohibitively time consuming.</p>
<p>What I've tried so far:</p>
<p>I have tried to apply the transformations from the FITS header data based on <a href="http://montage.ipac.caltech.edu/docs/headers.html" rel="nofollow noreferrer">this Caltech guide</a>, but the algorithm here does not match my data. In particular, the RA/Dec I calculate for a pixel's position using this algorithm is offset by up to an arcminute (in a ~1* FOV frame) from the RA/Dec I get using my cursor in AIJ to find the RA/Dec of the same spot, and I haven't managed to figure out why these transformations aren't working.</p>
<p>I know that these transformations assume a linear model and don't directly include distortion, so I included a correction for the <a href="https://irsa.ipac.caltech.edu/data/SPITZER/docs/files/spitzer/shupeADASS.pdf" rel="nofollow noreferrer">SIP distortion</a> but found that the distortions were sub-pixel in my images and couldn't account for the offset.</p>
<p>I then found <a href="https://photo.stackexchange.com/questions/6111/pixel-to-ra-dec-mapping-in-digitized-astrophotography">this StackExchange post</a> of a user with a similar problem and attempted to apply a correction for the gnomonic projection based on <a href="https://irsa.ipac.caltech.edu/IRASdocs/surveys/coordproj.html" rel="nofollow noreferrer">this algorithm</a>, and while that seemed to work, I'm still off by a few arcminutes.</p>
<p>So, after using the FITS algorithm, SIP distortion correction, and correcting for the gnomonic projection, I've accounted for everything I can think of but I still can't match AstroImageJ's numbers for a plate solved image. Can anyone help me find the missing piece here?</p>
<p>I would also be open to suggestions on alternate methods to get satellite streak RA/Dec from plate solved images, in case there's a simpler approach I'm missing.</p>

CC BY-SA 4.0
2022-07-21T01:03:35.900Z
NASA's FITS WCS page is a good place to start.
https://fits.gsfc.nasa.gov/fits_wcs.html
Greg Miller
CC BY-SA 4.0
2022-07-22T04:20:20.420Z
Thanks! I'd already checked the proper FITS documentation, but I hadn't seen that FITS software link before. I'll check that out.
Conor W
<artificial-satellite><astrometry><algorithm><fits-header>
How can I replicate AstroImageJ's pixel to RA/Dec algorithm in my own code?
Conor W
# How can I replicate AstroImageJ's pixel to RA/Dec algorithm in my own code?

(Asked by: username_0 on 2022-07-20)

**EDIT: solved, thanks to username_2's suggestion (in a comment on his answer) that I include the correction to my right ascension value as a function of my declination. All other comments were helpful but included information I'd already referenced, it was the the correction to RA that I was missing.**

As the title states, I'm trying to use the WCS data added to a FITS header by a plate solve to calculate RA and Dec at any given pixel for an image *in decimal form*. I am aware that programs like AstroImageJ can display this sort of info, but only in HHMSS/DDMMSS format and only by using the cursor. I'd like to automate this with my own code.

Problem summary:

I have a number of plate solved images of star fields with satellite streaks, and I'd like to automate the extraction of satellite RA/Dec from the images based on their positions within the FITS pixel space with a function that allows me to input FITS header data and pixel location to generate RA/Dec. I know that this is a solved problem because AstroImageJ can give me truth data from a plate solved image when I mouse over certain pixels, but I have so many satellite streaks to measure that using the cursor for RA/Dec readout for each point along them would be prohibitively time consuming.

What I've tried so far:

I have tried to apply the transformations from the FITS header data based on [this Caltech guide](http://montage.ipac.caltech.edu/docs/headers.html), but the algorithm here does not match my data. In particular, the RA/Dec I calculate for a pixel's position using this algorithm is offset by up to an arcminute (in a ~1\* FOV frame) from the RA/Dec I get using my cursor in AIJ to find the RA/Dec of the same spot, and I haven't managed to figure out why these transformations aren't working.

I know that these transformations assume a linear model and don't directly include distortion, so I included a correction for the [SIP distortion](https://irsa.ipac.caltech.edu/data/SPITZER/docs/files/spitzer/shupeADASS.pdf) but found that the distortions were sub-pixel in my images and couldn't account for the offset.

I then found [this StackExchange post](https://photo.stackexchange.com/questions/6111/pixel-to-ra-dec-mapping-in-digitized-astrophotography) of a user with a similar problem and attempted to apply a correction for the gnomonic projection based on [this algorithm](https://irsa.ipac.caltech.edu/IRASdocs/surveys/coordproj.html), and while that seemed to work, I'm still off by a few arcminutes.

So, after using the FITS algorithm, SIP distortion correction, and correcting for the gnomonic projection, I've accounted for everything I can think of but I still can't match AstroImageJ's numbers for a plate solved image. Can anyone help me find the missing piece here?

I would also be open to suggestions on alternate methods to get satellite streak RA/Dec from plate solved images, in case there's a simpler approach I'm missing.

**Tags:** artificial-satellite, astrometry, algorithm, fits-header

> **Comments:**
> **username_1**: NASA's FITS WCS page is a good place to start. https://fits.gsfc.nasa.gov/fits_wcs.html  
> **username_0**: Thanks! I'd already checked the proper FITS documentation, but I hadn't seen that FITS software link before. I'll check that out.  

---

## 1 answers

### Answer by username_2 on 2022-07-21. Score: 1 | Accepted answer

You don’t say what language your code is in, but there are Python functions to do this in the [WCS module in astropy](https://docs.astropy.org/en/stable/wcs/index.html).

In particular, look at the `pixel_to_world` function.

If you want to look at AstroImageJ’s implementation, the code is open source, [here](https://github.com/AstroImageJ/astroimagej). The file [WCS.java](https://github.com/AstroImageJ/astroimagej/blob/master/Astronomy_/src/main/java/astroj/WCS.java) has the main world coordinate routines.

> **Comments:**
> **username_0**: Thanks for the links, looks like good info! I'm working in Matlab right now, so I was askng to see if there were any algorithms I was missing beyond those I'd already implemented. I didn't think programming language would matter, but I can make a change if it's well solved in another language. I'll check out the links tomorrow.  
> **username_2**: In principle you can do it in any language and write your own implementation, but there are enough subtleties and potential gotchas in things like this that it's always better to use a well-tested library if one already exists.  
> **username_2**: Thinking more about this, I’d say that the easiest mistake to make in writing your own code for this is to forget the cos(Dec) factor that converts distance on the sky in the RA direction into actual change in RA to apply to your coordinates.  The linear FITS transformations will give you the angular distance components from the reference pixel to your desired pixel.  The change in RA is then that true angular distance (in the RA direction) divided by cos(Dec), accounting for convergence of lines of constant RA as you go from equator to poles.  
> **username_0**: I'll admit I haven't had time to check the links, but that sounds useful. While I was comparing the solutions I got from the linear FITS transformation I noticed that my RA would be off from AIJ's value much more than my Dec for a given point, and that the offset differed between images of different star fields. I looked into correcting for the gnomonic projection because I figured I had an issue with correcting between linear/nonlinear spaces, so this advice looks very promising. Thanks!  
> **username_0**: Yup, I think that's what I was looking for, thanks! My error in RA is now at least an order of magnitude smaller than it was without that correction. Let me do a bit more playing around with my data to verify it's repeatable, if it is, I'll edit my question to take this answer into account.  

---

