Skip to content

Commit

Permalink
Read/write for complex arrays, compat bounds
Browse files Browse the repository at this point in the history
Complex arrays may be reinterpreted as float
and written to FITS. Analogously floating point
arrays may be read from FITS by reinterpreting
them as (re,im) pairs stored contiguously. The
first dimension necessarily needs to be even
for this to work.

Shift from inequalities to caret bounds for
compatibilities of dependencies
  • Loading branch information
jishnub committed May 3, 2020
1 parent b59401a commit 803b883
Show file tree
Hide file tree
Showing 5 changed files with 103 additions and 16 deletions.
4 changes: 2 additions & 2 deletions Project.toml
Expand Up @@ -8,8 +8,8 @@ Libdl = "8f399da3-3557-5675-b5ff-fb832c97cbdb"
Printf = "de0858da-6303-5e67-8744-51eddeeeb8d7"

[compat]
BinaryProvider = "0.4.0"
julia = "≥ 1.0.0"
BinaryProvider = "0.4, 0.5"
julia = "1"

[extras]
Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
Expand Down
4 changes: 3 additions & 1 deletion src/FITSIO.jl
Expand Up @@ -43,7 +43,9 @@ using .Libcfitsio
import .Libcfitsio: libcfitsio,
fits_assert_ok,
fits_assert_isascii,
TYPE_FROM_BITPIX
TYPE_FROM_BITPIX,
ContiguousAbstractArray


# HDU Types
abstract type HDU end
Expand Down
12 changes: 7 additions & 5 deletions src/image.jl
Expand Up @@ -243,10 +243,11 @@ following array element types are supported: `UInt8`, `Int8`,
`Float64`. If a `FITSHeader` object is passed as the `header` keyword
argument, the header will also be added to the new HDU.
"""
function write(f::FITS, data::Array{T};
header::Union{Nothing, FITSHeader}=nothing,
name::Union{Nothing, String}=nothing,
ver::Union{Nothing, Integer}=nothing) where T
function write(f::FITS, data::ContiguousAbstractArray{T};
header::Union{Nothing, FITSHeader}=nothing,
name::Union{Nothing, String}=nothing,
ver::Union{Nothing, Integer}=nothing) where {T<:Real}

fits_assert_open(f.fitsfile)
s = size(data)
fits_create_img(f.fitsfile, T, [s...])
Expand All @@ -268,7 +269,8 @@ end
Write data to an existing image HDU.
"""
function write(hdu::ImageHDU, data::Array{T}) where T
function write(hdu::ImageHDU, data::ContiguousAbstractArray{T}) where {T<:Real}

fits_assert_open(hdu.fitsfile)
fits_movabs_hdu(hdu.fitsfile, hdu.ext)

Expand Down
20 changes: 13 additions & 7 deletions src/libcfitsio.jl
Expand Up @@ -120,8 +120,12 @@ export FITSFile,
using Libdl

import Base: FastContiguousSubArray
const ArrayOrFastContiguousSubArray{T} =
Union{Array{T},FastContiguousSubArray{T,N,<:Array{T}} where N}

const ArrayOrFastContiguousSubArray{T,N} = Union{Array{T,N},
FastContiguousSubArray{T,N,<:Array{T}}}

const ContiguousAbstractArray{T,N} = Union{ArrayOrFastContiguousSubArray{T,N},
Base.ReinterpretArray{T,N,Complex{T},<:ArrayOrFastContiguousSubArray{Complex{T},N}}}

const depsjl_path = joinpath(@__DIR__, "..", "deps", "deps.jl")
if !isfile(depsjl_path)
Expand Down Expand Up @@ -751,7 +755,9 @@ end
Write pixels from `data` into the FITS file.
"""
function fits_write_pix(f::FITSFile, fpixel::Vector{S},
nelements::Integer, data::Array{T}) where {S<:Integer,T}
nelements::Integer,
data::Union{Array{T},Base.ReinterpretArray{T}}) where {S<:Integer,T}

status = Ref{Cint}(0)
ccall((:ffppxll, libcfitsio), Cint,
(Ptr{Cvoid}, Cint, Ptr{Int64}, Int64, Ptr{Cvoid}, Ref{Cint}),
Expand All @@ -766,7 +772,7 @@ end

function fits_read_pix(f::FITSFile, fpixel::Vector{S},
nelements::Int, nullval::T,
data::ArrayOrFastContiguousSubArray{T}) where {S<:Integer,T}
data::ContiguousAbstractArray{T}) where {S<:Integer,T}
anynull = Ref{Cint}(0)
status = Ref{Cint}(0)
ccall((:ffgpxvll, libcfitsio), Cint,
Expand All @@ -785,7 +791,7 @@ Read pixels from the FITS file into `data`.
"""
function fits_read_pix(f::FITSFile, fpixel::Vector{S},
nelements::Int,
data::ArrayOrFastContiguousSubArray{T}) where {S<:Integer,T}
data::ContiguousAbstractArray{T}) where {S<:Integer,T}
anynull = Ref{Cint}(0)
status = Ref{Cint}(0)
ccall((:ffgpxvll, libcfitsio), Cint,
Expand All @@ -797,14 +803,14 @@ function fits_read_pix(f::FITSFile, fpixel::Vector{S},
anynull[]
end

function fits_read_pix(f::FITSFile, data::ArrayOrFastContiguousSubArray)
function fits_read_pix(f::FITSFile, data::ContiguousAbstractArray)
fits_read_pix(f, ones(Int64,length(size(data))), length(data), data)
end

function fits_read_subset(
f::FITSFile, fpixel::Vector{S1}, lpixel::Vector{S2},
inc::Vector{S3},
data::ArrayOrFastContiguousSubArray{T}) where {S1<:Integer,S2<:Integer,S3<:Integer,T}
data::ContiguousAbstractArray{T}) where {S1<:Integer,S2<:Integer,S3<:Integer,T}

anynull = Ref{Cint}(0)
status = Ref{Cint}(0)
Expand Down
79 changes: 78 additions & 1 deletion test/runtests.jl
Expand Up @@ -155,7 +155,84 @@ using Random # for `randstring`
finally
rm(fname1, force=true)
end
rm(fname1, force=true)
end

@testset "reinterpreted complex" begin
@testset "read" begin
fname = tempname() * ".fits"

function _testreadcomplex(arr::AbstractArray{Complex{T}}) where T
FITS(fname,"w") do f
write(f, collect(reinterpret(T,arr)) )
end
b = similar(arr)
FITS(fname,"r") do f
read!(f[1], reinterpret(T,b))
end
@test b == arr
c = FITS(fname,"r") do f
read(f[1])
end
@test c == reinterpret(T,arr)
end

function testreadcomplex(arr::Array{Complex{T}}) where T
# Given a complex array, reinterpret it as float and write to file
# Read it back and compare
_testreadcomplex(arr)

# Write a contiguous subsection instead of the entire array
region = 1:size(arr,1)
arrv = @view arr[region]
_testreadcomplex(arrv)
end

try
testreadcomplex(ones(ComplexF64,3))
testreadcomplex(ones(ComplexF64,3,4))
testreadcomplex(ones(ComplexF64,3,4,5))
finally
rm(fname,force=true)
end
end

@testset "write" begin
fname = tempname() * ".fits"

function _testwritecomplex(arr::AbstractArray{Complex{T}}) where T
FITS(fname,"w") do f
write(f, reinterpret(T,arr) )
end
a = FITS(fname,"r") do f
read(f[1])
end
b = collect(reinterpret(T,arr))
@test b == a
c = FITS(fname,"r") do f
reinterpret(Complex{T},read(f[1]))
end
@test c == arr
end

function testwritecomplex(arr::Array{Complex{T}}) where T
# Given a complex array, reinterpret it as float and write to file
# Read it back and compare
_testwritecomplex(arr)

# Write a subsection instead of the entire array
region = 1:size(arr,1)
arrv = @view arr[region]
_testwritecomplex(arrv)
end

try
testwritecomplex(ones(ComplexF64,3))
testwritecomplex(ones(ComplexF64,3,4))
testwritecomplex(ones(ComplexF64,3,4,5))
finally
rm(fname,force=true)
end
end
end
end

Expand Down

0 comments on commit 803b883

Please sign in to comment.