2024-11-01 21:12:33 +01:00
|
|
|
import SimpleITK as sitk
|
2024-11-02 21:14:27 +01:00
|
|
|
import numpy as np
|
2024-11-02 22:35:31 +01:00
|
|
|
import matplotlib.pyplot as plt
|
2024-10-30 17:50:07 +01:00
|
|
|
|
2024-11-01 21:12:33 +01:00
|
|
|
img = sitk.ReadImage("Thorax1.0B40f_Downsampled_by_2.nrrd")
|
2024-10-30 17:50:07 +01:00
|
|
|
|
|
|
|
xSp, ySp, zSp = img.GetSpacing()
|
2024-11-02 21:14:27 +01:00
|
|
|
volumePerVoxel = xSp * ySp * zSp
|
2024-10-30 17:50:07 +01:00
|
|
|
|
2024-10-30 20:13:38 +01:00
|
|
|
#print(f"Spacing: \n{xSp = }\n{ySp = }\n{zSp = }\n")
|
|
|
|
print(f"1) Die Voxel sind Quader mit den verschiedenen Seitenlängen: X: {xSp}, Y: {ySp}, Z: {zSp} \n")
|
2024-10-30 17:50:07 +01:00
|
|
|
|
|
|
|
sizeVector = img.GetSize()
|
|
|
|
xSi, ySi, zSi = sizeVector
|
|
|
|
|
|
|
|
print(f"Size: \n{xSi = }\n{ySi = }\n{zSi = }\n")
|
|
|
|
|
2024-10-30 20:13:38 +01:00
|
|
|
print(f"2) Das Bild hat in den drei Raumrichtungen {xSi * ySi * zSi} Voxel. \n")
|
|
|
|
print(f"3) Die Methode GetSize() gibt den Typ tuple zurück. \n")
|
2024-10-30 17:53:35 +01:00
|
|
|
|
2024-11-01 21:14:19 +01:00
|
|
|
#Sollte man kontrollieren
|
2024-11-02 21:14:27 +01:00
|
|
|
print(f"4) Der Bildausschnitt in die drei Raumrichtungen ist X: {xSp*xSi} mm, Y: {ySp*xSi} mm, Z: {zSp*xSi} mm\n")
|
2024-11-01 21:12:33 +01:00
|
|
|
|
2024-11-02 22:35:31 +01:00
|
|
|
# CT-Bild der Lunge erstellen
|
2024-11-01 21:12:33 +01:00
|
|
|
seedList = [(120,140,100)]
|
|
|
|
lungs = sitk.ConnectedThreshold(img, seedList , -1100, -200)
|
|
|
|
sitk.WriteImage(lungs, "lungs.nrrd")
|
|
|
|
|
2024-11-02 22:35:31 +01:00
|
|
|
# Volumen der Lunge berechnen
|
2024-11-02 21:14:27 +01:00
|
|
|
lungsArray = np.array(sitk.GetArrayFromImage(lungs))
|
|
|
|
voxelInLungs = lungsArray[lungsArray > 0].size * volumePerVoxel
|
|
|
|
lungsInLitre = voxelInLungs / 1000000
|
|
|
|
#print(lungsInLitre)
|
2024-11-01 15:11:23 +01:00
|
|
|
|
2024-11-02 22:35:31 +01:00
|
|
|
# Histogramm des CT-Bilds
|
|
|
|
thoraxArray = np.array(sitk.GetArrayFromImage(img))
|
|
|
|
thoraxArray1D = thoraxArray.flatten()
|
|
|
|
plt.hist(thoraxArray1D, 100)
|
|
|
|
plt.title("Histogramm eines Thorax")
|
|
|
|
plt.xlabel("Helligkeitswert in Hounsfield-Einheiten (HU)")
|
|
|
|
plt.ylabel("Häufigkeit in Millionen (mio)")
|
|
|
|
plt.show()
|
|
|
|
plt.clf()
|
|
|
|
|
|
|
|
|