全景图拼接——opencv之stitching
1 """ 2 Stitching sample (advanced) 3 =========================== 4 5 Show how to use Stitcher API from python. 6 """ 7 8 # Python 2/3 compatibility 9 from __future__ import print_function 10 11 import numpy as np 12 import cv2 as cv 13 14 import sys 15 import argparse 16 #创建parser实例 17 parser = argparse.ArgumentParser(prog='stitching_detailed.py', description='Rotation model images stitcher') 18 #继续增加参数 19 parser.add_argument('img_names', nargs='+',help='files to stitch',type=str) 20 #两个必须的入参,执行文件名,图片名 21 parser.add_argument('--preview',help='Run stitching in the preview mode. Works faster than usual mode but output image will have lower resolution.',type=bool,dest = 'preview' ) 22 parser.add_argument('--try_cuda',action = 'store', default = False,help='Try to use CUDA. The default value is no. All default values are for CPU mode.',type=bool,dest = 'try_cuda' ) 23 #默认不用cuda 24 #存储值到变量,action = store 25 parser.add_argument('--work_megapix',action = 'store', default = 0.6,help=' Resolution for image registration step. The default is 0.6 Mpx',type=float,dest = 'work_megapix' ) 26 parser.add_argument('--features',action = 'store', default = 'orb',help='Type of features used for images matching. The default is orb.',type=str,dest = 'features' ) 27 parser.add_argument('--matcher',action = 'store', default = 'homography',help='Matcher used for pairwise image matching.',type=str,dest = 'matcher' ) 28 parser.add_argument('--estimator',action = 'store', default = 'homography',help='Type of estimator used for transformation estimation.',type=str,dest = 'estimator' ) 29 parser.add_argument('--match_conf',action = 'store', default = 0.3,help='Confidence for feature matching step. The default is 0.65 for surf and 0.3 for orb.',type=float,dest = 'match_conf' ) 30 parser.add_argument('--conf_thresh',action = 'store', default = 1.0,help='Threshold for two images are from the same panorama confidence.The default is 1.0.',type=float,dest = 'conf_thresh' ) 31 parser.add_argument('--ba',action = 'store', default = 'ray',help='Bundle adjustment cost function. The default is ray.',type=str,dest = 'ba' ) 32 parser.add_argument('--ba_refine_mask',action = 'store', default = 'xxxxx',help='Set refinement mask for bundle adjustment. mask is "xxxxx"',type=str,dest = 'ba_refine_mask' ) 33 parser.add_argument('--wave_correct',action = 'store', default = 'horiz',help='Perform wave effect correction. The default is "horiz"',type=str,dest = 'wave_correct' ) 34 parser.add_argument('--save_graph',action = 'store', default = None,help='Save matches graph represented in DOT language to <file_name> file.',type=str,dest = 'save_graph' ) 35 parser.add_argument('--warp',action = 'store', default = 'plane',help='Warp surface type. The default is "spherical".',type=str,dest = 'warp' ) 36 parser.add_argument('--seam_megapix',action = 'store', default = 0.1,help=' Resolution for seam estimation step. The default is 0.1 Mpx.',type=float,dest = 'seam_megapix' ) 37 #缝合线尺度0.1M大概10万像素,320*320 38 parser.add_argument('--seam',action = 'store', default = 'no',help='Seam estimation method. The default is "gc_color".',type=str,dest = 'seam' ) 39 parser.add_argument('--compose_megapix',action = 'store', default = -1,help='Resolution for compositing step. Use -1 for original resolution.',type=float,dest = 'compose_megapix' ) 40 parser.add_argument('--expos_comp',action = 'store', default = 'no',help='Exposure compensation method. The default is "gain_blocks".',type=str,dest = 'expos_comp' ) 41 parser.add_argument('--expos_comp_nr_feeds',action = 'store', default = 1,help='Number of exposure compensation feed.',type=np.int32,dest = 'expos_comp_nr_feeds' ) 42 parser.add_argument('--expos_comp_nr_filtering',action = 'store', default = 2,help='Number of filtering iterations of the exposure compensation gains',type=float,dest = 'expos_comp_nr_filtering' ) 43 parser.add_argument('--expos_comp_block_size',action = 'store', default = 32,help='BLock size in pixels used by the exposure compensator.',type=np.int32,dest = 'expos_comp_block_size' ) 44 parser.add_argument('--blend',action = 'store', default = 'multiband',help='Blending method. The default is "multiband".',type=str,dest = 'blend' ) 45 parser.add_argument('--blend_strength',action = 'store', default = 5,help='Blending strength from [0,100] range.',type=np.int32,dest = 'blend_strength' ) 46 parser.add_argument('--output',action = 'store', default = 'result.jpg',help='The default is "result.jpg"',type=str,dest = 'output' ) 47 parser.add_argument('--timelapse',action = 'store', default = None,help='Output warped images separately as frames of a time lapse movie, with "fixed_" prepended to input file names.',type=str,dest = 'timelapse' ) 48 parser.add_argument('--rangewidth',action = 'store', default = -1,help='uses range_width to limit number of images to match with.',type=int,dest = 'rangewidth' ) 49 50 __doc__ += '\n' + parser.format_help() 51 # __doc__可以用函数名调用,会输出函数下的"""xxxxxx """三引号注释中的内容 52 53 def main(): 54 args = parser.parse_args() 55 #参数解析 56 #参数应用 57 img_names=args.img_names 58 #这种传参方法才是合理的,因为服务器传来的params是一个整体的json文件。因此参数解析要用类字典的形式。来的是整体,即便只有一个参数也不能直接类似params = 1如此传参 59 print(img_names) 60 preview = args.preview 61 try_cuda = args.try_cuda 62 work_megapix = args.work_megapix 63 seam_megapix = args.seam_megapix 64 compose_megapix = args.compose_megapix 65 conf_thresh = args.conf_thresh 66 features_type = args.features 67 matcher_type = args.matcher 68 estimator_type = args.estimator 69 ba_cost_func = args.ba 70 ba_refine_mask = args.ba_refine_mask 71 wave_correct = args.wave_correct 72 if wave_correct=='no': 73 do_wave_correct= False 74 else: 75 do_wave_correct=True 76 if args.save_graph is None: 77 save_graph = False 78 else: 79 save_graph =True 80 save_graph_to = args.save_graph 81 warp_type = args.warp 82 #曝光补偿供4种方式:增量补偿,块增量补偿,通道补偿,块通道补偿 83 if args.expos_comp=='no': 84 expos_comp_type = cv.detail.ExposureCompensator_NO 85 elif args.expos_comp=='gain': 86 expos_comp_type = cv.detail.ExposureCompensator_GAIN 87 elif args.expos_comp=='gain_blocks': 88 expos_comp_type = cv.detail.ExposureCompensator_GAIN_BLOCKS 89 elif args.expos_comp=='channel': 90 expos_comp_type = cv.detail.ExposureCompensator_CHANNELS 91 elif args.expos_comp=='channel_blocks': 92 expos_comp_type = cv.detail.ExposureCompensator_CHANNELS_BLOCKS 93 else: 94 print("Bad exposure compensation method") 95 exit() 96 expos_comp_nr_feeds = args.expos_comp_nr_feeds 97 expos_comp_nr_filtering = args.expos_comp_nr_filtering 98 expos_comp_block_size = args.expos_comp_block_size 99 match_conf = args.match_conf 100 seam_find_type = args.seam 101 blend_type = args.blend 102 blend_strength = args.blend_strength 103 result_name = args.output 104 if args.timelapse is not None: 105 timelapse = True 106 if args.timelapse=="as_is": 107 timelapse_type = cv.detail.Timelapser_AS_IS 108 elif args.timelapse=="crop": 109 timelapse_type = cv.detail.Timelapser_CROP 110 else: 111 print("Bad timelapse method") 112 exit() 113 else: 114 timelapse= False 115 range_width = args.rangewidth 116 #特征描述方法,三种方式 orb surf sift, 参数解析默认值是orb 117 if features_type=='orb': 118 finder= cv.ORB.create() 119 elif features_type=='surf': 120 finder= cv.xfeatures2d_SURF.create() 121 elif features_type=='sift': 122 finder= cv.xfeatures2d_SIFT.create() 123 else: 124 print ("Unknown descriptor type") 125 exit() 126 seam_work_aspect = 1 127 full_img_sizes=[] 128 features=[] 129 images=[] 130 is_work_scale_set = False 131 is_seam_scale_set = False 132 is_compose_scale_set = False; 133 for name in img_names:#循环 读图放到list中 134 full_img = cv.imread(cv.samples.findFile(name)) 135 136 if full_img is None: 137 print("Cannot read image ", name) 138 exit() 139 full_img_sizes.append((full_img.shape[1],full_img.shape[0])) 140 if work_megapix < 0: 141 img = full_img 142 work_scale = 1 143 is_work_scale_set = True 144 else: 145 if is_work_scale_set is False: 146 work_scale = min(1.0, np.sqrt(work_megapix * 1e6 / (full_img.shape[0]*full_img.shape[1]))) 147 is_work_scale_set = True 148 img = cv.resize(src=full_img, dsize=None, fx=work_scale, fy=work_scale, interpolation=cv.INTER_LINEAR_EXACT) 149 if is_seam_scale_set is False: 150 seam_scale = min(1.0, np.sqrt(seam_megapix * 1e6 / (full_img.shape[0]*full_img.shape[1]))) 151 seam_work_aspect = seam_scale / work_scale 152 #缝合区域大小/配准 153 is_seam_scale_set = True 154 imgFea= cv.detail.computeImageFeatures2(finder,img) 155 #单张图片的特征提取,finder是特征提取器对象 156 features.append(imgFea) 157 #根据以上来调整图片大小,保存到images图片列表中 158 img = cv.resize(src=full_img, dsize=None, fx=seam_scale, fy=seam_scale, interpolation=cv.INTER_LINEAR_EXACT) 159 images.append(img) 160 if matcher_type== "affine": 161 #匹配类型 仿射 162 matcher = cv.detail_AffineBestOf2NearestMatcher(False, try_cuda, match_conf) 163 elif range_width==-1: 164 matcher = cv.detail.BestOf2NearestMatcher_create(try_cuda, match_conf) 165 else: 166 matcher = cv.detail.BestOf2NearestRangeMatcher_create(range_width, try_cuda, match_conf) 167 #此时的features是所有入参图片的特征集合 168 p=matcher.apply2(features) 169 matcher.collectGarbage() 170 if save_graph: 171 f = open(save_graph_to,"w") 172 f.write(cv.detail.matchesGraphAsString(img_names, p, conf_thresh)) 173 f.close() 174 indices=cv.detail.leaveBiggestComponent(features,p,0.3) 175 img_subset =[] 176 img_names_subset=[] 177 full_img_sizes_subset=[] 178 num_images=len(indices) 179 for i in range(len(indices)): 180 img_names_subset.append(img_names[indices[i,0]]) 181 img_subset.append(images[indices[i,0]]) 182 full_img_sizes_subset.append(full_img_sizes[indices[i,0]]) 183 images = img_subset; 184 img_names = img_names_subset; 185 full_img_sizes = full_img_sizes_subset; 186 num_images = len(img_names) 187 if num_images < 2: 188 print("Need more images") 189 exit() 190 191 if estimator_type == "affine": 192 estimator = cv.detail_AffineBasedEstimator()#基于仿射变换特征 193 else: 194 estimator = cv.detail_HomographyBasedEstimator()#基于透视变换的单应阵 195 b, cameras =estimator.apply(features,p,None)#变换估计出变换矩阵和相机参数 196 if not b: 197 print("Homography estimation failed.") 198 exit() 199 for cam in cameras: 200 #相机的旋转矩阵 201 cam.R=cam.R.astype(np.float32) 202 203 if ba_cost_func == "reproj": 204 adjuster = cv.detail_BundleAdjusterReproj() 205 elif ba_cost_func == "ray":#默认,光束平差的损失函数类型 206 adjuster = cv.detail_BundleAdjusterRay() 207 elif ba_cost_func == "affine": 208 adjuster = cv.detail_BundleAdjusterAffinePartial() 209 elif ba_cost_func == "no": 210 adjuster = cv.detail_NoBundleAdjuster() 211 else: 212 print( "Unknown bundle adjustment cost function: ", ba_cost_func ) 213 exit() 214 adjuster.setConfThresh(1)#确定两幅图片属于同一张全景图的阈值 215 refine_mask=np.zeros((3,3),np.uint8) 216 #掩膜,3*3的0阵,如果要优化则将对应位置赋成1 217 """数字图像处理中,掩模为二维矩阵数组,有时也用多值图像。数字图像处理中,图像掩模主要用于: 218 219 ①提取感兴趣区,用预先制作的感兴趣区掩模与待处理图像相乘,得到感兴趣区图像,感兴趣区内图像值保持不变,而区外图像值都为0。 220 ②屏蔽作用,用掩模对图像上某些区域作屏蔽,使其不参加处理或不参加处理参数的计算,或仅对屏蔽区作处理或统计。 221 ③结构特征提取,用相似性变量或图像匹配方法检测和提取图像中与掩模相似的结构特征。 222 """ 223 224 if ba_refine_mask[0] == 'x': 225 refine_mask[0,0] = 1 226 if ba_refine_mask[1] == 'x': 227 refine_mask[0,1] = 1 228 if ba_refine_mask[2] == 'x': 229 refine_mask[0,2] = 1 230 if ba_refine_mask[3] == 'x': 231 refine_mask[1,1] = 1 232 if ba_refine_mask[4] == 'x': 233 refine_mask[1,2] = 1 234 adjuster.setRefinementMask(refine_mask) 235 b,cameras = adjuster.apply(features,p,cameras) 236 #光束平差,精修相机姿态等 237 #光束是指相机光心发出的射线。损失函数一般有重映射误差和射线距离误差两种 238 if not b: 239 print("Camera parameters adjusting failed.") 240 exit() 241 focals=[] 242 for cam in cameras: 243 focals.append(cam.focal) #焦距 244 sorted(focals) 245 #可以通过两张图的匹配点或单应阵求出焦距,用所有焦距的平均作为全景图焦距的近似。用这个焦距作为柱面的半径,先将图像投射到柱面,在柱面上区拼接 246 if len(focals)%2==1: 247 #根据焦距确定 图片投影的尺度 warped_image_scale 248 warped_image_scale = focals[len(focals) // 2] 249 else: 250 warped_image_scale = (focals[len(focals) // 2]+focals[len(focals) // 2-1])/2 251 if do_wave_correct: 252 #波形修正得到修正后的旋转矩阵 253 rmats=[] 254 for cam in cameras: 255 rmats.append(np.copy(cam.R)) 256 rmats = cv.detail.waveCorrect( rmats, cv.detail.WAVE_CORRECT_HORIZ) 257 for idx,cam in enumerate(cameras): 258 cam.R = rmats[idx] 259 corners=[] 260 mask=[] 261 masks_warped=[] 262 images_warped=[] 263 sizes=[] 264 masks=[] 265 for i in range(0,num_images): 266 um=cv.UMat(255*np.ones((images[i].shape[0],images[i].shape[1]),np.uint8)) 267 masks.append(um) 268 #处理相机旋转造成的扭曲 269 warper = cv.PyRotationWarper(warp_type,warped_image_scale*seam_work_aspect) # warper peut etre nullptr? 270 for idx in range(0,num_images): 271 K = cameras[idx].K().astype(np.float32)#相机内参 272 swa = seam_work_aspect 273 K[0,0] *= swa 274 K[0,2] *= swa 275 K[1,1] *= swa 276 K[1,2] *= swa 277 corner,image_wp =warper.warp(images[idx],K,cameras[idx].R,cv.INTER_LINEAR, cv.BORDER_REFLECT) 278 corners.append(corner) 279 sizes.append((image_wp.shape[1],image_wp.shape[0])) 280 images_warped.append(image_wp) 281 282 p,mask_wp =warper.warp(masks[idx],K,cameras[idx].R,cv.INTER_NEAREST, cv.BORDER_CONSTANT)#线性差值,连续边界处理 283 masks_warped.append(mask_wp.get()) 284 images_warped_f=[] 285 for img in images_warped:#对投影后的图片,单精度存储 286 imgf=img.astype(np.float32) 287 images_warped_f.append(imgf) 288 if cv.detail.ExposureCompensator_CHANNELS == expos_comp_type: 289 compensator = cv.detail_ChannelsCompensator(expos_comp_nr_feeds) 290 # compensator.setNrGainsFilteringIterations(expos_comp_nr_filtering) 291 elif cv.detail.ExposureCompensator_CHANNELS_BLOCKS == expos_comp_type: 292 compensator=cv.detail_BlocksChannelsCompensator(expos_comp_block_size, expos_comp_block_size,expos_comp_nr_feeds) 293 # compensator.setNrGainsFilteringIterations(expos_comp_nr_filtering) 294 else: 295 compensator=cv.detail.ExposureCompensator_createDefault(expos_comp_type) 296 #曝光补偿需要 分块 投影图片 和 掩膜 297 #对于分块补偿,得到分块曝光补偿系数后,此时直接补偿会出现明显的视觉分"块"效果还要对系数进行分段线性滤波处理,得到全局最终的补偿系数 298 compensator.feed(corners=corners, images=images_warped, masks=masks_warped) 299 #缝合线(两图最佳的一条融合曲线)查找描述子 300 #一般有三种方法:1基于距离的逐点法 2.动态规划查找法(快,省内存) 3.最大流图割法 301 if seam_find_type == "no": 302 seam_finder = cv.detail.SeamFinder_createDefault(cv.detail.SeamFinder_NO) 303 elif seam_find_type == "voronoi": 304 seam_finder = cv.detail.SeamFinder_createDefault(cv.detail.SeamFinder_VORONOI_SEAM); 305 elif seam_find_type == "gc_color": 306 seam_finder = cv.detail_GraphCutSeamFinder("COST_COLOR") 307 elif seam_find_type == "gc_colorgrad": 308 seam_finder = cv.detail_GraphCutSeamFinder("COST_COLOR_GRAD") 309 elif seam_find_type == "dp_color": 310 seam_finder = cv.detail_DpSeamFinder("COLOR") 311 elif seam_find_type == "dp_colorgrad": 312 seam_finder = cv.detail_DpSeamFinder("COLOR_GRAD") 313 if seam_finder is None: 314 print("Can't create the following seam finder ",seam_find_type) 315 exit() 316 seam_finder.find(images_warped_f, corners,masks_warped ) 317 imgListe=[] 318 compose_scale=1 319 corners=[] 320 sizes=[] 321 images_warped=[] 322 images_warped_f=[] 323 masks=[] 324 blender= None 325 timelapser=None 326 compose_work_aspect=1 327 for idx,name in enumerate(img_names): # https://github.com/opencv/opencv/blob/master/samples/cpp/stitching_detailed.cpp#L725 ? 328 full_img = cv.imread(name) 329 if not is_compose_scale_set: 330 if compose_megapix > 0: 331 #exp是2.718; e是10^ 332 compose_scale = min(1.0, np.sqrt(compose_megapix * 1e6 / (full_img.shape[0]*full_img.shape[1]))) 333 is_compose_scale_set = True; 334 compose_work_aspect = compose_scale / work_scale; 335 warped_image_scale *= compose_work_aspect 336 warper = cv.PyRotationWarper(warp_type,warped_image_scale) 337 for i in range(0,len(img_names)): 338 cameras[i].focal *= compose_work_aspect 339 cameras[i].ppx *= compose_work_aspect 340 cameras[i].ppy *= compose_work_aspect 341 sz = (full_img_sizes[i][0] * compose_scale,full_img_sizes[i][1]* compose_scale) 342 K = cameras[i].K().astype(np.float32) 343 roi = warper.warpRoi(sz, K, cameras[i].R); 344 #缝合区获得 345 corners.append(roi[0:2]) 346 sizes.append(roi[2:4]) 347 if abs(compose_scale - 1) > 1e-1: 348 img =cv.resize(src=full_img, dsize=None, fx=compose_scale, fy=compose_scale, interpolation=cv.INTER_LINEAR_EXACT) 349 else: 350 img = full_img; 351 img_size = (img.shape[1],img.shape[0]); 352 K=cameras[idx].K().astype(np.float32) 353 corner,image_warped =warper.warp(img,K,cameras[idx].R,cv.INTER_LINEAR, cv.BORDER_REFLECT) 354 mask =255*np.ones((img.shape[0],img.shape[1]),np.uint8) 355 p,mask_warped =warper.warp(mask,K,cameras[idx].R,cv.INTER_NEAREST, cv.BORDER_CONSTANT) 356 #块增益补偿 357 compensator.apply(idx,corners[idx],image_warped,mask_warped) 358 image_warped_s = image_warped.astype(np.int16) 359 image_warped=[] 360 #膨胀 361 dilated_mask = cv.dilate(masks_warped[idx],None) 362 seam_mask = cv.resize(dilated_mask,(mask_warped.shape[1],mask_warped.shape[0]),0,0,cv.INTER_LINEAR_EXACT) 363 mask_warped = cv.bitwise_and(seam_mask,mask_warped)#二进制的与运算,1&1 == 1 ,1&0 == 0 364 if blender==None and not timelapse: 365 blender = cv.detail.Blender_createDefault(cv.detail.Blender_NO) 366 dst_sz = cv.detail.resultRoi(corners=corners,sizes=sizes) 367 blend_width = np.sqrt(dst_sz[2]*dst_sz[3]) * blend_strength / 100 368 if blend_width < 1: 369 blender = cv.detail.Blender_createDefault(cv.detail.Blender_NO) 370 elif blend_type == "multiband": 371 blender = cv.detail_MultiBandBlender() 372 blender.setNumBands((np.log(blend_width)/np.log(2.) - 1.).astype(np.int)) 373 elif blend_type == "feather":#羽化融合 374 blender = cv.detail_FeatherBlender() 375 blender.setSharpness(1./blend_width) 376 blender.prepare(dst_sz) 377 elif timelapser==None and timelapse: 378 timelapser = cv.detail.Timelapser_createDefault(timelapse_type) 379 timelapser.initialize(corners, sizes) 380 if timelapse: 381 #延时处理 382 matones=np.ones((image_warped_s.shape[0],image_warped_s.shape[1]), np.uint8) 383 timelapser.process(image_warped_s, matones, corners[idx]) 384 pos_s = img_names[idx].rfind("/"); 385 if pos_s == -1: 386 fixedFileName = "fixed_" + img_names[idx]; 387 else: 388 fixedFileName = img_names[idx][:pos_s + 1 ]+"fixed_" + img_names[idx][pos_s + 1: ] 389 cv.imwrite(fixedFileName, timelapser.getDst()) 390 else: 391 blender.feed(cv.UMat(image_warped_s), mask_warped, corners[idx]) 392 if not timelapse: 393 result=None 394 result_mask=None 395 result,result_mask = blender.blend(result,result_mask) 396 cv.imwrite(result_name,result) 397 zoomx = 600.0 / result.shape[1] 398 dst=cv.normalize(src=result,dst=None,alpha=255.,norm_type=cv.NORM_MINMAX,dtype=cv.CV_8U) 399 dst=cv.resize(dst,dsize=None,fx=zoomx,fy=zoomx) 400 cv.imshow(result_name,dst) 401 cv.waitKey() 402 403 print('Done') 404 405 406 if __name__ == '__main__': 407 print(__doc__) 408 main() 409 cv.destroyAllWindows()