BrainWeb: Simulated Brain Database使用（原图以及GroundTruth切片可视化）和矢状冠状横断面

注：此篇文章参考脑图像的数据预处理
作者：凯鲁嘎吉 - 博客园 http://www.cnblogs.com/kailugaji/
由于我在做实验时需要脑部横断面，矢状，冠状的原始图片，以及ground truth，
重新梳理BrainWeb数据集的使用。

一、BrainWeb数据集简介？

BrainWeb: Simulated Brain Database
http://brainweb.bic.mni.mcgill.ca/brainweb/
数据集选自McGill大学Montreal神经所大脑成像中心的Brain Web反震脑部MR图像数据库。该数据库包含基于两种解剖模型的模拟脑MRI数据：正常和多发性硬化（MS），并且提供了人工合成三种模态(T1-, T2-, PD-)下的三维脑MR图像。图像中包含不同的扫描厚度、噪声以及偏移场，最重要的是且该数据库对于所有脑MR图像均提供了标准分割结果。
用户可对所需图片自定义进行下载：

1. 模态 T1,T2,PD

2. 协议 icmb

3. 脑图像名字：正常

4. 切片厚度 1mm,3mm,5mm,7mm,9mm

5. 噪声水平 pn0=0%;1%,3%,5%,7%,9%

6. 灰度不均匀水平 rf0=0%;20%,40%

二、使用步骤

1.文件下载

原始图片以ms病变脑数据库为例，模态= T1，协议= ICBM，幻像名称=正常，切片厚度= 1mm，噪声= 0％，INU = 0％：

groundtruth同样以以ms病变脑数据库为例

下载好之后就可以处理数据了

2.数据处理

得到第90切片

function g = readrawb(filename, num)
%函数readrawb(filename, num)中的第一个参数filename是欲读取的rawb文件的文件名，第二个参数num就是第多少张。
fid = fopen(filename);
%连续读取181*217*181个数据，这时候temp是一个长度为181*217*181的向量。
%先将rawb中的所有数据传递给temp数组
temp = fread(fid, 181 * 217 * 181);
%然后把它变成了一个181*217*181的数组
images = reshape(temp, 181 , 217, 181);  
%获取图片大小[xrange, yrange, zrange] = size(images);
%   不同维度（切面）的图的第num张切片data_volume(1:xrange,1:yrange) = images(:,:,num);
%     data_volume(1:xrange,1:zrange) = images(:,num ,:);
%     data_volume(1:yrange,1:zrange) = images(num,:,:);
g = data_volume;
fclose(fid);
end

function init_image(filename,num)
pic_type='.png';
save_path='E:\BrainImages\Origin_MS\';
% 函数init_image(filename,num)中的第一个参数filename是欲读取的rawb文件的文件名，第二个参数num就是第多少张。输出为原始图像，未处理
%例如：init_image('t1_icbm_normal_1mm_pn0_rf0.rawb','train.txt',90)， init_image('phantom_1.0mm_normal_csf.rawb','train.txt',90)
read=readrawb(filename, num);
% 旋转90°并显示出来
read=imrotate(read, 90);                                       
imshow(mat2gray(read));
imwrite(mat2gray(read),strcat(save_path,'origin_1_t3_90',pic_type));
end

命令行输入init_image('t1_ai_msles2_1mm_pn0_rf0.rawb', 90)即可获得结果
想要其他维度切片时只需将readrawb中的注释更改即可

GroundTruth

T1模态、icmb协议下，切片厚度为1mm，噪声水平为0，灰度不均匀水平为0的ms病变脑图像，第90层
只选取0、1、2、3、10类的数据，各类别官网有介绍
分割之后各类用以下像素值进行渲染：
0：0
1：50
2：150
3：255
10: 230

function Ground_truth(name, num)
%标准分割结果
%例如：Ground_truth('Ground truth.txt',90)
mark=Mark('phantom_1.0mm_msles2_crisp.rawb',num);
for i=1:181   %行for j=1:217    %列if mark(i,j)==1read_new(i,j)=50;elseif mark(i, j)==2read_new(i,j)=150;elseif mark(i, j)==3read_new(i,j)=255;elseif mark(i,j)==10;read_new(i,j)=230elseread_new(i,j)=0;endend
end
% 旋转90°并显示出来
read_new=imrotate(read_new, 90);
%mark=mark';
imshow(mat2gray(read_new));
pic_type='.png';
save_path='E:\BrainImages\GroundTruth\';
imwrite(mat2gray(read_new),strcat(save_path,'gt_1_t1_90',pic_type));

function mark=Mark(filename,num)
%将标签为1、2、3、10类分出来，其余为0，mark取值：0、1、2、3、10
%[mark_new,mark]=Mark('phantom_1.0mm_normal_crisp.rawb',90);
fp=fopen(filename);
temp=fread(fp, 181 * 217 * 181);
images=reshape(temp, 181 ,217, 181);  [xrange, yrange, zrange] = size(images);%不同方向的ground truthdata_volume(1:xrange,1:yrange) = images(:,:,num); %横截面
%     data_volume(1:xrange,1:zrange) = images(:,num ,:);%冠状
%     data_volume(1:yrange,1:zrange) = images(num,:,:);%矢状
mark_data=data_volume;
fclose(fp);
%将第0、1、2、3、10类标签所在的坐标点拿出来，其余置0
for i=1:181for j=1:217if (mark_data(i,j)==1)||(mark_data(i,j)==2)||(mark_data(i,j)==3)||(mark_data(i,j)==10)mark(i,j)=mark_data(i,j);elsemark(i,j)=0;endend
end

运行Ground_truth('Ground truth.txt',90)
即可得到结果，同样，需要其他维度切片，更改Mask.m中的注释即可。