cascade MRI reconstruction: dataset.py--1

最新推荐文章于 2022-04-07 11:48:43 发布
最新推荐文章于 2022-04-07 11:48:43 发布
阅读量436 收藏
点赞数
分类专栏： MR图像重建
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
本文链接：https://blog.csdn.net/xuru_0927/article/details/118654663
版权

.mat转成传入cascade网络的形式


首先配置文件不能少config.yaml 

# Model Parameters
network:
        num_cascades: 6
        num_layers: 5 # Number of layers in the CNN per cascade
        num_filters: 64
        kernel_size: 3
        stride: 1
        padding: 1 #A padding of 1 is needed to keep the image in the same size
        noise: null #Noise in the measurements. To be used in the data consistency step

#Dataset parameters
dataset:
        data_path: 'data1/'
        acceleration_factor: 4.0
        fraction: 0.8 #train set size
        shuffle: 3 #Seed for numpy random generator
        sample_n: 10
        acq_noise: 0 #acquisation noise
        centred: False
        norm: 'ortho'  #norm: 'ortho' or null. if 'ortho', performs unitary transform, otherwise normal dft

# Training parameters
train:
        batch_size: 1
        num_epochs: 5
        early_stop: 100

        # Adam Optimizer Parameters
        learning_rate: 0.001
        b_1: 0.9
        b_2: 0.999
        l2: 0.0000001

        # Miscellaneous
        output_path: 'logs'
        cuda: False

单步执行

import os
import torch
import numpy as np
from math import ceil
from helpers_1 import *
from scipy.io import loadmat
from numpy.lib.stride_tricks import as_strided

import yaml
args = yaml.load(open('config.yaml', 'r'), Loader=yaml.FullLoader)

在这里插入图片描述
dataset = OCMRDataset(fold=‘train’, **args[‘dataset’])

        self.evalset = evalset
        self.data_path = data_path
        self.acc = acceleration_factor
        self.sample_n = sample_n
        self.noise = acq_noise
        self.centred = centred
        self.norm = norm
        self.files = os.listdir(self.data_path)
        if shuffle: # shuffle=None
            np.random.seed(shuffle)
            np.random.shuffle(self.files) 
        if fold == 'train':
            self.files = self.files[:int(len(self.files) * fraction)]

在这里插入图片描述

    def __getitem__(self, idx):
        if self.evalset and idx == 0: # evalset=False
            np.random.seed(9001)
        data = loadmat(os.path.join(self.data_path, self.files[idx]))['xn'] * 1e3

在这里插入图片描述
plt.imshow(np.abs(data))
在这里插入图片描述

import scipy.io as scio
import matplotlib.pyplot as plt
path_mat = 'data1/fs_0050_1_5T_0_0.mat'
data_mat = scio.loadmat(path_mat)
plt.imshow(np.abs(data_mat['xn']))
plt.show()

在这里插入图片描述

‘full’: data_gnd

        data = np.expand_dims(data, 0) #扩展维度，np.expand_dims(a,axis=)即在 a 的相应的axis轴上扩展维度

在这里插入图片描述

        data_gnd = format_data(data)

def format_data(data, mask=False):
    if mask: 
        data = data * (1+1j)
    data = complex2real(data)

def complex2real(x):
	x_real = np.real(x) #np.real取实部数值

在这里插入图片描述
plt.imshow(x_real[0,:,:],cmap=‘Greys’)
在这里插入图片描述

x_imag = np.imag(x) # np.imag取虚步数值

在这里插入图片描述
plt.imshow(x_imag[0,:,:],cmap=‘Greys’)
在这里插入图片描述

y = np.array([x_real, x_imag]).astype(np.float)

在这里插入图片描述

    if x.ndim >= 3: # x.ndim=3
        y = y.swapaxes(0, 1)
    return y

在这里插入图片描述

def format_data(data, mask=False):

    data = complex2real(data)
    return data.squeeze(0) #将数组中维度0移除

在这里插入图片描述
data_gnd = data.squeeze(0) #2代表频域数据的实部和虚部, 图片的宽高（256 x 256)
在这里插入图片描述

‘mask’: mask.transpose(1,2,0)

        mask = self.cartesian_mask(data.shape)

    def cartesian_mask(self, shape):
        N, Nx, Ny = int(np.prod(shape[:-2])), shape[-2], shape[-1]
        pdf_x = normal_pdf(Nx, 0.5/(Nx/10.)**2)

shape[:-2]为(1, 256)，shape[:-2]为(1,)
np.prod为sum操作
在这里插入图片描述
pdf_x:
在这里插入图片描述

def normal_pdf(length, sensitivity):
    return np.exp(-sensitivity * (np.arange(length) - length / 2)**2)

在这里插入图片描述
在这里插入图片描述

        lmda = Nx/(2.*self.acc) # self.acc = acceleration_factor=4.0, lmda=32
        n_lines = int(Nx / self.acc) # 64

        # add uniform distribution
        pdf_x += lmda * 1./Nx

pdf_x:
在这里插入图片描述

        if self.sample_n: # sample_n=10
            pdf_x[Nx//2-self.sample_n//2:Nx//2+self.sample_n//2] = 0 #pdf_x[123:133]

在这里插入图片描述

            pdf_x /= np.sum(pdf_x)
            n_lines -= self.sample_n #54

pdf_x:
在这里插入图片描述

在这里插入图片描述
在这里插入图片描述

        mask = np.zeros((N, Nx))
        for i in range(N):
#numpy.random.choice(a, size=None, replace=True, p=None)
#从a(只要是ndarray都可以，但必须是一维的)中随机抽取数字，并组成指定大小(size)的数组
#replace:True表示可以取相同数字，False表示不可以取相同数字
#数组p：与数组a相对应，表示取数组a中每个元素的概率，默认为选取每个元素的概率相同。
            idx = np.random.choice(Nx, n_lines, False, pdf_x)
            mask[i, idx] = 1

在这里插入图片描述
在这里插入图片描述
在这里插入图片描述
在这里插入图片描述
在这里插入图片描述

        if self.sample_n:
            mask[:, Nx//2-self.sample_n//2:Nx//2+self.sample_n//2] = 1 # [123:133] shape(1, 256)

未执行=1时：
在这里插入图片描述
执行=1后：
在这里插入图片描述

        size = mask.itemsize # 8 float64
        mask = as_strided(mask, (N, Nx, Ny), (size * Nx, size, 0)) # (1, 256, 256)

在这里插入图片描述
plt.imshow(mask[0,:,:],cmap=‘Greys’)
在这里插入图片描述

mask = mask.reshape(shape) #(1, 256, 256)

reshape之后可视化条形图没变化
在这里插入图片描述

        if not self.centred:
            mask = ifftshift(mask, axes=(-1, -2))

        return mask

在这里插入图片描述
plt.imshow(mask[0,:,:],cmap=‘Greys’)
在这里插入图片描述
在这里插入图片描述

#这句要放到最后，否则后面代码assert x.shape == mask.shape会报错
        mask = format_data(mask, mask=True) 

def format_data(data, mask=False):
    if mask: 
        data = data * (1+1j)

    data = complex2real(data)
    return data.squeeze(0)

def complex2real(x):
    x_real = np.real(x)
    x_imag = np.imag(x)
    y = np.array([x_real, x_imag]).astype(np.float)
    # re-order in convenient order
    if x.ndim >= 3:
        y = y.swapaxes(0, 1)
    return y

在这里插入图片描述
在这里插入图片描述
在这里插入图片描述
data = complex2real(data)即等于y
在这里插入图片描述
mask = data.squeeze(0)
在这里插入图片描述
图是一样的, 因为mask * (1+1j)，这样实部和虚部数值一样
在这里插入图片描述

‘image’&‘k’:首先准备data_und, k_und

        data_und, k_und = self.undersample(data, mask)

在这里插入图片描述

        assert x.shape == mask.shape
        # zero mean complex Gaussian noise
        noise_power = self.noise # 0
        nz = np.sqrt(.5)*(np.random.normal(0, 1, x.shape) + 1j * np.random.normal(0, 1, x.shape)) # nz.shape(1, 256, 256),x.shape(1, 256, 256)

在这里插入图片描述

 nz = nz * np.sqrt(noise_power)

nz:
在这里插入图片描述
nz的可视化图是空白的，因为noise_power=0，所以min、max都是0j
在这里插入图片描述

        if self.norm == 'ortho':
            # multiplicative factor
            nz = nz * np.sqrt(np.prod(mask.shape[-2:])) #mask.shape[-2:]:(256, 256), np.prod():65536, np.sqrt():256.0

nz可视化图是空白的，因为上一个nz=0，相乘=0
在这里插入图片描述

        if self.centred: # centred=False
            x_f = fft2c(x, norm=self.norm)
            x_fu = mask * (x_f + nz)
            x_u = ifft2c(x_fu, norm=self.norm)
            return x_u, x_fu
        else:
            x_f = fft2(x, norm=self.norm) # norm='ortho'
            x_fu = mask * (x_f + nz)
            x_u = ifft2(x_fu, norm=self.norm)
            return x_u, x_fu

plt.imshow(np.log(np.abs(x_f[0,:,:])+1e-7),cmap=‘Greys’)
在这里插入图片描述
plt.imshow(np.log(np.abs(x_fu[0,:,:])+1e-7),cmap=‘Greys’)
在这里插入图片描述

在这里插入图片描述
在这里插入图片描述
data_und, k_und = x_u, x_fu
在这里插入图片描述
在这里插入图片描述
plt.imshow(np.abs(x_u[0,:,:]),cmap=‘Greys’)
在这里插入图片描述
plt.imshow(np.abs(data[0,:,:]),cmap=‘Greys’)
在这里插入图片描述
上面的流程总结如下：
在这里插入图片描述

‘image’: data_und

        data_und = format_data(data_und)

def format_data(data, mask=False):

    data = complex2real(data)
    return data.squeeze(0)

def complex2real(x):
    x_real = np.real(x)
    x_imag = np.imag(x)
    y = np.array([x_real, x_imag]).astype(np.float)
    # re-order in convenient order
    if x.ndim >= 3:
        y = y.swapaxes(0, 1)  #由(2,1,256,256)变(1,2,256,256)
    return y

data = complex2real(data)即等于y
在这里插入图片描述
在这里插入图片描述
plt.imshow(x_real[0,:,:],cmap=‘Greys’)
在这里插入图片描述
在这里插入图片描述
plt.imshow(x_imag[0,:,:],cmap=‘Greys’)
在这里插入图片描述

在这里插入图片描述
data_und = data.squeeze(0)
在这里插入图片描述

‘k’: k_und.transpose(1,2,0)

        k_und = format_data(k_und)

def format_data(data, mask=False):

    data = complex2real(data)
    return data.squeeze(0)

def complex2real(x):
    x_real = np.real(x)
    x_imag = np.imag(x)
    y = np.array([x_real, x_imag]).astype(np.float)
    # re-order in convenient order
    if x.ndim >= 3:
        y = y.swapaxes(0, 1)
    return y

在这里插入图片描述
在这里插入图片描述
plt.imshow(x_real[0,:,:],cmap=‘Greys’)
在这里插入图片描述
在这里插入图片描述
plt.imshow(x_imag[0,:,:],cmap=‘Greys’)
在这里插入图片描述

在这里插入图片描述
k_und = data.squeeze(0)
在这里插入图片描述

        return {
            'image': data_und,
            'k': k_und.transpose(1,2,0),
            'mask': mask.transpose(1,2,0),
            'full': data_gnd
        }

sample = dataset[0]即4个tensor
第一个tensor。‘image’: data_und
第二个tensor。‘k’: k_und.transpose(1,2,0)
在这里插入图片描述
第三个tensor。‘mask’: mask.transpose(1,2,0)
在这里插入图片描述
第四个tensor。‘full’: data_gnd

输出：
Sample image shape: (2, 256, 256)
Sample full shape: (2, 256, 256)
在这里插入图片描述

敲代码的妹子
关注
专栏目录
CT重建软件包：执行ct重建任务的功能集-matlab开发
06-01
执行 ct 重建任务的函数集，例如氡变换、简单的反投影、使用空间域中的卷积的滤波反投影、使用 2D 傅立叶变换的滤波反投影以及使用中心切片定理的滤波反投影。 它包括功能myCtReconstruction 这应该是开箱即用的，并使用 Matlab 的 Shepp Logan Phantom 进行演示。 也可以使用参数运行 myCtReconstruction myCtRecontruction(anyCtDataSet) 在自定义数据集上执行图像重建。 参考： [1] Justin Romberg 教授的一组 CT 重建函数@ http://www.clear.rice.edu/elec431/projects96/DSP/ [2] W. Schlegel, J. Bille: Medizinische Physik 2: Medizinische Strahlenphysik h
open-cv，haarcascade_frontalface_default.xml人脸检测模型包
09-23
1. **加载模型**：首先，你需要创建一个`cv::CascadeClassifier`对象，并使用`load`函数加载XML文件。 ```cpp cv::CascadeClassifier face_cascade; face_cascade.load("path/to/haarcascade_frontalface_default....
cascade MRI reconstruction:evaluate
xuru_0927的博客
08-02 350
if draw: pred = pred[0] grnd = grnd[0] err = np.abs(pred - grnd) ** 2 err -= err.min() err /= err.max() err *= 255 err = err.astype(np.uint8) err = cv2.applyColorMap(err, 2) plt.imsave(os.path.join(kwargs['output_path'],
计算机视觉-三维重建等几何视觉PPT--网上谣传的python实现SFM代码整理-数据集ReconstructionDataSet
baidu_40840693的博客
02-11 2263
http://leohope.com/%E8%A7%A3%E9%97%AE%E9%A2%98/2018/03/06/compare-re3d-system/ 参考视频 先说一些问题，一个完整的三维重建流程需要多个工具组合： 三维重建数据集： https://github.com/rperrot/ReconstructionDataSet 网上实现SFM代码整理： 1 ...
论文阅读Reconstructing CT from Biplanar X-Rays with Generative Adversarial Networks
qq_42042757的博客
04-06 3442
Reconstructing CT from Biplanar X-Rays with Generative Adversarial Networks 基于生成对抗网络的双平面x线重建CT 来源：CVPR2019 作者单位：腾讯，浙江大学，上海交通大学 摘要 计算机断层扫描(CT)可以提供患者内部器官的3D视图，有助于疾病诊断，但它会给患者带来更多的辐射剂量，而且CT扫描仪的成本也比x光机高得多。传统的CT重建方法需要通过对身体的全面旋转扫描进行数百个x射线投影，这在典型的x射线机器上是无法完成的。传统的C
CT迭代重建ART matlab论坛,CT重建包
weixin_39792519的博客
04-09 305
CT reconstruction packageSet of functions performing ct reconstruction tasks like a radon transformation, a simple backprojection, a filtered backprojection using a convolution in the spatial domain, ...
Opencv实现人脸识别文件：haarcascade-frontalface-alt.xml
09-15
haarcascade_frontalface_alt.xml是OpenCV中的一个经典的人脸检测器模型文件。它基于Haar特征分类器，用于检测图像或视频中的人脸。 具体来说，该模型文件使用AdaBoost算法训练得到，其中包含了一系列的弱分类器，...
camera.py基于opencv的摄像头人脸检测代码
12-09
基于opencv的摄像头人脸检测代码，基于haarcascade_eye.xml和haarcascade_smile.xml增加了人眼检测和笑脸检测
人脸harr级联分类器保存的.xml文件，包括：haarcascade-frontalface-default.xml等文件
06-30
在这个场景下，我们关注的是OpenCV中的人脸检测技术，尤其是基于Haar特征的级联分类器，这正是`haarcascade-frontalface-default.xml`文件所代表的。 Haar特征是一种用于图像分析的简单而强大的工具，由Paul Viola...
haarcascade-frontalface-alt.xml.zip
最新发布
12-15
标题中的"haarcascade-frontalface-alt.xml.zip"是一个压缩包文件，其中包含了OpenCV库的一个面部检测模型。OpenCV（开源计算机视觉库）是一个强大的工具，广泛应用于图像处理和计算机视觉领域。这个模型，"haar...
进行遥感图像语义分割时无法正确识别标签的问题
qq_44178263的博客
04-07 2563
在利用语义分割进行遥感图像水体提取时发现的问题
KBQA-Bert学习记录-CRF模型
Swayzzu的博客
12-15 3445
CRF模型的实现
cascade MRI reconstruction: cartesian_mask--2
xuru_0927的博客
07-12 355
def cartesian_mask(self, shape): N, Nx, Ny = int(np.prod(shape[:-2])), shape[-2], shape[-1] pdf_x = normal_pdf(Nx, 0.5/(Nx/10.)**2) # shape：256 pdf_x.reshape(16,16) lmda = Nx/(2.*self.acc) # 32 n_lines = int(Nx / self.acc) # 64 # ad...
numpy.transpose难点解析
一呆飞仙的博客
06-27 551
先来看一个例子 import numpy as np arr1 = np.array(range(1, 37)) print 'arr1:','\n',arr1 arr2 = arr1.reshape(4, 3, 3) print 'arr2:','\n',arr2 arr3 = arr2.transpose(1，0，2) print 'arr3:','\n',arr3 arr4 = arr...
Polar Mask：如何使用极坐标表示方法对单镜头实例进行分割？
Spuer_Tiger
03-06 1686
PolarMask：如何使用极坐标表示对单镜头实例进行分割? ​ PolarMask: Single Shot Instance Segmentation with Polar Representation (作者自译，转载需本人经许可) ​ _arXiv:1909.13226v2 [cs.CV] 10.Oct.2019 Abstract： ​ 在本文中，我们介绍了...
sql语句中出现笛卡尔乘积 SQL查询入门篇（转）
mypanlong的专栏
12-23 2076
本篇文章中，主要说明SQL中的各种连接以及使用范围，以及更进一步的解释关系代数法和关系演算法对在同一条查询的不同思路。   多表连接简介 在关系数据库中，一个查询往往会涉及多个表，因为很少有数据库只有一个表，而如果大多查询只涉及到一个表的，那么那个表也往往低于第三范式，存在大量冗余和异常。 因此，连接(Join)就是一种把多个表连接成一个表的重要手段. 比如简单两个表连接学生表(Stu
关于deeplabv3+ 使用自定义数据集训练出现问题2
u010358801的博客
01-04 777
关于deeplabv3+ 使用自定义数据集训练出现问题2 前言 见 ： 原创 关于deeplabv3+ 使用自定义数据集训练出现的问题1 1 新出现的问题 1.1 问题简述 ​ 首先由于我使用的数据处理是从dataloaders/datasets/pascal.py 复制过来的，所以很多代码都是从这里变化过来的。正因为如此，而我们的数据集瓷瓦裂痕的数据集，MASK 图片颜色就两种，背景是RGB：(0,0,0)，识别项是RGB(255,255,255) ​ 从vo.
美国高校开源迄今为止最大新冠肺炎CT数据集！
zandaoguang的博客
04-12 7212
点击上方“视学算法”，选择“星标”快速获得最新干货 本文转载自：新智元编辑：元子、白峰【导读】胸部计算机断层扫描（CT）图像在对新冠肺炎（COVID-19）提供准确、快速、廉价的...
2021-3-22科研进展
xuru_0927的博客
03-21 292
对于FastMRI数据集Knee path = r'/Volumes/exterD/codeTemp/pg_mri/datasets_knee/singlecoil_train/file1000026.h5' target = h5py.File(path, 'r')['reconstruction_esc'] Output： <HDF5 dataset "reconstruction_esc": shape (30, 320, 320), type "<f4"> The invers
usage: fadiga.py [-h] -c CASCADE -p SHAPE_PREDICTOR [-a ALARM] fadiga.py: error: the following arguments are required: -c/--cascade
05-28
这是一个Python脚本的命令行参数错误提示信息。 根据提示信息，该脚本需要两个必需参数：`-c`或`--cascade`和`-p`或`--shape_predictor`。其中，`-c`或`--cascade`指定级联分类器的路径，`-p`或`--shape_predictor`指定人脸关键点检测器的路径。此外，该脚本还支持一个可选参数`-a`或`--alarm`，用于指定疲劳检测的警报阈值。 因此，您需要在命令行中正确指定这些参数，例如： ``` python fadiga.py -c path/to/cascade.xml -p path/to/shape_predictor.dat -a 20 ``` 其中，`path/to/cascade.xml`和`path/to/shape_predictor.dat`分别是级联分类器和人脸关键点检测器的路径，`20`是警报阈值。请根据实际情况修改这些路径和阈值。
评论
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付：10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明：

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、付费专栏及课程。

余额充值