MATLAB中处理任意大图像的函数blockproc
MATLAB分块处理矩阵的函数很早就有了——blkproc，但缺点是blkproc必须一次性把要处理的矩阵全部导入内存中，这样大大限制了其应用范围，对一些超大型的图像就无能为力了。幸运的是随着MATLAB使用范围越来越广，Mathworks也紧跟用户需求，新版本的MATLAB推出了可以处理任意大图像的函数blockproc，其用法如下：[*]B = blockproc(A,,fun)
[*]B = blockproc(src_filename,,fun)
[*]B = blockproc(adapter,,fun)
[*]blockproc(...,param,val,...)
A是要处理的图像矩阵，如果图像太大不能完全导入内存，也可以用图像文件名src_filename来表示。是希望每次分块处理的矩阵大小，fun是函数句柄，即对每块矩阵的处理函数。
需要说明的是blockproc默认支持tiff/tif和jpeg2000格式的任意大图像处理，如果要读取其他格式的大图像需要针对该图像格式再写一个继承自MATLAB中ImageAdapter这个抽象类的子类adapter，来满足blockproc的输入要求。MATLAB帮助文档中有一个读取lan格式的LanAdapter示例类，大家可以参照那个格式来构造任意图像格式的Adapter类来实现blockproc函数对任意大图像文件的支持。
最后一种调用格式可以实现读取大图像文件，分块处理后再在指定路径写入处理后的图像文件，这个非常有用。下面给一个简单的例子,更多的用法希望大家讨论完善。
[*]
[*]fun = @(block_struct) repmat(block_struct.data,5,5) ;%注意处理函数的输入必须是结构体，其data域是我们的矩阵数据，这是由blockproc分块后的机制决定
[*]blockproc('moon.tif',,fun,'Destination','C:\moonmoon.tif');
[*]imshow('C:\moonmoon.tif')
可以看到分块处理后的效果，当然这是简单的把原图像分块（每一子块大小16*16）复制了25倍后的效果。
大家有没有试验一下，随便找个2G以上的tiff文件（手头没有可以利用这个函数把MATLAB自带的tiff文件repmat成2G左右）操作一下。我现在的疑问是，对于“操作本身”和“图像分块”彼此是独立的情形可以很方便利用blockproc，但是反之则不大方便，譬如对图像整体转置旋转等操作，如果简单用blockproc处理后得到的是每个块的转置和旋转。大家有好的办法没有？
问题难就难在一些大的图像（遥感、卫星影像可能达几G几十G）无法全部导入内存，因此先转置就比较困难了。在Mathworks新闻组我也问了这个问题，如下：
http://www.mathworks.com/matlabc ... /view_thread/298278
一个叫brendan的老外（不知道是不是Mathworks内部工程师）给出了他的办法，结果非常给力。简单说一下他的思路：
先构造一个“辅助图像矩阵”充当块计数器，利用blockproc遍历这个辅助图像的像素(其实对应的就是实际图像的每一个块)。充分利用用户自定义的那个函数，那个函数可以随心所欲，实现最后大图像相应的块处显示哪个图像。他给的那个自定义函数blockTranspose思路如下：
blockproc的分块机制会产生一个结构体，记录当前块的数据以及位置等信息，因此先构造一个100*100的辅助矩阵，分成100*100块，实际就是一个像素分一块，对应10000*10000大图像的100*100小块。blockTranspose接收辅助矩阵的每一像素，这个像素包含了其代表的原图像相应块的位置信息。然后调用imread函数实现读取大图像中相应位置的图像块实现转置。
当然Brendan也承认这种用法不同于blockproc的常规用法，甚至是剑走偏锋了，仔细看了他的回复后，不得不感叹思路之巧妙，自定义函数的用法之灵活，让人十分佩服！更加体会了思想有多远，MATLAB就能玩得有多神！下面给出我对moon.tif图像转置的示例（当然这个图像不大，可以全部读进内存再转置，这里是拿它来做例子，几个G的图像也可以按这种方法处理）
[*]function GlobalBlocDemo
[*]function output_block= myfun(block_struct)
[*]start_index = (block_struct.location - 1) * p + 1;
[*]end_index = start_index + p-1;
[*]output_rows = ;
[*]output_cols = ;
[*]input_block = imread('moon.tif','PixelRegion',{output_cols,output_rows});
[*]output_block = input_block';
[*]end
[*]info = imfinfo('moon.tif');
[*]p = gcd(info.Width,info.Height);
[*]blockproc(ones(info.Width/p,info.Height/p),,@myfun,'Destination','C:\moon2.tif')
[*]end
运行完毕后可以看下结果：
[*]imshow('moon.tif');
[*]figure;imshow('C:\moon2.tif')
由于分块矩阵的转置，在处理子矩阵的转置的同时，还是考虑到块本身的转置，
所以很难实现，的确brendan的程序完美的解答了这一问题，
而且构造辅助矩阵的想法确实体现了对matlab函数的深刻的理解，
将图像的处理放在自定义的函数中，确实给人一种高屋建瓴之感。
感谢roc对这个问题的深入讨论，使得我们在分块矩阵的转置有了新的认识。
-----------------------------------------------------------------------------------------------------------
我看了下blockproc，对没有小块的处理顺序，刚开始的时候感觉有点乱
好在blockproc的一些代码是可见的，才发现
= meshgrid(2:mblocks,2:nblocks);
rr = r(:);
cc = c(:);
= meshgrid(1,2:end_col);
rr = ;
cc = ;
= meshgrid(2:end_row,1);
rr = ;
cc = ;这些代码产生了执行小块的先后顺序，之所以这样，可以有通用性的考虑吧
为了能够理解4#roc给出的程序，我写了一个矩阵转置的小例子，供大家参考
主程序：
clear;clc;close all
m=2;% 子矩阵的行数
n=4;% 子矩阵的列数
M=3;% 大矩阵的块的行数
N=2;% 大矩阵的块的列数
maxM=3; % 子矩阵的最大数
A=randi(maxM,);% 生成子矩阵
%为了索引子矩阵处理方便，我们建立一个辅助的单元数组
Bcell=cell(M,N);
Bcell(:)=arrayfun(@(x)x*A,1:6,'UniformOutput',false);
Bnum=cell2mat(Bcell);
C=largeMatrixTransposeDemo(Bcell,M,N);
disp('原矩阵');
disp(Bnum);
disp('用blockproc函数转置后的矩阵');
disp(C);
disp('用matlab的transpose转置的结果');
disp(Bnum');调用的函数
function C=largeMatrixTransposeDemo(B,M,N)
function output_a= myfun(input_a)
pos=input_a.location;
output_a=B{pos(2),pos(1)}';
end
C=blockproc(zeros(N,M)+1,,@myfun);
end运行的结果：
原矩阵
3 3 3 2 12 12 12 8
1 2 1 1 4 8 4 4
6 6 6 4 15 15 15 10
2 4 2 2 5 10 5 5
9 9 9 6 18 18 18 12
3 6 3 3 6 12 6 6
用blockproc函数转置后的矩阵
3 1 6 2 9 3
3 2 6 4 9 6
3 1 6 2 9 3
2 1 4 2 6 3
12 4 15 5 18 6
12 8 15 10 18 12
12 4 15 5 18 6
8 4 10 5 12 6
用matlab的transpose转置的结果
3 1 6 2 9 3
3 2 6 4 9 6
3 1 6 2 9 3
2 1 4 2 6 3
12 4 15 5 18 6
12 8 15 10 18 12
12 4 15 5 18 6
8 4 10 5 12 6
页:
[1]