JuliaAstro / FITSIO.jl

Flexible Image Transport System (FITS) file support for Julia
http://juliaastro.org/FITSIO.jl/
MIT License
55 stars 29 forks source link

read and write StridedArrays as images #134

Closed jishnub closed 4 years ago

jishnub commented 4 years ago

Following the discussion in #131 , this is an initial attempt to replace the messy unions of contiguous types by StridedArrays. The cfitsio library seems to require a contiguous memory layout for the arrays to be written and read into, however AFAIK there isn't an obvious contiguous StridedArray type to dispatch on. To get around this currently we carry out a runtime check to ensure that the strides along the various axes are compatible with a contiguous layout. For the most part this seems to work, however further tests might be warranted.

julia> fitsfile=tempname();

julia> b = reshape([1:8;],2,2,2);

julia> b_view = @view b[:,:,1]
2×2 view(::Array{Int64,3}, :, :, 1) with eltype Int64:
 1  3
 2  4

julia> FITS(fitsfile,"w") do f
       write(f, b_view)
       end

julia> FITS(fitsfile,"r") do f
       read(f[1])
       end
2×2 Array{Int64,2}:
 1  3
 2  4

# Non-contiguous views throw an error
julia> b_view = @view b[:,1,:]
2×2 view(::Array{Int64,3}, :, 1, :) with eltype Int64:
 1  5
 2  6

julia> FITS(fitsfile,"w") do f
       write(f, b_view)
       end
ArgumentError: data to be written out needs to be contiguously stored

# Complex arrays may be reinterpreted as float as in #131 
julia> FITS(fitsfile,"w") do f
       write(f,reinterpret(Float64,zeros(ComplexF64,1,2)))
       end

julia> FITS(fitsfile,"r") do f
       read(f[1])
       end
2×2 Array{Float64,2}:
 0.0  0.0
 0.0  0.0
giordano commented 4 years ago

Have you tried to do some benchmarks? I'd expect the overhead of the runtime check to be negligible compared to the rest