Skip to content

fix!: Fix header value precision and HDU vector dispatch, add Base.size(::HDU) - #45

Draft
icweaver wants to merge 2 commits into
mainfrom
hdu-size
Draft

fix!: Fix header value precision and HDU vector dispatch, add Base.size(::HDU)#45
icweaver wants to merge 2 commits into
mainfrom
hdu-size

Conversation

@icweaver

@icweaver icweaver commented Jul 23, 2026

Copy link
Copy Markdown
Member

Hi Paul, these are three fixes/additions that surfaced while I was porting Reproject.jl over to FITSFiles.jl. The first changes the types of some parsed card values, so I've marked this as a breaking release: 0.3.2 --> 0.4.0 for your review.

1. Header real values are no longer silently corrupted (parse_number)

parse_number picked Float32 vs. Float64 by counting fractional digits (10 or more --> Float64), which seemed to unintentionally introduce a bug. Any value whose decimal representation didn't fit Float32 was silently corrupted:

Before:

hdu.cards["CRVAL1"] # Header says 266.400000
266.39999f0 # Float32, off by ~6e-6 deg (~0.02 arcsec)

After:

hdu.cards["CRVAL1"]
266.4 # Float64, exact

The parser now reads every real at Float64 precision and narrows to Float32 only when the narrowed value represents the written decimal exactly (Float64(Float32(x)) == x), preserving the smallest-type behavior for values like 1.0 or 0.5 while never losing information.

D-exponent values always parse as Float64 per the FITS convention, and E-notation values get the same round-trip treatment (previously, e.g., 1E30 also parsed lossily as Float32). The now-unused overflow helper is removed.

Downstream consequence: Scale keywords like BSCALE = 0.1 are now Float64, so scaled integer data promotes to Float64 (e.g., Int32 data with BSCALE = 0.1, BZERO = 1.0 now reads as exactly 1.1 rather than 1.1f0), matching what astropy and FITSIO produce. Test expectations that encoded the old lossy values are updated accordingly.

2. HDU vector methods now accept concrete vectors

All the Vector{HDU} signatures in fitscore.jl (write, info, show, haskey, getindex, get, findfirst) only matched vectors with the abstract element type, because parametric invariance means Vector{HDU{Primary}} is not a Vector{HDU}. Vectors returned by fits() happened to dispatch fine, but user-constructed ones did not:

write(fname, [HDU(data, cards)])
# ERROR: `write` is not supported on non-isbits arrays  (fell through to Base.write)

All signatures are now widened to AbstractVector{<:HDU}, with an additional Vector{<:HDU} method for write so it takes precedence over Base.write(::IO, ::Array) in dispatch.

3. New: Base.size(hdu)

Base.size(hdu::HDU{<:Union{Primary, Image}}) returns the data dimensions from the NAXISn cards without reading the data, so callers holding a lazy-loaded HDU can get its shape for free (Reproject.jl uses this to size output projections).

@codecov

codecov Bot commented Jul 23, 2026

Copy link
Copy Markdown

Codecov Report

❌ Patch coverage is 57.14286% with 9 lines in your changes missing coverage. Please review.
✅ Project coverage is 46.94%. Comparing base (870e4ee) to head (fbd4a23).

Files with missing lines Patch % Lines
src/fitscore.jl 35.71% 9 Missing ⚠️
Additional details and impacted files
@@            Coverage Diff             @@
##             main      #45      +/-   ##
==========================================
+ Coverage   46.35%   46.94%   +0.58%     
==========================================
  Files          17       17              
  Lines        2550     2550              
==========================================
+ Hits         1182     1197      +15     
+ Misses       1368     1353      -15     

☔ View full report in Codecov by Harness.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:
  • ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

@barrettp

Copy link
Copy Markdown
Member
  1. This fix presents a dilemma, because now most image arrays will be expanded to 64-bit floats which negates the goal of minimizing array sizes. A better solution needs to be found.

  2. I don't understand this fix, because there were no problems with the code working properly before.

  3. I don't understand this fix, because info returned the correct dimensions before. What is Reproject.jl?

These PRs do not provide sufficient descriptions to review them properly.

@icweaver

icweaver commented Jul 23, 2026

Copy link
Copy Markdown
Member Author

Hi Paul, sure, here is more detail based on my current understanding. Please correct me if I have misunderstood something. I've also updated the PR description in the Reproject.jl repo to include more detail as well:

  1. The header-value fix doesn't actually change image array sizes in the common case. Array eltypes are governed by BITPIX via BITS2TYPE, so a BITPIX = -32 image still loads as Float32 after this PR. Header card values are scalars, so widening some of them to Float64 costs essentially nothing in memory. The reason the fix matters is astrometric: CRVAL1 = 266.400000 parsed to 266.39999389... (Float32(266.4)), a ~0.02 arcsec error that made every WCS solution built from FITSFiles headers disagree with astropy/wcslib, which (like cfitsio and FITSIO.jl) parse header reals at double precision. And values that are exactly representable still narrow: parse_number("0.5") isa Float32, parse_number("1.0E2") isa Float32. Only values that Float32 cannot represent keep full precision.

    The one place arrays do widen is scaled integer data whose BSCALE/BZERO don't round-trip (e.g. 0.1). If that's the concern, I think there's a solution that serves both goals: compute the scaling at Float64 precision but store the result at the type declared in OUTTYPE (Int16 --> Float32, Int32 --> Float64, etc.). That keeps scaled arrays exactly the sizes this package already documents as intended. Actually, it seems that the previous behavior violated OUTTYPE for scaled Int32 data since it produced Float32 whereas the table says Float64. Happy to add that as a follow-up commit here if you'd like.

  2. It worked for the vectors returned by fits(), which happen to have the abstract element type HDU{<:AbstractHDU}, but not for vectors a user constructs. It looks like this is because of parametric invariance (Vector{HDU{Primary}} is not a subtype of Vector{HDU}).

    Before this PR:

    julia> write("out.fits", [HDU(data, cards)])
    ERROR: `write` is not supported on non-isbits arrays

    The call silently fell through to Base.write(::IO, ::Array) instead of this package's method. So writing a file you'd just read worked, but writing an HDU you'd just constructed did not. The new round-trip test in fits_tests.jl pins exactly this case.

  3. No fix to info intended, its Shape column was and is correct. info prints a human readable summary, Base.size(hdu) returns the dimensions as a Tuple so programs can use them, e.g., sizing an output buffer for a lazily-loaded HDU without forcing the data to be read. I think it complements info the same way Base.size(::Array) complements display.

    As for Reproject.jl: it's our pure-Julia implementation of astropy's reproject package, which transforms/resamples astronomical images between coordinate frames. Here are some usage examples: https://juliaastro.org/Reproject/stable/

Does this seem reasonable to you? Happy to adjust/provide additional examples if it is useful here

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

Successfully merging this pull request may close these issues.

2 participants