如何创建 julia 配色方案以显示 Ct 扫描 Makie.jl

问题描述

我使用带有 slicesNumb 的 makie.jl 来可视化 PET/CT 扫描,我有衰减值的 3d 数组,我使用滑块显示热图,并使用滑块更改切片 - 这很好用,我有两个问题

  1. 我不知道如何能够定义自定义颜色图(基本上我需要能够指定所有高于某个阈值的将是黑色的,所有低于白色的都将具有与衰减值成正比的灰度值)。

2)我希望能够显示在我的图像上显示(技术上的热图)另一个我可以控制透明度 - 像素的 alpha 值 - 以便显示一些注释/PET ...

可以运行但没有这两个功能代码以及它的外观

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

感谢您的帮助!

解决方法

您可以使用 ColorsColorSchemeTools,但您需要根据阈值添加方案的顶部和底部。

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,length = min_shown_white - max_shown_black - 1));
               fill(colorant"white",truemax - min_shown_white + 1)]

 

为了控制 alpha,我会添加一个带有 alpha 滑块的弹出窗口。查看一些可分发的 DICOM 工具以获取示例。

,

我最终管理它基本上我加载了存储在 hdf5 中的 3 维数据(我使用 python 从原始加载到 hdf5)

它可以查看横向切片并在将显示在主图像上的蒙版中注释 3d 路径

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)

现在定义所需的模块

```@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},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,:].= 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,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,0.01),stop=RGBA(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[:,$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,sliderXVal)

#displaying
colorB = layout[1,2]= Colorbar(scene,hm)
GLMakie.translate!(hmB,Vec3f0(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,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,sliceNumb) 
  maskArr[] = calculateMouseAndSetmask(maskArr,sliceNumb)
end


end #module

和辅助方法

```@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,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",truemax - min_shown_white + 1)]

end





end #module

enter image description here