regi.rar_NMI 配准_NMI matlab_nmi_nmi 配准_互信息

clear; clc tic % 图像路径 ReferPath = 'F:\手工分割心脏\配准前\f1\5(0.75)\'; FloatPath = 'F:\手工分割心脏\配准前\f2\5(0.75)\'; %读取三维图像 %Refer为参考图像,Floating为浮动图像 %参考图像 tempPath = ls(ReferPath); tempPath(1:2, :) = []; Refer = []; for i = 1:7:size(tempPath, 1) temp = dicomread([ReferPath tempPath(i, :)]); %合并为三维图像 Refer = cat(3, Refer, temp); end Refer = double(Refer); MinRefer = min(min(min(Refer))); % figure; % imshow(Refer(:,:,1),[]); %获取参考图像的DICOM信息:层厚、像素大小 tempinfor = dicominfo([ReferPath tempPath(1, :)]); ReferPixelSpacing = tempinfor.PixelSpacing; if size(tempPath, 1) == 1 ReferSliceThickness = tempinfor.SliceThickness; else tempinfor2 = dicominfo([ReferPath tempPath(2, :)]); ReferSliceThickness = abs(tempinfor.SliceLocation - tempinfor2.SliceLocation); end %位深度(假设两个序列的位深度一样) BitDepth = double(tempinfor.BitDepth); % 浮动图像 tempPath = ls(FloatPath); tempPath(1:2, :) = []; Floating = []; for i = 1:7:size(tempPath, 1) temp = dicomread([FloatPath tempPath(i, :)]); %合并为三维图像 Floating = cat(3, Floating, temp); end Floating = double(Floating); MinFloat = min(min(min(Floating))); % figure; % imshow(Floating(:,:,1),[]); %获取浮动图像的DICOM信息:层厚、像素大小 tempinfof = dicominfo([FloatPath tempPath(1, :)]); FloatPixelSpacing = tempinfof.PixelSpacing; if size(tempPath, 1) == 1 FloatSliceThickness = tempinfof.SliceThickness; else tempinfof2 = dicominfo([FloatPath tempPath(2, :)]); FloatSliceThickness = abs(tempinfof.SliceLocation - tempinfof2.SliceLocation); end clear i temp tempPath FloatPath ReferPath; %重采样 %找出像素代表的最小实际大小 PixelSize = min(min(min(FloatPixelSpacing), FloatSliceThickness),... min(min(ReferPixelSpacing), ReferSliceThickness)); %参考图像x方向 Refer = ImageResample(Refer, ReferPixelSpacing(1), PixelSize); %参考图像y方向 Refer = permute(Refer, [2 1 3]); Refer = ImageResample(Refer, ReferPixelSpacing(2), PixelSize); Refer = permute(Refer, [2 1 3]); %参考图像z方向 Refer = permute(Refer, [3 2 1]); Refer = ImageResample(Refer, ReferSliceThickness, PixelSize); Refer = permute(Refer, [3 2 1]); %浮动图像x方向 Floating = ImageResample(Floating, FloatPixelSpacing(1), PixelSize); %浮动图像y方向 Floating = permute(Floating, [2 1 3]); Floating = ImageResample(Floating, FloatPixelSpacing(2), PixelSize); Floating = permute(Floating, [2 1 3]); %浮动图像z方向 Floating = permute(Floating, [3 2 1]); Floating = ImageResample(Floating, FloatSliceThickness, PixelSize); Floating = permute(Floating, [3 2 1]); clear FloatPixelSpacing FloatSliceThickness ReferPixelSpacing ... ReferSliceThickness; % 重采样后把超出深度最大值的点约束到深度范围内 % 先考虑数值都为非负数 % 参考图像 Refer = round(Refer); id = find(Refer > 2^BitDepth-1); Refer(id) = 2^BitDepth-1; id = find(Refer < MinRefer); Refer(id) = MinRefer; %浮动图像 Floating = round(Floating); id = find(Floating > 2^BitDepth-1); Floating(id) = 2^BitDepth-1; id = find(Floating < MinFloat); Floating(id) = MinFloat; clear id MinRefer MinFloat; %加载数据 BitDepth = 12; load('F:\手工分割心脏\原m\f2.mat'); load('F:\手工分割心脏\原m\r.mat'); figure('NumberTitle', 'off', 'Name', 'Floating'); imshow(Floating(:,:,1),[]); figure('NumberTitle', 'off', 'Name', 'Reference'); imshow(Refer(:,:,1),[]); NMI = NormalisedMutualInfomation(Refer, Floating, BitDepth); %进行仿射变换 %NewFloat为仿射变换后的浮动图像 %仿射变换参数 %t = [tx, ty, tz]是平移参数 %s是缩放参数:x,y,z方向上的缩放系数均相同,以保持不变性 %theata = [thx, thy, thz]是旋转参数 %以下是初始化参数 t = [0 0 0]; s = 0.85; theata = [0 0 0]; Floating = AffineTransform( Floating, t, s, theata); figure; imshow(Floating(:,:,10),[]); clear t s theata; NMI1 = NormalisedMutualInfomation(Refer, Floating, BitDepth); load('F:\手工分割心脏\原m\f1.mat'); load('F:\手工分割心脏\原m\r.mat'); load('F:\手工分割心脏\原m\infor.mat'); load('F:\手工分割心脏\原m\infof.mat'); load('F:\手工分割心脏\原m\PixelSize.mat'); % 修改图像像素间距信息 tempinfof.PixelSpacing = PixelSize; tempinfor.PixelSpacing = PixelSize; % 将矩阵由double转成uint16型 Floating = uint16(Floating); Refer = uint16(Refer); %输出dicom图像 for i = 1:size(Refer,3) dicomwrite(Refer(:,:,i),[num2str(i),'r.dcm'],tempinfor); tempinfor2.SliceLocation = tempinfor2.SliceLocation + PixelSize; end for i = 1:size(Floating,3) dicomwrite(Floating(:,:,i),[num2str(i),'f.dcm'],tempinfof); tempinfof2.SliceLocation = tempinfof2.SliceLocation + PixelSize; end toc