贴张演示图
1、读取图片,我这里是dcm
dcm_path = ""
dcm = pydicom.read_file(dcm_path)
gt = dcm.pixel_array
2、对每一个slice进行归一化,加金属伪影,反归一化
withMA = np.zeros_like(gt)
for i in range(gt.shape[0]):
print(i)
ori = normalized(gt[i,:,:])
with_mask=Ma(ori) * (4000) - 1000
withMA[i] =with_mask
3、简单的归一化
def normalized(X):
maxX = 3000
minX = -1000
X=np.clip(X, minX, maxX)
X = (X - minX) / (maxX - minX)
return X
4、金属伪影,需要调整的:
vec =(-shape/2,shape/2),
x1,y1,为伪影位置,r1为伪影半径
a,c为伪影强度
def Ma(gt):
# 建立金属掩码
metalMask=np.zeros_like(gt)
vec=np.arange(-135, 135)
xx, yy=np.meshgrid(vec, vec)
x1=random.randint(-20, 20);
y1=-82 + np.abs(x1);
r1=4
metalMask[(xx - x1) ** 2 + (yy - y1) ** 2 <= r1 ** 2]=1
# 金属和非金属部分的Radon变换
radon_theta=np.arange(0, 180)
R_metal=radon(metalMask, theta=radon_theta, circle=True)
R_withMetal=radon(gt, theta=radon_theta, circle=True) + R_metal
# test_image=iradon(R_withMetal, theta=radon_theta, circle=True)
# 金属投影区域进行模拟变换
prj_metalSim=R_withMetal.copy()
a=6
c=0.02
# 假设 x 和 h 是根据某种方式计算得到的
x=np.linspace(a - 5, a + 5, 10) # 举例的 x 值范围
h=x - (2 * c * (x - a) + 1 - np.sqrt(np.maximum(4 * c * (x - a) + 1, 0))) / (2 * c) * np.sign(x - a)
# 参数 c 控制变换曲线的宽度,即金属伪影的过渡区域的宽度。
# 参数 a 是变换曲线中心的位置。通过调整这些参数的值,可以改变 h 函数的形状,进而影响金属伪影的外观。
# 创建插值函数,使用 fill_value='extrapolate' 来允许超出范围的插值
h_func=interp1d(x, h, kind='cubic', fill_value="extrapolate")
for i in range(prj_metalSim.shape[0]):
prj_metalSim[i, :]=h_func(prj_metalSim[i, :])
metal_image=iradon(prj_metalSim, theta=radon_theta, circle=True)
return metal_image
5、展示
plt.figure();plt.imshow(withMA[0], cmap='gray');plt.title('模拟金属伪影的CT图')
plt.show()