For use in creating the image transform needed to orthorectify overhead photographs taken in the Delta Basins, St. Anthony Falls Lab.
You will need:
- An overhead photograph (.jpg) of a series of crosshairs (at equal elevation) in (ideally) a regularly spaced grid in the area of the basin viewed by the camera, along with a list of their x,y basin coordinates.
- A copy of MATLAB
The procedure is as follows (the actual code is below the step, as needed)
1) Open MATLAB and load the image as a variable. Name this variable “skewed_img”.
2) Convert “skewed_img” to greyscale; name the new image “skewed_gray”
skewed_gray=rgb2gray(skewed_img);
3) Enter the list of coordinates of the crosshairs. Call this variable “base_pts”
4) Plot base_pts so that the plot consists of a series of crosshairs, and arrange the axes such that there is at least a 500 mm buffer between the axes (xmin will be the minimum value shown on the x axis, etc.) . Set the plot so that 0,0 is the upper left hand corner (as opposed to the lower left hand corner).
h=figure(1); clf, a=axes(‘Parent’,h);
plot(a,base_pts(:,1),base_pts(:,2),’+k’,’MarkerSize’,8);
set(a,’XLim’[xmin,xmax],’YLim’,[ymin,ymax]);
axis ij
5) Maximize the plot area and save as a .jpg image. Load the image back in, as “control_img”. Convert the image to greyscale as above; call the greyscale image “control_gray”.
6) Use the MATLAB built-in image registration tool cpselect to select pairs of points between the control image (the one you just made) and the skewed image (the one you took with the camera). The image registration tool will allow you to zoom and pan as needed to be precise about point selection. Select the points in the same order as they are listed in base_pts. Export the control points (both movingPoints and fixedPoints).
cpselect(skewed_gray,control_gray);
7) (Optional) if you are not able to select a certain crosshair as a control point (if, for example, you can’t actually see the crosshair), remove the coordinates of that crosshair from your list of positions.
8) Create a transformation matrix using movingPoints and base_pts:
transform=cp2tform(movingPoints,base_pts,’projective’);
9) Apply the transformation to your skewed image to check accuracy:
test_img=imtransform(skewed_img,transform);
imshow(test_img);
10) Save the transformation matrix for later use in an image correction script. Call it something descriptive and put it in a safe place. Or several safe places. Keep in mind that this transformation is only good for the camera, lens, and camera position you used to create the transform. If you modify any of these in any way, the transform won’t work right.
A note on imtransform (the MATALB command used here) vs. imwarp (the preferred, as of R2014a, command for image transformation calculations): they do not take the same inputs. The command imtransform looks for a 5-field struct, created by cp2tform(). It returns an image that has not been resized from the original, skewed, version.
Imwarp looks for a particular type of geometric transform object, which can be created using fitgeotrans(). (In using that function, you will need to specify the type as projective, as for cp2tform(). Check the documentation online for more details—it may be possible to feed it base_pts instead of fixedPoints, but this is untested.) The transformation created by fitgeotrans() does not preserve the size of the image. Instead, the corrected image will have the same dimensions as the image of the plot of the control points—often, this will be much smaller than the original image. Either because of this or because imwarp is simply more efficient, it runs much, much faster.
A simplified sample image correction script is below. Descriptive commentary is in green and prefaced by a % symbol; code is in black, blue, and purple.
%Asks user to select a folder containing distorted images, gets the directory
Images_Dir = uigetdir('','Select the Folder containing distorted overhead images');
%Creates a directory of all files with a .jpg format
Img_List = dir(strcat(Images_Dir,'\*.jpg'));
%Creates a directory of all files with a .jpg format
save_dir = strcat(Images_Dir, '\CorrectedImages\');
%Creates a directory of files with a jpg format for the new folder.
mkdir(save_dir);
% Loads the transformation information—note that the transform file must be saved to the same directory as this code
load transform
%Loops through the list of images, by determining the number of images in
%the directory. i acts as the index placemarker.
%defines the total number of images to be transformed
Img_tot=length(Img_List);
for i = 1:Img_tot
%Get the name of the image to be loaded
Img_Name = Img_List(i).name;
%Load the image into MATLAB
distorted_image = imread(strcat(Images_Dir, '\', Img_Name));
%Apply the transformation
corrected_image = imtransform(distorted_image,transform);
%Append _corrected to the old image name
Img_Name_corrected = strcat(Img_Name(1:(end-4)), '_corrected.jpg');
%Save the corrected image as a .jpg, with the new name
imwrite(corrected_image, strcat(save_dir, '\', ...
Img_Name_corrected), 'jpg')
end
