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/* diff --git a/.travis.yml b/.travis.yml index 3235644..fd5adb3 100644 --- a/.travis.yml +++ b/.travis.yml @@ -4,5 +4,8 @@ os: - osx julia: - 0.6 + - 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 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 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` + diff --git a/src/DICOM.jl b/src/DICOM.jl index ec00ebb..c9fb528 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,21 @@ const dcm_dict = dcm_init() const fieldname_dict = fieldname_init() _dcmdict_data_ = 0 +# These "empty" values are used internally. They are returned if a search fails. +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) + +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 +62,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] -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 + r = get(dcm_dict, gelt, emptyVR_lookup) + return(r[2]) end -function lookup(d, fieldnameString::String) - return(lookup(d, fieldname_dict[fieldnameString])) +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|| 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 +93,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 +137,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},Any}() 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] = data 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},Any}) + 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},Any},1}() while sz > 0 grp = read(st, UInt16) elt = read(st, UInt16) @@ -164,48 +175,69 @@ 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},Any},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) + # (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 + # (0x0028,0x0100) defines Bits Allocated + bitType = 16 + f = get(dcm, (0x0028,0x0100), nothing) + if f !== nothing + 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 + dtype = isSigned ? Int16 : UInt16 + end + yr=1 zr=1 - if vr=="OB" - xr = sz - dtype = UInt8 - else - 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), nothing) + if f !== nothing + xr = Int(f) + end + # (0028,0011) defines number of columns + f = get(dcm, (0x0028,0x0011), nothing) + if f !== nothing + yr = Int(f) + end + # (0028,0012) defines number of planes + f = get(dcm, (0x0028,0x0012), nothing) + if f !== nothing + 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 = 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 - 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) @@ -223,29 +255,28 @@ 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 -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 @@ -256,7 +287,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 @@ -273,22 +304,22 @@ 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))] -# 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, dVR=Dict{Tuple{UInt16,UInt16},String}()) 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 +329,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 +338,28 @@ function element(st, evr, dcm) # Assume private vr = "UN" end - if vr === false - error("dicom: unknown tag ", gelt) + if haskey(dVR, gelt) + vr = dVR[gelt] + end + if vr === emptyVR + if haskey(dVR, (0x0000,0x0000)) + vr = dVR[(0x0000,0x0000)] + elseif !haskey(dVR, gelt) + error("dicom: unknown tag ", gelt) + end 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)) : + vr=="ST" || vr=="LT" || vr=="UT" ? String(read(st, UInt8, sz)) : sz==0 || vr=="XX" ? Any[] : @@ -357,69 +402,96 @@ 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] + if length(data) == 1 + data = data[1] + end + end + + # Return vr by default + return(gelt, data, vr) end # todo: support maxlen -string_write(vals, maxlen) = join(vals, '\\') - -function element_write(st, evr, el::DcmElt) - gelt = el.tag - if el.vr != "" - vr = el.vr - else - vr = lookup_vr(el.tag) - if vr === false +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, '\\') + +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 el.tag == (0x7FE0, 0x0010) - return pixeldata_write(st, evr, el.data) + if gelt == (0x7FE0, 0x0010) + return pixeldata_write(st, evr, 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, data), vr) + end + + # Pack data into array container. This is to undo "data = data[1]" from element(). + if !isa(data,Array) && vr in ("FL","FD","SL","SS","UL","US") + data = [data] end + data = - isempty(el) ? UInt8[] : - vr in ("OB","OF","OW","ST","LT","UT") ? el[1] : + 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(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(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 + end + + dcm_store(st, gelt, s->write(s, data), vr) 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; header=true, maxGrp=0xffff, dVR=Dict{Tuple{UInt16,UInt16},String}()) st = open(fn) - evr = false - skip(st, 128) - 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 @@ -427,45 +499,72 @@ 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},Any}() + if giveVR + dcmVR = Dict{Tuple{UInt16,UInt16},String}() + end + while true - fld = element(st, evr, data) - if fld === false - return data + (gelt, data, vr) = element(st, evr, dcm, dVR) + if gelt === emptyTag || gelt[1] > maxGrp + break else - push!(data, fld) + dcm[gelt] = data + if giveVR + dcmVR[gelt] = vr + end 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 close(st) - return data + if giveVR + return(dcm,dcmVR) + else + return(dcm) + end end -function dcm_write(st, d) +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 - 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 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), "1.2.840.10008.1.2.1", "UI") + else + element_write(st, evr, (0x0002,0x0010), "1.2.840.10008.1.2", "UI") + end else - element_write(st, evr, DcmElt((0x0002,0x0010),[ "1.2.840.10008.1.2" ])) - end - for el in d - element_write(st, evr, el) + # Otherwise, use existing transfer UID, and overwrite evr accordingly + metaInfo = get(meta_uids, d[(0x0002,0x0010)], (false, true)) + evr = metaInfo[2] + end + # dVR is only used if it isn't empty and evr=true + if evr && !isempty(dVR) + for gelt in sort(collect(keys(d))) + # 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))) + element_write(st, evr, gelt, d[gelt]) + end end + close(st) end end 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" ], diff --git a/test/runtests.jl b/test/runtests.jl index 7be689d..f11d6a9 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -1,10 +1,150 @@ +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) +# 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") + +# 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-brain.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_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" + + @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") + @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") +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,vrMG); close(outIO) + dcm_write(outMG1b,dcmMG,vrMG) +end + +# Load saved DICOM data +dcmMR1 = dcm_parse(outMR1) +dcmCT1 = dcm_parse(outCT1) +(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, 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)] == "MR" + @test dcmCT1[(0x0008,0x0060)] == "CT" + @test dcmMG1[(0x0008,0x0060)] == "MG" + + @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)] == 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-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")`) + +# 2a. With user-supplied VRs +dVR_CTa = Dict( + (0x0008,0x0010) => "SH", + (0x0008,0x0040) => "US", + (0x0008,0x0041) => "LO", + (0x0018,0x1170) => "DS", + (0x0020,0x0030) => "DS", + (0x0020,0x0035) => "DS", + (0x0020,0x0050) => "DS", + (0x0020,0x0070) => "LO", + (0x0028,0x0005) => "US", + (0x0028,0x0040) => "CS", + (0x0028,0x0200) => "US") +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 + +fileMR_multiframe = joinpath(testdir, "MR-Explicit_Little_MultiFrame") +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 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