1

I use makie.jl with slicesNumb for visualization of PET/CT scans, I have 3d array of attenuation values and I display heatmap with changing slices using slider - this works well I have two problems

  1. I do not know how to be able to define custom colormaps (basically I need to be able to specify that all above some threshold value will be black and all below white and values between will have grey values proportional to attenuation value).

2)I would like to be able to display to display over my image (tachnically heatmap) another ones where I would be able to controll transparency - alpha value of pixels - in order to display some annotations/ PET ...

code that works but without those 2 functionalities and how it looks

using GLMakie

```@doc
simple display of single image - only in transverse plane
```
function singleCtScanDisplay(arr ::Array{Number, 3}) 
  

  fig = Figure()
  sl_x = Slider(fig[2, 1], range = 1:1:size(arr)[3], startvalue = 40)
  ax = Axis(fig[1, 1])
  hm = heatmap!(ax, lift(idx-> arr[:,:, floor(idx)], sl_x.value) ,colormap = :grays)
  Colorbar(fig[1, 2], hm)

  
  fig

end

enter image description here

Thanks for help !

Jakub Mitura
  • 159
  • 1
  • 14

2 Answers2

1

You can use Colors and ColorSchemeTools, but you will need to add the top and bottom of the scheme according to your thresholds.

using Colors, ColorSchemeTools

truemin = 0
truemax = 600
max_shown_black = 20
min_shown_white = 500

data = rand(truemin:truemax, (500, 500, 20))

grayscheme =  [fill(colorant"black", max_shown_black - truemin + 1);
               collect(make_colorscheme(identity, identity, identity,
                   length = min_shown_white - max_shown_black - 1));
               fill(colorant"white", truemax - min_shown_white + 1)]

 

For controlling alpha, I would add a popup window with an alpha slider. Take a look at some of the distributable DICOM tools for examples.

Bill
  • 5,600
  • 15
  • 27
0

I finally managed it basically I load 3 dimensional data stored in hdf5 (I loaded it into hdf5 from raw using python)

It enables viewing transverse slices and annotate 3d pathes in a mask that will be displayed over main image

exmpleH = @spawnat persistenceWorker Main.h5manag.getExample()
minimumm = -1000
maximumm = 2000
arrr= fetch(exmpleH)
imageDim = size(arrr)
using GLMakie
maskArr = Observable(BitArray(undef, imageDim))
MyImgeViewer.singleCtScanDisplay(arrr, maskArr,minimumm, maximumm)

Now definition of the required modules

```@doc
functions responsible for displaying medical image Data
```
using DrWatson
@quickactivate "Probabilistic medical segmentation"


module MyImgeViewer

using GLMakie
using Makie
#using GeometryBasics
using GeometricalPredicates
using ColorTypes
using Distributed
using GLMakie
using Main.imageViewerHelper
using Main.workerNumbers
## getting id of workers 

```@doc
simple display of single image - only in transverse plane we are adding also a mask  that 
arrr - main 3 dimensional data representing medical image for example in case of CT each voxel represents value of X ray attenuation
minimumm, maximumm - approximately minimum and maximum values we can have in our image
```
function singleCtScanDisplay(arrr ::Array{Number, 3}, maskArr , minimumm, maximumm) 
#we modify 2 pixels just in order to make the color range constant so slices will be displayed in the same windows
arrr[1,1,:].= minimumm 
arrr[2,1,:].= maximumm 
  

imageDim = size(arrr) # dimenstion of the primary image for example CT scan
slicesNumb =imageDim[3] # number of slices 

#defining layout variables
scene, layout = GLMakie.layoutscene(resolution = (600, 400))
ax1 = layout[1, 1] = GLMakie.Axis(scene, backgroundcolor = :transparent)
ax2 = layout[1, 1] = GLMakie.Axis(scene, backgroundcolor = :transparent)

#control widgets
sl_x =layout[2, 1]= GLMakie.Slider(scene, range = 1:1: slicesNumb , startvalue = slicesNumb/2 )
sliderXVal = sl_x.value


#color maps
cmwhite = cgrad(range(RGBA(10,10,10,0.01), stop=RGBA(0,0,255,0.4), length=10000));
greyss = createMedicalImageColorSchemeB(200,-200,maximumm, minimumm )
####heatmaps

#main heatmap that holds for example Ct scan
currentSliceMain = GLMakie.@lift(arrr[:,:, convert(Int32,$sliderXVal)])
hm = GLMakie.heatmap!(ax1, currentSliceMain ,colormap = greyss) 

#helper heatmap designed to respond to both changes in slider and changes in the bit matrix
currentSliceMask = GLMakie.@lift($maskArr[:,:, convert(Int32,$sliderXVal)])
hmB = GLMakie.heatmap!(ax1, currentSliceMask ,colormap = cmwhite) 

#adding ability to be able to add information to mask  where we clicked so in casse of mit matrix we will set the point where we clicked to 1 
indicatorC(ax1,imageDim,scene,maskArr,sliderXVal)

#displaying
colorB = layout[1,2]= Colorbar(scene, hm)
GLMakie.translate!(hmB, Vec3f0(0,0,5))  

scene

end

```@doc
inspired by   https://github.com/JuliaPlots/Makie.jl/issues/810
Generaly thanks to this function  the viewer is able to respond to clicking on the slices and records it in the supplied 3 dimensional AbstractArray
  ax - Axis which store our heatmap slices which we want to observe wheather user clicked on them and where
  dims - dimensions of  main image for example CT
  sc - Scene where our axis is
  maskArr - the 3 dimensional bit array  that has exactly the same dimensions as main Array storing image 
  sliceNumb - represents on what slide we are on currently on - ussually it just give information from slider 
```
function indicatorC(ax::Axis,dims::Tuple{Int64, Int64, Int64},sc::Scene,maskArr,sliceNumb::Observable{Any})
  register_interaction!(ax, :indicator) do event::GLMakie.MouseEvent, axis
  if event.type === MouseEventTypes.leftclick
    println("clicked")
    #@async begin
      #appropriately modyfing wanted pixels in mask array
   @async calculateMouseAndSetmaskWrap(maskArr, event,sc,dims,sliceNumb)            
    #  
    #  
    #  println("fetched" + fetch(maskA))

    #  finalize(maskA)
    #end
    return true
    #print("xMouse: $(xMouse)  yMouse: $(yMouse)   compBoxWidth: $(compBoxWidth)  compBoxHeight: $(compBoxHeight)   calculatedXpixel: $(calculatedXpixel)  calculatedYpixel: $(calculatedYpixel)      pixelsNumbInX  $(pixelsNumbInX)         ") 
  end
  
end
end
```@doc
wrapper for calculateMouseAndSetmask  - from imageViewerHelper module
  given mouse event modifies mask accordingly
  maskArr - the 3 dimensional bit array  that has exactly the same dimensions as main Array storing image 
  event - mouse event passed from Makie
  sc - scene we are using in Makie
  ```
function calculateMouseAndSetmaskWrap(maskArr, event,sc,dims,sliceNumb) 
  maskArr[] = calculateMouseAndSetmask(maskArr, event,sc,dims,sliceNumb)
end


end #module

and helper methods

```@doc
functions responsible for helping in image viewer - those functions are  meant to be invoked on separate process
- in parallel
```
using DrWatson
@quickactivate "Probabilistic medical segmentation"

module imageViewerHelper
using Documenter
using ColorTypes
using Colors, ColorSchemeTools
using Makie

export calculateMouseAndSetmask
export createMedicalImageColorSchemeB
# using AbstractPlotting
```@doc
  given mouse event modifies mask accordingly
  maskArr - the 3 dimensional bit array  that has exactly the same dimensions as main Array storing image 
  event - mouse event passed from Makie
  sc - scene we are using in Makie
  ```
function calculateMouseAndSetmask(maskArr, event,sc,dims,sliceNumb) 
  #position from top left corner 
  xMouse= Makie.to_world(sc,event.data)[1]
  yMouse= Makie.to_world(sc,event.data)[2]
  #data about height and width in layout
  compBoxWidth = 510 
  compBoxHeight = 510 
  #image dimensions - number of pixels  from medical image for example ct scan
  pixelsNumbInX =dims[1]
  pixelsNumbInY =dims[2]
  #calculating over which image pixel we are
  calculatedXpixel =convert(Int32, round( (xMouse/compBoxWidth)*pixelsNumbInX) )
  calculatedYpixel =  convert(Int32,round( (yMouse/compBoxHeight)*pixelsNumbInY ))
  sliceNumbConv =convert(Int32,round( sliceNumb[] ))
  #appropriately modyfing wanted pixels in mask array
  return markMaskArrayPatch( maskArr ,CartesianIndex(calculatedXpixel, calculatedYpixel, sliceNumbConv  ),2)
end



```@doc
  maskArr - the 3 dimensional bit array  that has exactly the same dimensions as main Array storing image 
  point - cartesian coordinates of point around which we want to modify the 3 dimensional array from 0 to 1

```
function markMaskArrayPatch(maskArr, pointCart::CartesianIndex{3}, patchSize ::Int64)

  ones = CartesianIndex(patchSize,patchSize,patchSize) # cartesian 3 dimensional index used for calculations to get range of the cartesian indicis to analyze
  maskArrB = maskArr[]
  for J in (pointCart-ones):(pointCart+ones)
    diff = J - pointCart # diffrence between dimensions relative to point of origin
      if cartesianTolinear(diff) <= patchSize
        maskArrB[J]=1
      end
      end
return maskArrB
end
```@doc
works only for 3d cartesian coordinates
  cart - cartesian coordinates of point where we will add the dimensions ...
```
function cartesianTolinear(pointCart::CartesianIndex{3}) :: Int16
   abs(pointCart[1])+ abs(pointCart[2])+abs(pointCart[3])
end

```@doc
creating grey scheme colors for proper display of medical image mainly CT scan
min_shown_white - max_shown_black range over which the gradint of greys will be shown
truemax - truemin the range of values in the image for which we are creating the scale
```
#taken from https://stackoverflow.com/questions/67727977/how-to-create-julia-color-scheme-for-displaying-ct-scan-makie-jl/67756158#67756158

function createMedicalImageColorSchemeB(min_shown_white,max_shown_black,truemax,truemin ) ::Vector{Any} 
# println("max_shown_black - truemin + 1")
# println(max_shown_black - truemin + 1)
# println(" min_shown_white - max_shown_black - 1")
# println( min_shown_white - max_shown_black - 1)
# println("truemax - min_shown_white + 1")
# println(truemax - min_shown_white + 1)

return  [fill(colorant"black", max_shown_black - truemin + 1);
               collect(make_colorscheme(identity, identity, identity,
                   length = min_shown_white - max_shown_black - 1));
               fill(colorant"white", truemax - min_shown_white + 1)]

end





end #module

enter image description here

Jakub Mitura
  • 159
  • 1
  • 14