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.

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, 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 but found that the distortions were sub-pixel in my images and couldn't account for the offset.

I then found this StackExchange post of a user with a similar problem and attempted to apply a correction for the gnomonic projection based on this algorithm, 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.

  • 1
    $\begingroup$ NASA's FITS WCS page is a good place to start. fits.gsfc.nasa.gov/fits_wcs.html $\endgroup$ Commented Jul 21, 2022 at 1:03
  • $\begingroup$ Thanks! I'd already checked the proper FITS documentation, but I hadn't seen that FITS software link before. I'll check that out. $\endgroup$
    – Conor W
    Commented Jul 22, 2022 at 4:20

1 Answer 1


You don’t say what language your code is in, but there are Python functions to do this in the WCS module in astropy.

In particular, look at the pixel_to_world function.

If you want to look at AstroImageJ’s implementation, the code is open source, here. The file WCS.java has the main world coordinate routines.

  • $\begingroup$ 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. $\endgroup$
    – Conor W
    Commented Jul 22, 2022 at 4:17
  • $\begingroup$ 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. $\endgroup$ Commented Jul 22, 2022 at 15:20
  • $\begingroup$ 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. $\endgroup$ Commented Jul 23, 2022 at 12:28
  • $\begingroup$ 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! $\endgroup$
    – Conor W
    Commented Jul 24, 2022 at 3:33
  • 1
    $\begingroup$ 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. $\endgroup$
    – Conor W
    Commented Jul 24, 2022 at 20:23

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .