Rubyを使ってDICOM画像で遊んでみる

RubyでDICOM画像を扱う情報がほとんどないので少しずつメモしてみる。

ruby-dicomhttps://github.com/dicom/ruby-dicom

インストール

gem install dicom

1	gem install dicom

DICOM画像を読み込む

require 'dicom'
dcm = DICOM::DObject.read("test.dcm")

1 2	require 'dicom' dcm = DICOM::DObject.read("test.dcm")

メタデータを取り出す

dcm.modality　# {"Modality"=>"CT"}

1	dcm.modality　# {"Modality"=>"CT"}

画像データを取り出す

pixels = dcm.pixels
pixels_narray = dcm.narray # 旧NArray ：画像によっては正しく動かないこともあるので注意

1 2	pixels = dcm.pixels pixels_narray = dcm.narray # 旧NArray ：画像によっては正しく動かないこともあるので注意

IRuby Notebook上にDICOM画像を表示する

imageメソッドで、ImageMagickオブジェクトに変換できる。
内部ではrmagickもしくはmini_magickが動いている。
rmagickは更新が滞っているので、mini_magickがおすすめ。

image = dcm.image # Magick::Image
# normalizeするともっと見栄えがよくなる
image = image.normalize

image = dcm.image # Magick::Image

# normalizeするともっと見栄えがよくなる

image = image.normalize

※画像はDSB2017に使用されたsampleの1つです

ヒストグラムを表示する

Guido ZuidhofFull Preprocessing Tutorial
の序盤をRubyに書き換えてみた。とりあえず動いてくれればOKという発想なので、例によって間違いあるかも。

require 'numo/narray'
require 'numo/gnuplot'
require 'dicom'
require 'histogram/array'

# Some contents
INPUT_FOLDER = 'stage1'
patients = Dir.entries(INPUT_FOLDER)
patients = patients[2..-1] # remove './', '../'

require 'numo/narray'

require 'numo/gnuplot'

require 'dicom'

require 'histogram/array'

# Some contents

INPUT_FOLDER = 'stage1'

patients = Dir.entries(INPUT_FOLDER)

patients = patients[2..-1] # remove './', '../'

class DICOM::Element
  def values
    value.split("\\").map(&:to_f)
  end
end

def load_scan(path)
  file_paths = Dir.glob(File.join(INPUT_FOLDER, path, '/*'))
  slices = file_paths.map{ |fp| DICOM::DObject.read(fp)}
  slices.sort_by!{ |s| s.image_position_patient.values[2]}
  slice_thickness = (
    slices[0].image_position_patient.values[2] -
      slices[1].image_position_patient.values[2] ||
      slices[0].slice_location.values[2] -
      slices[1].slice_location.values[2]
    ).abs
  slices.each do |s|
    s.slice_thickness = slice_thickness
  end
end

class DICOM::Element

def values

value.split("\\").map(&:to_f)

end

def load_scan(path)

file_paths = Dir.glob(File.join(INPUT_FOLDER, path, '/*'))

slices = file_paths.map{ |fp| DICOM::DObject.read(fp)}

slices.sort_by!{ |s| s.image_position_patient.values[2]}

slice_thickness = (

slices[0].image_position_patient.values[2] -

slices[1].image_position_patient.values[2] ||

slices[0].slice_location.values[2] -

slices[1].slice_location.values[2]

).abs

slices.each do |s|

s.slice_thickness = slice_thickness

end

def get_pixels_hu(slices)
  image = Numo::Int16[*slices.map(&:pixels)]
  slices.each_with_index do |slice, i|
    intercept = slice.rescale_intercept.value.to_f
    slope = slice.rescale_slope.value.to_f
    if slope != 1.0
      image[i,true] = slope * image[i,true] # 自動でfloatにキャストされ、なおかつInt16に戻るみたい
      image[i,true] += intercept
    end
  end
  image
end

def get_pixels_hu(slices)

image = Numo::Int16[*slices.map(&:pixels)]

slices.each_with_index do |slice, i|

intercept = slice.rescale_intercept.value.to_f

slope = slice.rescale_slope.value.to_f

if slope != 1.0

image[i,true] = slope * image[i,true] # 自動でfloatにキャストされ、なおかつInt16に戻るみたい

image[i,true] += intercept

end

image

end

first_10_patients = Array.new(10) { |pt| load_scan(patients[pt]) }
first_10_patients_pixels = first_10_patients.map! { |pt| get_pixels_hu(pt) }

1 2	first_10_patients = Array.new(10) { \|pt\| load_scan(patients[pt]) } first_10_patients_pixels = first_10_patients.map! { \|pt\| get_pixels_hu(pt) }

# histogram 旧NArray対応のみ...なのでRuby配列に一旦変換
first_10_patients_pixels.each do |ptp|
  first_patient_pixels_array = ptp.flatten.to_a
  bins,freq = first_patient_pixels_array.histogram(50)
  bins.map!{|b| b.round.to_s}

  noteplot = Numo.noteplot do
    set title:"ヒストグラム"
    set style:[:fill, :solid, {border:-1}]
    set xtics:["rotate by": -45]
    plot bins, freq, u:"2:xtic(1)", w: :boxes, lc:{rgb:"orange"}, t:"分布"
  end
  
  IRuby.display noteplot
end

# histogram 旧NArray対応のみ...なのでRuby配列に一旦変換

first_10_patients_pixels.each do |ptp|

first_patient_pixels_array = ptp.flatten.to_a

bins,freq = first_patient_pixels_array.histogram(50)

bins.map!{|b| b.round.to_s}

noteplot = Numo.noteplot do

set title:"ヒストグラム"

set style:[:fill, :solid, {border:-1}]

set xtics:["rotate by": -45]

plot bins, freq, u:"2:xtic(1)", w: :boxes, lc:{rgb:"orange"}, t:"分布"

end

IRuby.display noteplot

end

IRuby Notebook上にDICOM画像を表示する

ヒストグラムを表示する

コメントを残す コメントをキャンセル

コメントを残すコメントをキャンセル