Repository navigation
Slicing FITS_rec objects produces unexpected dtypes. #5490
Description
Consider the following snippet:
import numpy as np
from astropy.io import fits
def create_data():
tbhdu = fits.BinTableHDU.from_columns([
fits.Column(name='a', format='I', array=np.arange(2**16, dtype=np.int16)),
fits.Column(name='b', format='I', bscale=1, bzero=2**15, array=np.arange(2**16, dtype=np.uint16))])
tbhdu.writeto('test_scaled_slicing.fits')
def load_data(nonnegative=False):
with fits.open('test_scaled_slicing.fits') as hdulist:
data = hdulist[1].data
if nonnegative:
return data[data['a'] >= 0]
else:
return data
create_data()
t1 = load_data()
t2 = load_data(True)
assert t1['b'].dtype is np.dtype(np.uint16)
# assert t2['b'].dtype is np.dtype(np.float64)
assert t2['b'].dtype is np.dtype(np.uint16)The second assert fails because the slicing the FITS_rec object does not trigger the mechanism that handles the scaling of the column 'b'. In other words, some data access methods on FITS_rec trigger the scaling but some don't.
If however, I use a method that forces the scaling, then this problem goes away. Here's an example:
def load_data(nonnegative=False):
with fits.open('test_scaled_slicing.fits') as hdulist:
data = hdulist[1].data
#
# WORKAROUND
#
if any([c.bzero is not None for c in data.columns]):
foo = [data[k].dtype for k in data.dtype.names]
if nonnegative:
return data[data['a'] >= 0]
else:
return dataThat is, merely accessing the dtype of each column triggers the scaling. Looking into the code itself, it looks like the ._converted attribute needs to be filled in for slicing to work properly.
Activity
Fix uint conversion after slicing
Fix uint conversion after slicing
peppedilillo commentedon Sep 8, 2024
This issue was closed but I can still get the same bug with astropy 6.1.3.
neutrinoceros commentedon Sep 9, 2024
@peppedilillo This is surprising given that the resolving PR contained a regression test.
I also couldn't reproduce the problem with either astropy 6.1.3 or 7.0.0dev.
Can you please open a new issue ?
peppedilillo commentedon Sep 9, 2024
Yes sure i will. Are you sure you removed the comment to the second assert before trying it? I get:
neutrinoceros commentedon Sep 9, 2024
Oh, I didn't get that this was the problematic assert. Indeed I can reproduce it then, and I would recommend re-opening this issue instead of creating a new one.
peppedilillo commentedon Sep 9, 2024
@neutrinoceros I cannot re-open the issue. Likely no permissions.
neutrinoceros commentedon Sep 9, 2024
I know, this was meant for @weaverba137 or @saimn :)
weaverba137 commentedon Sep 9, 2024
I also do not see an obvious way to reopen the issue. I also note that #6837 is still open. Is that still genuinely open?
pllim commentedon Sep 9, 2024
I re-opened this issue as requested but I will defer to @astrofrog and @saimn for #6837
pllim commentedon Sep 9, 2024
Though I wonder if #5490 (comment) is from numpy 2.0 🤔 and also why didn't our regression test suite catch it?
pllim commentedon Sep 9, 2024
@peppedilillo , can you please give us a minimally reproducible code snippet that you encountered? What numpy version? Python version, OS, distro, etc?
neutrinoceros commentedon Sep 10, 2024
It's not, the same exception is raised with numpy 1.26
This is indeed confusing. Simon's test from #8982 makes sense to me. I keep wondering if Benjamin meant to refer to the commented
assertstatement inBecause it doesn't make much sense to me to expect
float64where only16bits data types were used originally.@weaverba137 Can you please clarify ?
peppedilillo commentedon Sep 10, 2024
This outputs:
The first four prints are garbage. The last four are correct, and take a long time in each statement. You can download the input file and the script here. Note that the problem occurs also if
fieldsis set to[PHA], which is signed. Some of these columns are variable length.This occurred to me because I am working on a FITS viewer application for the terminal. The application show FITS tables in short "pages". Each page consists of a slice of
FITS_rec. I could not get the right output for certain input files because of this problem. I'm trying to work on slices because, if some scaling or other slow astropy operation are carried out over the whole table, that would result in annoying delays on the readout.Just tested on Debian 12 bookworm, but observed the same behaviour also on Windows 11 and MacOS Ventura.
Yes i'm on numpy 2 and Python 3.12.5. The whole environment is the following, I'm on anaconda presently:
weaverba137 commentedon Sep 10, 2024
This issue was created in 2016(!) and I do not recall why this line was there:
# assert t2['b'].dtype is np.dtype(np.float64)I would not expect this to work, even after the bug was fixed, and I think that is reasonable.
When I used the term "The second
assertfails...", I currently believe that refers to this line:saimn commentedon Sep 10, 2024
I don't think the data type from the original issue should be float, this line is commented so maybe it was just some trial?
Ah thanks @weaverba137 for commenting :), I was just looking at the issue and writing a comment.
@peppedilillo - your issue is a different one, it has to do with VLA columns. Interestingly if you first access the full columns and then slice it, it works. Please open a new issue.
pllim commentedon Sep 10, 2024
Thanks, all! I opened #16990 for you. Please continue there.
peppedilillo commentedon Sep 10, 2024
I realize that the commented assert should not work indeed. Thank you for having opened the new issue.