FITS: BZERO/BSCALE are ignored and 16-bit data is read as unsigned #10309

Activity

1. sunxiuguo commented on Oct 25, 2026

   

   sunxiuguo

   on Oct 25, 2026

   More actions

   One constraint worth keeping explicit when choosing between physical and array values: scaled 32-bit data cannot always be represented losslessly by Pillow's F mode.

   For example, BITPIX=32, BSCALE=1, BZERO=0.5, and stored samples [16777216, 16777217, 16777218] have distinct physical values [16777216.5, 16777217.5, 16777218.5]. Putting those values into an F image produces [16777216.0, 16777218.0, 16777218.0]:

   from PIL import Image
   im = Image.new("F", (3, 1))
   im.putdata([16777216.5, 16777217.5, 16777218.5])
   print([im.getpixel((x, 0)) for x in range(3)])

   I ran this with the official Pillow 13.0.0 wheel on Linux/Python 3.12.14. This tests the proposed representation, not a patch to FITS decoding. Astropy's dtype selection selects float64 for this scaled 32-bit case, whereas Pillow's F mode is float32.

   So option 1 would need a documented precision boundary rather than a blanket match to Astropy; exposing scaling metadata with correctly signed array values would retain the caller's choice of representation. I have not implemented either option or tested a full main-branch build. Prepared with AI assistance; the example was executed.
2. radarhere commented on Oct 26, 2026

   

   radarhere

   on Oct 26, 2026

   Member

   More actions

   Unless you tell me that there are times when someone would want the original array values, Option 1 sounds like the most natural to me. Changing the mode of existing images is not a deal breaker, although I will consider whether data that happens to fit within the I;16 range can be left as I;16. Adding a new unpacker is not a huge deal.
3. radarhere commented on Oct 30, 2026

   

   radarhere

   on Oct 30, 2026

   Member

   More actions

   @bacescuandrei do you approve of the approach of #10312?
4. bacescuandrei commented on Nov 1, 2026

   

   bacescuandrei

   on Nov 1, 2026

   Author

   More actions

   Partly. The decoding in #10312 looks right and I'd keep it, but I wouldn't merge the approach as it stands. It is option 2, and with the results below I agree with you that option 1 is the better default.

   Tested at 0af16d7 on x86_64 with Python 3.12, astropy 8.0.2 and numpy 2.5.4. Tests/test_file_fits.py passes (12) and so does selftest.py (59). Big-endian wasn't run outside CI, but the s390x job passes at this commit.

   What works. These all match astropy's unscaled arrays: uncompressed BITPIX 8/16/32/-32/-64, and GZIP_1 8/16/32/-32 with consecutive samples. For 16-bit gzip that includes a 2x3 image with one tile and with row tiles, flipped vertically the way Pillow reads FITS. BZERO/BSCALE parse correctly as integers, decimals, D exponents and with a trailing comment, and 9007199254740993 stays exact.

   What I'd change:

   1. Saving as PNG regresses. Every BITPIX=16 file now opens as I, including positive data that main already reads correctly. On main that image saves as an I;16 PNG. On the PR, im.save(fp, "PNG") raises OSError: cannot write mode I as PNG.
   2. Camera frames still aren't usable without extra work. A frame with true values [0, 1000, 20000, 32768, 50000, 65535] stored with BZERO=32768 opens as [-32768, -31768, -12768, 0, 17232, 32767], and convert("L") gives [0, 0, 0, 0, 255, 255]. That is the most common case in the issue, and the caller now has to know to apply im.info["BZERO"].
   3. Truncated 16-bit gzip data is accepted. On a hand-built 5x1 GZIP_1 file, an 8-byte payload gives [258, 772, 1286, 1800, 0], because int.from_bytes(b"", ...) returns 0. A 16-byte payload is read as consecutive samples, [258, 772, 1286, 1800, 2314], ignoring the last 6 bytes. Main raises ValueError: not enough image data for both, and 8-bit and 32-bit still do on the PR. Accepting only the two expected payload lengths would fix it.
   4. There is no release note for the mode change and the new im.info keys, and codecov/patch is still red.

   Your earlier question: would anyone want the original array values? Yes, in the two cases where the physical values don't fit any Pillow mode:

   • Unsigned 32-bit (BZERO=2147483648). I is signed, so im.point() with that offset gives [0, 2147483647, -2147483648, -2147483648] where astropy gives [0, 2147483647, 2147483648, 4294967295] as uint32.
   • Scaled 32-bit. F is float32, so @sunxiuguo's [16777216.5, 16777217.5, 16777218.5] becomes [16777216.0, 16777218.0, 16777218.0]. astropy uses float64 here.

   What I'd suggest is option 1 where the result is exact, with this PR's behaviour as the fallback:

   header	mode	values
   BITPIX=16, BZERO=32768, BSCALE=1	I;16	physical
   BITPIX=16, no scaling	I	as stored, signed
   BITPIX=8 or 16, integer scaling that fits	I	physical
   BITPIX=8 or 16, non-integer scaling	F	physical (astropy also uses float32 for BITPIX=16)
   BITPIX=32, no scaling	I	as stored
   BITPIX=32, any scaling	I	as stored, with BZERO/BSCALE in im.info

   The first row covers the camera case and keeps the mode those files have today. It only needs the sign bit flipped in the raw bytes. A throwaway prototype of that (XOR on the bytes, then I;16B) matched astropy on [0, 1, 1000, 20000, 32767, 32768, 50000, 65535] and round-tripped through PNG. It isn't a patch against the plugin, so a proper unpacker is still to be written.

   One limit remains. The mode has to be chosen at open(), before the pixels are read, so BITPIX=16 without BZERO has to be I, and saving those as PNG needs im.convert("I;16"). That works for positive data, but it clips negatives to 0 ([-5, 0, 7] becomes [0, 0, 7]).

   If you'd rather merge option 2 as a first step, I'd want at least 3 and the release note fixed first, and 1 called out in that note.