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

LinAlgError: Singular matrix #7

Open
TrystanScottLambert opened this issue May 19, 2023 · 3 comments
Open

LinAlgError: Singular matrix #7

TrystanScottLambert opened this issue May 19, 2023 · 3 comments

Comments

@TrystanScottLambert
Copy link

Hi.

I'm getting the following error when trying to use sip_tpv.sip_to_tpv()

`LinAlgError Traceback (most recent call last)
in
----> 1 sip_tpv.sip_to_pv(header)

~/.local/lib/python3.10/site-packages/sip_tpv/sip_to_pv.py in sip_to_pv(header, tpv_format, preserve)
56 pvrange, tpvx, tpvy = sym_tpvexprs()
57 cd, ac, bc = get_sip_keywords(header)
---> 58 sipx, sipy = real_sipexprs(cd, ac, bc)
59 add_pv_keywords(header, sipx, sipy, pvrange, tpvx, tpvy,
60 int(header['B_ORDER']))

~/.local/lib/python3.10/site-packages/sip_tpv/pvsiputils.py in real_sipexprs(cd, ac, bc)
262 """
263 x, y = symbols("x y")
--> 264 cdinverse = cd**-1
265 uprime, vprime = cdinverse*Matrix([x, y])
266 usum = uprime

~/.local/lib/python3.10/site-packages/numpy/matrixlib/defmatrix.py in pow(self, other)
229
230 def pow(self, other):
--> 231 return matrix_power(self, other)
232
233 def ipow(self, other):

~/.local/lib/python3.10/site-packages/numpy/core/overrides.py in matrix_power(*args, **kwargs)

~/.local/lib/python3.10/site-packages/numpy/linalg/linalg.py in matrix_power(a, n)
648
649 elif n < 0:
--> 650 a = inv(a)
651 n = abs(n)
652

~/.local/lib/python3.10/site-packages/numpy/core/overrides.py in inv(*args, **kwargs)

~/.local/lib/python3.10/site-packages/numpy/linalg/linalg.py in inv(a)
550 signature = 'D->D' if isComplexType(t) else 'd->d'
551 extobj = get_linalg_error_extobj(_raise_linalgerror_singular)
--> 552 ainv = _umath_linalg.inv(a, signature=signature, extobj=extobj)
553 return wrap(ainv.astype(result_t, copy=False))
554

~/.local/lib/python3.10/site-packages/numpy/linalg/linalg.py in _raise_linalgerror_singular(err, flag)
87
88 def _raise_linalgerror_singular(err, flag):
---> 89 raise LinAlgError("Singular matrix")
90
91 def _raise_linalgerror_nonposdef(err, flag):

LinAlgError: Singular matrix
`

This arrises using the attached header.
Steps to reproduce:
header.txt

from sip_tpv import sip_to_tpv
from astropy.io import fits 
header  = fits.Header.fromtextfile('header.txt')
sip_to_pv(header)

The header was generated using the wcs.utils.fit_wcs_from_points function in astropy. Then writing that wcs object to header using the "relax=True" parameter (in order to actually write the keywords to the header)

Any ideas as to what is causing the crash? Am I just not running it correctly?
Thanks in advanced.

@TrystanScottLambert
Copy link
Author

I've realized what the issue is, namely that sip_pv looks for CD matrix keywords. I just set them equal to their PC counterparts and it runs but the solution isn't correct.

Is there any ideas to properly convert from PC matrix representations to CD matrix representations?

@stargaser
Copy link
Owner

@TrystanScottLambert thank you for the issue and for supplying a complete example. The CD matrix is obtained by multiplying the PC matrix values with the CDELT values. In the example header, the CDELT values are 1.0. So you can set the CD matrix values to the PC matrix values, and delete the PC matrix values and the CDELT values. I tried this in the attached text file header_cd.txt.

I made WCS objects out of the headers and I get the same transformation values out of the pixel_to_world method. I admit the SIP headers are not roundtripping completely in the world_to_pixel direction.

@TrystanScottLambert
Copy link
Author

@stargaser Great. Thanks for the reply. It runs, but now I'm getting some weird distortions after using swarp to Mosaic them together. Thanks for the help.

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

No branches or pull requests

2 participants