% imOut = AffineTransformImage_fa( im, T, szOut, METHOD ) % ------------------------------------------------------------------------- % % Applies a uniform scaling, oversampling or decimating, to an image, % using flat aperture interpolated values. % % The origin is preserved at the output image center. The resampling % grid is also centered on the input image. If an image dimension is even % a sampling point will occur at the center of the center pixel, and if % it is even sampling points will straddle the boundary between the % central pixels. % % Note: Best for decimation or flat aperture model oversampling. % See TransformImage.m for discrete sampling model. % % USAGE: imOut = AffineTransformImage_fa( im, [ .7 pi/6 ], size(im), 'linear' ); % imOut = AffineTransformImage_fa( im, [ .7 3/2 ], size(im) ); % imOut = AffineTransformImage_fa( im, [ .7 pi/6 ] ); % % INPUT: % % im - The image to be transformed. % T - The 3x3 homogeneous transformation matrix. Must be % affine. % METHOD - Interpolation method, same options as interp2_fa. % Default is 'quadratic'. % 'linear' linear interpolation % 'quadratic' piecewise cubic interpolation % szOut - Desired size of the transformed image. % Output pixels that are not mapped to the input will % be filled with NaN values. % Default is the size that includes all valid pixels % 0: Size of the input image, which will % result in cropping or a blank surround. % 1: Size that includes all valid pixels, the % default. % [ny mx]: Custom size. Output pixels that are not % mapped to the input will be filled with % NaN values. % % OUTPUT: % % imOut - The transformed image. % newT - The transformation matrix that relates transformed image % coordinates to the reference coordinates for use in a % function such as DIGIPLANE. % % NOTE: Input and output image origins are assumed to be at the geometric % centers of the full images (not the input region). % % NOTE: The region argument is used when one is inverting a perspective % transformation of a plane and the vanishing line of the plane lies % within the image. Attempts to transform any part of the vanishing line % will position you at infinity. Accordingly one should specify a region % that excludes any part of the vanishing line. % % NOTE: The szOut parameter is optionally used to control the size of the % output image. When inverting a perpective or affine transformation % the scale parameter is unknown/arbitrary, and without specifying % it explicitly the transformed image can end up being very small % or very large. % % NOTE: If your transformed image ends up as being two small bits of % image separated by a large black area then the chances are that you have % included the vanishing line of the plane within the specified region to % transform. If your image degenerates to a very thin triangular shape % part of your region is probably very close to the vanishing line of the % plane. % % FUNCTION CALLS: % % interp2_fa.m % -> interp1_fa.m % % % Mark Dow, http://lcni.uoregon.edu/~mark, dow@uoregon.edu % Mark Dow, April 10, 2007, derived from ScaleImage_fa % function imOut = AffineTransformImage_fa( im, T, szOut, METHOD ); if isa( im, 'uint8' ) im = double(im); % Make sure image is double end if nargin == 2 szOut = 1; METHOD = 'quadratic'; end if nargin == 3 METHOD = 'quadratic'; end % To Do: szOut sanity check. if szOut == 1 % size set based on resampling grid end if szOut == 0 szOut = size(im); end % To Do: transform sanity check % To Do: two element scale rotation if length(T) == 2 flScale = T(1); flAngle = T(2); T(1,1) = flScale*cos(flAngle); T(1,2) = -flScale*sin(flAngle); T(1,3) = 0; T(2,1) = flScale*sin(flAngle); T(2,2) = flScale*cos(flAngle); T(2,3) = 0; T(3,1) = 0; T(3,2) = 0; T(3,3) = 1; end % To Do: METHOD sanity check % To Do: image dimension sanity check bThreeD = ( ndims(im) == 3 ); % A colour image if bThreeD % Transform red, green, blue components separately. im = double(im)/255; r = DoTransformImage( im(:,:,1), T, szOut, METHOD ); g = DoTransformImage( im(:,:,2), T, szOut, METHOD ); b = DoTransformImage( im(:,:,3), T, szOut, METHOD ); imOut = repmat( uint8(0), [size(r),3] ); imOut(:,:,1) = uint8( round(r*255) ); imOut(:,:,2) = uint8( round(g*255) ); imOut(:,:,3) = uint8( round(b*255) ); else % Assume the image is greyscale. imOut = DoTransformImage( im, T, szOut, METHOD ); end %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % imOut = DoTransformImage( im, T, szOut, METHOD ) % The internal function that does all the work. function imOut = DoTransformImage( im, T, szOut, METHOD ); [ rows, cols ] = size(im); if szOut == 1 % Find the bounds of the new image. rcOut = bounds( T, size(im) ); szOut = [ rcOut(2) - rcOut(1) rcOut(4) - rcOut(3)]; end imOut( 1:szOut(2), 1:szOut(1) ) = NaN; imN( 1:szOut(2), 1:szOut(1) ) = 0; % To Do: specify center based on odd even % nOversampling = 8; nOversamplingMin = 4; TScaling = sqrt( T( 1, 1 )*T( 1, 1 ) + T( 2, 1 )*T( 2, 1 ) ); % To Do: should this be integral? nOversampling = min( nOversamplingMin, ceil( nOversamplingMin/TScaling ) ); TScalingOversampled = nOversampling*TScaling; %%%%%%%Note: To Do: remove this memory problem fudge if TScalingOversampled >= 2.2 TScalingOversampled = 2.2; nOversampling = TScalingOversampled/TScaling; end if TScalingOversampled ~= 1 if TScalingOversampled > 1 imOverSampled = ScaleImage_fa( im', TScalingOversampled, ceil(sqrt(2)*nOversampling*szOut), METHOD ); else imOverSampled = ScaleImage_fa( im', TScalingOversampled, 1, METHOD ); end else imOverSampled = im'; end TO = (1/TScalingOversampled)*T; TO(3,3) = 1; for i = 1 : size( imOverSampled , 1 ) for j = 1 : size( imOverSampled , 2 ) ii = i - ( size( imOverSampled , 1 ) + 1 )/2; jj = j - ( size( imOverSampled , 2 ) + 1 )/2; % Transform these xy coords to determine position on transformed image. px = floor( TO(1,1)*ii + TO(1,2)*jj + TO(1,3) + ( size( imOut , 1 ) + 1 )/2 + .5 ); py = floor( TO(2,1)*ii + TO(2,2)*jj + TO(2,3) + ( size( imOut , 2 ) + 1 )/2 + .5 ); if px > 0 & px <= size(imOut, 1) & py > 0 & py <= size(imOut, 2) if isnan( imOut( px, py ) ) imOut( px, py ) = 0; end imOut( px, py ) = imOut( px, py ) + imOverSampled( i, j ); imN( px, py ) = imN( px, py ) + 1; end end end for i = 1:size( imOut , 1 ) for j = 1:size( imOut , 2 ) % if imN(i,j) > 0 % if imN(i,j) > nOversampling*nOversampling/2 if imN(i,j) >= 1 imOut(i,j) = imOut(i,j)/imN(i,j); else imOut(i,j) = NaN; end end end imOut = imOut'; % % imOut = transformImage( im, T, region, szOut ) %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % B = bounds( T, szOut ) % % Internal function to find where the corners of a region, R % defined by [minrow maxrow mincol maxcol] are transformed to % by transform T and returns the bounds, B in the form % [minrow maxrow mincol maxcol] function B = bounds( T, szOut ) P = [ 1 szOut(2) szOut(2) 1 % homogeneous coords of region corners 1 1 szOut(1) szOut(1) 1 1 1 1 ]; % P = [ R(3) R(4) R(4) R(3) % homogeneous coords of region corners % R(1) R(1) R(2) R(2) % 1 1 1 1 ]; PT = round( homoTrans( T, P ) ); B = [ min( PT(2,:) ) max( PT(2,:) ) min( PT(1,:) ) max( PT(1,:) ) ]; % minrow maxrow mincol maxcol % % B = bounds( T, szOut ) %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%