From b47f80c93467521f7e47b2b34909b21f123ef828 Mon Sep 17 00:00:00 2001 From: Zaki A Date: Mon, 22 Jan 2018 12:31:29 -0500 Subject: [PATCH 01/16] Add de-identification fields (18,98-99) --- src/dcm_dict.jl | 2 ++ 1 file changed, 2 insertions(+) diff --git a/src/dcm_dict.jl b/src/dcm_dict.jl index 2fc6a72..a89c8d7 100644 --- a/src/dcm_dict.jl +++ b/src/dcm_dict.jl @@ -174,6 +174,8 @@ Any[ Any[ 18, 66 ], "Clinical Trial Subject Reading ID", "LO", "1" ], Any[ Any[ 18, 80 ], "Clinical Trial Time Point ID", "LO", "1" ], Any[ Any[ 18, 81 ], "Clinical Trial Time Point Description", "ST", "1" ], Any[ Any[ 18, 96 ], "Clinical Trial Coordinating Center Name", "LO", "1" ], +Any[ Any[ 18, 98 ], "Patient Identity Removed", "CS", "1" ], +Any[ Any[ 18, 99 ], "De-identification Method", "LO", "1-n" ], Any[ Any[ 24, 16 ], "Contrast/Bolus Agent", "LO", "1" ], Any[ Any[ 24, 18 ], "Contrast/Bolus Agent Sequence", "SQ", "1" ], Any[ Any[ 24, 20 ], "Contrast/Bolus Administration Route Sequence", "SQ", "1" ], From 3c106844f52bac4d4be161ba3a7a5b179ef92872 Mon Sep 17 00:00:00 2001 From: Zaki A Date: Mon, 22 Jan 2018 12:33:09 -0500 Subject: [PATCH 02/16] Switch from using Array{DcmElt} to Dict --- src/DICOM.jl | 333 +++++++++++++++++++++++++++++---------------------- 1 file changed, 187 insertions(+), 146 deletions(-) diff --git a/src/DICOM.jl b/src/DICOM.jl index ec00ebb..9423836 100644 --- a/src/DICOM.jl +++ b/src/DICOM.jl @@ -1,18 +1,29 @@ module DICOM -include("dcm_dict.jl") - export dcm_parse, dcm_write, lookup, lookup_vr +include("dcm_dict.jl") + +# Create dicom dictionary - used for reading/writing DICOM files +# Keys are tuple containing hex Group and Element of DICOM entry +# Stored value is String array: [Field Name, VR, Data Length], e.g: +# Julia> DICOM.dcm_dict[(0x0008,0x0005)] +# 3-element Array{String,1}: +# "Specific Character Set" +# "CS" +# "1-n" function dcm_init() - dcm_dict = Dict() + dcm_dict = Dict{Tuple{UInt16,UInt16},Array{String,1}}() for d in (_dcmdict_data_::Array{Any,1}) dcm_dict[(UInt16(d[1][1]),UInt16(d[1][2]))] = d[2:end] end - dcm_dict + return(dcm_dict) end -# Dict for looking up DICOM Tag (Tuple) from the fieldname (String) + +# For convenience, dictionary to get hex tag from field name, e.g: +# Julia> DICOM.fieldname_dict["Specific Character Set"] +# (0x0008, 0x0005) function fieldname_init() fieldname_dict = Dict{AbstractString, Tuple{UInt16,UInt16}}() for d in (_dcmdict_data_::Array{Any,1}) @@ -25,6 +36,30 @@ const dcm_dict = dcm_init() const fieldname_dict = fieldname_init() _dcmdict_data_ = 0 +# Composite type for storing DICOM entries/data in Julia +type DcmElt + data::Any + vr::String # "" except for non-standard VR + DcmElt(data,vr) = new(data,vr) + DcmElt(data) = new(data, "") +end + +# These "empty" values are used internally. They are returned if a search fails. +const emptyVR = "ZZ" # Can be any VR that doesn't exist +const emptyVR_dict = ["", emptyVR, ""] # Used in lookup_vr as failure state +const emptyTag = (0x0000,0x0000) +const emptyDcmElt = DcmElt([0], "") +const emptyDcmDict = Dict(DICOM.emptyTag => DICOM.emptyDcmElt) + +const VR_names = [ "AE","AS","AT","CS","DA","DS","DT","FL","FD","IS","LO","LT","OB","OF", + "OW","PN","SH","SL","SQ","SS","ST","TM","UI","UL","UN","US","UT" ] + +# mapping UID => bigendian? explicitvr? +const meta_uids = Dict([("1.2.840.10008.1.2", (false, false)), + ("1.2.840.10008.1.2.1", (false, true)), + ("1.2.840.10008.1.2.1.99", (false, true)), + ("1.2.840.10008.1.2.2", (true, true))]); + """ lookup_vr(tag::Tuple{Integer,Integer}) @@ -36,55 +71,30 @@ julia> lookup_vr((0x0028,0x0004)) "CS" ``` """ -function lookup_vr(gelt::Tuple{Integer,Integer}) +function lookup_vr(gelt::Tuple{UInt16,UInt16}) if gelt[1]&0xff00 == 0x5000 gelt = (0x5000,gelt[2]) elseif gelt[1]&0xff00 == 0x6000 gelt = (0x6000,gelt[2]) end - r = get(dcm_dict, gelt, false) - r !== false && r[2] + r = get(dcm_dict, gelt, emptyVR_dict) + return(r[2]) end -type DcmElt - tag::(Tuple{UInt16,UInt16}) - data::Array{Any,1} - vr::String # "" except for non-standard VR - DcmElt(tag, data) = new(tag,data,"") -end - -# takes dcm and tag (specify type?) -function lookup(d, t::Tuple) - for el in d - if isequal(el.tag,t) - return el - end - end - return false -end - -function lookup(d, fieldnameString::String) - return(lookup(d, fieldname_dict[fieldnameString])) +function lookup(d::Dict{Tuple{UInt16,UInt16},DICOM.DcmElt}, fieldnameString::String) + return(get(d, fieldname_dict[fieldnameString], emptyDcmElt)) end always_implicit(grp, elt) = (grp == 0xFFFE && (elt == 0xE0DD||elt == 0xE000|| elt == 0xE00D)) -VR_names = [ "AE","AS","AT","CS","DA","DS","DT","FL","FD","IS","LO","LT","OB","OF", - "OW","PN","SH","SL","SQ","SS","ST","TM","UI","UL","UN","US","UT" ] -# mapping UID => bigendian? explicitvr? -meta_uids = Dict([("1.2.840.10008.1.2", (false, false)), - ("1.2.840.10008.1.2.1", (false, true)), - ("1.2.840.10008.1.2.1.99", (false, true)), - ("1.2.840.10008.1.2.2", (true, true))]) - -dcm_store(st, grp, elt, writef) = dcm_store(st, grp, elt, writef, false) -function dcm_store(st, grp, elt, writef, vr) +dcm_store(st::IOStream, gelt::Tuple{UInt16,UInt16}, writef::Function) = dcm_store(st, gelt, writef, emptyVR) +function dcm_store(st::IOStream, gelt::Tuple{UInt16,UInt16}, writef::Function, vr::String) lentype = UInt32 - write(st, UInt16(grp)) - write(st, UInt16(elt)) - if vr !== false + write(st, UInt16(gelt[1])) # Grp + write(st, UInt16(gelt[2])) # Elt + if vr !== emptyVR write(st, vr) if vr in ("OB", "OW", "OF", "SQ", "UT", "UN") write(st, UInt16(0)) @@ -92,15 +102,25 @@ function dcm_store(st, grp, elt, writef, vr) lentype = UInt16 end end + # Write data first, then calculate length, then go back to write length p = position(st) - write(st, zero(lentype)) + write(st, zero(lentype)) # Placeholder for the data length writef(st) endp = position(st) - sz = endp-p-4 + # Remove placeholder's length, either 2 (UInt16) or 4 (UInt32) steps (1 step = 8 bits) + if lentype == UInt32 + sz = endp-p-4 + else + sz = endp-p-2 + end + szWasOdd = isodd(sz) # If length is odd, round up - UInt8(0) will be written at end + if szWasOdd + sz+=1 + end seek(st, p) - write(st, convert(lentype, sz)) + write(st, convert(lentype, max(0,sz))) seek(st, endp) - if isodd(sz) + if szWasOdd write(st, UInt8(0)) end end @@ -126,28 +146,28 @@ function undefined_length(st, vr) takebuf_array(data) end -function sequence_item(st, evr, sz) - item = Any[] +function sequence_item(st::IOStream, evr, sz) + item = Dict{Tuple{UInt16,UInt16},DcmElt}() endpos = position(st) + sz while position(st) < endpos - elt = element(st, evr) - if isequal(elt.tag, (0xFFFE,0xE00D)) + (gelt, data, vr) = element(st, evr) + if isequal(gelt, (0xFFFE,0xE00D)) break end - push!(item, elt) + item[gelt] = DcmElt(data,vr) end return item end -function sequence_item_write(st, evr, item) - for el in item - element_write(st, evr, el) +function sequence_item_write(st::IOStream, evr::Bool, items::Dict{Tuple{UInt16,UInt16},DICOM.DcmElt}) + for gelt in sort(collect(keys(items))) + element_write(st, evr, gelt, items[gelt]) end write(st, UInt16[0xFFFE, 0xE00D, 0x0000, 0x0000]) end function sequence_parse(st, evr, sz) - sq = Any[] + sq = Array{Dict{Tuple{UInt16,UInt16},DcmElt},1}() while sz > 0 grp = read(st, UInt16) elt = read(st, UInt16) @@ -164,16 +184,17 @@ function sequence_parse(st, evr, sz) return sq end -function sequence_write(st, evr, item) - for el in item - dcm_store(st, 0xFFFE, 0xE000, s->sequence_item_write(s, evr, el)) +function sequence_write(st::IOStream, evr::Bool, items::Array{Dict{Tuple{UInt16,UInt16},DICOM.DcmElt},1}) + for subitem in items + if length(subitem) > 0 + dcm_store(st, (0xFFFE,0xE000), s->sequence_item_write(s, evr, subitem)) + end end write(st, UInt16[0xFFFE, 0xE0DD, 0x0000, 0x0000]) end # always little-endian, "encapsulated" iff sz==0xffffffff -pixeldata_parse(st, sz, vr) = pixeldata_parse(st, sz, vr, false) -function pixeldata_parse(st, sz, vr, dcm) +function pixeldata_parse(st::IOStream, sz, vr::String, dcm=emptyDcmDict) yr=1 zr=1 if vr=="OB" @@ -183,22 +204,20 @@ function pixeldata_parse(st, sz, vr, dcm) xr = div(sz,2) dtype = UInt16 end - if dcm !== false - # (0028,0010) defines number of rows - f = lookup(dcm, (0x0028,0x0010)) - if f !== false - xr = f.data[1][1] - end - # (0028,0011) defines number of columns - f = lookup(dcm, (0x0028,0x0011)) - if f !== false - yr = f.data[1][1] - end - # (0028,0012) defines number of planes - f = lookup(dcm, (0x0028,0x0012)) - if f !== false - zr = f.data[1][1] - end + # (0028,0010) defines number of rows + f = get(dcm, (0x0028,0x0010), emptyDcmElt) + if f !== emptyDcmElt + xr = f.data[1] + end + # (0028,0011) defines number of columns + f = get(dcm, (0x0028,0x0011), emptyDcmElt) + if f !== emptyDcmElt + yr = f.data[1] + end + # (0028,0012) defines number of planes + f = get(dcm, (0x0028,0x0012), emptyDcmElt) + if f !== emptyDcmElt + zr = f.data[1] end if sz != 0xffffffff data = Array{dtype}(xr, yr, zr) @@ -223,22 +242,21 @@ function pixeldata_parse(st, sz, vr, dcm) return data end -function pixeldata_write(st, evr, el) - if length(el) > 1 - error("dicom: compression not supported") - end - d = el[1] +function pixeldata_write(st, evr, d) + # if length(el) > 1 + # error("dicom: compression not supported") + # end nt = eltype(d) vr = nt === UInt8 || nt === Int8 ? "OB" : nt === UInt16 || nt === Int16 ? "OW" : nt === Float32 ? "OF" : error("dicom: unsupported pixel format") if evr !== false - dcm_store(st, 0x7FE0, 0x0010, s->write(s,d), vr) + dcm_store(st, (0x7FE0,0x0010), s->write(s,d), vr) elseif vr != "OW" error("dicom: implicit VR only supports 16-bit pixels") else - dcm_store(st, 0x7FE0, 0x0010, s->write(s,d)) + dcm_store(st, (0x7FE0,0x0010), s->write(s,d)) end end @@ -275,20 +293,20 @@ end numeric_parse(st, T, sz) = [read(st, T) for i=1:div(sz,sizeof(T))] -# used internally to take a stream st that is reading a DICOM header -# returns DcmElt -element(st, evr) = element(st, evr, false) -function element(st, evr, dcm) +function element(st::IOStream, evr::Bool, dcm=emptyDcmDict) lentype = UInt32 diffvr = false local grp try grp = read(st, UInt16) catch - return false + return(emptyTag,0,emptyVR) end elt = read(st, UInt16) gelt = (grp,elt) + if grp <= 0x0002 + evr = true + end if evr && !always_implicit(grp,elt) vr = String(read(st, UInt8, 2)) if vr in ("OB", "OW", "OF", "SQ", "UT", "UN") @@ -298,7 +316,7 @@ function element(st, evr, dcm) end diffvr = !isequal(vr, lookup_vr(gelt)) else - vr = lookup_vr(gelt) + vr = elt == 0x0000 ? "UL" : lookup_vr(gelt) end if isodd(grp) && grp > 0x0008 && 0x0010 <= elt <+ 0x00FF # Private creator @@ -307,14 +325,14 @@ function element(st, evr, dcm) # Assume private vr = "UN" end - if vr === false + if vr === emptyVR error("dicom: unknown tag ", gelt) end sz = read(st,lentype) data = - vr=="ST" || vr=="LT" || vr=="UT" ? string(read(st, UInt8, sz)) : + vr=="ST" || vr=="LT" || vr=="UT" ? String(read(st, UInt8, sz)) : sz==0 || vr=="XX" ? Any[] : @@ -339,8 +357,8 @@ function element(st, evr, dcm) vr == "AS" ? String(read(st,UInt8,4)) : - vr == "DS" ? map(x->parse(Float64,x), string_parse(st, sz, 16, false)) : - vr == "IS" ? map(x->parse(Int,x), string_parse(st, sz, 12, false)) : + vr == "DS" ? [map(x->parse(Float64,x), string_parse(st, sz, 16, false))] : + vr == "IS" ? [map(x->parse(Int,x), string_parse(st, sz, 12, false))] : vr == "AE" ? string_parse(st, sz, 16, false) : vr == "CS" ? string_parse(st, sz, 16, false) : @@ -357,54 +375,71 @@ function element(st, evr, dcm) if isodd(sz) && sz != 0xffffffff skip(st, 1) end - delt = DcmElt(gelt, isa(data,Vector{Any}) ? data : Any[ data ]) - if diffvr - # record non-standard VR - delt.vr = vr - end - return delt + + # For convenience, get rid of array if it is just acting as a container + # Exception is "SQ", where array is part of structure + if length(data) == 1 && vr != "SQ" + data = data[1] + end + + # if diffvr + # # record non-standard VR + # return(gelt, data, vr) + # end + + # Saving vr by default + return(gelt, data, vr) end # todo: support maxlen -string_write(vals, maxlen) = join(vals, '\\') +string_write(vals::String, maxlen) = string_write([vals], maxlen) +string_write(vals::Array{String,1}, maxlen) = join(vals, '\\') -function element_write(st, evr, el::DcmElt) - gelt = el.tag +function element_write(st::IOStream, evr::Bool, gelt::Tuple{UInt16,UInt16}, el::DcmElt) if el.vr != "" vr = el.vr else - vr = lookup_vr(el.tag) - if vr === false + vr = lookup_vr(gelt) + if vr === emptyVR error("dicom: unknown tag ", gelt) end end - if el.tag == (0x7FE0, 0x0010) + if gelt == (0x7FE0, 0x0010) return pixeldata_write(st, evr, el.data) end - if evr !== false - evr = vr - end - el = el.data + if vr == "SQ" - return dcm_store(st, gelt[1], gelt[2], - s->sequence_write(s, evr, el), evr) + vr = evr ? vr : emptyVR + return dcm_store(st, gelt, + s->sequence_write(s, evr, el.data), vr) + end + + rawdata = el.data + # Pack data into array container. This is to undo "data = data[1]" from element(). + if length(rawdata) == 1 && vr in ("FL","FD","SL","SS","UL","US", "SQ") + rawdata =[rawdata] end + data = - isempty(el) ? UInt8[] : - vr in ("OB","OF","OW","ST","LT","UT") ? el[1] : + isempty(rawdata) ? UInt8[] : + vr in ("OB","OF","OW","ST","LT","UT") ? rawdata : vr in ("AE", "CS", "SH", "LO", "UI", "PN", "DA", "DT", "TM") ? - string_write(el, 0) : - vr == "FL" ? convert(Array{Float32,1}, el) : - vr == "FD" ? convert(Array{Float64,1}, el) : - vr == "SL" ? convert(Array{Int32,1}, el) : - vr == "SS" ? convert(Array{Int16,1}, el) : - vr == "UL" ? convert(Array{UInt32,1}, el) : - vr == "US" ? convert(Array{UInt16,1}, el) : - vr == "AT" ? [el...] : - vr in ("DS","IS") ? string_write(map(string,el), 0) : - el[1] - - dcm_store(st, gelt[1], gelt[2], s->write(s, data), evr) + string_write(rawdata, 0) : + vr == "FL" ? convert(Array{Float32,1}, rawdata) : + vr == "FD" ? convert(Array{Float64,1}, rawdata) : + vr == "SL" ? convert(Array{Int32,1}, rawdata) : + vr == "SS" ? convert(Array{Int16,1}, rawdata) : + vr == "UL" ? convert(Array{UInt32,1}, rawdata) : + vr == "US" ? convert(Array{UInt16,1}, rawdata) : + vr == "AT" ? [rawdata...] : + vr in ("DS","IS") ? string_write(map(string,rawdata), 0) : + rawdata + + if evr === false && gelt[1]>0x0002 + vr = emptyVR + end + + dcm_store(st, gelt, s->write(s, data), vr) end """ @@ -415,7 +450,9 @@ Reads file fn and returns an Array of DcmElt function dcm_parse(fn::AbstractString) st = open(fn) evr = false + # First 128 bytes are preamble - should be skipped skip(st, 128) + # "DICM" identifier must be after preamble sig = String(read(st,UInt8,4)) if sig != "DICM" error("dicom: invalid file header") @@ -427,24 +464,23 @@ function dcm_parse(fn::AbstractString) sig = String(read(st,UInt8,2)) evr = sig in VR_names skip(st, -6) - data = Any[] + dcm = Dict{Tuple{UInt16,UInt16},DcmElt}() while true - fld = element(st, evr, data) - if fld === false - return data + (gelt, data, vr) = element(st, evr, dcm) + if gelt === emptyTag + return dcm else - push!(data, fld) + dcm[gelt] = DcmElt(data,vr) end # look for transfer syntax UID - if fld.tag == (0x0002,0x0010) - fld = get(meta_uids, fld.data[1], false) - if fld !== false - evr = fld[2] - if fld[1] - # todo: set byte order to big - else - # todo: set byte order to little - end + if gelt == (0x0002,0x0010) + # Default is endian=little, explicitVR=true + metaInfo = get(meta_uids, data, (false, true)) + evr = metaInfo[2] + if metaInfo[1] + # todo: set byte order to big + else + # todo: set byte order to little end end end @@ -452,19 +488,24 @@ function dcm_parse(fn::AbstractString) return data end -function dcm_write(st, d) +function dcm_write(st::IOStream, d::Dict{Tuple{UInt16,UInt16},DICOM.DcmElt}) write(st, zeros(UInt8, 128)) write(st, "DICM") - # if any elements specify a VR then use explicit VR syntax - evr = anyp(x->x.vr!="", d) - # insert UID for our transfer syntax - if evr - element_write(st, evr, DcmElt((0x0002,0x0010),[ "1.2.840.10008.1.2.1" ])) + # If any elements specify a VR then use explicit VR syntax by default + evr = any(x->x.vr!="", values(d)) + # Insert UID for our transfer syntax, if it doesn't exist. Otherwise, use existing UID + if !haskey(d,(0x0002,0x0010)) + if evr + element_write(st, evr, (0x0002,0x0010), DcmElt("1.2.840.10008.1.2.1", "UI")) + else + element_write(st, evr, (0x0002,0x0010), DcmElt("1.2.840.10008.1.2", "UI")) + end else - element_write(st, evr, DcmElt((0x0002,0x0010),[ "1.2.840.10008.1.2" ])) + metaInfo = get(meta_uids, d[(0x0002,0x0010)].data, (false, true)) + evr = metaInfo[2] end - for el in d - element_write(st, evr, el) + for gelt in sort(collect(keys(d))) + element_write(st, evr, gelt, d[gelt]) end end From b171ffc4ff5c0b3ca95f4dadfad23822b56f8fb1 Mon Sep 17 00:00:00 2001 From: Zaki A Date: Sat, 27 Jan 2018 04:13:15 -0500 Subject: [PATCH 03/16] Ignore testdata folder --- .gitignore | 2 ++ 1 file changed, 2 insertions(+) create mode 100644 .gitignore diff --git a/.gitignore b/.gitignore new file mode 100644 index 0000000..6438ede --- /dev/null +++ b/.gitignore @@ -0,0 +1,2 @@ +# Ignore test data +test/testdata/* From b903a49648db18759224c026667eb0e7aa677fb2 Mon Sep 17 00:00:00 2001 From: Zaki A Date: Sat, 27 Jan 2018 06:15:13 -0500 Subject: [PATCH 04/16] Prevent skip_spaces() from exceeding length; Avoid single-element arrays. --- src/DICOM.jl | 16 ++++++++++------ 1 file changed, 10 insertions(+), 6 deletions(-) diff --git a/src/DICOM.jl b/src/DICOM.jl index 9423836..40947ba 100644 --- a/src/DICOM.jl +++ b/src/DICOM.jl @@ -260,10 +260,10 @@ function pixeldata_write(st, evr, d) end end -function skip_spaces(st) +function skip_spaces(st, endpos) while true c = read(st,Char) - if c != ' ' + if c != ' ' || position(st) == endpos return c end end @@ -274,7 +274,7 @@ function string_parse(st, sz, maxlen, spaces) data = [ "" ] first = true while position(st) < endpos - c = !first||spaces ? read(st,Char) : skip_spaces(st) + c = !first||spaces ? read(st,Char) : skip_spaces(st, endpos) if c == '\\' push!(data, "") first = true @@ -291,7 +291,7 @@ function string_parse(st, sz, maxlen, spaces) return data end -numeric_parse(st, T, sz) = [read(st, T) for i=1:div(sz,sizeof(T))] +numeric_parse(st::IOStream, T::DataType, sz) = T[read(st, T) for i=1:div(sz,sizeof(T))] function element(st::IOStream, evr::Bool, dcm=emptyDcmDict) lentype = UInt32 @@ -357,8 +357,8 @@ function element(st::IOStream, evr::Bool, dcm=emptyDcmDict) vr == "AS" ? String(read(st,UInt8,4)) : - vr == "DS" ? [map(x->parse(Float64,x), string_parse(st, sz, 16, false))] : - vr == "IS" ? [map(x->parse(Int,x), string_parse(st, sz, 12, false))] : + vr == "DS" ? map(x->parse(Float64,x), string_parse(st, sz, 16, false)): + vr == "IS" ? map(x->parse(Int,x), string_parse(st, sz, 12, false)): vr == "AE" ? string_parse(st, sz, 16, false) : vr == "CS" ? string_parse(st, sz, 16, false) : @@ -380,6 +380,9 @@ function element(st::IOStream, evr::Bool, dcm=emptyDcmDict) # Exception is "SQ", where array is part of structure if length(data) == 1 && vr != "SQ" data = data[1] + if length(data) == 1 + data = data[1] + end end # if diffvr @@ -392,6 +395,7 @@ function element(st::IOStream, evr::Bool, dcm=emptyDcmDict) end # todo: support maxlen +string_write(vals::Char, maxlen) = string_write(string(vals), maxlen) string_write(vals::String, maxlen) = string_write([vals], maxlen) string_write(vals::Array{String,1}, maxlen) = join(vals, '\\') From e6120a55418c3a9375759726f21388e8931cc4c4 Mon Sep 17 00:00:00 2001 From: Zaki A Date: Sat, 27 Jan 2018 08:05:15 -0500 Subject: [PATCH 05/16] Tests now include implicit and explicit VR types --- test/runtests.jl | 90 ++++++++++++++++++++++++++++++++++++++++++++---- 1 file changed, 84 insertions(+), 6 deletions(-) diff --git a/test/runtests.jl b/test/runtests.jl index 7be689d..2faf5f2 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -1,10 +1,88 @@ +using Base.Test using DICOM -dir = mktempdir() -filename = joinpath(dir,"dicomTestData.zip") +testdir = joinpath(Pkg.dir("DICOM"),"test","testdata") -download("http://www.dclunie.com/images/pixelspacingtestimages.zip", filename) -run(`unzip $filename -d $dir`) -dcm_parse(joinpath(dir, "DISCIMG/IMAGES/MGIMAGEA")) +if !isdir(testdir) + mkdir(testdir) +end -rm(filename) +fileMR = joinpath(testdir, "MR_Implicit_Little") +fileCT = joinpath(testdir, "CT_Explicit_Little") +fileMG = joinpath(testdir, "MG_Explicit_Little.zip") + +# Don't download files if they already exist +if !isfile(fileMR) && !isfile(fileMR) && !isfile(fileMR) + download("http://www.barre.nom.fr/medical/samples/files/MR-MONO2-16-head.gz", fileMR*".gz") + download("http://www.barre.nom.fr/medical/samples/files/CT-MONO2-16-ankle.gz", fileCT*".gz") + download("http://www.dclunie.com/images/pixelspacingtestimages.zip", fileMG) + + run(`gunzip -f $(fileMR*".gz")`) + run(`gunzip -f $(fileCT*".gz")`) + run(`unzip -o $fileMG -d $testdir`) +end + +# Load dicom data +dcmMR = dcm_parse(fileMR) +dcmCT = dcm_parse(fileCT) +dcmMG = dcm_parse(joinpath(testdir, "DISCIMG/IMAGES/MGIMAGEA")) + +@testset "Loading DICOM data" begin + @test dcmMR[(0x0008,0x0060)].data == "MR" + @test dcmCT[(0x0008,0x0060)].data == "CT" + @test dcmMG[(0x0008,0x0060)].data == "MG" + + @test length(dcmMR[(0x7FE0,0x0010)].data) == 65536 + @test length(dcmCT[(0x7FE0,0x0010)].data) == 262144 + @test length(dcmMG[(0x7FE0,0x0010)].data) == 262144 + + # Test lookup-by-fieldname + @test dcmMR[(0x0008,0x0060)] == lookup(dcmMR, "Modality") + @test dcmMR[(0x7FE0,0x0010)] == lookup(dcmMR, "Pixel Data") +end + +# Define two output files for each dcm - data will be saved, reloaded, then saved again +outMR1 = joinpath(testdir,"outMR1.dcm") +outMR2 = joinpath(testdir,"outMR2.dcm") + +outCT1 = joinpath(testdir,"outCT1.dcm") +outCT2 = joinpath(testdir,"outCT2.dcm") + +outMG1 = joinpath(testdir,"outMG1.dcm") +outMG2 = joinpath(testdir,"outMG2.dcm") + +@testset "Writing DICOM data" begin + # No specific test, just test if file is saved without errors + outIO = open(outMR1, "w+"); dcm_write(outIO,dcmMR); close(outIO) + outIO = open(outCT1, "w+"); dcm_write(outIO,dcmCT); close(outIO) + outIO = open(outMG1, "w+"); dcm_write(outIO,dcmMG); close(outIO) +end + +# Load saved DICOM data +dcmMR1 = dcm_parse(outMR1) +dcmCT1 = dcm_parse(outCT1) +dcmMG1 = dcm_parse(outMG1) + +@testset "Consistence of Reading/Writing data" begin + # Confirm that re-loading/saving leads to same file + outIO = open(outMR2, "w+"); dcm_write(outIO,dcmMR1); close(outIO) + outIO = open(outCT2, "w+"); dcm_write(outIO,dcmCT1); close(outIO) + outIO = open(outMG2, "w+"); dcm_write(outIO,dcmMG1); close(outIO) + + @test read(outMR1)==read(outMR2) + @test read(outCT1)==read(outCT2) + @test read(outMG1)==read(outMG2) + + # Repeat first tests on reloaded data + @test dcmMR1[(0x0008,0x0060)].data == "MR" + @test dcmCT1[(0x0008,0x0060)].data == "CT" + @test dcmMG1[(0x0008,0x0060)].data == "MG" + + @test length(dcmMR1[(0x7FE0,0x0010)].data) == 65536 + @test length(dcmCT1[(0x7FE0,0x0010)].data) == 262144 + @test length(dcmMG1[(0x7FE0,0x0010)].data) == 262144 + + # Test lookup-by-fieldname; cross-compare dcmMR with dcmMR1 + @test dcmMR1[(0x0008,0x0060)].data == lookup(dcmMR, "Modality").data + @test dcmMR1[(0x7FE0,0x0010)].data == lookup(dcmMR, "Pixel Data").data +end From 2b4e9ef113b3cd798937be72128193b3b4ca15ff Mon Sep 17 00:00:00 2001 From: Zaki A Date: Sat, 27 Jan 2018 09:26:56 -0500 Subject: [PATCH 06/16] Update travis to 0.6 and nightly --- .travis.yml | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/.travis.yml b/.travis.yml index 2df2416..2e02b46 100644 --- a/.travis.yml +++ b/.travis.yml @@ -3,6 +3,7 @@ os: - linux - osx julia: - - 0.5 + - 0.6 + - nightly notifications: email: false From 5168feef4e2865951d9a97ee03458431df55defd Mon Sep 17 00:00:00 2001 From: Zaki A Date: Sun, 28 Jan 2018 02:29:25 -0500 Subject: [PATCH 07/16] Uncouple data and VR in loaded Dict --- src/DICOM.jl | 158 +++++++++++++++++++++++++++-------------------- test/runtests.jl | 38 ++++++------ 2 files changed, 111 insertions(+), 85 deletions(-) diff --git a/src/DICOM.jl b/src/DICOM.jl index 40947ba..e6e2407 100644 --- a/src/DICOM.jl +++ b/src/DICOM.jl @@ -36,20 +36,11 @@ const dcm_dict = dcm_init() const fieldname_dict = fieldname_init() _dcmdict_data_ = 0 -# Composite type for storing DICOM entries/data in Julia -type DcmElt - data::Any - vr::String # "" except for non-standard VR - DcmElt(data,vr) = new(data,vr) - DcmElt(data) = new(data, "") -end - # These "empty" values are used internally. They are returned if a search fails. const emptyVR = "ZZ" # Can be any VR that doesn't exist -const emptyVR_dict = ["", emptyVR, ""] # Used in lookup_vr as failure state +const emptyVR_lookup = ["", emptyVR, ""] # Used in lookup_vr as failure state const emptyTag = (0x0000,0x0000) -const emptyDcmElt = DcmElt([0], "") -const emptyDcmDict = Dict(DICOM.emptyTag => DICOM.emptyDcmElt) +const emptyDcmDict = Dict(DICOM.emptyTag => nothing) const VR_names = [ "AE","AS","AT","CS","DA","DS","DT","FL","FD","IS","LO","LT","OB","OF", "OW","PN","SH","SL","SQ","SS","ST","TM","UI","UL","UN","US","UT" ] @@ -77,12 +68,12 @@ function lookup_vr(gelt::Tuple{UInt16,UInt16}) elseif gelt[1]&0xff00 == 0x6000 gelt = (0x6000,gelt[2]) end - r = get(dcm_dict, gelt, emptyVR_dict) + r = get(dcm_dict, gelt, emptyVR_lookup) return(r[2]) end -function lookup(d::Dict{Tuple{UInt16,UInt16},DICOM.DcmElt}, fieldnameString::String) - return(get(d, fieldname_dict[fieldnameString], emptyDcmElt)) +function lookup(d::Dict{Tuple{UInt16,UInt16},Any}, fieldnameString::String) + return(get(d, fieldname_dict[fieldnameString], nothing)) end always_implicit(grp, elt) = (grp == 0xFFFE && (elt == 0xE0DD||elt == 0xE000|| @@ -147,19 +138,19 @@ function undefined_length(st, vr) end function sequence_item(st::IOStream, evr, sz) - item = Dict{Tuple{UInt16,UInt16},DcmElt}() + item = Dict{Tuple{UInt16,UInt16},Any}() endpos = position(st) + sz while position(st) < endpos (gelt, data, vr) = element(st, evr) if isequal(gelt, (0xFFFE,0xE00D)) break end - item[gelt] = DcmElt(data,vr) + item[gelt] = data end return item end -function sequence_item_write(st::IOStream, evr::Bool, items::Dict{Tuple{UInt16,UInt16},DICOM.DcmElt}) +function sequence_item_write(st::IOStream, evr::Bool, items::Dict{Tuple{UInt16,UInt16},Any}) for gelt in sort(collect(keys(items))) element_write(st, evr, gelt, items[gelt]) end @@ -167,7 +158,7 @@ function sequence_item_write(st::IOStream, evr::Bool, items::Dict{Tuple{UInt16,U end function sequence_parse(st, evr, sz) - sq = Array{Dict{Tuple{UInt16,UInt16},DcmElt},1}() + sq = Array{Dict{Tuple{UInt16,UInt16},Any},1}() while sz > 0 grp = read(st, UInt16) elt = read(st, UInt16) @@ -184,7 +175,7 @@ function sequence_parse(st, evr, sz) return sq end -function sequence_write(st::IOStream, evr::Bool, items::Array{Dict{Tuple{UInt16,UInt16},DICOM.DcmElt},1}) +function sequence_write(st::IOStream, evr::Bool, items::Array{Dict{Tuple{UInt16,UInt16},Any},1}) for subitem in items if length(subitem) > 0 dcm_store(st, (0xFFFE,0xE000), s->sequence_item_write(s, evr, subitem)) @@ -205,22 +196,23 @@ function pixeldata_parse(st::IOStream, sz, vr::String, dcm=emptyDcmDict) dtype = UInt16 end # (0028,0010) defines number of rows - f = get(dcm, (0x0028,0x0010), emptyDcmElt) - if f !== emptyDcmElt - xr = f.data[1] + f = get(dcm, (0x0028,0x0010), nothing) + if f !== nothing + xr = f end # (0028,0011) defines number of columns - f = get(dcm, (0x0028,0x0011), emptyDcmElt) - if f !== emptyDcmElt - yr = f.data[1] + f = get(dcm, (0x0028,0x0011), nothing) + if f !== nothing + yr = f end # (0028,0012) defines number of planes - f = get(dcm, (0x0028,0x0012), emptyDcmElt) - if f !== emptyDcmElt - zr = f.data[1] + f = get(dcm, (0x0028,0x0012), nothing) + if f !== nothing + zr = f end if sz != 0xffffffff - data = Array{dtype}(xr, yr, zr) + data = + zr > 1 ? Array{dtype}(xr, yr, zr) : Array{dtype}(xr, yr) read!(st, data) else # start with Basic Offset Table Item @@ -399,45 +391,51 @@ string_write(vals::Char, maxlen) = string_write(string(vals), maxlen) string_write(vals::String, maxlen) = string_write([vals], maxlen) string_write(vals::Array{String,1}, maxlen) = join(vals, '\\') -function element_write(st::IOStream, evr::Bool, gelt::Tuple{UInt16,UInt16}, el::DcmElt) - if el.vr != "" - vr = el.vr - else - vr = lookup_vr(gelt) - if vr === emptyVR +element_write(st::IOStream, evr::Bool, gelt::Tuple{UInt16,UInt16}, data::Any) = element_write(st,evr,gelt,data,lookup_vr(gelt)) +function element_write(st::IOStream, evr::Bool, gelt::Tuple{UInt16,UInt16}, data::Any, vr::String) + if vr === emptyVR + # Element tags ending in 0x0000 are not included in dcm_dicm.jl, their vr is UL + if gelt[2] == 0x0000 + vr = "UL" + elseif isodd(gelt[1]) && gelt[1] > 0x0008 && 0x0010 <= gelt[2] <+ 0x00FF + # Private creator + vr = "LO" + elseif isodd(gelt[1]) && gelt[1] > 0x0008 + # Assume private + vr = "UN" + else error("dicom: unknown tag ", gelt) end end if gelt == (0x7FE0, 0x0010) - return pixeldata_write(st, evr, el.data) + return pixeldata_write(st, evr, data) end if vr == "SQ" vr = evr ? vr : emptyVR return dcm_store(st, gelt, - s->sequence_write(s, evr, el.data), vr) + s->sequence_write(s, evr, data), vr) end - rawdata = el.data # Pack data into array container. This is to undo "data = data[1]" from element(). - if length(rawdata) == 1 && vr in ("FL","FD","SL","SS","UL","US", "SQ") - rawdata =[rawdata] + if !isa(data,Array) && vr in ("FL","FD","SL","SS","UL","US") + data = [data] end data = - isempty(rawdata) ? UInt8[] : - vr in ("OB","OF","OW","ST","LT","UT") ? rawdata : + isempty(data) ? UInt8[] : + vr in ("OB","OF","OW","ST","LT","UT") ? data : vr in ("AE", "CS", "SH", "LO", "UI", "PN", "DA", "DT", "TM") ? - string_write(rawdata, 0) : - vr == "FL" ? convert(Array{Float32,1}, rawdata) : - vr == "FD" ? convert(Array{Float64,1}, rawdata) : - vr == "SL" ? convert(Array{Int32,1}, rawdata) : - vr == "SS" ? convert(Array{Int16,1}, rawdata) : - vr == "UL" ? convert(Array{UInt32,1}, rawdata) : - vr == "US" ? convert(Array{UInt16,1}, rawdata) : - vr == "AT" ? [rawdata...] : - vr in ("DS","IS") ? string_write(map(string,rawdata), 0) : - rawdata + string_write(data, 0) : + vr == "FL" ? convert(Array{Float32,1}, data) : + vr == "FD" ? convert(Array{Float64,1}, data) : + vr == "SL" ? convert(Array{Int32,1}, data) : + vr == "SS" ? convert(Array{Int16,1}, data) : + vr == "UL" ? convert(Array{UInt32,1}, data) : + vr == "US" ? convert(Array{UInt16,1}, data) : + vr == "AT" ? [data...] : + vr in ("DS","IS") ? string_write(map(string,data), 0) : + data if evr === false && gelt[1]>0x0002 vr = emptyVR @@ -449,9 +447,9 @@ end """ dcm_parse(fn::AbstractString) -Reads file fn and returns an Array of DcmElt +Reads file fn and returns a Dict """ -function dcm_parse(fn::AbstractString) +function dcm_parse(fn::AbstractString, giveVR=false) st = open(fn) evr = false # First 128 bytes are preamble - should be skipped @@ -468,13 +466,24 @@ function dcm_parse(fn::AbstractString) sig = String(read(st,UInt8,2)) evr = sig in VR_names skip(st, -6) - dcm = Dict{Tuple{UInt16,UInt16},DcmElt}() + dcm = Dict{Tuple{UInt16,UInt16},Any}() + if giveVR + dcmVR = Dict{Tuple{UInt16,UInt16},String}() + end + while true (gelt, data, vr) = element(st, evr, dcm) if gelt === emptyTag - return dcm + if giveVR + return(dcm,dcmVR) + else + return(dcm) + end else - dcm[gelt] = DcmElt(data,vr) + dcm[gelt] = data + if giveVR + dcmVR[gelt] = vr + end end # look for transfer syntax UID if gelt == (0x0002,0x0010) @@ -489,28 +498,43 @@ function dcm_parse(fn::AbstractString) end end close(st) - return data + if giveVR + return(dcm,dcmVR) + else + return(dcm) + end end -function dcm_write(st::IOStream, d::Dict{Tuple{UInt16,UInt16},DICOM.DcmElt}) +dcm_write(fn::String, d::Dict{Tuple{UInt16,UInt16},Any}, dVR=Dict{Tuple{UInt16,UInt16},String}()) = dcm_write(open(fn,"w+"),d,dVR) +function dcm_write(st::IOStream, d::Dict{Tuple{UInt16,UInt16},Any}, dVR=Dict{Tuple{UInt16,UInt16},String}()) write(st, zeros(UInt8, 128)) write(st, "DICM") - # If any elements specify a VR then use explicit VR syntax by default - evr = any(x->x.vr!="", values(d)) - # Insert UID for our transfer syntax, if it doesn't exist. Otherwise, use existing UID + # If no dictionary containing VRs is provided, then assume implicit VR - at first + evr = !isempty(dVR) if !haskey(d,(0x0002,0x0010)) + # Insert UID for our transfer syntax, if it doesn't exist if evr - element_write(st, evr, (0x0002,0x0010), DcmElt("1.2.840.10008.1.2.1", "UI")) + element_write(st, evr, (0x0002,0x0010), "1.2.840.10008.1.2.1", "UI") else - element_write(st, evr, (0x0002,0x0010), DcmElt("1.2.840.10008.1.2", "UI")) + element_write(st, evr, (0x0002,0x0010), "1.2.840.10008.1.2", "UI") end else - metaInfo = get(meta_uids, d[(0x0002,0x0010)].data, (false, true)) + # Otherwise, use existing transfer UID, and overwrite evr accordingly + metaInfo = get(meta_uids, d[(0x0002,0x0010)], (false, true)) evr = metaInfo[2] end - for gelt in sort(collect(keys(d))) - element_write(st, evr, gelt, d[gelt]) + # Possible cases: evr=true/false, dVR=empty/notEmpty. + # Only time dVR is used is if evr=true & dVR is not empty + if evr && !isempty(dVR) + for gelt in sort(collect(keys(d))) + element_write(st, evr, gelt, d[gelt], dVR[gelt]) + end + else + for gelt in sort(collect(keys(d))) + element_write(st, evr, gelt, d[gelt]) + end end + close(st) end end diff --git a/test/runtests.jl b/test/runtests.jl index 2faf5f2..5232105 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -25,16 +25,16 @@ end # Load dicom data dcmMR = dcm_parse(fileMR) dcmCT = dcm_parse(fileCT) -dcmMG = dcm_parse(joinpath(testdir, "DISCIMG/IMAGES/MGIMAGEA")) +(dcmMG, vrMG) = dcm_parse(joinpath(testdir, "DISCIMG/IMAGES/MGIMAGEA"), true) @testset "Loading DICOM data" begin - @test dcmMR[(0x0008,0x0060)].data == "MR" - @test dcmCT[(0x0008,0x0060)].data == "CT" - @test dcmMG[(0x0008,0x0060)].data == "MG" + @test dcmMR[(0x0008,0x0060)] == "MR" + @test dcmCT[(0x0008,0x0060)] == "CT" + @test dcmMG[(0x0008,0x0060)] == "MG" - @test length(dcmMR[(0x7FE0,0x0010)].data) == 65536 - @test length(dcmCT[(0x7FE0,0x0010)].data) == 262144 - @test length(dcmMG[(0x7FE0,0x0010)].data) == 262144 + @test length(dcmMR[(0x7FE0,0x0010)]) == 65536 + @test length(dcmCT[(0x7FE0,0x0010)]) == 262144 + @test length(dcmMG[(0x7FE0,0x0010)]) == 262144 # Test lookup-by-fieldname @test dcmMR[(0x0008,0x0060)] == lookup(dcmMR, "Modality") @@ -49,40 +49,42 @@ outCT1 = joinpath(testdir,"outCT1.dcm") outCT2 = joinpath(testdir,"outCT2.dcm") outMG1 = joinpath(testdir,"outMG1.dcm") +outMG1b = joinpath(testdir,"outMG1b.dcm") outMG2 = joinpath(testdir,"outMG2.dcm") @testset "Writing DICOM data" begin # No specific test, just test if file is saved without errors outIO = open(outMR1, "w+"); dcm_write(outIO,dcmMR); close(outIO) outIO = open(outCT1, "w+"); dcm_write(outIO,dcmCT); close(outIO) - outIO = open(outMG1, "w+"); dcm_write(outIO,dcmMG); close(outIO) + outIO = open(outMG1, "w+"); dcm_write(outIO,dcmMG,vrMG); close(outIO) + dcm_write(outMG1b,dcmMG,vrMG) end # Load saved DICOM data dcmMR1 = dcm_parse(outMR1) dcmCT1 = dcm_parse(outCT1) -dcmMG1 = dcm_parse(outMG1) +(dcmMG1, vrMG1) = dcm_parse(outMG1, true) @testset "Consistence of Reading/Writing data" begin # Confirm that re-loading/saving leads to same file outIO = open(outMR2, "w+"); dcm_write(outIO,dcmMR1); close(outIO) outIO = open(outCT2, "w+"); dcm_write(outIO,dcmCT1); close(outIO) - outIO = open(outMG2, "w+"); dcm_write(outIO,dcmMG1); close(outIO) + outIO = open(outMG2, "w+"); dcm_write(outIO,dcmMG1, vrMG1); close(outIO) @test read(outMR1)==read(outMR2) @test read(outCT1)==read(outCT2) @test read(outMG1)==read(outMG2) # Repeat first tests on reloaded data - @test dcmMR1[(0x0008,0x0060)].data == "MR" - @test dcmCT1[(0x0008,0x0060)].data == "CT" - @test dcmMG1[(0x0008,0x0060)].data == "MG" + @test dcmMR1[(0x0008,0x0060)] == "MR" + @test dcmCT1[(0x0008,0x0060)] == "CT" + @test dcmMG1[(0x0008,0x0060)] == "MG" - @test length(dcmMR1[(0x7FE0,0x0010)].data) == 65536 - @test length(dcmCT1[(0x7FE0,0x0010)].data) == 262144 - @test length(dcmMG1[(0x7FE0,0x0010)].data) == 262144 + @test length(dcmMR1[(0x7FE0,0x0010)]) == 65536 + @test length(dcmCT1[(0x7FE0,0x0010)]) == 262144 + @test length(dcmMG1[(0x7FE0,0x0010)]) == 262144 # Test lookup-by-fieldname; cross-compare dcmMR with dcmMR1 - @test dcmMR1[(0x0008,0x0060)].data == lookup(dcmMR, "Modality").data - @test dcmMR1[(0x7FE0,0x0010)].data == lookup(dcmMR, "Pixel Data").data + @test dcmMR1[(0x0008,0x0060)] == lookup(dcmMR, "Modality") + @test dcmMR1[(0x7FE0,0x0010)] == lookup(dcmMR, "Pixel Data") end From b26523d110a69137fc99c911db469fb18c308088 Mon Sep 17 00:00:00 2001 From: Zaki A Date: Sun, 28 Jan 2018 05:55:56 -0500 Subject: [PATCH 08/16] Allow partial, headless, and unsigned reading --- src/DICOM.jl | 47 ++++++++++++++++++++++++++--------------------- test/runtests.jl | 4 ++++ 2 files changed, 30 insertions(+), 21 deletions(-) diff --git a/src/DICOM.jl b/src/DICOM.jl index e6e2407..bdb0db2 100644 --- a/src/DICOM.jl +++ b/src/DICOM.jl @@ -37,7 +37,7 @@ const fieldname_dict = fieldname_init() _dcmdict_data_ = 0 # These "empty" values are used internally. They are returned if a search fails. -const emptyVR = "ZZ" # Can be any VR that doesn't exist +const emptyVR = "" # Can be any VR that doesn't exist const emptyVR_lookup = ["", emptyVR, ""] # Used in lookup_vr as failure state const emptyTag = (0x0000,0x0000) const emptyDcmDict = Dict(DICOM.emptyTag => nothing) @@ -186,14 +186,21 @@ end # always little-endian, "encapsulated" iff sz==0xffffffff function pixeldata_parse(st::IOStream, sz, vr::String, dcm=emptyDcmDict) + # (0x0028,0x0103) defined Pixel Representation + isSigned = false + f = get(dcm, (0x0028,0x0103), nothing) + if f !== nothing + # Data is signed if f==1 + isSigned = f == 1 + end yr=1 zr=1 if vr=="OB" xr = sz - dtype = UInt8 + dtype = isSigned ? Int8 : UInt8 else xr = div(sz,2) - dtype = UInt16 + dtype = isSigned ? Int16 : UInt16 end # (0028,0010) defines number of rows f = get(dcm, (0x0028,0x0010), nothing) @@ -449,16 +456,17 @@ end Reads file fn and returns a Dict """ -function dcm_parse(fn::AbstractString, giveVR=false) +function dcm_parse(fn::AbstractString, giveVR=false; header=true, maxGrp=0xffff) st = open(fn) - evr = false - # First 128 bytes are preamble - should be skipped - skip(st, 128) - # "DICM" identifier must be after preamble - sig = String(read(st,UInt8,4)) - if sig != "DICM" - error("dicom: invalid file header") - # seek(st, 0) + if header + # First 128 bytes are preamble - should be skipped + skip(st, 128) + # "DICM" identifier must be after preamble + sig = String(read(st,UInt8,4)) + if sig != "DICM" + error("dicom: invalid file header") + # seek(st, 0) + end end # a bit of a hack to detect explicit VR. seek past the first tag, # and check to see if a valid VR name is there @@ -473,12 +481,8 @@ function dcm_parse(fn::AbstractString, giveVR=false) while true (gelt, data, vr) = element(st, evr, dcm) - if gelt === emptyTag - if giveVR - return(dcm,dcmVR) - else - return(dcm) - end + if gelt === emptyTag || gelt[1] > maxGrp + break else dcm[gelt] = data if giveVR @@ -523,11 +527,12 @@ function dcm_write(st::IOStream, d::Dict{Tuple{UInt16,UInt16},Any}, dVR=Dict{Tup metaInfo = get(meta_uids, d[(0x0002,0x0010)], (false, true)) evr = metaInfo[2] end - # Possible cases: evr=true/false, dVR=empty/notEmpty. - # Only time dVR is used is if evr=true & dVR is not empty + # dVR is only used if it isn't empty and evr=true if evr && !isempty(dVR) for gelt in sort(collect(keys(d))) - element_write(st, evr, gelt, d[gelt], dVR[gelt]) + # dVR only needs to contain keys for cases where lookup_vr() fails + haskey(dVR, gelt) ? element_write(st, evr, gelt, d[gelt], dVR[gelt]) : + element_write(st, evr, gelt, d[gelt]) end else for gelt in sort(collect(keys(d))) diff --git a/test/runtests.jl b/test/runtests.jl index 5232105..4b978f6 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -23,11 +23,15 @@ if !isfile(fileMR) && !isfile(fileMR) && !isfile(fileMR) end # Load dicom data +dcmMR_partial = dcm_parse(fileMR, maxGrp=0x0008) dcmMR = dcm_parse(fileMR) dcmCT = dcm_parse(fileCT) (dcmMG, vrMG) = dcm_parse(joinpath(testdir, "DISCIMG/IMAGES/MGIMAGEA"), true) @testset "Loading DICOM data" begin + @test dcmMR_partial[(0x0008,0x0060)] == "MR" + @test haskey(dcmMR_partial, (0x7FE0,0x0010)) == false + @test dcmMR[(0x0008,0x0060)] == "MR" @test dcmCT[(0x0008,0x0060)] == "CT" @test dcmMG[(0x0008,0x0060)] == "MG" From a950a013359d0cf26d2f79d73a39324b9d089a72 Mon Sep 17 00:00:00 2001 From: Zaki A Date: Sun, 28 Jan 2018 22:43:59 -0500 Subject: [PATCH 09/16] Allow user-supplied VRs & Recognize pixel bits from dcm hdr --- src/DICOM.jl | 58 +++++++++++++++++++++++++++++++++++----------------- 1 file changed, 39 insertions(+), 19 deletions(-) diff --git a/src/DICOM.jl b/src/DICOM.jl index bdb0db2..b9d9c04 100644 --- a/src/DICOM.jl +++ b/src/DICOM.jl @@ -186,36 +186,51 @@ end # always little-endian, "encapsulated" iff sz==0xffffffff function pixeldata_parse(st::IOStream, sz, vr::String, dcm=emptyDcmDict) - # (0x0028,0x0103) defined Pixel Representation + # (0x0028,0x0103) defines Pixel Representation isSigned = false f = get(dcm, (0x0028,0x0103), nothing) if f !== nothing # Data is signed if f==1 isSigned = f == 1 end - yr=1 - zr=1 - if vr=="OB" - xr = sz + # (0x0028,0x0100) defines Pixel Representation + bitType = 16 + f = get(dcm, (0x0028,0x0100), nothing) + if f !== nothing + # Data is signed if f==1 + bitType = Int(f) + else + f = get(dcm, (0x0028,0x0101), nothing) + bitType = f !== nothing ? Int(f) : + vr == "OB" ? 8 : 16 + end + if bitType == 8 dtype = isSigned ? Int8 : UInt8 else - xr = div(sz,2) dtype = isSigned ? Int16 : UInt16 end + + yr=1 + zr=1 # (0028,0010) defines number of rows f = get(dcm, (0x0028,0x0010), nothing) if f !== nothing - xr = f + xr = Int(f) end # (0028,0011) defines number of columns f = get(dcm, (0x0028,0x0011), nothing) if f !== nothing - yr = f + yr = Int(f) end # (0028,0012) defines number of planes f = get(dcm, (0x0028,0x0012), nothing) if f !== nothing - zr = f + zr = Int(f) + end + # (0028,0008) defines number of frames + f = get(dcm, (0x0028,0x0008), nothing) + if f !== nothing + zr *= Int(f) end if sz != 0xffffffff data = @@ -223,7 +238,7 @@ function pixeldata_parse(st::IOStream, sz, vr::String, dcm=emptyDcmDict) read!(st, data) else # start with Basic Offset Table Item - data = Array{Any,1}(element(st, false)) + data = Array{Any,1}(element(st, false)[2]) while true grp = read(st, UInt16) elt = read(st, UInt16) @@ -292,7 +307,7 @@ end numeric_parse(st::IOStream, T::DataType, sz) = T[read(st, T) for i=1:div(sz,sizeof(T))] -function element(st::IOStream, evr::Bool, dcm=emptyDcmDict) +function element(st::IOStream, evr::Bool, dcm=emptyDcmDict, dVR=Dict{Tuple{UInt16,UInt16},String}()) lentype = UInt32 diffvr = false local grp @@ -324,12 +339,22 @@ function element(st::IOStream, evr::Bool, dcm=emptyDcmDict) # Assume private vr = "UN" end + if haskey(dVR, gelt) + vr = dVR[gelt] + end if vr === emptyVR error("dicom: unknown tag ", gelt) end sz = read(st,lentype) + # Empty VR can be supplied in dVR to skip an element + if vr == "" + sz = isodd(sz) ? sz+1 : sz + skip(st,sz) + return(element(st,evr,dcm,dVR)) + end + data = vr=="ST" || vr=="LT" || vr=="UT" ? String(read(st, UInt8, sz)) : @@ -383,13 +408,8 @@ function element(st::IOStream, evr::Bool, dcm=emptyDcmDict) data = data[1] end end - - # if diffvr - # # record non-standard VR - # return(gelt, data, vr) - # end - # Saving vr by default + # Return vr by default return(gelt, data, vr) end @@ -456,7 +476,7 @@ end Reads file fn and returns a Dict """ -function dcm_parse(fn::AbstractString, giveVR=false; header=true, maxGrp=0xffff) +function dcm_parse(fn::AbstractString, giveVR=false; header=true, maxGrp=0xffff, dVR=Dict{Tuple{UInt16,UInt16},String}()) st = open(fn) if header # First 128 bytes are preamble - should be skipped @@ -480,7 +500,7 @@ function dcm_parse(fn::AbstractString, giveVR=false; header=true, maxGrp=0xffff) end while true - (gelt, data, vr) = element(st, evr, dcm) + (gelt, data, vr) = element(st, evr, dcm, dVR) if gelt === emptyTag || gelt[1] > maxGrp break else From c643029c36c037603c28662280c0fc6846d5237b Mon Sep 17 00:00:00 2001 From: Zaki A Date: Sun, 28 Jan 2018 22:59:48 -0500 Subject: [PATCH 10/16] Add test for DICOM files without header or retired elements --- test/runtests.jl | 38 +++++++++++++++++++++++++++++++++++++- 1 file changed, 37 insertions(+), 1 deletion(-) diff --git a/test/runtests.jl b/test/runtests.jl index 4b978f6..d204a43 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -7,6 +7,8 @@ if !isdir(testdir) mkdir(testdir) end +# TEST SET 1: Simple Reading/Writing + fileMR = joinpath(testdir, "MR_Implicit_Little") fileCT = joinpath(testdir, "CT_Explicit_Little") fileMG = joinpath(testdir, "MG_Explicit_Little.zip") @@ -14,7 +16,7 @@ fileMG = joinpath(testdir, "MG_Explicit_Little.zip") # Don't download files if they already exist if !isfile(fileMR) && !isfile(fileMR) && !isfile(fileMR) download("http://www.barre.nom.fr/medical/samples/files/MR-MONO2-16-head.gz", fileMR*".gz") - download("http://www.barre.nom.fr/medical/samples/files/CT-MONO2-16-ankle.gz", fileCT*".gz") + download("http://www.barre.nom.fr/medical/samples/files/CT-MONO2-16-brain.gz", fileCT*".gz") download("http://www.dclunie.com/images/pixelspacingtestimages.zip", fileMG) run(`gunzip -f $(fileMR*".gz")`) @@ -92,3 +94,37 @@ dcmCT1 = dcm_parse(outCT1) @test dcmMR1[(0x0008,0x0060)] == lookup(dcmMR, "Modality") @test dcmMR1[(0x7FE0,0x0010)] == lookup(dcmMR, "Pixel Data") end + +# TEST SET 2: Reading uncommon datasets + +# 1. Loading DICOM file with missing header +fileOT = joinpath(testdir, "OT_Implicit_Little_Headless") +download("http://www.barre.nom.fr/medical/samples/files/OT-MONO2-8-a7.gz", fileOT*".gz") +run(`gunzip -f $(fileOT*".gz")`) + +dcmOT = dcm_parse(fileOT, header=false) + +# 2. Loading DICOM file with missing header and retired DICOM elements +fileCT = joinpath(testdir, "CT-Implicity_Little_Headless_Retired") +download("http://www.barre.nom.fr/medical/samples/files/CT-MONO2-12-lomb-an2.gz", fileCT*".gz") +run(`gunzip -f $(fileCT*".gz")`) + +dVR_CT = Dict( + (0x0008,0x0010) => "SH", + (0x0008,0x0040) => "US", + (0x0008,0x0041) => "LO", + (0x0018,0x1170) => "", + (0x0018,0x1190) => "", + (0x0020,0x0030) => "", + (0x0020,0x0035) => "DS", + (0x0020,0x0050) => "DS", + (0x0020,0x0070) => "LO", + (0x0028,0x0005) => "US", + (0x0028,0x0040) => "CS", + (0x0028,0x0200) => "US") +dcmCT = dcm_parse(fileCT, header=false, dVR=dVR_CT); + +@testset "Loading uncommon DICOM data" begin + @test dcmOT[(0x0008,0x0060)] == "OT" + @test dcmCT[(0x0008,0x0060)] == "CT" +end \ No newline at end of file From 9d56dcea52b1b765cef185c402d4ef115dff6043 Mon Sep 17 00:00:00 2001 From: Zaki A Date: Sun, 28 Jan 2018 23:02:54 -0500 Subject: [PATCH 11/16] Add code coverage --- .travis.yml | 2 ++ 1 file changed, 2 insertions(+) diff --git a/.travis.yml b/.travis.yml index 2e02b46..fd5adb3 100644 --- a/.travis.yml +++ b/.travis.yml @@ -7,3 +7,5 @@ julia: - nightly notifications: email: false +after_success: + - julia -e 'cd(Pkg.dir("DICOM")); Pkg.add("Coverage"); using Coverage; Codecov.submit(Codecov.process_folder())' \ No newline at end of file From 2956cfbcd30460ffb3b9aa00489cbdd965c95411 Mon Sep 17 00:00:00 2001 From: Zaki A Date: Sun, 28 Jan 2018 23:16:23 -0500 Subject: [PATCH 12/16] Add multi-frame test. Fix minor typo --- src/DICOM.jl | 3 +-- test/runtests.jl | 12 +++++++++++- 2 files changed, 12 insertions(+), 3 deletions(-) diff --git a/src/DICOM.jl b/src/DICOM.jl index b9d9c04..4b3fd5e 100644 --- a/src/DICOM.jl +++ b/src/DICOM.jl @@ -193,11 +193,10 @@ function pixeldata_parse(st::IOStream, sz, vr::String, dcm=emptyDcmDict) # Data is signed if f==1 isSigned = f == 1 end - # (0x0028,0x0100) defines Pixel Representation + # (0x0028,0x0100) defines Bits Allocated bitType = 16 f = get(dcm, (0x0028,0x0100), nothing) if f !== nothing - # Data is signed if f==1 bitType = Int(f) else f = get(dcm, (0x0028,0x0101), nothing) diff --git a/test/runtests.jl b/test/runtests.jl index d204a43..cc61a5d 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -105,7 +105,7 @@ run(`gunzip -f $(fileOT*".gz")`) dcmOT = dcm_parse(fileOT, header=false) # 2. Loading DICOM file with missing header and retired DICOM elements -fileCT = joinpath(testdir, "CT-Implicity_Little_Headless_Retired") +fileCT = joinpath(testdir, "CT-Implicit_Little_Headless_Retired") download("http://www.barre.nom.fr/medical/samples/files/CT-MONO2-12-lomb-an2.gz", fileCT*".gz") run(`gunzip -f $(fileCT*".gz")`) @@ -124,7 +124,17 @@ dVR_CT = Dict( (0x0028,0x0200) => "US") dcmCT = dcm_parse(fileCT, header=false, dVR=dVR_CT); +# 3. Loading DICOM file containing multiple frames + +fileMR_multiframe = joinpath(testdir, "MR-Explicit_Little_MultiFrame") +dlFile = "MR-heart.gz" +download("http://www.barre.nom.fr/medical/samples/files/MR-MONO2-8-16x-heart.gz", fileMR_multiframe*".gz") +run(`gunzip -f $(fileMR_multiframe*".gz")`) + +dcmMR_multiframe = dcm_parse(fileMR_multiframe) + @testset "Loading uncommon DICOM data" begin @test dcmOT[(0x0008,0x0060)] == "OT" @test dcmCT[(0x0008,0x0060)] == "CT" + @test dcmMR_multiframe[(0x0008,0x0060)] == "MR" end \ No newline at end of file From 29838ead8a8c2d1a5cfcda6f62fd762996c6f076 Mon Sep 17 00:00:00 2001 From: Zaki A Date: Mon, 29 Jan 2018 07:58:05 -0500 Subject: [PATCH 13/16] Allow master VR to handle all missing VR cases --- src/DICOM.jl | 6 +++++- test/runtests.jl | 23 +++++++++++++++++------ 2 files changed, 22 insertions(+), 7 deletions(-) diff --git a/src/DICOM.jl b/src/DICOM.jl index 4b3fd5e..337df22 100644 --- a/src/DICOM.jl +++ b/src/DICOM.jl @@ -342,7 +342,11 @@ function element(st::IOStream, evr::Bool, dcm=emptyDcmDict, dVR=Dict{Tuple{UInt1 vr = dVR[gelt] end if vr === emptyVR - error("dicom: unknown tag ", gelt) + if haskey(dVR, (0x0000,0x0000)) + vr = dVR[(0x0000,0x0000)] + elseif !haskey(dVR, gelt) + error("dicom: unknown tag ", gelt) + end end sz = read(st,lentype) diff --git a/test/runtests.jl b/test/runtests.jl index cc61a5d..8334f43 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -109,20 +109,26 @@ fileCT = joinpath(testdir, "CT-Implicit_Little_Headless_Retired") download("http://www.barre.nom.fr/medical/samples/files/CT-MONO2-12-lomb-an2.gz", fileCT*".gz") run(`gunzip -f $(fileCT*".gz")`) -dVR_CT = Dict( +# 2a. With user-supplied VRs +dVR_CTa = Dict( (0x0008,0x0010) => "SH", (0x0008,0x0040) => "US", (0x0008,0x0041) => "LO", - (0x0018,0x1170) => "", - (0x0018,0x1190) => "", - (0x0020,0x0030) => "", + (0x0018,0x1170) => "DS", + (0x0020,0x0030) => "DS", (0x0020,0x0035) => "DS", (0x0020,0x0050) => "DS", (0x0020,0x0070) => "LO", (0x0028,0x0005) => "US", (0x0028,0x0040) => "CS", (0x0028,0x0200) => "US") -dcmCT = dcm_parse(fileCT, header=false, dVR=dVR_CT); +dcmCTa = dcm_parse(fileCT, header=false, dVR=dVR_CTa); + +# 2b. With a master VR which skips elements +# Here we skip any element where lookup_vr() fails +# And we also force (0x0018,0x1170) to be read as float instead of integer +dVR_CTb = Dict( (0x0000,0x0000) => "", (0x0018,0x1170) => "DS") +dcmCTb = dcm_parse(fileCT, header=false, dVR=dVR_CTb); # 3. Loading DICOM file containing multiple frames @@ -135,6 +141,11 @@ dcmMR_multiframe = dcm_parse(fileMR_multiframe) @testset "Loading uncommon DICOM data" begin @test dcmOT[(0x0008,0x0060)] == "OT" - @test dcmCT[(0x0008,0x0060)] == "CT" + + @test dcmCTa[(0x0008,0x0060)] == "CT" + @test dcmCTb[(0x0008,0x0060)] == "CT" + @test haskey(dcmCTa, (0x0028,0x0040)) # dcmCTa should contain retired element + @test !haskey(dcmCTb, (0x0028,0x0040)) # dcmCTb skips retired elements + @test dcmMR_multiframe[(0x0008,0x0060)] == "MR" end \ No newline at end of file From e96ac557787f4c9183643cf412d94278544ccf30 Mon Sep 17 00:00:00 2001 From: Zaki A Date: Mon, 29 Jan 2018 07:59:16 -0500 Subject: [PATCH 14/16] Update license to 2018 --- LICENSE.txt | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/LICENSE.txt b/LICENSE.txt index fdc2848..e657c71 100644 --- a/LICENSE.txt +++ b/LICENSE.txt @@ -1,4 +1,4 @@ -Copyright (c) 2012-2013: Jeff Bezanson and contributors. +Copyright (c) 2012-2018: Jeff Bezanson and contributors. Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the From 36e255e0dd8d4b39d415b9010f84cc7ae119c4dc Mon Sep 17 00:00:00 2001 From: Zaki A Date: Mon, 29 Jan 2018 08:38:56 -0500 Subject: [PATCH 15/16] Add basic details to readme --- README.md | 74 ++++++++++++++++++++++++++++++++++++++++++++++++++++++- 1 file changed, 73 insertions(+), 1 deletion(-) diff --git a/README.md b/README.md index 56f1059..666fd75 100644 --- a/README.md +++ b/README.md @@ -1 +1,73 @@ -DICOM interface for the Julia language. +# DICOM.jl + +Julia interface for parsing/writing DICOM files + +## Usage + +**Installation** + +To install the package: +``` +julia> Pkg.add("DICOM") +``` + +Load the package by +``` +julia> using DICOM +``` + +**Reading Data** + +Read a DICOM file by +``` +julia> dcmData = dcm_parse("path/to/dicom/file") +``` +The data in `dcmData` is structured as a dictionary, and individual DICOM elements can be accessed by their hex tag. +For example, the hex tag of "Pixel Data" is `7FE0,0010`, and it can be accessed in Julia by `dcmData[(0x7FE0,0x0010)]`. + +**Writing Data** + +Data can be written to a DICOM file by +``` +julia> dcm_write("path/to/output/file", dcmData) +``` + +**Additional Notes** + +DICOM files use either explicit or implicit value representation (VR). For implicit files, `DICOM.jl` will use a lookup table to guess the VR from the DICOM element's hex tag. For explicit files, `DICOM.jl` will read the VRs from the file. + +- A user-defined dictionary can be supplied to override the default lookup table + For example, the "Instance Number" - tag `(0x0020,0x0013)` - is an integer (default VR = "IS"). We can read this as a float by setting the VR to "DS" by: + ``` + myVR = Dict( (0x0020,0x0013) => "DS" ) + dcmData = dcm_parse("path/to/dicom/file", dVR = myVR) + ``` + Now `dcmData[(0x0020,0x0013)]` will return a float instead of an integer. + +- It is possible to skip an element by setting its VR to `""`. + For example, we can skip reading the Instance Number by + ``` + myVR = Dict( (0x0020,0x0013) => "" ) + dcmData = dcm_parse("path/to/dicom/file", dVR = myVR) + ``` + and now `dcmData[(0x0020,0x0013)]` will return an error because the key `(0x0020,0x0013)` doesn't exist - it was skipped during reading. + +- The user-supplied VR can contain a master VR with the tag `(0x0000,0x0000)` which will be used whenever `DICOM.jl` is unable to guess the VR on its own. This is convenient for reading older dicom files and skipping retired elements - i.e. where the VR lookup fails - by: + ``` + myVR = Dict( (0x0000,0x0000) => "" ) + dcmData = dcm_parse("path/to/dicom/file", dVR = myVR) + ``` + +- A user-supplied VR can also be supplied during writing, e.g.: + ``` + # Note that dcm_write doesn't use a named input, unlike dcm_parse with "dVR =" + julia> dcm_write("path/to/output/file", dcmData, dcmVR) + ``` + where `dcmVR` is a dictionary which maps the hex tag to the VR. + +- A dictionary of VRs can be obtained by passing `true` as a 2nd argument to `dcm_parse()`, e.g.: + ``` + julia> (dcmData, dcmVR) = dcm_parse("path/to/dicom/file", true) + ``` + and `dcmVR` will contain a dictionary of VRs for all of the elements in `dcmData` + From ba26bd1807c2b063fc44b337f057f01072eee4a8 Mon Sep 17 00:00:00 2001 From: Zaki A Date: Tue, 30 Jan 2018 10:35:32 -0500 Subject: [PATCH 16/16] Handle substrings in string_write(). Fix typos. --- src/DICOM.jl | 6 ++++-- test/runtests.jl | 1 - 2 files changed, 4 insertions(+), 3 deletions(-) diff --git a/src/DICOM.jl b/src/DICOM.jl index 337df22..c9fb528 100644 --- a/src/DICOM.jl +++ b/src/DICOM.jl @@ -384,8 +384,8 @@ function element(st::IOStream, evr::Bool, dcm=emptyDcmDict, dVR=Dict{Tuple{UInt1 vr == "AS" ? String(read(st,UInt8,4)) : - vr == "DS" ? map(x->parse(Float64,x), string_parse(st, sz, 16, false)): - vr == "IS" ? map(x->parse(Int,x), string_parse(st, sz, 12, false)): + vr == "DS" ? map(x->parse(Float64,x), string_parse(st, sz, 16, false)) : + vr == "IS" ? map(x->parse(Int,x), string_parse(st, sz, 12, false)) : vr == "AE" ? string_parse(st, sz, 16, false) : vr == "CS" ? string_parse(st, sz, 16, false) : @@ -417,6 +417,8 @@ function element(st::IOStream, evr::Bool, dcm=emptyDcmDict, dVR=Dict{Tuple{UInt1 end # todo: support maxlen +string_write(vals::Array{SubString{String}}, maxlen) = string_write(convert(Array{String}, vals), maxlen) +string_write(vals::SubString{String}, maxlen) = string_write(convert(String, vals), maxlen) string_write(vals::Char, maxlen) = string_write(string(vals), maxlen) string_write(vals::String, maxlen) = string_write([vals], maxlen) string_write(vals::Array{String,1}, maxlen) = join(vals, '\\') diff --git a/test/runtests.jl b/test/runtests.jl index 8334f43..f11d6a9 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -133,7 +133,6 @@ dcmCTb = dcm_parse(fileCT, header=false, dVR=dVR_CTb); # 3. Loading DICOM file containing multiple frames fileMR_multiframe = joinpath(testdir, "MR-Explicit_Little_MultiFrame") -dlFile = "MR-heart.gz" download("http://www.barre.nom.fr/medical/samples/files/MR-MONO2-8-16x-heart.gz", fileMR_multiframe*".gz") run(`gunzip -f $(fileMR_multiframe*".gz")`)