This post contains full Python source code for Photo Merge. Source by Michal Lohnický.
<pre>#!/usr/bin/python import cv ####################################### # Class to hold helper functions ####################################### class Helper(): image_set = 0 max_distance = 10000000 r_plane = None g_plane = None b_plane = None hist1 = None hist2 = None is_init = False color_bins = 2 color_range = [0, 255] h_plane = None s_plane = None h_bins = 4#30 s_bins = 4#32 hist1_hsv = None hist2_hsv = None @staticmethod def init(grid_size): if Helper.is_init == True: return Helper.is_init = True Helper.r_plane = cv.CreateMat(grid_size, grid_size, cv.CV_8UC1) Helper.g_plane = cv.CreateMat(grid_size, grid_size, cv.CV_8UC1) Helper.b_plane = cv.CreateMat(grid_size, grid_size, cv.CV_8UC1) color_bins = Helper.color_bins color_range = Helper.color_range Helper.hist1 = cv.CreateHist([color_bins, color_bins, color_bins], cv.CV_HIST_SPARSE, [color_range, color_range, color_range], 1) Helper.hist2 = cv.CreateHist([color_bins, color_bins, color_bins], cv.CV_HIST_SPARSE, [color_range, color_range, color_range], 1) Helper.h_plane = cv.CreateMat(grid_size, grid_size, cv.CV_8UC1) Helper.s_plane = cv.CreateMat(grid_size, grid_size, cv.CV_8UC1) h_bins = Helper.h_bins s_bins = Helper.s_bins ranges = [[0, 180],[0, 255]] Helper.hist1_hsv = cv.CreateHist([h_bins, s_bins], cv.CV_HIST_ARRAY, ranges, 1) Helper.hist2_hsv = cv.CreateHist([h_bins, s_bins], cv.CV_HIST_ARRAY, ranges, 1) Helper.hsv_img = cv.CreateImage((grid_size, grid_size), 8, 3) @staticmethod def resize(im, new_width): size = (new_width, int(im.height/float(im.width)*new_width) ) src_img_1_resized = cv.CreateImage(size, im.depth, im.nChannels) cv.Resize(im, src_img_1_resized) return src_img_1_resized @staticmethod def get_histogram_hsv(src, hist_used_ID=1): # Convert to HSV cv.CvtColor(src, Helper.hsv_img, cv.CV_BGR2HSV) # Extract the H and S planes cv.Split(Helper.hsv_img, Helper.h_plane, Helper.s_plane, None, None) planes = [Helper.h_plane, Helper.s_plane] if hist_used_ID==1: cv.ClearHist(Helper.hist1_hsv) cv.CalcHist([cv.GetImage(i) for i in planes], Helper.hist1_hsv) return Helper.hist1_hsv cv.ClearHist(Helper.hist2_hsv) cv.CalcHist([cv.GetImage(i) for i in planes], Helper.hist2_hsv) return Helper.hist2_hsv @staticmethod def get_histogram(src, hist_used_ID=1): cv.Split(src, Helper.r_plane, Helper.g_plane, Helper.b_plane, None) planes = [Helper.r_plane, Helper.g_plane, Helper.b_plane] if hist_used_ID==1: cv.ClearHist(Helper.hist1) cv.CalcHist([cv.GetImage(i) for i in planes], Helper.hist1) return Helper.hist1 cv.ClearHist(Helper.hist2) cv.CalcHist([cv.GetImage(i) for i in planes], Helper.hist2) return Helper.hist2 @staticmethod def print_histogram(hist): sum = 0 print "Printing histogram" for r in range(0,Helper.color_bins): for g in range(0,Helper.color_bins): for b in range(0,Helper.color_bins): intensity = cv.QueryHistValue_3D(hist,r,g,b) sum+=intensity if intensity>0.0000001: print r,g,b,":", intensity print "Overal histogram sum:",sum @staticmethod def print_histogram_double(hist1, hist2): sum1 = 0 sum2 = 0 print "Printing histogram" for r in range(0,Helper.color_bins): for g in range(0,Helper.color_bins): for b in range(0,Helper.color_bins): intensity1 = cv.QueryHistValue_3D(hist1,r,g,b) intensity2 = cv.QueryHistValue_3D(hist2,r,g,b) sum1+=intensity1 sum2+=intensity2 if intensity1>0.0000001 or intensity2>0.0000001: print r,g,b,":", intensity1, intensity2 print "Overal histogram sum:",sum1, sum2 mode_list = { (cv.IPL_DEPTH_8U, 3) : "RGBA", (cv.IPL_DEPTH_8U, 3) : "RGB", (cv.IPL_DEPTH_8U, 1) : "L", (cv.IPL_DEPTH_32F, 1) : "F" } mode_list_r = { "RGBA" : (cv.IPL_DEPTH_8U, 3), "RGB" : (cv.IPL_DEPTH_8U, 3), "L" : (cv.IPL_DEPTH_8U, 1), "F" : (cv.IPL_DEPTH_32F, 1) } @staticmethod def convertToPIL(cv_im): from PIL import Image mode = Helper.mode_list[(cv_im.depth, cv_im.nChannels)] return Image.fromstring(mode, cv.GetSize(cv_im), cv_im.tostring()), mode @staticmethod def convertFromPIL(pi): mode = Helper.mode_list_r[pi.mode] cv_im = cv.CreateImageHeader(pi.size, mode[0], mode[1]) cv.SetData(cv_im, pi.tostring(), pi.size[0]*mode[1]) return cv_im @staticmethod def reduce_colors(cv_im, colors_count): from PIL import Image pi, mode = Helper.convertToPIL(cv_im) pi = pi.convert("RGB").convert("P", palette=Image.ADAPTIVE, colors=colors_count) pi = pi.convert(mode) return Helper.convertFromPIL(pi) @staticmethod def pasteMask(im1,im2,mask): pi1, mode = Helper.convertToPIL(im1) pi2, mode = Helper.convertToPIL(im2) pimask, mode = Helper.convertToPIL(mask) pi1.paste(pi2,(0,0),pimask) return Helper.convertFromPIL(pi1) ####################################### # Structure to hold matched values ####################################### class SURFFeatureMatch(): keypoint_img_1 = None descriptor_img_1 = None keypoint_img_2 = None descriptor_img_2 = None distance = Helper.max_distance ####################################### # Match feature points of 2 images ####################################### class SURFFeaturesMatcher(): #========================= # Match feature points of 2 images (main function) #========================= def findPairs(self, img_proc_1, img_proc_2): matches = [] for i in range(0,len(img_proc_1.keypoints)): featureMatch = SURFFeatureMatch() featureMatch.keypoint_img_1 = img_proc_1.keypoints[i] featureMatch.descriptor_img_1 = img_proc_1.descriptors[i] #find nearest feature point for j in range(0,len(img_proc_2.keypoints)): if featureMatch.keypoint_img_1[1] == img_proc_2.keypoints[j][1]: #laplacians has to be same if abs(featureMatch.keypoint_img_1[3] - img_proc_2.keypoints[j][3])<45: #similar orientation if self.realSqrDistance(featureMatch.keypoint_img_1[0], img_proc_2.keypoints[j][0])<4000: dist = self.descriptorsDistance(featureMatch.descriptor_img_1, img_proc_2.descriptors[j], featureMatch.distance) if dist current_best: break return total_cost #========================= # Equlide distance TODO consider image size #========================= def realSqrDistance(self, p1, p2): return (p1[0]-p2[0])*(p1[0]-p2[0])+(p1[1]-p2[1])*(p1[1]-p2[1]) ####################################### # Image processing class ####################################### class ImageProcessor: keypoints = None descriptors = None src_img_orig = None src_img = None src_img_colored = None src_img_path = "" matched_features = None binded_img_proc = None top_feature_points_count = 500 grid_size = 10 grid_tamplateMatching_offset = 8 def __init__(self, path): self.homography = None self.src_img_path = path self.src_img = Helper.resize(cv.LoadImage(self.src_img_path, cv.CV_LOAD_IMAGE_GRAYSCALE), 800) self.src_img_colored = Helper.resize(cv.LoadImage(self.src_img_path, cv.CV_LOAD_IMAGE_COLOR), 800) self.src_img_colored = Helper.reduce_colors(self.src_img_colored,255) self.src_img_orig = Helper.resize(cv.LoadImage(self.src_img_path, cv.CV_LOAD_IMAGE_COLOR), 800) Helper.init(self.grid_size) def extractSURFFeatures(self): try: self.border_w = int(self.src_img.width * 0.1) self.border_h = int(self.src_img.height * 0.1) cv.SetImageROI(self.src_img, ( self.border_w, self.border_h, self.src_img.width-2*self.border_w, self.src_img.height-2*self.border_h ) ) (self.keypoints, self.descriptors) = cv.ExtractSURF(self.src_img, None, cv.CreateMemStorage(), (1, 500, 3, 4)) print self.src_img_path+":\n"+" keypoints: "+str(len(self.keypoints))+"\n descriptors: "+str(len(self.descriptors)) cv.ResetImageROI(self.src_img) except Exception, e: print e def findSURFFeaturesPairs(self,img_proc): if self.keypoints == None: self.extractSURFFeatures() if img_proc.keypoints == None: img_proc.extractSURFFeatures() self.binded_img_proc = img_proc featureMatcher = SURFFeaturesMatcher() self.matched_features = featureMatcher.findPairs(self,img_proc) def getImportantFeaturePoints(self, top_points_count): ret_1 = [] ret_2 = [] for featureMatch in self.matched_features[:top_points_count]: point_1 = (int(featureMatch.keypoint_img_1[0][0]+self.border_w), int(featureMatch.keypoint_img_1[0][1]+self.border_h)) point_2 = (int(featureMatch.keypoint_img_2[0][0]+self.border_w), int(featureMatch.keypoint_img_2[0][1]+self.border_h)) ret_1.append(point_1) ret_2.append(point_2) return ret_1, ret_2 def getHomography(self): if self.homography!=None: return self.homography srcKeypoints, dstKeypoints = self.getImportantFeaturePoints(self.top_feature_points_count) dim=len(srcKeypoints) srcMat=cv.CreateMat(dim,2,cv.CV_32FC1) dstMat=cv.CreateMat(dim,2,cv.CV_32FC1) for i in range(1,dim): srcMat[i,0]=srcKeypoints[i][0] srcMat[i,1]=srcKeypoints[i][1] dstMat[i,0]=dstKeypoints[i][0] dstMat[i,1]=dstKeypoints[i][1] h=cv.CreateMat(3,3,cv.CV_64F) cv.FindHomography(dstMat,srcMat, h,cv.CV_RANSAC,5) self.homography = h return h def warpImage(self): src_img_warped = cv.CloneImage(self.binded_img_proc.src_img_colored) src_img_orig_warped = cv.CloneImage(self.binded_img_proc.src_img_orig) homo = self.getHomography() cv.WarpPerspective(self.binded_img_proc.src_img_colored, src_img_warped, homo) cv.WarpPerspective(self.binded_img_proc.src_img_orig, src_img_orig_warped, homo) return src_img_warped, src_img_orig_warped def drawSURFFeatures(self): window_name = "SURF "+self.src_img_path cv.NamedWindow(window_name, 1) for ((x, y), laplacian, size, dir, hessian) in self.keypoints: radio = size*1.2/9.*2 color = (255, 0, 0) if radio < 3: radio = 2 color = (0, 255, 0) cv.Circle(self.src_img_colored, (int(self.border_w+x),int(self.border_h+y)), int(radio), (0,255,0)) cv.ShowImage(window_name, self.src_img_colored) def drawMatchedPoints(self): #merge photos merged_size = (self.src_img.width, self.src_img.height + self.binded_img_proc.src_img.height) correspond = cv.CreateImage( merged_size, self.src_img.depth, self.src_img.nChannels) cv.SetImageROI( correspond, ( 0, 0, self.src_img.width, self.src_img.height ) ) cv.Copy( self.src_img, correspond ) cv.SetImageROI( correspond, ( 0, self.src_img.height, correspond.width, correspond.height ) ) cv.Copy( self.binded_img_proc.src_img, correspond ) cv.ResetImageROI( correspond ) #create window window_name = "Matched features" cv.NamedWindow(window_name, 1) #draw lines points_1, points_2 = self.getImportantFeaturePoints(self.top_feature_points_count) for i in range(0,len(points_1)): point_2 = (int(points_2[i][0]), int(points_2[i][1]+ self.src_img.height)) cv.Line( correspond, points_1[i], point_2, 128 ) #show image cv.ShowImage(window_name, correspond) def findMid(self, histogram): #normalize offset = len(histogram)/4 hist = [] for i in range(offset,len(histogram)-offset): hist.append((histogram[i-2]*0.3+histogram[i-1]*0.7+histogram[i]+histogram[i+1]*0.7+histogram[i+2]*0.3)/3) maxi = max(hist) for i in range(0,len(hist)): hist[i]/=maxi #find maxes max1_val = 0 max1_ind = 0 max2_val = 0 max2_ind = 0 ind1 = 0 ind2 = len(hist)-1 while ind1max1_val: max1_val = hist[ind1] max1_ind = ind1 if hist[ind2]>max2_val: max2_val = hist[ind2] max2_ind = ind2 ind1+=1 ind2-=1 print max1_ind, max2_ind #find min mid = len(hist)/2 min_ind = 0 min_val = 10000000 for i in range(max1_ind,max2_ind): new_val = hist[i]#*((abs(mid-i)/float(mid))*0.3+0.7) print i, ":", hist[i] if new_val=mid_ind: cv.Rectangle(mask2, (0, 0), (grid_size, grid_size), 0, thickness=-1) cv.ResetImageROI(mask1) cv.ResetImageROI(mask2) element_width_half = grid_size*2 element = cv.CreateStructuringElementEx(element_width_half*2+1, element_width_half*2+1, element_width_half, element_width_half, cv.CV_SHAPE_RECT) cv.Erode(mask2, mask2, element, 1) cv.Dilate(mask2, mask2, element, 1) cv.Erode(mask2, mask2, element, 1) cv.Erode(mask1, mask1, element, 1) cv.Dilate(mask1, mask1, element, 1) cv.Erode(mask1, mask1, element, 1) #finalize cv.Smooth(mask1, mask1, cv.CV_BLUR, grid_size*2, grid_size*2) cv.Smooth(mask2, mask2, cv.CV_BLUR, grid_size*2, grid_size*2) mask1 = Helper.resize(mask1, img_orig1.width) mask2 = Helper.resize(mask2, img_orig2.width) img_res1 = Helper.pasteMask(img_orig1,img_orig2,mask1) img_res2 = Helper.pasteMask(img_orig1,img_orig2,mask2) img_proc_1.crop(img_res1) cv.SaveImage("res1"+str(Helper.image_set)+".png",img_res1) cv.SaveImage("res2"+str(Helper.image_set)+".png",img_res2) def crop(self,im): roz = [[0, im.width, 0, im.width],[0, 0, im.height, im.height]] p = [] for i in range(0,4): p.append(cv.CreateMat(1,3,cv.CV_64FC1)) p[i][0,0]=roz[0][i] p[i][0,1]=roz[1][i] p[i][0,2]=1 cv.MatMul(p[i],self.homography,p[i]) p[i][0,0]/=p[i][0,2] p[i][0,1]/=p[i][0,2] cropped = cv.CreateImage( (im.width, im.height), im.depth, im.nChannel) src_region = cvGetSubRect(image, opencv.cvRect(left, top, new_width, new_height) ) cvCopy(src_region, cropped) import sys if __name__ == "__main__": Helper.image_set = str(sys.argv[1]) image_set = Helper.image_set img_proc_1 = ImageProcessor(".\\test_imgs\\"+str(image_set).zfill(2)+"_01.JPG") img_proc_2 = ImageProcessor(".\\test_imgs\\"+str(image_set).zfill(2)+"_02.JPG") img_final = cv.CreateImage( (img_proc_1.src_img_colored.width, img_proc_1.src_img_colored.height), img_proc_1.src_img_colored.depth, img_proc_1.src_img_colored.nChannels) img_histograms = cv.CreateImage( (img_proc_1.src_img_colored.width, img_proc_1.src_img_colored.height), 8, 1) img_proc_1.extractSURFFeatures() img_proc_2.extractSURFFeatures() img_proc_1.findSURFFeaturesPairs(img_proc_2) img_proc_1.drawMatchedPoints() src_img_warped, src_img_orig_warped = img_proc_1.warpImage() histograms_diffs = [] templateMatching_offset = 8 templateMatching_result = cv.CreateImage((templateMatching_offset+1, templateMatching_offset+1), cv.IPL_DEPTH_32F, 1) grid_size = ImageProcessor.grid_size act_width = templateMatching_offset/2 while act_width+grid_size(maxi-mini)*0.5: grid_mask_cols.append(1) frame_thickness = 1 if first and act_width<400 and act_height<230 and act_width>300 and act_height>180: hist1 = Helper.get_histogram(img_proc_1.src_img_colored, hist_used_ID=1) hist2 = Helper.get_histogram(src_img_warped, hist_used_ID=2) frame_thickness = -1 first = False cv.Copy( img_proc_1.src_img_colored, img_final ) frame_color = 128 cv.Rectangle(img_histograms, (1, 1), (grid_size-1, grid_size-1), 0) else: grid_mask_cols.append(0) cv.Copy( src_img_warped, img_final ) grid_mask.append(grid_mask_cols) cv.ResetImageROI(img_proc_1.src_img_colored) cv.ResetImageROI(src_img_warped) cv.ResetImageROI(img_final) cv.ResetImageROI(img_histograms) img_proc_1.combinePictures(img_proc_1.src_img_orig, src_img_orig_warped, img_proc_1.src_img_colored, src_img_warped, grid_mask) cv.ShowImage("tmp",img_proc_1.src_img_colored) cv.ShowImage("tmp2",src_img_warped) cv.ShowImage("tmp1",img_histograms) cv.ShowImage("FINAL",img_final) cv.WaitKey() cv.DestroyAllWindows()</pre>