Main Content

このページの翻訳は最新ではありません。ここをクリックして、英語の最新版を参照してください。

大きなイメージのワープ

この例では、大きなイメージに幾何学的変換 (ワーピング) を適用する方法を説明します。

イメージへの幾何学的変換の適用は、イメージ レジストレーションなどのような多数のイメージ処理アプリケーションの重要な手順です。メモリに収まる粗いイメージをワープさせるには、imwarp を使用します。メモリに収まらない大きな高解像度イメージの場合、イメージのブロックをワープさせます。ワープされたイメージの空間参照を、ピクセル範囲などのイメージの特性が維持されるように設定します。

CAMELYON16 データセットのイメージ "tumor_091.tif" の変更したバージョンを使用して bigimage を作成します。元のイメージは、腫瘍組織が含まれるリンパ節の学習イメージです。元のイメージには 8 つの解像度レベルがあり、最も細かいレベルの解像度は 53760 x 61440 です。変更したイメージには、3 つの粗い解像度レベルのみが含まれています。変更したイメージの空間参照は、縦横比が一定に維持され、各レベルで特徴がレジストレーションされるように調整されています。

bim = bigimage('tumor_091R.tif');

粗いイメージへの幾何学的変換の適用

アフィン幾何学的変換に関する情報を格納する affine2d オブジェクトを作成します。この変換は平行移動とせん断を適用します。

tform = affine2d([
    0.99 0.01 0
    0.17 0.98 0
    120 -30 1]);

最も粗い解像度レベルのイメージを取得します。

imCoarse = getFullLevel(bim);

imwarp を使用して粗いイメージをワープさせます。元のイメージとワープされたイメージをモンタージュに表示します。

imCoarseWarped = imwarp(imCoarse,tform);
montage({imCoarse,imCoarseWarped})

ワープされた細かいイメージの空間参照の作成

細かい解像度レベルのイメージに幾何学的変換を適用する前に、ワーピング後の大きなイメージの空間参照を計算します。ブロックを変換するときに、この空間参照を使用します。

その空間参照情報から元のイメージのピクセル範囲を取得します。

srcLevel = bim.FinestResolutionLevel;
inRef = bim.SpatialReferencing(srcLevel);
inPixelExtent = [inRef.PixelExtentInWorldX,inRef.PixelExtentInWorldY];

変換を適用するときの出力の水平方向と垂直方向の空間範囲を計算します。

[xout,yout] = outputLimits(tform,inRef.XWorldLimits,inRef.YWorldLimits);

ピクセル範囲を維持する出力イメージのサイズを計算します。[numrows, numcols] の形式でイメージ サイズを指定します。

outImgSize = [ceil(diff(yout)/inPixelExtent(2)),ceil(diff(xout)/inPixelExtent(1))];

ワープされたイメージの空間参照情報を格納する imref2d オブジェクトを作成します。ワープされたイメージのワールド座標範囲とイメージ サイズを設定します。

outRef = imref2d(outImgSize,xout,yout)
outRef = 
  imref2d with properties:

           XWorldLimits: [120.5800 6275]
           YWorldLimits: [-29.5050 4.9241e+03]
              ImageSize: [4954 6155]
    PixelExtentInWorldX: 0.9999
    PixelExtentInWorldY: 0.9999
    ImageExtentInWorldX: 6.1544e+03
    ImageExtentInWorldY: 4.9536e+03
       XIntrinsicLimits: [0.5000 6.1555e+03]
       YIntrinsicLimits: [0.5000 4.9545e+03]

出力の原点を左上のピクセルとして指定します。

outOrigin = [xout(1),yout(1)];

対応する出力ピクセルの次元を計算します。

outPixelExtent = [outRef.PixelExtentInWorldX,outRef.PixelExtentInWorldY];
halfPixWidth = outPixelExtent/2;

細かいイメージへのブロック単位のワーピングの適用

出力の空間参照情報を指定して、書き換え可能な bigimage を作成します。メモリを効率的に使用できる十分な大きさのブロック サイズを指定します。

outBlockSize = [1024 1024];
bwarped = bigimage(outRef,bim.Channels,bim.ClassUnderlying, ...
    'BlockSize',outBlockSize);

出力イメージを一度に 1 ブロックずつループ処理します。各出力ブロックに対して以下を行います。

  1. 出力ブロックの 4 つの隅の座標を求めます。

  2. これらの座標を入力に逆マッピングして、入力 (ソース) 領域を取得します。

  3. 入力領域のコンテンツを読み取ります。

  4. 入力領域を記述する空間参照を作成します。

  5. imwarp を使用して出力ブロックのコンテンツを計算します。

  6. 関数 setBlock を使用して出力ブロックを出力イメージに書き込みます。

Parallel Computing Toolbox™ がある場合、外側の for ステートメントを parfor ステートメントに置き換えて、このループを並列実行することができます。

colStarts = 1:outBlockSize(2):outRef.ImageSize(2);
for cInd = 1:numel(colStarts)
    colStart = colStarts(cInd);
    for rStart = 1:outBlockSize(1):outRef.ImageSize(1)
        
        % Center of top left pixel of this block in world units
        xyStart = [colStart,rStart].* outPixelExtent;
        xyStart = xyStart + outOrigin;
        
        % Center of bottom right pixel of this block in world units
        bsize = fliplr(bwarped.BlockSize); % (r,c) -> (x,y)
        xyEnd = ([colStart,rStart] + (bsize-1)).*outPixelExtent;
        xyEnd = xyEnd + outOrigin;
        
        % Output bounding box
        outbbox = [xyStart
            xyStart(1) xyEnd(2)
            xyEnd(1) xyStart(2)
            xyEnd];
        
        % Corresponding spatial reference which describes this rectangle
        outRegionRef = imref2d(bsize);
        outRegionRef.XWorldLimits = [xyStart(1)-halfPixWidth(1),xyEnd(1)+halfPixWidth(1)];
        outRegionRef.YWorldLimits = [xyStart(2)-halfPixWidth(2),xyEnd(2)+halfPixWidth(2)];
        
        % Get corresponding input region. The region is not rectangular if the transformation
        % includes shear
        inRegion = transformPointsInverse(tform,outbbox);   
        
        % Clamp to image extents (working with pixel centers)
        xcenterLims = inRef.XWorldLimits + [halfPixWidth(1)-halfPixWidth(1)];
        ycenterLims = inRef.YWorldLimits + [halfPixWidth(2)-halfPixWidth(2)];
        inRegion(:,1) = max(inRegion(:,1),xcenterLims(1));
        inRegion(:,1) = min(inRegion(:,1),xcenterLims(2));
        inRegion(:,2) = max(inRegion(:,2),ycenterLims(1));
        inRegion(:,2) = min(inRegion(:,2),ycenterLims(2));
        
        % Find the input bounding box (still using pixel centers)
        inbboxStart = [min(inRegion(:,1)) min(inRegion(:,2))];
        inbboxEnd = [max(inRegion(:,1)) max(inRegion(:,2))];
        
        % Read corresponding input region 
        inputRegion = getRegion(bim,srcLevel,inbboxStart,inbboxEnd);
        
        % Input region's spatial referencing
        inRegionRef = imref2d(size(inputRegion));
        % Convert to pixel edges from pixel centers
        inRegionRef.XWorldLimits = [inbboxStart(1)-halfPixWidth(1),inbboxEnd(1)+halfPixWidth(2)];
        inRegionRef.YWorldLimits = [inbboxStart(2)-halfPixWidth(1),inbboxEnd(2)+halfPixWidth(2)];
        
        warpedBlock = imwarp(inputRegion,inRegionRef,tform,'OutputView',outRegionRef);
        
        % Set the block data in the output bigimage
        setBlock(bwarped,1,xyStart,warpedBlock);
    end
end

ワープされたイメージを表示します。

bigimageshow(bwarped)

参考

| | | | | |