Reputation: 2399
I have an astronomical image fits
format with WCS (World Coordinate System) in its header.
I want to do some morphological operations on the image and understandably I want to correct the WCS.
Say I crop the image. Applying to the header would be:
from astropy.io import fits
from astropy.wcs import WCS
header: Header = fits.getheader("sample.fits")
x, y, w, h = 100, 100, 200, 100
wcs = WCS(header)
new_wcs = wcs[y:y + h, x:x + w]
I can confirm new_wcs is correct by checking the new pixel to sky coordinates:
print(wcs.pixel_to_world(2, 2))
print(new_wcs.pixel_to_world(2, 2))
The only problem here would be that I lost the information in the original header. So I tried to update the old header.
added_header: Header = header.copy()
added_header.extend(new_wcs.to_header(), unique=True, update=True)
added_header_wcs = WCS(added_header)
The whole code would be:
from astropy.io.fits import Header
from astropy.io import fits
from astropy.wcs import WCS
header: Header = fits.getheader("sample.fits")
x, y, w, h = 100, 100, 200, 100
wcs = WCS(header)
new_wcs = wcs[y:y + h, x:x + w]
added_header: Header = header.copy()
added_header.extend(new_wcs.to_header(), unique=True, update=True)
added_header_wcs = WCS(added_header) # boom
print(wcs.pixel_to_world(2, 2))
print(new_wcs.pixel_to_world(2, 2))
print(added_header_wcs.pixel_to_world(2, 2))
However at the added_header_wcs = WCS(added_header)
line python crashes:
Process finished with exit code -1073741819 (0xC0000005)
I tried to enable logging for astropy, but I was not able to do so.
Before creating an issue I wanted to ask more knowledgeable people. Where do I do wrong?
PS: copying a new header is not the problem, as I tried to modify the original header.
PS: the sample data: http://www.astropy.org/astropy-data/tutorials/FITS-images/HorseHead.fits
Executing the same code with a none-conda interpreter threw an exception:
WARNING: FITSFixedWarning: 'celfix' made the change 'Invalid parameter value'. [astropy.wcs.wcs]
Traceback (most recent call last):
File "[PATH]\main.py", line 15, in <module>
print(WCS(new_header).pixel_to_world(2, 2))
^^^^^^^^^^^^^^^
File "[PATH]\local-packages\Python311\site-packages\astropy\wcs\wcs.py", line 600, in __init__
self.wcs.set()
ValueError: ERROR 5 in wcsset() at line 2808 of file cextern\wcslib\C\wcs.c:
Invalid parameter value.
ERROR 4 in linset() at line 737 of file cextern\wcslib\C\lin.c:
Failed to initialize distortion functions.
ERROR 3 in tpdset() at line 1945 of file cextern\wcslib\C\dis.c:
Unrecognized field name for TPD on axis 1: DQ1.DSS.AMD.1.
Upvotes: 0
Views: 46