emmt / EasyFITS.jl

Using FITS files made easier for Julia
Other
8 stars 2 forks source link

Reading a fits file already in memory #9

Open bensetterholm opened 3 months ago

bensetterholm commented 3 months ago

I occasionally find it useful to read in a FITS file into Julia as a byte array first, do some manipulation, and then have the resulting array parsed into a FITS object.

One example of such a work around is for reading in raw instrument reads from CHARA/MIRC-X, which puts a 4D data-cube into a FITS file and then compresses it with fpack. Owing to a bug in the underlying CFITSIO library, which does not natively uncompress data-cubes with >3 dimensions, my work-around with the FITSIO.jl library was to read in the .fits.fz file directly into a Vector{UInt8} array, manipulate the header in such a way so as to trick CFITSIO into uncompressing the data anyway, and then reformatting the output.

using FITSIO

function read_compressed(path::AbstractString)
    str = read(path, String)
    idx_compim = findfirst("EXTNAME = 'COMPRESSED_IMAGE'", str) |> last
    idx_znaxis = findnext("ZNAXIS  = ", str, idx_compim) |> last
    znaxis = parse(Int64, str[idx_znaxis+1:idx_znaxis+20])
    idx_zna = [last(findnext(Regex("ZNAXIS$(n)[ ]*= "), str, idx_znaxis)) for n in 1:znaxis]
    zna = [parse(Int64, str[i+1:i+20]) for i in idx_zna]
    bytes = unsafe_wrap(Vector{UInt8}, str)
    bytes[idx_znaxis+1:idx_znaxis+20] = Vector{UInt8}(lpad("1", 20))
    bytes[idx_zna[1]+1:idx_zna[1]+20] = Vector{UInt8}(lpad(string(prod(zna)), 20))
    return FITS(bytes) do f
        reshape(read(f["COMPRESSED_IMAGE"]), zna...)
    end
end

Another example is shown in this OIFITS.jl issue, where the end goal is to be able to read in OIFITS files that have a mixture of version 1 and version 2 tables with a little more finesse than what is available with the hack_revn keyword argument.

I would like to start using EasyFITS.jl as a replacement for FITSIO.jl. Since the former does not support reading of raw bytes, I tried to mimic the FITSIO.jl solution by adding the following internal constructor to the end of mutable struct FitsFile <: AbstractVector{FitsHDU} in types.jl

    function FitsFile(mem::Vector{UInt8}; name::AbstractString = "")
        handle = Ref{Ptr{CFITSIO.fitsfile}}(0)
        status = Ref{Status}(0)
        num = Ref{Cint}(0)

        CFITSIO.fits_open_memfile(
            handle,
            name,
            zero(Cint),
            Ptr{Ptr{Cvoid}}(pointer(mem)),
            Ptr{Csize_t}(length(mem)),
            2880,
            C_NULL,
            status
        ) |> check
        if !iszero(CFITSIO.fits_get_num_hdus(handle[], num, status))
            errcode = status[]
            status[] = 0
            CFITSIO.fits_close_file(ptr, status) # ptr is undefined?
            throw(FitsError(errcode))
        end
        return finalizer(close_handle, new(handle[], :r, name, zero(Cint)))
    end

However, every time I try to use this code (at least in Windows, I have not tried elsewhere), julia crashes with EXCEPTION_ACCESS_VIOLATION. Admittedly, all of this ccall stuff is a bit beyond my experience, so it is quite probable I am making a simple mistake. Based on the error message, I think it might be due to all of the pointer of pointer stuff going on in the code above.

Windows error message ``` Please submit a bug report with steps to reproduce this fault, and any error messages that follow (in their entirety). Thanks. Exception: EXCEPTION_ACCESS_VIOLATION at 0x64751533 -- mem_openmem at /workspace/srcdir/cfitsio-4.3.1\drvrmem.c:204 in expression starting at REPL[4]:1 mem_openmem at /workspace/srcdir/cfitsio-4.3.1\drvrmem.c:204 ffomem at /workspace/srcdir/cfitsio-4.3.1\cfileio.c:173 ffomem at c:\Users\setterholm\Documents\EasyFITS.jl\deps\deps.jl:159 _#2 at c:\Users\setterholm\Documents\EasyFITS.jl\src\types.jl:346 FitsFile at c:\Users\setterholm\Documents\EasyFITS.jl\src\types.jl:341 unknown function (ip: 000001c635cca0ab) jl_apply at C:/workdir/src\julia.h:2154 [inlined] do_call at C:/workdir/src\interpreter.c:126 eval_value at C:/workdir/src\interpreter.c:223 eval_stmt_value at C:/workdir/src\interpreter.c:174 [inlined] eval_body at C:/workdir/src\interpreter.c:677 jl_interpret_toplevel_thunk at C:/workdir/src\interpreter.c:817 jl_toplevel_eval_flex at C:/workdir/src\toplevel.c:943 jl_toplevel_eval_flex at C:/workdir/src\toplevel.c:886 ijl_toplevel_eval at C:/workdir/src\toplevel.c:952 [inlined] ijl_toplevel_eval_in at C:/workdir/src\toplevel.c:994 eval at .\boot.jl:428 [inlined] eval_user_input at C:\workdir\usr\share\julia\stdlib\v1.11\REPL\src\REPL.jl:222 repl_backend_loop at C:\workdir\usr\share\julia\stdlib\v1.11\REPL\src\REPL.jl:319 #start_repl_backend#59 at C:\workdir\usr\share\julia\stdlib\v1.11\REPL\src\REPL.jl:304 start_repl_backend at C:\workdir\usr\share\julia\stdlib\v1.11\REPL\src\REPL.jl:301 #run_repl#72 at C:\workdir\usr\share\julia\stdlib\v1.11\REPL\src\REPL.jl:460 run_repl at C:\workdir\usr\share\julia\stdlib\v1.11\REPL\src\REPL.jl:446 jfptr_run_repl_12691 at C:\Users\setterholm\.julia\juliaup\julia-1.11.0-beta1+0.x64.w64.mingw32\share\julia\compiled\v1.11\REPL\u0gqU_e6ieL.dll (unknown line) #1131 at .\client.jl:441 jfptr_YY.1131_17495 at C:\Users\setterholm\.julia\juliaup\julia-1.11.0-beta1+0.x64.w64.mingw32\share\julia\compiled\v1.11\REPL\u0gqU_e6ieL.dll (unknown line) jl_apply at C:/workdir/src\julia.h:2154 [inlined] jl_f__call_latest at C:/workdir/src\builtins.c:875 #invokelatest#2 at .\essentials.jl:1030 [inlined] invokelatest at .\essentials.jl:1027 [inlined] run_main_repl at .\client.jl:425 repl_main at .\client.jl:562 [inlined] _start at .\client.jl:536 jfptr__start_72187.1 at C:\Users\setterholm\.julia\juliaup\julia-1.11.0-beta1+0.x64.w64.mingw32\lib\julia\sys.dll (unknown line) jl_apply at C:/workdir/src\julia.h:2154 [inlined] true_main at C:/workdir/src\jlapi.c:900 jl_repl_entrypoint at C:/workdir/src\jlapi.c:1059 mainCRTStartup at C:/workdir/cli\loader_exe.c:58 BaseThreadInitThunk at C:\Windows\System32\KERNEL32.DLL (unknown line) RtlUserThreadStart at C:\Windows\SYSTEM32\ntdll.dll (unknown line) Allocations: 3905264 (Pool: 3905146; Big: 118); GC: 6 ```