taigw / GeodisTK

geodesic distance transform of 2d/3d images
MIT License
130 stars 35 forks source link

Failed to calculate geodesic distance for three dimension image. #10

Open zhang-qiang-github opened 3 years ago

zhang-qiang-github commented 3 years ago

I have problem to calculate geodesic distance for three dimension image. Please have a look for the following code:

import GeodisTK
import SimpleITK as sitk
import os
import numpy as np
import matplotlib.pyplot as plt

path = '../data/arota'
img = sitk.GetArrayFromImage(sitk.ReadImage(os.path.join(path, 'input.mhd'))).astype(np.float32)
fore = sitk.GetArrayFromImage(sitk.ReadImage(os.path.join(path, 'foregroundMask.mhd'))) / 255
# for 2D
twoD_dis = GeodisTK.geodesic2d_raster_scan(img[:, 167, :], fore[:, 167, :].astype(np.uint8), 1, 2)

threeD_dis = GeodisTK.geodesic3d_raster_scan(img, fore.astype(np.uint8), 1, 1, 50)

plt.figure()
plt.subplot(1, 2, 1)
plt.imshow(twoD_dis)
plt.title('two dimension')
plt.subplot(1,2,2)
plt.imshow(threeD_dis[:, 167, :])
plt.title('three dimension')
plt.show()

The result is:

image

The left picture show the geodesic distance calculate from 2D image, and the right picture show the geodesic distance calculate from 3D image. The red arrow indicate the aorta, and it is where seed point locate at.

The left geodesic distance is what I expect, but the right distance is wrong.

Is there anything wrong?

taigw commented 3 years ago

Hi @zhang-qiang-github , please make sure you have read the demo for 3D images. https://github.com/taigw/GeodisTK/blob/master/demo3d.py If you can run that demo, please use similar hyper-parameters as that demo in your case.

zhang-qiang-github commented 3 years ago

@taigw, I am sorry that the wrong result is caused by wrong input. The input spacing should be [1, 1, 1], but I just input one value.

In addition, is it possible to speed the code up by multi-thread/OMP/cuda?

taigw commented 3 years ago

@zhang-qiang-github , Thanks for your suggestion. I think in a single forward pass, the value of each pixel depends on its left, up and left up neighbors, which makes it hard to use parallel computation. But computation for the left up and right lower regions of a single seed may be parallelly computed.