MATLAB Image Segmentation Algorithm Gives Empty Images in Julia?

matividalg · November 1, 2024, 9:39pm

Hi there! I’m using the DICOM, FileIO, ImageMagick and other packages to load DICOM images for segmentation, but I just seem unable to. I’ve checked every post and even the repo. This is my code:

using DICOM, Gtk4, ImageView, ImageMagick, Images, FileIO

filename = open_dialog("Choose a valid DICOM folder", select_folder = true)
data_array = dcmdir_parse(filename)
pixel_array = []

for i in range(1, size(data_array, 1))
    push!(pixel_array, (data_array[i])[(0x7FE0, 0x0010)][2])
end

for i in range(1, size(data_array, 1))
    filenamesave = raw"filedir" * string(i) *".png"
    FileIO.save(filenamesave, ImageMagick.load(IOBuffer(data_array[i][(0x7FE0, 0x0010)][2])))
end

Yes, this generates about 600 images in a folder, all black. Which is weird for me, as I have opened the dataset in MATLAB and it worked just fine. Does anyone have a clue on what could be happening? The JPEG2000 conversion is handled by ImageMagick.

ImageMagick.load(IOBuffer(data_array[i][(0x7FE0, 0x0010)][2]

The [2] is because the pixel values are a Vector(0, N). Basically, I’m flattening the vector so I can convert the image.

If anyone’s got a clue, I’d thank your help.
Thanks dearly,
Mati.

TheCedarPrince · November 1, 2024, 10:09pm

Hey @matividalg !

Oh no! Let’s see if we can get this sorted! Three quick questions:

Do you have an example or dummy DICOM image data set we could use that most closely matches the types of images you are working with?
What versions of these packages and Julia are you using? (i.e. the output of ] st in your Julia environment)
Could you share the MATLAB code here for reference?

I am going to CC a few folks who may be able to help more closely as I do not live and breathe DICOM images every day (but have worked with them): @Jakub_Mitura @cncastillo and @dilumaluthge

Additionally, just to let you know, you might want to check out the JuliaHealth organization. We are dedicated to bringing Julia to applications across the health research continuum. We live within the #health-and-medicine channel within the Julia Slack (if you haven’t joined, we would love to see you there: The Julia Language Slack). I wanted to mention this as there are a ton of clinical and health informaticist experts in the JuliaHealth community and you may be able to get more help also asking there.

Thanks!

matividalg · November 1, 2024, 10:32pm

Thanks for your quick and kind reply!
Sadly, I don’t have a dummy DICOM image for now, but I could provide the data array in a separate .txt file.

tatus `C:\Users\Matias Vidal\.julia\environments\v1.10\Project.toml`
  [a26e6606] DICOM v0.11.0
  [31c24e10] Distributions v0.25.112
  [7a1cc6ca] FFTW v1.8.0
  [5789e2e9] FileIO v1.16.4
⌅ [9db2cae5] Gtk4 v0.6.9
⌃ [6218d12a] ImageMagick v1.2.1
  [02fcd773] ImageTransformations v0.10.1
⌃ [86fae568] ImageView v0.12.1
  [916415d5] Images v0.26.1
⌃ [6a340f8b] KomaMRI v0.8.2
⌅ [d0bc0b20] KomaMRIBase v0.8.5
⌅ [4baa4f4d] KomaMRICore v0.8.3
⌅ [76db0263] KomaMRIPlots v0.8.3
  [23992714] MAT v0.10.7
  [efe261a4] NFFT v0.13.5
⌃ [f0f68f2c] PlotlyJS v0.18.13
  [91a5bcdd] Plots v1.40.8
⌃ [c3e4b0f8] Pluto v0.20.0
  [92933f4c] ProgressMeter v1.10.2

And the MATLAB code, I was dealing with an application, so as expected, used the Matlab GUI.

selpath = uigetdir;
app.CarpetaSeleccionadaEditField.Value = selpath;
[volumeData, volumeInfo] = dicomreadVolume(selpath);
imshow(volumeData(slice_selected,:,:), 'Parent', app.UIAxes);
figure(app.UIFigure);

Another thing I was thinking is that maybe my PC or Julia can’t access .JP2 data extensions… That’s quite weird, but when loading with ImageMagick, this appeared:

ERROR: MethodError: no method matching load(::Type{DataFormat{:JP2}}, ::IOBuffer)

But probably that’s only a sintax error.

I’m short on time right now so whenever I can I’ll provide a txt file with sample data.

TheCedarPrince · November 2, 2024, 12:27am

Hey @matividalg,

No worries! I actually googled around and found a test DICOM dataset and created a little toy repository to play around with this problem you encountered. Here’s a link to it:

It’s a publicly available DICOM mammogram dataset. I was able to read it effectively with DICOM as follows:

data = dcmdir_parse("dataset-uta4-dicom/dataset/mm/p1_d045f565-e0be-432c-aaab-8296b169c529")

fig = Figure(; size = (1200, 800));

ax = CairoMakie.Axis(fig[1, 1])
Makie.heatmap!(ax, data[1].PixelData; colormap = :magma)

Which gave the following figure:

Example Mammogram Loaded with DICOM

Because of this, I had a few observations:

First, the line in your example:

    push!(pixel_array, (data_array[i])[(0x7FE0, 0x0010)][2])

Only grabs one pixel from every single image. Is that what you wanted? Based on what you said here:

I don’t think so. You may want to have another look at that. If you wanted to get all the rows from the second column of each image, here is what I would do instead:

push!(pixel_array, (data[i])[(0x7FE0, 0x0010)][:, 2])

Another observation is that are you trying to downsample images (i.e. obtain an image with every other pixel dropped)? If you are, then that would look something like this (very naive but should get the idea across):

fig = Figure(; size = (1200, 800));

ax1 = CairoMakie.Axis(fig[1, 1])
ax2 = CairoMakie.Axis(fig[1, 2])
ds_img = []
img = (data[1])[(0x7FE0, 0x0010)]
for i in 1:size(img)[1]
    push!(ds_img, img[i, 1:2:end])
end

Makie.heatmap!(ax1, img; colormap = :magma)
Makie.heatmap!(ax2, reduce(hcat, ds_img)'; colormap = :magma)

Side by Side Downsampling

Finally, if you saved the images with your approach, you would get a vector of pixels from only one column of each image.

Hope some of those thoughts/observations were helpful – it seems like something might be going awry with your methods. If you are coming from MATLAB, be sure to give this resource from @cncastillo a look: Advanced Scientific Computing using Julia for the Uninitiated - HackMD

Let us know if you have any questions and would be happy to see you around the Slack!

P.S. Fun Friday evening hack! I’d be curious to hear what you are working on. We are always looking to showcase projects on the JuliaHealth Blog

P.P.S. Just a heads-up – I renamed your Discourse post title to perhaps reflect the problem here a little more precisely as I think this is not a DICOM-specific problem but something strange with the image segmentation algorithm you are trying to reimplement from MATLAB.

Jakub_Mitura · November 2, 2024, 8:39am

Hello @matividalg I do not work with JPEG - but It is not yet fully clear for me what do you want to achieve ? We are finilizing the updates of MedEye 1) that is used for visualization , and Medimages 2) for data loading into standardized struct, although the latter is sth like a month before it can be used , the first is already working and @Divital can guide you how to load data for visualization, if it is your goal; Although we are working generally on 3D data .

GitHub - JuliaHealth/MedEye3d.jl: Julia library for visualization and annotation medical images, specialized particularly for rapid development segmentation of 3 dimensional images like CT or PET/CT scans. Has full support of nuclear medicine Data.
GitHub - JuliaHealth/MedImages.jl

matividalg · November 2, 2024, 12:55pm

Hi, and yet again, thank you for your quick and kind responses. I have joined the Julia Slack and I’ll check the channel.

I apologize for not being clear. I’m loading a DICOM dataset that has been encoded and encapsulated with a JPEG format. I know this because of the field TransferSyntaxUID which is, 1.2.840.10008.1.2.4.70, which I thought it was JPEG2000 but it was JPEG Lossless, Nonhierarchical, First- Order Prediction. I wanted to extract the information in a 3d array, like I would do in MATLAB or Python, for implementing an algorithm of k-means segmentation for each image.

This is a project for University studies but I wanted to focus on a language that wasn’t MATLAB, and that was faster than Python.

I also noticed the Pixel Data of my dataset, adding to the fact that the size changes within every file, is in an encapsulated format, or so do I think. My own dataset gave the following variable:

julia> data[1].PixelData
2-element Vector{Any}:
 UInt8[]
 UInt8[0xff, 0xd8, 0xff, 0xe0, 0x00, 0x0f, 0x4c, 0x4a, 0x49, 0x46  …  0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x01, 0xff, 0xd9]

Whilst the dataset you kindly have provided is as following:

julia> data[1].PixelData
352×352 Matrix{UInt16}:
 0x0000  0x0000  0x0000  0x0000  0x0000  0x0000  0x0000  0x0000  …  0x0000  0x0000  0x0000  0x0000  0x0000  0x0000  0x0000
 0x0000  0x0000  0x0000  0x0000  0x0000  0x0000  0x0000  0x0000     0x0000  0x0000  0x0000  0x0000  0x0000  0x0000  0x0000
 0x0000  0x0000  0x0000  0x0000  0x0000  0x0000  0x0000  0x0000     0x0000  0x0000  0x0000  0x0000  0x0000  0x0000  0x0000   
 0x0000  0x0000  0x0000  0x0000  0x0000  0x0000  0x0000  0x0000     0x0000  0x0000  0x0000  0x0000  0x0000  0x0000  0x0000
      ⋮                                       ⋮                  ⋱       ⋮                                       ⋮  
 0x0000  0x0000  0x0000  0x0000  0x0000  0x0000  0x0000  0x0000     0x0000  0x0000  0x0000  0x0000  0x0000  0x0000  0x0000
 0x0000  0x0000  0x0000  0x0000  0x0000  0x0000  0x0000  0x0000     0x0000  0x0000  0x0000  0x0000  0x0000  0x0000  0x0000   
 0x0000  0x0000  0x0000  0x0000  0x0000  0x0000  0x0000  0x0000  …  0x0000  0x0000  0x0000  0x0000  0x0000  0x0000  0x0000
 0x0000  0x0000  0x0000  0x0000  0x0000  0x0000  0x0000  0x0000     0x0000  0x0000  0x0000  0x0000  0x0000  0x0000  0x0000

That’s why I was separating the pixel values and then transforming then with ImageMagick. My final goal here is obtain a matrix of that type that I know I can simply display and work with.

If any more information is needed, I’ll be happy of sharing it. I think I can manage to anonymize the dataset and provide you with just the DICOM file, but as it is I can’t share it, if that plays by the rules.

Also, thanks Jakub, although my current problem isn’t visualization but data manipulation. I’ll be checking the libraries you posted as they show promising for segmentation.

Thanks in advance!
Kindly, Matías.

matividalg · November 2, 2024, 12:59pm

I also failed to mention, sadly the code provided didn’t work with my files.
This is the error:

ERROR: ArgumentError:     Conversion failed for Heatmap (With conversion trait CellGrid()) with args: Tuple{Vector{Vector{UInt8}}} .
    Heatmap requires to convert to argument types Tuple{Union{AbstractVector{T} where T<:Real, MakieCore.EndPoints, AbstractMatrix{T} where T<:Real}, Union{AbstractVector{T} where T<:Real, MakieCore.EndPoints, AbstractMatrix{T} where T<:Real}, AbstractMatrix{<:Union{Float32, Float64, Colorant}}}, which convert_arguments didn't succeed in.
    To fix this overload convert_arguments(P, args...) for Heatmap or CellGrid() and return an object of type Tuple{Union{AbstractVector{T} where T<:Real, MakieCore.EndPoints, AbstractMatrix{T} where T<:Real}, Union{AbstractVector{T} where T<:Real, MakieCore.EndPoints, AbstractMatrix{T} where T<:Real}, AbstractMatrix{<:Union{Float32, Float64, Colorant}}}.`

TheCedarPrince · November 2, 2024, 4:36pm

Super rad that you are seeking out the challenge and opportunity here! I also appreciate that you can see the way of doing something like this in MATLAB is not exactly the same within Julia.

AH ok! I think I may still be a little unclear on what you are trying to do but if I can try understanding you are trying to:

Load DICOM images with Julia
Convert them using ImageMagick
Create a 600xMxN 3D matrix
Run a K-Means clustering algorithm across the 3D matrix

Is that what you are imagining? So at the end, you want a 3D matrix where each “slice” of the 3D matrix is an image loaded with ImageMagick – is that right? Maybe:

Could you describe more but what you mean by “simply display” and work with? My intuition feels like the MATLAB GUI you work with is probably hiding some complexity which is what is making the conversion process a bit tricky…

Oh I also meant to say if your PixelData looks like this:

julia> data[1].PixelData
2-element Vector{Any}:
 UInt8[]
 UInt8[0xff, 0xd8, 0xff, 0xe0, 0x00, 0x0f, 0x4c, 0x4a, 0x49, 0x46  …  0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x01, 0xff, 0xd9]

You would want to try this instead:

fig = Figure(; size = (1200, 800));

ax = CairoMakie.Axis(fig[1, 1])
Makie.heatmap!(ax, data[1].PixelData[2]; colormap = :magma)

fig # To show the figure

Yea, I was going to say the good news is is that it doesn’t seem like Julia or packages itself is the culprit here, but rather converting the algorithm or approach used in Python or MATLAB to Julia is tricky. If you have more code from Python or MATLAB that you can share, that would be great.