A few issues with `WCS_AST.c`
First of all, thank you for making your modifications to SoFiA publicly available. In an attempt to make SoFiA easier to install under Windows, I have copied across and updated your AST version of WCS conversion which is now incorporated into the SoFiA master branch (https://gitlab.com/SoFiA-Admin/SoFiA-2/-/blob/master/src/WCS_AST.c).
While assessing the changes that you had made, I found a few minor issues with your current implementation which I have fixed in the regular SoFiA repository:
- In lines 126-128 of
WCS_AST.cyou have removed the destructor calls from WCSLIB, but without adding the respective AST destructor calls. If I am not mistaken then these should beastAnnul(self->fitschan)andastAnnul(self->wcsinfo). - In
WCS_convertToWorld()there is currently a problem with the units. The issue arises because, unlike WCSLIB, AST ignores the units for celestial coordinates defined in the FITS header and always returns celestial coordinated in radians. SoFiA, on the other hand, expects the coordinates to be in the units defined in the header (usually degrees). Hence, with the current implementation SoFiA reports the wrong units in the output catalogue (deg) while the coordinates are actually in radians. I have solved this problem by converting the coordinates from AST to degrees by default and keeping track of which library (WCSLIB or AST) is actually supplying the coordinates so the correct units can be reported in the source catalogue. - My version of GCC refuses to compile line 264 in
WCS_AST.c, because the input coordinate array supplied toastTranP8()is expected to be of typeconst double *, while it is currently declared asdouble *in line 241. Adding aconstqualifier to the array declaration solved this problem and the code compiled fine after that. The same applies to line 305. - The code in lines 238-240 of
WCS_AST.cis problematic. The expression1.0 + n_axes > 0 ? x : 0does not behave as expected, because the+operator has a higher precedence than the ternary operator?:; hence the expression checks the value of1.0 + n_axesagainst zero and returns eitherxor0rather thanx + 1or1. This means that all coordinate transformations should currently be off by 1 pixel or channel. The correct expression would ben_axes > 0 ? x + 1.0 : 1.0(and likewise for the other two axes).
Thank you once more for sharing the AST implementation. This has been very helpful in getting AST support implemented in SoFiA. I have also implemented buffered storage in SoFiA’s dynamic array container which should address the issue of array handling being slow in Windows. This should hopefully resolve the most significant problems that you had encountered while trying to get SoFiA installed in Windows, but please let me know if there are other improvements that we can add to the standard version of SoFiA.